- coverity fix 10184
[u/mrichter/AliRoot.git] / PWG2 / 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"
18#include "AliLog.h"
19#include "AliESDEvent.h"
20#include "AliMultiplicity.h"
21#include "AliAnalysisManager.h"
22#include "AliInputEventHandler.h"
23#include "AliTriggerAnalysis.h"
24#include "AliPhysicsSelection.h"
25#include "AliAODForwardMult.h"
0bd4b00f 26#include "AliForwardUtil.h"
5e4d905e 27#include "AliCentrality.h"
8565b10b 28#include <TH1.h>
29#include <TList.h>
30#include <TDirectory.h>
0bd4b00f 31#include <TROOT.h>
32#include <iostream>
33#include <iomanip>
e58000b7 34#include "AliMCEvent.h"
35#include "AliGenPythiaEventHeader.h"
36#include "AliGenDPMjetEventHeader.h"
37#include "AliHeader.h"
38#include "AliMCEventHandler.h"
8565b10b 39//====================================================================
40AliFMDEventInspector::AliFMDEventInspector()
41 : TNamed(),
42 fHEventsTr(0),
43 fHEventsTrVtx(0),
44 fHTriggers(0),
0bd4b00f 45 fHType(0),
fe52e455 46 fHWords(0),
5e4d905e 47 fHCent(0),
8565b10b 48 fLowFluxCut(1000),
9d05ffeb 49 fMaxVzErr(0.2),
8565b10b 50 fList(0),
0bd4b00f 51 fEnergy(0),
52 fField(999),
53 fCollisionSystem(kUnknown),
8565b10b 54 fDebug(0)
55{
7984e5f7 56 //
57 // Constructor
58 //
8565b10b 59}
60
61//____________________________________________________________________
62AliFMDEventInspector::AliFMDEventInspector(const char* name)
63 : TNamed("fmdEventInspector", name),
64 fHEventsTr(0),
65 fHEventsTrVtx(0),
66 fHTriggers(0),
0bd4b00f 67 fHType(0),
fe52e455 68 fHWords(0),
5e4d905e 69 fHCent(0),
8565b10b 70 fLowFluxCut(1000),
9d05ffeb 71 fMaxVzErr(0.2),
8565b10b 72 fList(0),
0bd4b00f 73 fEnergy(0),
74 fField(999),
75 fCollisionSystem(kUnknown),
8565b10b 76 fDebug(0)
77{
7984e5f7 78 //
79 // Constructor
80 //
81 // Parameters:
82 // name Name of object
83 //
8565b10b 84}
85
86//____________________________________________________________________
87AliFMDEventInspector::AliFMDEventInspector(const AliFMDEventInspector& o)
88 : TNamed(o),
89 fHEventsTr(o.fHEventsTr),
90 fHEventsTrVtx(o.fHEventsTrVtx),
91 fHTriggers(o.fHTriggers),
0bd4b00f 92 fHType(o.fHType),
fe52e455 93 fHWords(o.fHWords),
5e4d905e 94 fHCent(o.fHCent),
6feacd76 95 fLowFluxCut(o.fLowFluxCut),
8565b10b 96 fMaxVzErr(o.fMaxVzErr),
97 fList(o.fList),
0bd4b00f 98 fEnergy(o.fEnergy),
99 fField(o.fField),
100 fCollisionSystem(o.fCollisionSystem),
8565b10b 101 fDebug(0)
102{
7984e5f7 103 //
104 // Copy constructor
105 //
106 // Parameters:
107 // o Object to copy from
108 //
8565b10b 109}
110
111//____________________________________________________________________
112AliFMDEventInspector::~AliFMDEventInspector()
113{
7984e5f7 114 //
115 // Destructor
116 //
8565b10b 117 if (fHEventsTr) delete fHEventsTr;
118 if (fHEventsTrVtx) delete fHEventsTrVtx;
119 if (fHTriggers) delete fHTriggers;
0bd4b00f 120 if (fHType) delete fHType;
fe52e455 121 if (fHWords) delete fHWords;
5e4d905e 122 if (fHCent) delete fHCent;
8565b10b 123 if (fList) delete fList;
124}
125//____________________________________________________________________
126AliFMDEventInspector&
127AliFMDEventInspector::operator=(const AliFMDEventInspector& o)
128{
7984e5f7 129 //
130 // Assignement operator
131 //
132 // Parameters:
133 // o Object to assign from
134 //
135 // Return:
136 // Reference to this object
137 //
8565b10b 138 TNamed::operator=(o);
139 fHEventsTr = o.fHEventsTr;
140 fHEventsTrVtx = o.fHEventsTrVtx;
141 fHTriggers = o.fHTriggers;
0bd4b00f 142 fHType = o.fHType;
fe52e455 143 fHWords = o.fHWords;
5e4d905e 144 fHCent = o.fHCent;
8565b10b 145 fLowFluxCut = o.fLowFluxCut;
146 fMaxVzErr = o.fMaxVzErr;
147 fDebug = o.fDebug;
148 fList = (o.fList ? new TList : 0);
0bd4b00f 149 fEnergy = o.fEnergy;
150 fField = o.fField;
151 fCollisionSystem = o.fCollisionSystem;
8565b10b 152 if (fList) {
153 fList->SetName(GetName());
154 if (fHEventsTr) fList->Add(fHEventsTr);
155 if (fHEventsTrVtx) fList->Add(fHEventsTrVtx);
156 if (fHTriggers) fList->Add(fHTriggers);
0bd4b00f 157 if (fHType) fList->Add(fHType);
fe52e455 158 if (fHWords) fList->Add(fHWords);
5e4d905e 159 if (fHCent) fList->Add(fHCent);
8565b10b 160 }
161 return *this;
162}
163
164//____________________________________________________________________
165Bool_t
fb3430ac 166AliFMDEventInspector::FetchHistograms(const TList* d,
8565b10b 167 TH1I*& hEventsTr,
168 TH1I*& hEventsTrVtx,
169 TH1I*& hTriggers) const
170{
7984e5f7 171 //
172 // Fetch our histograms from the passed list
173 //
174 // Parameters:
175 // d Input
176 // hEventsTr On return, pointer to histogram, or null
177 // hEventsTrVtx On return, pointer to histogram, or null
178 // hTriggers On return, pointer to histogram, or null
179 //
180 // Return:
181 // true on success, false otherwise
182 //
8565b10b 183 hEventsTr = 0;
184 hEventsTrVtx = 0;
185 hTriggers = 0;
186 TList* dd = dynamic_cast<TList*>(d->FindObject(GetName()));
187 if (!dd) return kFALSE;
188
189 hEventsTr = dynamic_cast<TH1I*>(dd->FindObject("nEventsTr"));
190 hEventsTrVtx = dynamic_cast<TH1I*>(dd->FindObject("nEventsTrVtx"));
191 hTriggers = dynamic_cast<TH1I*>(dd->FindObject("triggers"));
192
193 if (!hEventsTr || !hEventsTrVtx || !hTriggers) return kFALSE;
194 return kTRUE;
195}
196//____________________________________________________________________
197void
198AliFMDEventInspector::Init(const TAxis& vtxAxis)
199{
7984e5f7 200 //
201 // Initialize the object
202 //
203 // Parameters:
204 // vtxAxis Vertex axis in use
205 //
8565b10b 206 fHEventsTr = new TH1I("nEventsTr", "Number of events w/trigger",
207 vtxAxis.GetNbins(),
208 vtxAxis.GetXmin(),
209 vtxAxis.GetXmax());
210 fHEventsTr->SetXTitle("v_{z} [cm]");
211 fHEventsTr->SetYTitle("# of events");
212 fHEventsTr->SetFillColor(kRed+1);
213 fHEventsTr->SetFillStyle(3001);
214 fHEventsTr->SetDirectory(0);
215 // fHEventsTr->Sumw2();
216 fList->Add(fHEventsTr);
217
218 fHEventsTrVtx = new TH1I("nEventsTrVtx",
219 "Number of events w/trigger and vertex",
220 vtxAxis.GetNbins(),
221 vtxAxis.GetXmin(),
222 vtxAxis.GetXmax());
223 fHEventsTrVtx->SetXTitle("v_{z} [cm]");
224 fHEventsTrVtx->SetYTitle("# of events");
225 fHEventsTrVtx->SetFillColor(kBlue+1);
226 fHEventsTrVtx->SetFillStyle(3001);
227 fHEventsTrVtx->SetDirectory(0);
228 // fHEventsTrVtx->Sumw2();
229 fList->Add(fHEventsTrVtx);
230
231
232 fHTriggers = new TH1I("triggers", "Triggers", 10, 0, 10);
233 fHTriggers->SetFillColor(kRed+1);
234 fHTriggers->SetFillStyle(3001);
235 fHTriggers->SetStats(0);
236 fHTriggers->SetDirectory(0);
237 fHTriggers->GetXaxis()->SetBinLabel(kInel +1,"INEL");
238 fHTriggers->GetXaxis()->SetBinLabel(kInelGt0+1,"INEL>0");
239 fHTriggers->GetXaxis()->SetBinLabel(kNSD +1,"NSD");
240 fHTriggers->GetXaxis()->SetBinLabel(kEmpty +1,"Empty");
241 fHTriggers->GetXaxis()->SetBinLabel(kA +1,"A");
242 fHTriggers->GetXaxis()->SetBinLabel(kB +1,"B");
243 fHTriggers->GetXaxis()->SetBinLabel(kC +1,"C");
244 fHTriggers->GetXaxis()->SetBinLabel(kE +1,"E");
e58000b7 245 fHTriggers->GetXaxis()->SetBinLabel(kPileUp +1,"Pileup");
246 fHTriggers->GetXaxis()->SetBinLabel(kMCNSD +1,"nsd");
8565b10b 247 fList->Add(fHTriggers);
0bd4b00f 248
249 fHType = new TH1I("type", Form("Event type (cut: SPD mult>%d)",
250 fLowFluxCut), 2, -.5, 1.5);
251 fHType->SetFillColor(kRed+1);
252 fHType->SetFillStyle(3001);
253 fHType->SetStats(0);
254 fHType->SetDirectory(0);
255 fHType->GetXaxis()->SetBinLabel(1,"Low-flux");
256 fHType->GetXaxis()->SetBinLabel(2,"High-flux");
257 fList->Add(fHType);
fe52e455 258
259
260 fHWords = new TH1I("words", "Trigger words seen", 1, 0, 0);
261 fHWords->SetFillColor(kBlue+1);
262 fHWords->SetFillStyle(3001);
263 fHWords->SetStats(0);
264 fHWords->SetDirectory(0);
265 fHWords->SetBit(TH1::kCanRebin);
266 fList->Add(fHWords);
5e4d905e 267
268 fHCent = new TH1F("cent", "Centrality", 101, -1.5, 100.5);
269 fHCent->SetFillColor(kBlue+1);
270 fHCent->SetFillStyle(3001);
271 fHCent->SetStats(0);
272 fHCent->SetDirectory(0);
273 fHCent->SetXTitle("Centrality [%]");
274 fHCent->SetYTitle("Events");
275 fList->Add(fHCent);
8565b10b 276}
277
278//____________________________________________________________________
279void
280AliFMDEventInspector::DefineOutput(TList* dir)
281{
7984e5f7 282 //
283 // Define the output histograms. These are put in a sub list of the
284 // passed list. The histograms are merged before the parent task calls
285 // AliAnalysisTaskSE::Terminate
286 //
287 // dir Directory to add to
288 //
8565b10b 289 fList = new TList;
290 fList->SetName(GetName());
291 dir->Add(fList);
292}
293
294//____________________________________________________________________
295UInt_t
296AliFMDEventInspector::Process(const AliESDEvent* event,
297 UInt_t& triggers,
298 Bool_t& lowFlux,
0bd4b00f 299 UShort_t& ivz,
5e4d905e 300 Double_t& vz,
301 Double_t& cent)
8565b10b 302{
7984e5f7 303 //
304 // Process the event
305 //
306 // Parameters:
307 // event Input event
308 // triggers On return, the triggers fired
309 // lowFlux On return, true if the event is considered a low-flux
310 // event (according to the setting of fLowFluxCut)
311 // ivz On return, the found vertex bin (1-based). A zero
312 // means outside of the defined vertex range
313 // vz On return, the z position of the interaction
5e4d905e 314 // cent On return, the centrality - if not available < 0
7984e5f7 315 //
316 // Return:
317 // 0 (or kOk) on success, otherwise a bit mask of error codes
318 //
319
8565b10b 320 // Check that we have an event
321 if (!event) {
322 AliWarning("No ESD event found for input event");
323 return kNoEvent;
324 }
325
326 // Read trigger information from the ESD and store in AOD object
327 if (!ReadTriggers(event, triggers)) {
0bd4b00f 328 if (fDebug > 2) {
329 AliWarning("Failed to read triggers from ESD"); }
8565b10b 330 return kNoTriggers;
331 }
332
333 // Check if this is a high-flux event
334 const AliMultiplicity* testmult = event->GetMultiplicity();
335 if (!testmult) {
0bd4b00f 336 if (fDebug > 3) {
337 AliWarning("No central multiplicity object found"); }
8565b10b 338 return kNoSPD;
339 }
340 lowFlux = testmult->GetNumberOfTracklets() < fLowFluxCut;
5e4d905e 341
0bd4b00f 342 fHType->Fill(lowFlux ? 0 : 1);
5e4d905e 343
e58000b7 344 cent = -10;
5e4d905e 345 AliCentrality* centObj = const_cast<AliESDEvent*>(event)->GetCentrality();
e58000b7 346 if (centObj) {
9d05ffeb 347 // AliInfo(Form("Got centrality object %p with quality %d",
348 // centObj, centObj->GetQuality()));
349 // centObj->Print();
e58000b7 350 cent = centObj->GetCentralityPercentileUnchecked("V0M");
351 }
9d05ffeb 352 // AliInfo(Form("Centrality is %f", cent));
5e4d905e 353 fHCent->Fill(cent);
8565b10b 354
355 // Check the FMD ESD data
356 if (!event->GetFMDData()) {
0bd4b00f 357 if (fDebug > 3) {
358 AliWarning("No FMD data found in ESD"); }
8565b10b 359 return kNoFMD;
360 }
361
362 // Get the vertex information
363 vz = 0;
364 Bool_t vzOk = ReadVertex(event, vz);
365
366 fHEventsTr->Fill(vz);
367 if (!vzOk) {
0bd4b00f 368 if (fDebug > 3) {
369 AliWarning("Failed to read vertex from ESD"); }
8565b10b 370 return kNoVertex;
371 }
372 fHEventsTrVtx->Fill(vz);
373
374 // Get the vertex bin
0bd4b00f 375 ivz = fHEventsTr->GetXaxis()->FindBin(vz);
376 if (ivz <= 0 || ivz > fHEventsTr->GetXaxis()->GetNbins()) {
377 if (fDebug > 3) {
8565b10b 378 AliWarning(Form("Vertex @ %f outside of range [%f,%f]",
379 vz, fHEventsTr->GetXaxis()->GetXmin(),
0bd4b00f 380 fHEventsTr->GetXaxis()->GetXmax())); }
381 ivz = 0;
8565b10b 382 return kBadVertex;
383 }
e58000b7 384
385
8565b10b 386 return kOk;
387}
388
389//____________________________________________________________________
390Bool_t
391AliFMDEventInspector::ReadTriggers(const AliESDEvent* esd, UInt_t& triggers)
392{
7984e5f7 393 //
394 // Read the trigger information from the ESD event
395 //
396 // Parameters:
397 // esd ESD event
398 // triggers On return, contains the trigger bits
399 //
400 // Return:
401 // @c true on success, @c false otherwise
402 //
8565b10b 403 triggers = 0;
404
405 // Get the analysis manager - should always be there
406 AliAnalysisManager* am = AliAnalysisManager::GetAnalysisManager();
407 if (!am) {
408 AliWarning("No analysis manager defined!");
409 return kFALSE;
410 }
411
412 // Get the input handler - should always be there
413 AliInputEventHandler* ih =
414 static_cast<AliInputEventHandler*>(am->GetInputEventHandler());
415 if (!ih) {
416 AliWarning("No input handler");
417 return kFALSE;
418 }
e58000b7 419
420 AliMCEventHandler* mcHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
421 AliMCEvent* mcEvent = 0;
422 if(mcHandler)
423 mcEvent = mcHandler->MCEvent();
424
425 if(mcEvent) {
426 //Assign MC only triggers : True NSD etc.
427 AliHeader* header = mcEvent->Header();
428 AliGenEventHeader* genHeader = header->GenEventHeader();
429
430 AliGenPythiaEventHeader* pythiaGenHeader = dynamic_cast<AliGenPythiaEventHeader*>(genHeader);
431 AliGenDPMjetEventHeader* dpmHeader = dynamic_cast<AliGenDPMjetEventHeader*>(header->GenEventHeader());
432 Bool_t sd = kFALSE;
433 if(pythiaGenHeader) {
434 Int_t pythiaType = pythiaGenHeader->ProcessType();
435 if(pythiaType==92 || pythiaType==93)
436 sd = kTRUE;
437 }
438 if(dpmHeader) {
439 Int_t processType = dpmHeader->ProcessType();
440 if(processType == 5 || processType == 6)
441 sd = kTRUE;
442
443 }
444 if(!sd) {
445 triggers |= AliAODForwardMult::kMCNSD;
446 fHTriggers->Fill(kMCNSD+0.5);
447 }
448
449 }
450
e333578d 451 // Check if this is a collision candidate (INEL)
452 // Note, that we should use the value cached in the input
453 // handler rather than calling IsCollisionCandiate directly
454 // on the AliPhysicsSelection obejct. If we called the latter
455 // then the AliPhysicsSelection object would overcount by a
456 // factor of 2! :-(
457 Bool_t inel = ih->IsEventSelected();
8565b10b 458 if (inel) {
459 triggers |= AliAODForwardMult::kInel;
460 fHTriggers->Fill(kInel+0.5);
461 }
462
e333578d 463 // If this is inel, see if we have a tracklet
8565b10b 464 if (inel) {
465 const AliMultiplicity* spdmult = esd->GetMultiplicity();
466 if (!spdmult) {
467 AliWarning("No SPD multiplicity");
468 }
469 else {
e333578d 470 // Check if we have one or more tracklets
471 // in the range -1 < eta < 1 to set the INEL>0
472 // trigger flag.
8565b10b 473 Int_t n = spdmult->GetNumberOfTracklets();
474 for (Int_t j = 0; j < n; j++) {
475 if(TMath::Abs(spdmult->GetEta(j)) < 1) {
476 triggers |= AliAODForwardMult::kInelGt0;
477 fHTriggers->Fill(kInelGt0+.5);
478 break;
479 }
480 }
481 }
482 }
483
484 // Analyse some trigger stuff
485 AliTriggerAnalysis ta;
486 if (ta.IsOfflineTriggerFired(esd, AliTriggerAnalysis::kNSD1)) {
487 triggers |= AliAODForwardMult::kNSD;
488 fHTriggers->Fill(kNSD+.5);
489 }
e58000b7 490 //Check pileup
491 Bool_t pileup = esd->IsPileupFromSPD(3,0.8);
492 if (pileup) {
493 triggers |= AliAODForwardMult::kPileUp;
494 fHTriggers->Fill(kPileUp+.5);
495 }
496
8565b10b 497 // Get trigger stuff
498 TString trigStr = esd->GetFiredTriggerClasses();
fe52e455 499 // AliWarning(Form("Fired trigger classes: %s", trigStr.Data()));
500 fHWords->Fill(trigStr.Data(), 1);
501#if 0
502 if (trigStr.Contains("MB1") || trigStr.Contains("MBBG3"))
503 triggers |= AliAOODForwardMult::kB;
504 if (trigStr.Contains("COTA"))
505 triggers |= AliAODForwardMult::kA;
506 if (trigStr.Contains("COTC"))
507 triggers |= AliAODForwardMult::kC;
508#endif
8565b10b 509 if (trigStr.Contains("CBEAMB-ABCE-NOPF-ALL")) {
510 triggers |= AliAODForwardMult::kEmpty;
511 fHTriggers->Fill(kEmpty+.5);
512 }
513
514 if (trigStr.Contains("CINT1A-ABCE-NOPF-ALL")) {
515 triggers |= AliAODForwardMult::kA;
516 fHTriggers->Fill(kA+.5);
517 }
518
519 if (trigStr.Contains("CINT1B-ABCE-NOPF-ALL")) {
520 triggers |= AliAODForwardMult::kB;
521 fHTriggers->Fill(kB+.5);
522 }
523
524
525 if (trigStr.Contains("CINT1C-ABCE-NOPF-ALL")) {
526 triggers |= AliAODForwardMult::kC;
527 fHTriggers->Fill(kC+.5);
528 }
529
530 if (trigStr.Contains("CINT1-E-NOPF-ALL")) {
531 triggers |= AliAODForwardMult::kE;
532 fHTriggers->Fill(kE+.5);
533 }
534
535 return kTRUE;
536}
537//____________________________________________________________________
538Bool_t
539AliFMDEventInspector::ReadVertex(const AliESDEvent* esd, Double_t& vz)
540{
7984e5f7 541 //
542 // Read the vertex information from the ESD event
543 //
544 // Parameters:
545 // esd ESD event
546 // vz On return, the vertex Z position
547 //
548 // Return:
549 // @c true on success, @c false otherwise
550 //
8565b10b 551 vz = 0;
552 // Get the vertex
553 const AliESDVertex* vertex = esd->GetPrimaryVertexSPD();
554 if (!vertex) {
0bd4b00f 555 if (fDebug > 2) {
556 AliWarning("No SPD vertex found in ESD"); }
8565b10b 557 return kFALSE;
558 }
559
560 // Check that enough tracklets contributed
561 if(vertex->GetNContributors() <= 0) {
0bd4b00f 562 if (fDebug > 2) {
8565b10b 563 AliWarning(Form("Number of contributors to vertex is %d<=0",
0bd4b00f 564 vertex->GetNContributors())); }
8565b10b 565 vz = 0;
566 return kFALSE;
567 }
568
569 // Check that the uncertainty isn't too large
570 if (vertex->GetZRes() > fMaxVzErr) {
0bd4b00f 571 if (fDebug > 2) {
8565b10b 572 AliWarning(Form("Uncertaintity in Z of vertex is too large %f > %f",
0bd4b00f 573 vertex->GetZRes(), fMaxVzErr)); }
8565b10b 574 return kFALSE;
575 }
576
e58000b7 577 // HHD test
578 /*
579 if (vertex->IsFromVertexerZ()) {
580 return kFALSE;
581 }
582 if (TMath::Sqrt(TMath::Power(vertex->GetX(),2) + TMath::Power(vertex->GetY(),2)) > 3 ) {
583 return kFALSE;
584 }
585 */
586
8565b10b 587 // Get the z coordiante
588 vz = vertex->GetZ();
589 return kTRUE;
590}
0bd4b00f 591
592//____________________________________________________________________
593Bool_t
594AliFMDEventInspector::ReadRunDetails(const AliESDEvent* esd)
595{
7984e5f7 596 //
597 // Read the collision system, collision energy, and L3 field setting
598 // from the ESD
599 //
600 // Parameters:
601 // esd ESD to get information from
602 //
603 // Return:
604 // true on success, false
605 //
cc83fca2 606 // AliInfo(Form("Parameters from 1st ESD event: cms=%s, sNN=%f, field=%f",
607 // esd->GetBeamType(), 2*esd->GetBeamEnergy(),
608 // esd->GetMagneticField()));
0bd4b00f 609 fCollisionSystem =
610 AliForwardUtil::ParseCollisionSystem(esd->GetBeamType());
611 fEnergy =
612 AliForwardUtil::ParseCenterOfMassEnergy(fCollisionSystem,
613 2 * esd->GetBeamEnergy());
614 fField =
615 AliForwardUtil::ParseMagneticField(esd->GetMagneticField());
616
617 if (fCollisionSystem == AliForwardUtil::kUnknown ||
618 fEnergy <= 0 ||
619 TMath::Abs(fField) > 10)
620 return kFALSE;
621
622 return kTRUE;
623}
624
625//____________________________________________________________________
626void
627AliFMDEventInspector::Print(Option_t*) const
628{
7984e5f7 629 //
630 // Print information
631 //
632 // option Not used
633 //
0bd4b00f 634 char ind[gROOT->GetDirLevel()+1];
635 for (Int_t i = 0; i < gROOT->GetDirLevel(); i++) ind[i] = ' ';
636 ind[gROOT->GetDirLevel()] = '\0';
637 TString sNN(AliForwardUtil::CenterOfMassEnergyString(fEnergy));
638 sNN.Strip(TString::kBoth, '0');
639 sNN.ReplaceAll("GeV", " GeV");
640 TString field(AliForwardUtil::MagneticFieldString(fField));
641 field.ReplaceAll("p", "+");
642 field.ReplaceAll("m", "-");
643 field.ReplaceAll("kG", " kG");
644
645 std::cout << ind << "AliFMDEventInspector: " << GetName() << '\n'
646 << ind << " Low flux cut: " << fLowFluxCut << '\n'
647 << ind << " Max(delta v_z): " << fMaxVzErr << " cm\n"
648 << ind << " System: "
649 << AliForwardUtil::CollisionSystemString(fCollisionSystem) << '\n'
650 << ind << " CMS energy per nucleon: " << sNN << '\n'
651 << ind << " Field: " << field << std::endl;
652}
653
654
8565b10b 655//
656// EOF
657//
658