FMD/scripts/AliFMDAnaBg.C

   1 #include <AliFMDHit.h>
   2 #include <TH1D.h>
   3 #include <TH2D.h>
   4 #include <TParticle.h>
   5 #include <TBrowser.h>
   6 #include <TArrayF.h>
   7 #include <TCanvas.h>
   8 #include <TStyle.h>
   9 #include <TLegend.h>
  10 #include <TProfile.h>
  11 #include <TProfile2D.h>
  12 #include <TROOT.h>
  13 #include <AliFMDInput.h>
  14 #include <AliFMDAnaRing.h>
  15 #include <AliFMDGeometry.h>
  16 #include <AliCDBManager.h>
  17 #include <iostream>
  18
  19 //====================================================================
  20 class AliFMDAnaBgRing : public AliFMDAnaRing
  21 {
  22 public:
  23   AliFMDAnaBgRing(UShort_t det, Char_t ring,
  24                   Float_t deta, Float_t etaMin, Float_t etaMax)
  25     : AliFMDAnaRing(det, ring, 0, 0, 0),
  26       fTotal(Form("total_%s", fName), "Total number of hits",
  27              Int_t((etaMax-etaMin)/deta+.5), etaMin, etaMax,
  28              fNSeq, 0, 2*TMath::Pi()),
  29       fPrimary(Form("primary_%s", fName), "Number of hits from primary",
  30                Int_t((etaMax-etaMin)/deta+.5), etaMin, etaMax,
  31                fNSeq, 0, 2*TMath::Pi()),
  32       fRatio(Form("ratio_%s", fName), "Correction factors",
  33                Int_t((etaMax-etaMin)/deta+.5), etaMin, etaMax,
  34              fNSeq, 0, 2*TMath::Pi()),
  35       fInverse(Form("inverse_%s", fName), "Inverse Correction factors",
  36                Int_t((etaMax-etaMin)/deta+.5), etaMin, etaMax,
  37                fNSeq, 0, 2*TMath::Pi())
  38   {
  39     fTotal.SetDirectory(0);
  40     fTotal.Sumw2();
  41     fTotal.SetXTitle("#eta");
  42     fTotal.SetYTitle("#varphi");
  43     fTotal.SetZTitle("# of hits");
  44     fPrimary.SetDirectory(0);
  45     fPrimary.Sumw2();
  46     fPrimary.SetXTitle("#eta");
  47     fPrimary.SetYTitle("#varphi");
  48     fPrimary.SetZTitle("# of primary");
  49     fRatio.SetDirectory(0);
  50     fRatio.Sumw2();
  51     fRatio.SetXTitle("#eta");
  52     fRatio.SetYTitle("#varphi");
  53     fRatio.SetZTitle("# of primary/# of hits");
  54     fInverse.SetDirectory(0);
  55     fInverse.Sumw2();
  56     fInverse.SetXTitle("#eta");
  57     fInverse.SetYTitle("#varphi");
  58     fInverse.SetZTitle("# of hits/# of primary");
  59   }
  60   void Fill(Float_t phi, Float_t eta, Float_t mult)
  61   {
  62     fTotal.Fill(eta, phi);
  63     if (mult >= 1) fPrimary.Fill(eta,phi);
  64   }
  65   void Finish()
  66   {
  67     fRatio.Add(&fPrimary);
  68     fRatio.Divide(&fTotal);
  69     fInverse.Add(&fTotal);
  70     fInverse.Divide(&fPrimary);
  71   }
  72   void Browse(TBrowser* b)
  73   {
  74     // AliFMDAnaRing::Browse(b);
  75     b->Add(&fTotal);
  76     b->Add(&fPrimary);
  77     b->Add(&fRatio);
  78     b->Add(&fInverse);
  79   }
  80   TH2D fTotal;
  81   TH2D fPrimary;
  82   TH2D fRatio;
  83   TH2D fInverse;
  84 };
  85
  86 //====================================================================
  87 class AliFMDAnaBg : public AliFMDInput
  88 {
  89 public:
  90   //__________________________________________________________________
  91   /** Constructor */
  92   AliFMDAnaBg(Float_t deta=.1)
  93     : fEtaMin(-3.7),
  94       fEtaMax(5.3),
  95       fTotal("total", "Total number of hits",
  96              Int_t((fEtaMax-fEtaMin)/deta+.5), fEtaMin, fEtaMax,
  97              40, 0, 2*TMath::Pi()),
  98       fPrimary("primary", "Number of hits from primary",
  99                Int_t((fEtaMax-fEtaMin)/deta+.5), fEtaMin, fEtaMax,
 100                40, 0, 2*TMath::Pi()),
 101       fRatio("ratio", "Correction factors",
 102              Int_t((fEtaMax-fEtaMin)/deta+.5), fEtaMin, fEtaMax,
 103              40, 0, 2*TMath::Pi())
 104   {
 105     fTotal.SetDirectory(0);
 106     fTotal.Sumw2();
 107     fTotal.SetXTitle("#eta");
 108     fTotal.SetYTitle("#varphi");
 109     fTotal.SetZTitle("# of hits");
 110     fPrimary.SetDirectory(0);
 111     fPrimary.Sumw2();
 112     fPrimary.SetXTitle("#eta");
 113     fPrimary.SetYTitle("#varphi");
 114     fPrimary.SetZTitle("# of primary");
 115     fRatio.SetDirectory(0);
 116     fRatio.Sumw2();
 117     fRatio.SetXTitle("#eta");
 118     fRatio.SetYTitle("#varphi");
 119     fRatio.SetZTitle("# of primary/# of hits");
 120
 121     AliCDBManager::Instance()->SetRun(0);
 122
 123     AddLoad(kHits);
 124     AddLoad(kKinematics);
 125     AddLoad(kGeometry);
 126     for (int i = 0; i < 5; i++) {
 127       UShort_t det = (i == 0 ? 1 : (i <= 2 ? 2 : 3));
 128       Char_t   rng = (i == 0 ? 'I' :  (i % 2 == 1 ? 'I' : 'O'));
 129       fRing[i] = new AliFMDAnaBgRing(det,rng,deta,EtaMin(i)-.2,EtaMax(i)+.12);
 130     }
 131   }
 132   //__________________________________________________________________
 133   Int_t FindRing(UShort_t det, Char_t ring) const
 134   {
 135     Int_t idx = -1;
 136     switch (det) {
 137     case 1: idx = 0; break;
 138     case 2: idx = 1 + (ring == 'O' || ring == 'o' ? 1 : 0); break;
 139     case 3: idx = 3 + (ring == 'O' || ring == 'o' ? 1 : 0); break;
 140     }
 141     return idx;
 142   }
 143   //__________________________________________________________________
 144   /** Initialize */
 145   Bool_t Init()
 146   {
 147     Bool_t ret = AliFMDInput::Init();
 148     for (size_t i = 0; i < 5; i++) fRing[i]->Init();
 149     AliFMDGeometry::Instance()->Init();
 150     AliFMDGeometry::Instance()->InitTransformations();
 151     return ret;
 152   }
 153   //__________________________________________________________________
 154   /** Called at the beginning of an event */
 155   Bool_t Begin(Int_t ev)
 156   {
 157     for (size_t i = 0; i < 5; i++) fRing[i]->Begin();
 158     return AliFMDInput::Begin(ev);
 159   }
 160   //__________________________________________________________________
 161   /** Called at the beginning of an event */
 162   Bool_t ProcessHit(AliFMDHit* hit, TParticle* p)
 163   {
 164     Float_t          mult = p->IsPrimary() ? 1 : 0;
 165     UShort_t         det  = hit->Detector();
 166     Char_t           rng  = hit->Ring();
 167     UShort_t         sec  = hit->Sector();
 168     UShort_t         str  = hit->Strip();
 169     AliFMDAnaBgRing* frn  = fRing[FindRing(det, rng)];
 170     Double_t x, y, z;
 171     AliFMDGeometry::Instance()->Detector2XYZ(det,rng,sec,str,x,y,z);
 172     Double_t         phi  = TMath::ATan2(y,x);
 173     Double_t         r    = TMath::Sqrt(x * x + y * y);
 174     Double_t         the  = TMath::ATan2(r, z);
 175     Double_t         eta  = -TMath::Log(TMath::Tan(the/2));
 176     if (phi < 0) phi += 2*TMath::Pi();
 177     frn->Fill(phi, eta, mult);
 178     Fill(phi, eta, mult);
 179     return kTRUE;
 180   }
 181
 182   //__________________________________________________________________
 183   /** Fill in one strip */
 184   void Fill(Float_t phi, Float_t eta, Float_t mult)
 185   {
 186     fTotal.Fill(eta, phi);
 187     if (mult >= 1) fPrimary.Fill(eta,phi);
 188   }
 189   //__________________________________________________________________
 190   /** Called at the end of an event */
 191   Bool_t End()
 192   {
 193     for (UShort_t i = 0; i < 5; i++) fRing[i]->End();
 194     return AliFMDInput::End();
 195   }
 196   //__________________________________________________________________
 197   /** called at the end of a run */
 198   Bool_t Finish()
 199   {
 200     gStyle->SetPalette(1);
 201     gStyle->SetOptTitle(0);
 202     gStyle->SetOptStat(0);
 203
 204     fRatio.Add(&fPrimary);
 205     fRatio.Divide(&fTotal);
 206
 207     for (UShort_t i = 0; i < 5; i++) fRing[i]->Finish();
 208
 209     TCanvas* c = new TCanvas("c1", "Result", 800, 600);
 210     c->cd();
 211     c->SetGridy();
 212     c->SetTopMargin(0.02);
 213     c->SetRightMargin(0.02);
 214     c->SetFillColor(0);
 215     c->SetBorderSize(0);
 216     c->SetBorderMode(0);
 217     c->SetHighLightColor(0);
 218     c->Divide(5,1,0,0);
 219     for (UShort_t i = 0; i < 5; i++) {
 220       c->cd(i+1);
 221       fRing[i]->fRatio.Draw("colz");
 222     }
 223
 224     c = new TCanvas("c2", "All", 800, 600);
 225     c->cd();
 226     c->SetGridy();
 227     c->SetTopMargin(0.02);
 228     c->SetRightMargin(0.02);
 229     c->SetFillColor(0);
 230     c->SetBorderSize(0);
 231     c->SetBorderMode(0);
 232     c->SetHighLightColor(0);
 233     c->Divide(5,1,0,0);
 234     fRatio.Draw("colz");
 235
 236     // gROOT->GetListOfSpecials()->Add(this);
 237     return kTRUE;
 238   }
 239   //__________________________________________________________________
 240   Bool_t IsFolder() const { return kTRUE; }
 241   //__________________________________________________________________
 242   void Browse(TBrowser* b)
 243   {
 244     for (size_t i = 0; i < 5; i++)
 245       b->Add(fRing[i], fRing[i]->Name());
 246   }
 247   //__________________________________________________________________
 248   /** Get min pseudo rapidity of ring idx */
 249   Float_t EtaMin(Int_t idx) const {
 250     switch (idx) {
 251     case 0: return  3.61728; break; // FMD1I
 252     case 1: return  2.18122; break; // FMD2I
 253     case 2: return  1.79821; break; // FMD2O
 254     case 3: return -3.39913; break; // FMD3I
 255     case 4: return -2.28925; break; // FMD3O
 256     }
 257     return 0;
 258   }
 259   //__________________________________________________________________
 260   /** Get max pseudo rapidity of ring idx */
 261   Float_t EtaMax(Int_t idx) const {
 262     switch (idx) {
 263     case 0: return  5.02643; break; // FMD1I
 264     case 1: return  3.57899; break; // FMD2I
 265     case 2: return  2.39085; break; // FMD2O
 266     case 3: return -2.00644; break; // FMD3I
 267     case 4: return -1.70080; break; // FMD3O
 268     }
 269     return 0;
 270   }
 271   /** Flow rings */
 272   AliFMDAnaBgRing* fRing[5];
 273   Float_t fEtaMin;
 274   Float_t fEtaMax;
 275   TH2D fTotal;
 276   TH2D fPrimary;
 277   TH2D fRatio;
 278   ClassDef(AliFMDAnaBg,1)
 279 };
 280
 281 //____________________________________________________________________
 282 //
 283 // EOF
 284 //