1 //------------------------------------------------------------------------------
2 // Implementation of AliPerformancePtCalibMC class. It compares ESD, TPC track
3 // momenta with MC information.
4 // The output can be analysed with AliPerfAnalyzeInvPt* via AliPerformancePtCalibMC::Analyse():
5 // Projection of 1/pt vs theta and vs phi resp. histoprams will be fitted with either
6 // polynomial or gaussian fit function to extract minimum position of 1/pt.
7 // Fit options and theta, phi bins can be set by user.
8 // Attention: use the Set* functions of AliPerformancePtCalibMC when running
9 // AliPerformancePtCalibMC::Analyse()
10 // The result of the analysis (histograms/graphs) are stored in the folder which is
11 // a data member of AliPerformancePtCalib*.
13 // Author: S.Schuchmann 11/13/2009
14 //------------------------------------------------------------------------------
18 // after running comparison task, read the file, and get component
19 gROOT->LoadMacro("$ALICE_ROOT/PWG1/Macros/LoadMyLibs.C");
22 TFile f("Output.root");
23 AliPerformancePtCalib * compObj = (AliPerformancePtCalibMC*)coutput->FindObject("AliPerformancePtCalibMC");
25 // analyse comparison data
28 // the output histograms/graphs will be stored in the folder "folderRes"
29 compObj->GetAnalysisFolder()->ls("*");
31 // user can save whole comparison object (or only folder with anlysed histograms)
32 // in the seperate output file (e.g.)
33 TFile fout("Analysed_InvPt.root","recreate");
34 compObj->Write(); // compObj->GetAnalysisFolder()->Write();
46 #include "AliESDEvent.h"
47 #include "AliESDtrack.h"
48 #include "AliESDtrackCuts.h"
49 #include "AliMCEvent.h"
51 #include "AliESDfriendTrack.h"
52 #include "AliESDfriend.h"
54 #include "AliPerformancePtCalibMC.h"
55 #include "AliPerfAnalyzeInvPt.h"
60 ClassImp(AliPerformancePtCalibMC)
62 //________________________________________________________________________
63 AliPerformancePtCalibMC::AliPerformancePtCalibMC() :
64 AliPerformanceObject("AliPerformancePtCalibMC"),
65 // option parameter for AliPerformancePtCalibMC::Analyse()
72 // option parameter for user defined charge/pt shift
85 fMaxChi2PerClusterTPC(0),
88 fAcceptKinkDaughters(0),
89 fRequireSigmaToVertex(0),
102 fHistPrimaryVertexPosX(0),
103 fHistPrimaryVertexPosY(0),
104 fHistPrimaryVertexPosZ(0),
105 fHistTrackMultiplicity(0),
106 fHistTrackMultiplicityCuts(0),
107 fHistTPCMomentaPosP(0),
108 fHistTPCMomentaNegP(0),
109 fHistTPCMomentaPosPt(0),
110 fHistTPCMomentaNegPt(0),
111 fHistInvPtThetaMC(0),
121 fHistTPCMomentaPosInvPtMC(0),
122 fHistTPCMomentaNegInvPtMC(0),
123 fHistTPCMomentaPosPtMC(0),
124 fHistTPCMomentaNegPtMC(0),
125 fHistESDMomentaPosInvPtMC(0),
126 fHistESDMomentaNegInvPtMC(0),
127 fHistESDMomentaPosPtMC(0),
128 fHistESDMomentaNegPtMC(0),
137 fShift = kFALSE; // shift in charge/pt yes/no
138 fDeltaInvP = 0.00; // shift value
140 fOptTPC = kTRUE; // read TPC tracks yes/no
141 fESDcuts = kTRUE; // read ESD track cuts
142 fRefitTPC = kFALSE; // require TPC refit
143 fRefitITS = kFALSE; // require ITS refit
144 fDCAcut = kTRUE; // apply DCA cuts
148 fEtaAcceptance = 0.8;
149 fMinPt=0.15; // GeV/c
150 fMaxPt=1.e10; // GeV/c
151 fMinNClustersTPC = 50;
152 fMaxChi2PerClusterTPC = 4.0;
153 fMaxDCAtoVertexXY = 2.4; // cm
154 fMaxDCAtoVertexZ = 3.0; // cm
155 fAcceptKinkDaughters = kFALSE;
156 fRequireSigmaToVertex = kFALSE;
157 fDCAToVertex2D = kTRUE;
159 // options for function AliPerformancePtCalibMC::Analyse()
160 fFitGaus = kFALSE;// use gaussian function for fitting charge/pt yes/no
161 fNThetaBins = 0; //number of theta bins
162 fNPhiBins = 0; //number of phi bins
163 fRange = 0; //fit range around 0
164 fExclRange =0; //range of rejection of points around 0
165 fAnaMC = kTRUE; // analyse MC tracks yes/no
170 //________________________________________________________________________
171 AliPerformancePtCalibMC::AliPerformancePtCalibMC(const char *name= "AliPerformancePtCalibMC", const char *title="AliPerformancePtCalibMC")://,Int_t analysisMode=0,Bool_t hptGenerator=kFALSE) :
172 AliPerformanceObject(name,title),
173 // option parameter for AliPerformancePtCalibMC::Analyse()
180 // option parameter for user defined 1/pt shift
193 fMaxChi2PerClusterTPC(0),
194 fMaxDCAtoVertexXY(0),
196 fAcceptKinkDaughters(0),
197 fRequireSigmaToVertex(0),
210 fHistPrimaryVertexPosX(0),
211 fHistPrimaryVertexPosY(0),
212 fHistPrimaryVertexPosZ(0),
213 fHistTrackMultiplicity(0),
214 fHistTrackMultiplicityCuts(0),
215 fHistTPCMomentaPosP(0),
216 fHistTPCMomentaNegP(0),
217 fHistTPCMomentaPosPt(0),
218 fHistTPCMomentaNegPt(0),
219 fHistInvPtThetaMC(0),
229 fHistTPCMomentaPosInvPtMC(0),
230 fHistTPCMomentaNegInvPtMC(0),
231 fHistTPCMomentaPosPtMC(0),
232 fHistTPCMomentaNegPtMC(0),
233 fHistESDMomentaPosInvPtMC(0),
234 fHistESDMomentaNegInvPtMC(0),
235 fHistESDMomentaPosPtMC(0),
236 fHistESDMomentaNegPtMC(0),
246 fShift = kFALSE; // shift in charge/pt yes/no
247 fDeltaInvP = 0.00; // shift value
249 fOptTPC = kTRUE; // read TPC tracks yes/no
250 fESDcuts = kTRUE; // read ESD track cuts
251 fRefitTPC = kFALSE; // require TPC refit
252 fRefitITS = kFALSE; // require ITS refit
253 fDCAcut = kTRUE; // apply DCA cuts
257 fEtaAcceptance = 0.8;
258 fMinPt=0.15; // GeV/c
259 fMaxPt=1.e10; // GeV/c
260 fMinNClustersTPC = 50;
261 fMaxChi2PerClusterTPC = 4.0;
262 fMaxDCAtoVertexXY = 2.4; // cm
263 fMaxDCAtoVertexZ = 3.0; // cm
264 fAcceptKinkDaughters = kFALSE;
265 fRequireSigmaToVertex = kFALSE;
266 fDCAToVertex2D = kTRUE;
268 // options for function AliPerformancePtCalibMC::Analyse()
269 fFitGaus = kFALSE;// use gaussian function for fitting charge/pt yes/no
270 fNThetaBins = 0; //number of theta bins
271 fNPhiBins = 0; //number of phi bins
272 fRange = 0; //fit range around 0
273 fExclRange =0; //range of rejection of points around 0
274 fAnaMC = kTRUE; // analyse MC tracks yes/no
279 //________________________________________________________________________
280 AliPerformancePtCalibMC::~AliPerformancePtCalibMC() {
284 if(fAnalysisFolder) delete fAnalysisFolder; fAnalysisFolder=0;
287 //________________________________________________________________________
288 void AliPerformancePtCalibMC::Init()
296 fAnalysisFolder = CreateFolder("folderPt_TPC","Analysis Pt Resolution Folder");
297 fList->Add(fAnalysisFolder);
299 fHistPrimaryVertexPosX = new TH1F("fHistPrimaryVertexPosX", "Primary Vertex Position X;Primary Vertex Position X (cm);Events",100,-0.5,0.5);
300 fList->Add(fHistPrimaryVertexPosX);
301 fHistPrimaryVertexPosY = new TH1F("fHistPrimaryVertexPosY", "Primary Vertex Position Y;Primary Vertex Position Y (cm);Events",100,-0.5,0.5);
302 fList->Add(fHistPrimaryVertexPosY);
303 fHistPrimaryVertexPosZ = new TH1F("fHistPrimaryVertexPosZ", "Primary Vertex Position Z;Primary Vertex Position Z (cm);Events",200,-2.0,2.0);
304 fList->Add(fHistPrimaryVertexPosZ);
307 fHistTrackMultiplicity = new TH1F("fHistTrackMultiplicity", "Multiplicity distribution;Number of tracks;Events", 250, 0, 250);
308 fList->Add(fHistTrackMultiplicity);
309 fHistTrackMultiplicityCuts = new TH1F("fHistTrackMultiplicityCuts", "Multiplicity distribution;Number of tracks after cuts;Events", 250, 0, 250);
310 fList->Add(fHistTrackMultiplicityCuts);
313 fHistPtShift0 = new TH1F("fHistPtShift0","1/pt dN/pt vs. pt of ESD track ",600,0.0,6.0);
314 fList->Add(fHistPtShift0);
315 fHistInvPtTheta = new TH2F("fHistInvPtTheta","#theta vs 1/pt ",900, -4.5, 4.5,300,0.0,3.0);
316 fList->Add(fHistInvPtTheta);
317 fHistInvPtPhi = new TH2F("fHistInvPtPhi","#phi vs 1/pt",900, -4.5, 4.5,325,0.0,6.5);
318 fList->Add(fHistInvPtPhi);
319 fHistPtTheta = new TH2F("fHistPtTheta"," #theta vs pt ",300, 0.0, 15.0,300,0.0,3.0);
320 fList->Add(fHistPtTheta);
321 fHistPtPhi = new TH2F("fHistPtPhi"," #phi vs pt ",300, 0.0,15.0,325,0.0,6.5);
322 fList->Add(fHistPtPhi);
325 fHistTPCMomentaPosP = new TH2F("fHistTPCMomentaPosP","TPC p vs global esd track p pos",300,0.0,15.0,300,0.0,15.0);
326 fList->Add(fHistTPCMomentaPosP);
327 fHistTPCMomentaNegP = new TH2F("fHistTPCMomentaNegP","TPC p vs global esd track p neg",300,0.0,15.0,300,0.0,15.0);
328 fList->Add(fHistTPCMomentaNegP);
329 fHistTPCMomentaPosPt = new TH2F("fHistTPCMomentaPosPt","TPC pt vs global esd track pt pos",300,0.0,15.0,300,0.0,15.0);
330 fList->Add(fHistTPCMomentaPosPt);
331 fHistTPCMomentaNegPt = new TH2F("fHistTPCMomentaNegPt","TPC pt vs global esd track pt neg",300,0.0,15.0,300,0.0,15.0);
332 fList->Add(fHistTPCMomentaNegPt);
334 // mom test histos MC
335 fHistTPCMomentaPosInvPtMC = new TH2F("fHistTPCMomentaPosInvPtMC","TPC-MC of 1/pt vs global ESD-MC of 1/pt pos",500, -10.0, 10.0,500, -10.0,10.0);
336 fList->Add(fHistTPCMomentaPosInvPtMC);
337 fHistTPCMomentaNegInvPtMC = new TH2F("fHistTPCMomentaNegInvPtMC","TPC-MC of 1/pt vs global ESD-MC 1/pt neg",500, -10.0, 10.0,500, -10.0, 10.0);
338 fList->Add(fHistTPCMomentaNegInvPtMC);
339 fHistTPCMomentaPosPtMC = new TH2F("fHistTPCMomentaPosPtMC","TPC-MC of pt vs global ESD-MC of pt pos",600,-4.0,44.0,600,-4.0,44.0);
340 fList->Add(fHistTPCMomentaPosPtMC);
341 fHistTPCMomentaNegPtMC = new TH2F("fHistTPCMomentaNegPtMC","TPC-MC of pt vs global ESD-MC of pt neg",600,-4.0,44.0,600,-4.0,44.0);
342 fList->Add(fHistTPCMomentaNegPtMC);
343 fHistESDMomentaPosInvPtMC = new TH1F("fHistESDMomentaPosInvPtMC","ESD-MC of 1/pt ",500, -10.0, 10.0);
344 fList->Add(fHistESDMomentaPosInvPtMC);
345 fHistESDMomentaNegInvPtMC = new TH1F("fHistESDMomentaNegInvPtMC","ESD-MC of 1/pt",500, -10.0, 10.0);
346 fList->Add(fHistESDMomentaNegInvPtMC);
347 fHistESDMomentaPosPtMC = new TH1F("fHistESDMomentaPosPtMC","ESD-MC of pt ",600,-4.0,44.0);
348 fList->Add(fHistESDMomentaPosPtMC);
349 fHistESDMomentaNegPtMC = new TH1F("fHistESDMomentaNegPtMC","ESD-MC of pt ",600,-4.0,44.0);
350 fList->Add(fHistESDMomentaNegPtMC);
353 fHistInvPtThetaMC = new TH2F("fHistInvPtThetaMC","theta vs inv pt MC",900, -4.5, 4.5,300,0.0,3.0);
354 fList->Add(fHistInvPtThetaMC);
355 fHistInvPtPhiMC = new TH2F("fHistInvPtPhiMC","phi vs inv pt MC",900, -4.5, 4.5,325,0.0,6.5);
356 fList->Add(fHistInvPtPhiMC);
357 fHistPtThetaMC = new TH2F("fHistPtThetaMC","theta vs pt MC",300, 0.0, 15.0,300,0.0,3.0);
358 fList->Add(fHistPtThetaMC);
359 fHistPtPhiMC = new TH2F("fHistPtPhiMC"," phi vs pt MC",300, 0.0,15.0,325,0.0,6.5);
360 fList->Add(fHistPtPhiMC);
362 //correlation histos MC ESD or TPC
363 fHistInvPtMCESD = new TH2F("fHistInvPtMCESD","inv pt ESD vs MC",900, 0.0, 9.0,900, 0.0, 9.0);
364 fList->Add(fHistInvPtMCESD);
365 fHistPtMCESD = new TH2F("fHistPtMCESD"," pt ESD vs MC",300, 0.0, 15.0,300, 0.0, 15.0);
366 fList->Add(fHistPtMCESD);
367 fHistInvPtMCTPC = new TH2F("fHistInvPtMCTPC","inv pt TPC vs MC",900, 0.0, 9.0,900, 0.0, 9.0);
368 fList->Add(fHistInvPtMCTPC);
369 fHistPtMCTPC = new TH2F("fHistPtMCTPC"," pt TPC vs MC",300, 0.0, 15.0,300, 0.0, 15.0);
370 fList->Add(fHistPtMCTPC);
371 fHistMomresMCESD = new TH2F("fHistMomresMCESD"," (pt ESD - pt MC)/ptMC vs pt MC",300, 0.0, 15.0,400, -2.0, 2.0);
372 fList->Add(fHistMomresMCESD);
373 fHistMomresMCTPC = new TH2F("fHistMomresMCTPC"," (pt TPC - pt MC)/ptMC vs pt MC",300, 0.0, 15.0,400, -2.0, 2.0);
374 fList->Add(fHistMomresMCTPC);
377 //user pt shift check
378 fHistUserPtShift = new TH1F("fHistUserPtShift","user defined shift in 1/pt",100,-0.5,1.5);
379 fList->Add(fHistUserPtShift);
382 fESDTrackCuts = new AliESDtrackCuts("AliESDtrackCuts");
384 // //fESDTrackCuts->DefineHistoqgrams(1);
386 fESDTrackCuts->SetRequireSigmaToVertex(fRequireSigmaToVertex);
387 fESDTrackCuts->SetRequireTPCRefit(fRefitTPC);
388 fESDTrackCuts->SetAcceptKinkDaughters(fAcceptKinkDaughters);
389 fESDTrackCuts->SetMinNClustersTPC((Int_t)fMinNClustersTPC);
390 fESDTrackCuts->SetMaxChi2PerClusterTPC(fMaxChi2PerClusterTPC);
391 fESDTrackCuts->SetMaxDCAToVertexXY(fMaxDCAtoVertexXY);
392 fESDTrackCuts->SetMaxDCAToVertexZ(fMaxDCAtoVertexZ);
393 fESDTrackCuts->SetDCAToVertex2D(fDCAToVertex2D);
394 fESDTrackCuts->SetPtRange(fMinPt,fMaxPt);
399 //________________________________________________________________________
400 void AliPerformancePtCalibMC::SetPtShift(const Double_t shiftVal ) {
402 //set user defined shift in charge/pt
404 if(shiftVal) { fShift=kTRUE; fDeltaInvP = shiftVal; }
407 //________________________________________________________________________
408 void AliPerformancePtCalibMC::Exec(AliMCEvent* const mcEvent, AliESDEvent* const esdEvent , AliESDfriend*, Bool_t, Bool_t)
410 //exec: read MC and esd or tpc tracks
415 Printf("ERROR: Event not available");
419 if (!(esdEvent->GetNumberOfTracks())) {
420 Printf(" PtCalibMC task: There is no track in this event");
423 fHistTrackMultiplicity->Fill(esdEvent->GetNumberOfTracks());
426 Printf("ERROR: Could not retrieve MC event");
429 stack = mcEvent->Stack();
431 Printf("ERROR: Could not retrieve stack");
435 if(fShift) fHistUserPtShift->Fill(fDeltaInvP);
437 // read primary vertex info
438 Double_t tPrimaryVtxPosition[3];
439 // Double_t tPrimaryVtxCov[3];
440 const AliESDVertex *primaryVtx = esdEvent->GetPrimaryVertexTPC();
442 tPrimaryVtxPosition[0] = primaryVtx->GetXv();
443 tPrimaryVtxPosition[1] = primaryVtx->GetYv();
444 tPrimaryVtxPosition[2] = primaryVtx->GetZv();
446 fHistPrimaryVertexPosX->Fill(tPrimaryVtxPosition[0]);
447 fHistPrimaryVertexPosY->Fill(tPrimaryVtxPosition[1]);
448 fHistPrimaryVertexPosZ->Fill(tPrimaryVtxPosition[2]);
451 //fill histos for pt spectra and shift of transverse momentum
454 for(Int_t j = 0;j<esdEvent->GetNumberOfTracks();j++){
455 AliESDtrack *esdTrack = esdEvent->GetTrack(j);
456 if(!esdTrack) continue;
459 if(fESDcuts){Printf("esd cuts aplied");
460 if(!fESDTrackCuts->AcceptTrack(esdTrack)) continue;
464 if(fRefitTPC) if(AddTPCcuts(esdTrack)) continue;
465 if(fRefitITS) if(AddITScuts(esdTrack)) continue;
466 if(fDCAcut) if(AddDCAcuts(esdTrack)) continue ;
469 Int_t esdLabel = esdTrack->GetLabel();
470 if(esdLabel<0) continue;
471 TParticle * partMC = stack->Particle(esdLabel);
472 if (!partMC) continue;
474 // fill correlation histos MC ESD
475 Double_t pESD = esdTrack->GetP();
476 Double_t ptESD = esdTrack->GetSignedPt();
478 if(!ptESD || !(partMC->Pt()) ) continue;
479 Double_t mcPt = partMC->Pt();
480 Double_t invPtMC = 1.0/mcPt;
481 Int_t signMC = partMC->GetPdgCode();
483 if(signMC>0) signMC = 1;
487 fHistInvPtThetaMC->Fill(signMC*fabs(invPtMC),partMC->Theta());
488 fHistInvPtPhiMC->Fill(signMC*fabs(invPtMC),partMC->Phi());
489 fHistPtThetaMC->Fill(fabs(mcPt),partMC->Theta(),fabs(invPtMC));
490 fHistPtPhiMC->Fill(fabs(mcPt),partMC->Phi(),fabs(invPtMC));
492 //correlation histos MC ESD
493 fHistInvPtMCESD->Fill(fabs(invPtMC),fabs(1.0/ptESD));
494 fHistPtMCESD->Fill(fabs(mcPt),fabs(ptESD));
499 //TPC tracks and MC tracks
500 const AliExternalTrackParam *tpcTrack = esdTrack->GetTPCInnerParam();
501 if(!tpcTrack) continue;
502 if(fabs(tpcTrack->Eta())> fEtaAcceptance) continue;
504 Double_t signedPt = tpcTrack->GetSignedPt();
505 Double_t invPt = 0.0;
507 invPt = 1.0/signedPt;
509 fHistPtShift0->Fill(fabs(signedPt));
512 invPt += fDeltaInvP; //shift momentum for tests
513 if(invPt) signedPt = 1.0/invPt;
517 fHistInvPtTheta->Fill(invPt,tpcTrack->Theta());
518 fHistInvPtPhi->Fill(invPt,tpcTrack->Phi());
519 fHistPtTheta->Fill(fabs(signedPt),tpcTrack->Theta());
520 fHistPtPhi->Fill(fabs(signedPt),tpcTrack->Phi());
523 //correlation histos MC TPC
524 fHistInvPtMCTPC->Fill(fabs(invPtMC),fabs(invPt));
525 fHistPtMCTPC->Fill(fabs(mcPt),fabs(signedPt));
528 Double_t ptDiffESD = (fabs(ptESD)-fabs(mcPt))/pow(mcPt,2);
529 Double_t ptDiffTPC = (fabs(signedPt)-fabs(mcPt))/pow(mcPt,2);
530 Double_t invPtDiffESD = fabs(1.0/ptESD)-1.0/fabs(mcPt);
531 Double_t invPtDiffTPC = fabs(invPt)-1.0/fabs(mcPt);
532 Double_t pTPC = tpcTrack->GetP();
534 if(esdTrack->GetSign()>0){//compare momenta ESD track and TPC track
535 fHistTPCMomentaPosP->Fill(fabs(pESD),fabs(pTPC));
536 fHistTPCMomentaPosPt->Fill(fabs(ptESD),fabs(signedPt));
537 fHistTPCMomentaPosInvPtMC->Fill(invPtDiffESD,invPtDiffTPC);
538 fHistTPCMomentaPosPtMC->Fill(ptDiffESD,ptDiffTPC);
541 fHistTPCMomentaNegP->Fill(fabs(pESD),fabs(pTPC));
542 fHistTPCMomentaNegPt->Fill(fabs(ptESD),fabs(signedPt));
543 fHistTPCMomentaNegInvPtMC->Fill(invPtDiffESD,invPtDiffTPC);
544 fHistTPCMomentaNegPtMC->Fill(ptDiffESD,ptDiffTPC);
546 fHistMomresMCESD->Fill((fabs(mcPt)-fabs(ptESD))/fabs(mcPt),fabs(mcPt));
547 fHistMomresMCTPC->Fill((fabs(mcPt)-fabs(signedPt))/fabs(mcPt),fabs(mcPt));
554 // ESD tracks and MC tracks
555 Double_t invPt = 0.0;
559 fHistPtShift0->Fill(fabs(ptESD));
562 invPt += fDeltaInvP; //shift momentum for tests
563 if(invPt) ptESD = 1.0/invPt;
566 fHistInvPtTheta->Fill(invPt,esdTrack->Theta());
567 fHistInvPtPhi->Fill(invPt,esdTrack->Phi());
568 fHistPtTheta->Fill(ptESD,esdTrack->Theta());
569 fHistPtPhi->Fill(ptESD,esdTrack->Phi());
571 //differences MC ESD tracks
572 Double_t ptDiffESD = (fabs(ptESD)-fabs(mcPt))/pow(mcPt,2);
573 Double_t invPtdiffESD = fabs(1.0/ptESD)-1.0/fabs(mcPt);
574 if(esdTrack->GetSign()>0){
575 fHistESDMomentaPosInvPtMC->Fill(invPtdiffESD);
576 fHistESDMomentaPosPtMC->Fill(ptDiffESD);
579 fHistESDMomentaNegInvPtMC->Fill(invPtdiffESD);
580 fHistESDMomentaNegPtMC->Fill(ptDiffESD);
583 fHistMomresMCESD->Fill((fabs(mcPt)-fabs(ptESD))/fabs(mcPt),fabs(mcPt));
589 fHistTrackMultiplicityCuts->Fill(count);
593 //______________________________________________________________________________________________________________________
594 Bool_t AliPerformancePtCalibMC::AddTPCcuts(const AliESDtrack *esdTrack){
599 if (!(esdTrack->GetStatus()&AliESDtrack::kTPCrefit)) cut=kTRUE; // TPC refit
600 if (esdTrack->GetTPCNcls()<fMinNClustersTPC) cut=kTRUE; // min. nb. TPC clusters
601 if(cut) return kTRUE;
604 //______________________________________________________________________________________________________________________
605 Bool_t AliPerformancePtCalibMC::AddDCAcuts(const AliESDtrack *esdTrack){
609 Float_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z and impact parameters:
610 esdTrack->GetImpactParameters(dca,cov);
611 if(TMath::Abs(dca[0])>fMaxDCAtoVertexXY || TMath::Abs(dca[1])>fMaxDCAtoVertexZ) cut=kTRUE;
612 if(esdTrack->GetKinkIndex(0)>0) cut=kTRUE;
613 if(cut) return kTRUE;
617 //______________________________________________________________________________________________________________________
618 Bool_t AliPerformancePtCalibMC::AddITScuts(const AliESDtrack *esdTrack){
622 if (!(esdTrack->GetStatus()&AliESDtrack::kITSrefit)) cut=kTRUE; // ITS refit
623 Int_t clusterITS[200];
624 if(esdTrack->GetITSclusters(clusterITS)<2) cut=kTRUE; // min. nb. ITS clusters //3
626 if(cut) return kTRUE;
630 //______________________________________________________________________________________________________________________
632 void AliPerformancePtCalibMC::Analyse()
635 // analyse charge/pt spectra in bins of theta and phi. Bins can be set by user
637 AliPerfAnalyzeInvPt *ana = new AliPerfAnalyzeInvPt("AliPerfAnalyzeInvPt","AliPerfAnalyzeInvPt");
639 TH1::AddDirectory(kFALSE);
641 ana->SetProjBinsTheta(fThetaBins,fNThetaBins);
642 ana->SetProjBinsPhi(fPhiBins,fNPhiBins);
643 ana->SetMakeFitOption(fFitGaus,fExclRange,fRange);
645 TObjArray *aFolderObj = new TObjArray;
647 Printf("AliPerformancePtCalibMC: analysing MC!");
648 ana->StartAnalysis(fHistInvPtThetaMC,fHistInvPtPhiMC, aFolderObj);
651 Printf("AliPerformancePtCalibMC: analysing data!");
652 ana->StartAnalysis(fHistInvPtTheta,fHistInvPtPhi, aFolderObj);
655 // export objects to analysis folder
656 fAnalysisFolder = ExportToFolder(aFolderObj);
658 // delete only TObjArray
659 if(aFolderObj) delete aFolderObj;
664 //______________________________________________________________________________________________________________________
665 TFolder* AliPerformancePtCalibMC::ExportToFolder(TObjArray * array)
667 // recreate folder avery time and export objects to new one
669 AliPerformancePtCalibMC * comp=this;
670 TFolder *folder = comp->GetAnalysisFolder();
673 TFolder *newFolder = 0;
675 Int_t size = array->GetSize();
678 // get name and title from old folder
679 name = folder->GetName();
680 title = folder->GetTitle();
686 newFolder = CreateFolder(name.Data(),title.Data());
687 newFolder->SetOwner();
689 // add objects to folder
691 newFolder->Add(array->At(i));
699 //______________________________________________________________________________________________________________________
700 Long64_t AliPerformancePtCalibMC::Merge(TCollection* const list)
702 // Merge list of objects (needed by PROOF)
710 TIterator* iter = list->MakeIterator();
713 // collection of generated histograms
715 while((obj = iter->Next()) != 0)
717 AliPerformancePtCalibMC* entry = dynamic_cast<AliPerformancePtCalibMC*>(obj);
718 if (!entry) continue;
720 fHistInvPtTheta->Add(entry->fHistInvPtTheta);
721 fHistInvPtPhi->Add(entry-> fHistInvPtPhi);
722 fHistPtTheta->Add(entry->fHistPtTheta);
723 fHistPtPhi->Add(entry->fHistPtPhi);
725 fHistInvPtThetaMC->Add(entry->fHistInvPtThetaMC);
726 fHistInvPtPhiMC->Add(entry->fHistInvPtPhiMC);
727 fHistPtThetaMC->Add(entry->fHistPtThetaMC);
728 fHistPtPhiMC->Add(entry->fHistPtPhiMC);
729 fHistInvPtMCESD->Add(entry->fHistInvPtMCESD);
730 fHistPtMCESD->Add(entry->fHistPtMCESD);
731 fHistInvPtMCTPC->Add(entry->fHistInvPtMCTPC);
732 fHistPtMCTPC->Add(entry->fHistPtMCTPC);
733 fHistMomresMCESD->Add(entry->fHistMomresMCESD);
734 fHistMomresMCTPC->Add(entry->fHistMomresMCTPC);
736 fHistPtShift0->Add(entry->fHistPtShift0);
737 fHistPrimaryVertexPosX->Add(entry->fHistPrimaryVertexPosX);
738 fHistPrimaryVertexPosY->Add(entry->fHistPrimaryVertexPosY);
739 fHistPrimaryVertexPosZ->Add(entry->fHistPrimaryVertexPosZ);
740 fHistTrackMultiplicity->Add(entry->fHistTrackMultiplicity);
741 fHistTrackMultiplicityCuts->Add(entry->fHistTrackMultiplicityCuts);
743 fHistTPCMomentaPosP->Add(entry->fHistTPCMomentaPosP);
744 fHistTPCMomentaNegP->Add(entry->fHistTPCMomentaNegP);
745 fHistTPCMomentaPosPt->Add(entry->fHistTPCMomentaPosPt);
746 fHistTPCMomentaNegPt->Add(entry->fHistTPCMomentaNegPt);
747 fHistTPCMomentaPosInvPtMC->Add(entry->fHistTPCMomentaPosInvPtMC);
748 fHistTPCMomentaNegInvPtMC->Add(entry->fHistTPCMomentaNegInvPtMC);
749 fHistTPCMomentaPosPtMC->Add(entry->fHistTPCMomentaPosPtMC);
750 fHistTPCMomentaNegPtMC->Add(entry->fHistTPCMomentaNegPtMC);
751 fHistESDMomentaPosInvPtMC->Add(entry->fHistESDMomentaPosInvPtMC);
752 fHistESDMomentaNegInvPtMC->Add(entry->fHistESDMomentaNegInvPtMC);
753 fHistESDMomentaPosPtMC->Add(entry->fHistESDMomentaPosPtMC);
754 fHistESDMomentaNegPtMC->Add(entry->fHistESDMomentaNegPtMC);
761 //______________________________________________________________________________________________________________________
762 TFolder* AliPerformancePtCalibMC::CreateFolder(TString name,TString title) {
763 // create folder for analysed histograms
766 folder = new TFolder(name.Data(),title.Data());
772 // set variables for Analyse()
774 //______________________________________________________________________________________________________________________
775 void AliPerformancePtCalibMC::SetProjBinsPhi(const Double_t *phiBinArray,const Int_t nphBins){
777 // set phi bins for Analyse()
778 //set phi bins as array and set number of this array which is equal to number of bins analysed
779 //the last analysed bin will always be the projection from first to last bin in the array
783 for(Int_t k = 0;k<fNPhiBins;k++){
784 fPhiBins[k] = phiBinArray[k];
786 Printf("AliPerformancePtCalibMC: number of bins in phi set to %i",fNPhiBins);
788 else Printf("Warning AliPerformancePtCalibMC::SetProjBinsPhi: number of bins in phi NOT set!!! Default values are taken.");
790 //____________________________________________________________________________________________________________________________________________
791 void AliPerformancePtCalibMC::SetProjBinsTheta(const Double_t *thetaBinArray, const Int_t nthBins){
792 // set theta bins for Analyse()
793 //set theta bins as array and set number of this array which is equal to number of bins analysed
794 //the last analysed bin will always be the projection from first to last bin in the array
796 fNThetaBins = nthBins;
797 for(Int_t k = 0;k<fNThetaBins;k++){
798 fThetaBins[k] = thetaBinArray[k];
800 Printf("AliPerformancePtCalibMC: number of bins in theta set to %i",fNThetaBins);
802 else Printf("Warning AliPerformancePtCalibMC::SetProjBinsTheta: number of bins in theta NOT set!!! Default values are taken.");
804 //____________________________________________________________________________________________________________________________________________
805 void AliPerformancePtCalibMC::SetMakeFitOption(const Bool_t setGausFit, const Double_t exclusionR,const Double_t fitR ){
807 //set the fit options:
808 //for usage of gaussian function instead of polynomial (default) set setGausFit=kTRUE
809 //set the range of rejection of points around 0 via exclusionR
810 //set the fit range around 0 with fitR
812 fFitGaus = setGausFit;
813 fExclRange = exclusionR;
816 if(fFitGaus) Printf("AliPerformancePtCalibMC:set MakeGausFit with fit range %2.3f and exclusion range in 1/pt: %2.3f",fRange,fExclRange);
817 else Printf("AliPerformancePtCalibMC: set standard polynomial fit with fit range %2.3f and exclusion range in 1/pt: %2.3f",fRange,fExclRange);
820 //____________________________________________________________________________________________________________________________________________
821 void AliPerformancePtCalibMC::SetESDcutValues(const Double_t * esdCutValues){
822 // set ESD cut values as an array of size 6
824 fMinPt = esdCutValues[0];
825 fMaxPt = esdCutValues[1];
826 fMinNClustersTPC = esdCutValues[2];
827 fMaxChi2PerClusterTPC = esdCutValues[3];
828 fMaxDCAtoVertexXY = esdCutValues[4];
829 fMaxDCAtoVertexZ = esdCutValues[5];