From 4060dd2d9777f53696a7db8a25c870905dbd839b Mon Sep 17 00:00:00 2001 From: jotwinow Date: Thu, 24 Jun 2010 09:41:22 +0000 Subject: [PATCH] Changes by Simone Schuchmann --- PWG1/TPC/AliPerformancePtCalib.cxx | 405 +++++++++--------- PWG1/TPC/AliPerformancePtCalib.h | 91 ++-- PWG1/TPC/AliPerformancePtCalibMC.cxx | 617 ++++++++++++++------------- PWG1/TPC/AliPerformancePtCalibMC.h | 109 ++--- 4 files changed, 614 insertions(+), 608 deletions(-) diff --git a/PWG1/TPC/AliPerformancePtCalib.cxx b/PWG1/TPC/AliPerformancePtCalib.cxx index f3bb6b2ea0c..cf23522d142 100755 --- a/PWG1/TPC/AliPerformancePtCalib.cxx +++ b/PWG1/TPC/AliPerformancePtCalib.cxx @@ -1,14 +1,31 @@ + +/************************************************************************** + * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * + * * + * Author: The ALICE Off-line Project. * + * Contributors are mentioned in the code where appropriate. * + * * + * Permission to use, copy, modify and distribute this software and its * + * documentation strictly for non-commercial purposes is hereby granted * + * without fee, provided that the above copyright notice appears in all * + * copies and that both the copyright notice and this permission notice * + * appear in the supporting documentation. The authors make no claims * + * about the suitability of this software for any purpose. It is * + * provided "as is" without express or implied warranty. * + **************************************************************************/ + //------------------------------------------------------------------------------ // Implementation of AliPerformancePtCalib class. It fills histograms with ESD or -// TPC track information to study 1.pt spectra. +// TPC track information to study charge/pt spectra. // To analyse the output of AliPerformancePtCalib use AliPerfAnalyzeInvPt class: -// Projection of 1/pt vs theta and vs phi resp. histoprams will be fitted with either -// polynomial or gaussian fit function to extract minimum position of 1/pt. +// Projection of charge/pt vs theta and vs phi resp. histoprams will be fitted +// with either polynomial or gaussian fit function to extract minimum position of +// charge/pt. // Fit options and theta, phi bins can be set by user. -// Attention: use the Set* functions of AliPerformancePtCalib.h when running +// Attention: use the Set functions of AliPerformancePtCalib when running // AliPerformancePtCalib::Analyse() // The result of the analysis (histograms/graphs) are stored in the folder which is -// a data member of AliPerformancePtCalib*. +// a data member of AliPerformancePtCalib. // // Author: S.Schuchmann 11/13/2009 // sschuchm@ikf.uni-frankfurt.de @@ -16,15 +33,18 @@ /* -// after running comparison task, read the file, and get component -gROOT->LoadMacro("$ALICE_ROOT/PWG1/Macros/LoadMyLibs.C"); -LoadMyLibs(); +// after running the performance task, read the file, and get component TFile f("Output.root"); AliPerformancePtCalib * compObj = (AliPerformancePtCalib*)coutput->FindObject("AliPerformancePtCalib"); -// analyse comparison data +// set phi and theta bins for fitting and analyse comparison data +compObj->SetProjBinsTheta(thetaBins,nThetaBins,minPhi, maxPhi); +compObj->SetProjBinsPhi(phiBins,nPhiBins,minTheta,maxTheta); +compObj->SetMakeFitOption(kFALSE,exclRange,fitRange); +compObj->SetDoRebin(rebin); compObj->Analyse(); +//details see functions of this class // the output histograms/graphs will be stored in the folder "folderRes" compObj->GetAnalysisFolder()->ls("*"); @@ -38,15 +58,21 @@ fout.Close(); */ #include "TH1F.h" #include "TH2F.h" +#include "THnSparse.h" #include "TList.h" #include "TMath.h" #include "TFolder.h" + #include "AliESDEvent.h" #include "AliESDtrack.h" #include "AliESDfriendTrack.h" #include "AliESDfriend.h" #include "AliESDtrackCuts.h" + +#include "AliPhysicsSelection.h" +#include "AliTriggerAnalysis.h" + #include "AliPerfAnalyzeInvPt.h" #include "AliPerformancePtCalib.h" @@ -61,6 +87,10 @@ ClassImp(AliPerformancePtCalib) // option parameter for AliPerformancePtCalib::Analyse() fNThetaBins(0), fNPhiBins(0), + fMaxPhi(0), + fMinPhi(0), + fMaxTheta(0), + fMinTheta(0), fRange(0), fExclRange(0), fFitGaus(0) , @@ -73,30 +103,16 @@ ClassImp(AliPerformancePtCalib) //options for cuts fOptTPC(0), fESDcuts(0), - fRefitTPC(0), - fRefitITS(0), - fDCAcut(0), + fPhysSel(0), fEtaAcceptance(0), - fMinPt(0), - fMaxPt(0), - fMinNClustersTPC(0), - fMaxChi2PerClusterTPC(0), - fMaxDCAtoVertexXY(0), - fMaxDCAtoVertexZ(0), - fAcceptKinkDaughters(0), - fRequireSigmaToVertex(0), - fDCAToVertex2D(0), - + fTrigger(AliTriggerAnalysis::kMB1), + fPhysicsSelection(0), fCutsRC(0), fCutsMC(0), - fList(0), // histograms - fHistInvPtTheta(0), - fHistInvPtPhi(0), - fHistPtTheta(0), - fHistPtPhi(0), - + fHistInvPtPtThetaPhi(0), + fHistPtShift0(0), fHistPrimaryVertexPosX(0), fHistPrimaryVertexPosY(0), @@ -112,7 +128,6 @@ ClassImp(AliPerformancePtCalib) //esd track cuts fESDTrackCuts(0), - // analysis folder fAnalysisFolder(0) { @@ -123,79 +138,61 @@ ClassImp(AliPerformancePtCalib) fDeltaInvP = 0.00; // shift value //options for cuts fOptTPC = kTRUE; // read TPC tracks yes/no - fESDcuts = kTRUE; // read ESD track cuts - fRefitTPC = kFALSE; // require TPC refit - fRefitITS = kFALSE; // require ITS refit - fDCAcut = kTRUE; // apply DCA cuts - + fESDcuts = kFALSE; + fPhysSel = kFALSE; fCutsRC = NULL; fCutsMC = NULL; //set options for esd track cuts fEtaAcceptance = 0.8; - fMinPt=0.15; // GeV/c - fMaxPt=1.e10; // GeV/c - fMinNClustersTPC = 50; - fMaxChi2PerClusterTPC = 4.0; - fMaxDCAtoVertexXY = 2.4; // cm - fMaxDCAtoVertexZ = 3.2; // cm - fAcceptKinkDaughters = kFALSE; - fRequireSigmaToVertex = kFALSE; - fDCAToVertex2D = kTRUE; - - // options for function AliPerformancePtCalibMC::Analyse() + + // options for function AliPerformancePtCalib::Analyse() fFitGaus = kFALSE;// use gaussian function for fitting charge/pt yes/no fNThetaBins = 0; //number of theta bins fNPhiBins = 0; //number of phi bins fRange = 0; //fit range around 0 + fMaxPhi = 6.5;// max phi for 2D projection on theta and charge/pt axis + fMinPhi = 0.0;// min phi for 2D projection on theta and charge/pt axis + fMaxTheta = 3.0;// max theta for 2D projection on phi and charge/pt axis + fMinTheta = 0.0;// min theta for 2D projection on phi and charge/pt axis fExclRange =0; //range of rejection of points around 0 fDoRebin = kFALSE; fRebin = 0; - + Init(); } //________________________________________________________________________ -AliPerformancePtCalib::AliPerformancePtCalib(Char_t * name="AliPerformancePtCalib",Char_t* title ="AliPerformancePtCalib")://,Int_t analysisMode=0,Bool_t hptGenerator=kFALSE) : +AliPerformancePtCalib::AliPerformancePtCalib(Char_t * name="AliPerformancePtCalib",Char_t* title ="AliPerformancePtCalib"): AliPerformanceObject(name,title), // option parameter for AliPerformancePtCalib::Analyse() fNThetaBins(0), fNPhiBins(0), + fMaxPhi(0), + fMinPhi(0), + fMaxTheta(0), + fMinTheta(0), fRange(0), fExclRange(0), fFitGaus(0) , fDoRebin(0), fRebin(0), + // option parameter for user defined charge/pt shift fShift(0), fDeltaInvP(0), - //options for cuts + //options for cuts fOptTPC(0), fESDcuts(0), - fRefitTPC(0), - fRefitITS(0), - fDCAcut(0), + fPhysSel(0), fEtaAcceptance(0), - fMinPt(0), - fMaxPt(0), - fMinNClustersTPC(0), - fMaxChi2PerClusterTPC(0), - fMaxDCAtoVertexXY(0), - fMaxDCAtoVertexZ(0), - fAcceptKinkDaughters(0), - fRequireSigmaToVertex(0), - fDCAToVertex2D(0), - + fTrigger(AliTriggerAnalysis::kMB1), + fPhysicsSelection(0), fCutsRC(0), fCutsMC(0), fList(0), - // histograms - fHistInvPtTheta(0), - fHistInvPtPhi(0), - fHistPtTheta(0), - fHistPtPhi(0), - + fHistInvPtPtThetaPhi(0), fHistPtShift0(0), fHistPrimaryVertexPosX(0), fHistPrimaryVertexPosY(0), @@ -210,10 +207,9 @@ AliPerformancePtCalib::AliPerformancePtCalib(Char_t * name="AliPerformancePtCali fHistUserPtShift(0), //esd track cuts fESDTrackCuts(0), - // analysis folder fAnalysisFolder(0) - + { // Constructor @@ -221,30 +217,22 @@ AliPerformancePtCalib::AliPerformancePtCalib(Char_t * name="AliPerformancePtCali fDeltaInvP = 0.00; //options for cuts fOptTPC = kTRUE; // read TPC tracks yes/no - fESDcuts = kTRUE; // read ESD track cuts - fRefitTPC = kFALSE; // require TPC refit - fRefitITS = kFALSE; // require ITS refit - fDCAcut = kTRUE; // apply DCA cuts - + fESDcuts = kFALSE; + fPhysSel = kFALSE; fCutsRC = NULL; fCutsMC = NULL; //esd track cut options fEtaAcceptance = 0.8; - fMinPt=0.15; // GeV/c - fMaxPt=1.e10; // GeV/c - fMinNClustersTPC = 50; - fMaxChi2PerClusterTPC = 4.0; - fMaxDCAtoVertexXY = 2.4; // cm - fMaxDCAtoVertexZ = 3.2; // cm - fAcceptKinkDaughters = kFALSE; - fRequireSigmaToVertex = kFALSE; - fDCAToVertex2D = kTRUE; // options for function AliPerformancePtCalibMC::Analyse() fFitGaus = kFALSE;// use gaussian function for fitting charge/pt yes/no fNThetaBins = 0; //number of theta bins fNPhiBins = 0; //number of phi bins + fMaxPhi = 6.5;// max phi for 2D projection on theta and charge/pt axis + fMinPhi = 0.0;// min phi for 2D projection on theta and charge/pt axis + fMaxTheta = 3.0;// max theta for 2D projection on phi and charge/pt axis + fMinTheta = 0.0;// min theta for 2D projection on phi and charge/pt axis fRange = 0; //fit range around 0 fExclRange =0; //range of rejection of points around 0 fDoRebin = kFALSE; @@ -252,11 +240,30 @@ AliPerformancePtCalib::AliPerformancePtCalib(Char_t * name="AliPerformancePtCali Init(); } - +//________________________________________________________________________ AliPerformancePtCalib::~AliPerformancePtCalib() { // // destructor // + + if(fList) delete fList; + // histograms + if( fHistInvPtPtThetaPhi) delete fHistInvPtPtThetaPhi;fHistInvPtPtThetaPhi=0; + if(fHistPtShift0) delete fHistPtShift0;fHistPtShift0=0; + if(fHistPrimaryVertexPosX) delete fHistPrimaryVertexPosX;fHistPrimaryVertexPosX=0; + if(fHistPrimaryVertexPosY) delete fHistPrimaryVertexPosY;fHistPrimaryVertexPosY=0; + if(fHistPrimaryVertexPosZ) delete fHistPrimaryVertexPosZ;fHistPrimaryVertexPosZ=0; + if(fHistTrackMultiplicity) delete fHistTrackMultiplicity;fHistTrackMultiplicity=0; + if(fHistTrackMultiplicityCuts) delete fHistTrackMultiplicityCuts;fHistTrackMultiplicityCuts=0; + + if(fHistTPCMomentaPosP) delete fHistTPCMomentaPosP;fHistTPCMomentaPosP=0; + if(fHistTPCMomentaNegP) delete fHistTPCMomentaNegP;fHistTPCMomentaNegP=0; + if(fHistTPCMomentaPosPt) delete fHistTPCMomentaPosPt;fHistTPCMomentaPosPt=0; + if(fHistTPCMomentaNegPt) delete fHistTPCMomentaNegPt ;fHistTPCMomentaNegPt=0; + if(fHistUserPtShift) delete fHistUserPtShift;fHistUserPtShift=0; + //esd track cuts + if(fESDTrackCuts) delete fESDTrackCuts; + //analysis folder if(fAnalysisFolder) delete fAnalysisFolder; fAnalysisFolder=0; } @@ -285,18 +292,25 @@ void AliPerformancePtCalib::Init() // momentum histos //pt shift 0 only needed if shift in 1/pt is applied - fHistPtShift0 = new TH1F("fHistPtShift0","1/pt dN/pt vs. pt of ESD track ",600,0.0,6.0); + fHistPtShift0 = new TH1F("fHistPtShift0","1/pt dN/pt vs. pt of ESD track ",800,-20.0,20.0); fList->Add(fHistPtShift0); - fHistInvPtTheta = new TH2F("fHistInvPtTheta","#theta vs 1/pt ",900, -4.5, 4.5,300,0.0,3.0); - fList->Add(fHistInvPtTheta); - fHistInvPtPhi = new TH2F("fHistInvPtPhi","#phi vs 1/pt",900, -4.5, 4.5,325,0.0,6.5); - fList->Add(fHistInvPtPhi); - fHistPtTheta = new TH2F("fHistPtTheta"," #theta vs pt ",300, 0.0, 15.0,300,0.0,3.0); - fList->Add(fHistPtTheta); - fHistPtPhi = new TH2F("fHistPtPhi"," #phi vs pt ",300, 0.0,15.0,325,0.0,6.5); - fList->Add(fHistPtPhi); - - // mom test histos + + // THnSparse for 1/pt and pt spectra vs angles + const Int_t invPtDims = 4; + fMaxPhi = 6.5; + fMinPhi = 0.0; + fMaxTheta = 3.0; + fMinTheta = 0.0; + + Double_t xminInvPt[invPtDims] = {-4.5,-20.0,fMinTheta,fMinPhi}; + Double_t xmaxInvPt[invPtDims] = {4.5,20.0,fMaxTheta,fMaxPhi}; + Int_t binsInvPt[invPtDims] = {900,800,300,325}; + + + fHistInvPtPtThetaPhi = new THnSparseF("fHistInvPtPtThetaPhi","1/pt vs pt vs #theta vs #phi ",invPtDims,binsInvPt,xminInvPt,xmaxInvPt); + fList->Add(fHistInvPtPtThetaPhi); + + // momentum test histos fHistTPCMomentaPosP = new TH2F("fHistTPCMomentaPosP","TPC p vs global esd track p pos",300,0.0,15.0,300,0.0,15.0); fList->Add(fHistTPCMomentaPosP); fHistTPCMomentaNegP = new TH2F("fHistTPCMomentaNegP","TPC p vs global esd track p neg",300,0.0,15.0,300,0.0,15.0); @@ -312,18 +326,7 @@ void AliPerformancePtCalib::Init() // esd track cuts - fESDTrackCuts = new AliESDtrackCuts("AliESDtrackCuts"); - - //fESDTrackCuts->DefineHistoqgrams(1); - fESDTrackCuts->SetRequireSigmaToVertex(fRequireSigmaToVertex); - fESDTrackCuts->SetRequireTPCRefit(fRefitTPC); - fESDTrackCuts->SetAcceptKinkDaughters(fAcceptKinkDaughters); - fESDTrackCuts->SetMinNClustersTPC((Int_t)fMinNClustersTPC); - fESDTrackCuts->SetMaxChi2PerClusterTPC(fMaxChi2PerClusterTPC); - fESDTrackCuts->SetMaxDCAToVertexXY(fMaxDCAtoVertexXY); - fESDTrackCuts->SetMaxDCAToVertexZ(fMaxDCAtoVertexZ); - fESDTrackCuts->SetDCAToVertex2D(fDCAToVertex2D); - fESDTrackCuts->SetPtRange(fMinPt,fMaxPt); + fESDTrackCuts =NULL;//neu } //________________________________________________________________________ @@ -334,34 +337,62 @@ void AliPerformancePtCalib::SetPtShift(const Double_t shiftVal ) { } //________________________________________________________________________ -void AliPerformancePtCalib::Exec(AliMCEvent*, AliESDEvent* const esdEvent, AliESDfriend*, Bool_t, Bool_t) +void AliPerformancePtCalib::Exec(AliMCEvent* const /*mcEvent*/, AliESDEvent *const esdEvent, AliESDfriend * const /*esdFriend*/, const Bool_t /*bUseMC*/, const Bool_t /*bUseESDfriend*/) { + //exec: read esd or tpc tracksGetRunNumber - //exec: read esd or tpc tracks + if(!fESDTrackCuts) Printf("no esd track cut"); if (!esdEvent) { Printf("ERROR: Event not available"); return; } - if (!(esdEvent->GetNumberOfTracks())) { - Printf(" PtCalib task: There is no track in this event"); - return; - } + if (!(esdEvent->GetNumberOfTracks())) return; if (GetTriggerClass()){ - Bool_t isEventTriggered = esdEvent->IsTriggerClassFired(GetTriggerClass()); - if(!isEventTriggered) return; + Bool_t isEventTriggered = esdEvent->IsTriggerClassFired(GetTriggerClass()); + if(!isEventTriggered) return; + } + + + // trigger selection + Bool_t isEventTriggered = kTRUE; + AliPhysicsSelection *trigSel = NULL; + AliTriggerAnalysis *trigAna = NULL; + + if(fPhysSel){ trigSel = GetPhysicsTriggerSelection(); + if(!trigSel) { + Printf("cannot get trigSel \n"); + return; + } } + //if(esdEvent) isEventTriggered = trigSel->IsCollisionCandidate(esdEvent);//does not work, crashes + if(fTrigger == AliTriggerAnalysis::kV0AND) + { + trigAna = trigSel->GetTriggerAnalysis(); + if(!trigAna) + return; + isEventTriggered = trigAna->IsOfflineTriggerFired(esdEvent,GetTrigger()); + } + + if(fPhysSel && !isEventTriggered) return; + + + //vertex info for cut + const AliESDVertex *vtx = esdEvent->GetPrimaryVertex(); + if (!vtx->GetStatus()) return ; + + //histo fo user defined shift in charge/pt if(fShift) fHistUserPtShift->Fill(fDeltaInvP); - + + //trakc multiplicity fHistTrackMultiplicity->Fill(esdEvent->GetNumberOfTracks()); // read primary vertex info Double_t tPrimaryVtxPosition[3]; - // Double_t tPrimaryVtxCov[3]; const AliESDVertex *primaryVtx = esdEvent->GetPrimaryVertexTPC(); tPrimaryVtxPosition[0] = primaryVtx->GetXv(); @@ -385,10 +416,6 @@ void AliPerformancePtCalib::Exec(AliMCEvent*, AliESDEvent* const esdEvent, AliES if(!fESDTrackCuts->AcceptTrack(esdTrack))continue; } - //track cuts - if(fRefitTPC) if(AddTPCcuts(esdTrack)) continue; - if(fRefitITS) if(AddITScuts(esdTrack)) continue; - if(fDCAcut) if(AddDCAcuts(esdTrack)) continue ; // fill histos if(fOptTPC){ //TPC tracks @@ -401,19 +428,18 @@ void AliPerformancePtCalib::Exec(AliMCEvent*, AliESDEvent* const esdEvent, AliES if(signedPt) { invPt = 1.0/signedPt; - fHistPtShift0->Fill(fabs(signedPt)); + fHistPtShift0->Fill(signedPt); - if(fShift){ + if(fShift){Printf("user shift of momentum SET to non zero value!"); invPt += fDeltaInvP; //shift momentum for tests if(invPt) signedPt = 1.0/invPt; else continue; } - - fHistInvPtTheta->Fill(invPt,tpcTrack->Theta()); - fHistInvPtPhi->Fill(invPt,tpcTrack->Phi()); - fHistPtTheta->Fill(fabs(signedPt),tpcTrack->Theta()); - fHistPtPhi->Fill(fabs(signedPt),tpcTrack->Phi()); - + Double_t theta = tpcTrack->Theta(); + Double_t phi = tpcTrack->Phi(); + + Double_t momAng[4] = {invPt,signedPt,theta,phi}; + fHistInvPtPtThetaPhi->Fill(momAng); Double_t pTPC = tpcTrack->GetP(); Double_t pESD = esdTrack->GetP(); @@ -434,78 +460,44 @@ void AliPerformancePtCalib::Exec(AliMCEvent*, AliESDEvent* const esdEvent, AliES } else{// ESD tracks + if(fabs(esdTrack->Eta())> fEtaAcceptance) continue; Double_t invPt = 0.0; Double_t signedPt = esdTrack->GetSignedPt(); if(signedPt){ invPt = 1.0/signedPt; - fHistPtShift0->Fill(fabs(signedPt)); + fHistPtShift0->Fill(signedPt); - if(fShift){ + if(fShift){Printf("user shift of momentum SET to non zero value!"); invPt += fDeltaInvP;//shift momentum for tests if(invPt) signedPt = 1.0/invPt; else continue; } - fHistInvPtTheta->Fill(invPt,esdTrack->Theta()); - fHistInvPtPhi->Fill(invPt,esdTrack->Phi()); - fHistPtTheta->Fill(signedPt,esdTrack->Theta()); - fHistPtPhi->Fill(signedPt,esdTrack->Phi()); - + Double_t theta = esdTrack->Theta(); + Double_t phi = esdTrack->Phi(); + Double_t momAng[4] = {invPt,signedPt,theta,phi}; + fHistInvPtPtThetaPhi->Fill(momAng); count++; } } } fHistTrackMultiplicityCuts->Fill(count); - } - - - -//______________________________________________________________________________________________________________________ -Bool_t AliPerformancePtCalib::AddTPCcuts(const AliESDtrack *esdTrack){ - // apply TPC cuts - - Bool_t cut = kFALSE; - - if (!(esdTrack->GetStatus()&AliESDtrack::kTPCrefit)) cut=kTRUE; // TPC refit - if (esdTrack->GetTPCNcls()GetImpactParameters(dca,cov); - if(TMath::Abs(dca[0])>fMaxDCAtoVertexXY || TMath::Abs(dca[1])>fMaxDCAtoVertexZ) cut=kTRUE; - if(esdTrack->GetKinkIndex(0)>0) cut=kTRUE; - if(cut) return kTRUE; - return kFALSE; -} - -//______________________________________________________________________________________________________________________ -Bool_t AliPerformancePtCalib::AddITScuts(const AliESDtrack *esdTrack){ - //apply ITS cuts - Bool_t cut = kFALSE; - - if (!(esdTrack->GetStatus()&AliESDtrack::kITSrefit)) cut=kTRUE; // ITS refit - Int_t clusterITS[200]; - if(esdTrack->GetITSclusters(clusterITS)<2) cut=kTRUE; // min. nb. ITS clusters //3 - - if(cut) return kTRUE; - return kFALSE; -} - //______________________________________________________________________________________________________________________ void AliPerformancePtCalib::Analyse() { // analyse charge/pt spectra in bins of theta and phi. Bins can be set by user - + THnSparseF *copyTHnSparseTheta = (THnSparseF*)fHistInvPtPtThetaPhi->Clone("copyTHnSparseTheta"); + copyTHnSparseTheta->GetAxis(3)->SetRangeUser(fMinPhi,fMaxPhi); + TH2F *histInvPtTheta = (TH2F*)copyTHnSparseTheta->Projection(2,0); + + THnSparseF *copyTHnSparsePhi = (THnSparseF*)fHistInvPtPtThetaPhi->Clone("copyTHnSparsePhi"); + copyTHnSparsePhi->GetAxis(2)->SetRangeUser(fMinTheta,fMaxTheta); + TH2F *histInvPtPhi = (TH2F*)copyTHnSparsePhi->Projection(3,0); + AliPerfAnalyzeInvPt *ana = new AliPerfAnalyzeInvPt("AliPerfAnalyzeInvPt","AliPerfAnalyzeInvPt"); @@ -514,9 +506,10 @@ void AliPerformancePtCalib::Analyse() ana->SetProjBinsTheta(fThetaBins,fNThetaBins); ana->SetProjBinsPhi(fPhiBins,fNPhiBins); ana->SetMakeFitOption(fFitGaus,fExclRange,fRange); - ana->SetDoRebin(fRebin); + if(fDoRebin) ana->SetDoRebin(fRebin); TObjArray *aFolderObj = new TObjArray; - ana->StartAnalysis(fHistInvPtTheta,fHistInvPtPhi, aFolderObj); + + ana->StartAnalysis(histInvPtTheta,histInvPtPhi,aFolderObj); // export objects to analysis folder fAnalysisFolder = ExportToFolder(aFolderObj); @@ -524,7 +517,7 @@ void AliPerformancePtCalib::Analyse() // delete only TObjArray if(aFolderObj) delete aFolderObj; if(ana) delete ana; - + } //______________________________________________________________________________________________________________________ @@ -582,11 +575,7 @@ Long64_t AliPerformancePtCalib::Merge(TCollection* const list) { AliPerformancePtCalib* entry = dynamic_cast(obj); if (entry == 0) continue; - - fHistInvPtTheta->Add(entry->fHistInvPtTheta); - fHistInvPtPhi->Add(entry-> fHistInvPtPhi); - fHistPtTheta->Add(entry->fHistPtTheta); - fHistPtPhi->Add(entry->fHistPtPhi); + fHistInvPtPtThetaPhi->Add(entry->fHistInvPtPtThetaPhi); fHistPtShift0->Add(entry->fHistPtShift0); fHistPrimaryVertexPosX->Add(entry->fHistPrimaryVertexPosX); @@ -616,7 +605,7 @@ TFolder* AliPerformancePtCalib::CreateFolder(TString name,TString title) { return folder; } //______________________________________________________________________________________________________________________ -void AliPerformancePtCalib::SetProjBinsPhi(const Double_t *phiBinArray,const Int_t nphBins){ +void AliPerformancePtCalib::SetProjBinsPhi(const Double_t *phiBinArray,const Int_t nphBins, const Double_t minTheta, const Double_t maxTheta){ // set phi bins for Analyse() //set phi bins as array and set number of this array which is equal to number of bins analysed //the last analysed bin will always be the projection from first to last bin in the array @@ -628,10 +617,20 @@ void AliPerformancePtCalib::SetProjBinsPhi(const Double_t *phiBinArray,const Int } Printf("AliPerformancePtCalib::SetProjBinsPhi: number of bins in phi set to %i",fNPhiBins); } - else Printf("Warning AliPerformancePtCalib::SetProjBinsPhi: number of bins in phi NOT set!!! Default values are taken."); + else Printf("Warning AliPerformancePtCalib::SetProjBinsPhi: number of bins in phi NOT set!!! Default values are taken."); + + if(fabs(minTheta-maxTheta)<0.001){ + Printf("AliPerformancePtCalib::SetProjBinsPhi: theta range for projection on phi and charge/pt is too small. whole range of theta selected."); + } + else{ + Printf("AliPerformancePtCalib::SetProjBinsPhi: theta range for projection on phi and charge/pt is selected by user: %1.3f - %1.3f rad",minTheta,maxTheta); + fMinTheta = minTheta; + fMaxTheta = maxTheta; + } } + //____________________________________________________________________________________________________________________________________________ -void AliPerformancePtCalib::SetProjBinsTheta(const Double_t *thetaBinArray, const Int_t nthBins){ +void AliPerformancePtCalib::SetProjBinsTheta(const Double_t *thetaBinArray, const Int_t nthBins, const Double_t minPhi, const Double_t maxPhi){ // set theta bins for Analyse() //set theta bins as array and set number of this array which is equal to number of bins analysed //the last analysed bin will always be the projection from first to last bin in the array @@ -640,9 +639,18 @@ void AliPerformancePtCalib::SetProjBinsTheta(const Double_t *thetaBinArray, cons for(Int_t k = 0;kLoadMacro("$ALICE_ROOT/PWG1/Macros/LoadMyLibs.C"); -LoadMyLibs(); +// after running the performance task, read the file, and get component TFile f("Output.root"); -AliPerformancePtCalib * compObj = (AliPerformancePtCalibMC*)coutput->FindObject("AliPerformancePtCalibMC"); +AliPerformancePtCalibMC *compObj = (AliPerformancePtCalibMC*)coutput->FindObject("AliPerformancePtCalibMC"); -// analyse comparison data +// set phi and theta bins for fitting and analyse comparison data +compObj->SetProjBinsTheta(thetaBins,nThetaBins,minPhi, maxPhi); +compObj->SetProjBinsPhi(phiBins,nPhiBins,minTheta,maxTheta); +compObj->SetMakeFitOption(kFALSE,exclRange,fitRange); +compObj->SetDoRebin(rebin); +//compObj->SetAnaMCOff(); compObj->Analyse(); +//details see functions of this class // the output histograms/graphs will be stored in the folder "folderRes" compObj->GetAnalysisFolder()->ls("*"); @@ -39,6 +57,7 @@ fout.Close(); #include "TH1F.h" #include "TH2F.h" +#include "THnSparse.h" #include "TList.h" #include "TMath.h" #include "TFolder.h" @@ -51,6 +70,9 @@ fout.Close(); #include "AliESDfriendTrack.h" #include "AliESDfriend.h" +#include "AliPhysicsSelection.h" +#include "AliTriggerAnalysis.h" + #include "AliPerformancePtCalibMC.h" #include "AliPerfAnalyzeInvPt.h" @@ -65,9 +87,15 @@ ClassImp(AliPerformancePtCalibMC) // option parameter for AliPerformancePtCalibMC::Analyse() fNThetaBins(0), fNPhiBins(0), + fMaxPhi(0), + fMinPhi(0), + fMaxTheta(0), + fMinTheta(0), fRange(0), fExclRange(0), fFitGaus(0) , + fDoRebin(0), + fRebin(0), fAnaMC(0), // option parameter for user defined charge/pt shift fShift(0), @@ -75,29 +103,16 @@ ClassImp(AliPerformancePtCalibMC) //options for cuts fOptTPC(0), fESDcuts(0), - fRefitTPC(0), - fRefitITS(0), - fDCAcut(0), + fPhysSel(0), fEtaAcceptance(0), - fMinPt(0), - fMaxPt(0), - fMinNClustersTPC(0), - fMaxChi2PerClusterTPC(0), - fMaxDCAtoVertexXY(0), - fMaxDCAtoVertexZ(0), - fAcceptKinkDaughters(0), - fRequireSigmaToVertex(0), - fDCAToVertex2D(0), - + fTrigger(AliTriggerAnalysis::kMB1), + fPhysicsSelection(0), fCutsRC(0), fCutsMC(0), - fList(0), + // histograms - fHistInvPtTheta(0), - fHistInvPtPhi(0), - fHistPtTheta(0), - fHistPtPhi(0), + fHistInvPtPtThetaPhi(0), fHistPtShift0(0), fHistPrimaryVertexPosX(0), fHistPrimaryVertexPosY(0), @@ -108,10 +123,7 @@ ClassImp(AliPerformancePtCalibMC) fHistTPCMomentaNegP(0), fHistTPCMomentaPosPt(0), fHistTPCMomentaNegPt(0), - fHistInvPtThetaMC(0), - fHistInvPtPhiMC(0), - fHistPtThetaMC(0), - fHistPtPhiMC(0), + fHistInvPtPtThetaPhiMC(0), fHistInvPtMCESD(0), fHistInvPtMCTPC(0), fHistPtMCESD(0), @@ -127,6 +139,8 @@ ClassImp(AliPerformancePtCalibMC) fHistESDMomentaPosPtMC(0), fHistESDMomentaNegPtMC(0), fHistUserPtShift(0), + + //esd track cuts fESDTrackCuts(0), // analysis folder fAnalysisFolder(0) @@ -137,142 +151,125 @@ ClassImp(AliPerformancePtCalibMC) fShift = kFALSE; // shift in charge/pt yes/no fDeltaInvP = 0.00; // shift value //options for cuts + fOptTPC = kTRUE; // read TPC tracks yes/no - fESDcuts = kTRUE; // read ESD track cuts - fRefitTPC = kFALSE; // require TPC refit - fRefitITS = kFALSE; // require ITS refit - fDCAcut = kTRUE; // apply DCA cuts - + fESDcuts = kFALSE; + fPhysSel = kFALSE; fCutsRC = NULL; fCutsMC = NULL; + fEtaAcceptance = 0.8; - fMinPt=0.15; // GeV/c - fMaxPt=1.e10; // GeV/c - fMinNClustersTPC = 50; - fMaxChi2PerClusterTPC = 4.0; - fMaxDCAtoVertexXY = 2.4; // cm - fMaxDCAtoVertexZ = 3.0; // cm - fAcceptKinkDaughters = kFALSE; - fRequireSigmaToVertex = kFALSE; - fDCAToVertex2D = kTRUE; // options for function AliPerformancePtCalibMC::Analyse() fFitGaus = kFALSE;// use gaussian function for fitting charge/pt yes/no fNThetaBins = 0; //number of theta bins fNPhiBins = 0; //number of phi bins + fMaxPhi = 6.5;// max phi for 2D projection on theta and charge/pt axis + fMinPhi = 0.0;// min phi for 2D projection on theta and charge/pt axis + fMaxTheta = 3.0;// max theta for 2D projection on phi and charge/pt axis + fMinTheta = 0.0;// min theta for 2D projection on phi and charge/pt axis fRange = 0; //fit range around 0 fExclRange =0; //range of rejection of points around 0 fAnaMC = kTRUE; // analyse MC tracks yes/no + fDoRebin = kFALSE; + fRebin = 0; Init(); } //________________________________________________________________________ -AliPerformancePtCalibMC::AliPerformancePtCalibMC(const char *name= "AliPerformancePtCalibMC", const char *title="AliPerformancePtCalibMC")://,Int_t analysisMode=0,Bool_t hptGenerator=kFALSE) : +AliPerformancePtCalibMC::AliPerformancePtCalibMC(const char *name= "AliPerformancePtCalibMC", const char *title="AliPerformancePtCalibMC"): AliPerformanceObject(name,title), // option parameter for AliPerformancePtCalibMC::Analyse() fNThetaBins(0), - fNPhiBins(0), - fRange(0), - fExclRange(0), - fFitGaus(0) , - fAnaMC(0), - // option parameter for user defined 1/pt shift - fShift(0), - fDeltaInvP (0), - //options for cuts - fOptTPC(0), - fESDcuts(0), - fRefitTPC(0), - fRefitITS(0), - fDCAcut(0), - fEtaAcceptance(0), - fMinPt(0), - fMaxPt(0), - fMinNClustersTPC(0), - fMaxChi2PerClusterTPC(0), - fMaxDCAtoVertexXY(0), - fMaxDCAtoVertexZ(0), - fAcceptKinkDaughters(0), - fRequireSigmaToVertex(0), - fDCAToVertex2D(0), - - fCutsRC(0), - fCutsMC(0), - - fList(0), - // histograms - fHistInvPtTheta(0), - fHistInvPtPhi(0), - fHistPtTheta(0), - fHistPtPhi(0), - fHistPtShift0(0), - fHistPrimaryVertexPosX(0), - fHistPrimaryVertexPosY(0), - fHistPrimaryVertexPosZ(0), - fHistTrackMultiplicity(0), - fHistTrackMultiplicityCuts(0), - fHistTPCMomentaPosP(0), - fHistTPCMomentaNegP(0), - fHistTPCMomentaPosPt(0), - fHistTPCMomentaNegPt(0), - fHistInvPtThetaMC(0), - fHistInvPtPhiMC(0), - fHistPtThetaMC(0), - fHistPtPhiMC(0), - fHistInvPtMCESD(0), - fHistInvPtMCTPC(0), - fHistPtMCESD(0), - fHistPtMCTPC(0), - fHistMomresMCESD(0), - fHistMomresMCTPC(0), - fHistTPCMomentaPosInvPtMC(0), - fHistTPCMomentaNegInvPtMC(0), - fHistTPCMomentaPosPtMC(0), - fHistTPCMomentaNegPtMC(0), - fHistESDMomentaPosInvPtMC(0), - fHistESDMomentaNegInvPtMC(0), - fHistESDMomentaPosPtMC(0), - fHistESDMomentaNegPtMC(0), - fHistUserPtShift(0), - - fESDTrackCuts(0), - // analysis folder - fAnalysisFolder(0) - + fNPhiBins(0), + fMaxPhi(0), + fMinPhi(0), + fMaxTheta(0), + fMinTheta(0), + fRange(0), + fExclRange(0), + fFitGaus(0) , + fDoRebin(0), + fRebin(0), + fAnaMC(0), + // option parameter for user defined charge/pt shift + fShift(0), + fDeltaInvP(0), + //options for cuts + fOptTPC(0), + fESDcuts(0), + fPhysSel(0), + fEtaAcceptance(0), + fTrigger(AliTriggerAnalysis::kMB1), + fPhysicsSelection(0), + fCutsRC(0), + fCutsMC(0), + fList(0), + + // histograms + fHistInvPtPtThetaPhi(0), + fHistPtShift0(0), + fHistPrimaryVertexPosX(0), + fHistPrimaryVertexPosY(0), + fHistPrimaryVertexPosZ(0), + fHistTrackMultiplicity(0), + fHistTrackMultiplicityCuts(0), + fHistTPCMomentaPosP(0), + fHistTPCMomentaNegP(0), + fHistTPCMomentaPosPt(0), + fHistTPCMomentaNegPt(0), + fHistInvPtPtThetaPhiMC(0), + fHistInvPtMCESD(0), + fHistInvPtMCTPC(0), + fHistPtMCESD(0), + fHistPtMCTPC(0), + fHistMomresMCESD(0), + fHistMomresMCTPC(0), + fHistTPCMomentaPosInvPtMC(0), + fHistTPCMomentaNegInvPtMC(0), + fHistTPCMomentaPosPtMC(0), + fHistTPCMomentaNegPtMC(0), + fHistESDMomentaPosInvPtMC(0), + fHistESDMomentaNegInvPtMC(0), + fHistESDMomentaPosPtMC(0), + fHistESDMomentaNegPtMC(0), + fHistUserPtShift(0), + + //esd track cuts + fESDTrackCuts(0), + // analysis folder + fAnalysisFolder(0) { - // Constructor - + // Dummy constructor + + fShift = kFALSE; // shift in charge/pt yes/no fDeltaInvP = 0.00; // shift value //options for cuts + fOptTPC = kTRUE; // read TPC tracks yes/no - fESDcuts = kTRUE; // read ESD track cuts - fRefitTPC = kFALSE; // require TPC refit - fRefitITS = kFALSE; // require ITS refit - fDCAcut = kTRUE; // apply DCA cuts - + fESDcuts = kFALSE; + fPhysSel = kFALSE; fCutsRC = NULL; fCutsMC = NULL; + fEtaAcceptance = 0.8; - fMinPt=0.15; // GeV/c - fMaxPt=1.e10; // GeV/c - fMinNClustersTPC = 50; - fMaxChi2PerClusterTPC = 4.0; - fMaxDCAtoVertexXY = 2.4; // cm - fMaxDCAtoVertexZ = 3.0; // cm - fAcceptKinkDaughters = kFALSE; - fRequireSigmaToVertex = kFALSE; - fDCAToVertex2D = kTRUE; // options for function AliPerformancePtCalibMC::Analyse() fFitGaus = kFALSE;// use gaussian function for fitting charge/pt yes/no fNThetaBins = 0; //number of theta bins fNPhiBins = 0; //number of phi bins + fMaxPhi = 6.5;// max phi for 2D projection on theta and charge/pt axis + fMinPhi = 0.0;// min phi for 2D projection on theta and charge/pt axis + fMaxTheta = 3.0;// max theta for 2D projection on phi and charge/pt axis + fMinTheta = 0.0;// min theta for 2D projection on phi and charge/pt axis fRange = 0; //fit range around 0 fExclRange =0; //range of rejection of points around 0 fAnaMC = kTRUE; // analyse MC tracks yes/no - + fDoRebin = kFALSE;// flag for rebin + fRebin = 0;// bins for rebin + Init(); } @@ -281,7 +278,46 @@ AliPerformancePtCalibMC::~AliPerformancePtCalibMC() { // // destructor // - if(fAnalysisFolder) delete fAnalysisFolder; fAnalysisFolder=0; + + if(fList) delete fList; + // histograms + if(fHistInvPtPtThetaPhi) delete fHistInvPtPtThetaPhi;fHistInvPtPtThetaPhi=0; + if(fHistPtShift0) delete fHistPtShift0;fHistPtShift0=0; + if(fHistPrimaryVertexPosX) delete fHistPrimaryVertexPosX;fHistPrimaryVertexPosX=0; + if(fHistPrimaryVertexPosY) delete fHistPrimaryVertexPosY;fHistPrimaryVertexPosY=0; + if(fHistPrimaryVertexPosZ) delete fHistPrimaryVertexPosZ;fHistPrimaryVertexPosZ=0; + if(fHistTrackMultiplicity) delete fHistTrackMultiplicity;fHistTrackMultiplicity=0; + if(fHistTrackMultiplicityCuts) delete fHistTrackMultiplicityCuts;fHistTrackMultiplicityCuts=0; + + if(fHistTPCMomentaPosP) delete fHistTPCMomentaPosP;fHistTPCMomentaPosP=0; + if(fHistTPCMomentaNegP) delete fHistTPCMomentaNegP;fHistTPCMomentaNegP=0; + if(fHistTPCMomentaPosPt) delete fHistTPCMomentaPosPt;fHistTPCMomentaPosPt=0; + if(fHistTPCMomentaNegPt) delete fHistTPCMomentaNegPt ;fHistTPCMomentaNegPt=0; + if(fHistUserPtShift) delete fHistUserPtShift;fHistUserPtShift=0; + if(fHistInvPtPtThetaPhiMC) delete fHistInvPtPtThetaPhiMC;fHistInvPtPtThetaPhiMC=0; + if(fHistInvPtMCESD) delete fHistInvPtMCESD;fHistInvPtMCESD=0; + if(fHistInvPtMCTPC) delete fHistInvPtMCTPC;fHistInvPtMCTPC=0; + if(fHistPtMCESD) delete fHistPtMCESD;fHistPtMCESD=0; + if(fHistPtMCTPC) delete fHistPtMCTPC;fHistPtMCTPC=0; + if(fHistMomresMCESD) delete fHistMomresMCESD;fHistMomresMCESD=0; + if(fHistMomresMCTPC) delete fHistMomresMCTPC;fHistMomresMCTPC=0; + if(fHistTPCMomentaPosInvPtMC) delete fHistTPCMomentaPosInvPtMC;fHistTPCMomentaPosInvPtMC=0; + if(fHistTPCMomentaNegInvPtMC) delete fHistTPCMomentaNegInvPtMC;fHistTPCMomentaNegInvPtMC=0; + if(fHistTPCMomentaPosPtMC) delete fHistTPCMomentaPosPtMC;fHistTPCMomentaPosPtMC=0; + if(fHistTPCMomentaNegPtMC) delete fHistTPCMomentaNegPtMC;fHistTPCMomentaNegPtMC=0; + if(fHistESDMomentaPosInvPtMC) delete fHistESDMomentaPosInvPtMC;fHistESDMomentaPosInvPtMC=0; + if(fHistESDMomentaNegInvPtMC) delete fHistESDMomentaNegInvPtMC;fHistESDMomentaNegInvPtMC=0; + if(fHistESDMomentaPosPtMC) delete fHistESDMomentaPosPtMC;fHistESDMomentaPosPtMC=0; + if(fHistESDMomentaNegPtMC) delete fHistESDMomentaNegPtMC;fHistESDMomentaNegPtMC=0; + + + + //esd track cuts + if(fESDTrackCuts) delete fESDTrackCuts; + if(fAnalysisFolder) delete fAnalysisFolder; fAnalysisFolder=0; + + + } //________________________________________________________________________ @@ -310,18 +346,21 @@ void AliPerformancePtCalibMC::Init() fList->Add(fHistTrackMultiplicityCuts); // momentum histos - fHistPtShift0 = new TH1F("fHistPtShift0","1/pt dN/pt vs. pt of ESD track ",600,0.0,6.0); + fHistPtShift0 = new TH1F("fHistPtShift0","1/pt dN/pt vs. pt of ESD track ",800,-20.0,20.0); fList->Add(fHistPtShift0); - fHistInvPtTheta = new TH2F("fHistInvPtTheta","#theta vs 1/pt ",900, -4.5, 4.5,300,0.0,3.0); - fList->Add(fHistInvPtTheta); - fHistInvPtPhi = new TH2F("fHistInvPtPhi","#phi vs 1/pt",900, -4.5, 4.5,325,0.0,6.5); - fList->Add(fHistInvPtPhi); - fHistPtTheta = new TH2F("fHistPtTheta"," #theta vs pt ",300, 0.0, 15.0,300,0.0,3.0); - fList->Add(fHistPtTheta); - fHistPtPhi = new TH2F("fHistPtPhi"," #phi vs pt ",300, 0.0,15.0,325,0.0,6.5); - fList->Add(fHistPtPhi); - - // mom test histos + const Int_t invPtDims = 4; + fMaxPhi=6.5; + fMinPhi=0.0; + fMaxTheta=3.0; + fMinTheta=0.0; + Double_t xminInvPt[invPtDims] = {-4.5,-20.0,fMinTheta,fMinPhi}; + Double_t xmaxInvPt[invPtDims] = {4.5,20.0,fMaxTheta,fMaxPhi}; + Int_t binsInvPt[invPtDims] = {900,800,300,325}; + + fHistInvPtPtThetaPhi = new THnSparseF("fHistInvPtPtThetaPhi","1/pt vs pt vs #theta vs #phi ",invPtDims,binsInvPt,xminInvPt,xmaxInvPt); + fList->Add(fHistInvPtPtThetaPhi); + + // momentum test histos fHistTPCMomentaPosP = new TH2F("fHistTPCMomentaPosP","TPC p vs global esd track p pos",300,0.0,15.0,300,0.0,15.0); fList->Add(fHistTPCMomentaPosP); fHistTPCMomentaNegP = new TH2F("fHistTPCMomentaNegP","TPC p vs global esd track p neg",300,0.0,15.0,300,0.0,15.0); @@ -331,7 +370,7 @@ void AliPerformancePtCalibMC::Init() fHistTPCMomentaNegPt = new TH2F("fHistTPCMomentaNegPt","TPC pt vs global esd track pt neg",300,0.0,15.0,300,0.0,15.0); fList->Add(fHistTPCMomentaNegPt); - // mom test histos MC + // momentum test histos MC 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); fList->Add(fHistTPCMomentaPosInvPtMC); 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); @@ -349,15 +388,10 @@ void AliPerformancePtCalibMC::Init() fHistESDMomentaNegPtMC = new TH1F("fHistESDMomentaNegPtMC","ESD-MC of pt ",600,-4.0,44.0); fList->Add(fHistESDMomentaNegPtMC); - // MC info - fHistInvPtThetaMC = new TH2F("fHistInvPtThetaMC","theta vs inv pt MC",900, -4.5, 4.5,300,0.0,3.0); - fList->Add(fHistInvPtThetaMC); - fHistInvPtPhiMC = new TH2F("fHistInvPtPhiMC","phi vs inv pt MC",900, -4.5, 4.5,325,0.0,6.5); - fList->Add(fHistInvPtPhiMC); - fHistPtThetaMC = new TH2F("fHistPtThetaMC","theta vs pt MC",300, 0.0, 15.0,300,0.0,3.0); - fList->Add(fHistPtThetaMC); - fHistPtPhiMC = new TH2F("fHistPtPhiMC"," phi vs pt MC",300, 0.0,15.0,325,0.0,6.5); - fList->Add(fHistPtPhiMC); + // MC only info + fHistInvPtPtThetaPhiMC = new THnSparseF("fHistInvPtPtThetaPhiMC","MC 1/pt vs pt vs #theta vs #phi ",invPtDims,binsInvPt,xminInvPt,xmaxInvPt); + fList->Add(fHistInvPtPtThetaPhiMC); + //correlation histos MC ESD or TPC fHistInvPtMCESD = new TH2F("fHistInvPtMCESD","inv pt ESD vs MC",900, 0.0, 9.0,900, 0.0, 9.0); @@ -380,19 +414,7 @@ void AliPerformancePtCalibMC::Init() // esd track cuts fESDTrackCuts = new AliESDtrackCuts("AliESDtrackCuts"); - - // //fESDTrackCuts->DefineHistoqgrams(1); - - fESDTrackCuts->SetRequireSigmaToVertex(fRequireSigmaToVertex); - fESDTrackCuts->SetRequireTPCRefit(fRefitTPC); - fESDTrackCuts->SetAcceptKinkDaughters(fAcceptKinkDaughters); - fESDTrackCuts->SetMinNClustersTPC((Int_t)fMinNClustersTPC); - fESDTrackCuts->SetMaxChi2PerClusterTPC(fMaxChi2PerClusterTPC); - fESDTrackCuts->SetMaxDCAToVertexXY(fMaxDCAtoVertexXY); - fESDTrackCuts->SetMaxDCAToVertexZ(fMaxDCAtoVertexZ); - fESDTrackCuts->SetDCAToVertex2D(fDCAToVertex2D); - fESDTrackCuts->SetPtRange(fMinPt,fMaxPt); - + } @@ -405,11 +427,11 @@ void AliPerformancePtCalibMC::SetPtShift(const Double_t shiftVal ) { } //________________________________________________________________________ -void AliPerformancePtCalibMC::Exec(AliMCEvent* const mcEvent, AliESDEvent* const esdEvent , AliESDfriend*, Bool_t, Bool_t) +void AliPerformancePtCalibMC::Exec(AliMCEvent* const mcEvent, AliESDEvent *const esdEvent, AliESDfriend *const /*esdFriend*/, const Bool_t /*bUseMC*/, const Bool_t /*bUseESDfriend*/) { //exec: read MC and esd or tpc tracks - AliStack* stack; + AliStack* stack = NULL; if (!esdEvent) { Printf("ERROR: Event not available"); @@ -432,6 +454,46 @@ void AliPerformancePtCalibMC::Exec(AliMCEvent* const mcEvent, AliESDEvent* const return; } + if (GetTriggerClass()){ + Bool_t isEventTriggered = esdEvent->IsTriggerClassFired(GetTriggerClass()); + if(!isEventTriggered) return; + } + + /* + // trigger selection + Bool_t isEventTriggered = kTRUE; + AliPhysicsSelection *trigSel = NULL; + AliTriggerAnalysis *trigAna = NULL; + + if(fPhysSel){ trigSel = GetPhysicsTriggerSelection(); + if(!trigSel) { + Printf("cannot get trigSel \n"); + return; + } + } + +// if(esdEvent) isEventTriggered = trigSel->IsCollisionCandidate(esdEvent);// crashes + + if(fTrigger == AliTriggerAnalysis::kV0AND) + {Printf("physics selection V0AND");//test + trigAna = trigSel->GetTriggerAnalysis(); + if(!trigAna) + return; + + isEventTriggered = trigAna->IsOfflineTriggerFired(esdEvent,GetTrigger()); + + } + + if(fPhysSel && !isEventTriggered) return; + + + //vertex info for cut + const AliESDVertex *vtx = esdEvent->GetPrimaryVertex(); + if (!vtx->GetStatus()) return ; + */ + + + if(fShift) fHistUserPtShift->Fill(fDeltaInvP); // read primary vertex info @@ -456,15 +518,10 @@ void AliPerformancePtCalibMC::Exec(AliMCEvent* const mcEvent, AliESDEvent* const if(!esdTrack) continue; //esd track cuts - if(fESDcuts){Printf("esd cuts aplied"); + if(fESDcuts){ if(!fESDTrackCuts->AcceptTrack(esdTrack)) continue; } - - //more track cuts - if(fRefitTPC) if(AddTPCcuts(esdTrack)) continue; - if(fRefitITS) if(AddITScuts(esdTrack)) continue; - if(fDCAcut) if(AddDCAcuts(esdTrack)) continue ; - + // get MC info Int_t esdLabel = esdTrack->GetLabel(); if(esdLabel<0) continue; @@ -483,18 +540,22 @@ void AliPerformancePtCalibMC::Exec(AliMCEvent* const mcEvent, AliESDEvent* const if(signMC>0) signMC = 1; else signMC = -1; - //fill MC histos - fHistInvPtThetaMC->Fill(signMC*fabs(invPtMC),partMC->Theta()); - fHistInvPtPhiMC->Fill(signMC*fabs(invPtMC),partMC->Phi()); - fHistPtThetaMC->Fill(fabs(mcPt),partMC->Theta(),fabs(invPtMC)); - fHistPtPhiMC->Fill(fabs(mcPt),partMC->Phi(),fabs(invPtMC)); - - //correlation histos MC ESD - fHistInvPtMCESD->Fill(fabs(invPtMC),fabs(1.0/ptESD)); - fHistPtMCESD->Fill(fabs(mcPt),fabs(ptESD)); - + //fill MC information + Double_t thetaMC = partMC->Theta(); + Double_t phiMC = partMC->Phi(); + + Double_t momAngMC[4] = {signMC*(fabs(invPtMC)),signMC*(fabs(mcPt)),thetaMC,phiMC}; + + // fill only if MC track is in eta acceptance of TPC in order to be compareable to TPC tracks + if(fabs( partMC->Eta())<= fEtaAcceptance) { + fHistInvPtPtThetaPhiMC->Fill(momAngMC); + + //correlation histos MC ESD + fHistInvPtMCESD->Fill(fabs(invPtMC),fabs(1.0/ptESD)); + fHistPtMCESD->Fill(fabs(mcPt),fabs(ptESD)); + } - // fill histos + // fill histos TPC or ESD if(fOptTPC){ //TPC tracks and MC tracks const AliExternalTrackParam *tpcTrack = esdTrack->GetTPCInnerParam(); @@ -506,19 +567,20 @@ void AliPerformancePtCalibMC::Exec(AliMCEvent* const mcEvent, AliESDEvent* const if(signedPt) { invPt = 1.0/signedPt; - fHistPtShift0->Fill(fabs(signedPt)); + fHistPtShift0->Fill(signedPt);//changed - if(fShift){ + if(fShift){Printf("user shift of momentum SET to non zero value!"); invPt += fDeltaInvP; //shift momentum for tests if(invPt) signedPt = 1.0/invPt; else continue; } - fHistInvPtTheta->Fill(invPt,tpcTrack->Theta()); - fHistInvPtPhi->Fill(invPt,tpcTrack->Phi()); - fHistPtTheta->Fill(fabs(signedPt),tpcTrack->Theta()); - fHistPtPhi->Fill(fabs(signedPt),tpcTrack->Phi()); + Double_t theta = tpcTrack->Theta(); + Double_t phi = tpcTrack->Phi(); + + Double_t momAng[4] = {invPt,signedPt,theta,phi}; + fHistInvPtPtThetaPhi->Fill(momAng); //correlation histos MC TPC fHistInvPtMCTPC->Fill(fabs(invPtMC),fabs(invPt)); @@ -529,7 +591,7 @@ void AliPerformancePtCalibMC::Exec(AliMCEvent* const mcEvent, AliESDEvent* const Double_t ptDiffTPC = (fabs(signedPt)-fabs(mcPt))/pow(mcPt,2); Double_t invPtDiffESD = fabs(1.0/ptESD)-1.0/fabs(mcPt); Double_t invPtDiffTPC = fabs(invPt)-1.0/fabs(mcPt); - Double_t pTPC = tpcTrack->GetP(); + Double_t pTPC = tpcTrack->GetP(); if(esdTrack->GetSign()>0){//compare momenta ESD track and TPC track fHistTPCMomentaPosP->Fill(fabs(pESD),fabs(pTPC)); @@ -552,21 +614,24 @@ void AliPerformancePtCalibMC::Exec(AliMCEvent* const mcEvent, AliESDEvent* const else{ // ESD tracks and MC tracks + if(fabs(esdTrack->Eta())> fEtaAcceptance) continue; Double_t invPt = 0.0; if(ptESD) { invPt = 1.0/ptESD; - fHistPtShift0->Fill(fabs(ptESD)); + fHistPtShift0->Fill(ptESD);//changed - if(fShift){ + if(fShift){Printf("user shift of momentum SET to non zero value!"); invPt += fDeltaInvP; //shift momentum for tests if(invPt) ptESD = 1.0/invPt; else continue; } - fHistInvPtTheta->Fill(invPt,esdTrack->Theta()); - fHistInvPtPhi->Fill(invPt,esdTrack->Phi()); - fHistPtTheta->Fill(ptESD,esdTrack->Theta()); - fHistPtPhi->Fill(ptESD,esdTrack->Phi()); + + Double_t theta = esdTrack->Theta(); + Double_t phi = esdTrack->Phi(); + + Double_t momAng[4] = {invPt,ptESD,theta,phi}; + fHistInvPtPtThetaPhi->Fill(momAng); //differences MC ESD tracks Double_t ptDiffESD = (fabs(ptESD)-fabs(mcPt))/pow(mcPt,2); @@ -590,43 +655,6 @@ void AliPerformancePtCalibMC::Exec(AliMCEvent* const mcEvent, AliESDEvent* const } -//______________________________________________________________________________________________________________________ -Bool_t AliPerformancePtCalibMC::AddTPCcuts(const AliESDtrack *esdTrack){ - // apply TPC cuts - - Bool_t cut = kFALSE; - - if (!(esdTrack->GetStatus()&AliESDtrack::kTPCrefit)) cut=kTRUE; // TPC refit - if (esdTrack->GetTPCNcls()GetImpactParameters(dca,cov); - if(TMath::Abs(dca[0])>fMaxDCAtoVertexXY || TMath::Abs(dca[1])>fMaxDCAtoVertexZ) cut=kTRUE; - if(esdTrack->GetKinkIndex(0)>0) cut=kTRUE; - if(cut) return kTRUE; - return kFALSE; -} - -//______________________________________________________________________________________________________________________ -Bool_t AliPerformancePtCalibMC::AddITScuts(const AliESDtrack *esdTrack){ - //apply ITS cuts - Bool_t cut = kFALSE; - - if (!(esdTrack->GetStatus()&AliESDtrack::kITSrefit)) cut=kTRUE; // ITS refit - Int_t clusterITS[200]; - if(esdTrack->GetITSclusters(clusterITS)<2) cut=kTRUE; // min. nb. ITS clusters //3 - - if(cut) return kTRUE; - return kFALSE; -} - //______________________________________________________________________________________________________________________ void AliPerformancePtCalibMC::Analyse() @@ -634,6 +662,27 @@ void AliPerformancePtCalibMC::Analyse() // analyse charge/pt spectra in bins of theta and phi. Bins can be set by user + + THnSparseF *copyTHnSparseTheta; + THnSparseF *copyTHnSparsePhi; + + if(fAnaMC){ + Printf("AliPerformancePtCalibMC::Analyse: analysing MC!"); + copyTHnSparseTheta = (THnSparseF*)fHistInvPtPtThetaPhiMC->Clone("copyTHnSparseThetaMC"); + copyTHnSparsePhi = (THnSparseF*)fHistInvPtPtThetaPhiMC->Clone("copyTHnSparsePhiMC"); + } + else { + Printf("AliPerformancePtCalibMC::Analyse: analysing ESD (reco)!"); + copyTHnSparseTheta = (THnSparseF*)fHistInvPtPtThetaPhi->Clone("copyTHnSparseTheta"); + copyTHnSparsePhi = (THnSparseF*)fHistInvPtPtThetaPhi->Clone("copyTHnSparsePhi"); + } + + copyTHnSparseTheta->GetAxis(3)->SetRangeUser(fMinPhi,fMaxPhi); + copyTHnSparsePhi->GetAxis(2)->SetRangeUser(fMinTheta,fMaxTheta); + + TH2F *histInvPtTheta = (TH2F*)copyTHnSparseTheta->Projection(2,0); + TH2F *histInvPtPhi = (TH2F*)copyTHnSparsePhi->Projection(3,0); + AliPerfAnalyzeInvPt *ana = new AliPerfAnalyzeInvPt("AliPerfAnalyzeInvPt","AliPerfAnalyzeInvPt"); TH1::AddDirectory(kFALSE); @@ -641,17 +690,11 @@ void AliPerformancePtCalibMC::Analyse() ana->SetProjBinsTheta(fThetaBins,fNThetaBins); ana->SetProjBinsPhi(fPhiBins,fNPhiBins); ana->SetMakeFitOption(fFitGaus,fExclRange,fRange); - + if(fDoRebin) ana->SetDoRebin(fRebin); TObjArray *aFolderObj = new TObjArray; - if(fAnaMC){ - Printf("AliPerformancePtCalibMC: analysing MC!"); - ana->StartAnalysis(fHistInvPtThetaMC,fHistInvPtPhiMC, aFolderObj); - } - else { - Printf("AliPerformancePtCalibMC: analysing data!"); - ana->StartAnalysis(fHistInvPtTheta,fHistInvPtPhi, aFolderObj); - } - + + ana->StartAnalysis(histInvPtTheta,histInvPtPhi, aFolderObj); + // export objects to analysis folder fAnalysisFolder = ExportToFolder(aFolderObj); @@ -716,16 +759,10 @@ Long64_t AliPerformancePtCalibMC::Merge(TCollection* const list) { AliPerformancePtCalibMC* entry = dynamic_cast(obj); if (!entry) continue; - - fHistInvPtTheta->Add(entry->fHistInvPtTheta); - fHistInvPtPhi->Add(entry-> fHistInvPtPhi); - fHistPtTheta->Add(entry->fHistPtTheta); - fHistPtPhi->Add(entry->fHistPtPhi); - - fHistInvPtThetaMC->Add(entry->fHistInvPtThetaMC); - fHistInvPtPhiMC->Add(entry->fHistInvPtPhiMC); - fHistPtThetaMC->Add(entry->fHistPtThetaMC); - fHistPtPhiMC->Add(entry->fHistPtPhiMC); + fHistInvPtPtThetaPhi->Add(entry->fHistInvPtPtThetaPhi); + + fHistInvPtPtThetaPhiMC->Add(entry->fHistInvPtPtThetaPhiMC); + fHistInvPtMCESD->Add(entry->fHistInvPtMCESD); fHistPtMCESD->Add(entry->fHistPtMCESD); fHistInvPtMCTPC->Add(entry->fHistInvPtMCTPC); @@ -770,61 +807,67 @@ TFolder* AliPerformancePtCalibMC::CreateFolder(TString name,TString title) { // set variables for Analyse() - //______________________________________________________________________________________________________________________ -void AliPerformancePtCalibMC::SetProjBinsPhi(const Double_t *phiBinArray,const Int_t nphBins){ - +void AliPerformancePtCalibMC::SetProjBinsPhi(const Double_t *phiBinArray,const Int_t nphBins, const Double_t minTheta, const Double_t maxTheta){ // set phi bins for Analyse() - //set phi bins as array and set number of this array which is equal to number of bins analysed - //the last analysed bin will always be the projection from first to last bin in the array + // set phi bins as array and set number of this array which is equal to number of bins analysed + // the last analysed bin will always be the projection from first to last bin in the array if(nphBins){ fNPhiBins = nphBins; for(Int_t k = 0;k