PWG2/FORWARD/analysis2/AliFMDMCCorrections.cxx

   1 #include "AliFMDMCCorrections.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 <TProfile2D.h>
  12 #include <iostream>
  13 #include <iomanip>
  14
  15 ClassImp(AliFMDMCCorrections)
  16 #if 0
  17 ; // For Emacs
  18 #endif
  19
  20
  21 //____________________________________________________________________
  22 AliFMDMCCorrections::~AliFMDMCCorrections()
  23 {
  24   if (fComps) fComps->Clear();
  25   if (fFMD1i) delete fFMD1i;
  26   if (fFMD2i) delete fFMD2i;
  27   if (fFMD2o) delete fFMD2o;
  28   if (fFMD3i) delete fFMD3i;
  29   if (fFMD3o) delete fFMD3o;
  30 }
  31
  32 //____________________________________________________________________
  33 AliFMDMCCorrections&
  34 AliFMDMCCorrections::operator=(const AliFMDMCCorrections& o)
  35 {
  36   AliFMDCorrections::operator=(o);
  37
  38   return *this;
  39 }
  40
  41 //____________________________________________________________________
  42 Bool_t
  43 AliFMDMCCorrections::CorrectMC(AliForwardUtil::Histos& hists,
  44                                UShort_t                vtxbin)
  45 {
  46   AliForwardCorrectionManager& fcm = AliForwardCorrectionManager::Instance();
  47
  48   UShort_t uvb = vtxbin;
  49   for (UShort_t d=1; d<=3; d++) {
  50     UShort_t nr = (d == 1 ? 1 : 2);
  51     for (UShort_t q=0; q<nr; q++) {
  52       Char_t      r  = (q == 0 ? 'I' : 'O');
  53       TH2D*       h  = hists.Get(d,r);
  54       TH2D*       bg = fcm.GetSecondaryMap()->GetCorrection(d,r,uvb);
  55       TH2D*       ef = fcm.GetVertexBias()->GetCorrection(r, uvb);
  56       if (!bg) {
  57         AliWarning(Form("No secondary correction for FMDM%d%c in vertex bin %d",
  58                         d, r, uvb));
  59         continue;
  60       }
  61       if (!ef) {
  62         AliWarning(Form("No event vertex bias correction in vertex bin %d",
  63                         uvb));
  64         continue;
  65       }
  66
  67       // Divide by primary/total ratio
  68       h->Divide(bg);
  69
  70       // Divide by the event selection efficiency
  71       h->Divide(ef);
  72     }
  73   }
  74
  75   return kTRUE;
  76 }
  77
  78 //____________________________________________________________________
  79 void
  80 AliFMDMCCorrections::Init(const TAxis& eAxis)
  81 {
  82   fFMD1i = Make(1,'I',eAxis);
  83   fFMD2i = Make(2,'I',eAxis);
  84   fFMD2o = Make(2,'O',eAxis);
  85   fFMD3i = Make(3,'I',eAxis);
  86   fFMD3o = Make(3,'O',eAxis);
  87
  88   fComps->Add(fFMD1i);
  89   fComps->Add(fFMD2i);
  90   fComps->Add(fFMD2o);
  91   fComps->Add(fFMD3i);
  92   fComps->Add(fFMD3o);
  93 }
  94
  95 //____________________________________________________________________
  96 TProfile2D*
  97 AliFMDMCCorrections::Make(UShort_t d, Char_t r,
  98                                 const TAxis& axis) const
  99 {
 100   TProfile2D* ret = new TProfile2D(Form("FMD%d%c_esd_vs_mc", d, r),
 101                                    Form("ESD/MC signal for FMD%d%c", d, r),
 102                                    axis.GetNbins(),
 103                                    axis.GetXmin(),
 104                                    axis.GetXmax(),
 105                                    (r == 'I' || r == 'i') ? 20 : 40,
 106                                    0, 2*TMath::Pi());
 107   ret->GetXaxis()->SetTitle("#eta");
 108   ret->GetYaxis()->SetTitle("#varphi [degrees]");
 109   ret->GetZaxis()->SetTitle("#LT primary density ESD/MC#GT");
 110   ret->SetDirectory(0);
 111   return ret;
 112 }
 113 //____________________________________________________________________
 114 void
 115 AliFMDMCCorrections::Fill(UShort_t d, Char_t r, TH2* esd, TH2* mc)
 116 {
 117   if (!esd || !mc) return;
 118   TProfile2D* p = 0;
 119   switch (d) {
 120   case 1:  p = fFMD1i;                                   break;
 121   case 2:  p = (r == 'I' || r == 'i' ? fFMD2i : fFMD2o); break;
 122   case 3:  p = (r == 'I' || r == 'i' ? fFMD3i : fFMD3o); break;
 123   }
 124   if (!p) return;
 125
 126   for (Int_t iEta = 1; iEta <= esd->GetNbinsX(); iEta++) {
 127     Double_t eta = esd->GetXaxis()->GetBinCenter(iEta);
 128     for (Int_t iPhi = 1; iPhi <= esd->GetNbinsY(); iPhi++) {
 129       Double_t phi  = esd->GetYaxis()->GetBinCenter(iPhi);
 130       Double_t mEsd = esd->GetBinContent(iEta,iPhi);
 131       Double_t mMc  = mc->GetBinContent(iEta,iPhi);
 132
 133       p->Fill(eta, phi, (mMc > 0 ? mEsd / mMc : 0));
 134     }
 135   }
 136 }
 137
 138 //____________________________________________________________________
 139 Bool_t
 140 AliFMDMCCorrections::CompareResults(AliForwardUtil::Histos& esd,
 141                                           AliForwardUtil::Histos& mc)
 142 {
 143   Fill(1, 'I', esd.Get(1,'I'), mc.Get(1,'I'));
 144   Fill(2, 'I', esd.Get(2,'I'), mc.Get(2,'I'));
 145   Fill(2, 'O', esd.Get(2,'O'), mc.Get(2,'O'));
 146   Fill(3, 'I', esd.Get(3,'I'), mc.Get(3,'I'));
 147   Fill(3, 'O', esd.Get(3,'O'), mc.Get(3,'O'));
 148
 149   return kTRUE;
 150 }
 151
 152 //____________________________________________________________________
 153 void
 154 AliFMDMCCorrections::DefineOutput(TList* dir)
 155 {
 156   AliFMDCorrections::DefineOutput(dir);
 157   TList* d = static_cast<TList*>(dir->FindObject(GetName()));
 158
 159   fComps = new TList;
 160   fComps->SetName("esd_mc_comparison");
 161   d->Add(fComps);
 162 }
 163
 164 //____________________________________________________________________
 165 //
 166 // EOF
 167 //
 168
 169
 170