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 "AliForwardCorrectionManager.h"
   7 // #include "AliFMDAnaParameters.h"
   8 #include "AliLog.h"
   9 #include <TH2D.h>
  10 #include <TROOT.h>
  11 #include <iostream>
  12 #include <iomanip>
  13
  14 ClassImp(AliFMDCorrections)
  15 #if 0
  16 ; // For Emacs
  17 #endif
  18
  19 //____________________________________________________________________
  20 AliFMDCorrections::AliFMDCorrections()
  21   : TNamed(),
  22     fRingHistos(),
  23     fDebug(0)
  24 {}
  25
  26 //____________________________________________________________________
  27 AliFMDCorrections::AliFMDCorrections(const char* title)
  28   : TNamed("fmdCorrections", title),
  29     fRingHistos(),
  30     fDebug(0)
  31 {
  32   fRingHistos.SetName(GetName());
  33   fRingHistos.Add(new RingHistos(1, 'I'));
  34   fRingHistos.Add(new RingHistos(2, 'I'));
  35   fRingHistos.Add(new RingHistos(2, 'O'));
  36   fRingHistos.Add(new RingHistos(3, 'I'));
  37   fRingHistos.Add(new RingHistos(3, 'O'));
  38 }
  39
  40 //____________________________________________________________________
  41 AliFMDCorrections::AliFMDCorrections(const AliFMDCorrections& o)
  42   : TNamed(o),
  43     fRingHistos(),
  44     fDebug(o.fDebug)
  45 {
  46   TIter    next(&o.fRingHistos);
  47   TObject* obj = 0;
  48   while ((obj = next())) fRingHistos.Add(obj);
  49 }
  50
  51 //____________________________________________________________________
  52 AliFMDCorrections::~AliFMDCorrections()
  53 {
  54   fRingHistos.Delete();
  55 }
  56
  57 //____________________________________________________________________
  58 AliFMDCorrections&
  59 AliFMDCorrections::operator=(const AliFMDCorrections& o)
  60 {
  61   TNamed::operator=(o);
  62
  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                            UShort_t                vtxbin)
  91 {
  92   AliForwardCorrectionManager& fcm = AliForwardCorrectionManager::Instance();
  93
  94   UShort_t uvb = vtxbin;
  95   for (UShort_t d=1; d<=3; d++) {
  96     UShort_t nr = (d == 1 ? 1 : 2);
  97     for (UShort_t q=0; q<nr; q++) {
  98       Char_t      r = (q == 0 ? 'I' : 'O');
  99       TH2D*       h = hists.Get(d,r);
 100       RingHistos* rh= GetRingHistos(d,r);
 101       TH2D* bg = fcm.GetSecondaryMap()->GetCorrection(d,r,uvb);
 102       TH2D* ef = fcm.GetVertexBias()->GetCorrection(r, uvb);
 103       if (!bg) {
 104         AliWarning(Form("No secondary correction for FMDM%d%c in vertex bin %d",
 105                         d, r, uvb));
 106         continue;
 107       }
 108       if (!ef) {
 109         AliWarning(Form("No event vertex bias correction in vertex bin %d",
 110                         uvb));
 111         continue;
 112       }
 113
 114       // Divide by primary/total ratio
 115       h->Divide(bg);
 116
 117       // Divide by the event selection efficiency
 118       h->Divide(ef);
 119
 120
 121       // if(!pars->SharingEffPresent()) {
 122       //   AliWarning("No sharing efficiencies");
 123       //   continue;
 124       // }
 125       // TH1F* sf = pars->GetSharingEfficiencyTrVtx(d,r,vtxbin);
 126       if (!fcm.GetMergingEfficiency()) {
 127         AliWarning("No merging efficiencies");
 128         continue;
 129       }
 130       TH1D* sf = fcm.GetMergingEfficiency()->GetCorrection(d,r,uvb);
 131       if (!fcm.GetMergingEfficiency()->GetCorrection(d,r,uvb)) {
 132         AliWarning(Form("No merging efficiency for FMD%d%c at vertex bin %d",
 133                         d, r, uvb));
 134         continue;
 135       }
 136
 137
 138       for (Int_t ieta = 1; ieta <= h->GetNbinsX(); ieta++) {
 139         Float_t c  = sf->GetBinContent(ieta);
 140         Float_t ec = sf->GetBinError(ieta);
 141
 142         if (c == 0) continue;
 143
 144         for (Int_t iphi = 1; iphi <= h->GetNbinsY(); iphi++) {
 145           Double_t m  = h->GetBinContent(ieta, iphi) / c;
 146           Double_t em = h->GetBinError(ieta, iphi);
 147
 148           Double_t e  = TMath::Sqrt(em * em + (m * ec) * (m * ec));
 149
 150           h->SetBinContent(ieta,iphi,m);
 151           h->SetBinError(ieta,iphi,e);
 152         }
 153       }
 154
 155       rh->fDensity->Add(h);
 156     }
 157   }
 158
 159   return kTRUE;
 160 }
 161
 162 //____________________________________________________________________
 163 void
 164 AliFMDCorrections::ScaleHistograms(TList* dir, Int_t nEvents)
 165 {
 166   if (nEvents <= 0) return;
 167   TList* d = static_cast<TList*>(dir->FindObject(GetName()));
 168   if (!d) return;
 169
 170   TIter    next(&fRingHistos);
 171   RingHistos* o = 0;
 172   while ((o = static_cast<RingHistos*>(next())))
 173     o->ScaleHistograms(d, nEvents);
 174 }
 175 //____________________________________________________________________
 176 void
 177 AliFMDCorrections::DefineOutput(TList* dir)
 178 {
 179   TList* d = new TList;
 180   d->SetName(GetName());
 181   dir->Add(d);
 182   TIter    next(&fRingHistos);
 183   RingHistos* o = 0;
 184   while ((o = static_cast<RingHistos*>(next()))) {
 185     o->Output(d);
 186   }
 187 }
 188
 189 //____________________________________________________________________
 190 void
 191 AliFMDCorrections::Print(Option_t* /* option */) const
 192 {
 193   char ind[gROOT->GetDirLevel()+1];
 194   for (Int_t i = 0; i < gROOT->GetDirLevel(); i++) ind[i] = ' ';
 195   ind[gROOT->GetDirLevel()] = '\0';
 196   std::cout << ind << "AliFMDCorrections: " << GetName() <<  std::endl;
 197 }
 198
 199 //====================================================================
 200 AliFMDCorrections::RingHistos::RingHistos()
 201   : AliForwardUtil::RingHistos(),
 202     fDensity(0)
 203 {}
 204
 205 //____________________________________________________________________
 206 AliFMDCorrections::RingHistos::RingHistos(UShort_t d, Char_t r)
 207   : AliForwardUtil::RingHistos(d,r),
 208     fDensity(0)
 209 {
 210   fDensity = new TH2D(Form("FMD%d%c_Primary_Density", d, r),
 211                       Form("in FMD%d%c", d, r),
 212                       200, -4, 6, (r == 'I' || r == 'i' ? 20 : 40),
 213                       0, 2*TMath::Pi());
 214   fDensity->SetDirectory(0);
 215   fDensity->SetXTitle("#eta");
 216   fDensity->SetYTitle("#phi [radians]");
 217   fDensity->SetZTitle("Primary N_{ch} density");
 218 }
 219 //____________________________________________________________________
 220 AliFMDCorrections::RingHistos::RingHistos(const RingHistos& o)
 221   : AliForwardUtil::RingHistos(o),
 222     fDensity(o.fDensity)
 223 {}
 224
 225 //____________________________________________________________________
 226 AliFMDCorrections::RingHistos&
 227 AliFMDCorrections::RingHistos::operator=(const RingHistos& o)
 228 {
 229   AliForwardUtil::RingHistos::operator=(o);
 230
 231   if (fDensity) delete fDensity;
 232
 233   fDensity = static_cast<TH2D*>(o.fDensity->Clone());
 234
 235   return *this;
 236 }
 237 //____________________________________________________________________
 238 AliFMDCorrections::RingHistos::~RingHistos()
 239 {
 240   if (fDensity) delete fDensity;
 241 }
 242
 243 //____________________________________________________________________
 244 void
 245 AliFMDCorrections::RingHistos::Output(TList* dir)
 246 {
 247   TList* d = DefineOutputList(dir);
 248   d->Add(fDensity);
 249 }
 250
 251 //____________________________________________________________________
 252 void
 253 AliFMDCorrections::RingHistos::ScaleHistograms(TList* dir, Int_t nEvents)
 254 {
 255   TList* l = GetOutputList(dir);
 256   if (!l) return;
 257
 258   TH1* density = GetOutputHist(l,Form("%s_Primary_Density", fName.Data()));
 259   if (density) density->Scale(1./nEvents);
 260 }
 261
 262 //____________________________________________________________________
 263 //
 264 // EOF
 265 //
 266
 267
 268