PWGLF/FORWARD/analysis2/AliMCTruthdNdetaTask.cxx

   1 //====================================================================
   2 #include "AliMCTruthdNdetaTask.h"
   3 #include <TMath.h>
   4 #include <TH2D.h>
   5 #include <TH1D.h>
   6 #include <THStack.h>
   7 #include <TList.h>
   8 #include <TFile.h>
   9 #include <AliAnalysisManager.h>
  10 #include <AliAODEvent.h>
  11 #include <AliAODHandler.h>
  12 #include <AliAODInputHandler.h>
  13 #include "AliForwardUtil.h"
  14 #include "AliAODForwardMult.h"
  15
  16 //____________________________________________________________________
  17 AliMCTruthdNdetaTask::AliMCTruthdNdetaTask()
  18   : AliBasedNdetaTask(),
  19     fHasData(true)
  20 {
  21   //
  22   // Constructor
  23   //
  24 }
  25
  26 //____________________________________________________________________
  27 AliMCTruthdNdetaTask::AliMCTruthdNdetaTask(const char* /* name */)
  28   : AliBasedNdetaTask("MCTruth"),
  29     fHasData(true)
  30 {
  31   //
  32   // Constructor
  33   //
  34   // Paramters
  35   //   name    Name of task
  36 }
  37
  38
  39 //____________________________________________________________________
  40 AliBasedNdetaTask::CentralityBin*
  41 AliMCTruthdNdetaTask::MakeCentralityBin(const char* name, Short_t l, Short_t h)
  42   const
  43 {
  44   //
  45   // Make a new centrality bin
  46   //
  47   // Parameters:
  48   //    name   Histogram names
  49   //    l      Lower cut
  50   //    h      Upper cut
  51   //
  52   // Return:
  53   //    Newly allocated object (of our type)
  54   //
  55   return new AliMCTruthdNdetaTask::CentralityBin(name, l, h);
  56 }
  57
  58 //____________________________________________________________________
  59 TH2D*
  60 AliMCTruthdNdetaTask::GetHistogram(const AliAODEvent& aod, Bool_t mc)
  61 {
  62   //
  63   // Retrieve the histogram
  64   //
  65   // Parameters:
  66   //    aod AOD event
  67   //    mc  Whether to get the MC histogram or not
  68   //
  69   // Return:
  70   //    Retrieved histogram or null
  71   //
  72   if (!fHasData) return 0;
  73   if (mc) return 0;
  74
  75   TH2D* ret = GetPrimary(aod);
  76   if (!ret) {
  77     fHasData = false;
  78     return 0;
  79   }
  80   Int_t nY  = ret->GetNbinsY();
  81   // Need to fill under-/overflow bin with 1's
  82   for (Int_t i = 1; i <= ret->GetNbinsX(); i++)  {
  83     ret->SetBinContent(i, 0,    1);
  84     ret->SetBinContent(i, nY+1, 1);
  85   }
  86   return ret;
  87 }
  88
  89 //________________________________________________________________________
  90 Bool_t
  91 AliMCTruthdNdetaTask::Finalize()
  92 {
  93   //
  94   // Called at end of event processing..
  95   //
  96   // This is called once in the master
  97   //
  98   // Parameters:
  99   //    option Not used
 100   if (!fHasData) {
 101     AliInfo("The MC truth dN/deta task didn't get any data");
 102     return false;
 103   }
 104   AliBasedNdetaTask::Finalize();
 105
 106   THStack* truth      = new THStack("dndetaTruth", "dN/d#eta MC Truth");
 107   THStack* truthRebin = new THStack(Form("dndetaTruth_rebin%02d", fRebin),
 108                                     "dN/d#eta MC Truth");
 109
 110   TIter next(fListOfCentralities);
 111   CentralityBin* bin = 0;
 112   while ((bin = static_cast<CentralityBin*>(next()))) {
 113     if (HasCentrality() && bin->IsAllBin()) continue;
 114
 115     TList* results = bin->GetResults();
 116     if (!results) continue;
 117
 118     TH1* dndeta      = static_cast<TH1*>(results->FindObject("dndetaTruth"));
 119     TH1* dndetaRebin =
 120       static_cast<TH1*>(results->FindObject(Form("dndetaTruth_rebin%02d",
 121                                                 fRebin)));
 122     if (dndeta)      truth->Add(dndeta);
 123     if (dndetaRebin) truthRebin->Add(dndetaRebin);
 124   }
 125   // If available output rebinned stack
 126   if (!truth->GetHists() ||
 127       truth->GetHists()->GetEntries() <= 0) {
 128     AliWarning("No MC truth histograms found");
 129     delete truth;
 130     truth = 0;
 131   }
 132   if (truth) fResults->Add(truth);
 133
 134   // If available output rebinned stack
 135   if (!truthRebin->GetHists() ||
 136       truthRebin->GetHists()->GetEntries() <= 0) {
 137     AliWarning("No rebinned MC truth histograms found");
 138     delete truthRebin;
 139     truthRebin = 0;
 140   }
 141   if (truthRebin) fResults->Add(truthRebin);
 142   return true;
 143 }
 144
 145 //========================================================================
 146 Bool_t
 147 AliMCTruthdNdetaTask::CentralityBin::ProcessEvent(const AliAODForwardMult*
 148                                                   forward,
 149                                                   Int_t triggerMask,
 150                                                   Bool_t isZero,
 151                                                   Double_t vzMin,
 152                                                   Double_t vzMax,
 153                                                   const TH2D* primary,
 154                                                   const TH2D*)
 155 {
 156   // Check the centrality class unless this is the 'all' bin
 157   if (!primary) return false;
 158
 159   if (!IsAllBin()) {
 160     Double_t centrality = forward->GetCentrality();
 161     if (centrality < fLow || centrality >= fHigh) return false;
 162   }
 163
 164   if (!fSum) CreateSums(primary, 0);
 165   if (!fSumTruth) {
 166     fSumTruth = static_cast<TH2D*>(primary->Clone("truth"));
 167     fSumTruth->SetDirectory(0);
 168     fSumTruth->Reset();
 169     fSums->Add(fSumTruth);
 170   }
 171
 172   // translate real trigger mask to MC trigger mask
 173   Int_t mask = AliAODForwardMult::kB;
 174   if (triggerMask == AliAODForwardMult::kNSD) {
 175     mask ^= AliAODForwardMult::kNSD;
 176     mask =  AliAODForwardMult::kMCNSD;
 177   }
 178
 179   // Now use our normal check, but with the new mask, except ignore vertex
 180   if (forward->CheckEvent(mask, -10000, -10000, 0, 0, 0)) {
 181     fSumTruth->Add(primary);
 182
 183     // Store event count in left-most under- underflow bin
 184     Int_t cnt = Int_t(fSumTruth->GetBinContent(0,0));
 185     fSumTruth->SetBinContent(0,0, ++cnt);
 186   }
 187
 188   // Now use our normal check with the full trigger mask and vertex
 189   if (CheckEvent(forward, triggerMask, vzMin, vzMax))
 190     fSum->Add(primary, isZero);
 191   return true;
 192 }
 193
 194 //________________________________________________________________________
 195 void
 196 AliMCTruthdNdetaTask::CentralityBin::End(TList*      sums,
 197                                          TList*      results,
 198                                          UShort_t    scheme,
 199                                          const TH2F* shapeCorr,
 200                                          Double_t    trigEff,
 201                                          Double_t    trigEff0,
 202                                          Bool_t      symmetrice,
 203                                          Int_t       rebin,
 204                                          Bool_t      rootProj,
 205                                          Bool_t      corrEmpty,
 206                                          Bool_t      cutEdges,
 207                                          Int_t       triggerMask,
 208                                          Int_t       marker,
 209                                          Int_t       color,
 210                                          TList*      mclist,
 211                                          TList*      truthlist)
 212 {
 213 #if 0
 214   AliInfo(Form("At end with sums=%p, results=%p, scheme=%d, "
 215                "shapeCorr=%p, trigEff=%f, symmetrice=%d, rebin=%d, "
 216                "rootProj=%d, corrEmpty=%d, cutEdges=%d, triggerMask=0x%08x, "
 217                "marker=%d (%d)",
 218                sums, results, scheme, shapeCorr, trigEff, trigEff0, symmetrice,
 219                rebin, rootProj, corrEmpty, cutEdges, triggerMask, marker,
 220                GetMarkerStyle(kStar)));
 221 #endif
 222
 223   AliBasedNdetaTask::CentralityBin::End(sums, results, scheme,
 224                                         shapeCorr, trigEff, trigEff0,
 225                                         symmetrice, rebin,
 226                                         rootProj, corrEmpty, cutEdges,
 227                                         triggerMask, marker, color, mclist,
 228                                         truthlist);
 229
 230   fSumTruth     = static_cast<TH2D*>(fSums->FindObject("truth"));
 231
 232
 233   if (fSumTruth) {
 234     Int_t n0 = Int_t(fSumTruth->GetBinContent(0,0));
 235     Int_t n  = (triggerMask == AliAODForwardMult::kNSD ?
 236                 Int_t(fTriggers->GetBinContent(AliAODForwardMult::kBinMCNSD)) :
 237                 Int_t(fTriggers->GetBinContent(AliAODForwardMult::kBinAll)));
 238     DMSG(fDebug,0,"Normalising MC truth to %d (%d additions)", n, n0);
 239
 240     TH1D* dndetaTruth = fSumTruth->ProjectionX("dndetaTruth",1,
 241                                                fSumTruth->GetNbinsY(),"e");
 242     dndetaTruth->SetDirectory(0);
 243     dndetaTruth->Scale(1./n, "width");
 244
 245     SetHistogramAttributes(dndetaTruth, GetColor(color)+1,
 246                            GetMarkerStyle(kCross),
 247                            "Monte-Carlo truth");
 248
 249     fOutput->Add(dndetaTruth);
 250     fOutput->Add(Rebin(dndetaTruth, rebin, cutEdges));
 251   }
 252   TH1* dndeta          =                       GetResult(0,     false, "");
 253   TH1* dndetaSym       = symmetrice          ? GetResult(0,     true,  "") : 0;
 254   TH1* dndeta_rebin    = rebin               ? GetResult(rebin, false, "") : 0;
 255   TH1* dndetaSym_rebin = symmetrice && rebin ? GetResult(rebin, true,  "") : 0;
 256   if (dndeta)
 257     dndeta->SetTitle("Monte-Carlo truth (selected)");
 258   if (dndetaSym)
 259     dndetaSym->SetTitle("Monte-Carlo truth (selected,mirrored)");
 260   if (dndeta_rebin)
 261     dndeta_rebin->SetTitle("Monte-Carlo truth (selected)");
 262   if (dndetaSym_rebin)
 263     dndetaSym_rebin->SetTitle("Monte-Carlo truth (selected,mirrored)");
 264
 265 }
 266
 267 //________________________________________________________________________
 268 //
 269 // EOF
 270 //