FMD/scripts/DrawESD.C

   1 //____________________________________________________________________
   2 //
   3 // $Id$
   4 //
   5 // Script that contains a class to draw eloss from hits, versus ADC
   6 // counts from digits, using the AliFMDInputHits class in the util library.
   7 //
   8 // It draws the energy loss versus the p/(mq^2).  It can be overlayed
   9 // with the Bethe-Bloc curve to show how the simulation behaves
  10 // relative to the expected.
  11 //
  12 // Use the script `Compile.C' to compile this class using ACLic.
  13 //
  14 #include <TH1D.h>
  15 #include <AliFMDHit.h>
  16 #include <AliFMDDigit.h>
  17 #include <AliFMDInput.h>
  18 #include <AliFMDUShortMap.h>
  19 #include <AliFMDFloatMap.h>
  20 #include <AliFMDRecPoint.h>
  21 #include <AliESDFMD.h>
  22 #include <AliLog.h>
  23 #include <iostream>
  24 #include <TStyle.h>
  25 #include <TArrayF.h>
  26 #include <TCanvas.h>
  27 #include <TMath.h>
  28 #include <TF1.h>
  29 #include <TSpectrum.h>
  30 #include <TLegend.h>
  31 #include <TLine.h>
  32
  33 /** @class DrawESD
  34     @brief Draw digit ADC versus Rec point mult
  35     @code
  36     Root> .L Compile.C
  37     Root> Compile("DrawESD.C")
  38     Root> DrawESD c
  39     Root> c.Run();
  40     @endcode
  41     @ingroup FMD_script
  42  */
  43 class DrawESD : public AliFMDInput
  44 {
  45 private:
  46   TH1D* fMult; // Histogram
  47   const Double_t fCorr;
  48 public:
  49   //__________________________________________________________________
  50   TArrayF MakeLogScale(Int_t n, Double_t min, Double_t max)
  51   {
  52     TArrayF bins(n+1);
  53     Float_t dp   = n / TMath::Log10(max / min);
  54     Float_t pmin = TMath::Log10(min);
  55     bins[0]      = min;
  56     for (Int_t i = 1; i < n+1; i++) {
  57       Float_t p = pmin + i / dp;
  58       bins[i]   = TMath::Power(10, p);
  59     }
  60     return bins;
  61   }
  62   //__________________________________________________________________
  63   DrawESD(Int_t n=420, Double_t mmin=-0.5, Double_t mmax=20.5)
  64     : fCorr(1) // 0.68377 / 1.1)
  65   {
  66     AddLoad(kESD);
  67     fMult = new TH1D("mult", " Multiplicity (strip)", n, mmin, mmax);
  68     fMult->SetXTitle("Strip Multiplicity");
  69   }
  70   //__________________________________________________________________
  71   /** Begining of event
  72       @param ev Event number
  73       @return @c false on error */
  74   Bool_t Begin(Int_t ev)
  75   {
  76     return AliFMDInput::Begin(ev);
  77   }
  78   //__________________________________________________________________
  79   Bool_t ProcessESD(UShort_t /* det */,
  80                     Char_t   /* ring */,
  81                     UShort_t /* sec */,
  82                     UShort_t /* strip */,
  83                     Float_t  /* eta */,
  84                     Float_t  mult)
  85   {
  86     // Cache the energy loss
  87     fMult->Fill(mult/fCorr);
  88     return kTRUE;
  89   }
  90   //__________________________________________________________________
  91   Bool_t Finish()
  92   {
  93     gStyle->SetPalette(1);
  94     gStyle->SetOptTitle(0);
  95     gStyle->SetOptFit(1111111);
  96     gStyle->SetCanvasColor(0);
  97     gStyle->SetCanvasBorderSize(0);
  98     gStyle->SetPadColor(0);
  99     gStyle->SetPadBorderSize(0);
 100
 101     TCanvas* c = new TCanvas("c", "C");
 102     c->cd();
 103     c->SetLogy();
 104     // fMult->GetXaxis()->SetRangeUser(0.2,4);
 105     // fMult->Sumw2();
 106     // fMult->Scale(1. / fMult->GetMaximum());
 107     fMult->SetStats(kFALSE);
 108     fMult->SetFillColor(2);
 109     fMult->SetFillStyle(3001);
 110     fMult->Draw();
 111     Double_t x1 = .75; // .8;  // .65 / fCorr;
 112     Double_t x2 = 1.3; // 1.7; // fCorr;
 113     Double_t x3 = 2.5; // 2.7; // fCorr;
 114     Double_t x4 = 3.7; // fCorr;
 115     TF1* l1 = new TF1("landau1", "landau", x1, x2);
 116     TF1* l2 = new TF1("landau2", "landau", x2, x3);
 117     // TF1* l3 = new TF1("landau3", "landau", x3, x4);
 118     TF1* f  = new TF1("user", "landau(0)+landau(3)", x1, x3);
 119
 120     fMult->SetStats(kTRUE);
 121     fMult->GetXaxis()->SetRangeUser(0, 4);
 122     fMult->Fit(l1, "E0L", "", x1, x2);
 123     fMult->Fit(l2, "E0L", "", x2, x3);
 124     // fMult->Fit(l3, "E0L", "", x3, x4);
 125     f->SetParNames("A_{1}", "Mpv_{1}", "#sigma_{1}",
 126                    "A_{2}", "Mpv_{2}", "#sigma_{2}",
 127                    "A_{3}", "Mpv_{3}", "#sigma_{3}");
 128     f->SetParameters(l1->GetParameter(0),
 129                      l1->GetParameter(1),
 130                      l1->GetParameter(2),
 131                      l2->GetParameter(0),
 132                      l2->GetParameter(1),
 133                      l2->GetParameter(2),
 134                      l1->GetParameter(0),
 135                      l1->GetParameter(1) * 3,
 136                      l1->GetParameter(2) * 3);
 137     fMult->Fit(f, "E0L", "", x1, x3);
 138     fMult->Fit(f, "MEL", "E1", x1, x3);
 139
 140     TH1* h = static_cast<TH1*>(fMult->DrawClone("same hist"));
 141
 142     // l1->SetLineStyle(2);
 143     l1->SetLineColor(3);
 144     l1->SetLineWidth(2);
 145     l1->SetRange(0, 4);
 146     l1->Draw("same");
 147     // l2->SetLineStyle(3);
 148     l2->SetLineColor(4);
 149     l2->SetLineWidth(2);
 150     l2->SetRange(0, 4);
 151     l2->Draw("same");
 152     // l3->SetLineStyle(3);
 153     // l3->SetLineWidth(2);
 154     // l3->SetRange(0, 5);
 155     // l3->Draw("same");
 156     f->SetLineWidth(2);
 157     f->SetRange(0, 4);
 158     f->Draw("same");
 159
 160     TLegend* l = new TLegend(0.2, 0.6, .6, .89);
 161     l->AddEntry(l1, "1 particle Landau", "l");
 162     l->AddEntry(l2, "2 particle Landau", "l");
 163     l->AddEntry(f,  "1+2 particle Landau", "l");
 164     l->SetFillColor(0);
 165     l->Draw("same");
 166
 167
 168     c = new TCanvas("c2", "Landaus");
 169     c->SetLogy();
 170     fMult->DrawClone("axis");
 171     f->Draw("same");
 172     TF1* ll1 = new TF1("ll1", "landau", 0, 4);
 173     ll1->SetParameters(f->GetParameter(0),
 174                        f->GetParameter(1),
 175                        f->GetParameter(2));
 176     ll1->SetLineColor(3);
 177     ll1->Draw("same");
 178     TF1* ll2 = new TF1("ll2", "landau", 0, 4);
 179     ll2->SetParameters(f->GetParameter(3),
 180                        f->GetParameter(4),
 181                        f->GetParameter(5));
 182     ll2->SetLineColor(4);
 183     ll2->Draw("same");
 184
 185     Double_t y1  = fMult->GetMinimum() * 1.1;
 186     Double_t y2  = fMult->GetMaximum() * .9;
 187     Double_t xc1 = f->GetParameter(1)-3*f->GetParameter(2);
 188     Double_t xc2 = f->GetParameter(4)-2*f->GetParameter(5);
 189     TLine* c1 = new TLine(xc1, y1, xc1, y2);
 190     c1->Draw("same");
 191     TLine* c2 = new TLine(xc2, y1, xc2, y2);
 192     c2->Draw("same");
 193
 194     l = new TLegend(0.2, 0.6, .6, .89);
 195     l->AddEntry(ll1, "1 particle Landau", "l");
 196     l->AddEntry(ll2, "2 particle Landau", "l");
 197     l->AddEntry(f,  "1+2 particle Landau", "l");
 198     l->SetFillColor(0);
 199     l->Draw("same");
 200
 201 #if 0
 202     c = new TCanvas("s", "Spectrum");
 203     TSpectrum* s = new TSpectrum(16);
 204     h->GetXaxis()->SetRangeUser(0.3, 20);
 205     s->Search(h, .15, "", 0.01);
 206     c->Update();
 207     TH1* b = s->Background(h);
 208     b->SetFillColor(4);
 209     b->SetFillStyle(3001);
 210     b->Draw("same");
 211 #endif
 212
 213     return kTRUE;
 214   }
 215
 216   ClassDef(DrawESD,0);
 217
 218 };
 219
 220 //____________________________________________________________________
 221 //
 222 // EOF
 223 //