2 // Calculate the multiplicity in the forward regions event-by-event
10 // - AliAODForwardMult
16 #include "AliForwardMCMultiplicityTask.h"
17 #include "AliTriggerAnalysis.h"
18 #include "AliPhysicsSelection.h"
20 #include "AliESDEvent.h"
21 #include "AliAODHandler.h"
22 #include "AliMultiplicity.h"
23 #include "AliInputEventHandler.h"
24 #include "AliForwardCorrectionManager.h"
25 #include "AliAnalysisManager.h"
27 #include <TDirectory.h>
31 //====================================================================
32 AliForwardMCMultiplicityTask::AliForwardMCMultiplicityTask()
33 : AliForwardMultiplicityBase(),
56 //____________________________________________________________________
57 AliForwardMCMultiplicityTask::AliForwardMCMultiplicityTask(const char* name)
58 : AliForwardMultiplicityBase(name),
69 fEventInspector("event"),
70 fSharingFilter("sharing"),
71 fDensityCalculator("density"),
72 fCorrections("corrections"),
73 fHistCollector("collector"),
82 DefineOutput(1, TList::Class());
85 //____________________________________________________________________
86 AliForwardMCMultiplicityTask::AliForwardMCMultiplicityTask(const AliForwardMCMultiplicityTask& o)
87 : AliForwardMultiplicityBase(o),
92 fMCESDFMD(o.fMCESDFMD),
93 fMCHistos(o.fMCHistos),
94 fMCAODFMD(o.fMCAODFMD),
95 fRingSums(o.fRingSums),
96 fMCRingSums(o.fMCRingSums),
98 fEventInspector(o.fEventInspector),
99 fSharingFilter(o.fSharingFilter),
100 fDensityCalculator(o.fDensityCalculator),
101 fCorrections(o.fCorrections),
102 fHistCollector(o.fHistCollector),
109 // o Object to copy from
111 DefineOutput(1, TList::Class());
114 //____________________________________________________________________
115 AliForwardMCMultiplicityTask&
116 AliForwardMCMultiplicityTask::operator=(const AliForwardMCMultiplicityTask& o)
119 // Assignment operator
122 // o Object to assign from
125 // Reference to this object
127 AliForwardMultiplicityBase::operator=(o);
130 fEventInspector = o.fEventInspector;
131 fSharingFilter = o.fSharingFilter;
132 fDensityCalculator = o.fDensityCalculator;
133 fCorrections = o.fCorrections;
134 fHistCollector = o.fHistCollector;
137 fMCHistos = o.fMCHistos;
138 fMCAODFMD = o.fMCAODFMD;
139 fRingSums = o.fRingSums;
140 fMCRingSums = o.fMCRingSums;
141 fPrimary = o.fPrimary;
147 //____________________________________________________________________
149 AliForwardMCMultiplicityTask::SetDebug(Int_t dbg)
157 fEventInspector.SetDebug(dbg);
158 fSharingFilter.SetDebug(dbg);
159 fDensityCalculator.SetDebug(dbg);
160 fCorrections.SetDebug(dbg);
161 fHistCollector.SetDebug(dbg);
163 //____________________________________________________________________
165 AliForwardMCMultiplicityTask::SetOnlyPrimary(Bool_t use)
167 fSharingFilter.GetTrackDensity().SetUseOnlyPrimary(use);
168 fCorrections.SetSecondaryForMC(!use);
171 //____________________________________________________________________
173 AliForwardMCMultiplicityTask::InitializeSubs()
176 // Initialise the sub objects and stuff. Called on first event
182 if (!ReadCorrections(pe,pv,true)) return;
189 fMCRingSums.Init(*pe);
191 fHData = static_cast<TH2D*>(fAODFMD.GetHistogram().Clone("d2Ndetadphi"));
193 fHData->SetDirectory(0);
197 TList* rings = new TList;
198 rings->SetName("ringSums");
202 rings->Add(fRingSums.Get(1, 'I'));
203 rings->Add(fRingSums.Get(2, 'I'));
204 rings->Add(fRingSums.Get(2, 'O'));
205 rings->Add(fRingSums.Get(3, 'I'));
206 rings->Add(fRingSums.Get(3, 'O'));
207 fRingSums.Get(1, 'I')->SetMarkerColor(AliForwardUtil::RingColor(1, 'I'));
208 fRingSums.Get(2, 'I')->SetMarkerColor(AliForwardUtil::RingColor(2, 'I'));
209 fRingSums.Get(2, 'O')->SetMarkerColor(AliForwardUtil::RingColor(2, 'O'));
210 fRingSums.Get(3, 'I')->SetMarkerColor(AliForwardUtil::RingColor(3, 'I'));
211 fRingSums.Get(3, 'O')->SetMarkerColor(AliForwardUtil::RingColor(3, 'O'));
213 TList* mcRings = new TList;
214 mcRings->SetName("mcRingSums");
218 mcRings->Add(fMCRingSums.Get(1, 'I'));
219 mcRings->Add(fMCRingSums.Get(2, 'I'));
220 mcRings->Add(fMCRingSums.Get(2, 'O'));
221 mcRings->Add(fMCRingSums.Get(3, 'I'));
222 mcRings->Add(fMCRingSums.Get(3, 'O'));
223 fMCRingSums.Get(1, 'I')->SetMarkerColor(AliForwardUtil::RingColor(1, 'I'));
224 fMCRingSums.Get(2, 'I')->SetMarkerColor(AliForwardUtil::RingColor(2, 'I'));
225 fMCRingSums.Get(2, 'O')->SetMarkerColor(AliForwardUtil::RingColor(2, 'O'));
226 fMCRingSums.Get(3, 'I')->SetMarkerColor(AliForwardUtil::RingColor(3, 'I'));
227 fMCRingSums.Get(3, 'O')->SetMarkerColor(AliForwardUtil::RingColor(3, 'O'));
231 fEventInspector.Init(*pv);
232 fSharingFilter.Init();
233 fDensityCalculator.Init(*pe);
234 fCorrections.Init(*pe);
235 fHistCollector.Init(*pv,*pe);
240 //____________________________________________________________________
242 AliForwardMCMultiplicityTask::UserCreateOutputObjects()
245 // Create output objects
251 AliAnalysisManager* am = AliAnalysisManager::GetAnalysisManager();
253 dynamic_cast<AliAODHandler*>(am->GetOutputEventHandler());
254 if (!ah) AliFatal("No AOD output handler set in analysis manager");
257 TObject* obj = &fAODFMD;
258 ah->AddBranch("AliAODForwardMult", &obj);
260 TObject* mcobj = &fMCAODFMD;
261 ah->AddBranch("AliAODForwardMult", &mcobj);
263 fPrimary = new TH2D("primary", "MC Primaries",
264 200, -4, 6, 20, 0, 2*TMath::Pi());
265 fPrimary->SetXTitle("#eta");
266 fPrimary->SetYTitle("#varphi [radians]");
267 fPrimary->SetZTitle("d^{2}N_{ch}/d#etad#phi");
269 fPrimary->SetStats(0);
270 fPrimary->SetDirectory(0);
271 ah->AddBranch("TH2D", &fPrimary);
273 fEventInspector.DefineOutput(fList);
274 fSharingFilter.DefineOutput(fList);
275 fDensityCalculator.DefineOutput(fList);
276 fCorrections.DefineOutput(fList);
277 fHistCollector.DefineOutput(fList);
281 //____________________________________________________________________
283 AliForwardMCMultiplicityTask::UserExec(Option_t*)
286 // Process each event
292 // Get the input data
293 AliESDEvent* esd = GetESDEvent();
294 AliMCEvent* mcEvent = MCEvent();
305 Bool_t lowFlux = kFALSE;
310 UShort_t nClusters = 0;
311 UInt_t found = fEventInspector.Process(esd, triggers, lowFlux,
312 ivz, vz, cent, nClusters);
320 fEventInspector.ProcessMC(mcEvent, triggers, ivzMC, vzMC, b,
322 fEventInspector.CompareResults(vz, vzMC, cent, b, npart, nbin);
327 Bool_t isAccepted = true;
328 if (found & AliFMDEventInspector::kNoEvent) isAccepted = false; // return;
329 if (found & AliFMDEventInspector::kNoTriggers) isAccepted = false; // return;
330 //MarkEventForStore();
331 // Always set the B trigger - each MC event _is_ a collision
332 triggers |= AliAODForwardMult::kB;
333 // Set trigger bits, and mark this event for storage
334 fAODFMD.SetTriggerBits(triggers);
335 fAODFMD.SetSNN(fEventInspector.GetEnergy());
336 fAODFMD.SetSystem(fEventInspector.GetCollisionSystem());
337 fAODFMD.SetCentrality(cent);
338 fAODFMD.SetNClusters(nClusters);
340 fMCAODFMD.SetTriggerBits(triggers);
341 fMCAODFMD.SetSNN(fEventInspector.GetEnergy());
342 fMCAODFMD.SetSystem(fEventInspector.GetCollisionSystem());
343 fMCAODFMD.SetCentrality(cent);
344 fMCAODFMD.SetNClusters(nClusters);
346 //All events should be stored - HHD
347 //if (isAccepted) MarkEventForStore();
349 // Disable this check on SPD - will bias data
350 // if (found & AliFMDEventInspector::kNoSPD) isAccepted = false; // return;
351 if (found & AliFMDEventInspector::kNoFMD) isAccepted = false; // return;
352 if (found & AliFMDEventInspector::kNoVertex) isAccepted = false; // return;
356 fMCAODFMD.SetIpZ(vz);
358 if (found & AliFMDEventInspector::kBadVertex) isAccepted = false; // return;
360 // We we do not want to use low flux specific code, we disable it here.
361 if (!fEnableLowFlux) lowFlux = false;
366 AliESDFMD* esdFMD = esd->GetFMDData();
368 // Apply the sharing filter (or hit merging or clustering if you like)
369 if (isAccepted && !fSharingFilter.Filter(*esdFMD, lowFlux, fESDFMD)) {
370 AliWarning("Sharing filter failed!");
373 if (!fSharingFilter.FilterMC(*esdFMD, *mcEvent, vz, fMCESDFMD, fPrimary)) {
374 AliWarning("MC Sharing filter failed!");
377 if (!isAccepted) return; // Exit on MC event w/o trigger, vertex, data
378 // HHD if (!isAccepted) return; // Exit on MC event w/o trigger, vertex, data
380 //MarkEventForStore();
381 fSharingFilter.CompareResults(fESDFMD, fMCESDFMD);
383 // Calculate the inclusive charged particle density
384 if (!fDensityCalculator.Calculate(fESDFMD, fHistos, ivz, lowFlux)) {
385 AliWarning("Density calculator failed!");
388 if (!fDensityCalculator.CalculateMC(fMCESDFMD, fMCHistos)) {
389 AliWarning("MC Density calculator failed!");
392 fDensityCalculator.CompareResults(fHistos, fMCHistos);
394 // Do the secondary and other corrections.
395 if (!fCorrections.Correct(fHistos, ivz)) {
396 AliWarning("Corrections failed");
399 if (!fCorrections.CorrectMC(fMCHistos, ivz)) {
400 AliWarning("MC Corrections failed");
403 fCorrections.CompareResults(fHistos, fMCHistos);
405 if (!fHistCollector.Collect(fHistos, fRingSums,
406 ivz, fAODFMD.GetHistogram())) {
407 AliWarning("Histogram collector failed");
410 if (!fHistCollector.Collect(fMCHistos, fMCRingSums,
411 ivz, fMCAODFMD.GetHistogram())) {
412 AliWarning("MC Histogram collector failed");
416 if (fAODFMD.IsTriggerBits(AliAODForwardMult::kInel))
417 fHData->Add(&(fAODFMD.GetHistogram()));
422 //____________________________________________________________________
424 AliForwardMCMultiplicityTask::Terminate(Option_t*)
432 TList* list = dynamic_cast<TList*>(GetOutputData(1));
434 AliError(Form("No output list defined (%p)", GetOutputData(1)));
435 if (GetOutputData(1)) GetOutputData(1)->Print();
439 // Get our histograms from the container
441 TH1I* hEventsTrVtx = 0;
443 if (!fEventInspector.FetchHistograms(list, hEventsTr,
444 hEventsTrVtx, hTriggers)) {
445 AliError(Form("Didn't get histograms from event selector "
446 "(hEventsTr=%p,hEventsTrVtx=%p)",
447 hEventsTr, hEventsTrVtx));
452 TH2D* hData = static_cast<TH2D*>(list->FindObject("d2Ndetadphi"));
454 AliError(Form("Couldn't get our summed histogram from output "
455 "list %s (d2Ndetadphi=%p)", list->GetName(), hData));
460 // TH1D* dNdeta = fHData->ProjectionX("dNdeta", 0, -1, "e");
461 TH1D* dNdeta = hData->ProjectionX("dNdeta", 1, -1, "e");
462 TH1D* norm = hData->ProjectionX("norm", 0, 1, "");
463 dNdeta->SetTitle("dN_{ch}/d#eta in the forward regions");
464 dNdeta->SetYTitle("#frac{1}{N}#frac{dN_{ch}}{d#eta}");
465 dNdeta->Divide(norm);
467 dNdeta->Scale(Double_t(hEventsTrVtx->GetEntries())/hEventsTr->GetEntries(),
472 MakeRingdNdeta(list, "ringSums", list, "ringResults");
473 MakeRingdNdeta(list, "mcRingSums", list, "mcRingResults", 24);
475 fSharingFilter.ScaleHistograms(list,Int_t(hEventsTr->Integral()));
476 fDensityCalculator.ScaleHistograms(list,Int_t(hEventsTrVtx->Integral()));
477 fCorrections.ScaleHistograms(list,Int_t(hEventsTrVtx->Integral()));