]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGLF/FORWARD/analysis2/AliFMDEventInspector.cxx
Refactoring for AliFMDEnergyFitter to allow sub-classing
[u/mrichter/AliRoot.git] / PWGLF / FORWARD / analysis2 / AliFMDEventInspector.cxx
CommitLineData
7984e5f7 1//
2// This class inspects the event
3//
4// Input:
5// - AliESDFMD object possibly corrected for sharing
6//
7// Output:
8// - A histogram of v_z of events with triggers.
9// - A histogram of v_z of events with vertex and triggers
10// - A histogram of trigger counters
11//
12// Note, that these are added to the master output list
13//
14// Corrections used:
15// - None
16//
8565b10b 17#include "AliFMDEventInspector.h"
c8b1a7db 18#include "AliProdInfo.h"
8565b10b 19#include "AliLog.h"
20#include "AliESDEvent.h"
21#include "AliMultiplicity.h"
22#include "AliAnalysisManager.h"
11d40ecb 23#include "AliMCEventHandler.h"
8565b10b 24#include "AliInputEventHandler.h"
25#include "AliTriggerAnalysis.h"
26#include "AliPhysicsSelection.h"
e85a76b7 27#include "AliOADBPhysicsSelection.h"
8565b10b 28#include "AliAODForwardMult.h"
0bd4b00f 29#include "AliForwardUtil.h"
5e4d905e 30#include "AliCentrality.h"
8565b10b 31#include <TH1.h>
32#include <TList.h>
33#include <TDirectory.h>
0bd4b00f 34#include <TROOT.h>
e6463868 35#include <TParameter.h>
dd556bcd 36#include <TMatrixDSym.h>
c8b1a7db 37#include <TPRegexp.h>
0bd4b00f 38#include <iostream>
39#include <iomanip>
65abd48b 40#include "AliMCEvent.h"
41#include "AliHeader.h"
42#include "AliGenEventHeader.h"
43#include "AliCollisionGeometry.h"
da70cd6a 44#include "AliVVZERO.h"
e1f47419 45
8565b10b 46//====================================================================
5ca83fee 47const char* AliFMDEventInspector::fgkFolderName = "fmdEventInspector";
48
49//____________________________________________________________________
8565b10b 50AliFMDEventInspector::AliFMDEventInspector()
51 : TNamed(),
52 fHEventsTr(0),
53 fHEventsTrVtx(0),
5bb5d1f6 54 fHEventsAccepted(0),
96110c91 55 fHEventsAcceptedXY(0),
8565b10b 56 fHTriggers(0),
66cf95f2 57 fHTriggerCorr(0),
0bd4b00f 58 fHType(0),
5ca83fee 59 fHWords(0),
60 fHCent(0),
61 fHCentVsQual(0),
62 fHStatus(0),
8449e3e0 63 fHVtxStatus(0),
e65b8b56 64 fHTrgStatus(0),
5ca83fee 65 fLowFluxCut(1000),
66 fMaxVzErr(0.2),
67 fList(0),
68 fEnergy(0),
0bd4b00f 69 fField(999),
5ca83fee 70 fCollisionSystem(kUnknown),
e308a636 71 fDebug(0),
5bb5d1f6 72 fCentAxis(0),
e83d0620 73 fVtxAxis(10,-10,10),
5ca83fee 74 fUseFirstPhysicsVertex(false),
75 fUseV0AND(false),
0b7de667 76 fMinPileupContrib(3),
77 fMinPileupDistance(0.8),
e65b8b56 78 fUseDisplacedVertices(false),
0b7de667 79 fDisplacedVertex(),
80 fCollWords(),
81 fBgWords(),
82 fCentMethod("V0M"),
83 fMinCent(-1.0),
84 fMaxCent(-1.0),
85 fUsepA2012Vertex(false),
86 fRunNumber(0),
87 fMC(false),
88 fProdYear(-1),
89 fProdLetter('?'),
90 fProdPass(-1),
91 fProdSVN(-1),
92 fProdMC(false)
8565b10b 93{
7984e5f7 94 //
95 // Constructor
96 //
241cca4d 97 DGUARD(fDebug,1,"Default CTOR of AliFMDEventInspector");
8565b10b 98}
99
100//____________________________________________________________________
101AliFMDEventInspector::AliFMDEventInspector(const char* name)
5ca83fee 102 : TNamed(fgkFolderName, name),
8565b10b 103 fHEventsTr(0),
104 fHEventsTrVtx(0),
5bb5d1f6 105 fHEventsAccepted(0),
96110c91 106 fHEventsAcceptedXY(0),
8565b10b 107 fHTriggers(0),
66cf95f2 108 fHTriggerCorr(0),
0bd4b00f 109 fHType(0),
fe52e455 110 fHWords(0),
5e4d905e 111 fHCent(0),
e308a636 112 fHCentVsQual(0),
241cca4d 113 fHStatus(0),
8449e3e0 114 fHVtxStatus(0),
115 fHTrgStatus(0),
8565b10b 116 fLowFluxCut(1000),
9d05ffeb 117 fMaxVzErr(0.2),
8565b10b 118 fList(0),
0bd4b00f 119 fEnergy(0),
120 fField(999),
121 fCollisionSystem(kUnknown),
e308a636 122 fDebug(0),
5bb5d1f6 123 fCentAxis(0),
e83d0620 124 fVtxAxis(10,-10,10),
5ca83fee 125 fUseFirstPhysicsVertex(false),
e6463868 126 fUseV0AND(false),
127 fMinPileupContrib(3),
65abd48b 128 fMinPileupDistance(0.8),
129 fUseDisplacedVertices(false),
241cca4d 130 fDisplacedVertex(),
131 fCollWords(),
5934a3e3 132 fBgWords(),
4f9319f3 133 fCentMethod("V0M"),
8449e3e0 134 fMinCent(-1.0),
135 fMaxCent(-1.0),
136 fUsepA2012Vertex(false),
e65b8b56 137 fRunNumber(0),
0b7de667 138 fMC(false),
139 fProdYear(-1),
140 fProdLetter('?'),
141 fProdPass(-1),
142 fProdSVN(-1),
143 fProdMC(false)
8565b10b 144{
7984e5f7 145 //
146 // Constructor
147 //
148 // Parameters:
149 // name Name of object
150 //
241cca4d 151 DGUARD(fDebug,1,"Named CTOR of AliFMDEventInspector: %s", name);
8565b10b 152}
153
8565b10b 154//____________________________________________________________________
155AliFMDEventInspector::~AliFMDEventInspector()
156{
7984e5f7 157 //
158 // Destructor
159 //
c929bc03 160 DGUARD(fDebug,1,"DTOR of AliFMDEventInspector");
161 // if (fList) delete fList;
8565b10b 162}
8565b10b 163
5934a3e3 164//____________________________________________________________________
165void
166AliFMDEventInspector::SetCentralityMethod(ECentMethod m)
167{
168 switch (m) {
169 case kV0Multiplicity: fCentMethod = "VOM"; break; // VZERO multiplicity
170 case kV0Amplitude: fCentMethod = "V0A"; break; // VZERO amplitude
171 case kV0Charge: fCentMethod = "V0C"; break; // VZERO charge
172 case kFMDRough: fCentMethod = "FMD"; break; // FMD scaled energy l
173 case kNTracks: fCentMethod = "TRK"; break; // Number of tracks
174 case kLTracks: fCentMethod = "TKL"; break; // Number of tracks
175 case kCL0: fCentMethod = "CL0"; break; //
176 case kCL1: fCentMethod = "CL1"; break; //
177 case kCND: fCentMethod = "CND"; break; //
178 case kNParticles: fCentMethod = "NPA"; break; // Neutral particles
179 case kNeutrons: fCentMethod = "ZNA"; break; // ZDC neutron amplitu
180 case kV0vsFMD: fCentMethod = "V0MvsFMD"; break; // VZERO versus FMD
181 case kV0vsNTracks: fCentMethod = "TKLvsVOM"; break; // Tracks versus VZERO
182 case kZEMvsZDC: fCentMethod = "ZEMvsZDC"; break; // ZDC
183 default: fCentMethod = "V0M"; break;
184 }
185}
186
8449e3e0 187//____________________________________________________________________
188void
189AliFMDEventInspector::SetMinCentrality(Double_t minCent)
190{
191 AliWarning("\n"
192 "*******************************************************\n"
193 "* Setting centrality cuts in this stage is deprecated *\n"
194 "*******************************************************");
195 fMinCent = minCent;
196}
197//____________________________________________________________________
198void
199AliFMDEventInspector::SetMaxCentrality(Double_t maxCent)
200{
201 AliWarning("\n"
202 "*******************************************************\n"
203 "* Setting centrality cuts in this stage is deprecated *\n"
204 "*******************************************************");
205 fMaxCent = maxCent;
206}
207
8565b10b 208//____________________________________________________________________
209Bool_t
fb3430ac 210AliFMDEventInspector::FetchHistograms(const TList* d,
8565b10b 211 TH1I*& hEventsTr,
212 TH1I*& hEventsTrVtx,
5ca83fee 213 TH1I*& hEventsAcc,
8565b10b 214 TH1I*& hTriggers) const
215{
7984e5f7 216 //
217 // Fetch our histograms from the passed list
218 //
219 // Parameters:
220 // d Input
221 // hEventsTr On return, pointer to histogram, or null
222 // hEventsTrVtx On return, pointer to histogram, or null
223 // hTriggers On return, pointer to histogram, or null
224 //
225 // Return:
226 // true on success, false otherwise
227 //
6ab100ec 228 DGUARD(fDebug,3,"Fetch histograms in AliFMDEventInspector");
8565b10b 229 hEventsTr = 0;
230 hEventsTrVtx = 0;
5ca83fee 231 hEventsAcc = 0;
232 hTriggers = 0;
8565b10b 233 TList* dd = dynamic_cast<TList*>(d->FindObject(GetName()));
234 if (!dd) return kFALSE;
235
236 hEventsTr = dynamic_cast<TH1I*>(dd->FindObject("nEventsTr"));
237 hEventsTrVtx = dynamic_cast<TH1I*>(dd->FindObject("nEventsTrVtx"));
5ca83fee 238 hEventsAcc = dynamic_cast<TH1I*>(dd->FindObject("nEventsAccepted"));
8565b10b 239 hTriggers = dynamic_cast<TH1I*>(dd->FindObject("triggers"));
240
5ca83fee 241 if (!hEventsTr ||
242 !hEventsTrVtx ||
243 !hEventsAcc ||
244 !hTriggers) return kFALSE;
8565b10b 245 return kTRUE;
246}
241cca4d 247//____________________________________________________________________
248void
249AliFMDEventInspector::CacheConfiguredTriggerClasses(TList& cache,
250 const TList* classes,
251 AliOADBPhysicsSelection* o)
252{
253 TIter nextClass(classes);
254 TObjString* trigClass = 0;
255 // Loop over all trigger classes. Trigger classes have the format
256 //
257 // class := positive_words SPACE(s) negative_words
258 // positive_words :=
259 // | '+' words
260 // negative_words :=
261 // | '-' words
262 // words := word
263 // | word ',' words
264 //
265 while ((trigClass = static_cast<TObjString*>(nextClass()))) {
266 // Tokenize on space to get positive and negative parts
267 TString side = o->GetBeamSide(trigClass->String());
268 TObjArray* parts = trigClass->String().Tokenize(" ");
269 TObjString* part = 0;
270 TIter nextPart(parts);
271 while ((part = static_cast<TObjString*>(nextPart()))) {
272 // We only care about the positive ones
273 if (part->GetName()[0] != '+') continue;
274 part->String().Remove(0,1);
275
276 // Tokenize on a comma to get the words
277 TObjArray* words = part->String().Tokenize(",");
278 TObjString* word = 0;
279 TIter nextWord(words);
280 while ((word = static_cast<TObjString*>(nextWord()))) {
281 TNamed* store = new TNamed(word->String(), side);
282 cache.Add(store);
283 DMSG(fDebug,3,"Caching %s trigger word %s",
284 store->GetTitle(), store->GetName());
285 } // while (word)
286 delete words;
287 }
288 delete parts;
289 }
290}
291
8565b10b 292//____________________________________________________________________
293void
5934a3e3 294AliFMDEventInspector::SetupForData(const TAxis& vtxAxis)
8565b10b 295{
7984e5f7 296 //
241cca4d 297 // Initialize the object - this is called on the first seen event.
7984e5f7 298 //
299 // Parameters:
300 // vtxAxis Vertex axis in use
301 //
6ab100ec 302 DGUARD(fDebug,1,"Initialize in AliFMDEventInspector");
241cca4d 303
304 AliAnalysisManager* am = AliAnalysisManager::GetAnalysisManager();
305
306 // Get the input handler - should always be there
307 AliInputEventHandler* ih =
308 static_cast<AliInputEventHandler*>(am->GetInputEventHandler());
309 if (!ih) {
310 AliWarning("No input handler");
311 return;
312 }
313 // Get the physics selection - should always be there
314 AliPhysicsSelection* ps =
315 static_cast<AliPhysicsSelection*>(ih->GetEventSelection());
316 if (!ps) {
317 AliWarning("No physics selection");
318 return;
319 }
320 // Get the configured triggers
321 AliOADBPhysicsSelection* oadb =
322 const_cast<AliOADBPhysicsSelection*>(ps->GetOADBPhysicsSelection());
323 if (!oadb) {
324 AliWarning("No OADB physics selection object");
325 return;
326 }
327 // Get the configured trigger words from the physics selection
328 const TList* collTriggClasses = ps->GetCollisionTriggerClasses();
329 const TList* bgTriggClasses = ps->GetBGTriggerClasses();
330 if (!collTriggClasses) {
331 AliWarning("No configured collision trigger classes");
332 return;
333 }
334 if (!bgTriggClasses) {
335 AliWarning("No configured background trigger classes");
336 return;
337 }
338 CacheConfiguredTriggerClasses(fCollWords, collTriggClasses, oadb);
339 CacheConfiguredTriggerClasses(fBgWords, bgTriggClasses, oadb);
340 // fCollWords.ls();
341 // fBgWords.ls();
342
e308a636 343
8449e3e0 344 TArrayD limits;
345 if ((fMinCent < 0 && fMaxCent < 0) || fMaxCent <= fMinCent) {
346 // -1.5 -0.5 0.5 1.5 ... 89.5 ... 100.5
347 // ----- 92 number --------- ---- 1 ---
348 limits.Set(93);
349 for (Int_t i = 0; i < 92; i++) limits[i] = -1.5 + i;
350 limits[92] = 100.5;
351 }
352 else {
353 Int_t n = fMaxCent-fMinCent+2;
354 limits.Set(n);
355 for (Int_t i = 0; i < n; i++) {
356 limits[i] = fMinCent + i - .5;
357 }
358 }
359
5bb5d1f6 360 fVtxAxis.Set(vtxAxis.GetNbins(), vtxAxis.GetXmin(), vtxAxis.GetXmax());
e308a636 361
362 fCentAxis = new TAxis(limits.GetSize()-1, limits.GetArray());
8565b10b 363 fHEventsTr = new TH1I("nEventsTr", "Number of events w/trigger",
5bb5d1f6 364 4*vtxAxis.GetNbins(),
365 2*vtxAxis.GetXmin(),
366 2*vtxAxis.GetXmax());
8565b10b 367 fHEventsTr->SetXTitle("v_{z} [cm]");
368 fHEventsTr->SetYTitle("# of events");
369 fHEventsTr->SetFillColor(kRed+1);
370 fHEventsTr->SetFillStyle(3001);
371 fHEventsTr->SetDirectory(0);
372 // fHEventsTr->Sumw2();
373 fList->Add(fHEventsTr);
374
5bb5d1f6 375 fHEventsTrVtx = static_cast<TH1I*>(fHEventsTr->Clone("nEventsTrVtx"));
376 fHEventsTrVtx->SetTitle("Number of events w/trigger and vertex");
8565b10b 377 fHEventsTrVtx->SetFillColor(kBlue+1);
8565b10b 378 fHEventsTrVtx->SetDirectory(0);
379 // fHEventsTrVtx->Sumw2();
380 fList->Add(fHEventsTrVtx);
381
5bb5d1f6 382 fHEventsAccepted = new TH1I("nEventsAccepted",
383 "Number of events w/trigger and vertex in range",
384 2*vtxAxis.GetNbins(),
385 2*vtxAxis.GetXmin(),
386 2*vtxAxis.GetXmax());
387 fHEventsAccepted->SetXTitle("v_{z} [cm]");
388 fHEventsAccepted->SetYTitle("# of events");
389 fHEventsAccepted->SetFillColor(kGreen+1);
390 fHEventsAccepted->SetFillStyle(3001);
391 fHEventsAccepted->SetDirectory(0);
392 // fHEventsAccepted->Sumw2();
393 fList->Add(fHEventsAccepted);
394
96110c91 395 fHEventsAcceptedXY = new TH2D("nEventsAcceptedXY",
396 "XY vertex w/trigger and Z vertex in range",
397 1000,-1,1,1000,-1,1);
398
399 fHEventsAcceptedXY->SetXTitle("v_{x} [cm]");
400 fHEventsAcceptedXY->SetYTitle("v_{y} [cm]");
401 fHEventsAcceptedXY->SetDirectory(0);
402 // fHEventsAccepted->Sumw2();
403 fList->Add(fHEventsAcceptedXY);
404
405
0be6c8cd 406 fHTriggers = new TH1I("triggers", "Triggers", kOffline+1, 0, kOffline+1);
8565b10b 407 fHTriggers->SetFillColor(kRed+1);
408 fHTriggers->SetFillStyle(3001);
409 fHTriggers->SetStats(0);
410 fHTriggers->SetDirectory(0);
66cf95f2 411
412 fHTriggerCorr = new TH2I("triggerCorr", "Trigger correlation",
413 kOffline+1, 0, kOffline+1,
414 kOffline+1, 0, kOffline+1);
415 fHTriggerCorr->SetStats(0);
416 fHTriggerCorr->SetDirectory(0);
cd51b0dd 417 fHTriggerCorr->SetXTitle("Requirement");
418 fHTriggerCorr->SetYTitle("Companion");
419
66cf95f2 420 Int_t binNum[] = { kInel +1,
421 kInelGt0+1,
422 kNSD +1,
423 kV0AND +1,
424 kEmpty +1,
425 kA +1,
426 kB +1,
427 kC +1,
428 kE +1,
429 kPileUp +1,
430 kMCNSD +1,
8449e3e0 431 kSatellite+1,
66cf95f2 432 kOffline+1 };
433 const char* binLbl[] = { "INEL",
434 "INEL>0",
435 "NSD",
436 "VOAND",
437 "Empty",
438 "A",
439 "B",
440 "C",
441 "E",
442 "Pileup",
443 "NSD_{MC}",
8449e3e0 444 "Satellite",
66cf95f2 445 "Offline" };
cd51b0dd 446 for (Int_t i = 0; i < kOffline+1; i++) {
66cf95f2 447 fHTriggers->GetXaxis()->SetBinLabel(binNum[i], binLbl[i]);
448 fHTriggerCorr->GetXaxis()->SetBinLabel(binNum[i], binLbl[i]);
449 fHTriggerCorr->GetYaxis()->SetBinLabel(binNum[i], binLbl[i]);
450 }
8565b10b 451 fList->Add(fHTriggers);
66cf95f2 452 fList->Add(fHTriggerCorr);
453
0bd4b00f 454
455 fHType = new TH1I("type", Form("Event type (cut: SPD mult>%d)",
456 fLowFluxCut), 2, -.5, 1.5);
457 fHType->SetFillColor(kRed+1);
458 fHType->SetFillStyle(3001);
459 fHType->SetStats(0);
460 fHType->SetDirectory(0);
461 fHType->GetXaxis()->SetBinLabel(1,"Low-flux");
462 fHType->GetXaxis()->SetBinLabel(2,"High-flux");
463 fList->Add(fHType);
fe52e455 464
c56f8519 465#if 0
466 // This histogram disabled as it causes problems in the merge
fe52e455 467 fHWords = new TH1I("words", "Trigger words seen", 1, 0, 0);
468 fHWords->SetFillColor(kBlue+1);
469 fHWords->SetFillStyle(3001);
470 fHWords->SetStats(0);
471 fHWords->SetDirectory(0);
472 fHWords->SetBit(TH1::kCanRebin);
473 fList->Add(fHWords);
c56f8519 474#endif
5e4d905e 475
e308a636 476 fHCent = new TH1F("cent", "Centrality", limits.GetSize()-1,limits.GetArray());
5e4d905e 477 fHCent->SetFillColor(kBlue+1);
478 fHCent->SetFillStyle(3001);
479 fHCent->SetStats(0);
480 fHCent->SetDirectory(0);
481 fHCent->SetXTitle("Centrality [%]");
482 fHCent->SetYTitle("Events");
483 fList->Add(fHCent);
e308a636 484
485 fHCentVsQual = new TH2F("centVsQuality", "Quality vs Centrality",
486 5, 0, 5, limits.GetSize()-1, limits.GetArray());
487 fHCentVsQual->SetXTitle("Quality");
488 fHCentVsQual->SetYTitle("Centrality [%]");
489 fHCentVsQual->SetZTitle("Events");
490 fHCentVsQual->GetXaxis()->SetBinLabel(1, "OK");
491 fHCentVsQual->GetXaxis()->SetBinLabel(2, "Outside v_{z} cut");
492 fHCentVsQual->GetXaxis()->SetBinLabel(3, "V0 vs SPD outlier");
493 fHCentVsQual->GetXaxis()->SetBinLabel(4, "V0 vs TPC outlier");
494 fHCentVsQual->GetXaxis()->SetBinLabel(5, "V0 vs ZDC outlier");
b7ab8a2c 495 fHCentVsQual->SetDirectory(0);
e308a636 496 fList->Add(fHCentVsQual);
241cca4d 497
498 fHStatus = new TH1I("status", "Status", 7, 1, 8);
8449e3e0 499 fHStatus->SetFillColor(kBlue+1);
241cca4d 500 fHStatus->SetFillStyle(3001);
501 fHStatus->SetStats(0);
502 fHStatus->SetDirectory(0);
8449e3e0 503 TAxis* xAxis = fHStatus->GetXaxis();
504 xAxis->SetBinLabel(1, "OK");
505 xAxis->SetBinLabel(2, "No event");
506 xAxis->SetBinLabel(3, "No triggers");
507 xAxis->SetBinLabel(4, "No SPD");
508 xAxis->SetBinLabel(5, "No FMD");
509 xAxis->SetBinLabel(6, "No vertex");
510 xAxis->SetBinLabel(7, "Bad vertex");
241cca4d 511 fList->Add(fHStatus);
8449e3e0 512
513 fHVtxStatus = new TH1I("vtxStatus","Vertex Status",
514 kNotVtxZ,kVtxOK,kNotVtxZ+1);
515 fHVtxStatus->SetFillColor(kGreen+1);
516 fHVtxStatus->SetFillStyle(3001);
517 fHVtxStatus->SetStats(0);
518 fHVtxStatus->SetDirectory(0);
519 xAxis = fHVtxStatus->GetXaxis();
520 xAxis->SetBinLabel(kVtxOK, "OK");
521 xAxis->SetBinLabel(kNoVtx, "None/bad status");
522 xAxis->SetBinLabel(kNoSPDVtx, "No SPD/bad status");
523 xAxis->SetBinLabel(kFewContrib, "N_{contrib} <= 0");
524 xAxis->SetBinLabel(kUncertain, Form("#delta z > %4.2f", fMaxVzErr));
525 xAxis->SetBinLabel(kNotVtxZ, "Not Z vertexer");
526 fList->Add(fHVtxStatus);
527
528 fHTrgStatus = new TH1I("trgStatus", "Trigger Status",
529 kOther, kNoTrgWords, kOther+1);
530 fHTrgStatus->SetFillColor(kMagenta+1);
531 fHTrgStatus->SetFillStyle(3001);
532 fHTrgStatus->SetStats(0);
533 fHTrgStatus->SetDirectory(0);
534 xAxis = fHTrgStatus->GetXaxis();
535 xAxis->SetBinLabel(kNoTrgWords, "No words");
536 xAxis->SetBinLabel(kPP2760Fast, "FAST in pp@#sqrt{s}=2.76TeV");
537 xAxis->SetBinLabel(kMUON, "Muon trigger");
538 xAxis->SetBinLabel(kTriggered, "Triggered");
539 xAxis->SetBinLabel(kMinBias, "CINT1 (V0A||V0C||FASTOR)");
540 xAxis->SetBinLabel(kMinBiasNoSPD, "CINT5 (V0A||V0C)");
541 xAxis->SetBinLabel(kV0AndTrg, "CINT7 (V0A&&V0C)");
542 xAxis->SetBinLabel(kHighMult, "N>>0");
543 xAxis->SetBinLabel(kCentral, "Central");
544 xAxis->SetBinLabel(kSemiCentral, "Semi-central");
545 xAxis->SetBinLabel(kDiffractive, "Diffractive");
546 xAxis->SetBinLabel(kUser, "User");
547 xAxis->SetBinLabel(kOther, "Other");
548 fList->Add(fHTrgStatus);
549
bfab35d9 550 if (fUseDisplacedVertices) fDisplacedVertex.SetupForData(fList, "", false);
8565b10b 551}
552
f7cfc454 553//____________________________________________________________________
554void
8449e3e0 555AliFMDEventInspector::StoreInformation()
f7cfc454 556{
557 // Write TNamed objects to output list containing information about
558 // the running conditions
6ab100ec 559 DGUARD(fDebug,2,"Store information from AliFMDEventInspector");
f7cfc454 560 if (!fList) return;
561
241cca4d 562
563 fList->Add(AliForwardUtil::MakeParameter("sys", fCollisionSystem));
564 fList->Add(AliForwardUtil::MakeParameter("sNN", fEnergy));
565 fList->Add(AliForwardUtil::MakeParameter("field", fField));
8449e3e0 566 fList->Add(AliForwardUtil::MakeParameter("runNo", fRunNumber));
241cca4d 567 fList->Add(AliForwardUtil::MakeParameter("lowFlux", fLowFluxCut));
568 fList->Add(AliForwardUtil::MakeParameter("fpVtx",fUseFirstPhysicsVertex));
569 fList->Add(AliForwardUtil::MakeParameter("v0and",fUseV0AND));
570 fList->Add(AliForwardUtil::MakeParameter("nPileUp", fMinPileupContrib));
571 fList->Add(AliForwardUtil::MakeParameter("dPileup", fMinPileupDistance));
8449e3e0 572 fList->Add(AliForwardUtil::MakeParameter("satellite", fUseDisplacedVertices));
1ff25622 573 fList->Add(AliForwardUtil::MakeParameter("alirootRev",
574 AliForwardUtil::AliROOTRevision()));
575 fList->Add(AliForwardUtil::MakeParameter("alirootBranch",
576 AliForwardUtil::AliROOTBranch()));
e65b8b56 577 fList->Add(AliForwardUtil::MakeParameter("mc", fMC));
1ff25622 578
f7cfc454 579}
580
c8b1a7db 581//____________________________________________________________________
582void
583AliFMDEventInspector::StoreProduction()
584{
585 if (!fList) return;
586
587 AliAnalysisManager *mgr=AliAnalysisManager::GetAnalysisManager();
588 AliVEventHandler *inputHandler=mgr->GetInputEventHandler();
589 if (!inputHandler) {
590 AliWarning("Got no input handler");
591 return;
592 }
593 TList *uiList = inputHandler->GetUserInfo();
594 if (!uiList) {
595 AliWarning("Got no user list from input tree");
596 return;
597 }
598
599 AliProdInfo p(uiList);
600 p.List();
601 if (p.GetAlirootSvnVersion() <= 0) return;
602
603 // Make our output list
604 TList* out = new TList;
605 out->SetOwner(true);
606 out->SetName("production");
607 // out->SetTitle("ESD production details");
608 fList->Add(out);
609
610 TString period = p.GetLHCPeriod();
611 // TString aliROOTVersion = p.GetAlirootVersion();
0b7de667 612 fProdSVN = p.GetAlirootSvnVersion();
c8b1a7db 613 // TString rootVersion = p.GetRootVersion();
614 Int_t rootSVN = p.GetRootSvnVersion();
0b7de667 615 fProdPass = p.GetRecoPass();
616 fProdMC = p.IsMC();
c8b1a7db 617
618 TObjArray* pp = TPRegexp("LHC([0-9]+)([a-z]+)").MatchS(period);
0b7de667 619 fProdYear = static_cast<TObjString*>(pp->At(1))->String().Atoi();
620 fProdLetter = static_cast<TObjString*>(pp->At(2))->String()[0];
c8b1a7db 621 pp->Delete();
622
0b7de667 623 out->Add(AliForwardUtil::MakeParameter("year", fProdYear));
624 out->Add(AliForwardUtil::MakeParameter("letter", Int_t(fProdLetter)));
625 out->Add(AliForwardUtil::MakeParameter("alirootSVN", fProdSVN));
c8b1a7db 626 out->Add(AliForwardUtil::MakeParameter("rootSVN", rootSVN));
0b7de667 627 out->Add(AliForwardUtil::MakeParameter("pass", fProdPass));
628 out->Add(AliForwardUtil::MakeParameter("mc", fProdMC));
c8b1a7db 629}
630
631
632
633
8565b10b 634//____________________________________________________________________
635void
5934a3e3 636AliFMDEventInspector::CreateOutputObjects(TList* dir)
8565b10b 637{
7984e5f7 638 //
639 // Define the output histograms. These are put in a sub list of the
640 // passed list. The histograms are merged before the parent task calls
641 // AliAnalysisTaskSE::Terminate
642 //
643 // dir Directory to add to
644 //
6ab100ec 645 DGUARD(fDebug,1,"Define output from AliFMDEventInspector");
8565b10b 646 fList = new TList;
647 fList->SetName(GetName());
c929bc03 648 fList->SetOwner();
8565b10b 649 dir->Add(fList);
650}
651
652//____________________________________________________________________
653UInt_t
654AliFMDEventInspector::Process(const AliESDEvent* event,
655 UInt_t& triggers,
656 Bool_t& lowFlux,
0bd4b00f 657 UShort_t& ivz,
5ca83fee 658 TVector3& ip,
5bb5d1f6 659 Double_t& cent,
660 UShort_t& nClusters)
8565b10b 661{
7984e5f7 662 //
663 // Process the event
664 //
665 // Parameters:
666 // event Input event
667 // triggers On return, the triggers fired
668 // lowFlux On return, true if the event is considered a low-flux
669 // event (according to the setting of fLowFluxCut)
670 // ivz On return, the found vertex bin (1-based). A zero
671 // means outside of the defined vertex range
672 // vz On return, the z position of the interaction
5e4d905e 673 // cent On return, the centrality - if not available < 0
7984e5f7 674 //
675 // Return:
676 // 0 (or kOk) on success, otherwise a bit mask of error codes
677 //
4f9319f3 678 DGUARD(fDebug,1,"Process event in AliFMDEventInspector");
e1f47419 679 // --- Check that we have an event ---------------------------------
8565b10b 680 if (!event) {
681 AliWarning("No ESD event found for input event");
241cca4d 682 fHStatus->Fill(2);
8565b10b 683 return kNoEvent;
684 }
685
bfab35d9 686 // --- Process satellite event information is requested ------------
687 if (fUseDisplacedVertices) {
688 if (!fDisplacedVertex.Process(event))
689 AliWarning("Failed to process satellite event");
690 }
691
e1f47419 692 // --- Read trigger information from the ESD and store in AOD object
241cca4d 693 if (!ReadTriggers(*event, triggers, nClusters)) {
0bd4b00f 694 if (fDebug > 2) {
695 AliWarning("Failed to read triggers from ESD"); }
241cca4d 696 fHStatus->Fill(3);
8565b10b 697 return kNoTriggers;
698 }
699
e1f47419 700 // --- Check if this is a high-flux event --------------------------
8565b10b 701 const AliMultiplicity* testmult = event->GetMultiplicity();
702 if (!testmult) {
0bd4b00f 703 if (fDebug > 3) {
704 AliWarning("No central multiplicity object found"); }
8565b10b 705 }
e1f47419 706 else
707 lowFlux = testmult->GetNumberOfTracklets() < fLowFluxCut;
5e4d905e 708
0bd4b00f 709 fHType->Fill(lowFlux ? 0 : 1);
314f6077 710
711 // --- Get the interaction point -----------------------------------
712 Bool_t vzOk = ReadVertex(*event, ip);
713 fHEventsTr->Fill(ip.Z());
714 if (!vzOk) {
715 if (fDebug > 3) {
716 AliWarning("Failed to read vertex from ESD"); }
717 fHStatus->Fill(6);
718 return kNoVertex;
719 }
720
721 // --- Read centrality information
e308a636 722 cent = -10;
723 UShort_t qual = 0;
241cca4d 724 if (!ReadCentrality(*event, cent, qual)) {
e1f47419 725 if (fDebug > 3)
726 AliWarning("Failed to get centrality");
8565b10b 727 }
314f6077 728 // --- check centrality cut
4f9319f3 729
8449e3e0 730 if(fMinCent > -0.0001 && cent < fMinCent) return kNoEvent;
731 if(fMaxCent > -0.0001 && cent > fMaxCent) return kNoEvent;
e308a636 732 fHCent->Fill(cent);
733 if (qual == 0) fHCentVsQual->Fill(0., cent);
734 else {
735 for (UShort_t i = 0; i < 4; i++)
736 if (qual & (1 << i)) fHCentVsQual->Fill(Double_t(i+1), cent);
737 }
314f6077 738
8565b10b 739
5ca83fee 740 fHEventsTrVtx->Fill(ip.Z());
96110c91 741
e1f47419 742 // --- Get the vertex bin ------------------------------------------
5ca83fee 743 ivz = fVtxAxis.FindBin(ip.Z());
5bb5d1f6 744 if (ivz <= 0 || ivz > fVtxAxis.GetNbins()) {
0bd4b00f 745 if (fDebug > 3) {
8565b10b 746 AliWarning(Form("Vertex @ %f outside of range [%f,%f]",
5ca83fee 747 ip.Z(), fVtxAxis.GetXmin(), fVtxAxis.GetXmax()));
5bb5d1f6 748 }
0bd4b00f 749 ivz = 0;
241cca4d 750 fHStatus->Fill(7);
8565b10b 751 return kBadVertex;
752 }
5ca83fee 753 fHEventsAccepted->Fill(ip.Z());
754 fHEventsAcceptedXY->Fill(ip.X(),ip.Y());
e58000b7 755
e1f47419 756 // --- Check the FMD ESD data --------------------------------------
757 if (!event->GetFMDData()) {
758 if (fDebug > 3) {
759 AliWarning("No FMD data found in ESD"); }
241cca4d 760 fHStatus->Fill(5);
e1f47419 761 return kNoFMD;
762 }
763
241cca4d 764 fHStatus->Fill(1);
8565b10b 765 return kOk;
766}
767
e1f47419 768//____________________________________________________________________
769Bool_t
241cca4d 770AliFMDEventInspector::ReadCentrality(const AliESDEvent& esd,
e308a636 771 Double_t& cent,
772 UShort_t& qual) const
e1f47419 773{
774 //
775 // Read centrality from event
776 //
777 // Parameters:
778 // esd Event
779 // cent On return, the centrality or negative if not found
780 //
781 // Return:
782 // False on error, true otherwise
783 //
6ab100ec 784 DGUARD(fDebug,2,"Read the centrality in AliFMDEventInspector");
241cca4d 785
241cca4d 786 cent = -1;
8449e3e0 787 qual = 1;
241cca4d 788 AliCentrality* centObj = const_cast<AliESDEvent&>(esd).GetCentrality();
8449e3e0 789 if (centObj) {
790 cent = centObj->GetCentralityPercentile(fCentMethod);
791 qual = centObj->GetQuality();
792 }
241cca4d 793
bc31b177 794 // We overwrite with satellite events, so we can be sure to get the
795 // centrality determination from the satellite vertex selection
796 if (fUseDisplacedVertices && fDisplacedVertex.IsSatellite()) {
8449e3e0 797 cent = fDisplacedVertex.GetCentralityPercentile();
798 qual = 0;
799 }
e1f47419 800
801 return true;
802}
803
da70cd6a 804//____________________________________________________________________
805Bool_t
806AliFMDEventInspector::CheckpAExtraV0(const AliESDEvent& esd) const
807{
808 if (fCollisionSystem != AliForwardUtil::kPPb) return true;
809
810 AliVVZERO* esdV0 = esd.GetVZEROData();
811 if ((esdV0->GetV0ADecision()!=1) || (esdV0->GetV0CDecision()!=1))
812 return false;
813 return true;
814}
815
8565b10b 816//____________________________________________________________________
817Bool_t
241cca4d 818AliFMDEventInspector::ReadTriggers(const AliESDEvent& esd, UInt_t& triggers,
5bb5d1f6 819 UShort_t& nClusters)
8565b10b 820{
7984e5f7 821 //
822 // Read the trigger information from the ESD event
823 //
824 // Parameters:
825 // esd ESD event
826 // triggers On return, contains the trigger bits
827 //
828 // Return:
829 // @c true on success, @c false otherwise
830 //
6ab100ec 831 DGUARD(fDebug,2,"Read the triggers in AliFMDEventInspector");
8565b10b 832 triggers = 0;
833
834 // Get the analysis manager - should always be there
835 AliAnalysisManager* am = AliAnalysisManager::GetAnalysisManager();
241cca4d 836 DMSG(fDebug,10,"Got analysis manager %p", am);
8565b10b 837 if (!am) {
838 AliWarning("No analysis manager defined!");
839 return kFALSE;
840 }
841
842 // Get the input handler - should always be there
843 AliInputEventHandler* ih =
844 static_cast<AliInputEventHandler*>(am->GetInputEventHandler());
241cca4d 845 DMSG(fDebug,10,"Got input handler %p", ih);
8565b10b 846 if (!ih) {
847 AliWarning("No input handler");
848 return kFALSE;
849 }
e85a76b7 850
e1f47419 851 // Check if this is a collision candidate (MB)
8449e3e0 852 ///
853 // Historic remark: Note, that we should use the value cached in the
854 // input handler rather than calling IsCollisionCandiate directly
855 // on the AliPhysicsSelection obejct. If we called the latter
856 // then the AliPhysicsSelection object would overcount by a factor
857 // of 2! :-(
858 UInt_t trgMask = ih->IsEventSelected();
859 Bool_t offline = trgMask;
860 Bool_t fastonly = (trgMask & AliVEvent::kFastOnly);
241cca4d 861 TString trigStr = esd.GetFiredTriggerClasses();
862
8449e3e0 863 if (trigStr.IsNull()) fHTrgStatus->Fill(kNoTrgWords);
241cca4d 864 if (fHWords) fHWords->Fill(trigStr.Data(), 1);
11d40ecb 865
241cca4d 866 if(fUseDisplacedVertices) {
867 DMSG(fDebug,3,"Using displaced vertex stuff");
8449e3e0 868 // if (TMath::Abs(fDisplacedVertex.GetVertexZ()) >= 999) offline = false;
869 if (fDisplacedVertex.IsSatellite())
870 triggers |= AliAODForwardMult::kSatellite;
871 }
872
873 if (CheckFastPartition(fastonly)) {
874 fHTrgStatus->Fill(kPP2760Fast);
875 offline = false;
241cca4d 876 }
65abd48b 877
8449e3e0 878 if (offline && CheckCosmics(trigStr)) {
879 fHTrgStatus->Fill(kMUON);
880 offline = false;
881 }
882 if (offline) fHTrgStatus->Fill(kTriggered);
883 Int_t f = 0;
884 if (trgMask & AliVEvent::kMB) f += fHTrgStatus->Fill(kMinBias);
885 if (trgMask & AliVEvent::kCINT5) f += fHTrgStatus->Fill(kMinBiasNoSPD);
886 if (trgMask & AliVEvent::kINT7) f += fHTrgStatus->Fill(kV0AndTrg);
887 if (trgMask & AliVEvent::kHighMult) f += fHTrgStatus->Fill(kHighMult);
888 if (trgMask & AliVEvent::kCentral) f += fHTrgStatus->Fill(kCentral);
889 if (trgMask & AliVEvent::kSemiCentral) f += fHTrgStatus->Fill(kSemiCentral);
890 if (trgMask & AliVEvent::kDG5) f += fHTrgStatus->Fill(kDiffractive);
891 if (trgMask & AliVEvent::kUserDefined) f += fHTrgStatus->Fill(kUser);
892 if (f <= 0) fHTrgStatus->Fill(kOther);
893
c929bc03 894 // if (!CheckpAExtraV0(esd)) offline = false;
241cca4d 895
896 DMSG(fDebug,2,"Event is %striggered by off-line", offline ? "" : "NOT ");
897
dd556bcd 898 if (offline) triggers |= AliAODForwardMult::kOffline;
899 // Only flag as in-elastic if a min-bias trigger is here
900 if (trgMask & AliVEvent::kMB) {
241cca4d 901 triggers |= AliAODForwardMult::kInel;
241cca4d 902 CheckINELGT0(esd, nClusters, triggers);
65abd48b 903 }
904
241cca4d 905 CheckNSD(esd,triggers);
cd51b0dd 906 CheckPileup(esd, triggers);
907 CheckEmpty(trigStr, triggers);
908 // if (CheckPileup(esd, triggers)) fHTriggers->Fill(kPileUp+.5);
909 // if (CheckEmpty(trigStr, triggers)) fHTriggers->Fill(kEmpty+.5);
241cca4d 910
8449e3e0 911 CheckWords(esd, triggers);
241cca4d 912
66cf95f2 913#define TEST_TRIG_BIN(RET,BIN,TRIGGERS) \
914 do { switch (BIN) { \
8449e3e0 915 case kInel: RET = triggers & AliAODForwardMult::kInel; break; \
916 case kInelGt0: RET = triggers & AliAODForwardMult::kInelGt0; break; \
917 case kNSD: RET = triggers & AliAODForwardMult::kNSD; break; \
918 case kV0AND: RET = triggers & AliAODForwardMult::kV0AND; break; \
919 case kEmpty: RET = triggers & AliAODForwardMult::kEmpty; break; \
920 case kA: RET = triggers & AliAODForwardMult::kA; break; \
921 case kB: RET = triggers & AliAODForwardMult::kB; break; \
922 case kC: RET = triggers & AliAODForwardMult::kC; break; \
923 case kE: RET = triggers & AliAODForwardMult::kE; break; \
924 case kPileUp: RET = triggers & AliAODForwardMult::kPileUp; break; \
925 case kMCNSD: RET = triggers & AliAODForwardMult::kMCNSD; break; \
926 case kSatellite:RET = triggers & AliAODForwardMult::kSatellite; break; \
927 case kOffline: RET = triggers & AliAODForwardMult::kOffline; break; \
928 default: RET = false; } } while(false)
66cf95f2 929
8449e3e0 930 if (!fHTriggers) {
931 AliWarning("Histogram of triggers not defined - has init been called");
932 return false;
933 }
241cca4d 934
66cf95f2 935 for (Int_t i = 0; i < kOffline+1; i++) {
936 Bool_t hasX = false;
937 TEST_TRIG_BIN(hasX, i, triggers);
938 if (!hasX) continue;
cd51b0dd 939 fHTriggers->Fill(i+.5);
940 for (Int_t j = 0; j < kOffline+1; j++) {
66cf95f2 941 Bool_t hasY = false;
942 TEST_TRIG_BIN(hasY, j, triggers);
943 if (!hasY) continue;
944
945 fHTriggerCorr->Fill(i+.5, j+.5);
946 }
947 }
241cca4d 948 return kTRUE;
949}
950
951//____________________________________________________________________
952Bool_t
953AliFMDEventInspector::CheckFastPartition(bool fastonly) const
954{
11d40ecb 955 // For the 2.76 TeV p+p run, the FMD ran in the slow partition
956 // so it received no triggers from the fast partition. Therefore
957 // the fast triggers are removed here but not for MC where all
958 // triggers are fast.
241cca4d 959 if (TMath::Abs(fEnergy - 2750.) > 20) return false;
960 if (fCollisionSystem != AliForwardUtil::kPP) return false;
e65b8b56 961 if (fMC) return false; // All MC events for pp @ 2.76TeV are `fast'
241cca4d 962 if (fastonly)
963 DMSG(fDebug,1,"Fast trigger in pp @ sqrt(s)=2.76TeV removed");
964
965 return fastonly;
966}
967
968//____________________________________________________________________
969Bool_t
970AliFMDEventInspector::CheckCosmics(const TString& trigStr) const
971{
972 // MUON triggers are not strictly minimum bias (MB) so they are
973 // removed (HHD)
974 if(trigStr.Contains("CMUS1")) {
975 DMSG(fDebug,1,"Cosmic trigger ins't min-bias, removed");
976 return true;
977 }
978 return false;
979}
980
981//____________________________________________________________________
982Bool_t
983AliFMDEventInspector::CheckINELGT0(const AliESDEvent& esd,
984 UShort_t& nClusters,
985 UInt_t& triggers) const
986{
5bb5d1f6 987 nClusters = 0;
8565b10b 988
241cca4d 989 // If this is inel, see if we have a tracklet
990 const AliMultiplicity* spdmult = esd.GetMultiplicity();
991 if (!spdmult) {
992 AliWarning("No SPD multiplicity");
993 return false;
994 }
995
996 // Check if we have one or more tracklets
997 // in the range -1 < eta < 1 to set the INEL>0
998 // trigger flag.
999 //
1000 // Also count tracklets as a single cluster
1001 Int_t n = spdmult->GetNumberOfTracklets();
1002 for (Int_t j = 0; j < n; j++) {
dd556bcd 1003 if(TMath::Abs(spdmult->GetEta(j)) < 1) nClusters++;
8565b10b 1004 }
dd556bcd 1005 // If we at this point had non-zero nClusters, it's INEL>0
1006 if (nClusters > 0) triggers |= AliAODForwardMult::kInelGt0;
1007
1008 // Loop over single clusters
241cca4d 1009 n = spdmult->GetNumberOfSingleClusters();
1010 for (Int_t j = 0; j < n; j++) {
1011 Double_t eta = -TMath::Log(TMath::Tan(spdmult->GetThetaSingle(j)/2.));
1012 if (TMath::Abs(eta) < 1) nClusters++;
1013 }
1014 if (nClusters > 0) triggers |= AliAODForwardMult::kNClusterGt0;
1015
1016 return triggers & AliAODForwardMult::kNClusterGt0;
1017}
1018
1019//____________________________________________________________________
1020Bool_t
1021AliFMDEventInspector::CheckNSD(const AliESDEvent& esd, UInt_t& triggers) const
1022{
8565b10b 1023 // Analyse some trigger stuff
1024 AliTriggerAnalysis ta;
241cca4d 1025 if (ta.IsOfflineTriggerFired(&esd, AliTriggerAnalysis::kV0AND)) {
e6463868 1026 triggers |= AliAODForwardMult::kV0AND;
1027 if (fUseV0AND)
31554871 1028 triggers |= AliAODForwardMult::kNSD;
1029 }
241cca4d 1030 if (ta.IsOfflineTriggerFired(&esd, AliTriggerAnalysis::kNSD1))
e6463868 1031 triggers |= AliAODForwardMult::kNSD;
241cca4d 1032 return triggers & AliAODForwardMult::kNSD;
1033}
1034//____________________________________________________________________
1035Bool_t
1036AliFMDEventInspector::CheckPileup(const AliESDEvent& esd,
1037 UInt_t& triggers) const
1038{
5bb5d1f6 1039 // Check for multiple vertices (pile-up) with at least 3
1040 // contributors and at least 0.8cm from the primary vertex
dd556bcd 1041 // if(fCollisionSystem != AliForwardUtil::kPP) return false;
1042
1043 // Check for standard SPD pile-up
1044 Bool_t spdPileup = esd.IsPileupFromSPD(fMinPileupContrib,fMinPileupDistance);
241cca4d 1045
dd556bcd 1046 // Check for multi-vertex pileup
1047 Bool_t mvPileup = false; // CheckMultiVertex(esd);
1048
1049 // Check for out-of-bunch pileup
1050 Bool_t outPileup = (esd.GetHeader()->GetIRInt2ClosestInteractionMap() != 0 ||
1051 esd.GetHeader()->GetIRInt1ClosestInteractionMap() != 0);
1052
1053 // Our result
1054 Bool_t pileup = spdPileup || mvPileup || outPileup;
241cca4d 1055 if (pileup) triggers |= AliAODForwardMult::kPileUp;
1056 return pileup;
1057}
1058
dd556bcd 1059//____________________________________________________________________
1060Bool_t
c8b1a7db 1061AliFMDEventInspector::CheckMultiVertex(const AliESDEvent& esd,
1062 Bool_t checkOtherBC) const
dd556bcd 1063{
1064 // Adapted from AliAnalysisUtils
1065 //
1066 // Parameters
1067 const Double_t maxChi2nu = 5;
dd556bcd 1068
1069 // Number of vertices
1070 Int_t n = esd.GetNumberOfPileupVerticesTracks();
1071 if (n < 1)
1072 // No additional vertices from tracks -> no pileup
1073 return false;
1074
1075 // Primary vertex
1076 const AliESDVertex* primary = esd.GetPrimaryVertexTracks();
1077 if (primary->GetStatus() != 1)
1078 // No good primary vertex, but many candidates -> pileup
1079 return true;
1080
1081 // Bunch crossing number
1082 Int_t bcPrimary = primary->GetBC();
1083
1084 for (Int_t i = 0; i < n; i++) {
1085 const AliESDVertex* other = esd.GetPileupVertexTracks(i);
1086
1087 if (other->GetNContributors() < fMinPileupContrib)
1088 // Not enough contributors to this vertex
1089 continue;
1090
1091 if (other->GetChi2perNDF() > maxChi2nu)
1092 // Poorly determined vertex
1093 continue;
1094 if (checkOtherBC) {
1095 Int_t bcOther = other->GetBC();
1096 if (bcOther != AliVTrack::kTOFBCNA && TMath::Abs(bcOther-bcPrimary) > 2)
1097 // Pile-up from other BC
1098 return true;
1099 }
1100 // Calculate the distance
1101 Double_t d = -1;
1102 Double_t dx = primary->GetX() - other->GetX();
1103 Double_t dy = primary->GetY() - other->GetY();
1104 Double_t dz = primary->GetZ() - other->GetZ();
1105 Double_t covPrimary[6], covOther[6];
1106 primary->GetCovarianceMatrix(covPrimary);
1107 other->GetCovarianceMatrix(covOther);
1108 TMatrixDSym v(3);
1109 v(0,0) = covPrimary[0] + covOther[0]; // diagonal
1110 v(1,1) = covPrimary[2] + covOther[2]; // diagonal
1111 v(2,2) = covPrimary[5] + covOther[5]; // diagonal
1112 v(1,0) = v(0,1) = covPrimary[1]+covOther[1]; // Band
1113 v(0,2) = v(1,2) = v(2,0) = v(2,1) = 0; // Off-diagonal+band
1114 // Try to invert
1115 v.InvertFast();
1116 if (!v.IsValid())
1117 // Question if kStatus bit is every set after InvertFast?
1118 continue;
1119 d = (v(0,0) * dx * dx + v(1,1) * dy * dy + v(2,2) * dz * dz +
1120 2 * (v(0,1) * dx * dy + v(0,2) * dx * dz + v(1,2) * dy * dz));
1121
1122 if (d < 0 || TMath::Sqrt(d) < fMinPileupDistance)
1123 // Bad distance, or not fare enough from each other
1124 continue;
1125
1126 // Well separated vertices -> pileup
1127 return true;
1128 }
1129 return false;
1130}
1131
241cca4d 1132//____________________________________________________________________
1133Bool_t
1134AliFMDEventInspector::CheckEmpty(const TString& trigStr, UInt_t& triggers) const
1135{
8565b10b 1136 if (trigStr.Contains("CBEAMB-ABCE-NOPF-ALL")) {
1137 triggers |= AliAODForwardMult::kEmpty;
241cca4d 1138 return true;
8565b10b 1139 }
241cca4d 1140 return false;
1141}
1142//____________________________________________________________________
1143Bool_t
1144AliFMDEventInspector::CheckWords(const AliESDEvent& esd, UInt_t& triggers) const
1145{
1146 TObject* word = 0;
1147 TIter nextColl(&fCollWords);
1148 while ((word = nextColl())) {
1149 DMSG(fDebug,10,"Checking if %s trigger %s is fired",
1150 word->GetTitle(), word->GetName());
1151 if (!esd.IsTriggerClassFired(word->GetName())) continue;
8565b10b 1152
241cca4d 1153 TString beamSide = word->GetTitle();
1154 DMSG(fDebug,10,"Found it - this is a %s trigger", beamSide.Data());
e85a76b7 1155
241cca4d 1156 if (!beamSide.EqualTo("B")) continue;
1157 triggers |= AliAODForwardMult::kB;
1158 break; // No more to do here
e85a76b7 1159 }
241cca4d 1160 TIter nextBg(&fBgWords);
1161 UInt_t all = (AliAODForwardMult::kA |
1162 AliAODForwardMult::kC |
1163 AliAODForwardMult::kE);
1164 while ((word = nextBg())) {
1165 DMSG(fDebug,10,"Checking if %s trigger %s is fired",
1166 word->GetTitle(), word->GetName());
1167 if (!esd.IsTriggerClassFired(word->GetName())) continue;
0be6c8cd 1168
241cca4d 1169 TString beamSide = word->GetTitle();
1170 DMSG(fDebug,10,"Found it - this is a %s trigger", beamSide.Data());
e6463868 1171
241cca4d 1172 if (beamSide.Contains("A")) triggers |= AliAODForwardMult::kA;
1173 if (beamSide.Contains("C")) triggers |= AliAODForwardMult::kC;
1174 if (beamSide.Contains("E")) triggers |= AliAODForwardMult::kE;
1175
1176 if ((triggers & all) == all) break; // No more to do
0be6c8cd 1177 }
241cca4d 1178 return true;
8565b10b 1179}
241cca4d 1180
1181
8565b10b 1182//____________________________________________________________________
1183Bool_t
5ca83fee 1184AliFMDEventInspector::ReadVertex(const AliESDEvent& esd, TVector3& ip)
8565b10b 1185{
7984e5f7 1186 //
1187 // Read the vertex information from the ESD event
1188 //
1189 // Parameters:
1190 // esd ESD event
1191 // vz On return, the vertex Z position
1192 //
1193 // Return:
1194 // @c true on success, @c false otherwise
1195 //
6ab100ec 1196 DGUARD(fDebug,2,"Read the vertex in AliFMDEventInspector");
5ca83fee 1197 ip.SetXYZ(1024, 1024, 0);
65abd48b 1198
8449e3e0 1199 EVtxStatus s = kNoVtx;
bfab35d9 1200 if (fUseDisplacedVertices && fDisplacedVertex.IsSatellite()) {
1201 s = kVtxOK;
1202 ip.SetZ(fDisplacedVertex.GetVertexZ());
1203 }
1204 else if (fUseFirstPhysicsVertex)
8449e3e0 1205 s = CheckPWGUDVertex(esd, ip);
1206 else if (fUsepA2012Vertex)
1207 s = CheckpA2012Vertex(esd,ip);
1208 else
1209 s = CheckVertex(esd, ip);
1210
8449e3e0 1211 fHVtxStatus->Fill(s);
1212
1213 return s == kVtxOK;
241cca4d 1214}
1215
1216//____________________________________________________________________
8449e3e0 1217AliFMDEventInspector::EVtxStatus
241cca4d 1218AliFMDEventInspector::CheckPWGUDVertex(const AliESDEvent& esd,
5ca83fee 1219 TVector3& ip) const
241cca4d 1220{
1221 // This is the code used by the 1st physics people
1222 const AliESDVertex* vertex = esd.GetPrimaryVertex();
1223 if (!vertex || !vertex->GetStatus()) {
1224 DMSG(fDebug,2,"No primary vertex (%p) or bad status %d",
1225 vertex, (vertex ? vertex->GetStatus() : -1));
8449e3e0 1226 return kNoVtx;
241cca4d 1227 }
1228 const AliESDVertex* vertexSPD = esd.GetPrimaryVertexSPD();
1229 if (!vertexSPD || !vertexSPD->GetStatus()) {
1230 DMSG(fDebug,2,"No primary SPD vertex (%p) or bad status %d",
1231 vertexSPD, (vertexSPD ? vertexSPD->GetStatus() : -1));
8449e3e0 1232 return kNoSPDVtx;
241cca4d 1233 }
96110c91 1234
241cca4d 1235 // if vertex is from SPD vertexZ, require more stringent cuts
1236 if (vertex->IsFromVertexerZ()) {
1237 if (vertex->GetDispersion() > fMaxVzErr ||
1238 vertex->GetZRes() > 1.25 * fMaxVzErr) {
1239 DMSG(fDebug,2,"Dispersion %f > %f or resolution %f > %f",
1240 vertex->GetDispersion(), fMaxVzErr,
1241 vertex->GetZRes(), 1.25 * fMaxVzErr);
8449e3e0 1242 return kUncertain;
96110c91 1243 }
96110c91 1244 }
5ca83fee 1245 ip.SetZ(vertex->GetZ());
241cca4d 1246
1247 if(!vertex->IsFromVertexerZ()) {
5ca83fee 1248 ip.SetX(vertex->GetX());
1249 ip.SetY(vertex->GetY());
241cca4d 1250 }
8449e3e0 1251 return kVtxOK;
241cca4d 1252}
8449e3e0 1253//____________________________________________________________________
1254AliFMDEventInspector::EVtxStatus
1255AliFMDEventInspector::CheckpA2012Vertex(const AliESDEvent& esd,
1256 TVector3& ip) const
ffd35d33 1257{
8449e3e0 1258 const AliESDVertex *vertex = esd.GetPrimaryVertexSPD();
1259 if (!vertex) return kNoSPDVtx;
1260 if (vertex->GetNContributors() <= 0) return kFewContrib;
1261
1262 TString vtxTyp = vertex->GetTitle();
1263 if (vtxTyp.Contains("vertexer: Z")) return kNotVtxZ;
1264
1265 if (vertex->GetDispersion() >= 0.04 || vertex->GetZRes()>=0.25)
1266 return kUncertain;
1267
1268 ip.SetX(vertex->GetX());
1269 ip.SetY(vertex->GetY());
1270 ip.SetZ(vertex->GetZ());
1271
1272 return kVtxOK;
ffd35d33 1273}
1274
241cca4d 1275//____________________________________________________________________
8449e3e0 1276AliFMDEventInspector::EVtxStatus
241cca4d 1277AliFMDEventInspector::CheckVertex(const AliESDEvent& esd,
5ca83fee 1278 TVector3& ip) const
241cca4d 1279{
1280 // Use standard SPD vertex (perhaps preferable for Pb+Pb)
1281 // Get the vertex
1282 const AliESDVertex* vertex = esd.GetPrimaryVertexSPD();
1283 if (!vertex) {
1284 if (fDebug > 2) {
1285 AliWarning("No SPD vertex found in ESD"); }
8449e3e0 1286 return kNoSPDVtx;
241cca4d 1287 }
e83d0620 1288
8449e3e0 1289 // #if 0 // Check disabled - seem to kill a lot of PbPb events
241cca4d 1290 // Check that enough tracklets contributed
1291 if(vertex->GetNContributors() <= 0) {
1292 DMSG(fDebug,2,"Number of contributors to vertex is %d<=0",
1293 vertex->GetNContributors());
5ca83fee 1294 ip.SetZ(0);
8449e3e0 1295 return kFewContrib;
241cca4d 1296 }
8449e3e0 1297 // #endif
1298
241cca4d 1299 // Check that the uncertainty isn't too large
1300 if (vertex->GetZRes() > fMaxVzErr) {
1301 DMSG(fDebug,2,"Uncertaintity in Z of vertex is too large %f > %f",
1302 vertex->GetZRes(), fMaxVzErr);
8449e3e0 1303 return kUncertain;
241cca4d 1304 }
e83d0620 1305
241cca4d 1306 // Get the z coordiante
5ca83fee 1307 ip.SetZ(vertex->GetZ());
241cca4d 1308 const AliESDVertex* vertexXY = esd.GetPrimaryVertex();
96110c91 1309
5ca83fee 1310
241cca4d 1311 if(!vertexXY->IsFromVertexerZ()) {
5ca83fee 1312 ip.SetX(vertexXY->GetX());
1313 ip.SetY(vertexXY->GetY());
241cca4d 1314 }
8449e3e0 1315 return kVtxOK;
8565b10b 1316}
241cca4d 1317
0bd4b00f 1318//____________________________________________________________________
1319Bool_t
1320AliFMDEventInspector::ReadRunDetails(const AliESDEvent* esd)
1321{
7984e5f7 1322 //
1323 // Read the collision system, collision energy, and L3 field setting
1324 // from the ESD
1325 //
1326 // Parameters:
1327 // esd ESD to get information from
1328 //
1329 // Return:
1330 // true on success, false
1331 //
cc83fca2 1332 // AliInfo(Form("Parameters from 1st ESD event: cms=%s, sNN=%f, field=%f",
1333 // esd->GetBeamType(), 2*esd->GetBeamEnergy(),
1334 // esd->GetMagneticField()));
6ab100ec 1335 DGUARD(fDebug,2,"Read the run details in AliFMDEventInspector");
d4d486f8 1336 const char* sys = esd->GetBeamType();
1337 Float_t cms = 2 * esd->GetBeamEnergy();
1338 Float_t fld = esd->GetMagneticField();
1339 fCollisionSystem = AliForwardUtil::ParseCollisionSystem(sys);
1340 fEnergy = AliForwardUtil::ParseCenterOfMassEnergy(fCollisionSystem,
1341 cms);
1342 fField = AliForwardUtil::ParseMagneticField(fld);
8449e3e0 1343 fRunNumber = esd->GetRunNumber();
c8b1a7db 1344 StoreProduction();
8449e3e0 1345 StoreInformation();
d4d486f8 1346 if (fCollisionSystem == AliForwardUtil::kUnknown) {
1347 AliWarningF("Unknown collision system: %s - please check", sys);
1348 return false;
1349 }
1350 if (fEnergy <= 0) {
1351 AliWarningF("Unknown CMS energy: %f (%d) - please check", cms, fEnergy);
1352 return false;
1353 }
1354 if (TMath::Abs(fField) > 10) {
1355 AliWarningF("Unknown L3 field setting: %f (%d) - please check", fld,fField);
1356 return false;
1357 }
0bd4b00f 1358
d4d486f8 1359 return true;
0bd4b00f 1360}
1361
241cca4d 1362
1363//____________________________________________________________________
1364const Char_t*
1365AliFMDEventInspector::CodeString(UInt_t code)
1366{
1367 static TString s;
1368 s = "";
1369 if (code & kNoEvent) s.Append("NOEVENT ");
1370 if (code & kNoTriggers) s.Append("NOTRIGGERS ");
1371 if (code & kNoSPD) s.Append("NOSPD ");
1372 if (code & kNoFMD) s.Append("NOFMD ");
1373 if (code & kNoVertex) s.Append("NOVERTEX ");
1374 if (code & kBadVertex) s.Append("BADVERTEX ");
1375 return s.Data();
1376}
c8b1a7db 1377#define PF(N,V,...) \
1378 AliForwardUtil::PrintField(N,V, ## __VA_ARGS__)
1379#define PFB(N,FLAG) \
1380 do { \
1381 AliForwardUtil::PrintName(N); \
1382 std::cout << std::boolalpha << (FLAG) << std::noboolalpha << std::endl; \
1383 } while(false)
1384#define PFV(N,VALUE) \
1385 do { \
1386 AliForwardUtil::PrintName(N); \
1387 std::cout << (VALUE) << std::endl; } while(false)
1388
0bd4b00f 1389//____________________________________________________________________
1390void
1391AliFMDEventInspector::Print(Option_t*) const
1392{
7984e5f7 1393 //
1394 // Print information
1395 //
1396 // option Not used
1397 //
c8b1a7db 1398 AliForwardUtil::PrintTask(*this);
0bd4b00f 1399 TString sNN(AliForwardUtil::CenterOfMassEnergyString(fEnergy));
1400 sNN.Strip(TString::kBoth, '0');
1401 sNN.ReplaceAll("GeV", " GeV");
1402 TString field(AliForwardUtil::MagneticFieldString(fField));
1403 field.ReplaceAll("p", "+");
1404 field.ReplaceAll("m", "-");
1405 field.ReplaceAll("kG", " kG");
c8b1a7db 1406
1407 gROOT->IncreaseDirLevel();
0b7de667 1408 PFV("Production", "");
1409 gROOT->IncreaseDirLevel();
1410 PF("Period", "LHC%02d%c", fProdYear, fProdLetter);
1411 PFB("MC anchor", fProdMC);
1412 PFV("AliROOT Revision", fProdSVN);
1413 gROOT->DecreaseDirLevel();
c8b1a7db 1414 PFV("Vertex bins", fVtxAxis.GetNbins());
1415 PF("Vertex Range", "[%+6.1f,%+6.1f]", fVtxAxis.GetXmin(), fVtxAxis.GetXmax());
1416 PFV("Low flux cut", fLowFluxCut );
1417 PFV("Max(delta v_z)", fMaxVzErr );
1418 PFV("Min(nContrib_pileup)", fMinPileupContrib );
1419 PFV("Min(v-pileup)", fMinPileupDistance );
1420 PFV("System", AliForwardUtil::CollisionSystemString(fCollisionSystem));
1421 PFV("CMS energy per nucleon", sNN);
1422 PFV("Field", field);
1423 PFB("Satellite events", fUseDisplacedVertices);
1424 PFB("Use 2012 pA vertex", fUsepA2012Vertex );
1425 PFB("Use PWG-UD vertex", fUseFirstPhysicsVertex);
1426 PFB("Simulation input", fMC );
1427 PFV("Centrality method", fCentMethod);
1428 PFV("Centrality axis", (!fCentAxis ? "none" : ""));
1429 if (!fCentAxis) {return; }
d8244e9e 1430 Int_t nBin = fCentAxis->GetNbins();
d8244e9e 1431 for (Int_t i = 0; i < nBin; i++) {
0b7de667 1432 if (i != 0 && (i % 10) == 0) std::cout << '\n';
1433 if ((i % 10) == 0) std::cout << " ";
d8244e9e 1434 std::cout << std::setw(5) << fCentAxis->GetBinLowEdge(i+1) << '-';
1435 }
1436 std::cout << std::setw(5) << fCentAxis->GetBinUpEdge(nBin) << std::endl;
c8b1a7db 1437 gROOT->DecreaseDirLevel();
0bd4b00f 1438}
1439
1440
8565b10b 1441//
1442// EOF
1443//
1444