2 // Calculate the multiplicity in the forward regions event-by-event
14 #include "AliForwardQATask.h"
15 #include "AliForwardUtil.h"
16 #include "AliTriggerAnalysis.h"
17 #include "AliPhysicsSelection.h"
19 #include "AliESDEvent.h"
20 #include "AliAODHandler.h"
21 #include "AliMultiplicity.h"
22 #include "AliInputEventHandler.h"
23 #include "AliForwardCorrectionManager.h"
24 #include "AliAnalysisManager.h"
25 #include "AliAODForwardMult.h"
27 #include <TDirectory.h>
30 #include <TStopwatch.h>
32 //====================================================================
33 AliForwardQATask::AliForwardQATask()
34 : AliAnalysisTaskSE(),
35 fEnableLowFlux(false),
51 //____________________________________________________________________
52 AliForwardQATask::AliForwardQATask(const char* name)
53 : AliAnalysisTaskSE(name),
54 fEnableLowFlux(false),
59 fEventInspector("event"),
60 fEnergyFitter("energy"),
61 fSharingFilter("sharing"),
62 fDensityCalculator("density"),
71 DefineOutput(1, TList::Class());
72 DefineOutput(2, TList::Class());
73 fCorrManager = &AliForwardCorrectionManager::Instance();
74 fEnergyFitter.SetNParticles(1); // Just find the 1st peak
75 fEnergyFitter.SetDoMakeObject(false);
76 fEnergyFitter.SetUseIncreasingBins(true);
77 fEnergyFitter.SetDoFits(kTRUE);
78 fEnergyFitter.SetLowCut(0.4);
79 fEnergyFitter.SetFitRangeBinWidth(4);
80 fEnergyFitter.SetMinEntries(1000);
83 //____________________________________________________________________
84 AliForwardQATask::AliForwardQATask(const AliForwardQATask& o)
85 : AliAnalysisTaskSE(o),
86 fEnableLowFlux(o.fEnableLowFlux),
87 fFirstEvent(o.fFirstEvent),
88 fCorrManager(o.fCorrManager),
91 fEventInspector(o.fEventInspector),
92 fEnergyFitter(o.fEnergyFitter),
93 fSharingFilter(o.fSharingFilter),
94 fDensityCalculator(o.fDensityCalculator),
101 // o Object to copy from
103 DefineOutput(1, TList::Class());
104 DefineOutput(2, TList::Class());
107 //____________________________________________________________________
109 AliForwardQATask::operator=(const AliForwardQATask& o)
112 // Assignment operator
115 // o Object to assign from
118 // Reference to this object
120 AliAnalysisTaskSE::operator=(o);
122 fEnableLowFlux = o.fEnableLowFlux;
123 fFirstEvent = o.fFirstEvent;
124 fCorrManager = o.fCorrManager;
125 fEventInspector = o.fEventInspector;
126 fEnergyFitter = o.fEnergyFitter;
127 fSharingFilter = o.fSharingFilter;
128 fDensityCalculator = o.fDensityCalculator;
135 //____________________________________________________________________
137 AliForwardQATask::SetDebug(Int_t dbg)
145 fEventInspector.SetDebug(dbg);
146 fEnergyFitter.SetDebug(dbg);
147 fSharingFilter.SetDebug(dbg);
148 fDensityCalculator.SetDebug(dbg);
151 //____________________________________________________________________
153 AliForwardQATask::CheckCorrections(UInt_t what) const
156 // Check if all needed corrections are there and accounted for. If not,
160 // what Which corrections is needed
163 // true if all present, false otherwise
166 AliForwardCorrectionManager& fcm = AliForwardCorrectionManager::Instance();
167 // Check that we have the energy loss fits, needed by
168 // AliFMDSharingFilter
169 // AliFMDDensityCalculator
170 if (what & AliForwardCorrectionManager::kELossFits && !fcm.GetELossFit()) {
171 AliWarning("No energy loss fits");
177 //____________________________________________________________________
179 AliForwardQATask::ReadCorrections(const TAxis*& pe,
187 UInt_t what = AliForwardCorrectionManager::kAll;
188 what ^= AliForwardCorrectionManager::kDoubleHit;
189 what ^= AliForwardCorrectionManager::kVertexBias;
190 what ^= AliForwardCorrectionManager::kAcceptance;
191 what ^= AliForwardCorrectionManager::kMergingEfficiency;
193 AliForwardCorrectionManager& fcm = AliForwardCorrectionManager::Instance();
194 if (!fcm.Init(GetEventInspector().GetCollisionSystem(),
195 GetEventInspector().GetEnergy(),
196 GetEventInspector().GetField(),
199 if (!CheckCorrections(what)) {
203 // Sett our persistency pointer
204 // fCorrManager = &fcm;
206 // Get the eta axis from the secondary maps - if read in
208 pe = fcm.GetEtaAxis();
209 if (!pe) AliFatal("No eta axis defined");
211 // Get the vertex axis from the secondary maps - if read in
213 pv = fcm.GetVertexAxis();
214 if (!pv) AliFatal("No vertex axis defined");
220 //____________________________________________________________________
222 AliForwardQATask::GetESDEvent()
225 // Get the ESD event. IF this is the first event, initialise
227 AliESDEvent* esd = dynamic_cast<AliESDEvent*>(InputEvent());
229 AliWarning("No ESD event found for input event");
233 // On the first event, initialize the parameters
234 if (fFirstEvent && esd->GetESDRun()) {
235 GetEventInspector().ReadRunDetails(esd);
237 AliInfoF("Initializing with parameters from the ESD:\n"
238 " AliESDEvent::GetBeamEnergy() ->%f\n"
239 " AliESDEvent::GetBeamType() ->%s\n"
240 " AliESDEvent::GetCurrentL3() ->%f\n"
241 " AliESDEvent::GetMagneticField()->%f\n"
242 " AliESDEvent::GetRunNumber() ->%d\n",
243 esd->GetBeamEnergy(),
246 esd->GetMagneticField(),
247 esd->GetRunNumber());
250 if (!InitializeSubs()) {
251 AliWarning("Initialisation of sub algorithms failed!");
254 AliInfoF("Clearing first event flag from %s to false",
255 fFirstEvent ? "true" : "false");
260 //____________________________________________________________________
262 AliForwardQATask::InitializeSubs()
265 // Initialise the sub objects and stuff. Called on first event
271 if (!ReadCorrections(pe,pv)) {
272 AliWarning("Failed to read corrections");
278 fEventInspector.Init(*pv);
279 fEnergyFitter.Init(*pe);
280 fSharingFilter.Init();
281 fDensityCalculator.Init(*pe);
288 //____________________________________________________________________
290 AliForwardQATask::UserCreateOutputObjects()
293 // Create output objects
299 fEventInspector.DefineOutput(fList);
300 fEnergyFitter.DefineOutput(fList);
301 fSharingFilter.DefineOutput(fList);
302 fDensityCalculator.DefineOutput(fList);
306 //____________________________________________________________________
308 AliForwardQATask::UserExec(Option_t*)
311 // Process each event
317 // static Int_t cnt = 0;
319 // Get the input data
320 AliESDEvent* esd = GetESDEvent();
322 AliWarning("Got no ESD event");
326 // If the first event flag wasn't cleared in the above call to
327 // GetESDEvent, we should not do anything, since nothing has been
328 // initialised yet, so we opt out here (with a warning)
329 AliWarning("Nothing has been initialized yet, opt'ing out");
337 Bool_t lowFlux = kFALSE;
342 UShort_t nClusters = 0;
343 UInt_t found = fEventInspector.Process(esd, triggers, lowFlux,
344 ivz, vz, cent, nClusters);
346 if (found & AliFMDEventInspector::kNoEvent) return;
347 if (found & AliFMDEventInspector::kNoTriggers) return;
348 if (found & AliFMDEventInspector::kNoSPD) return;
349 if (found & AliFMDEventInspector::kNoFMD) return;
350 if (found & AliFMDEventInspector::kNoVertex) return;
351 if (triggers & AliAODForwardMult::kPileUp) return;
352 if (found & AliFMDEventInspector::kBadVertex) return;
354 // We we do not want to use low flux specific code, we disable it here.
355 if (!fEnableLowFlux) lowFlux = false;
358 AliESDFMD* esdFMD = esd->GetFMDData();
360 // Run the energy loss fitter
361 if (!fEnergyFitter.Accumulate(*esdFMD, cent,
362 triggers & AliAODForwardMult::kEmpty)) {
363 AliWarning("Energy fitter failed");
367 // // Apply the sharing filter (or hit merging or clustering if you like)
368 if (!fSharingFilter.Filter(*esdFMD, lowFlux, fESDFMD)) {
369 AliWarning("Sharing filter failed!");
373 // Calculate the inclusive charged particle density
374 if (!fDensityCalculator.Calculate(fESDFMD, fHistos, ivz, lowFlux)) {
375 // if (!fDensityCalculator.Calculate(*esdFMD, fHistos, ivz, lowFlux)) {
376 AliWarning("Density calculator failed!");
382 //____________________________________________________________________
384 AliForwardQATask::Terminate(Option_t*)
395 TList* list = dynamic_cast<TList*>(GetOutputData(1));
397 AliError(Form("No output list defined (%p)", GetOutputData(1)));
398 if (GetOutputData(1)) GetOutputData(1)->Print();
402 // Get our histograms from the container
403 TH1I* hEventsTr = 0;//static_cast<TH1I*>(list->FindObject("nEventsTr"));
404 TH1I* hEventsTrVtx = 0;//static_cast<TH1I*>(list->FindObject("nEventsTrVtx"));
406 if (!fEventInspector.FetchHistograms(list, hEventsTr,
407 hEventsTrVtx, hTriggers)) {
408 AliError(Form("Didn't get histograms from event selector "
409 "(hEventsTr=%p,hEventsTrVtx=%p)",
410 hEventsTr, hEventsTrVtx));
416 fEnergyFitter.Fit(list);
418 AliInfoF("Fitting took %d real-time seconds, and %f CPU seconds",
419 Int_t(swf.RealTime()), swf.CpuTime());
421 fSharingFilter.ScaleHistograms(list,Int_t(hEventsTr->Integral()));
422 fDensityCalculator.ScaleHistograms(list,Int_t(hEventsTrVtx->Integral()));
424 // Make a deep copy and post that as output 2
425 TList* list2 = static_cast<TList*>(list->Clone(Form("%sResults",
431 AliInfoF("Terminate took %d real-time seconds, and %f CPU seconds",
432 Int_t(swt.RealTime()), swt.CpuTime());
436 //____________________________________________________________________
438 AliForwardQATask::Print(Option_t* option) const
447 std::cout << ClassName() << ": " << GetName() << "\n"
448 << " Enable low flux code: " << (fEnableLowFlux ? "yes" : "no")
450 << " Off-line trigger mask: 0x"
451 << std::hex << std::setfill('0')
452 << std::setw (8) << fOfflineTriggerMask
453 << std::dec << std::setfill (' ') << std::endl;
454 gROOT->IncreaseDirLevel();
455 if (fCorrManager) fCorrManager->Print();
457 std::cout << " Correction manager not set yet" << std::endl;
458 GetEventInspector() .Print(option);
459 GetEnergyFitter() .Print(option);
460 GetSharingFilter() .Print(option);
461 gROOT->DecreaseDirLevel();