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