1 /*************************************************************************
2 * Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
16 // -----------------------------------------------------------------------
17 // In this analysis task we compare multiple AOD filtermasks.
18 // -----------------------------------------------------------------------
19 // Author: Misha Veldhoen (misha.veldhoen@cern.ch)
34 #include "AliAnalysisManager.h"
35 #include "AliAnalysisTaskSE.h"
36 #include "AliInputEventHandler.h"
39 #include "AliPIDResponse.h"
42 #include "AliAODEvent.h"
43 #include "AliAODTrack.h"
44 #include "AliAODHandler.h"
45 #include "AliAODVertex.h"
46 #include "AliAODInputHandler.h"
47 #include "AliAODMCParticle.h"
48 #include "AliAODMCHeader.h"
50 // Includes from own library
51 #include "AliTrackDiHadronPID.h"
52 #include "AliAODTrackCutsDiHadronPID.h"
53 #include "AliAODEventCutsDiHadronPID.h"
55 // AnalysisTask Header
56 #include "AliAnalysisTaskCompareAODTrackCuts.h"
60 ClassImp(AliAnalysisTaskCompareAODTrackCuts);
62 // -----------------------------------------------------------------------
63 AliAnalysisTaskCompareAODTrackCuts::AliAnalysisTaskCompareAODTrackCuts():
69 fCalculateTOFMismatch(kFALSE),
70 fUseMismatchFileFromHomeDir(kTRUE),
71 fUseNSigmaOnPIDAxes(kFALSE),
77 fGlobalPtvsTPCPt(0x0),
78 fLvsEtaProjections(0x0),
79 fCurrentAODEvent(0x0),
80 fCurrentAODTrack(0x0),
81 fCurrentDiHadronPIDTrack(0x0),
82 fGlobalTracksArray(0x0)
86 if (fDebug > 0) {cout << Form("File: %s, Line: %i, Function: %s",__FILE__,__LINE__,__func__) << endl;}
90 // -----------------------------------------------------------------------
91 AliAnalysisTaskCompareAODTrackCuts::AliAnalysisTaskCompareAODTrackCuts(const char* name):
92 AliAnalysisTaskSE(name),
97 fCalculateTOFMismatch(kFALSE),
98 fUseMismatchFileFromHomeDir(kTRUE),
99 fUseNSigmaOnPIDAxes(kFALSE),
102 fInclusiveTimes(0x0),
105 fGlobalPtvsTPCPt(0x0),
106 fLvsEtaProjections(0x0),
107 fCurrentAODEvent(0x0),
108 fCurrentAODTrack(0x0),
109 fCurrentDiHadronPIDTrack(0x0),
110 fGlobalTracksArray(0x0)
114 if (fDebug > 0) {cout << Form("File: %s, Line: %i, Function: %s",__FILE__,__LINE__,__func__) << endl;}
116 // Named Constructor.
117 fTrackCuts = new TObjArray();
118 fTrackCuts->SetName("fTrackCuts");
120 DefineInput(0,TChain::Class());
121 DefineOutput(1, TList::Class());
125 // -----------------------------------------------------------------------
126 AliAnalysisTaskCompareAODTrackCuts::~AliAnalysisTaskCompareAODTrackCuts() {
128 if (fDebug > 0) {cout << Form("File: %s, Line: %i, Function: %s",__FILE__,__LINE__,__func__) << endl;}
132 // -----------------------------------------------------------------------
133 void AliAnalysisTaskCompareAODTrackCuts::UserCreateOutputObjects() {
135 if (fDebug > 0) {cout << Form("File: %s, Line: %i, Function: %s",__FILE__,__LINE__,__func__) << endl;}
137 // Getting a pointer to the analysis manager and input handler.
138 AliAnalysisManager* manager = AliAnalysisManager::GetAnalysisManager();
139 if (!manager) {AliFatal("Could not obtain analysis manager.");}
140 AliInputEventHandler* inputHandler = dynamic_cast<AliInputEventHandler*> (manager->GetInputEventHandler());
141 if (!inputHandler) {AliFatal("Could not obtain input handler."); return;}
143 // Getting the pointer to the PID response object.
144 fPIDResponse = inputHandler->GetPIDResponse();
145 if (!fPIDResponse) {AliFatal("Could not obtain PID response."); return;}
147 // Create the output list.
148 fOutputList = new TList();
149 fOutputList->SetOwner(kTRUE);
151 // Adding Event Cuts to the output.
152 fEventCuts->CreateHistos(); // Generating all requested histograms.
153 fOutputList->Add(fEventCuts);
155 // Adding Track Cuts to the output.
156 for (Int_t iCuts = 0; iCuts < fTrackCuts->GetEntries(); iCuts++) {
157 AliAODTrackCutsDiHadronPID* currentcuts = (AliAODTrackCutsDiHadronPID*)fTrackCuts->At(iCuts);
158 currentcuts->CreateHistos(); // Generating all requested histograms.
159 fOutputList->Add(currentcuts);
162 // Creating inclusive times histogram.
163 fInclusiveTimes = new TH2F("fInclusiveTimes","Inclusive Times;#eta;t (ps)",100,0,1.,1500,10000.,25000.);
164 fOutputList->Add(fInclusiveTimes);
166 // Create the diagram correlating TPC and Global transverse momenta.
167 fGlobalPtvsTPCPt = new TH2F("fGlobalPtvsTPCPt","Global p_{T} vs TPC p_{T}; Global p_{T}; TPC p_{T}",100,0,10,100,0,10);
168 fOutputList->Add(fGlobalPtvsTPCPt);
170 // Creating Global Tracks Array
171 fGlobalTracksArray = new TObjArray();
173 // Loading the appropriate external mismatch histograms.
174 if (fCalculateTOFMismatch) LoadExternalMismatchHistos();
176 PostData(1,fOutputList);
180 // -----------------------------------------------------------------------
181 void AliAnalysisTaskCompareAODTrackCuts::FillGlobalTracksArray() {
183 if (fDebug > 0) {cout << Form("File: %s, Line: %i, Function: %s",__FILE__,__LINE__,__func__) << endl;}
185 // Initialize the mapping for corresponding PID tracks.
186 // See AliAnalysisTaskESDFilter.cxx for explanation about the mapping
187 // between TPC-Only and Global tracks
189 // Create TObjArray for Global Tracks.
191 //cout<<"Global Tracks Array: "<<fGlobalTracksArray<<endl;
193 // Clear previous tracks...
194 fGlobalTracksArray->Clear();
196 //if (fGlobalTracksArray) delete fGlobalTracksArray;
197 //fGlobalTracksArray = new TObjArray();
199 AliAODTrack* track = 0x0;
201 for (Int_t iTrack = 0; iTrack < fCurrentAODEvent->GetNumberOfTracks(); iTrack++) {
203 track = fCurrentAODEvent->GetTrack(iTrack);
204 if (track->GetID()>-1) fGlobalTracksArray->AddAtAndExpand(track,track->GetID());
209 // -----------------------------------------------------------------------
210 AliAODTrack* AliAnalysisTaskCompareAODTrackCuts::GetGlobalTrack(AliAODTrack* track) {
212 // Returns Global Track corresponding to a TPC-Only Track.
213 if (fDebug > 0) {cout << Form("File: %s, Line: %i, Function: %s",__FILE__,__LINE__,__func__) << endl;}
215 if (track->GetID() >= 0) {
216 cout<<"AliAnalysisTaskCompareAODTrackCuts::GetGlobalTrack -> Input Track is not TPC-Only."<<endl;
220 if (!fGlobalTracksArray) {
221 cout<<"AliAnalysisTaskCompareAODTrackCuts::GetGlobalTrack -> Global Tracks Array Does not Exist."<<endl;
225 AliAODTrack* globaltrack = (AliAODTrack*)(fGlobalTracksArray->At(-track->GetID()-1));
227 cout<<"AliAnalysisTaskCompareAODTrackCuts::GetGlobalTrack -> No Global Track Found."<<endl;
235 // -----------------------------------------------------------------------
236 void AliAnalysisTaskCompareAODTrackCuts::UserExec(Option_t*) {
238 if (fDebug > 0) {cout << Form("File: %s, Line: %i, Function: %s",__FILE__,__LINE__,__func__) << endl;}
239 Int_t nmismatched = 0;
242 // Input Current Event.
243 fCurrentAODEvent = dynamic_cast<AliAODEvent*>(InputEvent());
244 if (!fCurrentAODEvent) {
245 cout<<"AliAnalysisTaskCompareAODTrackCuts::UserExec -> AOD Event not found."<<endl;
249 // Check Event Cuts Object.
250 if (!fEventCuts) AliFatal("No Event Cuts Object Found!");
252 // Perform Event Cuts.
253 if (!fEventCuts->IsSelected(fCurrentAODEvent)) {return;}
255 // Check Track Cuts Array.
256 if (!fTrackCuts) AliFatal("No Track Cuts Array Found!");
257 if (fTrackCuts->GetEntries() == 0) AliFatal("Track Cuts Array is Empty!");
259 // If MC, then fill MC reconstructed QA histograms.
260 TClonesArray* mcArray = 0x0;
261 TObjArray* mcArrayLabel = 0x0;
264 mcArray = dynamic_cast<TClonesArray*>(fCurrentAODEvent->FindListObject(AliAODMCParticle::StdBranchName()));
266 AliFatal("No MC array found in the AOD.");
270 mcArrayLabel = new TObjArray(10000);
272 // Loop over MC particles.
273 for (Int_t iParticle = 0; iParticle < mcArray->GetEntriesFast(); iParticle++) {
275 // Put the MC Particle in the Label array.
276 AliAODMCParticle* CurrentAODMCParticle = (AliAODMCParticle*) mcArray->At(iParticle);
277 //cout << "PR: " << CurrentAODMCParticle->IsPhysicalPrimary() << " SW: " << CurrentAODMCParticle->IsSecondaryFromWeakDecay() << " SM: " << CurrentAODMCParticle->IsSecondaryFromMaterial() << endl;
278 mcArrayLabel->AddAtAndExpand(CurrentAODMCParticle,CurrentAODMCParticle->Label());
279 //cout<<"Index: "<<iParticle<<" Label: "<<CurrentAODMCParticle->Label()<<endl;
281 //if (CurrentAODMCParticle->Label()!=iParticle) cout<<"Index unequal to particle's label!"<<endl;
283 // Loop over all Track Cuts.
284 for (Int_t iCuts = 0; iCuts < fTrackCuts->GetEntries(); iCuts++) {
285 ((AliAODTrackCutsDiHadronPID*)fTrackCuts->At(iCuts))->IsSelectedGeneratedMC(CurrentAODMCParticle);
290 //for (Int_t ii=0; ii<200; ii++) cout<<fLvsEtaProj[ii]<<endl;
292 // Create mapping to Global Tracks.
293 FillGlobalTracksArray();
295 // Tell the track cuts object that a new event has started.
296 for (Int_t iCuts = 0; iCuts < fTrackCuts->GetEntries(); iCuts++) {
297 AliAODTrackCutsDiHadronPID* currentcuts = (AliAODTrackCutsDiHadronPID*)fTrackCuts->At(iCuts);
298 if (fVerbose) AliInfo(Form("Starting new event for cuts: %s",currentcuts->GetName()));
299 currentcuts->StartNewEvent();
302 // Loop over all reconstructed tracks.
303 for (Int_t iTrack = 0; iTrack < fCurrentAODEvent->GetNumberOfTracks(); iTrack++) {
305 // Get the Current Track.
306 fCurrentAODTrack = (AliAODTrack*)fCurrentAODEvent->GetTrack(iTrack);
307 if (!fCurrentAODTrack) {
308 cout<<"AliAnalysisTaskCompareAODTrackCuts::UserExec -> AOD Track not found."<<endl;
312 // Ignore muon tracks.
313 if (fCurrentAODTrack->IsMuonTrack()) {continue;}
315 // Copy Track info into AliTrackDiHadronPID object.
316 fCurrentDiHadronPIDTrack = 0x0;
318 // Check if we can indeed find a MC particle corresponding to the reconstructed track.
320 AliAODMCParticle* MCparticleCheck = (AliAODMCParticle*)mcArrayLabel->At(TMath::Abs(fCurrentAODTrack->GetLabel()));
321 Double_t MCPt1 = -999.;
322 if (MCparticleCheck) MCPt1 = MCparticleCheck->Pt();
323 AliAODMCParticle* MCparticleCheck2 = (AliAODMCParticle*)mcArray->At(TMath::Abs(fCurrentAODTrack->GetLabel()));
324 Double_t MCPt2 = -999.;
325 if (MCparticleCheck2) MCPt2 = MCparticleCheck2->Pt();
327 cout<<"AOD track abs label: "<<TMath::Abs(fCurrentAODTrack->GetLabel())<<" MC particle with same label found: "<<(Bool_t)MCparticleCheck<<" and with same index: "<<(Bool_t)MCparticleCheck2<<endl;
328 cout<<"AOD track pt: "<<fCurrentAODTrack->Pt()<<" same label pt: "<<MCPt1<<" same index pt: "<<MCPt2<<endl;
329 //if (!MCparticleCheck) cout<<"MC Particle for a reconstructed track not found."<<endl;
331 // Check if it's a TPC-Only or Global Track.
332 if (fCurrentAODTrack->GetID() < 0) {
333 // Q: Do we really need to create the arraylabel object or is the original MC array already nicely ordered. (NOTE THAT WE TAKE THE MC PARTICLES FROM THE mcArray)
334 if (fIsMC) fCurrentDiHadronPIDTrack = new AliTrackDiHadronPID(fCurrentAODTrack,GetGlobalTrack(fCurrentAODTrack),(AliAODMCParticle*)mcArray->At(TMath::Abs(fCurrentAODTrack->GetLabel())),fPIDResponse);
335 else fCurrentDiHadronPIDTrack = new AliTrackDiHadronPID(fCurrentAODTrack,GetGlobalTrack(fCurrentAODTrack),0x0,fPIDResponse);
337 // Fill histogram of Global p_T vs TPC-only p_T
338 fGlobalPtvsTPCPt->Fill(GetGlobalTrack(fCurrentAODTrack)->Pt(), fCurrentAODTrack->Pt());
341 if (fIsMC) fCurrentDiHadronPIDTrack = new AliTrackDiHadronPID(fCurrentAODTrack,0x0,(AliAODMCParticle*)mcArray->At(TMath::Abs(fCurrentAODTrack->GetLabel())),fPIDResponse);
342 else fCurrentDiHadronPIDTrack = new AliTrackDiHadronPID(fCurrentAODTrack,0x0,0x0,fPIDResponse);
345 if (!fCurrentDiHadronPIDTrack) {
346 cout<<"AliAnalysisTaskCompareAODTrackCuts::UserExec -> Copying to DiHadronPIDTrack failed."<<endl;
350 // Filling random times histogram:
351 ULong_t requestedflags = (UInt_t)(AliAODTrack::kTOFout)|(UInt_t)(AliAODTrack::kTIME);
352 ULong_t trackflags = fCurrentDiHadronPIDTrack->GetFlags();
353 if (requestedflags&trackflags) {fInclusiveTimes->Fill(TMath::Abs(fCurrentDiHadronPIDTrack->Eta()),fCurrentDiHadronPIDTrack->GetTOFsignal());}
355 Double_t rndhittime = -1.e21;
356 if (fCalculateTOFMismatch) rndhittime = GenerateRandomHit(fCurrentDiHadronPIDTrack->Eta());
358 // Loop over all Track Cuts.
359 for (Int_t iCuts = 0; iCuts < fTrackCuts->GetEntries(); iCuts++) {
360 if (fIsMC) ((AliAODTrackCutsDiHadronPID*)fTrackCuts->At(iCuts))->IsSelectedReconstructedMC(fCurrentDiHadronPIDTrack);
361 else ((AliAODTrackCutsDiHadronPID*)fTrackCuts->At(iCuts))->IsSelectedData(fCurrentDiHadronPIDTrack,rndhittime);
364 Int_t tofmatchstat = fCurrentDiHadronPIDTrack->GetTOFMatchingStatus();
365 if (tofmatchstat==0) {nmatched++;}
366 if (tofmatchstat==1) {nmismatched++;}
367 if (tofmatchstat==2) {nnotof++;}
369 // Delete Current DiHadronPIDTrack.
370 delete fCurrentDiHadronPIDTrack;
371 fCurrentDiHadronPIDTrack = 0x0;
375 // Tell the track cuts object that the event has been processed.
376 for (Int_t iCuts = 0; iCuts < fTrackCuts->GetEntries(); iCuts++) {
377 ((AliAODTrackCutsDiHadronPID*)fTrackCuts->At(iCuts))->EventIsDone(fIsMC);
380 //cout << "Matched: "<<nmatched<<" Mismatched: "<<nmismatched<<" No TOF hit: "<<nnotof<<endl;
382 PostData(1,fOutputList);
386 // -----------------------------------------------------------------------
387 Bool_t AliAnalysisTaskCompareAODTrackCuts::LoadExternalMismatchHistos() {
390 // Attempting to load a root file containing information needed
391 // to generate random TOF hits.
394 if (fDebug > 0) {cout << Form("File: %s, Line: %i, Function: %s",__FILE__,__LINE__,__func__) << endl;}
396 // Opening external TOF file.
399 // The default is that the file TOFmismatchHistos.root is taken from the /rootfiles/ directory,
400 // and this works fine when running on the train. When the user submits the jobs himself, he can
401 // choose to not take the file from the home dir, but upload one along with the source code.
402 // If in this case the file is not found, the program tries to get the file from the root directory
404 if (fUseMismatchFileFromHomeDir == kFALSE) {
405 fin = TFile::Open("TOFmismatchHistos.root");
406 if (!fin) {AliWarning("Tried to open uploaded TOFmismatchHistos.root, but failed");}
409 if (!fin) {fin = TFile::Open("alien:///alice/cern.ch/user/m/mveldhoe/rootfiles/TOFmismatchHistos.root");}
412 AliWarning("Couln't open TOFmismatchHistos.root, will not calculate mismatches...");
413 fCalculateTOFMismatch = kFALSE;
416 AliInfo("Sucessfully loaded TOFmismatchHistos.root");
419 // Check if the required histograms are present.
420 TH1F* tmp1 = (TH1F*)fin->Get("hNewT0Fill");
422 AliWarning("Couln't find hNewT0Fill, will not calculate mismatches...");
423 fCalculateTOFMismatch = kFALSE;
426 TH2F* tmp2 = (TH2F*)fin->Get("hLvsEta");
428 AliWarning("Couln't find hLvsEta, will not calculate mismatches...");
429 fCalculateTOFMismatch = kFALSE;
433 // Make a deep copy of the files in the histogram.
434 fT0Fill = (TH1F*)tmp1->Clone("fT0Fill");
435 fLvsEta = (TH2F*)tmp2->Clone("fLvsEta");
437 // Close the external file.
438 AliInfo("Closing external file.");
441 // Creating a TObjArray for LvsEta projections.
442 const Int_t nbinseta = fLvsEta->GetNbinsX();
443 fLvsEtaProjections = new TObjArray(nbinseta);
444 fLvsEtaProjections->SetOwner(kTRUE);
446 // Making the projections needed (excluding underflow/ overflow).
447 for (Int_t iEtaBin = 1; iEtaBin < (nbinseta + 1); iEtaBin++) {
448 TH1F* tmp = (TH1F*)fLvsEta->ProjectionY(Form("LvsEtaProjection_%i",iEtaBin),iEtaBin,iEtaBin);
449 tmp->SetDirectory(0);
450 fLvsEtaProjections->AddAt(tmp,iEtaBin - 1);
457 // -----------------------------------------------------------------------
458 Double_t AliAnalysisTaskCompareAODTrackCuts::GenerateRandomHit(Double_t eta) {
461 // Returns a random TOF time.
464 if (fDebug > 0) {cout << Form("File: %s, Line: %i, Function: %s",__FILE__,__LINE__,__func__) << endl;}
466 // Default (error) value:
467 Double_t rndhittime = -1.e21;
469 // TOF mismatch flag is not turned on.
470 if (!fCalculateTOFMismatch) {
471 AliFatal("Called GenerateRandomHit() method, but flag fCalculateTOFMismatch not set.");
475 // TOF doesn't extend much further than 0.8.
476 if (TMath::Abs(eta) > 0.8) {
477 if (fDebug) {AliInfo("Tried to get a random hit for a track with eta > 0.8.");}
481 // Finding the bin of the eta.
482 TAxis* etaAxis = fLvsEta->GetXaxis();
483 Int_t etaBin = etaAxis->FindBin(eta);
484 if (etaBin == 0 || (etaBin == etaAxis->GetNbins() + 1)) {return rndhittime;}
486 const TH1F* lengthDistribution = (const TH1F*)fLvsEtaProjections->At(etaBin - 1);
488 if (!lengthDistribution) {
489 AliFatal("length Distribution not found.");
493 Double_t currentRndLength = lengthDistribution->GetRandom(); // in cm.
495 // Similar to Roberto's code.
496 Double_t currentRndTime = currentRndLength / (TMath::C() * 1.e2 / 1.e12);
497 Double_t t0fill = -1.26416e+04;
498 rndhittime = fT0Fill->GetRandom() - t0fill + currentRndTime;
504 // -----------------------------------------------------------------------
505 void AliAnalysisTaskCompareAODTrackCuts::SetUseNSigmaOnPIDAxes(Bool_t UseNSigma) {
507 // Will use NSigma on all PID axes. Will also change all track cuts objects
508 // owned by this task.
509 if (fDebug > 0) {cout << Form("File: %s, Line: %i, Function: %s",__FILE__,__LINE__,__func__) << endl;}
511 fUseNSigmaOnPIDAxes = UseNSigma;
514 for (Int_t iCut = 0; iCut < fTrackCuts->GetSize(); ++iCut) {
515 AliAODTrackCutsDiHadronPID* cutstmp = (AliAODTrackCutsDiHadronPID*)(fTrackCuts->At(iCut));
516 if (cutstmp) {cutstmp->SetUseNSigmaOnPIDAxes(UseNSigma);}
517 else {cout << Form("%s -> WARNING: Found an empty spot in the track cuts array...",__func__) << endl;}
522 // -----------------------------------------------------------------------
523 void AliAnalysisTaskCompareAODTrackCuts::SetEventCuts(AliAODEventCutsDiHadronPID* eventcuts) {
525 if (fDebug > 0) {cout << Form("File: %s, Line: %i, Function: %s",__FILE__,__LINE__,__func__) << endl;}
526 if (!eventcuts) {cout << Form("%s -> ERROR: No Event Cuts Object provided.",__func__) << endl; return;}
528 fEventCuts = eventcuts;
532 // -----------------------------------------------------------------------
533 void AliAnalysisTaskCompareAODTrackCuts::AddTrackCuts(AliAODTrackCutsDiHadronPID* trackcuts) {
535 if (fDebug > 0) {cout << Form("File: %s, Line: %i, Function: %s",__FILE__,__LINE__,__func__) << endl;}
536 if (!trackcuts) {cout << Form("%s -> ERROR: No Track Cuts Object provided.",__func__) << endl; return;}
537 if (!fTrackCuts) {cout << Form("%s -> ERROR: No Track Cuts array available.",__func__) << endl; return;}
539 // The setting of the task propagates to the imported track cuts object.
540 trackcuts->SetUseNSigmaOnPIDAxes(fUseNSigmaOnPIDAxes);
542 fTrackCuts->AddLast(trackcuts);
546 // -----------------------------------------------------------------------
547 void AliAnalysisTaskCompareAODTrackCuts::SetDebugLevel(Int_t debuglvl) {
549 // Sets debug level to a certain value, as well as the debug level of the
550 // track cuts objects and event cut object.
551 if (fDebug > 0) {cout << Form("File: %s, Line: %i, Function: %s",__FILE__,__LINE__,__func__) << endl;}
555 if (fEventCuts) {fEventCuts->SetDebugLevel(debuglvl);}
558 for (Int_t iTrackCutObj = 0; iTrackCutObj < fTrackCuts->GetEntriesFast(); ++iTrackCutObj) {
559 ((AliTrackDiHadronPID*)fTrackCuts->At(iTrackCutObj))->SetDebugLevel(debuglvl);
565 // -----------------------------------------------------------------------
566 void AliAnalysisTaskCompareAODTrackCuts::Terminate(Option_t*) {
568 if (fDebug > 0) {cout << Form("File: %s, Line: %i, Function: %s",__FILE__,__LINE__,__func__) << endl;}
570 if (fCalculateTOFMismatch) {
575 delete fLvsEtaProjections;
576 fLvsEtaProjections = 0x0;
581 // -----------------------------------------------------------------------