From a0b94e5c61c0e9d68eba33e7d01dbd367d345940 Mon Sep 17 00:00:00 2001 From: kaamodt Date: Tue, 15 Sep 2009 05:56:16 +0000 Subject: [PATCH] Added correction framework (Kathrin) Added new functions to calculate the xyz point of conversion (Markus) Added CORRFW to the makefile Changed the configuration macro for the new functionality. --- .../AliAnalysisTaskGammaConversion.cxx | 682 ++++---- .../AliAnalysisTaskGammaConversion.h | 106 +- PWG4/GammaConv/AliV0Reader.cxx | 1533 ++++++++++------- PWG4/GammaConv/AliV0Reader.h | 1092 ++++++------ PWG4/Makefile | 4 + PWG4/macros/ConfigGammaConversion.C | 304 ++-- 6 files changed, 2075 insertions(+), 1646 deletions(-) diff --git a/PWG4/GammaConv/AliAnalysisTaskGammaConversion.cxx b/PWG4/GammaConv/AliAnalysisTaskGammaConversion.cxx index 7d720968e56..66c8839c40c 100644 --- a/PWG4/GammaConv/AliAnalysisTaskGammaConversion.cxx +++ b/PWG4/GammaConv/AliAnalysisTaskGammaConversion.cxx @@ -28,6 +28,8 @@ #include "AliLog.h" #include "AliESDtrackCuts.h" #include "TNtuple.h" +#include "AliCFManager.h" // for CF +#include "AliCFContainer.h" // for CF class AliKFVertex; class AliAODHandler; @@ -50,8 +52,11 @@ AliAnalysisTaskGammaConversion::AliAnalysisTaskGammaConversion(): AliAnalysisTaskSE(), fV0Reader(NULL), fStack(NULL), + fMCTruth(NULL), // for CF + fMCEvent(NULL), // for CF fESDEvent(NULL), fOutputContainer(NULL), + fCFManager(0x0), // for CF fHistograms(NULL), fDoMCTruth(kFALSE), fDoNeutralMeson(kFALSE), @@ -96,7 +101,7 @@ AliAnalysisTaskSE(), fLeadingChargedIndex(-1), fAODBranch(NULL), fAODBranchName("GammaConv")//, -// fAODObjects(NULL) + // fAODObjects(NULL) { // Default constructor // Common I/O in slot 0 @@ -105,7 +110,7 @@ AliAnalysisTaskSE(), // Your private output DefineOutput(1, TList::Class()); - + // Define standard ESD track cuts for Gamma-hadron correlation SetESDtrackCuts(); } @@ -114,8 +119,11 @@ AliAnalysisTaskGammaConversion::AliAnalysisTaskGammaConversion(const char* name) AliAnalysisTaskSE(name), fV0Reader(NULL), fStack(NULL), + fMCTruth(NULL), // for CF + fMCEvent(NULL), // for CF fESDEvent(NULL), fOutputContainer(0x0), + fCFManager(0x0), // for CF fHistograms(NULL), fDoMCTruth(kFALSE), fDoNeutralMeson(kFALSE), @@ -168,7 +176,9 @@ AliAnalysisTaskGammaConversion::AliAnalysisTaskGammaConversion(const char* name) // Your private output DefineOutput(1, TList::Class()); - + DefineOutput(2, AliCFContainer::Class()); // for CF + + // Define standard ESD track cuts for Gamma-hadron correlation SetESDtrackCuts(); } @@ -187,6 +197,12 @@ AliAnalysisTaskGammaConversion::~AliAnalysisTaskGammaConversion() if(fV0Reader){ delete fV0Reader; } + + // for CF + if(fCFManager){ + delete fCFManager; + } + if (fAODBranch) { fAODBranch->Clear(); delete fAODBranch ; @@ -202,19 +218,19 @@ void AliAnalysisTaskGammaConversion::Init() void AliAnalysisTaskGammaConversion::SetESDtrackCuts() { // SetESDtrackCuts - + fEsdTrackCuts = new AliESDtrackCuts("AliESDtrackCuts"); -//standard cuts from: -//http://aliceinfo.cern.ch/alicvs/viewvc/PWG0/dNdEta/CreateCuts.C?revision=1.4&view=markup -//fEsdTrackCuts->SetMinNClustersTPC(50); -//fEsdTrackCuts->SetMaxChi2PerClusterTPC(3.5); -//fEsdTrackCuts->SetMaxCovDiagonalElements(2,2,0.5,0.5,2); + //standard cuts from: + //http://aliceinfo.cern.ch/alicvs/viewvc/PWG0/dNdEta/CreateCuts.C?revision=1.4&view=markup + //fEsdTrackCuts->SetMinNClustersTPC(50); + //fEsdTrackCuts->SetMaxChi2PerClusterTPC(3.5); + //fEsdTrackCuts->SetMaxCovDiagonalElements(2,2,0.5,0.5,2); fEsdTrackCuts->SetRequireTPCRefit(kTRUE); fEsdTrackCuts->SetRequireITSRefit(kTRUE); fEsdTrackCuts->SetMaxNsigmaToVertex(3); fEsdTrackCuts->SetRequireSigmaToVertex(kTRUE); // fEsdTrackCuts->SetAcceptKinkDaughters(kFALSE); - + } void AliAnalysisTaskGammaConversion::Exec(Option_t */*option*/) @@ -225,7 +241,7 @@ void AliAnalysisTaskGammaConversion::Exec(Option_t */*option*/) //Each event needs an empty branch fAODBranch->Clear(); - + if(fKFReconstructedGammasTClone == NULL){ fKFReconstructedGammasTClone = new TClonesArray("AliKFParticle",0); } @@ -247,7 +263,7 @@ void AliAnalysisTaskGammaConversion::Exec(Option_t */*option*/) if(fChargedParticles == NULL){ fChargedParticles = new TClonesArray("AliESDtrack",0); } - + //clear TClones fKFReconstructedGammasTClone->Clear(); fCurrentEventPosElectronTClone->Clear(); @@ -255,7 +271,7 @@ void AliAnalysisTaskGammaConversion::Exec(Option_t */*option*/) fKFReconstructedGammasCutTClone->Clear(); fPreviousEventTLVNegElectronTClone->Clear(); fPreviousEventTLVPosElectronTClone->Clear(); - + //clear vectors // fKFReconstructedGammas.clear(); fElectronv1.clear(); @@ -266,40 +282,40 @@ void AliAnalysisTaskGammaConversion::Exec(Option_t */*option*/) fChargedParticles->Clear(); fChargedParticlesId.clear(); - + //Clear the data in the v0Reader fV0Reader->UpdateEventByEventData(); - - + + // Process the MC information if(fDoMCTruth){ ProcessMCData(); } - + //Process the v0 information with no cuts ProcessV0sNoCut(); - + // Process the v0 information ProcessV0s(); - + //Fill Gamma AOD FillAODWithConversionGammas() ; - + //calculate background if flag is set if(fCalculateBackground){ CalculateBackground(); } - - - - + + + + // Process reconstructed gammas if(fDoNeutralMeson == kTRUE){ - ProcessGammasForNeutralMesonAnalysis(); + ProcessGammasForNeutralMesonAnalysis(); } - + CheckV0Efficiency(); - + //Process reconstructed gammas electrons for Chi_c Analysis if(fDoChic == kTRUE){ ProcessGammaElectronsForChicAnalysis(); @@ -308,8 +324,9 @@ void AliAnalysisTaskGammaConversion::Exec(Option_t */*option*/) if(fDoJet == kTRUE){ ProcessGammasForGammaJetAnalysis(); } - + PostData(1, fOutputContainer); + PostData(2, fCFManager->GetParticleContainer()); // for CF } @@ -328,11 +345,21 @@ void AliAnalysisTaskGammaConversion::ProcessMCData(){ // see header file for documentation fStack = fV0Reader->GetMCStack(); - + fMCTruth = fV0Reader->GetMCTruth(); // for CF + fMCEvent = fV0Reader->GetMCEvent(); // for CF + + + // for CF + if(!fMCEvent) cout << "NO MC INFO FOUND" << endl; + fCFManager->SetEventInfo(fMCEvent); + Double_t containerInput[3]; + // end for CF + + if(fV0Reader->CheckForPrimaryVertex() == kFALSE){ return; // aborts if the primary vertex does not have contributors. } - + for (Int_t iTracks = 0; iTracks < fStack->GetNtrack(); iTracks++) { TParticle* particle = (TParticle *)fStack->Particle(iTracks); @@ -340,12 +367,15 @@ void AliAnalysisTaskGammaConversion::ProcessMCData(){ //print warning here continue; } - + + + + ///////////////////////Begin Chic Analysis///////////////////////////// - - + + if(particle->GetPdgCode() == 443){//Is JPsi - + if(particle->GetNDaughters()==2){ if(TMath::Abs(fStack->Particle(particle->GetFirstDaughter())->GetPdgCode()) == 11 && TMath::Abs(fStack->Particle(particle->GetLastDaughter())->GetPdgCode()) == 11){ @@ -353,13 +383,13 @@ void AliAnalysisTaskGammaConversion::ProcessMCData(){ TParticle* daug1 = fStack->Particle(particle->GetLastDaughter()); if(TMath::Abs(daug0->Eta()) < 0.9 && TMath::Abs(daug1->Eta()) < 0.9) fHistograms->FillTable("Table_Electrons",3);//e+ e- from J/Psi inside acceptance - + if( TMath::Abs(daug0->Eta()) < 0.9){ if(daug0->GetPdgCode() == -11) fHistograms->FillTable("Table_Electrons",1);//e+ from J/Psi inside acceptance else fHistograms->FillTable("Table_Electrons",2);//e- from J/Psi inside acceptance - + } if(TMath::Abs(daug1->Eta()) < 0.9){ if(daug1->GetPdgCode() == -11) @@ -387,14 +417,14 @@ void AliAnalysisTaskGammaConversion::ProcessMCData(){ if( particle->GetNDaughters() == 2){ TParticle* daug0 = fStack->Particle(particle->GetFirstDaughter()); TParticle* daug1 = fStack->Particle(particle->GetLastDaughter()); - + if( (daug0->GetPdgCode() == 443 || daug0->GetPdgCode() == 22) && (daug1->GetPdgCode() == 443 || daug1->GetPdgCode() == 22) ){ if( daug0->GetPdgCode() == 443){ TParticle* daugE0 = fStack->Particle(daug0->GetFirstDaughter()); TParticle* daugE1 = fStack->Particle(daug0->GetLastDaughter()); if( TMath::Abs(daug1->Eta()) < 1.2 && TMath::Abs(daugE0->Eta()) < 0.9 && TMath::Abs(daugE1->Eta()) < 0.9 ) fHistograms->FillTable("Table_Electrons",18); - + }//if else if (daug1->GetPdgCode() == 443){ TParticle* daugE0 = fStack->Particle(daug1->GetFirstDaughter()); @@ -405,13 +435,13 @@ void AliAnalysisTaskGammaConversion::ProcessMCData(){ }//gamma o Jpsi }//GetNDaughters } - - + + /////////////////////End Chic Analysis//////////////////////////// - - + + if(TMath::Abs(particle->Eta())> fV0Reader->GetEtaCut() ) continue; - + if(particle->R()>fV0Reader->GetMaxRCut()) continue; // cuts on distance from collision point Double_t tmpPhi=particle->Phi(); @@ -434,17 +464,22 @@ void AliAnalysisTaskGammaConversion::ProcessMCData(){ if(particle->GetMother(0) >-1 && fStack->Particle(particle->GetMother(0))->GetPdgCode() == 22){ continue; // no photon as mothers! } - + if(particle->GetMother(0) >= fStack->GetNprimary()){ continue; // the gamma has a mother, and it is not a primary particle } - + fHistograms->FillHistogram("MC_allGamma_Energy", particle->Energy()); fHistograms->FillHistogram("MC_allGamma_Pt", particle->Pt()); fHistograms->FillHistogram("MC_allGamma_Eta", particle->Eta()); fHistograms->FillHistogram("MC_allGamma_Phi", tmpPhi); fHistograms->FillHistogram("MC_allGamma_Rapid", rapidity); + // for CF + containerInput[0] = particle->Pt(); + containerInput[1] = particle->Eta(); + containerInput[2] = fStack->Particle(particle->GetMother(0))->GetMass(); + fCFManager->GetParticleContainer()->Fill(containerInput,kStepGenerated); // generated gamma if(particle->GetMother(0) < 0){ // direct gamma fHistograms->FillHistogram("MC_allDirectGamma_Energy",particle->Energy()); @@ -454,7 +489,6 @@ void AliAnalysisTaskGammaConversion::ProcessMCData(){ fHistograms->FillHistogram("MC_allDirectGamma_Rapid", rapidity); } - // looking for conversion (electron + positron from pairbuilding (= 5) ) TParticle* ePos = NULL; TParticle* eNeg = NULL; @@ -489,23 +523,28 @@ void AliAnalysisTaskGammaConversion::ProcessMCData(){ if(ePos->Pt()GetPtCut() || eNeg->Pt()GetPtCut()){ continue; // no reconstruction below the Pt cut } - + if(TMath::Abs(ePos->Eta())> fV0Reader->GetEtaCut() || TMath::Abs(eNeg->Eta())> fV0Reader->GetEtaCut()){ continue; } - + if(ePos->R()>fV0Reader->GetMaxRCut()){ continue; // cuts on distance from collision point } - - if(TMath::Abs(ePos->Vz()) > 240){ + + if(TMath::Abs(ePos->Vz()) > fV0Reader->GetMaxZCut()){ continue; // outside material } + if((TMath::Abs(ePos->Vz()) * fV0Reader->GetLineCutZRSlope()) - fV0Reader->GetLineCutZValue() > ePos->R()){ continue; // line cut to exclude regions where we do not reconstruct } + + // for CF + fCFManager->GetParticleContainer()->Fill(containerInput,kStepReconstructable); // reconstructable gamma + fHistograms->FillHistogram("MC_ConvGamma_Energy", particle->Energy()); fHistograms->FillHistogram("MC_ConvGamma_Pt", particle->Pt()); fHistograms->FillHistogram("MC_ConvGamma_Eta", particle->Eta()); @@ -524,9 +563,6 @@ void AliAnalysisTaskGammaConversion::ProcessMCData(){ fHistograms->FillHistogram("MC_P_Phi", ePosPhi); - - //cout << "filled histos for converted gamma, ePos, eNeg" << endl; - // begin Mapping Int_t rBin = fHistograms->GetRBin(ePos->R()); Int_t phiBin = fHistograms->GetPhiBin(particle->Phi()); @@ -553,9 +589,6 @@ void AliAnalysisTaskGammaConversion::ProcessMCData(){ fHistograms->FillHistogram("MC_Conversion_XY",ePos->Vx(),ePos->Vy()); fHistograms->FillHistogram("MC_Conversion_OpeningAngle",GetMCOpeningAngle(ePos, eNeg)); - //cout << "mapping is done" << endl; - - if(particle->GetMother(0) < 0){ // no mother = direct gamma, still inside converted fHistograms->FillHistogram("MC_ConvDirectGamma_Energy",particle->Energy()); fHistograms->FillHistogram("MC_ConvDirectGamma_Pt", particle->Pt()); @@ -566,11 +599,11 @@ void AliAnalysisTaskGammaConversion::ProcessMCData(){ } // end direct gamma else{ // mother exits /* if( fStack->Particle(particle->GetMother(0))->GetPdgCode()==10441 ||//chic0 - fStack->Particle(particle->GetMother(0))->GetPdgCode()==20443 ||//psi2S - fStack->Particle(particle->GetMother(0))->GetPdgCode()==445 //chic2 - ){ - fMCGammaChic.push_back(particle); - } + fStack->Particle(particle->GetMother(0))->GetPdgCode()==20443 ||//psi2S + fStack->Particle(particle->GetMother(0))->GetPdgCode()==445 //chic2 + ){ + fMCGammaChic.push_back(particle); + } */ } // end if mother exits } // end if particle is a photon @@ -661,7 +694,7 @@ void AliAnalysisTaskGammaConversion::ProcessMCData(){ } } - + if(particle->GetPdgCode()==111){ //Pi0 @@ -696,7 +729,7 @@ void AliAnalysisTaskGammaConversion::ProcessMCData(){ fHistograms->FillHistogram("MC_Pi0_ZR", particle->Vz(),particle->R()); fHistograms->FillHistogram("MC_Pi0_GammaDaughter_OpeningAngle", GetMCOpeningAngle(daughter0,daughter1)); fHistograms->FillHistogram("MC_Pi0_XY", particle->Vx(), particle->Vy());//only fill from one daughter to avoid multiple filling - + if(gammaEtaCut && gammaRCut){ // if(TMath::Abs(daughter0->Eta()) <= fV0Reader->GetEtaCut() && TMath::Abs(daughter1->Eta()) <= fV0Reader->GetEtaCut() ){ fHistograms->FillHistogram("MC_Pi0_Pt_Eta_withinAcceptance", particle->Pt(),particle->Eta()); @@ -722,7 +755,7 @@ void AliAnalysisTaskGammaConversion::ProcessMCData(){ fHistograms->FillHistogram("MC_Eta_XY", particle->Vx(), particle->Vy());//only fill from one daughter to avoid multiple filling if(gammaEtaCut && gammaRCut){ - // if(TMath::Abs(daughter0->Eta()) <= fV0Reader->GetEtaCut() && TMath::Abs(daughter1->Eta()) <= fV0Reader->GetEtaCut() ){ + // if(TMath::Abs(daughter0->Eta()) <= fV0Reader->GetEtaCut() && TMath::Abs(daughter1->Eta()) <= fV0Reader->GetEtaCut() ){ fHistograms->FillHistogram("MC_Eta_Pt_Eta_withinAcceptance", particle->Pt(),particle->Eta()); fHistograms->FillHistogram("MC_Eta_Pt_Rapid_withinAcceptance", particle->Pt(),rapidity); if(daughter0Electron && daughter0Positron && daughter1Electron && daughter1Positron){ @@ -747,9 +780,9 @@ void AliAnalysisTaskGammaConversion::ProcessMCData(){ fHistograms->FillHistogram("MC_Mother_Rapid", rapidity); fHistograms->FillHistogram("MC_Mother_Phi",tmpPhi); fHistograms->FillHistogram("MC_Mother_InvMass_vs_Pt",particle->GetMass(),particle->Pt()); - + if(gammaEtaCut && gammaRCut){ - // if(TMath::Abs(daughter0->Eta()) <= fV0Reader->GetEtaCut() && TMath::Abs(daughter1->Eta()) <= fV0Reader->GetEtaCut() ){ + // if(TMath::Abs(daughter0->Eta()) <= fV0Reader->GetEtaCut() && TMath::Abs(daughter1->Eta()) <= fV0Reader->GetEtaCut() ){ fHistograms->FillHistogram("MC_Mother_Pt_Eta_withinAcceptance", particle->Pt(),particle->Eta()); fHistograms->FillHistogram("MC_Mother_Pt_Rapid_withinAcceptance", particle->Pt(),rapidity); fHistograms->FillHistogram("MC_Mother_InvMass_vs_Pt_withinAcceptance",particle->GetMass(),particle->Pt()); @@ -757,25 +790,23 @@ void AliAnalysisTaskGammaConversion::ProcessMCData(){ fHistograms->FillHistogram("MC_Mother_Pt_Eta_ConvGamma_withinAcceptance", particle->Pt(),particle->Eta()); fHistograms->FillHistogram("MC_Mother_Pt_Rapid_ConvGamma_withinAcceptance", particle->Pt(),rapidity); fHistograms->FillHistogram("MC_Mother_InvMass_vs_Pt_ConvGamma_withinAcceptance",particle->GetMass(),particle->Pt()); - + } } // end passed R and eta cut - + } // end if(particle->GetNDaughters() == 2) }// end for (Int_t iTracks = 0; iTracks < fStack->GetNtrack(); iTracks++) - - //cout << "right before the end of processMCdata" << endl; - + } // end ProcessMCData void AliAnalysisTaskGammaConversion::FillNtuple(){ //Fills the ntuple with the different values - + if(fGammaNtuple == NULL){ return; } @@ -825,33 +856,33 @@ void AliAnalysisTaskGammaConversion::FillNtuple(){ void AliAnalysisTaskGammaConversion::ProcessV0sNoCut(){ // Process all the V0's without applying any cuts to it - + Int_t numberOfV0s = fV0Reader->GetNumberOfV0s(); for(Int_t i=0;iGetV0(i); - + if(fV0Reader->CheckForPrimaryVertex() == kFALSE){ return; } - + if(fDoMCTruth){ - + if(fV0Reader->HasSameMCMother() == kFALSE){ continue; } - + TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle(); TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle(); - + if(TMath::Abs(negativeMC->GetPdgCode())!=11 || TMath::Abs(positiveMC->GetPdgCode())!=11){ continue; } if(negativeMC->GetPdgCode()==positiveMC->GetPdgCode()){ continue; } - + if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22){ - + fHistograms->FillHistogram("ESD_NoCutConvGamma_Pt", fV0Reader->GetMotherCandidatePt()); fHistograms->FillHistogram("ESD_NoCutConvGamma_Energy", fV0Reader->GetMotherCandidateEnergy()); fHistograms->FillHistogram("ESD_NoCutConvGamma_Eta", fV0Reader->GetMotherCandidateEta()); @@ -862,7 +893,7 @@ void AliAnalysisTaskGammaConversion::ProcessV0sNoCut(){ fHistograms->FillHistogram("ESD_NoCutConvGamma_NDF", fV0Reader->GetMotherCandidateNDF()); fHistograms->FillHistogram("ESD_NoCutConvGamma_Rapid", fV0Reader->GetMotherCandidateRapidity()); fHistograms->FillHistogram("ESD_NoCutConvGamma_Pt_Eta", fV0Reader->GetMotherCandidatePt(),fV0Reader->GetMotherCandidateEta()); - + fHistograms->FillHistogram("ESD_NoCutConvGamma_Pt_Chi2", fV0Reader->GetMotherCandidatePt(), fV0Reader->GetMotherCandidateChi2()); fHistograms->FillHistogram("ESD_NoCutConvGamma_Eta_Chi2", fV0Reader->GetMotherCandidateEta(), fV0Reader->GetMotherCandidateChi2()); @@ -875,8 +906,6 @@ void AliAnalysisTaskGammaConversion::ProcessV0sNoCut(){ fHistograms->FillHistogram("ESD_NoCutConvGamma_MC_Pt_Eta", fV0Reader->GetMotherMCParticle()->Pt(),fV0Reader->GetMotherMCParticle()->Eta()); fHistograms->FillHistogram("ESD_NoCutConversion_MC_ZR", negativeMC->Vz(),negativeMC->R()); fHistograms->FillHistogram("ESD_NoCutConversion_MC_XY", negativeMC->Vx(),negativeMC->Vy()); - - } } } @@ -926,7 +955,6 @@ void AliAnalysisTaskGammaConversion::ProcessV0s(){ fHistograms->FillHistogram("ESD_ConvGamma_Eta_Chi2", fV0Reader->GetMotherCandidateEta(), fV0Reader->GetMotherCandidateChi2()); - // begin mapping Int_t rBin = fHistograms->GetRBin(fV0Reader->GetXYRadius()); Int_t phiBin = fHistograms->GetPhiBin(fV0Reader->GetNegativeTrackPhi()); @@ -950,18 +978,18 @@ void AliAnalysisTaskGammaConversion::ProcessV0s(){ // end mapping new((*fKFReconstructedGammasTClone)[fKFReconstructedGammasTClone->GetEntriesFast()]) AliKFParticle(*fV0Reader->GetMotherCandidateKFCombination()); - + // fKFReconstructedGammas.push_back(*fV0Reader->GetMotherCandidateKFCombination()); fElectronv1.push_back(fV0Reader->GetCurrentV0()->GetPindex()); fElectronv2.push_back(fV0Reader->GetCurrentV0()->GetNindex()); - + //----------------------------------- checking for "real" conversions (MC match) -------------------------------------- if(fDoMCTruth){ TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle(); // TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle(); // not used any longer - + if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22){ fHistograms->FillHistogram("ESD_TrueConvGamma_Pt", fV0Reader->GetMotherCandidatePt()); @@ -978,21 +1006,21 @@ void AliAnalysisTaskGammaConversion::ProcessV0s(){ fHistograms->FillHistogram("ESD_TrueConvGamma_TrackLength", /*fV0Reader->GetPositiveTrackLength()*/fV0Reader->GetPositiveNTPCClusters()); fHistograms->FillHistogram("ESD_TrueConvGamma_TrackLengthVSInvMass",/*fV0Reader->GetNegativeTrackLength()*/fV0Reader->GetNegativeNTPCClusters(),fV0Reader->GetMotherCandidateMass()); fHistograms->FillHistogram("ESD_TrueConvGamma_TrackLengthVSInvMass",/*fV0Reader->GetPositiveTrackLength()*/fV0Reader->GetPositiveNTPCClusters(),fV0Reader->GetMotherCandidateMass()); - + fHistograms->FillHistogram("ESD_TrueConvGamma_Pt_Chi2", fV0Reader->GetMotherCandidatePt(), fV0Reader->GetMotherCandidateChi2()); fHistograms->FillHistogram("ESD_TrueConvGamma_Eta_Chi2", fV0Reader->GetMotherCandidateEta(), fV0Reader->GetMotherCandidateChi2()); - + fHistograms->FillHistogram("ESD_TrueConversion_XY", fV0Reader->GetX(),fV0Reader->GetY()); fHistograms->FillHistogram("ESD_TrueConversion_R", fV0Reader->GetXYRadius()); fHistograms->FillHistogram("ESD_TrueConversion_ZR", fV0Reader->GetZ(),fV0Reader->GetXYRadius()); fHistograms->FillHistogram("ESD_TrueConversion_OpeningAngle", fV0Reader->GetOpeningAngle()); - + //store MCTruth properties fHistograms->FillHistogram("ESD_TrueConvGamma_MC_Pt_Eta", fV0Reader->GetMotherMCParticle()->Pt(),fV0Reader->GetMotherMCParticle()->Eta()); fHistograms->FillHistogram("ESD_TrueConversion_MC_ZR", negativeMC->Vz(),negativeMC->R()); fHistograms->FillHistogram("ESD_TrueConversion_MC_XY", negativeMC->Vx(),negativeMC->Vy()); - + //resolution Double_t mcpt = fV0Reader->GetMotherMCParticle()->Pt(); Double_t esdpt = fV0Reader->GetMotherCandidatePt(); @@ -1018,7 +1046,7 @@ void AliAnalysisTaskGammaConversion::ProcessV0s(){ if(fV0Reader->GetNegativeMCParticle()->R() != 0){ resdR = ((fV0Reader->GetXYRadius() - fV0Reader->GetNegativeMCParticle()->R())/fV0Reader->GetNegativeMCParticle()->R())*100; } - + fHistograms->FillHistogram("Resolution_dR", fV0Reader->GetNegativeMCParticle()->R(), resdR); fHistograms->FillHistogram("Resolution_MC_R", fV0Reader->GetNegativeMCParticle()->R()); fHistograms->FillHistogram("Resolution_ESD_R", fV0Reader->GetXYRadius()); @@ -1028,32 +1056,29 @@ void AliAnalysisTaskGammaConversion::ProcessV0s(){ }//while(fV0Reader->NextV0) fHistograms->FillHistogram("ESD_NumberOfSurvivingV0s", nSurvivingV0s); fHistograms->FillHistogram("ESD_NumberOfV0s", fV0Reader->GetNumberOfV0s()); - - //cout << "nearly at the end of doMCTruth" << endl; - } void AliAnalysisTaskGammaConversion::FillAODWithConversionGammas(){ // Fill AOD with reconstructed Gamma - + for(Int_t gammaIndex=0;gammaIndexGetEntriesFast();gammaIndex++){ // for(UInt_t gammaIndex=0;gammaIndexPx(),gammakf->Py(),gammakf->Pz(), gammakf->E()); gamma.SetLabel(-1);//How to get the MC label of the reconstructed gamma? - gamma.SetCaloLabel(-1,-1); //How to get the MC label of the 2 electrons that form the gamma? - gamma.SetDetector("CTS"); //tag the gamma as reconstructed in the central barrel - gamma.SetPdg(AliCaloPID::kPhotonConv); //photon id - gamma.SetTag(-1); //Here I usually put a flag saying that montecarlo says it is prompt, decay fragmentation photon, or hadrons or whatever - - //Add it to the aod list - Int_t i = fAODBranch->GetEntriesFast(); - new((*fAODBranch)[i]) AliAODPWG4Particle(gamma); + gamma.SetCaloLabel(-1,-1); //How to get the MC label of the 2 electrons that form the gamma? + gamma.SetDetector("CTS"); //tag the gamma as reconstructed in the central barrel + gamma.SetPdg(AliCaloPID::kPhotonConv); //photon id + gamma.SetTag(-1); //Here I usually put a flag saying that montecarlo says it is prompt, decay fragmentation photon, or hadrons or whatever + + //Add it to the aod list + Int_t i = fAODBranch->GetEntriesFast(); + new((*fAODBranch)[i]) AliAODPWG4Particle(gamma); */ // AliKFParticle * gammakf = &fKFReconstructedGammas[gammaIndex]; AliKFParticle * gammakf = (AliKFParticle *)fKFReconstructedGammasTClone->At(gammaIndex); @@ -1065,8 +1090,8 @@ void AliAnalysisTaskGammaConversion::FillAODWithConversionGammas(){ aodObject.SetLabel2(fElectronv2[gammaIndex]); Int_t i = fAODBranch->GetEntriesFast(); new((*fAODBranch)[i]) AliGammaConversionAODObject(aodObject); - } - + } + } @@ -1080,17 +1105,17 @@ void AliAnalysisTaskGammaConversion::ProcessGammasForNeutralMesonAnalysis(){ // AliKFParticle * twoGammaDecayCandidateDaughter0 = &fKFReconstructedGammas[firstGammaIndex]; // AliKFParticle * twoGammaDecayCandidateDaughter1 = &fKFReconstructedGammas[secondGammaIndex]; - + AliKFParticle * twoGammaDecayCandidateDaughter0 = (AliKFParticle *)fKFReconstructedGammasTClone->At(firstGammaIndex); AliKFParticle * twoGammaDecayCandidateDaughter1 = (AliKFParticle *)fKFReconstructedGammasTClone->At(secondGammaIndex); - + if(fElectronv1[firstGammaIndex]==fElectronv1[secondGammaIndex] || fElectronv1[firstGammaIndex]==fElectronv2[secondGammaIndex]){ continue; } if(fElectronv2[firstGammaIndex]==fElectronv1[secondGammaIndex] || fElectronv2[firstGammaIndex]==fElectronv2[secondGammaIndex]){ continue; } - + AliKFParticle *twoGammaCandidate = new AliKFParticle(*twoGammaDecayCandidateDaughter0,*twoGammaDecayCandidateDaughter1); Double_t massTwoGammaCandidate = 0.; @@ -1115,7 +1140,7 @@ void AliAnalysisTaskGammaConversion::ProcessGammasForNeutralMesonAnalysis(){ } if(openingAngleTwoGammaCandidate < fMinOpeningAngleGhostCut) continue; // minimum opening angle to avoid using ghosttracks - + fHistograms->FillHistogram("ESD_Mother_GammaDaughter_OpeningAngle", openingAngleTwoGammaCandidate); fHistograms->FillHistogram("ESD_Mother_Energy", twoGammaCandidate->GetE()); fHistograms->FillHistogram("ESD_Mother_Pt", momentumVectorTwoGammaCandidate.Pt()); @@ -1131,10 +1156,6 @@ void AliAnalysisTaskGammaConversion::ProcessGammasForNeutralMesonAnalysis(){ } } delete twoGammaCandidate; - - //cout << "nearly at the end of processgamma for neutral meson ..." << endl; - - } } } @@ -1144,12 +1165,12 @@ void AliAnalysisTaskGammaConversion::CalculateBackground(){ vector vectorCurrentEventGoodV0s = fV0Reader->GetCurrentEventGoodV0s(); vector vectorPreviousEventGoodV0s = fV0Reader->GetPreviousEventGoodV0s(); - + for(UInt_t iCurrent=0;iCurrentGetX(),backgroundCandidate->GetY(),backgroundCandidate->GetZ()); Double_t openingAngleBG = currentEventGoodV0->GetAngle(*previousGoodV0); - + Double_t rapidity; if(backgroundCandidate->GetE() - backgroundCandidate->GetPz() == 0 || backgroundCandidate->GetE() + backgroundCandidate->GetPz() == 0) rapidity=0; else rapidity = 0.5*(TMath::Log((backgroundCandidate->GetE() +backgroundCandidate->GetPz()) / (backgroundCandidate->GetE()-backgroundCandidate->GetPz()))); - + if(openingAngleBG < fMinOpeningAngleGhostCut ) continue; // minimum opening angle to avoid using ghosttracks - + fHistograms->FillHistogram("ESD_Background_GammaDaughter_OpeningAngle", openingAngleBG); fHistograms->FillHistogram("ESD_Background_Energy", backgroundCandidate->GetE()); @@ -1190,8 +1211,6 @@ void AliAnalysisTaskGammaConversion::CalculateBackground(){ } } delete backgroundCandidate; - //cout << "nearly at the end of background" << endl; - } } } @@ -1200,20 +1219,20 @@ void AliAnalysisTaskGammaConversion::CalculateBackground(){ void AliAnalysisTaskGammaConversion::ProcessGammasForGammaJetAnalysis(){ //ProcessGammasForGammaJetAnalysis - + Double_t distIsoMin; - + CreateListOfChargedParticles(); - - + + // for(UInt_t gammaIndex=0;gammaIndexGetEntriesFast();gammaIndex++){ AliKFParticle * currentGamma = (AliKFParticle*)fKFReconstructedGammasTClone->At(gammaIndex); TVector3 momentumVectorCurrentGamma(currentGamma->GetPx(),currentGamma->GetPy(),currentGamma->GetPz()); - + if( momentumVectorCurrentGamma.Pt()> fMinPtForGammaJet){ distIsoMin=GetMinimumDistanceToCharge(gammaIndex); - + if (distIsoMin > fMinIsoConeSize && fLeadingChargedIndex>=0){ CalculateJetCone(gammaIndex); } @@ -1223,15 +1242,15 @@ void AliAnalysisTaskGammaConversion::ProcessGammasForGammaJetAnalysis(){ void AliAnalysisTaskGammaConversion::CreateListOfChargedParticles(){ // CreateListOfChargedParticles - + fESDEvent = fV0Reader->GetESDEvent(); for(Int_t iTracks = 0; iTracks < fESDEvent->GetNumberOfTracks(); iTracks++){ AliESDtrack* curTrack = fESDEvent->GetTrack(iTracks); - + if(!curTrack){ continue; } - + if(fEsdTrackCuts->AcceptTrack(curTrack) ){ new((*fChargedParticles)[fChargedParticles->GetEntriesFast()]) AliESDtrack(*curTrack); // fChargedParticles.push_back(curTrack); @@ -1241,25 +1260,25 @@ void AliAnalysisTaskGammaConversion::CreateListOfChargedParticles(){ } void AliAnalysisTaskGammaConversion::CalculateJetCone(Int_t gammaIndex){ // CaculateJetCone - + Double_t cone; Double_t coneSize=0.3; Double_t ptJet=0; - + // AliKFParticle * currentGamma = &fKFReconstructedGammas[gammaIndex]; AliKFParticle * currentGamma = (AliKFParticle*)fKFReconstructedGammasTClone->At(gammaIndex); TVector3 momentumVectorCurrentGamma(currentGamma->GetPx(),currentGamma->GetPy(),currentGamma->GetPz()); - + AliESDtrack* leadingCharged = (AliESDtrack*)(fChargedParticles->At(fLeadingChargedIndex)); Double_t momLeadingCharged[3]; leadingCharged->GetConstrainedPxPyPz(momLeadingCharged); - + TVector3 momentumVectorLeadingCharged(momLeadingCharged[0],momLeadingCharged[1],momLeadingCharged[2]); - + Double_t phi1=momentumVectorLeadingCharged.Phi(); Double_t eta1=momentumVectorLeadingCharged.Eta(); Double_t phi3=momentumVectorCurrentGamma.Phi(); - + for(Int_t iCh=0;iChGetEntriesFast();iCh++){ AliESDtrack* curTrack = (AliESDtrack*)(fChargedParticles->At(iCh)); Int_t chId = fChargedParticlesId[iCh]; @@ -1269,8 +1288,8 @@ void AliAnalysisTaskGammaConversion::CalculateJetCone(Int_t gammaIndex){ TVector3 momentumVectorChargedParticle(mom[0],mom[1],mom[2]); Double_t phi2=momentumVectorChargedParticle.Phi(); Double_t eta2=momentumVectorChargedParticle.Eta(); - - + + cone=100.; if( TMath::Abs(phi2 - phi1) <= ( TMath::TwoPi()-coneSize) ){ cone = TMath::Sqrt( TMath::Power((eta2-eta1),2)+ TMath::Power((phi2-phi1),2) ); @@ -1282,45 +1301,45 @@ void AliAnalysisTaskGammaConversion::CalculateJetCone(Int_t gammaIndex){ cone = TMath::Sqrt( TMath::Power((eta2-eta1),2)+ TMath::Power((phi2+TMath::TwoPi()-phi1),2) ); } } - + if(cone fMinPtJetCone ){ ptJet+= momentumVectorChargedParticle.Pt(); Double_t ffzHdrGam = momentumVectorChargedParticle.Pt()/momentumVectorCurrentGamma.Pt(); Double_t imbalanceHdrGam=-momentumVectorChargedParticle.Dot(momentumVectorCurrentGamma)/momentumVectorCurrentGamma.Mag2(); fHistograms->FillHistogram("ESD_FFzHdrGam",ffzHdrGam); fHistograms->FillHistogram("ESD_ImbalanceHdrGam",imbalanceHdrGam); - + } - + Double_t dphiHdrGam=phi3-phi2; if ( dphiHdrGam < (-TMath::PiOver2())){ dphiHdrGam+=(TMath::TwoPi()); } - + if ( dphiHdrGam > (3.*TMath::PiOver2()) ){ dphiHdrGam-=(TMath::TwoPi()); } - + if (momentumVectorChargedParticle.Pt()>fMinPtGamChargedCorr){ fHistograms->FillHistogram("ESD_dphiHdrGamIsolated",dphiHdrGam); } }//track loop - - + + } Double_t AliAnalysisTaskGammaConversion::GetMinimumDistanceToCharge(Int_t indexHighestPtGamma){ // GetMinimumDistanceToCharge - + Double_t fIsoMin=100.; Double_t ptLeadingCharged=-1.; - + AliKFParticle * gammaHighestPt = (AliKFParticle*)fKFReconstructedGammasTClone->At(indexHighestPtGamma); TVector3 momentumVectorgammaHighestPt(gammaHighestPt->GetPx(),gammaHighestPt->GetPy(),gammaHighestPt->GetPz()); - + Double_t phi1=momentumVectorgammaHighestPt.Phi(); Double_t eta1=momentumVectorgammaHighestPt.Eta(); - + for(Int_t iCh=0;iChGetEntriesFast();iCh++){ AliESDtrack* curTrack = (AliESDtrack*)(fChargedParticles->At(iCh)); Int_t chId = fChargedParticlesId[iCh]; @@ -1331,57 +1350,57 @@ Double_t AliAnalysisTaskGammaConversion::GetMinimumDistanceToCharge(Int_t indexH Double_t phi2=momentumVectorChargedParticle.Phi(); Double_t eta2=momentumVectorChargedParticle.Eta(); Double_t iso=pow( (pow( (eta1-eta2),2)+ pow((phi1-phi2),2)),0.5 ); - + if(momentumVectorChargedParticle.Pt()>fMinPtIsoCone ){ if (iso (3.*TMath::PiOver2()) ){ dphiHdrGam-=(TMath::TwoPi()); } if (momentumVectorChargedParticle.Pt()>fMinPtGamChargedCorr){ fHistograms->FillHistogram("ESD_dphiHdrGam",dphiHdrGam); } - + if (dphiHdrGam>0.9*TMath::Pi() && dphiHdrGam<1.1*TMath::Pi()){ - if (momentumVectorChargedParticle.Pt()> ptLeadingCharged && momentumVectorChargedParticle.Pt()>0.1*momentumVectorgammaHighestPt.Pt()){ + if (momentumVectorChargedParticle.Pt()> ptLeadingCharged && momentumVectorChargedParticle.Pt()>0.1*momentumVectorgammaHighestPt.Pt()){ ptLeadingCharged=momentumVectorChargedParticle.Pt(); fLeadingChargedIndex=iCh; } } - + }//track loop fHistograms->FillHistogram("ESD_MinimumIsoDistance",fIsoMin); return fIsoMin; - + } - Int_t AliAnalysisTaskGammaConversion::GetIndexHighestPtGamma(){ - //GetIndexHighestPtGamma - +Int_t AliAnalysisTaskGammaConversion::GetIndexHighestPtGamma(){ + //GetIndexHighestPtGamma + Int_t indexHighestPtGamma=-1; //Double_t fGammaPtHighest = -100.; - + for(Int_t firstGammaIndex=0;firstGammaIndexGetEntriesFast();firstGammaIndex++){ AliKFParticle * gammaHighestPtCandidate = (AliKFParticle*)fKFReconstructedGammasTClone->At(firstGammaIndex); - TVector3 momentumVectorgammaHighestPtCandidate(gammaHighestPtCandidate->GetPx(),gammaHighestPtCandidate->GetPy(),gammaHighestPtCandidate->GetPz()); - if (momentumVectorgammaHighestPtCandidate.Pt() > fGammaPtHighest){ - fGammaPtHighest=momentumVectorgammaHighestPtCandidate.Pt(); - //gammaHighestPt = gammaHighestPtCandidate; - indexHighestPtGamma=firstGammaIndex; - } + TVector3 momentumVectorgammaHighestPtCandidate(gammaHighestPtCandidate->GetPx(),gammaHighestPtCandidate->GetPy(),gammaHighestPtCandidate->GetPz()); + if (momentumVectorgammaHighestPtCandidate.Pt() > fGammaPtHighest){ + fGammaPtHighest=momentumVectorgammaHighestPtCandidate.Pt(); + //gammaHighestPt = gammaHighestPtCandidate; + indexHighestPtGamma=firstGammaIndex; + } } - + return indexHighestPtGamma; - + } @@ -1398,7 +1417,7 @@ void AliAnalysisTaskGammaConversion::UserCreateOutputObjects() fAODBranch = new TClonesArray("AliGammaConversionAODObject", 0); fAODBranch->SetName(fAODBranchName); AddAODBranch("TClonesArray", &fAODBranch); - + // Create the output container if(fOutputContainer != NULL){ delete fOutputContainer; @@ -1437,15 +1456,15 @@ Double_t AliAnalysisTaskGammaConversion::GetMCOpeningAngle(TParticle* const daug void AliAnalysisTaskGammaConversion::CheckV0Efficiency(){ // see header file for documentation - + vector indexOfGammaParticle; - + fStack = fV0Reader->GetMCStack(); - + if(fV0Reader->CheckForPrimaryVertex() == kFALSE){ return; // aborts if the primary vertex does not have contributors. } - + for (Int_t iTracks = 0; iTracks < fStack->GetNprimary(); iTracks++) { TParticle* particle = (TParticle *)fStack->Particle(iTracks); if(particle->GetPdgCode()==22){ //Gamma @@ -1471,28 +1490,28 @@ void AliAnalysisTaskGammaConversion::CheckV0Efficiency(){ } } } - + Int_t nFoundGammas=0; Int_t nNotFoundGammas=0; - + Int_t numberOfV0s = fV0Reader->GetNumberOfV0s(); for(Int_t i=0;iGetV0(i); - + if(fV0Reader->HasSameMCMother() == kFALSE){ continue; } - + TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle(); TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle(); - + if(TMath::Abs(negativeMC->GetPdgCode())!=11 || TMath::Abs(positiveMC->GetPdgCode())!=11){ continue; } if(negativeMC->GetPdgCode()==positiveMC->GetPdgCode()){ continue; } - + if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22){ //TParticle * v0Gamma = fV0Reader->GetMotherMCParticle(); for(UInt_t mcIndex=0;mcIndexGetESDEvent(); - - - TClonesArray * vESDeNegTemp = new TClonesArray("AliESDtrack",0); + + + TClonesArray * vESDeNegTemp = new TClonesArray("AliESDtrack",0); TClonesArray * vESDePosTemp = new TClonesArray("AliESDtrack",0); TClonesArray * vESDxNegTemp = new TClonesArray("AliESDtrack",0); TClonesArray * vESDxPosTemp = new TClonesArray("AliESDtrack",0); TClonesArray * vESDeNegNoJPsi = new TClonesArray("AliESDtrack",0); TClonesArray * vESDePosNoJPsi = new TClonesArray("AliESDtrack",0); - + /* vector vESDeNegTemp(0); vector vESDePosTemp(0); @@ -1530,113 +1548,113 @@ void AliAnalysisTaskGammaConversion::ProcessGammaElectronsForChicAnalysis(){ vector vESDeNegNoJPsi(0); vector vESDePosNoJPsi(0); */ - - + + fHistograms->FillTable("Table_Electrons",0);//Count number of Events - + for(Int_t iTracks = 0; iTracks < fESDEvent->GetNumberOfTracks(); iTracks++){ AliESDtrack* curTrack = fESDEvent->GetTrack(iTracks); - + if(!curTrack){ //print warning here continue; } - + double p[3];if(!curTrack->GetConstrainedPxPyPz(p))continue; double r[3];curTrack->GetConstrainedXYZ(r); - + TVector3 rXYZ(r); - + fHistograms->FillTable("Table_Electrons",4);//Count number of ESD tracks - + Bool_t flagKink = kTRUE; Bool_t flagTPCrefit = kTRUE; Bool_t flagTRDrefit = kTRUE; Bool_t flagITSrefit = kTRUE; Bool_t flagTRDout = kTRUE; Bool_t flagVertex = kTRUE; - - + + //Cuts --------------------------------------------------------------- - + if(curTrack->GetKinkIndex(0) > 0){ fHistograms->FillHistogram("Table_Electrons",5);//Count kink flagKink = kFALSE; } - + ULong_t trkStatus = curTrack->GetStatus(); - + ULong_t tpcRefit = (trkStatus & AliESDtrack::kTPCrefit); - + if(!tpcRefit){ fHistograms->FillHistogram("Table_Electrons",9);//Count not TPCrefit flagTPCrefit = kFALSE; } - + ULong_t itsRefit = (trkStatus & AliESDtrack::kITSrefit); if(!itsRefit){ fHistograms->FillHistogram("Table_Electrons",10);//Count not ITSrefit flagITSrefit = kFALSE; } - + ULong_t trdRefit = (trkStatus & AliESDtrack::kTRDrefit); - + if(!trdRefit){ fHistograms->FillHistogram("Table_Electrons",8); //Count not TRDrefit flagTRDrefit = kFALSE; } - + ULong_t trdOut = (trkStatus & AliESDtrack::kTRDout); - + if(!trdOut) { fHistograms->FillHistogram("Table_Electrons",7); //Count not TRDout flagTRDout = kFALSE; } - + double nsigmaToVxt = GetSigmaToVertex(curTrack); - + if(nsigmaToVxt > 3){ fHistograms->FillHistogram("Table_Electrons",6); //Count Tracks with number of sigmas > 3 flagVertex = kFALSE; } - + if(! (flagKink && flagTPCrefit && flagITSrefit && flagTRDrefit && flagTRDout && flagVertex ) ) continue; fHistograms->FillHistogram("Table_Electrons",11);//Count Tracks passed Cuts - - + + Stat_t pid, weight; GetPID(curTrack, pid, weight); - + if(pid!=0){ fHistograms->FillHistogram("Table_Electrons",12); //Count Tracks with pid != 0 } - + if(pid == 0){ fHistograms->FillHistogram("Table_Electrons",13); //Count Tracks with pid != 0 } - - - - + + + + Int_t labelMC = TMath::Abs(curTrack->GetLabel()); TParticle* curParticle = fStack->Particle(labelMC); - - - - + + + + TLorentzVector curElec; curElec.SetXYZM(p[0],p[1],p[2],fElectronMass); - - - - + + + + if(curTrack->GetSign() > 0){ - + // vESDxPosTemp.push_back(curTrack); new((*vESDxPosTemp)[vESDxPosTemp->GetEntriesFast()]) AliESDtrack(*curTrack); - + if( pid == 0){ - + fHistograms->FillHistogram("ESD_ElectronPosNegPt",curElec.Pt()); fHistograms->FillHistogram("ESD_ElectronPosPt",curElec.Pt()); fHistograms->FillHistogram("MC_ElectronPosNegPt",curParticle->Pt()); @@ -1644,23 +1662,23 @@ void AliAnalysisTaskGammaConversion::ProcessGammaElectronsForChicAnalysis(){ fHistograms->FillHistogram("MC_ElectronPosNegEta",curParticle->Eta()); // vESDePosTemp.push_back(curTrack); new((*vESDePosTemp)[vESDePosTemp->GetEntriesFast()]) AliESDtrack(*curTrack); - + } - + } else { // vESDxNegTemp.push_back(curTrack); - if(vESDxNegTemp == NULL){ - cout<<"TCloes is zero god damn it"<GetEntriesFast()]) AliESDtrack(*curTrack); - + if( pid == 0){ - + fHistograms->FillHistogram("ESD_ElectronPosNegPt",curElec.Pt()); fHistograms->FillHistogram("ESD_ElectronNegPt",curElec.Pt()); fHistograms->FillHistogram("MC_ElectronPosNegPt",curParticle->Pt()); @@ -1668,21 +1686,21 @@ void AliAnalysisTaskGammaConversion::ProcessGammaElectronsForChicAnalysis(){ fHistograms->FillHistogram("MC_ElectronPosNegEta",curParticle->Eta()); //vESDeNegTemp.push_back(curTrack); new((*vESDeNegTemp)[vESDeNegTemp->GetEntriesFast()]) AliESDtrack(*curTrack); - + } - + } - + } - - + + Bool_t ePosJPsi = kFALSE; Bool_t eNegJPsi = kFALSE; Bool_t ePosPi0 = kFALSE; Bool_t eNegPi0 = kFALSE; UInt_t iePosJPsi=0,ieNegJPsi=0,iePosPi0=0,ieNegPi0=0; - + for(Int_t iNeg=0; iNeg < vESDeNegTemp->GetEntriesFast(); iNeg++){ if(fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDeNegTemp->At(iNeg)))->GetLabel()))->GetPdgCode() == 11) if(fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDeNegTemp->At(iNeg)))->GetLabel()))->GetMother(0) > -1){ @@ -1702,10 +1720,10 @@ void AliAnalysisTaskGammaConversion::ProcessGammaElectronsForChicAnalysis(){ new((*vESDeNegNoJPsi)[vESDeNegNoJPsi->GetEntriesFast()]) AliESDtrack(*(AliESDtrack*)(vESDeNegTemp->At(iNeg))); // cout<<"ESD No Positivo JPsi "<GetEntriesFast(); iPos++){ if(fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDePosTemp->At(iPos)))->GetLabel()))->GetPdgCode() == -11) if(fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDePosTemp->At(iPos)))->GetLabel()))->GetMother(0) > -1){ @@ -1725,7 +1743,7 @@ void AliAnalysisTaskGammaConversion::ProcessGammaElectronsForChicAnalysis(){ new((*vESDePosNoJPsi)[vESDePosNoJPsi->GetEntriesFast()]) AliESDtrack(*(AliESDtrack*)(vESDePosTemp->At(iPos))); // cout<<"ESD No Negativo JPsi "<FillHistogram("MC_ElectronPosNegPi0Angle",GetMCOpeningAngle(fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDeNegTemp->At(ieNegPi0)))->GetLabel())), fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDePosTemp->At(iePosPi0)))->GetLabel())))); } - - + + FillAngle("ESD_eNegePosAngleBeforeCut",GetTLorentzVector(vESDeNegTemp),GetTLorentzVector(vESDePosTemp)); - + CleanWithAngleCuts(*vESDeNegTemp,*vESDePosTemp,*fKFReconstructedGammasTClone); // vector vCurrentTLVeNeg = GetTLorentzVector(fCurrentEventNegElectron); // vector vCurrentTLVePos = GetTLorentzVector(fCurrentEventPosElectron); - + TClonesArray vCurrentTLVeNeg = GetTLorentzVector(fCurrentEventNegElectronTClone); TClonesArray vCurrentTLVePos = GetTLorentzVector(fCurrentEventPosElectronTClone); - - + + FillAngle("ESD_eNegePosAngleAfterCut",vCurrentTLVeNeg,vCurrentTLVePos); - - - - + + + + //FillAngle("ESD_eNegePosAngleAfterCut",CurrentTLVeNeg,CurrentTLVePos); - - + + FillElectronInvMass("ESD_InvMass_ePluseMinus",vCurrentTLVeNeg,vCurrentTLVePos); FillElectronInvMass("ESD_InvMass_xPlusxMinus",GetTLorentzVector(vESDxNegTemp),GetTLorentzVector(vESDxPosTemp)); - - - + + + FillGammaElectronInvMass("ESD_InvMass_GammaePluseMinusChiC","ESD_InvMass_GammaePluseMinusChiCDiff",*fKFReconstructedGammasCutTClone,vCurrentTLVeNeg,vCurrentTLVePos); - + FillGammaElectronInvMass("ESD_InvMass_GammaePluseMinusPi0","ESD_InvMass_GammaePluseMinusPi0Diff", *fKFReconstructedGammasCutTClone,vCurrentTLVeNeg,vCurrentTLVePos); - + //BackGround - + //Like Sign e+e- ElectronBackground("ESD_ENegBackground",vCurrentTLVeNeg); ElectronBackground("ESD_EPosBackground",vCurrentTLVePos); ElectronBackground("ESD_EPosENegBackground",vCurrentTLVeNeg); ElectronBackground("ESD_EPosENegBackground",vCurrentTLVePos); - + // Like Sign e+e- no JPsi ElectronBackground("ESD_EPosENegNoJPsiBG",GetTLorentzVector(vESDeNegNoJPsi)); ElectronBackground("ESD_EPosENegNoJPsiBG",GetTLorentzVector(vESDePosNoJPsi)); - + //Mixed Event - + if( fCurrentEventPosElectronTClone->GetEntriesFast() > 0 && fCurrentEventNegElectronTClone->GetEntriesFast() > 0 && fKFReconstructedGammasCutTClone->GetEntriesFast() > 0 ){ FillGammaElectronInvMass("ESD_EPosENegGammaBackgroundMX","ESD_EPosENegGammaBackgroundMXDiff", *fKFReconstructedGammasCutTClone,*fPreviousEventTLVNegElectronTClone,*fPreviousEventTLVPosElectronTClone); *fPreviousEventTLVNegElectronTClone = vCurrentTLVeNeg; *fPreviousEventTLVPosElectronTClone = vCurrentTLVePos; - + } - + /* //Photons P Double_t vtx[3]; vtx[0]=0;vtx[1]=0;vtx[2]=0; for(UInt_t i=0;iGetMother(0) < 0) continue; - - - + + + Int_t tempLabel = fStack->Particle(fMCGammaChicTempCut[i]->GetMother(0))->GetPdgCode(); // cout<<" Label Pedro Gonzalez " <FillHistogram("ESD_PhotonsMomentum",fKFReconstructedGammasChic[i].GetMomentum()); - - + + } - - + + */ - - + + } /* @@ -1862,19 +1880,19 @@ void AliAnalysisTaskGammaConversion::FillGammaElectronInvMass(TString histoMass, TClonesArray const tlVeNeg,TClonesArray const tlVePos) { //see header file for documentation - + for( Int_t iNeg=0; iNeg < tlVeNeg.GetEntriesFast(); iNeg++ ){ - + for (Int_t iPos=0; iPos < tlVePos.GetEntriesFast(); iPos++){ - + TLorentzVector xy = *((TLorentzVector *)(tlVePos.At(iPos))) + *((TLorentzVector *)(tlVeNeg.At(iNeg))); - + for (Int_t iGam=0; iGam < fKFGammas.GetEntriesFast(); iGam++){ - + // AliKFParticle * gammaCandidate = &fKFGammas[iGam]; AliKFParticle * gammaCandidate = (AliKFParticle *)(fKFGammas.At(iGam)); TLorentzVector g; - + g.SetXYZM(gammaCandidate->GetPx(),gammaCandidate->GetPy(),gammaCandidate->GetPz(),fGammaMass); TLorentzVector xyg = xy + g; fHistograms->FillHistogram(histoMass.Data(),xyg.M()); @@ -1882,7 +1900,7 @@ void AliAnalysisTaskGammaConversion::FillGammaElectronInvMass(TString histoMass, } } } - + } void AliAnalysisTaskGammaConversion::ElectronBackground(TString hBg, TClonesArray e) { @@ -1892,7 +1910,7 @@ void AliAnalysisTaskGammaConversion::ElectronBackground(TString hBg, TClonesArra for (Int_t j=i+1; j < e.GetEntriesFast(); j++) { TLorentzVector ee = (*(TLorentzVector*)(e.At(i))) + (*(TLorentzVector*)(e.At(j))); - + fHistograms->FillHistogram(hBg.Data(),ee.M()); } } @@ -1903,41 +1921,41 @@ void AliAnalysisTaskGammaConversion::CleanWithAngleCuts(TClonesArray const negat TClonesArray const positiveElectrons, TClonesArray const gammas){ // see header file for documentation - + UInt_t sizeN = negativeElectrons.GetEntriesFast(); UInt_t sizeP = positiveElectrons.GetEntriesFast(); UInt_t sizeG = gammas.GetEntriesFast(); - - - + + + vector xNegBand(sizeN); vector xPosBand(sizeP); vector gammaBand(sizeG); - - + + for(UInt_t iNeg=0; iNeg < sizeN; iNeg++) xNegBand[iNeg]=kTRUE; for(UInt_t iPos=0; iPos < sizeP; iPos++) xPosBand[iPos]=kTRUE; for(UInt_t iGam=0; iGam < sizeG; iGam++) gammaBand[iGam]=kTRUE; - - - for(UInt_t iPos=0; iPos < sizeP; iPos++){ + + for(UInt_t iPos=0; iPos < sizeP; iPos++){ + Double_t aP[3]; ((AliESDtrack*)(positiveElectrons.At(iPos)))->GetConstrainedPxPyPz(aP); - + TVector3 ePosV(aP[0],aP[1],aP[2]); - + for(UInt_t iNeg=0; iNeg < sizeN; iNeg++){ - + Double_t aN[3]; ((AliESDtrack*)(negativeElectrons.At(iNeg)))->GetConstrainedPxPyPz(aN); TVector3 eNegV(aN[0],aN[1],aN[2]); - + if(ePosV.Angle(eNegV) < 0.05){ //e+e- from gamma xPosBand[iPos]=kFALSE; xNegBand[iNeg]=kFALSE; } - + for(UInt_t iGam=0; iGam < sizeG; iGam++){ AliKFParticle* gammaCandidate = (AliKFParticle*)gammas.At(iGam); TVector3 gammaCandidateVector(gammaCandidate->Px(),gammaCandidate->Py(),gammaCandidate->Pz()); @@ -1946,10 +1964,10 @@ void AliAnalysisTaskGammaConversion::CleanWithAngleCuts(TClonesArray const negat } } } - - - - + + + + for(UInt_t iPos=0; iPos < sizeP; iPos++){ if(xPosBand[iPos]){ new((*fCurrentEventPosElectronTClone)[fCurrentEventPosElectronTClone->GetEntriesFast()]) AliESDtrack((*(AliESDtrack*)(positiveElectrons.At(iPos)))); @@ -1976,33 +1994,33 @@ void AliAnalysisTaskGammaConversion::GetPID(AliESDtrack *track, Stat_t &pid, St // see header file for documentation pid = -1; weight = -1; - + double wpart[5]; double wpartbayes[5]; - + //get probability of the diffenrent particle types track->GetESDpid(wpart); - + // Tentative particle type "concentrations" double c[5]={0.01, 0.01, 0.85, 0.10, 0.05}; - + //Bayes' formula double rcc = 0.; for (int i = 0; i < 5; i++) { rcc+=(c[i] * wpart[i]); } - - - + + + for (int i=0; i<5; i++) { if( rcc!=0){ wpartbayes[i] = c[i] * wpart[i] / rcc; } } - - - + + + Float_t max=0.; int ipid=-1; //find most probable particle in ESD pid @@ -2011,18 +2029,18 @@ void AliAnalysisTaskGammaConversion::GetPID(AliESDtrack *track, Stat_t &pid, St { if (wpartbayes[i] > max) { - ipid = i; - max = wpartbayes[i]; + ipid = i; + max = wpartbayes[i]; } } - + pid = ipid; weight = max; } double AliAnalysisTaskGammaConversion::GetSigmaToVertex(AliESDtrack* t) { // Calculates the number of sigma to the vertex. - + Float_t b[2]; Float_t bRes[2]; Float_t bCov[3]; @@ -2033,7 +2051,7 @@ double AliAnalysisTaskGammaConversion::GetSigmaToVertex(AliESDtrack* t) } bRes[0] = TMath::Sqrt(bCov[0]); bRes[1] = TMath::Sqrt(bCov[2]); - + // ----------------------------------- // How to get to a n-sigma cut? // @@ -2044,18 +2062,18 @@ double AliAnalysisTaskGammaConversion::GetSigmaToVertex(AliESDtrack* t) // // It means that for a 2-dim gauss: n_sigma(d) = Sqrt(2)*ErfInv(1 - Exp((-x**2)/2) // Can this be expressed in a different way? - + if (bRes[0] == 0 || bRes[1] ==0) return -1; - + double d = TMath::Sqrt(TMath::Power(b[0]/bRes[0],2) + TMath::Power(b[1]/bRes[1],2)); - + // stupid rounding problem screws up everything: // if d is too big, TMath::Exp(...) gets 0, and TMath::ErfInverse(1) that should be infinite, gets 0 :( if (TMath::Exp(-d * d / 2) < 1e-10) return 1000; - - + + d = TMath::ErfInverse(1 - TMath::Exp(-d * d / 2)) * TMath::Sqrt(2); return d; } @@ -2065,7 +2083,7 @@ TClonesArray AliAnalysisTaskGammaConversion::GetTLorentzVector(TClonesArray *con //Return TLoresntz vector of track? // vector tlVtrack(0); TClonesArray array("TLorentzVector",0); - + for(Int_t itrack=0; itrack < esdTrack->GetEntriesFast(); itrack++){ double p[3]; //esdTrack[itrack]->GetConstrainedPxPyPz(p); @@ -2075,9 +2093,9 @@ TClonesArray AliAnalysisTaskGammaConversion::GetTLorentzVector(TClonesArray *con new((array)[array.GetEntriesFast()]) TLorentzVector(currentTrack); // tlVtrack.push_back(currentTrack); } - + return array; - + // return tlVtrack; } diff --git a/PWG4/GammaConv/AliAnalysisTaskGammaConversion.h b/PWG4/GammaConv/AliAnalysisTaskGammaConversion.h index 042a1ef52cf..9ff90c11a81 100644 --- a/PWG4/GammaConv/AliAnalysisTaskGammaConversion.h +++ b/PWG4/GammaConv/AliAnalysisTaskGammaConversion.h @@ -1,6 +1,6 @@ #ifndef ALIANALYSISTASKGAMMACONVERSION_H #define ALIANALYSISTASKGAMMACONVERSION_H - + /* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * * See cxx source for full Copyright notice */ @@ -9,7 +9,7 @@ // Class used to do analysis on conversion pairs //--------------------------------------------- //////////////////////////////////////////////// - + #include "AliAnalysisTaskSE.h" #include #include "AliV0Reader.h" @@ -22,18 +22,43 @@ class AliKFParticle; class AliESDInputHandler; class AliESDEvent; class AliAODEvent; +class AliMCEvent; class TList; class AliStack; class AliESDtrackCuts; - +class AliCFManager; // for CF +class AliCFContainer; // for CF class AliAnalysisTaskGammaConversion : public AliAnalysisTaskSE { + + // for CF + enum{ + kStepGenerated = 0, + kStepReconstructable = 1, + kStepLikeSign = 2, + kStepTPCRefit = 3, + kStepKinks = 4, + kStepGetOnFly = 5, + kStepNContributors = 6, + kStepTPCPID = 7, + kStepR = 8, + kStepLine = 9, + kStepZ = 10, + kStepNDF = 11, + kStepChi2 = 12, + kStepEta = 13, + kStepPt = 14 + }; + + + + public: AliAnalysisTaskGammaConversion(); AliAnalysisTaskGammaConversion(const char* name); virtual ~AliAnalysisTaskGammaConversion() ;// virtual destructor - + // Implementation of interface methods virtual void UserCreateOutputObjects(); virtual void Init(); @@ -41,20 +66,25 @@ class AliAnalysisTaskGammaConversion : public AliAnalysisTaskSE virtual void Exec(Option_t *option); virtual void Terminate(Option_t *option); virtual void ConnectInputData(Option_t *); - + void ProcessMCData(); void ProcessV0sNoCut(); void ProcessV0s(); void ProcessGammasForNeutralMesonAnalysis(); - + + // for CF + void SetCFManager(AliCFManager *io) {fCFManager = io;}; + AliCFManager *GetCFManager() const {return fCFManager;} + + // AOD TString GetAODBranchName() const {return fAODBranchName;} void SetAODBranchName(TString name) {fAODBranchName = name ;} void FillAODWithConversionGammas(); // end AOD - - - // for GammaJetAnalysis + + + // for GammaJetAnalysis void ProcessGammasForGammaJetAnalysis(); void CreateListOfChargedParticles(); Double_t GetMinimumDistanceToCharge(Int_t indexHighestPtGamma); @@ -62,19 +92,19 @@ class AliAnalysisTaskGammaConversion : public AliAnalysisTaskSE Int_t GetIndexHighestPtGamma(); void SetESDtrackCuts(); // end of Gamma Jet - + void SetMinPtForGammaJet(Double_t minPtForGammaJet){fMinPtForGammaJet=minPtForGammaJet;} void SetMinIsoConeSize(Double_t minIsoConeSize){fMinIsoConeSize=minIsoConeSize;} void SetMinPtIsoCone(Double_t minPtIsoCone){fMinPtIsoCone=minPtIsoCone;} void SetMinPtGamChargedCorr(Double_t minPtGamChargedCorr){fMinPtGamChargedCorr=minPtGamChargedCorr;} void SetMinPtJetCone(Double_t minPtJetCone){fMinPtJetCone=minPtJetCone;} - + void SetHistograms(AliGammaConversionHistograms *const histograms){fHistograms=histograms;} void SetDoMCTruth(Bool_t flag){fDoMCTruth=flag;} void SetDoNeutralMeson(Bool_t flag){fDoNeutralMeson=flag;} void SetDoJet(Bool_t flag){fDoJet=flag;} void SetDoChic(Bool_t flag){fDoChic=flag;} - + void SetElectronMass(Double_t electronMass){fElectronMass = electronMass;} void SetGammaMass(Double_t gammaMass){fGammaMass = gammaMass;} void SetGammaWidth(Double_t gammaWidth){fGammaWidth = gammaWidth;} @@ -90,8 +120,8 @@ class AliAnalysisTaskGammaConversion : public AliAnalysisTaskSE void FillNtuple(); Double_t GetMCOpeningAngle(TParticle* const daughter0, TParticle* const daughter1) const; void CheckV0Efficiency(); - - + + //////////////////Chi_c Analysis//////////////////////////// void GetPID(AliESDtrack *track, Stat_t &pid, Stat_t &weight); double GetSigmaToVertex(AliESDtrack* t); @@ -103,36 +133,42 @@ class AliAnalysisTaskGammaConversion : public AliAnalysisTaskSE TClonesArray GetTLorentzVector(TClonesArray* esdTrack); void ProcessGammaElectronsForChicAnalysis(); /////////////////////////////////////////////////////////////// - - + + private: AliAnalysisTaskGammaConversion(const AliAnalysisTaskGammaConversion&); // Not implemented AliAnalysisTaskGammaConversion& operator=(const AliAnalysisTaskGammaConversion&); // Not implemented - + AliV0Reader* fV0Reader; // The V0 reader object - + AliStack * fStack; // pointer to the MC particle stack + AliMCEventHandler *fMCTruth; // for CF pointer to MCTruth + AliMCEvent *fMCEvent; // for CF pointer to the MC Event AliESDEvent* fESDEvent; //pointer to the ESDEvent - TList * fOutputContainer ; // Histogram container + TList * fOutputContainer; // Histogram container + AliCFManager *fCFManager; // for CF + // AliCFContainer *container; // for CF + + AliGammaConversionHistograms *fHistograms; // Pointer to the histogram handling class - + Bool_t fDoMCTruth; // Flag to switch on/off MC truth Bool_t fDoNeutralMeson; // flag Bool_t fDoJet; // flag Bool_t fDoChic; // flag - + TClonesArray * fKFReconstructedGammasTClone; //! transient TClonesArray * fCurrentEventPosElectronTClone; //! transient TClonesArray * fCurrentEventNegElectronTClone; //! transient TClonesArray * fKFReconstructedGammasCutTClone; //! transient TClonesArray * fPreviousEventTLVNegElectronTClone; //! transient TClonesArray * fPreviousEventTLVPosElectronTClone; //! transient - + // vector fKFReconstructedGammas; // vector containing all reconstructed gammas vector fElectronv1; // vector containing index of electron 1 vector fElectronv2; // vector containing index of electron 2 - + ///////Chi_c Analysis/////////////////////////// // vector fCurrentEventPosElectron; // comment here // vector fCurrentEventNegElectron; // comment here @@ -140,32 +176,32 @@ class AliAnalysisTaskGammaConversion : public AliAnalysisTaskSE // vector fPreviousEventTLVNegElectron; // comment here // vector fPreviousEventTLVPosElectron; // comment here ////////////////////////////////////////////////// - + //mass defines Double_t fElectronMass; //electron mass Double_t fGammaMass; //gamma mass Double_t fPi0Mass; //pi0mass Double_t fEtaMass; //eta mass - + // width defines Double_t fGammaWidth; //gamma width cut Double_t fPi0Width; // pi0 width cut Double_t fEtaWidth; // eta width cut - + Double_t fMinOpeningAngleGhostCut; // minimum angle cut - + AliESDtrackCuts* fEsdTrackCuts; // Object containing the parameters of the esd track cuts - + Bool_t fCalculateBackground; //flag to set backgrount calculation on/off Bool_t fWriteNtuple; // flag to set if writing to ntuple on/off TNtuple *fGammaNtuple; // Ntuple for gamma values TNtuple *fNeutralMesonNtuple;// NTuple for mesons - + Int_t fTotalNumberOfAddedNtupleEntries; // number of added ntuple entries - + TClonesArray* fChargedParticles; //! transient vector fChargedParticlesId; //! transient - + Double_t fGammaPtHighest; //! transient Double_t fMinPtForGammaJet; //! transient Double_t fMinIsoConeSize; //! transient @@ -173,13 +209,13 @@ class AliAnalysisTaskGammaConversion : public AliAnalysisTaskSE Double_t fMinPtGamChargedCorr; //! transient Double_t fMinPtJetCone; //! transient Int_t fLeadingChargedIndex; //! transient - + TClonesArray* fAODBranch ; //! selected particles branch TString fAODBranchName; // New AOD branch name - + // TClonesArray *fAODObjects; - + ClassDef(AliAnalysisTaskGammaConversion, 4); // Analysis task for gamma conversions }; - + #endif //ALIANALYSISTASKGAMMA_H diff --git a/PWG4/GammaConv/AliV0Reader.cxx b/PWG4/GammaConv/AliV0Reader.cxx index 970e4e6f999..e9cd81b5485 100644 --- a/PWG4/GammaConv/AliV0Reader.cxx +++ b/PWG4/GammaConv/AliV0Reader.cxx @@ -1,664 +1,869 @@ -/************************************************************************** - * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * - * * - * Author: Ana Marin, Kathrin Koch, Kenneth Aamodt * - * Version 1.0 * - * * - * 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. * - **************************************************************************/ - -//////////////////////////////////////////////// -//--------------------------------------------- -// Class used to do analysis on conversion pairs -//--------------------------------------------- -//////////////////////////////////////////////// - -// --- ROOT system --- -#include - -//---- ANALYSIS system ---- -#include "AliV0Reader.h" -#include "AliAnalysisManager.h" -#include "AliESDInputHandler.h" -#include "AliESDtrack.h" -#include "AliMCEvent.h" -#include "AliKFVertex.h" - -#include "AliStack.h" -#include "AliMCEventHandler.h" - - -class iostream; -class AliESDv0; -class TFormula; - -using namespace std; - -ClassImp(AliV0Reader) - - - -AliV0Reader::AliV0Reader() : -TObject(), - fMCStack(NULL), - fMCTruth(NULL), - fChain(NULL), - fESDHandler(NULL), - fESDEvent(NULL), - fHistograms(NULL), - fCurrentV0IndexNumber(0), - fCurrentV0(NULL), - fCurrentNegativeKFParticle(NULL), - fCurrentPositiveKFParticle(NULL), - fCurrentMotherKFCandidate(NULL), - fCurrentNegativeESDTrack(NULL), - fCurrentPositiveESDTrack(NULL), - fNegativeTrackLorentzVector(NULL), - fPositiveTrackLorentzVector(NULL), - fMotherCandidateLorentzVector(NULL), - fCurrentXValue(0), - fCurrentYValue(0), - fCurrentZValue(0), - fPositiveTrackPID(0), - fNegativeTrackPID(0), - fNegativeMCParticle(NULL), - fPositiveMCParticle(NULL), - fMotherMCParticle(NULL), - fMotherCandidateKFMass(0), - fMotherCandidateKFWidth(0), - fUseKFParticle(kTRUE), - fUseESDTrack(kFALSE), - fDoMC(kFALSE), - fMaxR(10000),// 100 meter(outside of ALICE) - fEtaCut(0.), - fPtCut(0.), - fLineCutZRSlope(0.), - fLineCutZValue(0.), - fChi2CutConversion(0.), - fChi2CutMeson(0.), - fPIDProbabilityCutNegativeParticle(0), - fPIDProbabilityCutPositiveParticle(0), - fXVertexCut(0.), - fYVertexCut(0.), - fZVertexCut(0.), - fNSigmaMass(0.), - fUseImprovedVertex(kFALSE), - fCurrentEventGoodV0s(), - fPreviousEventGoodV0s() -{ - -} - - -AliV0Reader::AliV0Reader(const AliV0Reader & original) : - TObject(original), - fMCStack(original.fMCStack), - fMCTruth(original.fMCTruth), - fChain(original.fChain), - fESDHandler(original.fESDHandler), - fESDEvent(original.fESDEvent), - fHistograms(original.fHistograms), - fCurrentV0IndexNumber(original.fCurrentV0IndexNumber), - fCurrentV0(original.fCurrentV0), - fCurrentNegativeKFParticle(original.fCurrentNegativeKFParticle), - fCurrentPositiveKFParticle(original.fCurrentPositiveKFParticle), - fCurrentMotherKFCandidate(original.fCurrentMotherKFCandidate), - fCurrentNegativeESDTrack(original.fCurrentNegativeESDTrack), - fCurrentPositiveESDTrack(original.fCurrentPositiveESDTrack), - fNegativeTrackLorentzVector(original.fNegativeTrackLorentzVector), - fPositiveTrackLorentzVector(original.fPositiveTrackLorentzVector), - fMotherCandidateLorentzVector(original.fMotherCandidateLorentzVector), - fCurrentXValue(original.fCurrentXValue), - fCurrentYValue(original.fCurrentYValue), - fCurrentZValue(original.fCurrentZValue), - fPositiveTrackPID(original.fPositiveTrackPID), - fNegativeTrackPID(original.fNegativeTrackPID), - fNegativeMCParticle(original.fNegativeMCParticle), - fPositiveMCParticle(original.fPositiveMCParticle), - fMotherMCParticle(original.fMotherMCParticle), - fMotherCandidateKFMass(original.fMotherCandidateKFMass), - fMotherCandidateKFWidth(original.fMotherCandidateKFWidth), - fUseKFParticle(kTRUE), - fUseESDTrack(kFALSE), - fDoMC(kFALSE), - fMaxR(original.fMaxR), - fEtaCut(original.fEtaCut), - fPtCut(original.fPtCut), - fLineCutZRSlope(original.fLineCutZRSlope), - fLineCutZValue(original.fLineCutZValue), - fChi2CutConversion(original.fChi2CutConversion), - fChi2CutMeson(original.fChi2CutMeson), - fPIDProbabilityCutNegativeParticle(original.fPIDProbabilityCutNegativeParticle), - fPIDProbabilityCutPositiveParticle(original.fPIDProbabilityCutPositiveParticle), - fXVertexCut(original.fXVertexCut), - fYVertexCut(original.fYVertexCut), - fZVertexCut(original.fZVertexCut), - fNSigmaMass(original.fNSigmaMass), - fUseImprovedVertex(original.fUseImprovedVertex), - fCurrentEventGoodV0s(original.fCurrentEventGoodV0s), - fPreviousEventGoodV0s(original.fPreviousEventGoodV0s) -{ - -} - - -AliV0Reader & AliV0Reader::operator = (const AliV0Reader & /*source*/) -{ - // assignment operator - return *this; -} - -void AliV0Reader::Initialize(){ - //see header file for documentation - - // Get the input handler from the manager - fESDHandler = (AliESDInputHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()); - if(fESDHandler == NULL){ - //print warning here - } - - // Get pointer to esd event from input handler - fESDEvent = fESDHandler->GetEvent(); - if(fESDEvent == NULL){ - //print warning here - } - - //Get pointer to MCTruth - fMCTruth = (AliMCEventHandler*)((AliAnalysisManager::GetAnalysisManager())->GetMCtruthEventHandler()); - if(fMCTruth == NULL){ - //print warning here - } - - //Get pointer to the mc stack - fMCStack = fMCTruth->MCEvent()->Stack(); - if(fMCStack == NULL){ - //print warning here - } - - AliKFParticle::SetField(fESDEvent->GetMagneticField()); - -} - -AliESDv0* AliV0Reader::GetV0(Int_t index){ - //see header file for documentation - - fCurrentV0 = fESDEvent->GetV0(index); - UpdateV0Information(); - return fCurrentV0; -} - -Bool_t AliV0Reader::CheckForPrimaryVertex(){ - return fESDEvent->GetPrimaryVertex()->GetNContributors()>0; -} - - - -Bool_t AliV0Reader::NextV0(){ - //see header file for documentation - - Bool_t iResult=kFALSE; - while(fCurrentV0IndexNumberGetNumberOfV0s()){ - fCurrentV0 = fESDEvent->GetV0(fCurrentV0IndexNumber); - - //checks if on the fly mode is set - if ( !fCurrentV0->GetOnFlyStatus() ){ - if(fHistograms != NULL){ - fHistograms->FillHistogram("ESD_CutGetOnFly_InvMass",GetMotherCandidateMass()); - } - fCurrentV0IndexNumber++; - continue; - } - - //checks if we have a prim vertex - if(fESDEvent->GetPrimaryVertex()->GetNContributors()<=0) { - if(fHistograms != NULL){ - fHistograms->FillHistogram("ESD_CutNContributors_InvMass",GetMotherCandidateMass()); - } - fCurrentV0IndexNumber++; - continue; - } - - //Check the pid probability - if(CheckPIDProbability(fPIDProbabilityCutNegativeParticle,fPIDProbabilityCutPositiveParticle)==kFALSE){ - if(fHistograms != NULL){ - fHistograms->FillHistogram("ESD_CutPIDProb_InvMass",GetMotherCandidateMass()); - } - fCurrentV0IndexNumber++; - continue; - } - - - fCurrentV0->GetXYZ(fCurrentXValue,fCurrentYValue,fCurrentZValue); - - - if(GetXYRadius()>fMaxR){ // cuts on distance from collision point - if(fHistograms != NULL){ - fHistograms->FillHistogram("ESD_CutR_InvMass",GetMotherCandidateMass()); - } - fCurrentV0IndexNumber++; - continue; - } - - - if((TMath::Abs(fCurrentZValue)*fLineCutZRSlope)-fLineCutZValue > GetXYRadius() ){ // cuts out regions where we do not reconstruct - if(fHistograms != NULL){ - fHistograms->FillHistogram("ESD_CutLine_InvMass",GetMotherCandidateMass()); - } - fCurrentV0IndexNumber++; - continue; - } - - - if(UpdateV0Information() == kFALSE){ - fCurrentV0IndexNumber++; - continue; - } - - if(fUseKFParticle){ - if(fCurrentMotherKFCandidate->GetNDF()<=0){ - if(fHistograms != NULL){ - fHistograms->FillHistogram("ESD_CutNDF_InvMass",GetMotherCandidateMass()); - } - fCurrentV0IndexNumber++; - continue; - } - - - Double_t chi2V0 = fCurrentMotherKFCandidate->GetChi2()/fCurrentMotherKFCandidate->GetNDF(); - if(chi2V0 > fChi2CutConversion || chi2V0 <=0){ - if(fHistograms != NULL){ - fHistograms->FillHistogram("ESD_CutChi2_InvMass",GetMotherCandidateMass()); - } - fCurrentV0IndexNumber++; - continue; - } - - - if(TMath::Abs(fMotherCandidateLorentzVector->Eta())> fEtaCut){ - if(fHistograms != NULL){ - fHistograms->FillHistogram("ESD_CutEta_InvMass",GetMotherCandidateMass()); - } - fCurrentV0IndexNumber++; - continue; - } - - - if(fMotherCandidateLorentzVector->Pt()FillHistogram("ESD_CutPt_InvMass",GetMotherCandidateMass()); - } - fCurrentV0IndexNumber++; - continue; - } - - - } - else if(fUseESDTrack){ - //TODO - } - - fCurrentEventGoodV0s.push_back(*fCurrentMotherKFCandidate); - - iResult=kTRUE;//means we have a v0 who survived all the cuts applied - - fCurrentV0IndexNumber++; - - break; - } - return iResult; -} - -Bool_t AliV0Reader::UpdateV0Information(){ - //see header file for documentation - - Bool_t iResult=kTRUE; // for taking out not refitted, kinks and like sign tracks - - Bool_t switchTracks = kFALSE; - - fCurrentNegativeESDTrack = fESDEvent->GetTrack(fCurrentV0->GetNindex()); - fCurrentPositiveESDTrack = fESDEvent->GetTrack(fCurrentV0->GetPindex()); - - if(fCurrentNegativeESDTrack->GetSign() == fCurrentPositiveESDTrack->GetSign()){ // avoid like sign - iResult=kFALSE; - if(fHistograms != NULL){ - fHistograms->FillHistogram("ESD_CutLikeSign_InvMass",GetMotherCandidateMass()); - } - } - - if(fCurrentPositiveESDTrack->GetSign() == -1 && fCurrentNegativeESDTrack->GetSign() == 1){ // switch wrong signed tracks - fCurrentNegativeESDTrack = fESDEvent->GetTrack(fCurrentV0->GetPindex()); - fCurrentPositiveESDTrack = fESDEvent->GetTrack(fCurrentV0->GetNindex()); - switchTracks = kTRUE; - } - - if( !(fCurrentNegativeESDTrack->GetStatus() & AliESDtrack::kTPCrefit) || - !(fCurrentPositiveESDTrack->GetStatus() & AliESDtrack::kTPCrefit) ){ - // if( !(fCurrentNegativeESDTrack->GetStatus() & AliESDtrack::kITSrefit) || - // !(fCurrentPositiveESDTrack->GetStatus() & AliESDtrack::kITSrefit) ){ - - iResult=kFALSE; - if(fHistograms != NULL){ - fHistograms->FillHistogram("ESD_CutRefit_InvMass",GetMotherCandidateMass()); - } - } - - - if( fCurrentNegativeESDTrack->GetKinkIndex(0) > 0 || - fCurrentPositiveESDTrack->GetKinkIndex(0) > 0) { - - iResult=kFALSE; - if(fHistograms != NULL){ - fHistograms->FillHistogram("ESD_CutKink_InvMass",GetMotherCandidateMass()); - } - } - - - - if(fCurrentNegativeKFParticle != NULL){ - delete fCurrentNegativeKFParticle; - } - if(switchTracks == kFALSE){ - fCurrentNegativeKFParticle = new AliKFParticle(*(fCurrentV0->GetParamN()),fNegativeTrackPID); - } - else{ - fCurrentNegativeKFParticle = new AliKFParticle(*(fCurrentV0->GetParamP()),fNegativeTrackPID); - } - - if(fCurrentPositiveKFParticle != NULL){ - delete fCurrentPositiveKFParticle; - } - if(switchTracks == kFALSE){ - fCurrentPositiveKFParticle = new AliKFParticle(*(fCurrentV0->GetParamP()),fPositiveTrackPID); - } - else{ - fCurrentPositiveKFParticle = new AliKFParticle(*(fCurrentV0->GetParamN()),fPositiveTrackPID); - } - - if(fCurrentMotherKFCandidate != NULL){ - delete fCurrentMotherKFCandidate; - } - fCurrentMotherKFCandidate = new AliKFParticle(*fCurrentNegativeKFParticle,*fCurrentPositiveKFParticle); - - - if(fPositiveTrackPID==-11 && fNegativeTrackPID==11){ - fCurrentMotherKFCandidate->SetMassConstraint(0,fNSigmaMass); - } - - - - - if(fUseImprovedVertex == kTRUE){ - AliKFVertex primaryVertexImproved(*GetPrimaryVertex()); - primaryVertexImproved+=*fCurrentMotherKFCandidate; - fCurrentMotherKFCandidate->SetProductionVertex(primaryVertexImproved); - } - - fCurrentMotherKFCandidate->GetMass(fMotherCandidateKFMass,fMotherCandidateKFWidth); - - - if(fNegativeTrackLorentzVector != NULL){ - delete fNegativeTrackLorentzVector; - } - if(fUseKFParticle){ - fNegativeTrackLorentzVector = new TLorentzVector(fCurrentNegativeKFParticle->Px(),fCurrentNegativeKFParticle->Py(),fCurrentNegativeKFParticle->Pz()); - } - else if(fUseESDTrack){ - fNegativeTrackLorentzVector = new TLorentzVector(fCurrentNegativeESDTrack->Px(),fCurrentNegativeESDTrack->Py(),fCurrentNegativeESDTrack->Pz()); - } - - if(fPositiveTrackLorentzVector != NULL){ - delete fPositiveTrackLorentzVector; - } - if(fUseKFParticle){ - fPositiveTrackLorentzVector = new TLorentzVector(fCurrentPositiveKFParticle->Px(),fCurrentPositiveKFParticle->Py(),fCurrentPositiveKFParticle->Pz()); - } - else if(fUseESDTrack){ - fPositiveTrackLorentzVector = new TLorentzVector(fCurrentPositiveESDTrack->Px(),fCurrentPositiveESDTrack->Py(),fCurrentPositiveESDTrack->Pz()); - } - - if(fMotherCandidateLorentzVector != NULL){ - delete fMotherCandidateLorentzVector; - } - if(fUseKFParticle){ - fMotherCandidateLorentzVector = new TLorentzVector(*fNegativeTrackLorentzVector + *fPositiveTrackLorentzVector); - } - else if(fUseESDTrack){ - fMotherCandidateLorentzVector = new TLorentzVector(*fNegativeTrackLorentzVector + *fPositiveTrackLorentzVector); - } - - if(fPositiveTrackPID==-11 && fNegativeTrackPID==11){ - fMotherCandidateLorentzVector->SetXYZM(fMotherCandidateLorentzVector->Px() ,fMotherCandidateLorentzVector->Py(),fMotherCandidateLorentzVector->Pz(),0.); - } - - - if(fDoMC == kTRUE){ - fMotherMCParticle= NULL; - fNegativeMCParticle = fMCStack->Particle(TMath::Abs(fESDEvent->GetTrack(fCurrentV0->GetNindex())->GetLabel())); - fPositiveMCParticle = fMCStack->Particle(TMath::Abs(fESDEvent->GetTrack(fCurrentV0->GetPindex())->GetLabel())); - if(fPositiveMCParticle->GetMother(0)>-1){ - fMotherMCParticle = fMCStack->Particle(fPositiveMCParticle->GetMother(0)); - } - } - - // if(iResult==kTRUE){ - // fCurrentEventGoodV0s.push_back(*fCurrentMotherKFCandidate); // moved it to NextV0() after all the cuts are applied - // } - - return iResult; -} - - - -Bool_t AliV0Reader::HasSameMCMother(){ - //see header file for documentation - - Bool_t iResult = kFALSE; - if(fDoMC == kTRUE){ - if(fNegativeMCParticle != NULL && fPositiveMCParticle != NULL){ - if(fNegativeMCParticle->GetMother(0) == fPositiveMCParticle->GetMother(0)) - if(fMotherMCParticle){ - iResult = kTRUE; - } - } - } - return iResult; -} - -Bool_t AliV0Reader::CheckPIDProbability(Double_t negProbCut, Double_t posProbCut){ - //see header file for documentation - - Bool_t iResult=kFALSE; - - Double_t *posProbArray = new Double_t[10]; - Double_t *negProbArray = new Double_t[10]; - AliESDtrack* negTrack = fESDEvent->GetTrack(fCurrentV0->GetNindex()); - AliESDtrack* posTrack = fESDEvent->GetTrack(fCurrentV0->GetPindex()); - - negTrack->GetTPCpid(negProbArray); - posTrack->GetTPCpid(posProbArray); - - if(negProbArray!=NULL && posProbArray!=NULL){ - if(negProbArray[GetSpeciesIndex(-1)]>=negProbCut && posProbArray[GetSpeciesIndex(1)]>=posProbCut){ - iResult=kTRUE; - } - } - delete [] posProbArray; - delete [] negProbArray; - return iResult; -} - -void AliV0Reader::GetPIDProbability(Double_t &negPIDProb,Double_t & posPIDProb){ - - Double_t *posProbArray = new Double_t[10]; - Double_t *negProbArray = new Double_t[10]; - AliESDtrack* negTrack = fESDEvent->GetTrack(fCurrentV0->GetNindex()); - AliESDtrack* posTrack = fESDEvent->GetTrack(fCurrentV0->GetPindex()); - - negTrack->GetTPCpid(negProbArray); - posTrack->GetTPCpid(posProbArray); - - if(negProbArray!=NULL && posProbArray!=NULL){ - negPIDProb = negProbArray[GetSpeciesIndex(-1)]; - posPIDProb = posProbArray[GetSpeciesIndex(1)]; - } - delete [] posProbArray; - delete [] negProbArray; -} - -void AliV0Reader::UpdateEventByEventData(){ - //see header file for documentation - - if(fCurrentEventGoodV0s.size() >0 ){ - // fPreviousEventGoodV0s.clear(); - // fPreviousEventGoodV0s = fCurrentEventGoodV0s; - if(fPreviousEventGoodV0s.size()>19){ - for(UInt_t nCurrent=0;nCurrentPhi()> TMath::Pi()){ - offset = -2*TMath::Pi(); - } - return fNegativeTrackLorentzVector->Phi()+offset; -} - -Double_t AliV0Reader::GetPositiveTrackPhi() const{ - //see header file for documentation - - Double_t offset=0; - if(fPositiveTrackLorentzVector->Phi()> TMath::Pi()){ - offset = -2*TMath::Pi(); - } - return fPositiveTrackLorentzVector->Phi()+offset; -} - -Double_t AliV0Reader::GetMotherCandidatePhi() const{ - //see header file for documentation - - Double_t offset=0; - if(fMotherCandidateLorentzVector->Phi()> TMath::Pi()){ - offset = -2*TMath::Pi(); - } - return fMotherCandidateLorentzVector->Phi()+offset; -} - - -Double_t AliV0Reader::GetMotherCandidateRapidity() const{ - //see header file for documentation - - Double_t rapidity=0; - if(fMotherCandidateLorentzVector->Energy() - fMotherCandidateLorentzVector->Pz() == 0 || fMotherCandidateLorentzVector->Energy() + fMotherCandidateLorentzVector->Pz() == 0) rapidity=0; - else rapidity = 0.5*(TMath::Log((fMotherCandidateLorentzVector->Energy() + fMotherCandidateLorentzVector->Pz()) / (fMotherCandidateLorentzVector->Energy()-fMotherCandidateLorentzVector->Pz()))); - return rapidity; - -} - - - - - -Int_t AliV0Reader::GetSpeciesIndex(Int_t chargeOfTrack){ - //see header file for documentation - - Int_t iResult = 10; // Unknown particle - - if(chargeOfTrack==-1){ //negative track - switch(abs(fNegativeTrackPID)){ - case 11: //electron - iResult = 0; - break; - case 13: //muon - iResult = 1; - break; - case 211: //pion - iResult = 2; - break; - case 321: //kaon - iResult = 3; - break; - case 2212: //proton - iResult = 4; - break; - case 22: //photon - iResult = 5; - break; - case 111: //pi0 - iResult = 6; - break; - case 2112: //neutron - iResult = 7; - break; - case 311: //K0 - iResult = 8; - break; - - //Put in here for kSPECIES::kEleCon ???? - } - } - else if(chargeOfTrack==1){ //positive track - switch(abs(fPositiveTrackPID)){ - case 11: //electron - iResult = 0; - break; - case 13: //muon - iResult = 1; - break; - case 211: //pion - iResult = 2; - break; - case 321: //kaon - iResult = 3; - break; - case 2212: //proton - iResult = 4; - break; - case 22: //photon - iResult = 5; - break; - case 111: //pi0 - iResult = 6; - break; - case 2112: //neutron - iResult = 7; - break; - case 311: //K0 - iResult = 8; - break; - - //Put in here for kSPECIES::kEleCon ???? - } - } - else{ - //Wrong parameter.. Print warning - } - return iResult; -} +/************************************************************************** + * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * + * * + * Author: Ana Marin, Kathrin Koch, Kenneth Aamodt * + * Version 1.0 * + * * + * 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. * + **************************************************************************/ + +//////////////////////////////////////////////// +//--------------------------------------------- +// Class used to do analysis on conversion pairs +//--------------------------------------------- +//////////////////////////////////////////////// + +// --- ROOT system --- +#include + +//---- ANALYSIS system ---- +#include "AliV0Reader.h" +#include "AliAnalysisManager.h" +#include "AliESDInputHandler.h" +#include "AliESDtrack.h" +#include "AliMCEvent.h" +#include "AliKFVertex.h" + +#include "AliStack.h" +#include "AliMCEventHandler.h" + + +class iostream; +class AliESDv0; +class TFormula; + +using namespace std; + +ClassImp(AliV0Reader) + + + +AliV0Reader::AliV0Reader() : + TObject(), + fMCStack(NULL), + fMCTruth(NULL), + fMCEvent(NULL), // for CF + fChain(NULL), + fESDHandler(NULL), + fESDEvent(NULL), + fHistograms(NULL), + fCurrentV0IndexNumber(0), + fCurrentV0(NULL), + fCurrentNegativeKFParticle(NULL), + fCurrentPositiveKFParticle(NULL), + fCurrentMotherKFCandidate(NULL), + fCurrentNegativeESDTrack(NULL), + fCurrentPositiveESDTrack(NULL), + fNegativeTrackLorentzVector(NULL), + fPositiveTrackLorentzVector(NULL), + fMotherCandidateLorentzVector(NULL), + fCurrentXValue(0), + fCurrentYValue(0), + fCurrentZValue(0), + fPositiveTrackPID(0), + fNegativeTrackPID(0), + fNegativeMCParticle(NULL), + fPositiveMCParticle(NULL), + fMotherMCParticle(NULL), + fMotherCandidateKFMass(0), + fMotherCandidateKFWidth(0), + fUseKFParticle(kTRUE), + fUseESDTrack(kFALSE), + fDoMC(kFALSE), + fMaxR(10000),// 100 meter(outside of ALICE) + fEtaCut(0.), + fPtCut(0.), + fMaxZ(0.), + fLineCutZRSlope(0.), + fLineCutZValue(0.), + fChi2CutConversion(0.), + fChi2CutMeson(0.), + fPIDProbabilityCutNegativeParticle(0), + fPIDProbabilityCutPositiveParticle(0), + fXVertexCut(0.), + fYVertexCut(0.), + fZVertexCut(0.), + fNSigmaMass(0.), + fUseImprovedVertex(kFALSE), + fUseOwnXYZCalculation(kFALSE), + fCurrentEventGoodV0s(), + fPreviousEventGoodV0s() +{ + +} + + +AliV0Reader::AliV0Reader(const AliV0Reader & original) : + TObject(original), + fMCStack(original.fMCStack), + fMCTruth(original.fMCTruth), + fMCEvent(original.fMCEvent), // for CF + fChain(original.fChain), + fESDHandler(original.fESDHandler), + fESDEvent(original.fESDEvent), + fHistograms(original.fHistograms), + fCurrentV0IndexNumber(original.fCurrentV0IndexNumber), + fCurrentV0(original.fCurrentV0), + fCurrentNegativeKFParticle(original.fCurrentNegativeKFParticle), + fCurrentPositiveKFParticle(original.fCurrentPositiveKFParticle), + fCurrentMotherKFCandidate(original.fCurrentMotherKFCandidate), + fCurrentNegativeESDTrack(original.fCurrentNegativeESDTrack), + fCurrentPositiveESDTrack(original.fCurrentPositiveESDTrack), + fNegativeTrackLorentzVector(original.fNegativeTrackLorentzVector), + fPositiveTrackLorentzVector(original.fPositiveTrackLorentzVector), + fMotherCandidateLorentzVector(original.fMotherCandidateLorentzVector), + fCurrentXValue(original.fCurrentXValue), + fCurrentYValue(original.fCurrentYValue), + fCurrentZValue(original.fCurrentZValue), + fPositiveTrackPID(original.fPositiveTrackPID), + fNegativeTrackPID(original.fNegativeTrackPID), + fNegativeMCParticle(original.fNegativeMCParticle), + fPositiveMCParticle(original.fPositiveMCParticle), + fMotherMCParticle(original.fMotherMCParticle), + fMotherCandidateKFMass(original.fMotherCandidateKFMass), + fMotherCandidateKFWidth(original.fMotherCandidateKFWidth), + fUseKFParticle(kTRUE), + fUseESDTrack(kFALSE), + fDoMC(kFALSE), + fMaxR(original.fMaxR), + fEtaCut(original.fEtaCut), + fPtCut(original.fPtCut), + fMaxZ(original.fMaxZ), + fLineCutZRSlope(original.fLineCutZRSlope), + fLineCutZValue(original.fLineCutZValue), + fChi2CutConversion(original.fChi2CutConversion), + fChi2CutMeson(original.fChi2CutMeson), + fPIDProbabilityCutNegativeParticle(original.fPIDProbabilityCutNegativeParticle), + fPIDProbabilityCutPositiveParticle(original.fPIDProbabilityCutPositiveParticle), + fXVertexCut(original.fXVertexCut), + fYVertexCut(original.fYVertexCut), + fZVertexCut(original.fZVertexCut), + fNSigmaMass(original.fNSigmaMass), + fUseImprovedVertex(original.fUseImprovedVertex), + fUseOwnXYZCalculation(original.fUseOwnXYZCalculation), + fCurrentEventGoodV0s(original.fCurrentEventGoodV0s), + fPreviousEventGoodV0s(original.fPreviousEventGoodV0s) +{ + +} + + +AliV0Reader & AliV0Reader::operator = (const AliV0Reader & /*source*/) +{ + // assignment operator + return *this; +} + +void AliV0Reader::Initialize(){ + //see header file for documentation + + // Get the input handler from the manager + fESDHandler = (AliESDInputHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()); + if(fESDHandler == NULL){ + //print warning here + } + + // Get pointer to esd event from input handler + fESDEvent = fESDHandler->GetEvent(); + if(fESDEvent == NULL){ + //print warning here + } + + //Get pointer to MCTruth + fMCTruth = (AliMCEventHandler*)((AliAnalysisManager::GetAnalysisManager())->GetMCtruthEventHandler()); + if(fMCTruth == NULL){ + //print warning here + } + + //Get pointer to the mc stack + fMCStack = fMCTruth->MCEvent()->Stack(); + if(fMCStack == NULL){ + //print warning here + } + + + // for CF + //Get pointer to the mc event + fMCEvent = fMCTruth->MCEvent(); + if(fMCEvent == NULL){ + //print warning here + } + + + AliKFParticle::SetField(fESDEvent->GetMagneticField()); + +} + +AliESDv0* AliV0Reader::GetV0(Int_t index){ + //see header file for documentation + fCurrentV0 = fESDEvent->GetV0(index); + UpdateV0Information(); + return fCurrentV0; +} + +Bool_t AliV0Reader::CheckForPrimaryVertex(){ + return fESDEvent->GetPrimaryVertex()->GetNContributors()>0; +} + + + +Bool_t AliV0Reader::NextV0(){ + //see header file for documentation + + Bool_t iResult=kFALSE; + while(fCurrentV0IndexNumberGetNumberOfV0s()){ + fCurrentV0 = fESDEvent->GetV0(fCurrentV0IndexNumber); + + // moved it up here so that the correction framework can access pt and eta information + if(UpdateV0Information() == kFALSE){ + fCurrentV0IndexNumber++; + continue; + } + + Double_t containerInput[3]; + containerInput[0] = GetMotherCandidatePt(); + containerInput[1] = GetMotherCandidateEta(); + containerInput[2] = GetMotherCandidateMass(); + + + //checks if on the fly mode is set + if ( !fCurrentV0->GetOnFlyStatus() ){ + if(fHistograms != NULL){ + fHistograms->FillHistogram("ESD_CutGetOnFly_InvMass",GetMotherCandidateMass()); + } + fCurrentV0IndexNumber++; + continue; + } + fCFManager->GetParticleContainer()->Fill(containerInput,kStepGetOnFly); // for CF + + //checks if we have a prim vertex + if(fESDEvent->GetPrimaryVertex()->GetNContributors()<=0) { + if(fHistograms != NULL){ + fHistograms->FillHistogram("ESD_CutNContributors_InvMass",GetMotherCandidateMass()); + } + fCurrentV0IndexNumber++; + continue; + } + fCFManager->GetParticleContainer()->Fill(containerInput,kStepNContributors); // for CF + + + //Check the pid probability + if(CheckPIDProbability(fPIDProbabilityCutNegativeParticle,fPIDProbabilityCutPositiveParticle)==kFALSE){ + if(fHistograms != NULL){ + fHistograms->FillHistogram("ESD_CutPIDProb_InvMass",GetMotherCandidateMass()); + } + fCurrentV0IndexNumber++; + continue; + } + fCFManager->GetParticleContainer()->Fill(containerInput,kStepTPCPID); // for CF + + + + + fCurrentV0->GetXYZ(fCurrentXValue,fCurrentYValue,fCurrentZValue); + + + if(GetXYRadius()>fMaxR){ // cuts on distance from collision point + if(fHistograms != NULL){ + fHistograms->FillHistogram("ESD_CutR_InvMass",GetMotherCandidateMass()); + } + fCurrentV0IndexNumber++; + continue; + } + fCFManager->GetParticleContainer()->Fill(containerInput,kStepR); // for CF + + + + if((TMath::Abs(fCurrentZValue)*fLineCutZRSlope)-fLineCutZValue > GetXYRadius() ){ // cuts out regions where we do not reconstruct + if(fHistograms != NULL){ + fHistograms->FillHistogram("ESD_CutLine_InvMass",GetMotherCandidateMass()); + } + fCurrentV0IndexNumber++; + continue; + } + fCFManager->GetParticleContainer()->Fill(containerInput,kStepLine); // for CF + + + if(TMath::Abs(fCurrentZValue) > fMaxZ ){ // cuts out regions where we do not reconstruct + if(fHistograms != NULL){ + fHistograms->FillHistogram("ESD_CutZ_InvMass",GetMotherCandidateMass()); + } + fCurrentV0IndexNumber++; + continue; + } + fCFManager->GetParticleContainer()->Fill(containerInput,kStepZ); // for CF + + + /* Moved further up so corr framework can work + if(UpdateV0Information() == kFALSE){ + fCurrentV0IndexNumber++; + continue; + } + */ + + if(fUseKFParticle){ + if(fCurrentMotherKFCandidate->GetNDF()<=0){ + if(fHistograms != NULL){ + fHistograms->FillHistogram("ESD_CutNDF_InvMass",GetMotherCandidateMass()); + } + fCurrentV0IndexNumber++; + continue; + } + fCFManager->GetParticleContainer()->Fill(containerInput,kStepNDF); // for CF + + + Double_t chi2V0 = fCurrentMotherKFCandidate->GetChi2()/fCurrentMotherKFCandidate->GetNDF(); + if(chi2V0 > fChi2CutConversion || chi2V0 <=0){ + if(fHistograms != NULL){ + fHistograms->FillHistogram("ESD_CutChi2_InvMass",GetMotherCandidateMass()); + } + fCurrentV0IndexNumber++; + continue; + } + fCFManager->GetParticleContainer()->Fill(containerInput,kStepChi2); // for CF + + + if(TMath::Abs(fMotherCandidateLorentzVector->Eta())> fEtaCut){ + if(fHistograms != NULL){ + fHistograms->FillHistogram("ESD_CutEta_InvMass",GetMotherCandidateMass()); + } + fCurrentV0IndexNumber++; + continue; + } + fCFManager->GetParticleContainer()->Fill(containerInput,kStepEta); // for CF + + + if(fMotherCandidateLorentzVector->Pt()FillHistogram("ESD_CutPt_InvMass",GetMotherCandidateMass()); + } + fCurrentV0IndexNumber++; + continue; + } + fCFManager->GetParticleContainer()->Fill(containerInput,kStepPt); // for CF + + + } + else if(fUseESDTrack){ + //TODO + } + + fCurrentEventGoodV0s.push_back(*fCurrentMotherKFCandidate); + + iResult=kTRUE;//means we have a v0 who survived all the cuts applied + + fCurrentV0IndexNumber++; + + break; + } + return iResult; +} + +Bool_t AliV0Reader::UpdateV0Information(){ + //see header file for documentation + + Bool_t iResult=kTRUE; // for taking out not refitted, kinks and like sign tracks + + Bool_t switchTracks = kFALSE; + + fCurrentNegativeESDTrack = fESDEvent->GetTrack(fCurrentV0->GetNindex()); + fCurrentPositiveESDTrack = fESDEvent->GetTrack(fCurrentV0->GetPindex()); + + if(fCurrentNegativeESDTrack->GetSign() == fCurrentPositiveESDTrack->GetSign()){ // avoid like sign + iResult=kFALSE; + if(fHistograms != NULL){ + fHistograms->FillHistogram("ESD_CutLikeSign_InvMass",GetMotherCandidateMass()); + } + } + + if(fCurrentPositiveESDTrack->GetSign() == -1 && fCurrentNegativeESDTrack->GetSign() == 1){ // switch wrong signed tracks + fCurrentNegativeESDTrack = fESDEvent->GetTrack(fCurrentV0->GetPindex()); + fCurrentPositiveESDTrack = fESDEvent->GetTrack(fCurrentV0->GetNindex()); + switchTracks = kTRUE; + } + + if( !(fCurrentNegativeESDTrack->GetStatus() & AliESDtrack::kTPCrefit) || + !(fCurrentPositiveESDTrack->GetStatus() & AliESDtrack::kTPCrefit) ){ + // if( !(fCurrentNegativeESDTrack->GetStatus() & AliESDtrack::kITSrefit) || + // !(fCurrentPositiveESDTrack->GetStatus() & AliESDtrack::kITSrefit) ){ + + iResult=kFALSE; + if(fHistograms != NULL){ + fHistograms->FillHistogram("ESD_CutRefit_InvMass",GetMotherCandidateMass()); + } + } + + if( fCurrentNegativeESDTrack->GetKinkIndex(0) > 0 || + fCurrentPositiveESDTrack->GetKinkIndex(0) > 0) { + + iResult=kFALSE; + if(fHistograms != NULL){ + fHistograms->FillHistogram("ESD_CutKink_InvMass",GetMotherCandidateMass()); + } + } + + if(fCurrentNegativeKFParticle != NULL){ + delete fCurrentNegativeKFParticle; + } + if(switchTracks == kFALSE){ + fCurrentNegativeKFParticle = new AliKFParticle(*(fCurrentV0->GetParamN()),fNegativeTrackPID); + } + else{ + fCurrentNegativeKFParticle = new AliKFParticle(*(fCurrentV0->GetParamP()),fNegativeTrackPID); + } + + if(fCurrentPositiveKFParticle != NULL){ + delete fCurrentPositiveKFParticle; + } + if(switchTracks == kFALSE){ + fCurrentPositiveKFParticle = new AliKFParticle(*(fCurrentV0->GetParamP()),fPositiveTrackPID); + } + else{ + fCurrentPositiveKFParticle = new AliKFParticle(*(fCurrentV0->GetParamN()),fPositiveTrackPID); + } + + if(fCurrentMotherKFCandidate != NULL){ + delete fCurrentMotherKFCandidate; + } + fCurrentMotherKFCandidate = new AliKFParticle(*fCurrentNegativeKFParticle,*fCurrentPositiveKFParticle); + + + if(fPositiveTrackPID==-11 && fNegativeTrackPID==11){ + fCurrentMotherKFCandidate->SetMassConstraint(0,fNSigmaMass); + } + + if(fUseImprovedVertex == kTRUE){ + AliKFVertex primaryVertexImproved(*GetPrimaryVertex()); + primaryVertexImproved+=*fCurrentMotherKFCandidate; + fCurrentMotherKFCandidate->SetProductionVertex(primaryVertexImproved); + } + + fCurrentMotherKFCandidate->GetMass(fMotherCandidateKFMass,fMotherCandidateKFWidth); + + + if(fNegativeTrackLorentzVector != NULL){ + delete fNegativeTrackLorentzVector; + } + if(fUseKFParticle){ + fNegativeTrackLorentzVector = new TLorentzVector(fCurrentNegativeKFParticle->Px(),fCurrentNegativeKFParticle->Py(),fCurrentNegativeKFParticle->Pz()); + } + else if(fUseESDTrack){ + fNegativeTrackLorentzVector = new TLorentzVector(fCurrentNegativeESDTrack->Px(),fCurrentNegativeESDTrack->Py(),fCurrentNegativeESDTrack->Pz()); + } + + if(fPositiveTrackLorentzVector != NULL){ + delete fPositiveTrackLorentzVector; + } + if(fUseKFParticle){ + fPositiveTrackLorentzVector = new TLorentzVector(fCurrentPositiveKFParticle->Px(),fCurrentPositiveKFParticle->Py(),fCurrentPositiveKFParticle->Pz()); + } + else if(fUseESDTrack){ + fPositiveTrackLorentzVector = new TLorentzVector(fCurrentPositiveESDTrack->Px(),fCurrentPositiveESDTrack->Py(),fCurrentPositiveESDTrack->Pz()); + } + + if(fMotherCandidateLorentzVector != NULL){ + delete fMotherCandidateLorentzVector; + } + if(fUseKFParticle){ + fMotherCandidateLorentzVector = new TLorentzVector(*fNegativeTrackLorentzVector + *fPositiveTrackLorentzVector); + } + else if(fUseESDTrack){ + fMotherCandidateLorentzVector = new TLorentzVector(*fNegativeTrackLorentzVector + *fPositiveTrackLorentzVector); + } + + if(fPositiveTrackPID==-11 && fNegativeTrackPID==11){ + fMotherCandidateLorentzVector->SetXYZM(fMotherCandidateLorentzVector->Px() ,fMotherCandidateLorentzVector->Py(),fMotherCandidateLorentzVector->Pz(),0.); + } + + + if(fDoMC == kTRUE){ + fMotherMCParticle= NULL; + fNegativeMCParticle = fMCStack->Particle(TMath::Abs(fESDEvent->GetTrack(fCurrentV0->GetNindex())->GetLabel())); + fPositiveMCParticle = fMCStack->Particle(TMath::Abs(fESDEvent->GetTrack(fCurrentV0->GetPindex())->GetLabel())); + if(fPositiveMCParticle->GetMother(0)>-1){ + fMotherMCParticle = fMCStack->Particle(fPositiveMCParticle->GetMother(0)); + } + } + + // if(iResult==kTRUE){ + // fCurrentEventGoodV0s.push_back(*fCurrentMotherKFCandidate); // moved it to NextV0() after all the cuts are applied + // } + + + // for CF + Double_t containerInput[3]; + containerInput[0] = GetMotherCandidatePt(); + containerInput[1] = GetMotherCandidateEta(); + containerInput[2] = GetMotherCandidateMass(); + + fCFManager->GetParticleContainer()->Fill(containerInput,kStepLikeSign); // for CF + fCFManager->GetParticleContainer()->Fill(containerInput,kStepTPCRefit); // for CF + fCFManager->GetParticleContainer()->Fill(containerInput,kStepKinks); // for CF + + return iResult; +} + + + +Bool_t AliV0Reader::HasSameMCMother(){ + //see header file for documentation + + Bool_t iResult = kFALSE; + if(fDoMC == kTRUE){ + if(fNegativeMCParticle != NULL && fPositiveMCParticle != NULL){ + if(fNegativeMCParticle->GetMother(0) == fPositiveMCParticle->GetMother(0)) + if(fMotherMCParticle){ + iResult = kTRUE; + } + } + } + return iResult; +} + +Bool_t AliV0Reader::CheckPIDProbability(Double_t negProbCut, Double_t posProbCut){ + //see header file for documentation + + Bool_t iResult=kFALSE; + + Double_t *posProbArray = new Double_t[10]; + Double_t *negProbArray = new Double_t[10]; + AliESDtrack* negTrack = fESDEvent->GetTrack(fCurrentV0->GetNindex()); + AliESDtrack* posTrack = fESDEvent->GetTrack(fCurrentV0->GetPindex()); + + negTrack->GetTPCpid(negProbArray); + posTrack->GetTPCpid(posProbArray); + + if(negProbArray!=NULL && posProbArray!=NULL){ + if(negProbArray[GetSpeciesIndex(-1)]>=negProbCut && posProbArray[GetSpeciesIndex(1)]>=posProbCut){ + iResult=kTRUE; + } + } + delete [] posProbArray; + delete [] negProbArray; + return iResult; +} + +void AliV0Reader::GetPIDProbability(Double_t &negPIDProb,Double_t & posPIDProb){ + + Double_t *posProbArray = new Double_t[10]; + Double_t *negProbArray = new Double_t[10]; + AliESDtrack* negTrack = fESDEvent->GetTrack(fCurrentV0->GetNindex()); + AliESDtrack* posTrack = fESDEvent->GetTrack(fCurrentV0->GetPindex()); + + negTrack->GetTPCpid(negProbArray); + posTrack->GetTPCpid(posProbArray); + + if(negProbArray!=NULL && posProbArray!=NULL){ + negPIDProb = negProbArray[GetSpeciesIndex(-1)]; + posPIDProb = posProbArray[GetSpeciesIndex(1)]; + } + delete [] posProbArray; + delete [] negProbArray; +} + +void AliV0Reader::UpdateEventByEventData(){ + //see header file for documentation + + if(fCurrentEventGoodV0s.size() >0 ){ + // fPreviousEventGoodV0s.clear(); + // fPreviousEventGoodV0s = fCurrentEventGoodV0s; + if(fPreviousEventGoodV0s.size()>19){ + for(UInt_t nCurrent=0;nCurrentPhi()> TMath::Pi()){ + offset = -2*TMath::Pi(); + } + return fNegativeTrackLorentzVector->Phi()+offset; +} + +Double_t AliV0Reader::GetPositiveTrackPhi() const{ + //see header file for documentation + + Double_t offset=0; + if(fPositiveTrackLorentzVector->Phi()> TMath::Pi()){ + offset = -2*TMath::Pi(); + } + return fPositiveTrackLorentzVector->Phi()+offset; +} + +Double_t AliV0Reader::GetMotherCandidatePhi() const{ + //see header file for documentation + + Double_t offset=0; + if(fMotherCandidateLorentzVector->Phi()> TMath::Pi()){ + offset = -2*TMath::Pi(); + } + return fMotherCandidateLorentzVector->Phi()+offset; +} + + +Double_t AliV0Reader::GetMotherCandidateRapidity() const{ + //see header file for documentation + + Double_t rapidity=0; + if(fMotherCandidateLorentzVector->Energy() - fMotherCandidateLorentzVector->Pz() == 0 || fMotherCandidateLorentzVector->Energy() + fMotherCandidateLorentzVector->Pz() == 0) rapidity=0; + else rapidity = 0.5*(TMath::Log((fMotherCandidateLorentzVector->Energy() + fMotherCandidateLorentzVector->Pz()) / (fMotherCandidateLorentzVector->Energy()-fMotherCandidateLorentzVector->Pz()))); + return rapidity; + +} + + + + + +Int_t AliV0Reader::GetSpeciesIndex(Int_t chargeOfTrack){ + //see header file for documentation + + Int_t iResult = 10; // Unknown particle + + if(chargeOfTrack==-1){ //negative track + switch(abs(fNegativeTrackPID)){ + case 11: //electron + iResult = 0; + break; + case 13: //muon + iResult = 1; + break; + case 211: //pion + iResult = 2; + break; + case 321: //kaon + iResult = 3; + break; + case 2212: //proton + iResult = 4; + break; + case 22: //photon + iResult = 5; + break; + case 111: //pi0 + iResult = 6; + break; + case 2112: //neutron + iResult = 7; + break; + case 311: //K0 + iResult = 8; + break; + + //Put in here for kSPECIES::kEleCon ???? + } + } + else if(chargeOfTrack==1){ //positive track + switch(abs(fPositiveTrackPID)){ + case 11: //electron + iResult = 0; + break; + case 13: //muon + iResult = 1; + break; + case 211: //pion + iResult = 2; + break; + case 321: //kaon + iResult = 3; + break; + case 2212: //proton + iResult = 4; + break; + case 22: //photon + iResult = 5; + break; + case 111: //pi0 + iResult = 6; + break; + case 2112: //neutron + iResult = 7; + break; + case 311: //K0 + iResult = 8; + break; + + //Put in here for kSPECIES::kEleCon ???? + } + } + else{ + //Wrong parameter.. Print warning + } + return iResult; +} + +Bool_t GetHelixCenter(AliESDtrack* track, Double_t b,Int_t charge, Double_t center[2]){ + // see header file for documentation + + Double_t pi = 3.14159265358979323846; + + Double_t helix[6]; + track->GetHelixParameters(helix,b); + + Double_t xpos = helix[5]; + Double_t ypos = helix[0]; + Double_t radius = TMath::Abs(1./helix[4]); + Double_t phi = helix[2]; + + if(phi < 0){ + phi = phi + 2*pi; + } + + phi -= pi/2.; + Double_t xpoint = radius * TMath::Cos(phi); + Double_t ypoint = radius * TMath::Sin(phi); + + if(charge > 0){ + xpoint = - xpoint; + ypoint = - ypoint; + } + + if(charge < 0){ + xpoint = xpoint; + ypoint = ypoint; + } + center[0] = xpos + xpoint; + center[1] = ypos + ypoint; + + return 1; +} + +Bool_t GetConvPosXY(AliESDtrack* ptrack,AliESDtrack* ntrack, Double_t b, Double_t convpos[2]){ + //see header file for documentation + + Double_t helixcenterpos[2]; + GetHelixCenter(ptrack,b,ptrack->Charge(),helixcenterpos); + + Double_t helixcenterneg[2]; + GetHelixCenter(ntrack,b,ntrack->Charge(),helixcenterneg); + + Double_t poshelix[6]; + ptrack->GetHelixParameters(poshelix,b); + Double_t posradius = TMath::Abs(1./poshelix[4]); + + Double_t neghelix[6]; + ntrack->GetHelixParameters(neghelix,b); + Double_t negradius = TMath::Abs(1./neghelix[4]); + + Double_t xpos = helixcenterpos[0]; + Double_t ypos = helixcenterpos[1]; + Double_t xneg = helixcenterneg[0]; + Double_t yneg = helixcenterneg[1]; + + convpos[0] = (xpos*negradius + xneg*posradius)/(negradius+posradius); + convpos[1] = (ypos*negradius+ yneg*posradius)/(negradius+posradius); + + return 1; +} + + + +Double_t GetConvPosZ(AliESDtrack* ptrack,AliESDtrack* ntrack, Double_t b){ + //see header file for documentation + + Double_t helixpos[6]; + ptrack->GetHelixParameters(helixpos,b); + + Double_t helixneg[6]; + ntrack->GetHelixParameters(helixneg,b); + + Double_t negtrackradius = TMath::Abs(1./helixneg[4]); + Double_t postrackradius = TMath::Abs(1./helixpos[4]); + + Double_t pi = 3.14159265358979323846; + + Double_t convpos[2]; + GetConvPosXY(ptrack,ntrack,b,convpos); + + Double_t convposx = convpos[0]; + Double_t convposy = convpos[1]; + + Double_t helixcenterpos[2]; + GetHelixCenter(ptrack,b,ptrack->Charge(),helixcenterpos); + + Double_t helixcenterneg[2]; + GetHelixCenter(ntrack,b,ntrack->Charge(),helixcenterneg); + + Double_t xpos = helixcenterpos[0]; + Double_t ypos = helixcenterpos[1]; + Double_t xneg = helixcenterneg[0]; + Double_t yneg = helixcenterneg[1]; + + Double_t delta_x_pos = convposx - xpos; + Double_t delta_y_pos = convposy - ypos; + + Double_t delta_x_neg = convposx - xneg; + Double_t delta_y_neg = convposy - yneg; + + Double_t alpha_pos = pi + TMath::ATan2(-delta_y_pos,-delta_x_pos); + Double_t alpha_neg = pi + TMath::ATan2(-delta_y_neg,-delta_x_neg); + + Double_t vertex_x_neg = xneg + TMath::Abs(negtrackradius)* + TMath::Cos(alpha_neg); + Double_t vertex_y_neg = yneg + TMath::Abs(negtrackradius)* + TMath::Sin(alpha_neg); + + Double_t vertex_x_pos = xpos + TMath::Abs(postrackradius)* + TMath::Cos(alpha_pos); + Double_t vertex_y_pos = ypos + TMath::Abs(postrackradius)* + TMath::Sin(alpha_pos); + + Double_t x0neg = helixneg[5]; + Double_t y0neg = helixneg[0]; + + Double_t x0pos = helixpos[5]; + Double_t y0pos = helixpos[0]; + + Double_t d_neg = TMath::Sqrt((vertex_x_neg - x0neg)*(vertex_x_neg - x0neg) + +(vertex_y_neg - y0neg)*(vertex_y_neg - y0neg)); + + Double_t d_pos = TMath::Sqrt((vertex_x_pos - x0pos)*(vertex_x_pos - x0pos) + +(vertex_y_pos - y0pos)*(vertex_y_pos - y0pos)); + + Double_t r_neg = TMath::Sqrt(negtrackradius*negtrackradius - + d_neg*d_neg/4.); + + Double_t r_pos = TMath::Sqrt(postrackradius*postrackradius - + d_pos*d_pos/4.); + + Double_t deltabeta_neg = 2*(pi + TMath::ATan2(-d_neg/2.,-r_neg)); + Double_t deltabeta_pos = 2*(pi + TMath::ATan2(-d_pos/2.,-r_pos)); + + Double_t delta_U_neg = negtrackradius*deltabeta_neg; + Double_t delta_U_pos = postrackradius*deltabeta_pos; + + Double_t zphase_neg = ntrack->GetZ() + delta_U_neg * ntrack->GetTgl(); + Double_t zphase_pos = ptrack->GetZ() + delta_U_pos * ptrack->GetTgl(); + + Double_t convposz = + (zphase_pos*negtrackradius+zphase_neg*postrackradius)/(negtrackradius+postrackradius); + + return convposz; +} diff --git a/PWG4/GammaConv/AliV0Reader.h b/PWG4/GammaConv/AliV0Reader.h index 82c6d39c5f2..589b68bcc83 100644 --- a/PWG4/GammaConv/AliV0Reader.h +++ b/PWG4/GammaConv/AliV0Reader.h @@ -1,506 +1,586 @@ -#ifndef ALIV0READER_H -#define ALIV0READER_H -/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * - * See cxx source for full Copyright notice */ - -//////////////////////////////////////////////// -//--------------------------------------------- -// Class used to do analysis on conversion pairs -//--------------------------------------------- -//////////////////////////////////////////////// - -// --- ROOT system --- -#include "TObject.h" -#include "AliESDv0.h" -#include "AliESDEvent.h" -#include "AliKFParticle.h" -#include "TParticle.h" -#include "AliGammaConversionHistograms.h" -#include - -class TClonesArray; -class TFormula; -class Riostream; -class TChain; -//--- AliRoot system --- - -class AliStack; -class AliESDEvent; -class AliMCEventHandler; -class AliESDInputHandler; -class AliESDVertex; -class AliLog; -class TChain; -class TChain; - - - -class AliV0Reader : public TObject { - - public: - - AliV0Reader(); //constructor - AliV0Reader(const AliV0Reader & g); //copy constructor - AliV0Reader & operator = (const AliV0Reader & g); //assignment operator - virtual ~AliV0Reader() {;} //virtual destructor - /* - *Initialize the reader - */ - void Initialize(); - - /* - * Returns AliESDEvent - */ - AliESDEvent* GetESDEvent() const{return fESDEvent;} - - /* - *Returns the number of v0s in the event, no cuts applied. - */ - Int_t GetNumberOfV0s() const{return fESDEvent->GetNumberOfV0s();} - - /* - * Check if there are any more good v0s left in the v0 stack - * if so, fCurrent v0 is set to this v0 and can be retrieved - * by GetCurrentV0 function. - * returns kFALSE if there is no more good v0s in the v0 stack - */ - Bool_t NextV0(); - - /* - * Returns the v0 at the given index, no checks are done on the v0. - */ - AliESDv0* GetV0(Int_t index); - - /* - * Returns the current v0 - */ - AliESDv0* GetCurrentV0() const{return fCurrentV0;} - - /* - * Returns the negative ESD track which belongs to fCurrentV0 - */ - AliESDtrack* GetNegativeESDTrack(){return fESDEvent->GetTrack(fCurrentV0->GetNindex());} - - /* - * Returns the positive ESD track which belongs to fCurrentV0 - */ - AliESDtrack* GetPositiveESDTrack(){return fESDEvent->GetTrack(fCurrentV0->GetPindex());} - - /* - * Returns the negative KF particle which belongs to fCurrentV0 - */ - AliKFParticle* GetNegativeKFParticle() const{return fCurrentNegativeKFParticle;} - - /* - * Returns the positive KF particle which belongs to fCurrentV0 - */ - AliKFParticle* GetPositiveKFParticle() const{return fCurrentPositiveKFParticle;} - - /* - * Returns the KFParticle object of the 2 tracks. - */ - AliKFParticle* GetMotherCandidateKFCombination() const{return fCurrentMotherKFCandidate;} - - /* - * Checks the probablity that the PID of the particle is what we want it to be. - */ - Bool_t CheckPIDProbability(Double_t negProbCut, Double_t posProbCut); - - /* - * Checks if the PID of the two particles are within our cuts. - */ - void GetPIDProbability(Double_t &negPIDProb, Double_t &posPIDProb); - - /* - *Get the negative MC TParticle from the stack - */ - TParticle * GetNegativeMCParticle() const{return fNegativeMCParticle;} - - /* - *Get the positive MC TParticle from the stack - */ - TParticle * GetPositiveMCParticle() const{return fPositiveMCParticle;} - - /* - *Get the mother MC TParticle from the stack - */ - TParticle * GetMotherMCParticle() const{return fMotherMCParticle;} - - /* - * Flag to see if the v0 particles share the same mother - */ - Bool_t HasSameMCMother(); - - - /* - *Get the PID of the MC mother particle - */ - Int_t GetMotherMCParticlePDGCode() const{if(fMotherMCParticle != NULL){ cout<<"MCParticle exists"<GetPdgCode();} - - /* - *Get the MC stack - */ - AliStack* GetMCStack() const{return fMCStack;} - - /* - *Get the magnetic field from the ESD event - */ - Double_t GetMagneticField() const{return fESDEvent->GetMagneticField();} - - /* - *Get the primary vertex from the esd event - */ - const AliESDVertex *GetPrimaryVertex() const {return fESDEvent->GetPrimaryVertex();} - - /* - * Set the PID of the negative track - */ - void SetNegativeTrackPID(Int_t negTrackPID){fNegativeTrackPID=negTrackPID;} - - /* - * Set the PID of the positive track - */ - void SetPositiveTrackPID(Int_t posTrackPID){fPositiveTrackPID=posTrackPID;} - - /* - * Set the flag to use the kfparticle class. Will also disable the use of esd tracks - */ - void UseKFParticle(){fUseKFParticle = kTRUE; fUseESDTrack = kFALSE;} - - /* - * Set the flag to use the esd track class. Will also disable the use of kf particles - */ - void UseESDTrack(){fUseESDTrack = kTRUE; fUseKFParticle = kFALSE;} - - /* - * Set the flag to use improved vertex or not - */ - void SetUseImprovedVertex(Bool_t useImprovedVertex){fUseImprovedVertex=useImprovedVertex;} - - /* - * Return the number in the species array belonging to the negative or positive track pid. - */ - Int_t GetSpeciesIndex(Int_t chargeOfTrack); - - /* - * Return the x coordinate of the v0 - */ - Double_t GetX() const{return fCurrentXValue;} - - /* - * Return the y coordinate of the v0 - */ - Double_t GetY() const{return fCurrentYValue;} - - /* - * Return the Z coordinate of the v0 - */ - Double_t GetZ() const{return fCurrentZValue;} - - /* - * Return the radius of the v0 - */ - Double_t GetXYRadius() const{return sqrt((Double_t)(fCurrentXValue*fCurrentXValue + fCurrentYValue*fCurrentYValue));} - - /* - * Get the opening angle between the two tracks - */ - Double_t GetOpeningAngle(){return fNegativeTrackLorentzVector->Angle(fPositiveTrackLorentzVector->Vect());} - - /* - * Gets the Energy of the negative track. - */ - Double_t GetNegativeTrackEnergy() const{return fCurrentNegativeKFParticle->E();} - - /* - * Gets the Energy of the positive track. - */ - Double_t GetPositiveTrackEnergy() const{return fCurrentPositiveKFParticle->E();} - - /* - * Gets the Energy of the mother candidate. - */ - Double_t GetMotherCandidateEnergy() const{return fCurrentMotherKFCandidate->E();} - - /* - * Gets the Pt of the negative track. - */ - Double_t GetNegativeTrackPt() const{return fNegativeTrackLorentzVector->Pt();} - - /* - * Gets the Pt of the positive track. - */ - Double_t GetPositiveTrackPt() const{return fPositiveTrackLorentzVector->Pt();} - - /* - * Gets the Pt of the mother candidate. - */ - Double_t GetMotherCandidatePt() const{return fMotherCandidateLorentzVector->Pt();} - - /* - * Gets the Eta of the negative track. - */ - Double_t GetNegativeTrackEta() const{return fNegativeTrackLorentzVector->Eta();} - /* - * Gets the Eta of the positive track. - */ - Double_t GetPositiveTrackEta() const{return fPositiveTrackLorentzVector->Eta();} - /* - * Gets the Eta of the mother candidate. - */ - Double_t GetMotherCandidateEta() const{return fMotherCandidateLorentzVector->Eta();} - - /* - * Gets the NDF of the mother candidate. - */ - Double_t GetMotherCandidateNDF() const{return fCurrentMotherKFCandidate->GetNDF();} - - /* - * Gets the Chi2 of the mother candidate. - */ - Double_t GetMotherCandidateChi2() const{return fCurrentMotherKFCandidate->GetChi2();} - - /* - * Gets the Mass of the mother candidate. - */ - Double_t GetMotherCandidateMass() const{return fMotherCandidateKFMass;} - - /* - * Gets the Width of the mother candidate. - */ - Double_t GetMotherCandidateWidth() const{return fMotherCandidateKFWidth;} - - /* - * Gets the Phi of the negative track. - */ - Double_t GetNegativeTrackPhi() const; - - /* - * Gets the Phi of the positive track. - */ - Double_t GetPositiveTrackPhi() const; - - /* - * Gets the Phi of the mother candidate. - */ - Double_t GetMotherCandidatePhi() const; - - /* - * Gets the Rapidity of the mother candidate. - */ - Double_t GetMotherCandidateRapidity() const; - - - /* - * Update data which need to be updated every event. - */ - void UpdateEventByEventData(); - - /* - * Gets the MaxRCut value. - */ - Double_t GetMaxRCut() const{return fMaxR;} - - /* - * Gets the Eta cut value. - */ - Double_t GetEtaCut() const{return fEtaCut;} - - /* - * Gets the Pt cut value. - */ - Double_t GetPtCut() const{return fPtCut;} - - /* - * Gets the line cut values. - */ - Double_t GetLineCutZRSlope() const{return fLineCutZRSlope;} - Double_t GetLineCutZValue() const{return fLineCutZValue;} - - /* - * Gets the Chi2 cut value for the conversions. - */ - Double_t GetChi2CutConversion() const{return fChi2CutConversion;} - - /* - * Gets the Chi2 cut value for the mesons. - */ - Double_t GetChi2CutMeson() const{return fChi2CutMeson;} - - Double_t GetPositiveTrackLength() const{return fCurrentPositiveESDTrack->GetIntegratedLength();} - Double_t GetNegativeTrackLength() const{return fCurrentNegativeESDTrack->GetIntegratedLength();} - - Double_t GetPositiveNTPCClusters() const{return fCurrentPositiveESDTrack->GetTPCNcls();} - Double_t GetNegativeNTPCClusters() const{return fCurrentNegativeESDTrack->GetTPCNcls();} - - /* - * Sets the MaxRCut value. - */ - void SetMaxRCut(Double_t maxR){fMaxR=maxR;} - - /* - * Sets the EtaCut value. - */ - void SetEtaCut(Double_t etaCut){fEtaCut=etaCut;} - - /* - * Sets the PtCut value. - */ - void SetPtCut(Double_t ptCut){fPtCut=ptCut;} - - /* - * Sets the LineCut values. - */ - void SetLineCutZRSlope(Double_t LineCutZRSlope){fLineCutZRSlope=LineCutZRSlope;} - void SetLineCutZValue(Double_t LineCutZValue){fLineCutZValue=LineCutZValue;} - - /* - * Sets the Chi2Cut value for conversions. - */ - void SetChi2CutConversion(Double_t chi2){fChi2CutConversion=chi2;} - - /* - * Sets the Chi2Cut for the mesons. - */ - void SetChi2CutMeson(Double_t chi2){fChi2CutMeson=chi2;} - - /* - * Sets the XVertexCut value. - */ - void SetXVertexCut(Double_t xVtx){fCurrentXValue=xVtx;} - - /* - * Sets the YVertexCut value. - */ - void SetYVertexCut(Double_t yVtx){fCurrentYValue=yVtx;} - - /* - * Sets the ZVertexCut value. - */ - void SetZVertexCut(Double_t zVtx){fCurrentZValue=zVtx;} - - /* - * Sets the PIDProbabilityCut value for track particles. - */ - void SetPIDProbability(Double_t pidProb){fPIDProbabilityCutPositiveParticle=pidProb; fPIDProbabilityCutNegativeParticle=pidProb;} - - /* - * Sets the PIDProbability cut value for the negative track. - */ - void SetPIDProbabilityNegativeParticle(Double_t pidProb){fPIDProbabilityCutNegativeParticle=pidProb;} - - /* - * Sets the PIDProbability cut value for the positive track. - */ - void SetPIDProbabilityPositiveParticle(Double_t pidProb){fPIDProbabilityCutPositiveParticle=pidProb;} - - /* - * Sets the SigmaMassCut value. - */ - void SetSigmaMass(Double_t sigmaMass){fNSigmaMass=sigmaMass;} - - /* - * Sets the flag to enable/disable the usage of MC information. - */ - void SetDoMCTruth(Bool_t doMC){fDoMC = doMC;} - - /* - * Updates the V0 information of the current V0. - */ - Bool_t UpdateV0Information(); - - /* - * Resets the V0 index. - */ - void ResetV0IndexNumber(){fCurrentV0IndexNumber=0;} - - /* - * Sets the histograms. - */ - void SetHistograms(AliGammaConversionHistograms *histograms){fHistograms=histograms;} - - /* - * Check for primary vertex. - */ - Bool_t CheckForPrimaryVertex(); - - /* - * Gets a vector of good v0s. - */ - vector GetCurrentEventGoodV0s() const{return fCurrentEventGoodV0s;} - - /* - * Gets the vector of previous events v0s (for bacground analysis) - */ - vector GetPreviousEventGoodV0s() const{return fPreviousEventGoodV0s;} - - private: - AliStack * fMCStack; // pointer to MonteCarlo particle stack - AliMCEventHandler* fMCTruth; // pointer to the MC event handler - TChain * fChain; // pointer to the TChain - - AliESDInputHandler* fESDHandler; //! pointer to esd object - AliESDEvent *fESDEvent; //! pointer to esd object - - AliGammaConversionHistograms *fHistograms; //! pointer to histogram handling class - - Int_t fCurrentV0IndexNumber; - AliESDv0 * fCurrentV0; //! pointer to the current v0 - AliKFParticle * fCurrentNegativeKFParticle; //! pointer to the negative KF particle - AliKFParticle * fCurrentPositiveKFParticle; //! pointer to the positive KF particle - AliKFParticle * fCurrentMotherKFCandidate; //! pointer to the positive KF particle - - AliESDtrack * fCurrentNegativeESDTrack; //! pointer to the negative ESD track - AliESDtrack * fCurrentPositiveESDTrack; //! pointer to the positive ESD track - - TLorentzVector * fNegativeTrackLorentzVector; //! pointer to the negative Track Lorentz Vector - TLorentzVector * fPositiveTrackLorentzVector; //! pointer to the positive Track Lorentz Vector - TLorentzVector * fMotherCandidateLorentzVector; //! pointer to the mother candidate Track Lorentz Vector - - Double_t fCurrentXValue; // current x value - Double_t fCurrentYValue; // current y value - Double_t fCurrentZValue; // current z value - - Int_t fPositiveTrackPID; // positive track pid - Int_t fNegativeTrackPID; // negative track pid - - TParticle *fNegativeMCParticle; //! - TParticle *fPositiveMCParticle; //! - TParticle *fMotherMCParticle; //! - - Double_t fMotherCandidateKFMass; // mass of mother candidate KF particle - Double_t fMotherCandidateKFWidth; // width of mother candidate KF particle - - Bool_t fUseKFParticle; // flag - Bool_t fUseESDTrack; // flag - Bool_t fDoMC; // flag - - //cuts - Double_t fMaxR; //r cut - Double_t fEtaCut; //eta cut - Double_t fPtCut; // pt cut - Double_t fLineCutZRSlope; //linecut - Double_t fLineCutZValue; //linecut - Double_t fChi2CutConversion; //chi2cut - Double_t fChi2CutMeson; //chi2cut - Double_t fPIDProbabilityCutNegativeParticle; //pid cut - Double_t fPIDProbabilityCutPositiveParticle; //pid cut - Double_t fXVertexCut; //vertex cut - Double_t fYVertexCut; //vertex cut - Double_t fZVertexCut; // vertexcut - - Double_t fNSigmaMass; //nsigma cut - - Bool_t fUseImprovedVertex; //flag - - vector fCurrentEventGoodV0s; //vector of good v0s - vector fPreviousEventGoodV0s; // vector of good v0s from prevous events - - ClassDef(AliV0Reader,2) -}; - - -#endif - - - +#ifndef ALIV0READER_H +#define ALIV0READER_H +/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * + * See cxx source for full Copyright notice */ + +//////////////////////////////////////////////// +//--------------------------------------------- +// Class used to do analysis on conversion pairs +//--------------------------------------------- +//////////////////////////////////////////////// + +// --- ROOT system --- +#include "TObject.h" +#include "AliMCEvent.h" // for CF +#include "AliESDv0.h" +#include "AliESDEvent.h" +#include "AliKFParticle.h" +#include "TParticle.h" +#include "AliGammaConversionHistograms.h" +#include +#include "AliCFManager.h" + + +class TClonesArray; +class TFormula; +class Riostream; +class TChain; + +//--- AliRoot system --- + +class AliStack; +class AliMCEvent; // for CF +class AliESDEvent; +class AliMCEventHandler; +class AliESDInputHandler; +class AliESDVertex; +class AliLog; +class TChain; +class TChain; +class AliCFManager; // for CF +class AliCFContainer; // for CF + + + +class AliV0Reader : public TObject { + + public: + + + // for CF + enum{ + kStepGenerated = 0, + kStepReconstructable = 1, + kStepLikeSign = 2, + kStepTPCRefit = 3, + kStepKinks = 4, + kStepGetOnFly = 5, + kStepNContributors = 6, + kStepTPCPID = 7, + kStepR = 8, + kStepLine = 9, + kStepZ = 10, + kStepNDF = 11, + kStepChi2 = 12, + kStepEta = 13, + kStepPt = 14 + }; + + + + AliV0Reader(); //constructor + AliV0Reader(const AliV0Reader & g); //copy constructor + AliV0Reader & operator = (const AliV0Reader & g); //assignment operator + virtual ~AliV0Reader() {;} //virtual destructor + /* + *Initialize the reader + */ + void Initialize(); + + + // for CF + void SetCFManager(AliCFManager *io){fCFManager = io;}; + AliCFManager *GetCFManager() const {return fCFManager;} + + + + + /* + * Returns AliESDEvent + */ + AliESDEvent* GetESDEvent() const{return fESDEvent;} + + /* + *Returns the number of v0s in the event, no cuts applied. + */ + Int_t GetNumberOfV0s() const{return fESDEvent->GetNumberOfV0s();} + + /* + * Check if there are any more good v0s left in the v0 stack + * if so, fCurrent v0 is set to this v0 and can be retrieved + * by GetCurrentV0 function. + * returns kFALSE if there is no more good v0s in the v0 stack + */ + Bool_t NextV0(); + + /* + * Returns the v0 at the given index, no checks are done on the v0. + */ + AliESDv0* GetV0(Int_t index); + + /* + * Returns the current v0 + */ + AliESDv0* GetCurrentV0() const{return fCurrentV0;} + + /* + * Returns the negative ESD track which belongs to fCurrentV0 + */ + AliESDtrack* GetNegativeESDTrack(){return fESDEvent->GetTrack(fCurrentV0->GetNindex());} + + /* + * Returns the positive ESD track which belongs to fCurrentV0 + */ + AliESDtrack* GetPositiveESDTrack(){return fESDEvent->GetTrack(fCurrentV0->GetPindex());} + + /* + * Returns the negative KF particle which belongs to fCurrentV0 + */ + AliKFParticle* GetNegativeKFParticle() const{return fCurrentNegativeKFParticle;} + + /* + * Returns the positive KF particle which belongs to fCurrentV0 + */ + AliKFParticle* GetPositiveKFParticle() const{return fCurrentPositiveKFParticle;} + + /* + * Returns the KFParticle object of the 2 tracks. + */ + AliKFParticle* GetMotherCandidateKFCombination() const{return fCurrentMotherKFCandidate;} + + /* + * Checks the probablity that the PID of the particle is what we want it to be. + */ + Bool_t CheckPIDProbability(Double_t negProbCut, Double_t posProbCut); + + /* + * Checks if the PID of the two particles are within our cuts. + */ + void GetPIDProbability(Double_t &negPIDProb, Double_t &posPIDProb); + + /* + *Get the negative MC TParticle from the stack + */ + TParticle * GetNegativeMCParticle() const{return fNegativeMCParticle;} + + /* + *Get the positive MC TParticle from the stack + */ + TParticle * GetPositiveMCParticle() const{return fPositiveMCParticle;} + + /* + *Get the mother MC TParticle from the stack + */ + TParticle * GetMotherMCParticle() const{return fMotherMCParticle;} + + /* + * Flag to see if the v0 particles share the same mother + */ + Bool_t HasSameMCMother(); + + + /* + *Get the PID of the MC mother particle + */ + Int_t GetMotherMCParticlePDGCode() const{if(fMotherMCParticle != NULL){ cout<<"MCParticle exists"<GetPdgCode();} + + /* + *Get the MC stack + */ + AliStack* GetMCStack() const{return fMCStack;} + + + /* + * Setup AliMCEventHandler + */ + AliMCEventHandler* GetMCTruth() const{return fMCTruth;} // for CF + + + /* + *Get the MC stack + */ + AliMCEvent* GetMCEvent() const{return fMCEvent;} // for CF + + + /* + *Get the magnetic field from the ESD event + */ + Double_t GetMagneticField() const{return fESDEvent->GetMagneticField();} + + /* + *Get the primary vertex from the esd event + */ + const AliESDVertex *GetPrimaryVertex() const {return fESDEvent->GetPrimaryVertex();} + + /* + * Set the PID of the negative track + */ + void SetNegativeTrackPID(Int_t negTrackPID){fNegativeTrackPID=negTrackPID;} + + /* + * Set the PID of the positive track + */ + void SetPositiveTrackPID(Int_t posTrackPID){fPositiveTrackPID=posTrackPID;} + + /* + * Set the flag to use the kfparticle class. Will also disable the use of esd tracks + */ + void UseKFParticle(){fUseKFParticle = kTRUE; fUseESDTrack = kFALSE;} + + /* + * Set the flag to use the esd track class. Will also disable the use of kf particles + */ + void UseESDTrack(){fUseESDTrack = kTRUE; fUseKFParticle = kFALSE;} + + /* + * Set the flag to use improved vertex or not + */ + void SetUseImprovedVertex(Bool_t useImprovedVertex){fUseImprovedVertex=useImprovedVertex;} + + /* + * Return the number in the species array belonging to the negative or positive track pid. + */ + Int_t GetSpeciesIndex(Int_t chargeOfTrack); + + /* + * Return the x coordinate of the v0 + */ + Double_t GetX() const{return fCurrentXValue;} + + /* + * Return the y coordinate of the v0 + */ + Double_t GetY() const{return fCurrentYValue;} + + /* + * Return the Z coordinate of the v0 + */ + Double_t GetZ() const{return fCurrentZValue;} + + /* + * Return the radius of the v0 + */ + Double_t GetXYRadius() const{return sqrt((Double_t)(fCurrentXValue*fCurrentXValue + fCurrentYValue*fCurrentYValue));} + + /* + * Get the opening angle between the two tracks + */ + Double_t GetOpeningAngle(){return fNegativeTrackLorentzVector->Angle(fPositiveTrackLorentzVector->Vect());} + + /* + * Gets the Energy of the negative track. + */ + Double_t GetNegativeTrackEnergy() const{return fCurrentNegativeKFParticle->E();} + + /* + * Gets the Energy of the positive track. + */ + Double_t GetPositiveTrackEnergy() const{return fCurrentPositiveKFParticle->E();} + + /* + * Gets the Energy of the mother candidate. + */ + Double_t GetMotherCandidateEnergy() const{return fCurrentMotherKFCandidate->E();} + + /* + * Gets the Pt of the negative track. + */ + Double_t GetNegativeTrackPt() const{return fNegativeTrackLorentzVector->Pt();} + + /* + * Gets the Pt of the positive track. + */ + Double_t GetPositiveTrackPt() const{return fPositiveTrackLorentzVector->Pt();} + + /* + * Gets the Pt of the mother candidate. + */ + Double_t GetMotherCandidatePt() const{return fMotherCandidateLorentzVector->Pt();} + + /* + * Gets the Eta of the negative track. + */ + Double_t GetNegativeTrackEta() const{return fNegativeTrackLorentzVector->Eta();} + /* + * Gets the Eta of the positive track. + */ + Double_t GetPositiveTrackEta() const{return fPositiveTrackLorentzVector->Eta();} + /* + * Gets the Eta of the mother candidate. + */ + Double_t GetMotherCandidateEta() const{return fMotherCandidateLorentzVector->Eta();} + + /* + * Gets the NDF of the mother candidate. + */ + Double_t GetMotherCandidateNDF() const{return fCurrentMotherKFCandidate->GetNDF();} + + /* + * Gets the Chi2 of the mother candidate. + */ + Double_t GetMotherCandidateChi2() const{return fCurrentMotherKFCandidate->GetChi2();} + + /* + * Gets the Mass of the mother candidate. + */ + Double_t GetMotherCandidateMass() const{return fMotherCandidateKFMass;} + + /* + * Gets the Width of the mother candidate. + */ + Double_t GetMotherCandidateWidth() const{return fMotherCandidateKFWidth;} + + /* + * Gets the Phi of the negative track. + */ + Double_t GetNegativeTrackPhi() const; + + /* + * Gets the Phi of the positive track. + */ + Double_t GetPositiveTrackPhi() const; + + /* + * Gets the Phi of the mother candidate. + */ + Double_t GetMotherCandidatePhi() const; + + /* + * Gets the Rapidity of the mother candidate. + */ + Double_t GetMotherCandidateRapidity() const; + + + /* + * Update data which need to be updated every event. + */ + void UpdateEventByEventData(); + + /* + * Gets the MaxRCut value. + */ + Double_t GetMaxRCut() const{return fMaxR;} + + /* + * Gets the Eta cut value. + */ + Double_t GetEtaCut() const{return fEtaCut;} + + /* + * Gets the Pt cut value. + */ + Double_t GetPtCut() const{return fPtCut;} + + + /* + * Gets the MaxZCut value. + */ + Double_t GetMaxZCut() const{return fMaxZ;} + + + /* + * Gets the line cut values. + */ + Double_t GetLineCutZRSlope() const{return fLineCutZRSlope;} + Double_t GetLineCutZValue() const{return fLineCutZValue;} + + /* + * Gets the Chi2 cut value for the conversions. + */ + Double_t GetChi2CutConversion() const{return fChi2CutConversion;} + + /* + * Gets the Chi2 cut value for the mesons. + */ + Double_t GetChi2CutMeson() const{return fChi2CutMeson;} + + Double_t GetPositiveTrackLength() const{return fCurrentPositiveESDTrack->GetIntegratedLength();} + Double_t GetNegativeTrackLength() const{return fCurrentNegativeESDTrack->GetIntegratedLength();} + + Double_t GetPositiveNTPCClusters() const{return fCurrentPositiveESDTrack->GetTPCNcls();} + Double_t GetNegativeNTPCClusters() const{return fCurrentNegativeESDTrack->GetTPCNcls();} + + /* + * Sets the MaxRCut value. + */ + void SetMaxRCut(Double_t maxR){fMaxR=maxR;} + + /* + * Sets the EtaCut value. + */ + void SetEtaCut(Double_t etaCut){fEtaCut=etaCut;} + + /* + * Sets the PtCut value. + */ + void SetPtCut(Double_t ptCut){fPtCut=ptCut;} + + + /* + * Sets the MaxZCut value. + */ + void SetMaxZCut(Double_t maxZ){fMaxZ=maxZ;} + + + /* + * Sets the LineCut values. + */ + void SetLineCutZRSlope(Double_t LineCutZRSlope){fLineCutZRSlope=LineCutZRSlope;} + void SetLineCutZValue(Double_t LineCutZValue){fLineCutZValue=LineCutZValue;} + + /* + * Sets the Chi2Cut value for conversions. + */ + void SetChi2CutConversion(Double_t chi2){fChi2CutConversion=chi2;} + + /* + * Sets the Chi2Cut for the mesons. + */ + void SetChi2CutMeson(Double_t chi2){fChi2CutMeson=chi2;} + + /* + * Sets the XVertexCut value. + */ + void SetXVertexCut(Double_t xVtx){fCurrentXValue=xVtx;} + + /* + * Sets the YVertexCut value. + */ + void SetYVertexCut(Double_t yVtx){fCurrentYValue=yVtx;} + + /* + * Sets the ZVertexCut value. + */ + void SetZVertexCut(Double_t zVtx){fCurrentZValue=zVtx;} + + /* + * Sets the PIDProbabilityCut value for track particles. + */ + void SetPIDProbability(Double_t pidProb){fPIDProbabilityCutPositiveParticle=pidProb; fPIDProbabilityCutNegativeParticle=pidProb;} + + /* + * Sets the PIDProbability cut value for the negative track. + */ + void SetPIDProbabilityNegativeParticle(Double_t pidProb){fPIDProbabilityCutNegativeParticle=pidProb;} + + /* + * Sets the PIDProbability cut value for the positive track. + */ + void SetPIDProbabilityPositiveParticle(Double_t pidProb){fPIDProbabilityCutPositiveParticle=pidProb;} + + /* + * Sets the SigmaMassCut value. + */ + void SetSigmaMass(Double_t sigmaMass){fNSigmaMass=sigmaMass;} + + /* + * Sets the flag to enable/disable the usage of MC information. + */ + void SetDoMCTruth(Bool_t doMC){fDoMC = doMC;} + + /* + * Updates the V0 information of the current V0. + */ + Bool_t UpdateV0Information(); + + /* + * Resets the V0 index. + */ + void ResetV0IndexNumber(){fCurrentV0IndexNumber=0;} + + /* + * Sets the histograms. + */ + void SetHistograms(AliGammaConversionHistograms *histograms){fHistograms=histograms;} + + /* + * Check for primary vertex. + */ + Bool_t CheckForPrimaryVertex(); + + /* + * Gets a vector of good v0s. + */ + vector GetCurrentEventGoodV0s() const{return fCurrentEventGoodV0s;} + + /* + * Gets the vector of previous events v0s (for bacground analysis) + */ + vector GetPreviousEventGoodV0s() const{return fPreviousEventGoodV0s;} + + void SetUseOwnXYZCalculation(Bool_t flag){fUseOwnXYZCalculation=flag;} + + Bool_t GetHelixCenter(AliESDtrack* track, Double_t b,Int_t charge, Double_t center[2]); + + Bool_t GetConvPosXY(AliESDtrack* ptrack,AliESDtrack* ntrack, Double_t b, Double_t convpos[2]); + + Double_t GetConvPosZ(AliESDtrack* ptrack,AliESDtrack* ntrack, Double_t b); + + private: + AliStack * fMCStack; // pointer to MonteCarlo particle stack + AliMCEventHandler* fMCTruth; // for CF pointer to the MC object + AliMCEvent *fMCEvent; // for CF pointer to MC event + TChain * fChain; // pointer to the TChain + + AliESDInputHandler* fESDHandler; //! pointer to esd object + AliESDEvent *fESDEvent; //! pointer to esd object + + + // for CF + AliCFManager *fCFManager; + AliCFContainer *container; + + + AliGammaConversionHistograms *fHistograms; //! pointer to histogram handling class + + Int_t fCurrentV0IndexNumber; + AliESDv0 * fCurrentV0; //! pointer to the current v0 + AliKFParticle * fCurrentNegativeKFParticle; //! pointer to the negative KF particle + AliKFParticle * fCurrentPositiveKFParticle; //! pointer to the positive KF particle + AliKFParticle * fCurrentMotherKFCandidate; //! pointer to the positive KF particle + + AliESDtrack * fCurrentNegativeESDTrack; //! pointer to the negative ESD track + AliESDtrack * fCurrentPositiveESDTrack; //! pointer to the positive ESD track + + TLorentzVector * fNegativeTrackLorentzVector; //! pointer to the negative Track Lorentz Vector + TLorentzVector * fPositiveTrackLorentzVector; //! pointer to the positive Track Lorentz Vector + TLorentzVector * fMotherCandidateLorentzVector; //! pointer to the mother candidate Track Lorentz Vector + + Double_t fCurrentXValue; // current x value + Double_t fCurrentYValue; // current y value + Double_t fCurrentZValue; // current z value + + Int_t fPositiveTrackPID; // positive track pid + Int_t fNegativeTrackPID; // negative track pid + + TParticle *fNegativeMCParticle; //! + TParticle *fPositiveMCParticle; //! + TParticle *fMotherMCParticle; //! + + Double_t fMotherCandidateKFMass; // mass of mother candidate KF particle + Double_t fMotherCandidateKFWidth; // width of mother candidate KF particle + + Bool_t fUseKFParticle; // flag + Bool_t fUseESDTrack; // flag + Bool_t fDoMC; // flag + + //cuts + Double_t fMaxR; //r cut + Double_t fEtaCut; //eta cut + Double_t fPtCut; // pt cut + Double_t fMaxZ; //z cut + Double_t fLineCutZRSlope; //linecut + Double_t fLineCutZValue; //linecut + Double_t fChi2CutConversion; //chi2cut + Double_t fChi2CutMeson; //chi2cut + Double_t fPIDProbabilityCutNegativeParticle; //pid cut + Double_t fPIDProbabilityCutPositiveParticle; //pid cut + Double_t fXVertexCut; //vertex cut + Double_t fYVertexCut; //vertex cut + Double_t fZVertexCut; // vertexcut + + Double_t fNSigmaMass; //nsigma cut + + Bool_t fUseImprovedVertex; //flag + + Bool_t fUseOwnXYZCalculation; //flag that determines if we use our own calculation of xyz (markus) + + vector fCurrentEventGoodV0s; //vector of good v0s + vector fPreviousEventGoodV0s; // vector of good v0s from prevous events + + ClassDef(AliV0Reader,3) +}; +#endif + + + diff --git a/PWG4/Makefile b/PWG4/Makefile index dd4a853846c..aa8ff70183a 100644 --- a/PWG4/Makefile +++ b/PWG4/Makefile @@ -34,6 +34,10 @@ ifneq ($(JETAN_INCLUDE),) ALICEINC += -I../$(JETAN_INCLUDE) endif +ifneq ($(CORRFW_INCLUDE),) + ALICEINC += -I../$(CORRFW_INCLUDE) +endif + ifneq ($(PHOSUtils_INCLUDE),) ALICEINC += -I../$(PHOSUtils_INCLUDE) CXXFLAGS+=-D__PHOSUTIL__ diff --git a/PWG4/macros/ConfigGammaConversion.C b/PWG4/macros/ConfigGammaConversion.C index f6db6d11bf3..cfa9459250a 100644 --- a/PWG4/macros/ConfigGammaConversion.C +++ b/PWG4/macros/ConfigGammaConversion.C @@ -23,6 +23,7 @@ Double_t kGCLineCutZValue = 7.; Double_t kGCmaxRCut = 180.; Double_t kGCetaCut = 1.2; Double_t kGCptCut = 0.02; +Double_t kGCmaxZCut = 240.; Double_t kGCchi2CutConversion = 20.; Double_t kGCchi2CutMeson = 20.; @@ -54,8 +55,8 @@ Double_t kGCminOpeningAngleGhostCut = 0.01; /** ----------------------------------end define cuts here----------------------------------*/ /** -------------------------------- Phi/R Mapping ---------------------------------------*/ -Int_t kGCnPhiIndex = 18; -Int_t kGCnRIndex = 40; +Int_t kGCnPhiIndex = 8; +Int_t kGCnRIndex = 4; Double_t kGCminRadius = 0.; Double_t kGCmaxRadius = 200.; @@ -63,6 +64,10 @@ Double_t kGCminPhi = -TMath::Pi(); Double_t kGCmaxPhi = TMath::Pi(); /** ------------------------------- end Phi/R Mapping ------------------------------------*/ +Bool_t kGCdoOwnXYZCalculation = kFALSE; + +Bool_t fWriteStandardAOD =kFALSE; + /** ------------------- define which histograms to plot here --------------------------------*/ /** NB: to change the bin numbers, see below the histogram flags */ @@ -262,7 +267,7 @@ Bool_t kGCplotESDBackgroundZR = kTRUE; Bool_t kGCplotESDBackgroundXY = kTRUE; Bool_t kGCplotESDBackgroundRapid = kTRUE; -Bool_t kGCplotMapping = kFALSE; +Bool_t kGCplotMapping = kTRUE; Bool_t kGCplotResolutiondPt = kTRUE; Bool_t kGCplotResolutiondR = kTRUE; @@ -676,9 +681,13 @@ Bool_t scanArguments(TString arguments){ cout<<"Switching off kGCdoMCTruth"<=pTokens->GetEntries()))) break; - kGCoutputFileAppendix = "_"+((TObjString*)pTokens->At(i))->GetString(); + kGCoutputFileAppendix = TString("_")+((TObjString*)pTokens->At(i))->GetString(); if(kGCoutputFileAppendix.IsNull()){ cout<<"-appending-to-output-file is NULL"<Getenv("ALICE_ROOT"); file+="/PWG4/PartCorr/AliAnalysisTaskGammaConversion.cxx"; - + ifstream stream; stream.open(file.Data()); - + if(!stream){ kGCusePWG4PartCorr=kFALSE; } @@ -724,15 +733,15 @@ AliAnalysisTaskGammaConversion* ConfigGammaConversion(TString arguments,AliAnaly if(!scanArguments(arguments)){ break; } - + SetVersionLibrary(); // checks if PWG4GammaConv or PWG4PartCorr is used - + if(cin_esd == NULL && kGCrunOnTrain == kTRUE){ cout<<"Error: kGCrunOnTrain flag is set to true but the input AliAnalysisDataContainer is NULL"<LoadMacro("$ALICE_ROOT/PWG0/CreateESDChain.C"); // load the CreateChain macro - + AliLog::SetGlobalLogLevel(AliLog::kError); + + + // ------------------------------------------------------------------------ + + // for CF + + //Container def. + const Double_t ptmin = kGCfirstXBinPt; + const Double_t ptmax = kGClastXBinPt; + const Double_t etamin = kGCfirstXBinEta; + const Double_t etamax = kGClastXBinEta; + const Double_t massmin = kGCfirstXBinPi0Mass; + const Double_t massmax = kGClastXBinPi0Mass; + + + // sensitive variables + UInt_t ipt = 0; + UInt_t ieta = 1; + UInt_t imass = 2; + + //how many selection steps + UInt_t nstep = 15; + const Int_t nvar = 3; + const Int_t nbin0 = kGCnXBinsPt; + const Int_t nbin1 = kGCnXBinsEta; + const Int_t nbin2 = kGCnXBinsPi0Mass; + + //arrays for the number of bins in each dimension + Int_t iBin[nvar]; + iBin[0] = nbin0; + iBin[1] = nbin1; + iBin[2] = nbin2; + + //arrays for lower bounds + Double_t *binLim0 = new Double_t[nbin0+1]; + Double_t *binLim1 = new Double_t[nbin1+1]; + Double_t *binLim2 = new Double_t[nbin2+1]; + + // values for lower bounds + for(Int_t i = 0; i <= nbin0; i++) binLim0[i] = ptmin + (ptmax - ptmin)/nbin0*i; + for(Int_t i = 0; i <= nbin1; i++) binLim1[i] = etamin + (etamax - etamin)/nbin1*i; + for(Int_t i = 0; i <= nbin2; i++) binLim2[i] = massmin + (massmax - massmin)/nbin2*i; + + // create container + AliCFContainer *container = new AliCFContainer("container","container for gammaconversion", nstep,nvar,iBin); + container->SetBinLimits(ipt,binLim0); + container->SetBinLimits(ieta,binLim1); + container->SetBinLimits(imass,binLim2); + + AliCFManager *man = new AliCFManager(); + man->SetParticleContainer(container); + + // end --------------------------------------------------------------------------- + + } - + AliGammaConversionHistograms* histograms = new AliGammaConversionHistograms(); AddHistograms(histograms); // Create the Analysis manager AliAnalysisManager *mgr =NULL; if(kGCrunOnTrain == kFALSE){ - mgr = new AliAnalysisManager("My Manager", "My Analysis"); + mgr = new AliAnalysisManager("My Manager", "My Analysis"); } else{ mgr = AliAnalysisManager::GetAnalysisManager(); } - + if (!mgr) { ::Error("ConfigGammaConversion", "No analysis manager to connect to."); return NULL; @@ -787,7 +851,7 @@ AliAnalysisTaskGammaConversion* ConfigGammaConversion(TString arguments,AliAnaly } } AliESDInputHandler* inpHandler = NULL; - + if(kGCrunOnTrain == kFALSE){ // Define Input Event Handler inpHandler = new AliESDInputHandler(); @@ -806,22 +870,22 @@ AliAnalysisTaskGammaConversion* ConfigGammaConversion(TString arguments,AliAnaly return NULL; } } - - + // Define Output Event Handler and ad - if(kGCrunOnTrain == kFALSE){ + // if(kGCrunOnTrain == kFALSE){ + if(fWriteStandardAOD == kTRUE){ AliAODHandler* aodHandler = new AliAODHandler(); TString fileOutAOD = "AOD_"+ kGCoutputFileName + kGCoutputFileAppendix + ".root"; aodHandler->SetOutputFileName(fileOutAOD); mgr->SetOutputEventHandler (aodHandler); } - + if(kGCrunOnTrain == kFALSE){ mgr->SetInputEventHandler (inpHandler); mgr->SetMCtruthEventHandler(mcHandler); } // Be sure you are told what you are doing - // mgr->SetDebugLevel(10); + mgr->SetDebugLevel(10); // Declare Common Input Tchain AliAnalysisDataContainer *cinput1 = NULL; @@ -835,7 +899,7 @@ AliAnalysisTaskGammaConversion* ConfigGammaConversion(TString arguments,AliAnaly } else{ if(kGCrunOnTrain == kFALSE){ - cinput1 = mgr->GetCommonInputContainer(); + cinput1 = mgr->GetCommonInputContainer(); // added by kenneth to avoid writing the standard AOD } else{ // cinput = cin_esd; @@ -849,7 +913,9 @@ AliAnalysisTaskGammaConversion* ConfigGammaConversion(TString arguments,AliAnaly coutput1 = mgr->CreateContainer("tree",TTree::Class(),AliAnalysisManager::kOutputContainer, "default"); } else{ - coutput1 = mgr->GetCommonOutputContainer(); + if(fWriteStandardAOD){ + coutput1 = mgr->GetCommonOutputContainer(); + } } // Private output objects @@ -860,20 +926,24 @@ AliAnalysisTaskGammaConversion* ConfigGammaConversion(TString arguments,AliAnaly kGCoutputFileAppendix.ReplaceAll(".root",""); } TString fileOut = kGCoutputFileName + kGCoutputFileAppendix + ".root"; - + AliAnalysisDataContainer *coutput2 = mgr->CreateContainer("histogramsAliGammaConversion", TList::Class(),AliAnalysisManager::kOutputContainer, fileOut); + // for CF + AliAnalysisDataContainer *coutput3 = mgr->CreateContainer("ccontainer0",AliCFContainer::Class(),AliAnalysisManager::kOutputContainer,fileOut); + //------------------------ END: Define input/output handlers --------------------------------------------------- //check for errors in the specified data if(kGCuseKFParticle == kTRUE && kGCuseESDTrack == kTRUE){ //Print warning, cannot use both ::Error("ConfigGammaConversion","Both kGCuseKFParticle and kGCuseESDTracks can be true at the same time") - } + } if(kGCuseKFParticle == kFALSE && kGCuseESDTrack == kFALSE){ //Print warning, one have to be specified ::Error("ConfigGammaConversion","Both kGCuseKFParticle and kGCuseESDTracks can be false at the same time") - } - + } + + //Create the V0Reader AliV0Reader * v0Reader = new AliV0Reader(); if(kGCuseKFParticle){ @@ -889,6 +959,7 @@ AliAnalysisTaskGammaConversion* ConfigGammaConversion(TString arguments,AliAnaly v0Reader->SetPtCut(kGCptCut); v0Reader->SetLineCutZRSlope(kGCLineCutZRSlope); v0Reader->SetLineCutZValue(kGCLineCutZValue); + v0Reader->SetMaxZCut(kGCmaxZCut); v0Reader->SetChi2CutConversion(kGCchi2CutConversion); v0Reader->SetChi2CutMeson(kGCchi2CutMeson); v0Reader->SetPIDProbability(kGCprobElectron); @@ -898,6 +969,10 @@ AliAnalysisTaskGammaConversion* ConfigGammaConversion(TString arguments,AliAnaly v0Reader->SetSigmaMass(kGCsigmaCutGammaMass); v0Reader->SetUseImprovedVertex(kGCuseImprovedVertex); v0Reader->SetDoMCTruth(kGCdoMCTruth); + v0Reader->SetUseOwnXYZCalculation(kGCdoOwnXYZCalculation); + // for CF + v0Reader->SetCFManager(man); + // Create the GammaConversionTask AliAnalysisTaskGammaConversion *gammaconversion = new AliAnalysisTaskGammaConversion("GammaConversionTask"); @@ -917,7 +992,7 @@ AliAnalysisTaskGammaConversion* ConfigGammaConversion(TString arguments,AliAnaly gammaconversion->SetGammaWidth(kGCgammaWidth); gammaconversion->SetPi0Width(kGCpi0Width); gammaconversion->SetEtaWidth(kGCetaWidth); - + gammaconversion->SetMinOpeningAngleGhostCut(kGCminOpeningAngleGhostCut); // define the width constraint used by KF particle. @@ -929,35 +1004,39 @@ AliAnalysisTaskGammaConversion* ConfigGammaConversion(TString arguments,AliAnaly v0Reader->SetHistograms(histograms);// also give the pointer to the v0reader, for debugging cuts gammaconversion->SetDoMCTruth(kGCdoMCTruth); - + gammaconversion->SetDoNeutralMeson(kGCrunNeutralMeson); gammaconversion->SetDoJet(kGCrunJet); gammaconversion->SetDoChic(kGCrunChic); - + // for CF + gammaconversion->SetCFManager(man); + // Add task to the manager mgr->AddTask(gammaconversion); // Connect I/O to the task mgr->ConnectInput (gammaconversion, 0, cinput1); - - if(kGCrunOnTrain == kFALSE){ + + if(fWriteStandardAOD){ mgr->ConnectOutput(gammaconversion, 0, coutput1); - mgr->ConnectOutput(gammaconversion, 1, coutput2); } + mgr->ConnectOutput(gammaconversion, 1, coutput2); + mgr->ConnectOutput(gammaconversion, 2, coutput3); + if(kGCrunOnTrain == kFALSE){ if(kGCdataList.IsNull()){ cout<<"Data list is not set, aborting."<InitAnalysis(); - + mgr->PrintStatus(); - + mgr->StartAnalysis("local",chain); } } @@ -965,7 +1044,7 @@ AliAnalysisTaskGammaConversion* ConfigGammaConversion(TString arguments,AliAnaly } void build() { - + TStopwatch timer; timer.Start(); gSystem->Load("libTree.so"); @@ -1004,10 +1083,17 @@ void build() { //// //Setting up ANALYSISalice.par// //// - cout<<"compiling ANALUSISalice"<Load("libANALYSISalice.so"); + //// + //Setting up CORRFW.par// + //// + cout<<"compiling CORRFW"<Load("CORRFW.so"); + //// //Setting up PWG4GammaConv.par// //// @@ -1059,47 +1145,47 @@ void AddHistograms(AliGammaConversionHistograms *histograms){ if (kGCplotdPhiHdrGam == kTRUE){ histograms->AddHistogram("ESD_dphiHdrGam","ESD_dphiHdrGam", kGCnXBinsdphiHdrGam,kGCfirstXBindphiHdrGam,kGClastXBindphiHdrGam,"dphiHdrGam (rad)","Counts"); } - + if (kGCplotdPhiHdrGamIsolated == kTRUE){ histograms->AddHistogram("ESD_dphiHdrGamIsolated","ESD_dphiHdrGamIsolated", kGCnXBinsdphiHdrGam,kGCfirstXBindphiHdrGam,kGClastXBindphiHdrGam,"dphiHdrGamIsolated (rad)","Counts"); } - + if (kGCplotMinimumIsoDistance == kTRUE){ histograms->AddHistogram("ESD_MinimumIsoDistance","ESD_MinimumIsoDistance", kGCnXBinsMinimumIsoDistance,kGCfirstXBinMinimumIsoDistance,kGClastXBinMinimumIsoDistance,"Minimum Iso Distance (rad)","Counts"); } - + if (kGCplotFFzHdrGam == kTRUE){ histograms->AddHistogram("ESD_FFzHdrGam","ESD_FFzHdrGam", kGCnXBinsFFzHdrGam, kGCfirstXBinFFzHdrGam,kGClastXBinFFzHdrGam,"FFz Hdr Gam","Counts"); } - + if (kGCplotImbalanceHdrGam == kTRUE){ histograms->AddHistogram("ESD_ImbalanceHdrGam","ESD_ImbalanceHdrGam", kGCnXBinsImbalanceHdrGam, kGCfirstXBinImbalanceHdrGam,kGClastXBinImbalanceHdrGam,"Imbalance Hdr Gam","Counts"); } }//end if(kGCrunJet) - + //---------------------------------------------- Chi_c --------------------------------------------------------- if(kGCrunChic){ - + if(kGCplotESDInvMassePluseMinus == kTRUE){histograms->AddHistogram("ESD_InvMass_ePluseMinus","",kGCnXBinsJPsiMass, kGCfirstXBinJPsiMass, kGClastXBinJPsiMass, "", - "");} + "");} if(kGCplotESDInvMassePluseMinus == kTRUE){histograms->AddHistogram("ESD_InvMass_ePluseMinusTest","",kGCnXBinsJPsiMass, kGCfirstXBinJPsiMass, kGClastXBinJPsiMass, - "","");} + "","");} if(kGCplotESDInvMassePluseMinus == kTRUE){histograms->AddHistogram("ESD_InvMass_xPlusxMinus","",kGCnXBinsJPsiMass, kGCfirstXBinJPsiMass, kGClastXBinJPsiMass, "", - "");} + "");} if(kGCplotESDElectronPosNegPt == kTRUE){histograms->AddHistogram("ESD_ElectronPosNegPt","",kGCnXBinsEPosNegPt,kGCfirstXBinEPosNegPt,kGClastXBinEPosNegPt,"","");} if(kGCplotESDElectronPosNegEta == kTRUE){histograms->AddHistogram("ESD_ElectronPosNegEta","",kGCnXBinsEPosNegEta,kGCfirstXBinEPosNegEta,kGClastXBinEPosNegEta,""," -");} - + ");} + if(kGCplotESDElectronPosNegPt == kTRUE){histograms->AddHistogram("ESD_ElectronPosPt","",kGCnXBinsEPosNegPt,kGCfirstXBinEPosNegPt,kGClastXBinEPosNegPt,"","");} if(kGCplotESDElectronPosNegPt == kTRUE){histograms->AddHistogram("ESD_ElectronNegPt","",kGCnXBinsEPosNegPt,kGCfirstXBinEPosNegPt,kGClastXBinEPosNegPt,"","");} - + if(kGCplotESDElectronPosNegAngle == kTRUE){histograms->AddHistogram("ESD_ElectronPosNegJPsiAngle","",kGCnXBinsEPosNegAngle,kGCfirstXBinEPosNegAngle,kGClastXBinEPo - sNegAngle,"","");} + sNegAngle,"","");} if(kGCplotMCElectronPosNegPt == kTRUE){histograms->AddHistogram("MC_ElectronPosNegPt","",kGCnXBinsEPosNegPt,kGCfirstXBinEPosNegPt,kGClastXBinEPosNegPt,"","");} if(kGCplotMCElectronPosNegEta == kTRUE){histograms->AddHistogram("MC_ElectronPosNegEta","",kGCnXBinsEPosNegEta,kGCfirstXBinEPosNegEta,kGClastXBinEPosNegEta,"","") ;} if(kGCplotMCElectronPosNegJPsiAngle == kTRUE){histograms->AddHistogram("MC_ElectronPosNegJPsiAngle","",kGCnXBinsEPosNegAngle,kGCfirstXBinEPosNegAngle,kGClastXBinE - PosNegAngle,"","");} + PosNegAngle,"","");} if(kGCplotESDePoseNegAngle == kTRUE){histograms->AddHistogram("ESD_eNegePosAngleBeforeCut","",kGCnXBinsEPosNegAngle,kGCfirstXBinEPosNegAngle,kGClastXBinEPosNegAngle,"","");} if(kGCplotESDePoseNegAngle == kTRUE){histograms->AddHistogram("ESD_eNegePosAngleAfterCut","",kGCnXBinsEPosNegAngle,kGCfirstXBinEPosNegAngle,kGClastXBinEPosNegAngle,"","");} if(kGCplotESDInvMassGammaePluseMinusChiC == kTRUE) {histograms->AddHistogram("ESD_InvMass_GammaePluseMinusChiC","",kGCnXBinsChicMass,kGCfirstXBinChicMass,kGClastXBinChicMass,"","");} @@ -1107,64 +1193,64 @@ void AddHistograms(AliGammaConversionHistograms *histograms){ if(kGCplotESDInvMassGammaePluseMinusPi0 == kTRUE) {histograms->AddHistogram("ESD_InvMass_GammaePluseMinusPi0","",kGCnXBinsPi0Mass,kGCfirstXBinPi0Mass,kGClastXBinPi0Mass,"","");} if(kGCplotESDElectronPosNegPi0Angle == kTRUE){histograms->AddHistogram("ESD_ElectronPosNegPi0Angle","",kGCnXBinsEPosNegAngle,kGCfirstXBinEPosNegAngle,kGClastXBinEPosNegAngle,"","");} if(kGCplotMCElectronPosNegPi0Angle == kTRUE){histograms->AddHistogram("MC_ElectronPosNegPi0Angle","",kGCnXBinsEPosNegAngle,kGCfirstXBinEPosNegAngle,kGClastXBinEPosNegAngle,"","");} - + if(kGCplotESDEPosBackground == kTRUE){histograms->AddHistogram("ESD_EPosBackground","",kGCnXBinsEBackground,kGCfirstXBinEBackground,kGClastXBinEBackground,"","");} - + if(kGCplotESDEPosBackground == kTRUE){histograms->AddHistogram("ESD_EPosENegNoJPsiBG","",kGCnXBinsEBackground,kGCfirstXBinEBackground,kGClastXBinEBackground,"","");} - - + + if(kGCplotESDENegBackground == kTRUE){histograms->AddHistogram("ESD_ENegBackground","",kGCnXBinsEBackground,kGCfirstXBinEBackground,kGClastXBinEBackground,"","");} if(kGCplotESDEPosENegBackground == kTRUE){histograms->AddHistogram("ESD_EPosENegBackground","",kGCnXBinsEBackground,kGCfirstXBinEBackground,kGClastXBinEBackground,"","");} if(kGCplotESDEPosENegBackgroundCut == kTRUE){histograms->AddHistogram("ESD_EPosENegBackgroundCut","",kGCnXBinsEBackgroundCut,kGCfirstXBinEBackgroundCut,kGClastXBinEBackgroundCut,"","");} - + if(kGCplotESDEPosENegGammaBackgroundMX == kTRUE){histograms->AddHistogram("ESD_EPosENegGammaBackgroundMX","",kGCnXBinsEBackground,kGCfirstXBinEBackground,kGClastXBinEBackground,"","");} if(kGCplotESDEPosENegGammaBackgroundMX == kTRUE){histograms->AddHistogram("ESD_EPosENegGammaBackgroundMXDiff","",kGCnXBinsEBackground,kGCfirstXBinEBackground,kGClastXBinEBackground,"","");} - + if(kGCplotTableElectrons == kTRUE){ histograms->AddTable("Table_Electrons","",kGCnElementsElectronTable,kGCelectronTable);} }// end kGCrunChic - + //---------------------------------------------- Neutral Meson --------------------------------------------------------- if(kGCrunNeutralMeson){ if(kGCplotMCConversionR == kTRUE){ histograms->AddHistogram("MC_Conversion_R","Radius of gamma conversion points",kGCnXBinsR, kGCfirstXBinR, kGClastXBinR,"counts","cm");} if(kGCplotMCConversionZR == kTRUE){ histograms->AddHistogram("MC_Conversion_ZR","Radius of gamma conversion points vs Z",kGCnXBinsZR, kGCfirstXBinZR, kGClastXBinZR, kGCnYBinsZR, kGCfirstYBinZR, kGClastYBinZR, "cm", "cm");} if(kGCplotMCConversionXY == kTRUE){ histograms->AddHistogram("MC_Conversion_XY","Gamma XY converison point.",kGCnXBinsXY, kGCfirstXBinXY, kGClastXBinXY, kGCnYBinsXY, kGCfirstYBinXY, kGClastYBinXY, "cm", "cm");} if(kGCplotMCConversionOpeningAngle == kTRUE){ histograms->AddHistogram("MC_Conversion_OpeningAngle","Opening angle of e+e- pairs from gamma conversion",kGCnXBinsOpeningAngle, kGCfirstXBinOpeningAngle, kGClastXBinOpeningAngle, "counts", "cm");} - + if(kGCplotMCEEnergy == kTRUE){ histograms->AddHistogram("MC_E_Energy" ,"" , kGCnXBinsEnergy, kGCfirstXBinEnergy, kGClastXBinEnergy, "", "");} if(kGCplotMCEPt == kTRUE){ histograms->AddHistogram("MC_E_Pt" ,"" , kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt, "", "");} if(kGCplotMCEEta == kTRUE){ histograms->AddHistogram("MC_E_Eta" ,"" , kGCnXBinsEta, kGCfirstXBinEta, kGClastXBinEta, "", "");} if(kGCplotMCEPhi == kTRUE){ histograms->AddHistogram("MC_E_Phi" ,"" , kGCnXBinsPhi, kGCfirstXBinPhi, kGClastXBinPhi, "", "");} - + if(kGCplotMCPEnergy == kTRUE){ histograms->AddHistogram("MC_P_Energy" ,"" , kGCnXBinsEnergy, kGCfirstXBinEnergy, kGClastXBinEnergy, "", "");} if(kGCplotMCPPt == kTRUE){ histograms->AddHistogram("MC_P_Pt" ,"" , kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt, "", "");} if(kGCplotMCPEta == kTRUE){ histograms->AddHistogram("MC_P_Eta" ,"" , kGCnXBinsEta, kGCfirstXBinEta, kGClastXBinEta, "", "");} if(kGCplotMCPPhi == kTRUE){ histograms->AddHistogram("MC_P_Phi" ,"" , kGCnXBinsPhi, kGCfirstXBinPhi, kGClastXBinPhi, "", "");} - + if(kGCplotMCallGammaEnergy == kTRUE){ histograms->AddHistogram("MC_allGamma_Energy" ,"" , kGCnXBinsEnergy, kGCfirstXBinEnergy, kGClastXBinEnergy, "", "");} if(kGCplotMCallGammaPt == kTRUE){ histograms->AddHistogram("MC_allGamma_Pt" ,"" , kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt, "", "");} if(kGCplotMCallGammaEta == kTRUE){ histograms->AddHistogram("MC_allGamma_Eta" ,"" , kGCnXBinsEta, kGCfirstXBinEta, kGClastXBinEta, "", "");} if(kGCplotMCallGammaPhi == kTRUE){ histograms->AddHistogram("MC_allGamma_Phi" ,"" , kGCnXBinsPhi, kGCfirstXBinPhi, kGClastXBinPhi, "", "");} if(kGCplotMCallGammaRapid == kTRUE){ histograms->AddHistogram("MC_allGamma_Rapid" ,"" , kGCnXBinsRapid, kGCfirstXBinRapid, kGClastXBinRapid, "", "");} - + if(kGCplotMCConvGammaEnergy == kTRUE){ histograms->AddHistogram("MC_ConvGamma_Energy" ,"" , kGCnXBinsEnergy, kGCfirstXBinEnergy, kGClastXBinEnergy, "", "");} if(kGCplotMCConvGammaPt == kTRUE){ histograms->AddHistogram("MC_ConvGamma_Pt" ,"" , kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt, "", "");} if(kGCplotMCConvGammaEta == kTRUE){ histograms->AddHistogram("MC_ConvGamma_Eta" ,"" , kGCnXBinsEta, kGCfirstXBinEta, kGClastXBinEta, "", "");} if(kGCplotMCConvGammaPhi == kTRUE){ histograms->AddHistogram("MC_ConvGamma_Phi" ,"" , kGCnXBinsPhi, kGCfirstXBinPhi, kGClastXBinPhi, "", "");} if(kGCplotMCConvGammaRapid == kTRUE){ histograms->AddHistogram("MC_ConvGamma_Rapid" ,"" , kGCnXBinsRapid, kGCfirstXBinRapid, kGClastXBinRapid, "", "");} if(kGCplotMCConvGammaPtvsEta == kTRUE){ histograms->AddHistogram("MC_ConvGamma_Pt_Eta","", kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt,kGCnXBinsEta, kGCfirstXBinEta, kGClastXBinEta,"","");} - + if(kGCplotMCallDirectGammaEnergy == kTRUE){ histograms->AddHistogram("MC_allDirectGamma_Energy" ,"" , kGCnXBinsEnergy, kGCfirstXBinEnergy, kGClastXBinEnergy, "", "");} if(kGCplotMCallDirectGammaPt == kTRUE){ histograms->AddHistogram("MC_allDirectGamma_Pt" ,"" , kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt, "", "");} if(kGCplotMCallDirectGammaEta == kTRUE){ histograms->AddHistogram("MC_allDirectGamma_Eta" ,"" , kGCnXBinsEta, kGCfirstXBinEta, kGClastXBinEta, "", "");} if(kGCplotMCallDirectGammaPhi == kTRUE){ histograms->AddHistogram("MC_allDirectGamma_Phi" ,"" , kGCnXBinsPhi, kGCfirstXBinPhi, kGClastXBinPhi, "", "");} if(kGCplotMCallDirectGammaRapid == kTRUE){ histograms->AddHistogram("MC_allDirectGamma_Rapid" ,"" , kGCnXBinsRapid, kGCfirstXBinRapid, kGClastXBinRapid, "", "");} - + if(kGCplotMCConvDirectGammaEnergy == kTRUE){ histograms->AddHistogram("MC_ConvDirectGamma_Energy" ,"" , kGCnXBinsEnergy, kGCfirstXBinEnergy, kGClastXBinEnergy, "", "");} if(kGCplotMCConvDirectGammaPt == kTRUE){ histograms->AddHistogram("MC_ConvDirectGamma_Pt" ,"" , kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt, "", "");} if(kGCplotMCConvDirectGammaEta == kTRUE){ histograms->AddHistogram("MC_ConvDirectGamma_Eta" ,"" , kGCnXBinsEta, kGCfirstXBinEta, kGClastXBinEta, "", "");} if(kGCplotMCConvDirectGammaPhi == kTRUE){ histograms->AddHistogram("MC_ConvDirectGamma_Phi" ,"" , kGCnXBinsPhi, kGCfirstXBinPhi, kGClastXBinPhi, "", "");} if(kGCplotMCConvDirectGammaRapid == kTRUE){ histograms->AddHistogram("MC_ConvDirectGamma_Rapid" ,"" , kGCnXBinsRapid, kGCfirstXBinRapid, kGClastXBinRapid, "", "");} - + if(kGCplotMCMotherEta == kTRUE){ histograms->AddHistogram("MC_Mother_Eta" ,"" , kGCnXBinsEta, kGCfirstXBinEta, kGClastXBinEta, "", "");} if(kGCplotMCMotherPhi == kTRUE){ histograms->AddHistogram("MC_Mother_Phi" ,"" , kGCnXBinsPhi, kGCfirstXBinPhi, kGClastXBinPhi, "", "");} if(kGCplotMCMotherRapid == kTRUE){ histograms->AddHistogram("MC_Mother_Rapid" ,"" , kGCnXBinsRapid, kGCfirstXBinRapid, kGClastXBinRapid, "", "");} @@ -1179,14 +1265,14 @@ void AddHistograms(AliGammaConversionHistograms *histograms){ if(kGCplotMCMotherPtvsRapidWithinAcceptance == kTRUE){ histograms->AddHistogram("MC_Mother_Pt_Rapid_withinAcceptance" ,"" , kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt, kGCnXBinsRapid, kGCfirstXBinRapid, kGClastXBinRapid, "", "");} if(kGCplotMCMotherPtvsEtaConvGammaWithinAcceptance == kTRUE){ histograms->AddHistogram("MC_Mother_Pt_Eta_ConvGamma_withinAcceptance" ,"" , kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt, kGCnXBinsEta, kGCfirstXBinEta, kGClastXBinEta, "", "");} if(kGCplotMCMotherPtvsRapidConvGammaWithinAcceptance == kTRUE){ histograms->AddHistogram("MC_Mother_Pt_Rapid_ConvGamma_withinAcceptance" ,"" , kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt, kGCnXBinsRapid, kGCfirstXBinRapid, kGClastXBinRapid, "", "");} - + if(kGCplotMCMotherSpectra == kTRUE){ histograms->AddHistogram("MC_Mother_InvMass_vs_Pt" ,"" ,kGCnXBinsSpectra, kGCfirstXBinSpectra, kGClastXBinSpectra, kGCnYBinsSpectra, kGCfirstYBinSpectra, kGClastYBinSpectra, "", ""); histograms->AddHistogram("MC_Mother_InvMass_vs_Pt_withinAcceptance" ,"" ,kGCnXBinsSpectra, kGCfirstXBinSpectra, kGClastXBinSpectra, kGCnYBinsSpectra, kGCfirstYBinSpectra, kGClastYBinSpectra, "", ""); histograms->AddHistogram("MC_Mother_InvMass_vs_Pt_ConvGamma_withinAcceptance" ,"" ,kGCnXBinsSpectra, kGCfirstXBinSpectra, kGClastXBinSpectra, kGCnYBinsSpectra, kGCfirstYBinSpectra, kGClastYBinSpectra, "", ""); } - - + + if(kGCplotMCPi0Eta == kTRUE){ histograms->AddHistogram("MC_Pi0_Eta" ,"" , kGCnXBinsEta, kGCfirstXBinEta, kGClastXBinEta, "", "");} if(kGCplotMCPi0Rapid == kTRUE){ histograms->AddHistogram("MC_Pi0_Rapid" ,"" , kGCnXBinsRapid, kGCfirstXBinRapid, kGClastXBinRapid, "", "");} if(kGCplotMCPi0Phi == kTRUE){ histograms->AddHistogram("MC_Pi0_Phi" ,"" , kGCnXBinsPhi, kGCfirstXBinPhi, kGClastXBinPhi, "", "");} @@ -1202,8 +1288,8 @@ void AddHistograms(AliGammaConversionHistograms *histograms){ if(kGCplotMCPi0PtvsEtaConvGammaWithinAcceptance == kTRUE){ histograms->AddHistogram("MC_Pi0_Pt_Eta_ConvGamma_withinAcceptance" ,"" , kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt, kGCnXBinsEta, kGCfirstXBinEta, kGClastXBinEta, "", "");} if(kGCplotMCPi0PtvsRapidConvGammaWithinAcceptance == kTRUE){ histograms->AddHistogram("MC_Pi0_Pt_Rapid_ConvGamma_withinAcceptance" ,"" , kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt, kGCnXBinsRapid, kGCfirstXBinRapid, kGClastXBinRapid, "", "");} if(kGCplotMCPi0ZRConvGammaWithinAcceptance == kTRUE){ histograms->AddHistogram("MC_Pi0_ZR_ConvGamma_withinAcceptance" ,"" , kGCnXBinsZR, kGCfirstXBinZR, kGClastXBinZR, kGCnYBinsZR, kGCfirstYBinZR, kGClastYBinZR, "", "");} - - + + if(kGCplotMCPi0SecondaryEta == kTRUE){ histograms->AddHistogram("MC_Pi0_Secondaries_Eta" ,"" , kGCnXBinsEta, kGCfirstXBinEta, kGClastXBinEta, "", "");} if(kGCplotMCPi0SecondaryRapid == kTRUE){ histograms->AddHistogram("MC_Pi0_Secondaries_Rapid" ,"" , kGCnXBinsRapid, kGCfirstXBinRapid, kGClastXBinRapid, "", "");} if(kGCplotMCPi0SecondaryPhi == kTRUE){ histograms->AddHistogram("MC_Pi0_Secondaries_Phi" ,"" , kGCnXBinsPhi, kGCfirstXBinPhi, kGClastXBinPhi, "", "");} @@ -1218,9 +1304,9 @@ void AddHistograms(AliGammaConversionHistograms *histograms){ if(kGCplotMCPi0SecondaryPtvsRapidWithinAcceptance == kTRUE){ histograms->AddHistogram("MC_Pi0_Secondaries_Pt_Rapid_withinAcceptance" ,"" , kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt, kGCnXBinsRapid, kGCfirstXBinRapid, kGClastXBinRapid, "", "");} if(kGCplotMCPi0SecondaryPtvsEtaConvGammaWithinAcceptance == kTRUE){ histograms->AddHistogram("MC_Pi0_Secondaries_Pt_Eta_ConvGamma_withinAcceptance" ,"" , kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt, kGCnXBinsEta, kGCfirstXBinEta, kGClastXBinEta, "", "");} if(kGCplotMCPi0SecondaryPtvsRapidConvGammaWithinAcceptance == kTRUE){ histograms->AddHistogram("MC_Pi0_Secondaries_Pt_Rapid_ConvGamma_withinAcceptance" ,"" , kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt, kGCnXBinsRapid, kGCfirstXBinRapid, kGClastXBinRapid, "", "");} - - - + + + if(kGCplotMCEtaEta == kTRUE){ histograms->AddHistogram("MC_Eta_Eta" ,"" , kGCnXBinsEta, kGCfirstXBinEta, kGClastXBinEta, "", "");} if(kGCplotMCEtaRapid == kTRUE){ histograms->AddHistogram("MC_Eta_Rapid" ,"" , kGCnXBinsRapid, kGCfirstXBinRapid, kGClastXBinRapid, "", "");} if(kGCplotMCEtaPhi == kTRUE){ histograms->AddHistogram("MC_Eta_Phi" ,"" , kGCnXBinsPhi, kGCfirstXBinPhi, kGClastXBinPhi, "", "");} @@ -1236,19 +1322,19 @@ void AddHistograms(AliGammaConversionHistograms *histograms){ if(kGCplotMCEtaPtvsEtaConvGammaWithinAcceptance == kTRUE){ histograms->AddHistogram("MC_Eta_Pt_Eta_ConvGamma_withinAcceptance" ,"" , kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt, kGCnXBinsEta, kGCfirstXBinEta, kGClastXBinEta, "", "");} if(kGCplotMCEtaPtvsRapidConvGammaWithinAcceptance == kTRUE){ histograms->AddHistogram("MC_Eta_Pt_Rapid_ConvGamma_withinAcceptance" ,"" , kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt, kGCnXBinsRapid, kGCfirstXBinRapid, kGClastXBinRapid, "", "");} if(kGCplotMCEtaZRConvGammaWithinAcceptance == kTRUE){ histograms->AddHistogram("MC_Eta_ZR_ConvGamma_withinAcceptance" ,"" , kGCnXBinsZR, kGCfirstXBinZR, kGClastXBinZR, kGCnYBinsZR, kGCfirstYBinZR, kGClastYBinZR, "", "");} - - + + // Histograms from esd tracks if(kGCplotESDEEnergy == kTRUE){ histograms->AddHistogram("ESD_E_Energy" ,"" , kGCnXBinsEnergy, kGCfirstXBinEnergy, kGClastXBinEnergy, "", "");} if(kGCplotESDEPt == kTRUE){ histograms->AddHistogram("ESD_E_Pt" ,"" , kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt, "", "");} if(kGCplotESDEEta == kTRUE){ histograms->AddHistogram("ESD_E_Eta" ,"" , kGCnXBinsEta, kGCfirstXBinEta, kGClastXBinEta, "", "");} if(kGCplotESDEPhi == kTRUE){ histograms->AddHistogram("ESD_E_Phi" ,"" , kGCnXBinsPhi, kGCfirstXBinPhi, kGClastXBinPhi, "", "");} - + if(kGCplotESDPEnergy == kTRUE){ histograms->AddHistogram("ESD_P_Energy" ,"" , kGCnXBinsEnergy, kGCfirstXBinEnergy, kGClastXBinEnergy, "", "");} if(kGCplotESDPPt == kTRUE){ histograms->AddHistogram("ESD_P_Pt" ,"" , kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt, "", "");} if(kGCplotESDPEta == kTRUE){ histograms->AddHistogram("ESD_P_Eta" ,"" , kGCnXBinsEta, kGCfirstXBinEta, kGClastXBinEta, "", "");} if(kGCplotESDPPhi == kTRUE){ histograms->AddHistogram("ESD_P_Phi" ,"" , kGCnXBinsPhi, kGCfirstXBinPhi, kGClastXBinPhi, "", "");} - + if(kGCplotESDConvGammaEnergy == kTRUE){ histograms->AddHistogram("ESD_ConvGamma_Energy" ,"" , kGCnXBinsEnergy, kGCfirstXBinEnergy, kGClastXBinEnergy, "", "");} if(kGCplotESDConvGammaPt == kTRUE){ histograms->AddHistogram("ESD_ConvGamma_Pt" ,"" , kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt, "", "");} if(kGCplotESDConvGammaEta == kTRUE){ histograms->AddHistogram("ESD_ConvGamma_Eta" ,"" , kGCnXBinsEta, kGCfirstXBinEta, kGClastXBinEta, "", "");} @@ -1261,15 +1347,15 @@ void AddHistograms(AliGammaConversionHistograms *histograms){ if(kGCplotESDConvGammaPtvsEta == kTRUE){ histograms->AddHistogram("ESD_ConvGamma_Pt_Eta","", kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt,kGCnXBinsEta, kGCfirstXBinEta, kGClastXBinEta,"","" );} if(kGCplotESDConvGammaPtvsChi2 == kTRUE){ histograms->AddHistogram("ESD_ConvGamma_Pt_Chi2" ,"" ,kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt, kGCnXBinsGammaChi2, kGCfirstXBinGammaChi2, kGClastXBinGammaChi2, "", "");} if(kGCplotESDConvGammaEtavsChi2 == kTRUE){ histograms->AddHistogram("ESD_ConvGamma_Eta_Chi2" ,"" ,kGCnXBinsEta, kGCfirstXBinEta, kGClastXBinEta, kGCnXBinsGammaChi2, kGCfirstXBinGammaChi2, kGClastXBinGammaChi2, "", "");} - - - + + + if(kGCplotESDConversionR == kTRUE){ histograms->AddHistogram("ESD_Conversion_R" ,"" , kGCnXBinsR, kGCfirstXBinR, kGClastXBinR, "", "");} if(kGCplotESDConversionZR == kTRUE){ histograms->AddHistogram("ESD_Conversion_ZR" ,"" , kGCnXBinsZR, kGCfirstXBinZR, kGClastXBinZR, kGCnYBinsZR, kGCfirstYBinZR, kGClastYBinZR, "", "");} if(kGCplotESDConversionXY == kTRUE){ histograms->AddHistogram("ESD_Conversion_XY" ,"" , kGCnXBinsXY, kGCfirstXBinXY, kGClastXBinXY, kGCnYBinsXY, kGCfirstYBinXY, kGClastYBinXY, "", "");} if(kGCplotESDConversionOpeningAngle == kTRUE){ histograms->AddHistogram("ESD_Conversion_OpeningAngle" ,"" , kGCnXBinsOpeningAngle, kGCfirstXBinOpeningAngle, kGClastXBinOpeningAngle, "", "");} - - + + if(kGCplotESDTrueConvGammaEnergy == kTRUE){ histograms->AddHistogram("ESD_TrueConvGamma_Energy" ,"" , kGCnXBinsEnergy, kGCfirstXBinEnergy, kGClastXBinEnergy, "", "");} if(kGCplotESDTrueConvGammaPt == kTRUE){ histograms->AddHistogram("ESD_TrueConvGamma_Pt" ,"" , kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt, "", "");} if(kGCplotESDTrueConvGammaEta == kTRUE){ histograms->AddHistogram("ESD_TrueConvGamma_Eta" ,"" , kGCnXBinsEta, kGCfirstXBinEta, kGClastXBinEta, "", "");} @@ -1282,18 +1368,18 @@ void AddHistograms(AliGammaConversionHistograms *histograms){ if(kGCplotESDTrueConvGammaPtvsEta == kTRUE){ histograms->AddHistogram("ESD_TrueConvGamma_Pt_Eta" ,"" , kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt,kGCnXBinsEta, kGCfirstXBinEta, kGClastXBinEta, "", "");} if(kGCplotESDTrueConvGammaPtvsChi2 == kTRUE){ histograms->AddHistogram("ESD_TrueConvGamma_Pt_Chi2" ,"" ,kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt, kGCnXBinsGammaChi2, kGCfirstXBinGammaChi2, kGClastXBinGammaChi2, "", "");} if(kGCplotESDTrueConvGammaEtavsChi2 == kTRUE){ histograms->AddHistogram("ESD_TrueConvGamma_Eta_Chi2" ,"" ,kGCnXBinsEta, kGCfirstXBinEta, kGClastXBinEta, kGCnXBinsGammaChi2, kGCfirstXBinGammaChi2, kGClastXBinGammaChi2, "", "");} - + if(kGCplotESDTrueConversionR == kTRUE){ histograms->AddHistogram("ESD_TrueConversion_R" ,"" , kGCnXBinsR, kGCfirstXBinR, kGClastXBinR, "", "");} if(kGCplotESDTrueConversionZR == kTRUE){ histograms->AddHistogram("ESD_TrueConversion_ZR" ,"" , kGCnXBinsZR, kGCfirstXBinZR, kGClastXBinZR, kGCnYBinsZR, kGCfirstYBinZR, kGClastYBinZR, "", "");} if(kGCplotESDTrueConversionXY == kTRUE){ histograms->AddHistogram("ESD_TrueConversion_XY" ,"" , kGCnXBinsXY, kGCfirstXBinXY, kGClastXBinXY, kGCnYBinsXY, kGCfirstYBinXY, kGClastYBinXY, "", "");} if(kGCplotESDTrueConversionOpeningAngle == kTRUE){ histograms->AddHistogram("ESD_TrueConversion_OpeningAngle" ,"" , kGCnXBinsOpeningAngle, kGCfirstXBinOpeningAngle, kGClastXBinOpeningAngle, "", "");} - + if(kGCplotESDTrueConvGammaMCPtEta == kTRUE){ histograms->AddHistogram("ESD_TrueConvGamma_MC_Pt_Eta" ,"" , kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt, kGCnXBinsEta, kGCfirstXBinEta, kGClastXBinEta, "", "");} if(kGCplotESDTrueConversionMCZR == kTRUE){ histograms->AddHistogram("ESD_TrueConversion_MC_ZR" ,"" , kGCnXBinsZR, kGCfirstXBinZR, kGClastXBinZR, kGCnYBinsZR, kGCfirstYBinZR, kGClastYBinZR, "", "");} if(kGCplotESDTrueConversionMCXY == kTRUE){ histograms->AddHistogram("ESD_TrueConversion_MC_XY" ,"" , kGCnXBinsXY, kGCfirstXBinXY, kGClastXBinXY, kGCnYBinsXY, kGCfirstYBinXY, kGClastYBinXY, "", "");} - - - + + + if(kGCplotESDNoCutConvGammaEnergy == kTRUE){ histograms->AddHistogram("ESD_NoCutConvGamma_Energy" ,"" , kGCnXBinsEnergy, kGCfirstXBinEnergy, kGClastXBinEnergy, "", "");} if(kGCplotESDNoCutConvGammaPt == kTRUE){ histograms->AddHistogram("ESD_NoCutConvGamma_Pt" ,"" , kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt, "", "");} if(kGCplotESDNoCutConvGammaEta == kTRUE){ histograms->AddHistogram("ESD_NoCutConvGamma_Eta" ,"" , kGCnXBinsEta, kGCfirstXBinEta, kGClastXBinEta, "", "");} @@ -1306,18 +1392,18 @@ void AddHistograms(AliGammaConversionHistograms *histograms){ if(kGCplotESDNoCutConvGammaPtvsEta == kTRUE){ histograms->AddHistogram("ESD_NoCutConvGamma_Pt_Eta" ,"" , kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt,kGCnXBinsEta, kGCfirstXBinEta, kGClastXBinEta, "", "");} if(kGCplotESDNoCutConvGammaPtvsChi2 == kTRUE){ histograms->AddHistogram("ESD_NoCutConvGamma_Pt_Chi2" ,"" ,kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt, kGCnXBinsGammaChi2, kGCfirstXBinGammaChi2, kGClastXBinGammaChi2, "", "");} if(kGCplotESDNoCutConvGammaEtavsChi2 == kTRUE){ histograms->AddHistogram("ESD_NoCutConvGamma_Eta_Chi2" ,"" ,kGCnXBinsEta, kGCfirstXBinEta, kGClastXBinEta, kGCnXBinsGammaChi2, kGCfirstXBinGammaChi2, kGClastXBinGammaChi2, "", "");} - + if(kGCplotESDNoCutConversionR == kTRUE){ histograms->AddHistogram("ESD_NoCutConversion_R" ,"" , kGCnXBinsR, kGCfirstXBinR, kGClastXBinR, "", "");} if(kGCplotESDNoCutConversionZR == kTRUE){ histograms->AddHistogram("ESD_NoCutConversion_ZR" ,"" , kGCnXBinsZR, kGCfirstXBinZR, kGClastXBinZR, kGCnYBinsZR, kGCfirstYBinZR, kGClastYBinZR, "", "");} if(kGCplotESDNoCutConversionXY == kTRUE){ histograms->AddHistogram("ESD_NoCutConversion_XY" ,"" , kGCnXBinsXY, kGCfirstXBinXY, kGClastXBinXY, kGCnYBinsXY, kGCfirstYBinXY, kGClastYBinXY, "", "");} if(kGCplotESDNoCutConversionOpeningAngle == kTRUE){ histograms->AddHistogram("ESD_NoCutConversion_OpeningAngle" ,"" , kGCnXBinsOpeningAngle, kGCfirstXBinOpeningAngle, kGClastXBinOpeningAngle, "", "");} - + if(kGCplotESDNoCutConvGammaMCPtEta == kTRUE){ histograms->AddHistogram("ESD_NoCutConvGamma_MC_Pt_Eta" ,"" , kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt, kGCnXBinsEta, kGCfirstXBinEta, kGClastXBinEta, "", "");} if(kGCplotESDNoCutConversionMCZR == kTRUE){ histograms->AddHistogram("ESD_NoCutConversion_MC_ZR" ,"" , kGCnXBinsZR, kGCfirstXBinZR, kGClastXBinZR, kGCnYBinsZR, kGCfirstYBinZR, kGClastYBinZR, "", "");} if(kGCplotESDNoCutConversionMCXY == kTRUE){ histograms->AddHistogram("ESD_NoCutConversion_MC_XY" ,"" , kGCnXBinsXY, kGCfirstXBinXY, kGClastXBinXY, kGCnYBinsXY, kGCfirstYBinXY, kGClastYBinXY, "", "");} - - - + + + if(kGCplotESDMotherOpeningAngleGamma == kTRUE){ histograms->AddHistogram("ESD_Mother_GammaDaughter_OpeningAngle" ,"" , kGCnXBinsOpeningAngle, kGCfirstXBinOpeningAngle, kGClastXBinOpeningAngle, "", "");} if(kGCplotESDMotherEnergy == kTRUE){ histograms->AddHistogram("ESD_Mother_Energy" ,"" , kGCnXBinsEnergy, kGCfirstXBinEnergy, kGClastXBinEnergy, "", "");} if(kGCplotESDMotherPt == kTRUE){ histograms->AddHistogram("ESD_Mother_Pt" ,"" , kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt, "", "");} @@ -1328,8 +1414,8 @@ void AddHistograms(AliGammaConversionHistograms *histograms){ if(kGCplotESDMotherZR == kTRUE){ histograms->AddHistogram("ESD_Mother_ZR" ,"" , kGCnXBinsZR, kGCfirstXBinZR, kGClastXBinZR, kGCnYBinsZR, kGCfirstYBinZR, kGClastYBinZR, "", "");} if(kGCplotESDMotherXY == kTRUE){ histograms->AddHistogram("ESD_Mother_XY" ,"" , kGCnXBinsXY, kGCfirstXBinXY, kGClastXBinXY, kGCnYBinsXY, kGCfirstYBinXY, kGClastYBinXY, "", "");} if(kGCplotESDMotherRapid == kTRUE){ histograms->AddHistogram("ESD_Mother_Rapid" ,"" , kGCnXBinsRapid, kGCfirstXBinRapid, kGClastXBinRapid, "", "");} - - + + if(kGCplotESDBackgroundOpeningAngleGamma == kTRUE){ histograms->AddHistogram("ESD_Background_GammaDaughter_OpeningAngle" ,"" , kGCnXBinsOpeningAngle, kGCfirstXBinOpeningAngle, kGClastXBinOpeningAngle, "", "");} if(kGCplotESDBackgroundEnergy == kTRUE){ histograms->AddHistogram("ESD_Background_Energy" ,"" , kGCnXBinsEnergy, kGCfirstXBinEnergy, kGClastXBinEnergy, "", "");} if(kGCplotESDBackgroundPt == kTRUE){ histograms->AddHistogram("ESD_Background_Pt" ,"" , kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt, "", "");} @@ -1340,30 +1426,30 @@ void AddHistograms(AliGammaConversionHistograms *histograms){ if(kGCplotESDBackgroundZR == kTRUE){ histograms->AddHistogram("ESD_Background_ZR" ,"" , kGCnXBinsZR, kGCfirstXBinZR, kGClastXBinZR, kGCnYBinsZR, kGCfirstYBinZR, kGClastYBinZR, "", "");} if(kGCplotESDBackgroundXY == kTRUE){ histograms->AddHistogram("ESD_Background_XY" ,"" , kGCnXBinsXY, kGCfirstXBinXY, kGClastXBinXY, kGCnYBinsXY, kGCfirstYBinXY, kGClastYBinXY, "", "");} if(kGCplotESDBackgroundRapid == kTRUE){ histograms->AddHistogram("ESD_Background_Rapid" ,"" , kGCnXBinsRapid, kGCfirstXBinRapid, kGClastXBinRapid, "", "");} - - + + if(kGCplotMapping == kTRUE){ - histograms->InitializeMappingValues(nPhiIndex,nRIndex,kGCnXBinsMapping,minRadius,maxRadius,kGCnYBinsMapping,minPhi,maxPhi); - histograms->AddMappingHistograms(nPhiIndex,nRIndex,kGCnXBinsMapping,minRadius,maxRadius,kGCnYBinsMapping,minPhi,maxPhi); + histograms->InitializeMappingValues(kGCnPhiIndex,kGCnRIndex,kGCnXBinsMapping,kGCminRadius,kGCmaxRadius,kGCnYBinsMapping,kGCminPhi,kGCmaxPhi); + histograms->AddMappingHistograms(kGCnPhiIndex,kGCnRIndex,kGCnXBinsMapping,kGCminRadius,kGCmaxRadius,kGCnYBinsMapping,kGCminPhi,kGCmaxPhi); } - + if(kGCplotResolutiondPt == kTRUE){histograms->AddHistogram("Resolution_dPt" ,"" , kGCnXBinsResdPt, kGCfirstXBinResdPt, kGClastXBinResdPt, kGCnYBinsResdPt, kGCfirstYBinResdPt, kGClastYBinResdPt, "", "");} if(kGCplotResolutiondR == kTRUE){histograms->AddHistogram("Resolution_dR" ,"" , kGCnXBinsResdR, kGCfirstXBinResdR, kGClastXBinResdR, kGCnYBinsResdR, kGCfirstYBinResdR, kGClastYBinResdR, "", "");} if(kGCplotResolutiondZ == kTRUE){histograms->AddHistogram("Resolution_dZ" ,"" , kGCnXBinsResdZ, kGCfirstXBinResdZ, kGClastXBinResdZ, kGCnYBinsResdZ, kGCfirstYBinResdZ, kGClastYBinResdZ, "", "");} - + if(kGCplotResolutiondRdPt == kTRUE){histograms->AddHistogram("Resolution_dR_dPt" ,"" , kGCnXBinsResdRdPt, kGCfirstXBinResdRdPt, kGClastXBinResdRdPt, kGCnYBinsResdRdPt, kGCfirstYBinResdRdPt, kGClastYBinResdRdPt, "", "");} - + if(kGCplotResolutionMCPt == kTRUE){histograms->AddHistogram("Resolution_MC_Pt" ,"" , kGCnXBinsResPt, kGCfirstXBinResPt, kGClastXBinResPt,"","");} if(kGCplotResolutionMCR == kTRUE){histograms->AddHistogram("Resolution_MC_R" ,"" , kGCnXBinsResR, kGCfirstXBinResR, kGClastXBinResR,"","");} if(kGCplotResolutionMCZ == kTRUE){histograms->AddHistogram("Resolution_MC_Z" ,"" , kGCnXBinsResZ, kGCfirstXBinResZ, kGClastXBinResZ,"","");} - + if(kGCplotResolutionESDPt == kTRUE){histograms->AddHistogram("Resolution_ESD_Pt" ,"" , kGCnXBinsResPt, kGCfirstXBinResPt, kGClastXBinResPt,"","");} if(kGCplotResolutionESDR == kTRUE){histograms->AddHistogram("Resolution_ESD_R" ,"" , kGCnXBinsResR, kGCfirstXBinResR, kGClastXBinResR,"","");} if(kGCplotResolutionESDZ == kTRUE){histograms->AddHistogram("Resolution_ESD_Z" ,"" , kGCnXBinsResZ, kGCfirstXBinResZ, kGClastXBinResZ,"","");} - + if(kGCplotESDNumberOfV0s == kTRUE){histograms->AddHistogram("ESD_NumberOfV0s","Number of v0s",100, 0, 100,"","");} if(kGCplotESDNumberOfSurvivingV0s == kTRUE){histograms->AddHistogram("ESD_NumberOfSurvivingV0s","Number of surviving v0s",100, 0, 100,"","");} - + // debug histograms if(kGCplotESDCutGetOnFly == kTRUE){histograms->AddHistogram("ESD_CutGetOnFly_InvMass" ,"Not GetOnFly" , kGCnXBinsGammaMass, kGCfirstXBinGammaMass, kGClastXBinGammaMass,"","");} if(kGCplotESDCutNContributors == kTRUE){histograms->AddHistogram("ESD_CutNContributors_InvMass" ,"NContributors <= 0" , kGCnXBinsGammaMass, kGCfirstXBinGammaMass, kGClastXBinGammaMass,"","");} @@ -1379,8 +1465,8 @@ void AddHistograms(AliGammaConversionHistograms *histograms){ if(kGCplotESDCutLine == kTRUE){histograms->AddHistogram("ESD_CutLine_InvMass" ,"Out of reconstruction area" , kGCnXBinsGammaMass, kGCfirstXBinGammaMass, kGClastXBinGammaMass,"","");} if(kGCplotESDTrueConvGammaTrackLength == kTRUE){histograms->AddHistogram("ESD_TrueConvGamma_TrackLength","Track length of TrueConvGamma",kGCnXBinsTrackLength,kGCfirstXBinTrackLength,kGClastXBinTrackLength,"","");} if(kGCplotESDTrueConvGammaTrackLengthVSInvMass == kTRUE){histograms->AddHistogram("ESD_TrueConvGamma_TrackLengthVSInvMass","Track length of TrueConvGamma vs Inv mass",kGCnXBinsTrackLength,kGCfirstXBinTrackLength,kGClastXBinTrackLength,kGCnXBinsPt, kGCfirstXBinPt, kGClastXBinPt,"","");} - - + + if(kGCplotPi0Spectra == kTRUE){ histograms->AddHistogram("ESD_Mother_InvMass_vs_Pt" ,"Invariant Mass vs Pt" , kGCnXBinsSpectra, kGCfirstXBinSpectra, kGClastXBinSpectra,kGCnYBinsSpectra, kGCfirstYBinSpectra, kGClastYBinSpectra,"InvMass [GeV]","Pt [GeV]"); histograms->AddHistogram("ESD_Mother_InvMass","Invariant mass",kGCnXBinsSpectra,kGCfirstXBinSpectra, kGClastXBinSpectra,"InvMass [GeV]","Counts"); -- 2.43.0