PWG2/FORWARD/analysis2/AliFMDCorrections.cxx

   1 #include "AliFMDCorrections.h"
   2 #include <AliESDFMD.h>
   3 #include <TAxis.h>
   4 #include <TList.h>
   5 #include <TMath.h>
   6 #include "AliFMDAnaParameters.h"
   7 #include "AliLog.h"
   8 #include <TH2D.h>
   9
  10 ClassImp(AliFMDCorrections)
  11 #if 0
  12 ; // For Emacs
  13 #endif
  14
  15 //____________________________________________________________________
  16 AliFMDCorrections::AliFMDCorrections()
  17   : TNamed(),
  18     fRingHistos(),
  19     fMultCut(0.3)
  20 {}
  21
  22 //____________________________________________________________________
  23 AliFMDCorrections::AliFMDCorrections(const char* title)
  24   : TNamed("fmdCorrections", title),
  25     fRingHistos(),
  26     fMultCut(0.3)
  27 {
  28   fRingHistos.SetName(GetName());
  29   fRingHistos.Add(new RingHistos(1, 'I'));
  30   fRingHistos.Add(new RingHistos(2, 'I'));
  31   fRingHistos.Add(new RingHistos(2, 'O'));
  32   fRingHistos.Add(new RingHistos(3, 'I'));
  33   fRingHistos.Add(new RingHistos(3, 'O'));
  34 }
  35
  36 //____________________________________________________________________
  37 AliFMDCorrections::AliFMDCorrections(const AliFMDCorrections& o)
  38   : TNamed(o),
  39     fRingHistos(),
  40     fMultCut(o.fMultCut)
  41 {
  42   TIter    next(&o.fRingHistos);
  43   TObject* obj = 0;
  44   while ((obj = next())) fRingHistos.Add(obj);
  45 }
  46
  47 //____________________________________________________________________
  48 AliFMDCorrections::~AliFMDCorrections()
  49 {
  50   fRingHistos.Delete();
  51 }
  52
  53 //____________________________________________________________________
  54 AliFMDCorrections&
  55 AliFMDCorrections::operator=(const AliFMDCorrections& o)
  56 {
  57   SetName(o.GetName());
  58   SetTitle(o.GetTitle());
  59
  60   fMultCut = o.fMultCut;
  61
  62   fRingHistos.Delete();
  63   TIter    next(&o.fRingHistos);
  64   TObject* obj = 0;
  65   while ((obj = next())) fRingHistos.Add(obj);
  66
  67   return *this;
  68 }
  69
  70 //____________________________________________________________________
  71 AliFMDCorrections::RingHistos*
  72 AliFMDCorrections::GetRingHistos(UShort_t d, Char_t r) const
  73 {
  74   Int_t idx = -1;
  75   switch (d) {
  76   case 1: idx = 0; break;
  77   case 2: idx = 1 + (r == 'I' || r == 'i' ? 0 : 1); break;
  78   case 3: idx = 3 + (r == 'I' || r == 'i' ? 0 : 1); break;
  79   }
  80   if (idx < 0 || idx >= fRingHistos.GetEntries()) return 0;
  81
  82   return static_cast<RingHistos*>(fRingHistos.At(idx));
  83 }
  84
  85 //____________________________________________________________________
  86 Bool_t
  87 AliFMDCorrections::Correct(AliForwardUtil::Histos& hists,
  88                            Int_t                   vtxbin)
  89 {
  90   AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
  91
  92   for (UShort_t d=1; d<=3; d++) {
  93     UShort_t nr = (d == 1 ? 1 : 2);
  94     for (UShort_t q=0; q<nr; q++) {
  95       Char_t      r = (q == 0 ? 'I' : 'O');
  96       TH2D*       h = hists.Get(d,r);
  97       RingHistos* rh= GetRingHistos(d,r);
  98       TH2F*       bg= pars->GetBackgroundCorrection(d, r, vtxbin);
  99       TH2F*       ef= pars->GetEventSelectionEfficiency("INEL",vtxbin,r);
 100       if (!bg) {
 101         AliWarning(Form("No secondary correction for FMDM%d%c in vertex bin %d",
 102                         d, r, vtxbin));
 103         continue;
 104       }
 105       if (!ef) {
 106         AliWarning(Form("No event %s selection efficiency in vertex bin %d",
 107                         "INEL", vtxbin));
 108         continue;
 109       }
 110
 111       // Divide by primary/total ratio
 112       h->Divide(bg);
 113
 114       // Divide by the event selection efficiency
 115       h->Divide(ef);
 116
 117       if(!pars->SharingEffPresent()) {
 118         AliWarning("No sharing efficiencies");
 119         continue;
 120       }
 121
 122       TH1F* sf = pars->GetSharingEfficiencyTrVtx(d,r,vtxbin);
 123
 124       for (Int_t ieta = 1; ieta <= h->GetNbinsX(); ieta++) {
 125         Float_t c  = sf->GetBinContent(ieta);
 126         Float_t ec = sf->GetBinError(ieta);
 127
 128         if (c == 0) continue;
 129
 130         for (Int_t iphi = 1; iphi <= h->GetNbinsY(); iphi++) {
 131           Double_t m  = h->GetBinContent(ieta, iphi) / c;
 132           Double_t em = h->GetBinError(ieta, iphi);
 133
 134           Double_t e  = TMath::Sqrt(em * em + (m * ec) * (m * ec));
 135
 136           h->SetBinContent(ieta,iphi,m);
 137           h->SetBinError(ieta,iphi,e);
 138         }
 139       }
 140
 141       rh->fDensity->Add(h);
 142     }
 143   }
 144
 145   return kTRUE;
 146 }
 147
 148 //____________________________________________________________________
 149 void
 150 AliFMDCorrections::ScaleHistograms(TList* dir, Int_t nEvents)
 151 {
 152   if (nEvents <= 0) return;
 153   TList* d = static_cast<TList*>(dir->FindObject(GetName()));
 154   if (!d) return;
 155
 156   TIter    next(&fRingHistos);
 157   RingHistos* o = 0;
 158   while ((o = static_cast<RingHistos*>(next())))
 159     o->ScaleHistograms(d, nEvents);
 160 }
 161 //____________________________________________________________________
 162 void
 163 AliFMDCorrections::DefineOutput(TList* dir)
 164 {
 165   TList* d = new TList;
 166   d->SetName(GetName());
 167   dir->Add(d);
 168   TIter    next(&fRingHistos);
 169   RingHistos* o = 0;
 170   while ((o = static_cast<RingHistos*>(next()))) {
 171     o->Output(d);
 172   }
 173 }
 174
 175 //====================================================================
 176 AliFMDCorrections::RingHistos::RingHistos()
 177   : AliForwardUtil::RingHistos(),
 178     fDensity(0)
 179 {}
 180
 181 //____________________________________________________________________
 182 AliFMDCorrections::RingHistos::RingHistos(UShort_t d, Char_t r)
 183   : AliForwardUtil::RingHistos(d,r),
 184     fDensity(0)
 185 {
 186   fDensity = new TH2D(Form("FMD%d%c_Primary_Density", d, r),
 187                       Form("in FMD%d%c", d, r),
 188                       200, -4, 6, (r == 'I' || r == 'i' ? 20 : 40),
 189                       0, 2*TMath::Pi());
 190   fDensity->SetDirectory(0);
 191   fDensity->SetXTitle("#eta");
 192   fDensity->SetYTitle("#phi [radians]");
 193   fDensity->SetZTitle("Primary N_{ch} density");
 194 }
 195 //____________________________________________________________________
 196 AliFMDCorrections::RingHistos::RingHistos(const RingHistos& o)
 197   : AliForwardUtil::RingHistos(o),
 198     fDensity(o.fDensity)
 199 {}
 200
 201 //____________________________________________________________________
 202 AliFMDCorrections::RingHistos&
 203 AliFMDCorrections::RingHistos::operator=(const RingHistos& o)
 204 {
 205   AliForwardUtil::RingHistos::operator=(o);
 206
 207   if (fDensity) delete fDensity;
 208
 209   fDensity = static_cast<TH2D*>(o.fDensity->Clone());
 210
 211   return *this;
 212 }
 213 //____________________________________________________________________
 214 AliFMDCorrections::RingHistos::~RingHistos()
 215 {
 216   if (fDensity) delete fDensity;
 217 }
 218
 219 //____________________________________________________________________
 220 void
 221 AliFMDCorrections::RingHistos::Output(TList* dir)
 222 {
 223   TList* d = DefineOutputList(dir);
 224   d->Add(fDensity);
 225 }
 226
 227 //____________________________________________________________________
 228 void
 229 AliFMDCorrections::RingHistos::ScaleHistograms(TList* dir, Int_t nEvents)
 230 {
 231   TList* l = GetOutputList(dir);
 232   if (!l) return;
 233
 234   TH1* density = GetOutputHist(l,Form("%s_Primary_Density", fName.Data()));
 235   if (density) density->Scale(1./nEvents);
 236 }
 237
 238 //____________________________________________________________________
 239 //
 240 // EOF
 241 //
 242
 243
 244