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