PWG2/FORWARD/analysis2/AliFMDMCCorrector.cxx

   1 //
   2 // This class calculates the exclusive charged particle density
   3 // in each for the 5 FMD rings.
   4 //
   5 // Input:
   6 //   - 5 RingHistos objects - each with a number of vertex dependent
   7 //     2D histograms of the inclusive charge particle density
   8 //
   9 // Output:
  10 //   - 5 RingHistos objects - each with a number of vertex dependent
  11 //     2D histograms of the exclusive charge particle density
  12 //
  13 // Corrections used:
  14 //   - AliFMDCorrSecondaryMap;
  15 //   - AliFMDCorrVertexBias
  16 //   - AliFMDCorrMergingEfficiency
  17 //
  18 #include "AliFMDMCCorrector.h"
  19 #include <AliESDFMD.h>
  20 #include <TAxis.h>
  21 #include <TList.h>
  22 #include <TMath.h>
  23 #include "AliForwardCorrectionManager.h"
  24 #include "AliLog.h"
  25 #include <TH2D.h>
  26 #include <TROOT.h>
  27 #include <TProfile2D.h>
  28 #include <iostream>
  29 #include <iomanip>
  30
  31 ClassImp(AliFMDMCCorrector)
  32 #if 0
  33 ; // For Emacs
  34 #endif
  35
  36
  37 //____________________________________________________________________
  38 AliFMDMCCorrector::~AliFMDMCCorrector()
  39 {
  40   //
  41   // Destructor
  42   //
  43   if (fComps) fComps->Clear();
  44   if (fFMD1i) delete fFMD1i;
  45   if (fFMD2i) delete fFMD2i;
  46   if (fFMD2o) delete fFMD2o;
  47   if (fFMD3i) delete fFMD3i;
  48   if (fFMD3o) delete fFMD3o;
  49 }
  50
  51 //____________________________________________________________________
  52 AliFMDMCCorrector&
  53 AliFMDMCCorrector::operator=(const AliFMDMCCorrector& o)
  54 {
  55   //
  56   // Assignement operator
  57   //
  58   // Parameters:
  59   //    o Object to assign from
  60   //
  61   // Return:
  62   //    Reference to this object
  63   //
  64   AliFMDCorrector::operator=(o);
  65
  66   return *this;
  67 }
  68
  69 //____________________________________________________________________
  70 Bool_t
  71 AliFMDMCCorrector::CorrectMC(AliForwardUtil::Histos& hists,
  72                                UShort_t                vtxbin)
  73 {
  74   //
  75   // Do the calculations
  76   //
  77   // Parameters:
  78   //    hists    Cache of histograms
  79   //    vtxBin   Vertex bin
  80   //
  81   // Return:
  82   //    true on successs
  83   //
  84   AliForwardCorrectionManager& fcm = AliForwardCorrectionManager::Instance();
  85
  86   UShort_t uvb = vtxbin;
  87   for (UShort_t d=1; d<=3; d++) {
  88     UShort_t nr = (d == 1 ? 1 : 2);
  89     for (UShort_t q=0; q<nr; q++) {
  90       Char_t      r  = (q == 0 ? 'I' : 'O');
  91       TH2D*       h  = hists.Get(d,r);
  92       TH2D*       bg = fcm.GetSecondaryMap()->GetCorrection(d,r,uvb);
  93       TH2D*       ef = fcm.GetVertexBias()->GetCorrection(r, uvb);
  94       if (!bg) {
  95         AliWarning(Form("No secondary correction for FMDM%d%c in vertex bin %d",
  96                         d, r, uvb));
  97         continue;
  98       }
  99       if (!ef) {
 100         AliWarning(Form("No event vertex bias correction in vertex bin %d",
 101                         uvb));
 102         continue;
 103       }
 104
 105       // Divide by primary/total ratio
 106       h->Divide(bg);
 107
 108       // Divide by the event selection efficiency
 109       h->Divide(ef);
 110     }
 111   }
 112
 113   return kTRUE;
 114 }
 115
 116 //____________________________________________________________________
 117 void
 118 AliFMDMCCorrector::Init(const TAxis& eAxis)
 119 {
 120   //
 121   // Initialize this object
 122   //
 123   // Parameters:
 124   //    etaAxis Eta axis to use
 125   //
 126   fFMD1i = Make(1,'I',eAxis);
 127   fFMD2i = Make(2,'I',eAxis);
 128   fFMD2o = Make(2,'O',eAxis);
 129   fFMD3i = Make(3,'I',eAxis);
 130   fFMD3o = Make(3,'O',eAxis);
 131
 132   fComps->Add(fFMD1i);
 133   fComps->Add(fFMD2i);
 134   fComps->Add(fFMD2o);
 135   fComps->Add(fFMD3i);
 136   fComps->Add(fFMD3o);
 137 }
 138
 139 //____________________________________________________________________
 140 TProfile2D*
 141 AliFMDMCCorrector::Make(UShort_t d, Char_t r,
 142                                 const TAxis& axis) const
 143 {
 144   //
 145   // MAke comparison profiles
 146   //
 147   // Parameters:
 148   //    d     Detector
 149   //    r     Ring
 150   //    axis  Eta axis
 151   //
 152   // Return:
 153   //    Newly allocated profile object
 154   //
 155   TProfile2D* ret = new TProfile2D(Form("FMD%d%c_esd_vs_mc", d, r),
 156                                    Form("ESD/MC signal for FMD%d%c", d, r),
 157                                    axis.GetNbins(),
 158                                    axis.GetXmin(),
 159                                    axis.GetXmax(),
 160                                    (r == 'I' || r == 'i') ? 20 : 40,
 161                                    0, 2*TMath::Pi());
 162   ret->GetXaxis()->SetTitle("#eta");
 163   ret->GetYaxis()->SetTitle("#varphi [degrees]");
 164   ret->GetZaxis()->SetTitle("#LT primary density ESD/MC#GT");
 165   ret->SetDirectory(0);
 166   return ret;
 167 }
 168 //____________________________________________________________________
 169 void
 170 AliFMDMCCorrector::Fill(UShort_t d, Char_t r, TH2* esd, TH2* mc)
 171 {
 172   //
 173   // Fill comparison profiles
 174   //
 175   // Parameters:
 176   //    d    Detector
 177   //    r    Ring
 178   //    esd  ESD histogram
 179   //    mc   MC histogram
 180   //
 181   if (!esd || !mc) return;
 182   TProfile2D* p = 0;
 183   switch (d) {
 184   case 1:  p = fFMD1i;                                   break;
 185   case 2:  p = (r == 'I' || r == 'i' ? fFMD2i : fFMD2o); break;
 186   case 3:  p = (r == 'I' || r == 'i' ? fFMD3i : fFMD3o); break;
 187   }
 188   if (!p) return;
 189
 190   for (Int_t iEta = 1; iEta <= esd->GetNbinsX(); iEta++) {
 191     Double_t eta = esd->GetXaxis()->GetBinCenter(iEta);
 192     for (Int_t iPhi = 1; iPhi <= esd->GetNbinsY(); iPhi++) {
 193       Double_t phi  = esd->GetYaxis()->GetBinCenter(iPhi);
 194       Double_t mEsd = esd->GetBinContent(iEta,iPhi);
 195       Double_t mMc  = mc->GetBinContent(iEta,iPhi);
 196
 197       p->Fill(eta, phi, (mMc > 0 ? mEsd / mMc : 0));
 198     }
 199   }
 200 }
 201
 202 //____________________________________________________________________
 203 Bool_t
 204 AliFMDMCCorrector::CompareResults(AliForwardUtil::Histos& esd,
 205                                           AliForwardUtil::Histos& mc)
 206 {
 207   //
 208   // Compare the result of analysing the ESD for
 209   // the inclusive charged particle density to analysing
 210   // MC truth
 211   //
 212   // Parameters:
 213   //    esd
 214   //    mc
 215   //
 216   // Return:
 217   //   true
 218   //
 219   Fill(1, 'I', esd.Get(1,'I'), mc.Get(1,'I'));
 220   Fill(2, 'I', esd.Get(2,'I'), mc.Get(2,'I'));
 221   Fill(2, 'O', esd.Get(2,'O'), mc.Get(2,'O'));
 222   Fill(3, 'I', esd.Get(3,'I'), mc.Get(3,'I'));
 223   Fill(3, 'O', esd.Get(3,'O'), mc.Get(3,'O'));
 224
 225   return kTRUE;
 226 }
 227
 228 //____________________________________________________________________
 229 void
 230 AliFMDMCCorrector::DefineOutput(TList* dir)
 231 {
 232   //
 233   // Output diagnostic histograms to directory
 234   //
 235   // Parameters:
 236   //    dir List to write in
 237   //
 238   AliFMDCorrector::DefineOutput(dir);
 239   TList* d = static_cast<TList*>(dir->FindObject(GetName()));
 240
 241   fComps = new TList;
 242   fComps->SetName("esd_mc_comparison");
 243   d->Add(fComps);
 244 }
 245
 246 //____________________________________________________________________
 247 //
 248 // EOF
 249 //
 250
 251
 252