#include "TArrayF.h"
#include "TF1.h"
#include "TRandom.h"
+#include "TROOT.h"
#include "AliAnalysisTaskSE.h"
#include "AliAnalysisManager.h"
AliAnalysisTaskBFPsi::AliAnalysisTaskBFPsi(const char *name)
: AliAnalysisTaskSE(name),
fDebugLevel(kFALSE),
- fArrayMC(0), //+++++++++++++
+ fArrayMC(0),
fBalance(0),
fRunShuffling(kFALSE),
fShuffledBalance(0),
fHistProbTPCTOFvsPtbeforePID(NULL),
fHistNSigmaTPCvsPtbeforePID(NULL),
fHistNSigmaTOFvsPtbeforePID(NULL),
- fHistBetaVsdEdXbeforePID(NULL), //+++++++
- fHistNSigmaTPCTOFvsPtbeforePID(NULL), //++++++
+ fHistBetaVsdEdXbeforePID(NULL),
+ fHistNSigmaTPCTOFvsPtbeforePID(NULL),
+ fHistNSigmaTPCTOFPbefPID(NULL),
fHistdEdxVsPTPCafterPID(NULL),
fHistBetavsPTOFafterPID(NULL),
fHistProbTPCvsPtafterPID(NULL),
fHistProbTPCTOFvsPtafterPID(NULL),
fHistNSigmaTPCvsPtafterPID(NULL),
fHistNSigmaTOFvsPtafterPID(NULL),
- fHistBetaVsdEdXafterPID(NULL), //+++++++
- fHistNSigmaTPCTOFvsPtafterPID(NULL), //+++++++
- fHistdEdxVsPTPCbeforePIDelectron(NULL), //+++++++
- fHistNSigmaTPCvsPtbeforePIDelectron(NULL), //+++++++
- fHistdEdxVsPTPCafterPIDelectron(NULL), //+++++++
- fHistNSigmaTPCvsPtafterPIDelectron(NULL), //+++++++
+ fHistBetaVsdEdXafterPID(NULL),
+ fHistNSigmaTPCTOFvsPtafterPID(NULL),
+ fHistNSigmaTPCTOFPafterPID(NULL),
+ fHistdEdxVsPTPCbeforePIDelectron(NULL),
+ fHistNSigmaTPCvsPtbeforePIDelectron(NULL),
+ fHistdEdxVsPTPCafterPIDelectron(NULL),
+ fHistNSigmaTPCvsPtafterPIDelectron(NULL),
fCentralityArrayBinsForCorrections(kCENTRALITY),
fCentralityWeights(0x0),
fPIDResponse(0x0),
fnAODtrackCutBit(128),
fPtMin(0.3),
fPtMax(1.5),
- fPtMinForCorrections(0.3),//=================================correction
- fPtMaxForCorrections(1.5),//=================================correction
- fPtBinForCorrections(36), //=================================correction
fEtaMin(-0.8),
- fEtaMax(-0.8),
- fEtaMinForCorrections(-0.8),//=================================correction
- fEtaMaxForCorrections(0.8),//=================================correction
- fEtaBinForCorrections(16), //=================================correction
+ fEtaMax(0.8),
fPhiMin(0.),
fPhiMax(360.),
- fPhiMinForCorrections(0.),//=================================correction
- fPhiMaxForCorrections(360.),//=================================correction
- fPhiBinForCorrections(100), //=================================correction
fDCAxyCut(-1),
fDCAzCut(-1),
fTPCchi2Cut(-1),
fNClustersTPCCut(-1),
+ fTPCsharedCut(-1),
fAcceptanceParameterization(0),
fDifferentialV2(0),
fUseFlowAfterBurner(kFALSE),
fList->Add(fHistRapidity);
fHistPhi = new TH2F("fHistPhi","#phi distribution;#phi (rad);Centrality percentile",200,0.0,2.*TMath::Pi(),220,-5,105);
fList->Add(fHistPhi);
- fHistEtaPhiPos = new TH3F("fHistEtaPhiPos","#eta-#phi distribution (+);#eta;#phi (rad);Centrality percentile",80,-2,2,72,0.0,2.*TMath::Pi(),220,-5,105);
+ fHistEtaPhiPos = new TH3F("fHistEtaPhiPos","#eta-#phi distribution (+);#eta;#phi (rad);Centrality percentile",40,-1.6,1.6,72,0.,2.*TMath::Pi(),220,-5,105);
fList->Add(fHistEtaPhiPos);
- fHistEtaPhiNeg = new TH3F("fHistEtaPhiNeg","#eta-#phi distribution (-);#eta;#phi (rad);Centrality percentile",80,-2,2,72,0.0,2.*TMath::Pi(),220,-5,105);
+ fHistEtaPhiNeg = new TH3F("fHistEtaPhiNeg","#eta-#phi distribution (-);#eta;#phi (rad);Centrality percentile",40,-1.6,1.6,72,0.,2.*TMath::Pi(),220,-5,105);
fList->Add(fHistEtaPhiNeg);
fHistPhiBefore = new TH2F("fHistPhiBefore","#phi distribution;#phi;Centrality percentile",200,0.,2*TMath::Pi(),220,-5,105);
fList->Add(fHistPhiBefore);
Int_t nPsiBins;
psibins = fBalance->GetBinning(fBalance->GetBinningString(), "eventPlane", nPsiBins);
-
+
// run the event mixing also in bins of event plane (statistics!)
if(fRunMixingEventPlane){
if(fEventClass=="Multiplicity"){
}
}
}
+
+ if(centbins) delete [] centbins;
+ if(multbins) delete [] multbins;
+ if(vtxbins) delete [] vtxbins;
+ if(psibins) delete [] psibins;
+ // set minimum values for track depth, fraction, and number of events
+ fPoolMgr->SetTargetValues(fMixingTracks, 0.1, 5);
+
// check pool manager
if(!fPoolMgr){
AliError("Event Mixing required, but Pool Manager not initialized...");
fHistProbTPCTOFvsPtbeforePID =new TH2D ("ProbTPCTOFvsPtbefore","ProbTPCTOFvsPtbefore", 1000, -50, 50, 1000, 0, 2.0);
fHistListPIDQA->Add(fHistProbTPCTOFvsPtbeforePID);
- fHistNSigmaTPCvsPtbeforePID = new TH2D ("NSigmaTPCvsPtbefore","NSigmaTPCvsPtbefore", 1000, -10, 10, 1000, 0, 500);
+ fHistNSigmaTPCvsPtbeforePID = new TH2D ("NSigmaTPCvsPtbefore","NSigmaTPCvsPtbefore", 1000, -10, 10, 1000, -25, 25);
fHistListPIDQA->Add(fHistNSigmaTPCvsPtbeforePID);
- fHistNSigmaTOFvsPtbeforePID = new TH2D ("NSigmaTOFvsPtbefore","NSigmaTOFvsPtbefore", 1000, -10, 10, 1000, 0, 500);
+ fHistNSigmaTOFvsPtbeforePID = new TH2D ("NSigmaTOFvsPtbefore","NSigmaTOFvsPtbefore", 1000, -10, 10, 1000, -25, 25);
fHistListPIDQA->Add(fHistNSigmaTOFvsPtbeforePID);
fHistBetaVsdEdXbeforePID = new TH2D ("BetaVsdEdXbefore","BetaVsdEdXbefore", 1000, 0., 1000, 1000, 0, 1.2);
- fHistListPIDQA->Add(fHistBetaVsdEdXbeforePID); //++++++++
+ fHistListPIDQA->Add(fHistBetaVsdEdXbeforePID);
- fHistNSigmaTPCTOFvsPtbeforePID = new TH2D ("NSigmaTPCTOFvsPtbefore","NSigmaTPCTOFvsPtbefore", 1000, -10., 10., 1000, 0, 500);
- fHistListPIDQA->Add(fHistNSigmaTPCTOFvsPtbeforePID); //++++++++
+ fHistNSigmaTPCTOFvsPtbeforePID = new TH2D ("NSigmaTPCTOFvsPtbefore","NSigmaTPCTOFvsPtbefore", 1000, -10., 10., 1000, -25, 25);
+ fHistListPIDQA->Add(fHistNSigmaTPCTOFvsPtbeforePID);
+ //+++++++++++++++++//
+ //p array
+ const Int_t pBins = 36;
+ Double_t nArrayP[pBins+1]={0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0, 1.1, 1.2, 1.3, 1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2.0, 2.2, 2.4, 2.6, 2.8, 3.0, 3.5, 4.0, 4.5, 5.0, 5.5, 6.0, 7.0, 8.0, 9.0, 10.0, 15.0, 20.0};
+ //nSigma Array
+ const Int_t nSigmaBins = 250;
+ Double_t nArrayS[nSigmaBins+1];
+ for (Int_t i = 0; i <= nSigmaBins; i++){
+ nArrayS[i]=i-125; //i+1
+ //Printf("nS: %lf - i: %d", nSigmaArray[i], i);
+ }
+
+ fHistNSigmaTPCTOFPbefPID = new TH3D ("fHistNSigmaTPCTOFPbefPID","fHistNSigmaTPCTOFPbefPID;#sigma_{TPC};#sigma_{TOF};p_{T} (GeV/c)", nSigmaBins, nArrayS, nSigmaBins, nArrayS, pBins,nArrayP);
+ fHistListPIDQA->Add(fHistNSigmaTPCTOFPbefPID);
+ //++++++++++++++//
+
fHistdEdxVsPTPCafterPID = new TH2D ("dEdxVsPTPCafter","dEdxVsPTPCafter", 1000, -10, 10, 1000, 0, 1000);
fHistListPIDQA->Add(fHistdEdxVsPTPCafterPID);
fHistProbTPCTOFvsPtafterPID =new TH2D ("ProbTPCTOFvsPtafter","ProbTPCTOFvsPtafter", 1000, -50, 50, 1000, 0, 2.0);
fHistListPIDQA->Add(fHistProbTPCTOFvsPtafterPID);
- fHistNSigmaTPCvsPtafterPID = new TH2D ("NSigmaTPCvsPtafter","NSigmaTPCvsPtafter", 1000, -10, 10, 1000, 0, 500);
+ fHistNSigmaTPCvsPtafterPID = new TH2D ("NSigmaTPCvsPtafter","NSigmaTPCvsPtafter", 1000, -10, 10, 1000, -25, 25);
fHistListPIDQA->Add(fHistNSigmaTPCvsPtafterPID);
- fHistNSigmaTOFvsPtafterPID = new TH2D ("NSigmaTOFvsPtafter","NSigmaTOFvsPtafter", 1000, -10, 10, 1000, 0, 500);
+ fHistNSigmaTOFvsPtafterPID = new TH2D ("NSigmaTOFvsPtafter","NSigmaTOFvsPtafter", 1000, -10, 10, 1000, -25, 25);
fHistListPIDQA->Add(fHistNSigmaTOFvsPtafterPID);
fHistBetaVsdEdXafterPID = new TH2D ("BetaVsdEdXafter","BetaVsdEdXafter", 1000, 0., 1000, 1000, 0, 1.2);
- fHistListPIDQA->Add(fHistBetaVsdEdXafterPID); //++++++++
+ fHistListPIDQA->Add(fHistBetaVsdEdXafterPID);
+
+ fHistNSigmaTPCTOFvsPtafterPID = new TH2D ("NSigmaTPCTOFvsPtafter","NSigmaTPCTOFvsPtafter", 1000, -10., 10., 1000, -25, 25);
+ fHistListPIDQA->Add(fHistNSigmaTPCTOFvsPtafterPID);
- fHistNSigmaTPCTOFvsPtafterPID = new TH2D ("NSigmaTPCTOFvsPtafter","NSigmaTPCTOFvsPtafter", 1000, -10., 10., 1000, 0, 500);
- fHistListPIDQA->Add(fHistNSigmaTPCTOFvsPtafterPID); //++++++++
+ fHistNSigmaTPCTOFPafterPID = new TH3D ("fHistNSigmaTPCTOFPafterPID","fHistNSigmaTPCTOFPafterPID;#sigma_{TPC};#sigma_{TOF};p_{T} (GeV/c)", nSigmaBins, nArrayS, nSigmaBins, nArrayS, pBins,nArrayP);
+ fHistListPIDQA->Add(fHistNSigmaTPCTOFPafterPID); //++++++++++++++
}
// for electron rejection only TPC nsigma histograms
AliInfo("Finished setting up the Output");
TH1::AddDirectory(oldStatus);
+
}
return;
}
}//loop over centralities: ONLY the PbPb case is covered
-
- if(fHistCorrectionPlus[0]){
- fEtaMinForCorrections = fHistCorrectionPlus[0]->GetXaxis()->GetXmin();
- fEtaMaxForCorrections = fHistCorrectionPlus[0]->GetXaxis()->GetXmax();
- fEtaBinForCorrections = fHistCorrectionPlus[0]->GetNbinsX();
-
- fPtMinForCorrections = fHistCorrectionPlus[0]->GetYaxis()->GetXmin();
- fPtMaxForCorrections = fHistCorrectionPlus[0]->GetYaxis()->GetXmax();
- fPtBinForCorrections = fHistCorrectionPlus[0]->GetNbinsY();
-
- fPhiMinForCorrections = fHistCorrectionPlus[0]->GetZaxis()->GetXmin();
- fPhiMaxForCorrections = fHistCorrectionPlus[0]->GetZaxis()->GetXmax();
- fPhiBinForCorrections = fHistCorrectionPlus[0]->GetNbinsZ();
- }
}
//________________________________________________________________________
return;
}
// get the reaction plane
- if(fEventClass != "Multiplicity") {
+ if(fEventClass != "Multiplicity" && gAnalysisLevel!="AODnano") {
gReactionPlane = GetEventPlane(eventMain);
fHistEventPlane->Fill(gReactionPlane,lMultiplicityVar);
if(gReactionPlane < 0){
//pool->SetDebug(1);
- if (pool->IsReady() || pool->NTracksInPool() > fMixingTracks / 10 || pool->GetCurrentNEvents() >= 5){
+ if (pool->IsReady()){
Int_t nMix = pool->GetCurrentNEvents();
Double_t AliAnalysisTaskBFPsi::GetRefMultiOrCentrality(AliVEvent *event){
// Checks the Event cuts
// Fills Event statistics histograms
-
+
Float_t gCentrality = -1.;
Double_t gMultiplicity = -1.;
TString gAnalysisLevel = fBalance->GetAnalysisLevel();
}//AOD header
}//AOD
+
+ // calculate centrality always (not only in centrality mode)
+ else if(gAnalysisLevel == "AODnano" ) { //centrality via JF workaround
+
+ AliAODHeader *header = (AliAODHeader*) event->GetHeader();
+ if(header){
+ gCentrality = (Float_t) gROOT->ProcessLine(Form("100.0 + 100.0 * ((AliNanoAODHeader*) %p)->GetCentrality(\"%s\")", header,fCentralityEstimator.Data())) / 100 - 1.0;
+
+ // QA histogram
+ fHistCentStatsUsed->Fill(0.,gCentrality);
+
+ }//AOD header
+ }//AODnano
else if(gAnalysisLevel == "ESD" || gAnalysisLevel == "MCESD"){ // centrality class for ESDs or MC-ESDs
AliCentrality *centrality = event->GetCentrality();
// -- Get efficiency correction of particle dependent on (eta, phi, pt, charge, centrality)
Double_t correction = 1.;
- Int_t binEta = 0, binPt = 0, binPhi = 0;
-
- //Printf("EtaMAx: %lf - EtaMin: %lf - EtaBin: %lf", fEtaMaxForCorrections,fEtaMinForCorrections,fEtaBinForCorrections);
- if(fEtaBinForCorrections != 0) {
- Double_t widthEta = (fEtaMaxForCorrections - fEtaMinForCorrections)/fEtaBinForCorrections;
- if(fEtaMaxForCorrections != fEtaMinForCorrections)
- binEta = (Int_t)((vEta-fEtaMinForCorrections)/widthEta)+1;
- }
-
- if(fPtBinForCorrections != 0) {
- Double_t widthPt = (fPtMaxForCorrections - fPtMinForCorrections)/fPtBinForCorrections;
- if(fPtMaxForCorrections != fPtMinForCorrections)
- binPt = (Int_t)((vPt-fPtMinForCorrections)/widthPt) + 1;
- }
-
- if(fPhiBinForCorrections != 0) {
- Double_t widthPhi = (fPhiMaxForCorrections - fPhiMinForCorrections)/fPhiBinForCorrections;
- if(fPhiMaxForCorrections != fPhiMinForCorrections)
- binPhi = (Int_t)((vPhi-fPhiMinForCorrections)/widthPhi)+ 1;
- }
-
Int_t gCentralityInt = -1;
+
for (Int_t i=0; i<fCentralityArrayBinsForCorrections-1; i++){
if((fCentralityArrayForCorrections[i] <= gCentrality)&&(gCentrality <= fCentralityArrayForCorrections[i+1])){
gCentralityInt = i;
if(fHistCorrectionPlus[gCentralityInt]){
if (vCharge > 0) {
- correction = fHistCorrectionPlus[gCentralityInt]->GetBinContent(fHistCorrectionPlus[gCentralityInt]->GetBin(binEta, binPt, binPhi));
+ correction = fHistCorrectionPlus[gCentralityInt]->GetBinContent(fHistCorrectionPlus[gCentralityInt]->FindBin(vEta,vPt,vPhi));
//Printf("CORRECTIONplus: %.2f | Centrality %d",correction,gCentralityInt);
}
if (vCharge < 0) {
- correction = fHistCorrectionMinus[gCentralityInt]->GetBinContent(fHistCorrectionMinus[gCentralityInt]->GetBin(binEta, binPt, binPhi));
- //Printf("CORRECTIONminus: %.2f | Centrality %d",correction,gCentralityInt);
+ correction = fHistCorrectionPlus[gCentralityInt]->GetBinContent(fHistCorrectionMinus[gCentralityInt]->FindBin(vEta,vPt,vPhi));
+ //Printf("CORRECTIONminus: %.2f | Centrality %d",correction,gCentralityInt);
}
}
else {
Double_t probTPCTOF[AliPID::kSPECIES]={0.};
Double_t nSigma = 0.;
- Double_t nSigmaTPC = 0.; //++++
- Double_t nSigmaTOF = 0.; //++++
- Double_t nSigmaTPCTOF = 0.; //++++
+ Double_t nSigmaTPC = 0.;
+ Double_t nSigmaTOF = 0.;
+ Double_t nSigmaTPCTOF = 0.;
UInt_t detUsedTPC = 0;
UInt_t detUsedTOF = 0;
UInt_t detUsedTPCTOF = 0;
- nSigmaTPC = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(aodTrack,(AliPID::EParticleType)fParticleOfInterest));
- nSigmaTOF = TMath::Abs(fPIDResponse->NumberOfSigmasTOF(aodTrack,(AliPID::EParticleType)fParticleOfInterest));
+ nSigmaTPC = fPIDResponse->NumberOfSigmasTPC(aodTrack,(AliPID::EParticleType)fParticleOfInterest);
+ nSigmaTOF = fPIDResponse->NumberOfSigmasTOF(aodTrack,(AliPID::EParticleType)fParticleOfInterest);
nSigmaTPCTOF = TMath::Sqrt(nSigmaTPC*nSigmaTPC + nSigmaTOF*nSigmaTOF);
+ if (nSigmaTOF == 999 || nSigmaTOF == -999){
+ nSigmaTPCTOF = nSigmaTPC;
+ }
//Decide what detector configuration we want to use
switch(fPidDetectorConfig) {
Double_t tofTime = -999., length = 999., tof = -999.;
Double_t c = TMath::C()*1.E-9;// m/ns
Double_t beta = -999.;
- if ( (aodTrack->IsOn(AliAODTrack::kTOFin)) &&
- (aodTrack->IsOn(AliAODTrack::kTIME)) ) {
- tofTime = aodTrack->GetTOFsignal();//in ps
- length = aodTrack->GetIntegratedLength();
- tof = tofTime*1E-3; // ns
-
- if (tof <= 0) {
- //Printf("WARNING: track with negative TOF time found! Skipping this track for PID checks\n");
- continue;
- }
- if (length <= 0){
- //printf("WARNING: track with negative length found!Skipping this track for PID checks\n");
- continue;
- }
-
- length = length*0.01; // in meters
- tof = tof*c;
- beta = length/tof;
-
- fHistBetavsPTOFbeforePID ->Fill(aodTrack->P()*aodTrack->Charge(),beta);
- fHistProbTOFvsPtbeforePID ->Fill(aodTrack->Pt(),probTOF[fParticleOfInterest]);
- fHistNSigmaTOFvsPtbeforePID ->Fill(aodTrack->Pt(),nSigmaTOF);
- }//TOF signal
-
- fHistdEdxVsPTPCbeforePID -> Fill(aodTrack->P()*aodTrack->Charge(),aodTrack->GetTPCsignal());
- fHistProbTPCvsPtbeforePID -> Fill(aodTrack->Pt(),probTPC[fParticleOfInterest]);
- fHistNSigmaTPCvsPtbeforePID -> Fill(aodTrack->Pt(),nSigmaTPC);
-
- fHistProbTPCTOFvsPtbeforePID -> Fill(aodTrack->Pt(),probTPCTOF[fParticleOfInterest]);
- //combined TPC&TOF
- fHistBetaVsdEdXbeforePID->Fill(aodTrack->GetTPCsignal(),beta); //+++++++++
- fHistNSigmaTPCTOFvsPtbeforePID -> Fill(aodTrack->Pt(),nSigmaTPCTOF);
- Printf("NSIGMA: %lf - nSigmaTOF: %lf - nSigmaTPC: %lf - nSigmaTPCTOF: %lf",nSigma,nSigmaTOF,nSigmaTPC,nSigmaTPCTOF);
-
- //end of QA-before pid
-
- if ((detUsedTPC != 0)||(detUsedTOF != 0)||(detUsedTPCTOF != 0)) {
- //Make the decision based on the n-sigma
- if(fUsePIDnSigma) {
- if(nSigma > fPIDNSigma) continue;
+ Float_t probMis = fPIDResponse->GetTOFMismatchProbability(aodTrack);
+ if (probMis < 0.01) { //if u want to reduce mismatch using also TPC
+
+ if ((aodTrack->IsOn(AliAODTrack::kITSin)) && (aodTrack->IsOn(AliAODTrack::kTOFpid)) ) { //leonardo's analysis
- fHistNSigmaTOFvsPtafterPID ->Fill(aodTrack->Pt(),nSigmaTOF);
- fHistNSigmaTPCvsPtafterPID ->Fill(aodTrack->Pt(),nSigmaTPC);
- fHistNSigmaTPCTOFvsPtafterPID ->Fill(aodTrack->Pt(),nSigmaTPCTOF);
- }
- //Make the decision based on the bayesian
- else if(fUsePIDPropabilities) {
- if(fParticleOfInterest != TMath::LocMax(AliPID::kSPECIES,prob)) continue;
- if (prob[fParticleOfInterest] < fMinAcceptedPIDProbability) continue;
+ tofTime = aodTrack->GetTOFsignal();//in ps
+ length = aodTrack->GetIntegratedLength();
+ tof = tofTime*1E-3; // ns
+ if (tof <= 0) {
+ //Printf("WARNING: track with negative TOF time found! Skipping this track for PID checks\n");
+ continue;
+ }
+ if (length <= 0){
+ //printf("WARNING: track with negative length found!Skipping this track for PID checks\n");
+ continue;
+ }
+ length = length*0.01; // in meters
+ tof = tof*c;
+ beta = length/tof;
+
+ fHistBetavsPTOFbeforePID ->Fill(aodTrack->P()*aodTrack->Charge(),beta);
+ fHistProbTOFvsPtbeforePID ->Fill(aodTrack->Pt(),probTOF[fParticleOfInterest]);
+ fHistNSigmaTOFvsPtbeforePID ->Fill(aodTrack->Pt(),nSigmaTOF);
+ }//TOF signal
- fHistProbTOFvsPtafterPID ->Fill(aodTrack->Pt(),probTOF[fParticleOfInterest]);
- fHistProbTPCvsPtafterPID ->Fill(aodTrack->Pt(),probTPC[fParticleOfInterest]);
- fHistProbTPCTOFvsPtafterPID ->Fill(aodTrack->Pt(),probTPCTOF[fParticleOfInterest]);
+ fHistdEdxVsPTPCbeforePID -> Fill(aodTrack->GetTPCmomentum()*aodTrack->Charge(),aodTrack->GetTPCsignal()); //aodTrack->P()*aodTrack->Charge()
+ fHistProbTPCvsPtbeforePID -> Fill(aodTrack->Pt(),probTPC[fParticleOfInterest]);
+ fHistNSigmaTPCvsPtbeforePID -> Fill(aodTrack->Pt(),nSigmaTPC);
+ fHistProbTPCTOFvsPtbeforePID -> Fill(aodTrack->Pt(),probTPCTOF[fParticleOfInterest]);
- }
+ //combined TPC&TOF
+ fHistBetaVsdEdXbeforePID->Fill(aodTrack->GetTPCsignal(),beta);
+ fHistNSigmaTPCTOFvsPtbeforePID -> Fill(aodTrack->Pt(),nSigmaTPCTOF);
+ fHistNSigmaTPCTOFPbefPID ->Fill(nSigmaTPC,nSigmaTOF,aodTrack->P()); //++++++++++++++
+ //Printf("NSIGMA: %lf - nSigmaTOF: %lf - nSigmaTPC: %lf - nSigmaTPCTOF: %lf",nSigma,nSigmaTOF,nSigmaTPC,nSigmaTPCTOF);
+ //end of QA-before pid
+
+ if ((detUsedTPC != 0)||(detUsedTOF != 0)||(detUsedTPCTOF != 0)) {
+ //Make the decision based on the n-sigma
+ if(fUsePIDnSigma) {
+ if(nSigma > fPIDNSigma || nSigma < -fPIDNSigma) continue;
+
+ fHistNSigmaTOFvsPtafterPID ->Fill(aodTrack->Pt(),nSigmaTOF);
+ fHistNSigmaTPCvsPtafterPID ->Fill(aodTrack->Pt(),nSigmaTPC);
+ fHistNSigmaTPCTOFvsPtafterPID ->Fill(aodTrack->Pt(),nSigmaTPCTOF);
+ fHistNSigmaTPCTOFPafterPID ->Fill(nSigmaTPC,nSigmaTOF,aodTrack->P()); //++++++++++++++
+ }
+ //Make the decision based on the bayesian
+ else if(fUsePIDPropabilities) {
+ if(fParticleOfInterest != TMath::LocMax(AliPID::kSPECIES,prob)) continue;
+ if (prob[fParticleOfInterest] < fMinAcceptedPIDProbability) continue;
+
+ fHistProbTOFvsPtafterPID ->Fill(aodTrack->Pt(),probTOF[fParticleOfInterest]);
+ fHistProbTPCvsPtafterPID ->Fill(aodTrack->Pt(),probTPC[fParticleOfInterest]);
+ fHistProbTPCTOFvsPtafterPID ->Fill(aodTrack->Pt(),probTPCTOF[fParticleOfInterest]);
+ }
- //Fill QA after the PID
- fHistBetavsPTOFafterPID ->Fill(aodTrack->P()*aodTrack->Charge(),beta);
- fHistdEdxVsPTPCafterPID ->Fill(aodTrack->P()*aodTrack->Charge(),aodTrack->GetTPCsignal());
- fHistBetaVsdEdXafterPID->Fill(aodTrack->GetTPCsignal(),beta); //+++++++++
+ //Fill QA after the PID
+ fHistBetavsPTOFafterPID ->Fill(aodTrack->P()*aodTrack->Charge(),beta);
+ fHistdEdxVsPTPCafterPID ->Fill(aodTrack->P()*aodTrack->Charge(),aodTrack->GetTPCsignal());
+ fHistBetaVsdEdXafterPID ->Fill(aodTrack->GetTPCsignal(),beta);
+ }
}
}
//===========================PID===============================//
if( fNClustersTPCCut != -1 && aodTrack->GetTPCNcls() < fNClustersTPCCut){
continue;
}
+
+ // Extra cut on shared clusters
+ if( fTPCsharedCut != -1 && aodTrack->GetTPCnclsS() > fTPCsharedCut){
+ continue;
+ }
// fill QA histograms
fHistClus->Fill(aodTrack->GetITSNcls(),aodTrack->GetTPCNcls());
}//track loop
}// AOD analysis
+
+ // nano AODs
+ else if(gAnalysisLevel == "AODnano") { // not fully supported yet (PID missing)
+ // Loop over tracks in event
+
+ for (Int_t iTracks = 0; iTracks < event->GetNumberOfTracks(); iTracks++) {
+ AliVTrack* aodTrack = dynamic_cast<AliVTrack *>(event->GetTrack(iTracks));
+ if (!aodTrack) {
+ AliError(Form("Could not receive track %d", iTracks));
+ continue;
+ }
+
+ // AOD track cuts (not needed)
+ //if(!aodTrack->TestFilterBit(fnAODtrackCutBit)) continue;
+
+ vCharge = aodTrack->Charge();
+ vEta = aodTrack->Eta();
+ vY = -999.;
+ vPhi = aodTrack->Phi();// * TMath::RadToDeg();
+ vPt = aodTrack->Pt();
+
+
+ // Kinematics cuts from ESD track cuts
+ if( vPt < fPtMin || vPt > fPtMax) continue;
+ if( vEta < fEtaMin || vEta > fEtaMax) continue;
+
+
+ // fill QA histograms
+ fHistPt->Fill(vPt,gCentrality);
+ fHistEta->Fill(vEta,gCentrality);
+ fHistRapidity->Fill(vY,gCentrality);
+ if(vCharge > 0) fHistPhiPos->Fill(vPhi,gCentrality);
+ else if(vCharge < 0) fHistPhiNeg->Fill(vPhi,gCentrality);
+ fHistPhi->Fill(vPhi,gCentrality);
+ if(vCharge > 0) fHistEtaPhiPos->Fill(vEta,vPhi,gCentrality);
+ else if(vCharge < 0) fHistEtaPhiNeg->Fill(vEta,vPhi,gCentrality);
+
+ //=======================================correction
+ Double_t correction = GetTrackbyTrackCorrectionMatrix(vEta, vPhi, vPt, vCharge, gCentrality);
+ //Printf("CORRECTIONminus: %.2f | Centrality %lf",correction,gCentrality);
+
+ // add the track to the TObjArray
+ tracksAccepted->Add(new AliBFBasicParticle(vEta, vPhi, vPt, vCharge, correction));
+ }//track loop
+ }// AOD nano analysis
+
+
//==============================================================================================================
else if(gAnalysisLevel == "MCAOD") {
continue;
}
+ // Extra cut on shared clusters
+ if( fTPCsharedCut != -1 && aodTrack->GetTPCnclsS() > fTPCsharedCut){
+ continue;
+ }
+
//Exclude resonances
if(fExcludeResonancesInMC) {
Bool_t kExcludeParticle = kFALSE;
Int_t gMotherIndex = particle->GetFirstMother();
+
if(gMotherIndex != -1) {
- AliMCParticle* motherTrack = dynamic_cast<AliMCParticle *>(event->GetTrack(gMotherIndex));
- if(motherTrack) {
- TParticle *motherParticle = motherTrack->Particle();
- if(motherParticle) {
- Int_t pdgCodeOfMother = motherParticle->GetPdgCode();
- //if((pdgCodeOfMother == 113)||(pdgCodeOfMother == 213)||(pdgCodeOfMother == 221)||(pdgCodeOfMother == 223)||(pdgCodeOfMother == 331)||(pdgCodeOfMother == 333)) {
- if(pdgCodeOfMother == 113 // rho0
- || pdgCodeOfMother == 213 || pdgCodeOfMother == -213 // rho+
- // || pdgCodeOfMother == 221 // eta
- // || pdgCodeOfMother == 331 // eta'
- // || pdgCodeOfMother == 223 // omega
- // || pdgCodeOfMother == 333 // phi
- || pdgCodeOfMother == 311 || pdgCodeOfMother == -311 // K0
- // || pdgCodeOfMother == 313 || pdgCodeOfMother == -313 // K0*
- // || pdgCodeOfMother == 323 || pdgCodeOfMother == -323 // K+*
- || pdgCodeOfMother == 3122 || pdgCodeOfMother == -3122 // Lambda
- || pdgCodeOfMother == 111 // pi0 Dalitz
- ) {
- kExcludeParticle = kTRUE;
- }
- }
- }
+ // AliMCParticle* motherTrack = dynamic_cast<AliMCParticle *>(event->GetTrack(gMotherIndex));
+ // if(motherTrack) {
+ // TParticle *motherParticle = motherTrack->Particle();
+ // if(motherParticle) {
+ // Int_t pdgCodeOfMother = motherParticle->GetPdgCode();
+ // //if((pdgCodeOfMother == 113)||(pdgCodeOfMother == 213)||(pdgCodeOfMother == 221)||(pdgCodeOfMother == 223)||(pdgCodeOfMother == 331)||(pdgCodeOfMother == 333)) {
+ // if(pdgCodeOfMother == 113 // rho0
+ // || pdgCodeOfMother == 213 || pdgCodeOfMother == -213 // rho+
+ // // || pdgCodeOfMother == 221 // eta
+ // // || pdgCodeOfMother == 331 // eta'
+ // // || pdgCodeOfMother == 223 // omega
+ // // || pdgCodeOfMother == 333 // phi
+ // || pdgCodeOfMother == 311 || pdgCodeOfMother == -311 // K0
+ // // || pdgCodeOfMother == 313 || pdgCodeOfMother == -313 // K0*
+ // // || pdgCodeOfMother == 323 || pdgCodeOfMother == -323 // K+*
+ // || pdgCodeOfMother == 3122 || pdgCodeOfMother == -3122 // Lambda
+ // || pdgCodeOfMother == 111 // pi0 Dalitz
+ // ) {
+ kExcludeParticle = kTRUE;
+ // }
+ // }
+ //}
}
//Exclude from the analysis decay products of rho0, rho+, eta, eta' and phi
if(kExcludeParticle) continue;
+
}
//Exclude electrons with PDG