PWGLF/FORWARD/analysis2/AliForwardMultiplicityTask.cxx

   1 //
   2 // Calculate the multiplicity in the forward regions event-by-event
   3 //
   4 // Inputs:
   5 //   - AliESDEvent
   6 //
   7 // Outputs:
   8 //   - AliAODForwardMult
   9 //
  10 // Histograms
  11 //
  12 // Corrections used
  13 //
  14 #include "AliForwardMultiplicityTask.h"
  15 #include "AliTriggerAnalysis.h"
  16 #include "AliPhysicsSelection.h"
  17 #include "AliLog.h"
  18 #include "AliESDEvent.h"
  19 #include "AliAODHandler.h"
  20 #include "AliMultiplicity.h"
  21 #include "AliInputEventHandler.h"
  22 #include "AliForwardCorrectionManager.h"
  23 #include "AliAnalysisManager.h"
  24 #include <TH1.h>
  25 #include <TH3D.h>
  26 #include <TDirectory.h>
  27 #include <TTree.h>
  28 #include <TROOT.h>
  29 #include <TStopwatch.h>
  30 #include <TProfile.h>
  31 // #define ENABLE_TIMING
  32 #ifndef ENABLE_TIMING
  33 # define MAKE_SW(NAME) do {} while(false)
  34 # define START_SW(NAME) do {} while(false)
  35 # define FILL_SW(NAME,WHICH) do {} while(false)
  36 #else
  37 # define MAKE_SW(NAME) TStopwatch NAME
  38 # define START_SW(NAME) if (fDoTiming) NAME.Start(true)
  39 # define FILL_SW(NAME,WHICH)                            \
  40   if (fDoTiming) fHTiming->Fill(WHICH,NAME.CpuTime())
  41 #endif
  42
  43 //====================================================================
  44 AliForwardMultiplicityTask::AliForwardMultiplicityTask()
  45   : AliForwardMultiplicityBase(),
  46     fESDFMD(),
  47     fEventInspector(),
  48     fESDFixer(),
  49     fSharingFilter(),
  50     fDensityCalculator(),
  51     fCorrections(),
  52     fHistCollector(),
  53     fEventPlaneFinder()
  54 {
  55   //
  56   // Constructor
  57   //
  58   DGUARD(fDebug, 3,"Default CTOR of AliForwardMultiplicityTask");
  59 }
  60
  61 //____________________________________________________________________
  62 AliForwardMultiplicityTask::AliForwardMultiplicityTask(const char* name)
  63   : AliForwardMultiplicityBase(name),
  64     fESDFMD(),
  65     fEventInspector("event"),
  66     fESDFixer("esdFizer"),
  67     fSharingFilter("sharing"),
  68     fDensityCalculator("density"),
  69     fCorrections("corrections"),
  70     fHistCollector("collector"),
  71     fEventPlaneFinder("eventplane")
  72 {
  73   //
  74   // Constructor
  75   //
  76   // Parameters:
  77   //    name Name of task
  78   //
  79   DGUARD(fDebug, 3,"named CTOR of AliForwardMultiplicityTask: %s", name);
  80 }
  81
  82
  83 //____________________________________________________________________
  84 void
  85 AliForwardMultiplicityTask::SetDoTiming(Bool_t enable)
  86 {
  87 #ifndef ENABLE_TIMING
  88   if (enable)
  89     AliWarning("Timing of task explicitly disabled in compilation");
  90 #else
  91   fDoTiming = enable;
  92 #endif
  93 }
  94
  95 //____________________________________________________________________
  96 void
  97 AliForwardMultiplicityTask::PreCorrections(const AliESDEvent* esd)
  98 {
  99   if (!esd) return;
 100
 101   AliESDFMD* esdFMD = esd->GetFMDData();
 102   if (!esdFMD) return;
 103
 104   Int_t tgt = GetESDFixer().FindTargetNoiseFactor(*esdFMD, false);
 105   if (tgt <= 0) {
 106     // If the target noise factor is 0 or less, disable the noise/gain
 107     // correction.
 108     GetESDFixer().SetRecoNoiseFactor(4);
 109     fNeededCorrections ^= AliForwardCorrectionManager::kNoiseGain;
 110   }
 111   else
 112     AliWarning("The noise corrector has been enabled!");
 113 }
 114 //____________________________________________________________________
 115 Bool_t
 116 AliForwardMultiplicityTask::PreEvent()
 117 {
 118   // Clear stuff
 119   fHistos.Clear();
 120   fESDFMD.Clear();
 121   fAODFMD.Clear();
 122   fAODEP.Clear();
 123   return true;
 124 }
 125 //____________________________________________________________________
 126 Bool_t
 127 AliForwardMultiplicityTask::Event(AliESDEvent& esd)
 128 {
 129   //
 130   // Process each event
 131   //
 132   // Parameters:
 133   //    option Not used
 134   //
 135   MAKE_SW(total);
 136   MAKE_SW(individual);
 137   START_SW(total);
 138
 139   DGUARD(fDebug,1,"Process the input event");
 140
 141   // Inspect the event
 142   START_SW(individual);
 143   Bool_t   lowFlux   = kFALSE;
 144   UInt_t   triggers  = 0;
 145   UShort_t ivz       = 0;
 146   TVector3 ip;
 147   Double_t cent      = -1;
 148   UShort_t nClusters = 0;
 149   UInt_t   found     = fEventInspector.Process(&esd, triggers, lowFlux,
 150                                                ivz, ip, cent, nClusters);
 151   FILL_SW(individual,kTimingEventInspector);
 152
 153   if (found & AliFMDEventInspector::kNoEvent)    {
 154     fHStatus->Fill(1);
 155     return false;
 156   }
 157   if (found & AliFMDEventInspector::kNoTriggers) {
 158     fHStatus->Fill(2);
 159     return false;
 160   }
 161
 162   // Set trigger bits, and mark this event for storage
 163   fAODFMD.SetTriggerBits(triggers);
 164   fAODFMD.SetSNN(fEventInspector.GetEnergy());
 165   fAODFMD.SetSystem(fEventInspector.GetCollisionSystem());
 166   fAODFMD.SetCentrality(cent);
 167   fAODFMD.SetNClusters(nClusters);
 168   MarkEventForStore();
 169
 170   // Do not check if SPD data is there - potential bias
 171   // if (found & AliFMDEventInspector::kNoSPD) {
 172   //   fHStatus->Fill(3);
 173   //   return false;
 174   // }
 175   if (found    & AliFMDEventInspector::kNoFMD) {
 176     fHStatus->Fill(4);
 177     return false;
 178   }
 179   if (found    & AliFMDEventInspector::kNoVertex) {
 180     fHStatus->Fill(5);
 181     return false;
 182   }
 183   // Also analyse pile-up events - we'll remove them in later steps.
 184   if (triggers & AliAODForwardMult::kPileUp) {
 185     fHStatus->Fill(6);
 186     return false;
 187   }
 188   fAODFMD.SetIpZ(ip.Z());
 189   if (found & AliFMDEventInspector::kBadVertex) {
 190     fHStatus->Fill(7);
 191     return false;
 192   }
 193
 194   // We we do not want to use low flux specific code, we disable it here.
 195   if (!fEnableLowFlux) lowFlux = false;
 196
 197   // Get FMD data
 198   AliESDFMD* esdFMD = esd.GetFMDData();
 199
 200   // Fix up the the ESD
 201   GetESDFixer().Fix(*esdFMD, ip.Z());
 202
 203   // Apply the sharing filter (or hit merging or clustering if you like)
 204   START_SW(individual);
 205   if (!fSharingFilter.Filter(*esdFMD, lowFlux, fESDFMD, ip.Z())) {
 206     AliWarning("Sharing filter failed!");
 207     fHStatus->Fill(8);
 208     return false;
 209   }
 210   FILL_SW(individual,kTimingSharingFilter);
 211
 212   // Calculate the inclusive charged particle density
 213   START_SW(individual);
 214   if (!fDensityCalculator.Calculate(fESDFMD, fHistos, lowFlux, cent, ip)) {
 215     // if (!fDensityCalculator.Calculate(*esdFMD, fHistos, ivz, lowFlux)) {
 216     AliWarning("Density calculator failed!");
 217     fHStatus->Fill(9);
 218     return false;
 219   }
 220   FILL_SW(individual,kTimingDensityCalculator);
 221
 222   // Check if we should do the event plane finder
 223   if (fEventInspector.GetCollisionSystem() == AliFMDEventInspector::kPbPb) {
 224     START_SW(individual);
 225     if (!fEventPlaneFinder.FindEventplane(&esd, fAODEP,
 226                                           &(fAODFMD.GetHistogram()), &fHistos)){
 227       AliWarning("Eventplane finder failed!");
 228       fHStatus->Fill(10);
 229     }
 230     FILL_SW(individual,kTimingEventPlaneFinder);
 231   }
 232
 233   // Check how many rings have been marked for skipping
 234   Int_t nSkip = 0;
 235   for (UShort_t d=1; d<=3; d++) {
 236     for (UShort_t q=0; q<=(d/2); q++) {
 237       TH2D* h = fHistos.Get(d,q == 0 ? 'I' : 'O');
 238       if (h && h->TestBit(AliForwardUtil::kSkipRing)) nSkip++;
 239     }
 240   }
 241   if (nSkip > 0) {
 242     // Skip the rest if we have too many outliers
 243     fHStatus->Fill(11);
 244     return false;
 245   }
 246
 247   // Do the secondary and other corrections.
 248   START_SW(individual);
 249   if (!fCorrections.Correct(fHistos, ivz)) {
 250     AliWarning("Corrections failed");
 251     fHStatus->Fill(12);
 252     return false;
 253   }
 254   FILL_SW(individual,kTimingCorrections);
 255
 256   // Check if we should add to internal caches
 257   Bool_t add = (fAODFMD.IsTriggerBits(AliAODForwardMult::kInel) &&
 258                 !(triggers & AliAODForwardMult::kPileUp) && nSkip < 1);
 259
 260   // Collect our `super' histogram
 261   START_SW(individual);
 262   if (!fHistCollector.Collect(fHistos,
 263                               fRingSums,
 264                               ivz,
 265                               fAODFMD.GetHistogram(),
 266                               fAODFMD.GetCentrality(),
 267                               false,
 268                               add)) {
 269     AliWarning("Histogram collector failed");
 270     fHStatus->Fill(13);
 271     return false;
 272   }
 273   FILL_SW(individual,kTimingHistCollector);
 274
 275   if (!add) {
 276     fHStatus->Fill(14);
 277   }
 278   else {
 279     // Collect rough Min. Bias result
 280     fHData->Add(&(fAODFMD.GetHistogram()));
 281     fHStatus->Fill(15);
 282   }
 283   FILL_SW(total,kTimingTotal);
 284
 285   return true;
 286 }
 287
 288
 289 //
 290 // EOF
 291 //