X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=PWGDQ%2Fdielectron%2FAliDielectronVarManager.h;h=b50a28c3f67a99878eac8b6551dad0b8c3b7b15b;hb=5720c7652fa1ef2686d5d9e4d83d91a249654ca9;hp=752e181852ba491dd77ffc60a9f26253f3b62796;hpb=36cfef5c70284218a56e779fe849ceeffe2bbe53;p=u%2Fmrichter%2FAliRoot.git diff --git a/PWGDQ/dielectron/AliDielectronVarManager.h b/PWGDQ/dielectron/AliDielectronVarManager.h index 752e181852b..b50a28c3f67 100644 --- a/PWGDQ/dielectron/AliDielectronVarManager.h +++ b/PWGDQ/dielectron/AliDielectronVarManager.h @@ -20,8 +20,12 @@ //# # //############################################################# - #include +#include +#include +#include +#include +#include #include #include @@ -29,6 +33,10 @@ #include #include #include +#include + +#include +#include #include #include @@ -78,7 +86,10 @@ public: kM, // mass kCharge, // charge kNclsITS, // number of clusters assigned in the ITS + kITSchi2Cl, // chi2/cl in the ITS kNclsTPC, // number of clusters assigned in the TPC + kNclsSTPC, // number of shared clusters assigned in the TPC + kNclsSFracTPC, // fraction of shared clusters assigned in the TPC kNclsTPCiter1, // number of clusters assigned in the TPC after first iteration kNFclsTPC, // number of findable clusters in the TPC kNFclsTPCr, // number of findable clusters in the TPC with more robust definition @@ -86,6 +97,7 @@ public: kTPCsignalN, // number of points used for dEdx kTPCsignalNfrac, // fraction of points used for dEdx / cluster used for tracking kTPCchi2Cl, // chi2/cl in TPC + kTPCclsDiff, // TPC cluster difference kTrackStatus, // track status bits kNclsTRD, // number of clusters assigned in the TRD @@ -93,24 +105,30 @@ public: kTRDpidQuality, // number of TRD tracklets used for PID kTRDprobEle, // TRD electron pid probability kTRDprobPio, // TRD electron pid probability + kTRDphi, // Phi angle of the track at the entrance of the TRD + kTRDpidEffLeg, // TRD pid efficiency from conversion electrons kImpactParXY, // Impact parameter in XY plane kImpactParZ, // Impact parameter in Z kTrackLength, // Track length - kPdgCode, // PDG code + kPdgCode, // PDG code kPdgCodeMother, // PDG code of the mother - kPdgCodeGrandMother, // PDG code of the mother - + kPdgCodeGrandMother, // PDG code of the grandmother kNumberOfDaughters, // number of daughters kHaveSameMother, // check that particles have the same mother (MC) kIsJpsiPrimary, // check if the particle is primary (MC) - kITSsignal, // ITS dE/dx signal - kITSsignalSSD1, // SSD1 dE/dx signal - kITSsignalSSD2, // SSD2 dE/dx signal - kITSsignalSDD1, // SDD1 dE/dx signal - kITSsignalSDD2, // SDD2 dE/dx signal - kITSclusterMap, // ITS cluster map + kNumberOfJPsis, // number of generated inclusive jpsis per event (MC) + kNumberOfJPsisPrompt, // number of generated prompt jpsis per event (MC) + kNumberOfJPsisNPrompt, // number of generated non-prompt jpsis per event (MC) + + kITSsignal, // ITS dE/dx signal + kITSsignalSSD1, // SSD1 dE/dx signal + kITSsignalSSD2, // SSD2 dE/dx signal + kITSsignalSDD1, // SDD1 dE/dx signal + kITSsignalSDD2, // SDD2 dE/dx signal + kITSclusterMap, // ITS cluster map + kITSLayerFirstCls, // No of innermost ITS layer with a cluster of a track kITSnSigmaEle, // number of sigmas to the dE/dx electron line in the ITS kITSnSigmaPio, // number of sigmas to the dE/dx pion line in the ITS kITSnSigmaMuo, // number of sigmas to the dE/dx muon line in the ITS @@ -118,10 +136,13 @@ public: kITSnSigmaPro, // number of sigmas to the dE/dx proton line in the ITS kPIn, // momentum at inner wall of TPC (if available), used for PID + kPOut, // momentum at outer wall of TPC, used for TRD studies + kYsignedIn, // signed local y at inner wall of TPC kTPCsignal, // TPC dE/dx signal - kTOFsignal, // TOF signal - kTOFbeta, // TOF beta + kTOFsignal, // TOF signal + kTOFbeta, // TOF beta + kTOFPIDBit, // TOF PID bit (1:set, 0:TOF not available) kTPCnSigmaEle, // number of sigmas to the dE/dx electron line in the TPC kTPCnSigmaPio, // number of sigmas to the dE/dx pion line in the TPC @@ -129,17 +150,19 @@ public: kTPCnSigmaKao, // number of sigmas to the dE/dx kaon line in the TPC kTPCnSigmaPro, // number of sigmas to the dE/dx proton line in the TPC - kTOFnSigmaEle, // number of sigmas to the pion line in the TOF + kTOFnSigmaEle, // number of sigmas to the electron line in the TOF kTOFnSigmaPio, // number of sigmas to the pion line in the TOF kTOFnSigmaMuo, // number of sigmas to the muon line in the TOF kTOFnSigmaKao, // number of sigmas to the kaon line in the TOF kTOFnSigmaPro, // number of sigmas to the proton line in the TOF + kEMCALnSigmaEle, // number of sigmas to the proton line in the TOF + kKinkIndex0, // kink index 0 kParticleMax, // // TODO: kRNClusters ?? - // AliDielectronPair specific variables + // AliDielectronPair specific variables kChi2NDF = kParticleMax, // Chi^2/NDF kDecayLength, // decay length kR, // distance to the origin @@ -147,10 +170,18 @@ public: // helicity picture: Z-axis is considered the direction of the mother's 3-momentum vector kThetaHE, // theta in mother's rest frame in the helicity picture kPhiHE, // phi in mother's rest frame in the helicity picture + kCos2PhiHE, // Cosine of 2*phi in mother's rest frame in the helicity picture // Collins-Soper picture: Z-axis is considered the direction of the vectorial difference between // the 3-mom vectors of target and projectile beams kThetaCS, // theta in mother's rest frame in Collins-Soper picture kPhiCS, // phi in mother's rest frame in Collins-Soper picture + kCos2PhiCS, // Cosine of 2*phi in mother's rest frame in the Collins-Soper picture + kDeltaPhiV0ArpH2, // Delta phi of the pair with respect to the 2nd order harmonic reaction plane from V0-A + kDeltaPhiV0CrpH2, // Delta phi of the pair with respect to the 2nd order harmonic reaction plane from V0-C + kDeltaPhiV0ACrpH2, // Delta phi of the pair with respect to the 2nd order harmonic reaction plane from V0-A + V0-C + kV0ArpH2FlowV2, // v2 coefficient with respect to the 2nd order reaction plane from V0-A + kV0CrpH2FlowV2, // v2 coefficient with respect to the 2nd order reaction plane from V0-C + kV0ACrpH2FlowV2, // v2 coefficient with respect to the 2nd order reaction plane from V0-A + V0-C kLegDist, // distance of the legs kLegDistXY, // distance of the legs in XY kDeltaEta, // Absolute value of Delta Eta for the legs @@ -159,6 +190,9 @@ public: kDCA, // distance of closest approach TODO: not implemented yet kPairType, // type of the pair, like like sign ++ unlikesign ... kPseudoProperTime, // pseudo proper time + kPseudoProperTimeResolution, // resolution for pseudo proper decay time (reconstructed - MC truth) + kPseudoProperTimePull, // normalizd resolution for pseudo proper time = (reco - MC truth)/dReco + kTRDpidEffPair, // TRD pid efficieny from conversion electrons kPairMax, // // Event specific variables kXvPrim=kPairMax, // prim vertex @@ -167,13 +201,101 @@ public: kXRes, // primary vertex x-resolution kYRes, // primary vertex y-resolution kZRes, // primary vertex z-resolution + + // v0 reaction plane quantities from AliEPSelectionTaks + kv0ArpH2, // VZERO-A reaction plane of the Q vector for 2nd harmonic + kv0CrpH2, // reaction plane + kv0ACrpH2, // VZERO-AC reaction plane of the Q vector for 2nd harmonic + kv0ATPCDiffH2, // V0A-TPC reaction plane difference for 2nd harmonic + kv0CTPCDiffH2, // V0C-TPC reaction plane difference for 2nd harmonic + kv0Av0CDiffH2, // V0A-V0C reaction plane difference for 2nd harmonic + + kMultV0A, // VZERO multiplicity and ADC amplitudes + kMultV0C, + kMultV0, + kAdcV0A, + kAdcV0C, + kAdcV0, + kVZEROchMult, + // VZERO reaction plane quantities + kV0AxH2=kVZEROchMult+64, // VZERO-A x-component of the Q vector for 2nd harmonic + kV0AyH2, // VZERO-A y-component of the Q vector for 2nd harmonic + kV0ArpH2, // VZERO-A reaction plane of the Q vector for 2nd harmonic + kV0CxH2, // VZERO-C x-component of the Q vector for 2nd harmonic + kV0CyH2, // y-component + kV0CrpH2, // reaction plane + kV0ACxH2, // VZERO-AC x-component of the Q vector for 2nd harmonic + kV0ACyH2, // VZERO-AC y-component of the Q vector for 2nd harmonic + kV0ACrpH2, // VZERO-AC reaction plane of the Q vector for 2nd harmonic + kV0ArpResH2, // 2nd harmonic reaction plane resolution for V0A + kV0CrpResH2, // V0C + kV0ACrpResH2, // V0A+V0C + kV0XaXcH2, // Correlation quantities to check V0 reaction plane quality + kV0XaYaH2, + kV0XaYcH2, + kV0YaXcH2, + kV0YaYcH2, + kV0XcYcH2, + kV0ATPCDiffH2, // V0A-TPC reaction plane difference for 2nd harmonic + kV0CTPCDiffH2, // V0C-TPC reaction plane difference for 2nd harmonic + kV0AV0CDiffH2, // V0A-V0C reaction plane difference for 2nd harmonic + // TPC reaction plane quantities + kTPCxH2, // TPC x-component of the Q vector for 2nd harmonic (corrected) + kTPCyH2, // TPC y-component of the Q vector for 2nd harmonic (corrected) + kTPCrpH2, // TPC reaction plane of the Q vector for 2nd harmonic (corrected) + kTPCsub1xH2, // TPC x-component of the Q vector for 2nd harmonic (corrected, sub event 1) + kTPCsub1yH2, // TPC y-component of the Q vector for 2nd harmonic (corrected, sub event 1) + kTPCsub1rpH2, // TPC reaction plane of the Q vector for 2nd harmonic (corrected, sub event 1) + kTPCsub2xH2, // TPC x-component of the Q vector for 2nd harmonic (corrected, sub event 2) + kTPCsub2yH2, // TPC y-component of the Q vector for 2nd harmonic (corrected, sub event 2) + kTPCsub2rpH2, // TPC reaction plane of the Q vector for 2nd harmonic (corrected, sub event 2) + kTPCsub12DiffH2, // TPC reaction plane difference of sub event 1,2 for 2nd harmonic + kTPCsub12DiffH2Sin, // TPC reaction plane difference of sub event 1,2 for 2nd harmonic, sinus term + + kTPCxH2uc, // TPC x-component of the Q vector for 2nd harmonic (uncorrected) + kTPCyH2uc, // TPC y-component of the Q vector for 2nd harmonic (uncorrected) + kTPCrpH2uc, // TPC reaction plane of the Q vector for 2nd harmonic (uncorrected) + kTPCsub1xH2uc, // TPC x-component of the Q vector for 2nd harmonic (uncorrected, sub event 1) + kTPCsub1yH2uc, // TPC y-component of the Q vector for 2nd harmonic (uncorrected, sub event 1) + kTPCsub1rpH2uc, // TPC reaction plane of the Q vector for 2nd harmonic (uncorrected, sub event 1) + kTPCsub2xH2uc, // TPC x-component of the Q vector for 2nd harmonic (uncorrected, sub event 2) + kTPCsub2yH2uc, // TPC y-component of the Q vector for 2nd harmonic (uncorrected, sub event 2) + kTPCsub2rpH2uc, // TPC reaction plane of the Q vector for 2nd harmonic (uncorrected, sub event 2) + kTPCsub12DiffH2uc, // TPC reaction plane difference of sub event 1,2 for 2nd harmonic (uncorrected) + kNTrk, // number of tracks (or tracklets) TODO: ambiguous kTracks, // ESD tracks TODO: ambiguous kNacc, // Number of accepted tracks - kNaccTrcklts, // Number of accepted tracklets (MUON definition) - kNch, // Number of charged MC tracks + kNaccTrcklts, // number of accepted SPD tracklets in |eta|<1.6 + kNaccTrcklts0916, // number of accepted SPD tracklets in 0.9<|eta|<1.6 + + kNaccTrckltsEsd05, // number of accepted SPD tracklets in |eta|<0.5 (AliESDEvent::EstimateMultiplicity()) + kNaccTrckltsEsd10, // number of accepted SPD tracklets in |eta|<1.0 (AliESDEvent::EstimateMultiplicity()) + kNaccTrckltsEsd16, // number of accepted SPD tracklets in |eta|<1.6 (AliESDEvent::EstimateMultiplicity()) + kNaccTrckltsEsd05Corr, // + kNaccTrckltsEsd10Corr, // + kNaccTrckltsEsd16Corr, // + kNaccItsTpcEsd05, // ITS-TPC tracks + ITS SA complementary tracks + tracklets from unassigned tracklets in |eta|<0.5 (AliESDEvent::EstimateMultiplicity()) + kNaccItsTpcEsd10, // ITS-TPC tracks + ITS SA complementary tracks + tracklets from unassigned tracklets in |eta|<1.0 (AliESDEvent::EstimateMultiplicity()) + kNaccItsTpcEsd16, // ITS-TPC tracks + ITS SA complementary tracks + tracklets from unassigned tracklets in |eta|<1.6 (AliESDEvent::EstimateMultiplicity()) + kNaccItsTpcEsd05Corr, // + kNaccItsTpcEsd10Corr, // + kNaccItsTpcEsd16Corr, // + + kNaccItsPureEsd05, // ITS SA tracks + tracklets from unassigned tracklets in |eta|<0.5 (AliESDEvent::EstimateMultiplicity()) + kNaccItsPureEsd10, // ITS SA tracks + tracklets from unassigned tracklets in |eta|<1.0 (AliESDEvent::EstimateMultiplicity()) + kNaccItsPureEsd16, // ITS SA tracks + tracklets from unassigned tracklets in |eta|<1.6 (AliESDEvent::EstimateMultiplicity()) + kNaccItsPureEsd05Corr, // + kNaccItsPureEsd10Corr, // + kNaccItsPureEsd16Corr, // + + kNch, // MC true number of charged particles in |eta|<1.6 + kNch05, // MC true number of charged particles in |eta|<0.5 + kNch10, // MC true number of charged particles in |eta|<1.0 + kCentrality, // event centrality fraction kNevents, // event counter + kRunNumber, // run number kNMaxValues // // TODO: (for A+A) ZDCEnergy, impact parameter, Iflag?? }; @@ -187,12 +309,18 @@ public: static void InitESDpid(Int_t type=0); static void InitAODpidUtil(Int_t type=0); - static void SetESDpid(AliESDpid * const pid) {fgPIDResponse=pid;} - static AliESDpid* GetESDpid() {return (AliESDpid*)fgPIDResponse;} - static AliAODpidUtil* GetAODpidUtil() {return (AliAODpidUtil*)fgPIDResponse;} + static void InitEstimatorAvg(const Char_t* filename); + static void InitTRDpidEffHistograms(const Char_t* filename); static void SetPIDResponse(AliPIDResponse *pidResponse) {fgPIDResponse=pidResponse;} + static AliPIDResponse* GetPIDResponse() { return fgPIDResponse; } static void SetEvent(AliVEvent * const ev); + static void SetEventData(const Double_t data[AliDielectronVarManager::kNMaxValues]); static Bool_t GetDCA(const AliAODTrack *track, Double_t d0z0[2]); + static void SetTPCEventPlane(AliEventplane *const evplane); + static void GetVzeroRP(const AliVEvent* event, Double_t* qvec, Int_t sideOption); // 0- V0A; 1- V0C; 2- V0A+V0C + + static TProfile* GetEstimatorHistogram(Int_t period, Int_t type) {return fgMultEstimatorAvg[period][type];} + static Double_t GetTRDpidEfficiency(Int_t runNo, Double_t centrality, Double_t eta, Double_t trdPhi, Double_t pout, Double_t& effErr); static const AliKFVertex* GetKFVertex() {return fgKFVertex;} @@ -217,10 +345,15 @@ private: static void FillVarESDEvent(const AliESDEvent *event, Double_t * const values); static void FillVarAODEvent(const AliAODEvent *event, Double_t * const values); static void FillVarMCEvent(const AliMCEvent *event, Double_t * const values); - + static void FillVarTPCEventPlane(const AliEventplane *evplane, Double_t * const values); + static AliPIDResponse *fgPIDResponse; // PID response object static AliVEvent *fgEvent; // current event pointer + static AliEventplane *fgTPCEventPlane; // current event plane pointer static AliKFVertex *fgKFVertex; // kf vertex + static TProfile *fgMultEstimatorAvg[4][9]; // multiplicity estimator averages (4 periods x 9 estimators) + static Double_t fgTRDpidEffCentRanges[10][4]; // centrality ranges for the TRD pid efficiency histograms + static TH3D *fgTRDpidEff[10][4]; // TRD pid efficiencies from conversion electrons static Double_t fgData[kNMaxValues]; //! data @@ -237,6 +370,7 @@ inline void AliDielectronVarManager::Fill(const TObject* object, Double_t * cons // // Main function to fill all available variables according to the type of particle // +// if (!object) return; if (object->IsA() == AliESDtrack::Class()) FillVarESDtrack(static_cast(object), values); else if (object->IsA() == AliAODTrack::Class()) FillVarAODTrack(static_cast(object), values); else if (object->IsA() == AliMCParticle::Class()) FillVarMCParticle(static_cast(object), values); @@ -249,6 +383,7 @@ inline void AliDielectronVarManager::Fill(const TObject* object, Double_t * cons else if (object->IsA() == AliESDEvent::Class()) FillVarESDEvent(static_cast(object), values); else if (object->IsA() == AliAODEvent::Class()) FillVarAODEvent(static_cast(object), values); else if (object->IsA() == AliMCEvent::Class()) FillVarMCEvent(static_cast(object), values); + else if (object->IsA() == AliEventplane::Class()) FillVarTPCEventPlane(static_cast(object), values); // else printf(Form("AliDielectronVarManager::Fill: Type %s is not supported by AliDielectronVarManager!", object->ClassName())); //TODO: implement without object needed } @@ -293,12 +428,25 @@ inline void AliDielectronVarManager::FillVarESDtrack(const AliESDtrack *particle // Fill common AliVParticle interface information FillVarVParticle(particle, values); + AliESDtrack *esdTrack=0x0; + Double_t origdEdx=particle->GetTPCsignal(); + + // apply ETa correction, remove once this is in the tender + if( AliDielectronPID::GetEtaCorrFunction() ){ + esdTrack=const_cast(particle); + esdTrack->SetTPCsignal(origdEdx/AliDielectronPID::GetEtaCorr(esdTrack),esdTrack->GetTPCsignalSigma(),esdTrack->GetTPCsignalN()); + } + Double_t pidProbs[AliPID::kSPECIES]; // Fill AliESDtrack interface specific information Double_t tpcNcls=particle->GetTPCNcls(); + Double_t tpcNclsS = particle->GetTPCnclsS(); + Double_t itsNcls=particle->GetNcls(0); Double_t tpcSignalN=particle->GetTPCsignalN(); - values[AliDielectronVarManager::kNclsITS] = particle->GetNcls(0); // TODO: get rid of the plain numbers + values[AliDielectronVarManager::kNclsITS] = itsNcls; // TODO: get rid of the plain numbers values[AliDielectronVarManager::kNclsTPC] = tpcNcls; // TODO: get rid of the plain numbers + values[AliDielectronVarManager::kNclsSTPC] = tpcNclsS; + values[AliDielectronVarManager::kNclsSFracTPC] = tpcNcls>0?tpcNclsS/tpcNcls:0; values[AliDielectronVarManager::kNclsTPCiter1] = particle->GetTPCNclsIter1(); // TODO: get rid of the plain numbers values[AliDielectronVarManager::kNFclsTPC] = particle->GetTPCNclsF(); values[AliDielectronVarManager::kNFclsTPCr] = particle->GetTPCClusterInfo(2,1); @@ -308,11 +456,14 @@ inline void AliDielectronVarManager::FillVarESDtrack(const AliESDtrack *particle values[AliDielectronVarManager::kNclsTRD] = particle->GetNcls(2); // TODO: get rid of the plain numbers values[AliDielectronVarManager::kTRDntracklets] = particle->GetTRDntracklets(); // TODO: GetTRDtracklets/GetTRDntracklets? values[AliDielectronVarManager::kTRDpidQuality] = particle->GetTRDpidQuality(); + values[AliDielectronVarManager::kTPCclsDiff] = tpcSignalN-tpcNcls; values[AliDielectronVarManager::kTrackStatus] = (Double_t)particle->GetStatus(); values[AliDielectronVarManager::kTPCchi2Cl] = -1; if (tpcNcls>0) values[AliDielectronVarManager::kTPCchi2Cl] = particle->GetTPCchi2() / tpcNcls; + values[AliDielectronVarManager::kITSchi2Cl] = -1; + if (itsNcls>0) values[AliDielectronVarManager::kITSchi2Cl] = particle->GetITSchi2() / itsNcls; //TRD pidProbs particle->GetTRDpid(pidProbs); values[AliDielectronVarManager::kTRDprobEle] = pidProbs[AliPID::kElectron]; @@ -361,13 +512,45 @@ inline void AliDielectronVarManager::FillVarESDtrack(const AliESDtrack *particle values[AliDielectronVarManager::kITSsignalSDD1] = itsdEdx[2]; values[AliDielectronVarManager::kITSsignalSDD2] = itsdEdx[3]; values[AliDielectronVarManager::kITSclusterMap] = particle->GetITSClusterMap(); + values[AliDielectronVarManager::kITSLayerFirstCls] = -1.; + + for (Int_t iC=0; iC<6; iC++) { + if (((particle->GetITSClusterMap()) & (1<<(iC))) > 0) { + values[AliDielectronVarManager::kITSLayerFirstCls] = iC; + break; + } + } + values[AliDielectronVarManager::kTrackLength] = particle->GetIntegratedLength(); //dEdx information Double_t mom = particle->GetP(); const AliExternalTrackParam *in=particle->GetInnerParam(); - if (in) mom = in->GetP(); + Double_t ysignedIn=-100; + if (in) { + mom = in->GetP(); + ysignedIn=particle->Charge()*in->GetY(); + } values[AliDielectronVarManager::kPIn]=mom; + values[AliDielectronVarManager::kYsignedIn]=ysignedIn; + const AliExternalTrackParam *out=particle->GetOuterParam(); + if(out) values[AliDielectronVarManager::kPOut] = out->GetP(); + else values[AliDielectronVarManager::kPOut] = mom; + if(out && fgEvent) { + Double_t localCoord[3]={0.0}; + Bool_t localCoordGood = out->GetXYZAt(298.0, ((AliESDEvent*)fgEvent)->GetMagneticField(), localCoord); + values[AliDielectronVarManager::kTRDphi] = (localCoordGood && TMath::Abs(localCoord[0])>1.0e-6 && TMath::Abs(localCoord[1])>1.0e-6 ? TMath::ATan2(localCoord[1], localCoord[0]) : -999.); + } + if(mc->HasMC() && fgTRDpidEff[0][0]) { + Int_t runNo = (fgEvent ? fgEvent->GetRunNumber() : -1); + Float_t centrality=-1.0; + AliCentrality *esdCentrality = (fgEvent ? fgEvent->GetCentrality() : 0x0); + if(esdCentrality) centrality = esdCentrality->GetCentralityPercentile("V0M"); + Double_t effErr=0.0; + values[kTRDpidEffLeg] = GetTRDpidEfficiency(runNo, centrality, values[AliDielectronVarManager::kEta], + values[AliDielectronVarManager::kTRDphi], + values[AliDielectronVarManager::kPOut], effErr); + } values[AliDielectronVarManager::kTPCsignal]=particle->GetTPCsignal(); values[AliDielectronVarManager::kTOFsignal]=particle->GetTOFsignal(); @@ -388,6 +571,7 @@ inline void AliDielectronVarManager::FillVarESDtrack(const AliESDtrack *particle Float_t beta = v / TMath::C(); values[AliDielectronVarManager::kTOFbeta]=beta; } + values[AliDielectronVarManager::kTOFPIDBit]=(particle->GetStatus()&AliESDtrack::kTOFpid? 1: 0); // nsigma to Electron band // TODO: for the moment we set the bethe bloch parameters manually @@ -410,6 +594,11 @@ inline void AliDielectronVarManager::FillVarESDtrack(const AliESDtrack *particle values[AliDielectronVarManager::kTOFnSigmaMuo]=fgPIDResponse->NumberOfSigmasTOF(particle,AliPID::kMuon); values[AliDielectronVarManager::kTOFnSigmaKao]=fgPIDResponse->NumberOfSigmasTOF(particle,AliPID::kKaon); values[AliDielectronVarManager::kTOFnSigmaPro]=fgPIDResponse->NumberOfSigmasTOF(particle,AliPID::kProton); + + values[AliDielectronVarManager::kEMCALnSigmaEle]=fgPIDResponse->NumberOfSigmasEMCAL(particle,AliPID::kElectron); + + //restore TPC signal if it was changed + if (esdTrack) esdTrack->SetTPCsignal(origdEdx,esdTrack->GetTPCsignalSigma(),esdTrack->GetTPCsignalN()); } inline void AliDielectronVarManager::FillVarAODTrack(const AliAODTrack *particle, Double_t * const values) @@ -420,11 +609,18 @@ inline void AliDielectronVarManager::FillVarAODTrack(const AliAODTrack *particle // Fill common AliVParticle interface information FillVarVParticle(particle, values); - Double_t tpcNcls=particle->GetTPCNcls(); + + //GetNclsS not present in AODtrack + //Replace with method as soon as available + TBits tpcSharedMap = particle->GetTPCSharedMap(); + Double_t tpcNclsS= tpcSharedMap.CountBits(0)-tpcSharedMap.CountBits(159); // Reset AliESDtrack interface specific information values[AliDielectronVarManager::kNclsITS] = 0; + values[AliDielectronVarManager::kITSchi2Cl] = -1; values[AliDielectronVarManager::kNclsTPC] = tpcNcls; + values[AliDielectronVarManager::kNclsSTPC] = tpcNclsS; + values[AliDielectronVarManager::kNclsSFracTPC] = tpcNcls>0?tpcNclsS/tpcNcls:0; values[AliDielectronVarManager::kNclsTPCiter1] = tpcNcls; // not really available in AOD values[AliDielectronVarManager::kNFclsTPC] = 0; values[AliDielectronVarManager::kNFclsTPCr] = 0; @@ -471,6 +667,7 @@ inline void AliDielectronVarManager::FillVarAODTrack(const AliAODTrack *particle Double_t tpcSignalN=pid->GetTPCsignalN(); values[AliDielectronVarManager::kTPCsignalN] = tpcSignalN; values[AliDielectronVarManager::kTPCsignalN] = tpcNcls>0?tpcSignalN/tpcNcls:0; + values[AliDielectronVarManager::kTPCclsDiff] = tpcSignalN-tpcNcls; Double_t tpcNsigmaEle=fgPIDResponse->NumberOfSigmasTPC(particle,AliPID::kElectron); Double_t tpcNsigmaPio=fgPIDResponse->NumberOfSigmasTPC(particle,AliPID::kPion); Double_t tpcNsigmaMuo=fgPIDResponse->NumberOfSigmasTPC(particle,AliPID::kMuon); @@ -488,6 +685,18 @@ inline void AliDielectronVarManager::FillVarAODTrack(const AliAODTrack *particle values[AliDielectronVarManager::kTRDntracklets] = 0; values[AliDielectronVarManager::kTRDpidQuality] = 0; + + Double_t tofNsigmaEle=fgPIDResponse->NumberOfSigmasTOF(particle,AliPID::kElectron); + Double_t tofNsigmaPio=fgPIDResponse->NumberOfSigmasTOF(particle,AliPID::kPion); + Double_t tofNsigmaMuo=fgPIDResponse->NumberOfSigmasTOF(particle,AliPID::kMuon); + Double_t tofNsigmaKao=fgPIDResponse->NumberOfSigmasTOF(particle,AliPID::kKaon); + Double_t tofNsigmaPro=fgPIDResponse->NumberOfSigmasTOF(particle,AliPID::kProton); + + values[AliDielectronVarManager::kTOFnSigmaEle]=tofNsigmaEle; + values[AliDielectronVarManager::kTOFnSigmaPio]=tofNsigmaPio; + values[AliDielectronVarManager::kTOFnSigmaMuo]=tofNsigmaMuo; + values[AliDielectronVarManager::kTOFnSigmaKao]=tofNsigmaKao; + values[AliDielectronVarManager::kTOFnSigmaPro]=tofNsigmaPro; } @@ -514,6 +723,7 @@ inline void AliDielectronVarManager::FillVarAODTrack(const AliAODTrack *particle values[AliDielectronVarManager::kNumberOfDaughters]=mc->NumberOfDaughters(particle); } //if(mc->HasMC()) + values[AliDielectronVarManager::kTOFPIDBit]=(particle->GetStatus()&AliESDtrack::kTOFpid? 1: 0); } inline void AliDielectronVarManager::FillVarMCParticle(const AliMCParticle *particle, Double_t * const values) @@ -523,7 +733,10 @@ inline void AliDielectronVarManager::FillVarMCParticle(const AliMCParticle *part // values[AliDielectronVarManager::kNclsITS] = 0; + values[AliDielectronVarManager::kITSchi2Cl] = 0; values[AliDielectronVarManager::kNclsTPC] = 0; + values[AliDielectronVarManager::kNclsSTPC] = 0; + values[AliDielectronVarManager::kNclsSFracTPC] = 0; values[AliDielectronVarManager::kNclsTPCiter1] = 0; values[AliDielectronVarManager::kNFclsTPC] = 0; values[AliDielectronVarManager::kNFclsTPCr] = 0; @@ -536,10 +749,12 @@ inline void AliDielectronVarManager::FillVarMCParticle(const AliMCParticle *part values[AliDielectronVarManager::kTRDprobEle] = 0; values[AliDielectronVarManager::kTRDprobPio] = 0; values[AliDielectronVarManager::kTPCsignalN] = 0; + values[AliDielectronVarManager::kTPCclsDiff] = 0; values[AliDielectronVarManager::kTPCsignalNfrac] = 0; values[AliDielectronVarManager::kImpactParXY] = 0; values[AliDielectronVarManager::kImpactParZ] = 0; values[AliDielectronVarManager::kPIn] = 0; + values[AliDielectronVarManager::kYsignedIn] = 0; values[AliDielectronVarManager::kTPCsignal] = 0; values[AliDielectronVarManager::kTOFsignal] = 0; values[AliDielectronVarManager::kTOFbeta] = 0; @@ -580,7 +795,10 @@ inline void AliDielectronVarManager::FillVarMCParticle2(const AliVParticle *p1, // values[AliDielectronVarManager::kNclsITS] = 0; + values[AliDielectronVarManager::kITSchi2Cl] = -1; values[AliDielectronVarManager::kNclsTPC] = 0; + values[AliDielectronVarManager::kNclsSTPC] = 0; + values[AliDielectronVarManager::kNclsSFracTPC] = 0; values[AliDielectronVarManager::kNclsTPCiter1] = 0; values[AliDielectronVarManager::kNFclsTPC] = 0; values[AliDielectronVarManager::kNFclsTPCr] = 0; @@ -593,10 +811,12 @@ inline void AliDielectronVarManager::FillVarMCParticle2(const AliVParticle *p1, values[AliDielectronVarManager::kTRDprobEle] = 0; values[AliDielectronVarManager::kTRDprobPio] = 0; values[AliDielectronVarManager::kTPCsignalN] = 0; + values[AliDielectronVarManager::kTPCclsDiff] = 0; values[AliDielectronVarManager::kTPCsignalNfrac] = 0; values[AliDielectronVarManager::kImpactParXY] = 0; values[AliDielectronVarManager::kImpactParZ] = 0; values[AliDielectronVarManager::kPIn] = 0; + values[AliDielectronVarManager::kYsignedIn] = 0; values[AliDielectronVarManager::kTPCsignal] = 0; values[AliDielectronVarManager::kTPCnSigmaEle] = 0; values[AliDielectronVarManager::kTPCnSigmaPio] = 0; @@ -614,8 +834,14 @@ inline void AliDielectronVarManager::FillVarMCParticle2(const AliVParticle *p1, Int_t mLabel2 = mc->GetMothersLabel(p2->GetLabel()); if(mLabel1==mLabel2) mother = mc->GetMCTrackFromMCEvent(mLabel1); - if(mother) // same mother + + values[AliDielectronVarManager::kPseudoProperTime] = -2e10; + if(mother) { // same mother FillVarVParticle(mother, values); + Double_t lxy = ((mother->Xv()- values[AliDielectronVarManager::kXvPrim]) * mother->Px() + + (mother->Yv()- values[AliDielectronVarManager::kYvPrim]) * mother->Py() )/mother->Pt(); + values[AliDielectronVarManager::kPseudoProperTime] = lxy*(TDatabasePDG::Instance()->GetParticle(443)->Mass())/mother->Pt(); + } else { // the 2 particles come from different mothers so 2-particles quantities are calculated here // AliVParticle part values[AliDielectronVarManager::kPx] = p1->Px()+p2->Px(); @@ -642,11 +868,17 @@ inline void AliDielectronVarManager::FillVarMCParticle2(const AliVParticle *p1, values[AliDielectronVarManager::kY] = ((values[AliDielectronVarManager::kE]-values[AliDielectronVarManager::kPz])>1.0e-6 && (values[AliDielectronVarManager::kE]+values[AliDielectronVarManager::kPz])>1.0e-6 ? 0.5*TMath::Log((values[AliDielectronVarManager::kE]+values[AliDielectronVarManager::kPz])/(values[AliDielectronVarManager::kE]-values[AliDielectronVarManager::kPz])) : -9999.); values[AliDielectronVarManager::kCharge] = p1->Charge()+p2->Charge(); - values[AliDielectronVarManager::kM] = p1->M()*p1->M()+p2->M()*p2->M()+ - 2.0*(p1->E()*p2->E()-p1->Px()*p2->Px()-p1->Py()*p2->Py()-p1->Pz()*p2->Pz()); + values[AliDielectronVarManager::kM] = TMath::Sqrt(p1->M()*p1->M()+p2->M()*p2->M()+ + 2.0*(p1->E()*p2->E()-p1->Px()*p2->Px()-p1->Py()*p2->Py()-p1->Pz()*p2->Pz())); if ( fgEvent ) AliDielectronVarManager::Fill(fgEvent, values); } + values[AliDielectronVarManager::kThetaHE] = AliDielectronPair::ThetaPhiCM(p1,p2,kTRUE, kTRUE); + values[AliDielectronVarManager::kPhiHE] = AliDielectronPair::ThetaPhiCM(p1,p2,kTRUE, kFALSE); + values[AliDielectronVarManager::kCos2PhiHE] = TMath::Cos(2*values[AliDielectronVarManager::kPhiHE]); + values[AliDielectronVarManager::kThetaCS] = AliDielectronPair::ThetaPhiCM(p1,p2,kFALSE, kTRUE); + values[AliDielectronVarManager::kPhiCS] = AliDielectronPair::ThetaPhiCM(p1,p2,kFALSE, kFALSE); + values[AliDielectronVarManager::kCos2PhiCS] = TMath::Cos(2*values[AliDielectronVarManager::kPhiCS]); } @@ -657,7 +889,10 @@ inline void AliDielectronVarManager::FillVarAODMCParticle(const AliAODMCParticle // values[AliDielectronVarManager::kNclsITS] = 0; + values[AliDielectronVarManager::kITSchi2Cl] = -1; values[AliDielectronVarManager::kNclsTPC] = 0; + values[AliDielectronVarManager::kNclsSTPC] = 0; + values[AliDielectronVarManager::kNclsSFracTPC] = 0; values[AliDielectronVarManager::kNclsTPCiter1] = 0; values[AliDielectronVarManager::kNFclsTPC] = 0; values[AliDielectronVarManager::kNclsTRD] = 0; @@ -668,10 +903,12 @@ inline void AliDielectronVarManager::FillVarAODMCParticle(const AliAODMCParticle values[AliDielectronVarManager::kTRDprobEle] = 0; values[AliDielectronVarManager::kTRDprobPio] = 0; values[AliDielectronVarManager::kTPCsignalN] = 0; + values[AliDielectronVarManager::kTPCclsDiff] = 0; values[AliDielectronVarManager::kTPCsignalNfrac]= 0; values[AliDielectronVarManager::kImpactParXY] = 0; values[AliDielectronVarManager::kImpactParZ] = 0; values[AliDielectronVarManager::kPIn] = 0; + values[AliDielectronVarManager::kYsignedIn] = 0; values[AliDielectronVarManager::kTPCsignal] = 0; values[AliDielectronVarManager::kTPCnSigmaEle] = 0; values[AliDielectronVarManager::kTPCnSigmaPio] = 0; @@ -728,31 +965,83 @@ inline void AliDielectronVarManager::FillVarDielectronPair(const AliDielectronPa Double_t phiCS=0; pair->GetThetaPhiCM(thetaHE,phiHE,thetaCS,phiCS); - + values[AliDielectronVarManager::kChi2NDF] = kfPair.GetChi2()/kfPair.GetNDF(); values[AliDielectronVarManager::kDecayLength] = kfPair.GetDecayLength(); values[AliDielectronVarManager::kR] = kfPair.GetR(); values[AliDielectronVarManager::kOpeningAngle] = pair->OpeningAngle(); values[AliDielectronVarManager::kThetaHE] = thetaHE; values[AliDielectronVarManager::kPhiHE] = phiHE; + values[AliDielectronVarManager::kCos2PhiHE] = TMath::Cos(2.0*phiHE); values[AliDielectronVarManager::kThetaCS] = thetaCS; values[AliDielectronVarManager::kPhiCS] = phiCS; + values[AliDielectronVarManager::kCos2PhiCS] = TMath::Cos(2.0*phiCS); values[AliDielectronVarManager::kLegDist] = pair->DistanceDaughters(); values[AliDielectronVarManager::kLegDistXY] = pair->DistanceDaughtersXY(); values[AliDielectronVarManager::kDeltaEta] = pair->DeltaEta(); values[AliDielectronVarManager::kDeltaPhi] = pair->DeltaPhi(); values[AliDielectronVarManager::kMerr] = kfPair.GetErrMass()>1e-30&&kfPair.GetMass()>1e-30?kfPair.GetErrMass()/kfPair.GetMass():1000000; values[AliDielectronVarManager::kPairType] = pair->GetType(); - values[AliDielectronVarManager::kPseudoProperTime] = pair->GetPseudoProperTime(fgEvent->GetPrimaryVertex()); - + Double_t errPseudoProperTime2 = -1; + values[AliDielectronVarManager::kPseudoProperTime] = fgEvent ? pair->GetKFParticle().GetPseudoProperDecayTime(*(fgEvent->GetPrimaryVertex()), TDatabasePDG::Instance()->GetParticle(443)->Mass(), &errPseudoProperTime2 ) : -1e10; + // values[AliDielectronVarManager::kPseudoProperTime] = fgEvent ? pair->GetPseudoProperTime(fgEvent->GetPrimaryVertex()): -1e10; + + // Flow quantities + Double_t delta=0.0; + // v2 with respect to VZERO-A event plane + delta = values[AliDielectronVarManager::kPhi] - fgData[AliDielectronVarManager::kV0ArpH2]; + if(delta>TMath::Pi()) delta -= 2.0*TMath::Pi(); // keep the [-pi,+pi] interval + if(delta<-1.0*TMath::Pi()) delta += 2.0*TMath::Pi(); + values[AliDielectronVarManager::kV0ArpH2FlowV2] = TMath::Cos(2.0*delta); // 2nd harmonic flow coefficient + values[AliDielectronVarManager::kDeltaPhiV0ArpH2] = delta; + // v2 with respect to VZERO-C event plane + delta = values[AliDielectronVarManager::kPhi] - fgData[AliDielectronVarManager::kV0CrpH2]; + if(delta>TMath::Pi()) delta -= 2.0*TMath::Pi(); // keep the [-pi,+pi] interval + if(delta<-1.0*TMath::Pi()) delta += 2.0*TMath::Pi(); + values[AliDielectronVarManager::kV0CrpH2FlowV2] = TMath::Cos(2.0*delta); // 2nd harmonic flow coefficient + values[AliDielectronVarManager::kDeltaPhiV0CrpH2] = delta; + // v2 with respect to the combined VZERO-A and VZERO-C event plane + delta = values[AliDielectronVarManager::kPhi] - fgData[AliDielectronVarManager::kV0ACrpH2]; + if(delta>TMath::Pi()) delta -= 2.0*TMath::Pi(); // keep the [-pi,+pi] interval + if(delta<-1.0*TMath::Pi()) delta += 2.0*TMath::Pi(); + values[AliDielectronVarManager::kV0ACrpH2FlowV2] = TMath::Cos(2.0*delta); // 2nd harmonic flow coefficient + values[AliDielectronVarManager::kDeltaPhiV0ACrpH2] = delta; + + AliDielectronMC *mc=AliDielectronMC::Instance(); if (mc->HasMC()){ + values[AliDielectronVarManager::kPseudoProperTimeResolution] = -10.0e+10; Bool_t samemother = mc->HaveSameMother(pair); values[AliDielectronVarManager::kIsJpsiPrimary] = mc->IsJpsiPrimary(pair); values[AliDielectronVarManager::kHaveSameMother] = samemother ; - + + // fill kPseudoProperTimeResolution + values[AliDielectronVarManager::kPseudoProperTimeResolution] = -1e10; + // values[AliDielectronVarManager::kPseudoProperTimePull] = -1e10; + if(samemother && fgEvent) { + if(pair->GetFirstDaughter()->GetLabel() > 0) { + const AliVParticle* d1 = mc->GetMCTrackFromMCEvent(pair->GetFirstDaughter()->GetLabel()); + const AliVParticle* motherMC = mc->GetMCTrackFromMCEvent(((AliMCParticle*)d1)->GetMother()); + const AliVVertex* mcVtx = mc->GetMCEvent()->GetPrimaryVertex(); + if(motherMC && mcVtx) { + const Double_t lxyMC = ( (motherMC->Xv() - mcVtx->GetX()) * motherMC->Px() + + (motherMC->Yv() - mcVtx->GetY()) * motherMC->Py() ) / motherMC->Pt(); + const Double_t pseudoMC = lxyMC * (TDatabasePDG::Instance()->GetParticle(443)->Mass())/motherMC->Pt(); + values[AliDielectronVarManager::kPseudoProperTimeResolution] = values[AliDielectronVarManager::kPseudoProperTime] - pseudoMC; + if (errPseudoProperTime2 > 0) + values[AliDielectronVarManager::kPseudoProperTimePull] = values[AliDielectronVarManager::kPseudoProperTimeResolution]/sqrt(errPseudoProperTime2); + } + } + } + Double_t valuesLeg1[AliDielectronVarManager::kNMaxValues]; + Double_t valuesLeg2[AliDielectronVarManager::kNMaxValues]; + AliVParticle* leg1 = pair->GetFirstDaughter(); + AliVParticle* leg2 = pair->GetSecondDaughter(); + Fill(leg1, valuesLeg1); + Fill(leg2, valuesLeg2); + values[AliDielectronVarManager::kTRDpidEffPair] = valuesLeg1[AliDielectronVarManager::kTRDpidEffLeg]*valuesLeg2[AliDielectronVarManager::kTRDpidEffLeg]; }//if (mc->HasMC()) @@ -784,7 +1073,10 @@ inline void AliDielectronVarManager::FillVarKFParticle(const AliKFParticle *part values[AliDielectronVarManager::kCharge] = particle->GetQ(); values[AliDielectronVarManager::kNclsITS] = 0; + values[AliDielectronVarManager::kITSchi2Cl] = -1; values[AliDielectronVarManager::kNclsTPC] = 0; + values[AliDielectronVarManager::kNclsSTPC] = 0; + values[AliDielectronVarManager::kNclsSFracTPC] = 0; values[AliDielectronVarManager::kNclsTPCiter1] = 0; values[AliDielectronVarManager::kNFclsTPC] = 0; values[AliDielectronVarManager::kNclsTRD] = 0; @@ -795,10 +1087,12 @@ inline void AliDielectronVarManager::FillVarKFParticle(const AliKFParticle *part values[AliDielectronVarManager::kTRDprobEle] = 0; values[AliDielectronVarManager::kTRDprobPio] = 0; values[AliDielectronVarManager::kTPCsignalN] = 0; + values[AliDielectronVarManager::kTPCclsDiff] = 0; values[AliDielectronVarManager::kTPCsignalNfrac]= 0; values[AliDielectronVarManager::kImpactParXY] = 0; values[AliDielectronVarManager::kImpactParZ] = 0; values[AliDielectronVarManager::kPIn] = 0; + values[AliDielectronVarManager::kYsignedIn] = 0; values[AliDielectronVarManager::kTPCsignal] = 0; values[AliDielectronVarManager::kTOFsignal] = 0; values[AliDielectronVarManager::kTOFbeta] = 0; @@ -814,7 +1108,10 @@ inline void AliDielectronVarManager::FillVarKFParticle(const AliKFParticle *part values[AliDielectronVarManager::kPdgCodeGrandMother] = -1; - if ( fgEvent ) AliDielectronVarManager::Fill(fgEvent, values); +// if ( fgEvent ) AliDielectronVarManager::Fill(fgEvent, values); + for (Int_t i=AliDielectronVarManager::kPairMax; iGetRunNumber(); const AliVVertex *primVtx = event->GetPrimaryVertex(); values[AliDielectronVarManager::kXvPrim] = 0; values[AliDielectronVarManager::kYvPrim] = 0; values[AliDielectronVarManager::kZvPrim] = 0; // values[AliDielectronVarManager::kChi2NDF] = 0; //This is the pair value!!! - + values[AliDielectronVarManager::kNTrk] = 0; values[AliDielectronVarManager::kNacc] = 0; values[AliDielectronVarManager::kNaccTrcklts] = 0; + values[AliDielectronVarManager::kNaccTrcklts0916] = 0; values[AliDielectronVarManager::kNevents] = 0; //always fill bin 0; - + if (!primVtx) return; values[AliDielectronVarManager::kXvPrim] = primVtx->GetX(); @@ -841,9 +1140,107 @@ inline void AliDielectronVarManager::FillVarVEvent(const AliVEvent *event, Doubl values[AliDielectronVarManager::kZvPrim] = primVtx->GetZ(); // values[AliDielectronVarManager::kChi2NDF] = primVtx->GetChi2perNDF(); //this is the pair value - values[AliDielectronVarManager::kNTrk] = event->GetNumberOfTracks(); - values[AliDielectronVarManager::kNacc] = AliDielectronHelper::GetNacc(event); - values[AliDielectronVarManager::kNaccTrcklts] = AliDielectronHelper::GetNaccTrcklts(event); + values[AliDielectronVarManager::kNTrk] = event->GetNumberOfTracks(); + values[AliDielectronVarManager::kNacc] = AliDielectronHelper::GetNacc(event); + values[AliDielectronVarManager::kNaccTrcklts] = AliDielectronHelper::GetNaccTrcklts(event); // etaRange = 1.6 (default) + values[AliDielectronVarManager::kNaccTrcklts0916] = AliDielectronHelper::GetNaccTrcklts(event,1.6)-AliDielectronHelper::GetNaccTrcklts(event,.9); + // values[AliDielectronVarManager::kNaccTrcklts05] = AliDielectronHelper::GetNaccTrcklts(event, 0.5); + // values[AliDielectronVarManager::kNaccTrcklts10] = AliDielectronHelper::GetNaccTrcklts(event, 1.0); + // values[AliDielectronVarManager::kNaccTrckltsCorr] = AliDielectronHelper::GetNaccTrckltsCorrected(event, values[AliDielectronVarManager::kNaccTrcklts], values[AliDielectronVarManager::kZvPrim]); + + // TPC event plane (corrected) + if(fgTPCEventPlane) FillVarTPCEventPlane(fgTPCEventPlane, values); + + // VZERO event plane quantities from the AliEPSelectionTask + AliEventplane *ep = const_cast(event)->GetEventplane(); + values[AliDielectronVarManager::kv0ACrpH2] = TVector2::Phi_mpi_pi(ep->GetEventplane("V0", event, 2)); + values[AliDielectronVarManager::kv0ArpH2] = TVector2::Phi_mpi_pi(ep->GetEventplane("V0A",event, 2)); + values[AliDielectronVarManager::kv0CrpH2] = TVector2::Phi_mpi_pi(ep->GetEventplane("V0C",event, 2)); + + // ESD VZERO information + AliVVZERO* vzeroData = event->GetVZEROData(); + values[AliDielectronVarManager::kMultV0A] = 0.0; + values[AliDielectronVarManager::kMultV0C] = 0.0; + values[AliDielectronVarManager::kAdcV0A] = 0.0; + values[AliDielectronVarManager::kAdcV0C] = 0.0; + for(Int_t i=0; i<32; ++i) { + //values[AliDielectronVarManager::kVZEROchMult+i] = vzeroData->GetMultiplicity(i); + //values[AliDielectronVarManager::kVZEROchMult+32+i] = vzeroData->GetMultiplicity(i+32); + values[AliDielectronVarManager::kVZEROchMult+i] = event->GetVZEROEqMultiplicity(i); + values[AliDielectronVarManager::kVZEROchMult+32+i] = event->GetVZEROEqMultiplicity(i+32); + values[AliDielectronVarManager::kMultV0A] += vzeroData->GetMultiplicityV0A(i); + values[AliDielectronVarManager::kMultV0C] += vzeroData->GetMultiplicityV0C(i); + //values[AliDielectronVarManager::kAdcV0A] += vzeroData->GetAdcV0A(i); + //values[AliDielectronVarManager::kAdcV0C] += vzeroData->GetAdcV0C(i); + } + values[AliDielectronVarManager::kMultV0] = values[AliDielectronVarManager::kMultV0A] + values[AliDielectronVarManager::kMultV0C]; + values[AliDielectronVarManager::kAdcV0] = values[AliDielectronVarManager::kAdcV0A] + values[AliDielectronVarManager::kAdcV0C]; + // VZERO event plane quantities + Double_t qvec[3]={0.0}; + GetVzeroRP(event, qvec,0); // V0-A + values[AliDielectronVarManager::kV0AxH2] = qvec[0]; values[AliDielectronVarManager::kV0AyH2] = qvec[1]; + values[AliDielectronVarManager::kV0ArpH2] = qvec[2]; + qvec[0]=0.0; qvec[1]=0.0; qvec[2]=0.0; + GetVzeroRP(event, qvec,1); // V0-C + values[AliDielectronVarManager::kV0CxH2] = qvec[0]; values[AliDielectronVarManager::kV0CyH2] = qvec[1]; + values[AliDielectronVarManager::kV0CrpH2] = qvec[2]; + qvec[0]=0.0; qvec[1]=0.0; qvec[2]=0.0; + GetVzeroRP(event, qvec,2); // V0-A and V0-C combined + values[AliDielectronVarManager::kV0ACxH2] = qvec[0]; values[AliDielectronVarManager::kV0ACyH2] = qvec[1]; + values[AliDielectronVarManager::kV0ACrpH2] = qvec[2]; + // VZERO event plane resolution + values[AliDielectronVarManager::kV0ArpResH2] = 1.0; + values[AliDielectronVarManager::kV0CrpResH2] = 1.0; + values[AliDielectronVarManager::kV0ACrpResH2] = 1.0; + // Q vector components correlations + values[AliDielectronVarManager::kV0XaXcH2] = values[AliDielectronVarManager::kV0AxH2]*values[AliDielectronVarManager::kV0CxH2]; + values[AliDielectronVarManager::kV0XaYaH2] = values[AliDielectronVarManager::kV0AxH2]*values[AliDielectronVarManager::kV0AyH2]; + values[AliDielectronVarManager::kV0XaYcH2] = values[AliDielectronVarManager::kV0AxH2]*values[AliDielectronVarManager::kV0CyH2]; + values[AliDielectronVarManager::kV0YaXcH2] = values[AliDielectronVarManager::kV0AyH2]*values[AliDielectronVarManager::kV0CxH2]; + values[AliDielectronVarManager::kV0YaYcH2] = values[AliDielectronVarManager::kV0AyH2]*values[AliDielectronVarManager::kV0CyH2]; + values[AliDielectronVarManager::kV0XcYcH2] = values[AliDielectronVarManager::kV0CxH2]*values[AliDielectronVarManager::kV0CyH2]; + + // TPC event plane quantities (uncorrected) + AliEventplane *evplane = const_cast(event)->GetEventplane(); + TVector2 *qstd = evplane->GetQVector(); // This is the "standard" Q-Vector + TVector2 *qsub1 = evplane->GetQsub1(); + TVector2 *qsub2 = evplane->GetQsub2(); + if(!qstd || !qsub1 || !qsub2) return; + + values[AliDielectronVarManager::kTPCxH2uc] = qstd->X(); + values[AliDielectronVarManager::kTPCyH2uc] = qstd->Y(); + values[AliDielectronVarManager::kTPCrpH2uc] = ((TMath::Abs(qstd->X())>1.0e-10) ? TMath::ATan2(qstd->Y(),qstd->X())/2.0 : 0.0); + + values[AliDielectronVarManager::kTPCsub1xH2uc] = qsub1->X(); + values[AliDielectronVarManager::kTPCsub1yH2uc] = qsub1->Y(); + values[AliDielectronVarManager::kTPCsub1rpH2uc] = ((TMath::Abs(qsub1->X())>1.0e-10) ? TMath::ATan2(qsub1->Y(),qsub1->X())/2.0 : 0.0); + + values[AliDielectronVarManager::kTPCsub2xH2uc] = qsub2->X(); + values[AliDielectronVarManager::kTPCsub2yH2uc] = qsub2->Y(); + values[AliDielectronVarManager::kTPCsub2rpH2uc] = ((TMath::Abs(qsub2->X())>1.0e-10) ? TMath::ATan2(qsub2->Y(),qsub2->X())/2.0 : 0.0); + + values[AliDielectronVarManager::kTPCsub12DiffH2uc] = TMath::Cos( 2.*(values[AliDielectronVarManager::kTPCsub1rpH2uc] - + values[AliDielectronVarManager::kTPCsub2rpH2uc]) ); + + // using corrected tpc quantities + values[AliDielectronVarManager::kV0ATPCDiffH2] = TMath::Cos( 2.*(values[AliDielectronVarManager::kV0ArpH2] - + values[AliDielectronVarManager::kTPCrpH2]) ); + + values[AliDielectronVarManager::kV0CTPCDiffH2] = TMath::Cos( 2.*(values[AliDielectronVarManager::kV0CrpH2] - + values[AliDielectronVarManager::kTPCrpH2]) ); + + values[AliDielectronVarManager::kV0AV0CDiffH2] = TMath::Cos( 2.*(values[AliDielectronVarManager::kV0ArpH2] - + values[AliDielectronVarManager::kV0CrpH2]) ); + + values[AliDielectronVarManager::kv0ATPCDiffH2] = TMath::Cos( 2.*(values[AliDielectronVarManager::kv0ArpH2] - + values[AliDielectronVarManager::kTPCrpH2]) ); + + values[AliDielectronVarManager::kv0CTPCDiffH2] = TMath::Cos( 2.*(values[AliDielectronVarManager::kv0CrpH2] - + values[AliDielectronVarManager::kTPCrpH2]) ); + + values[AliDielectronVarManager::kv0Av0CDiffH2] = TMath::Cos( 2.*(values[AliDielectronVarManager::kv0ArpH2] - + values[AliDielectronVarManager::kv0CrpH2]) ); + } inline void AliDielectronVarManager::FillVarESDEvent(const AliESDEvent *event, Double_t * const values) @@ -865,6 +1262,44 @@ inline void AliDielectronVarManager::FillVarESDEvent(const AliESDEvent *event, D values[AliDielectronVarManager::kYRes] = primVtx->GetYRes(); values[AliDielectronVarManager::kZRes] = primVtx->GetZRes(); values[AliDielectronVarManager::kCentrality] = centralityF; + + // Event multiplicity estimators + Int_t nTrSPD05=0; Int_t nTrITSTPC05=0; Int_t nTrITSSA05=0; + event->EstimateMultiplicity(nTrSPD05, nTrITSTPC05, nTrITSSA05, 0.5); + values[AliDielectronVarManager::kNaccTrckltsEsd05] = nTrSPD05; + values[AliDielectronVarManager::kNaccItsTpcEsd05] = nTrITSTPC05; + values[AliDielectronVarManager::kNaccItsPureEsd05] = nTrITSSA05; + values[AliDielectronVarManager::kNaccTrckltsEsd05Corr] = + AliDielectronHelper::GetNaccTrckltsCorrected(event,Double_t(nTrSPD05),values[AliDielectronVarManager::kZvPrim],0); + values[AliDielectronVarManager::kNaccItsTpcEsd05Corr] = + AliDielectronHelper::GetNaccTrckltsCorrected(event,Double_t(nTrITSTPC05),values[AliDielectronVarManager::kZvPrim],3); + values[AliDielectronVarManager::kNaccItsPureEsd05Corr] = + AliDielectronHelper::GetNaccTrckltsCorrected(event,Double_t(nTrITSSA05),values[AliDielectronVarManager::kZvPrim],6); + + Int_t nTrSPD10=0; Int_t nTrITSTPC10=0; Int_t nTrITSSA10=0; + event->EstimateMultiplicity(nTrSPD10, nTrITSTPC10, nTrITSSA10, 1.0); + values[AliDielectronVarManager::kNaccTrckltsEsd10] = nTrSPD10; + values[AliDielectronVarManager::kNaccItsTpcEsd10] = nTrITSTPC10; + values[AliDielectronVarManager::kNaccItsPureEsd10] = nTrITSSA10; + values[AliDielectronVarManager::kNaccTrckltsEsd10Corr] = + AliDielectronHelper::GetNaccTrckltsCorrected(event,Double_t(nTrSPD10),values[AliDielectronVarManager::kZvPrim],1); + values[AliDielectronVarManager::kNaccItsTpcEsd10Corr] = + AliDielectronHelper::GetNaccTrckltsCorrected(event,Double_t(nTrITSTPC10),values[AliDielectronVarManager::kZvPrim],4); + values[AliDielectronVarManager::kNaccItsPureEsd10Corr] = + AliDielectronHelper::GetNaccTrckltsCorrected(event,Double_t(nTrITSSA10),values[AliDielectronVarManager::kZvPrim],7); + + Int_t nTrSPD16=0; Int_t nTrITSTPC16=0; Int_t nTrITSSA16=0; + event->EstimateMultiplicity(nTrSPD16, nTrITSTPC16, nTrITSSA16, 1.6); + values[AliDielectronVarManager::kNaccTrckltsEsd16] = nTrSPD16; + values[AliDielectronVarManager::kNaccItsTpcEsd16] = nTrITSTPC16; + values[AliDielectronVarManager::kNaccItsPureEsd16] = nTrITSSA16; + values[AliDielectronVarManager::kNaccTrckltsEsd16Corr] = + AliDielectronHelper::GetNaccTrckltsCorrected(event,Double_t(nTrSPD16),values[AliDielectronVarManager::kZvPrim],2); + values[AliDielectronVarManager::kNaccItsTpcEsd16Corr] = + AliDielectronHelper::GetNaccTrckltsCorrected(event,Double_t(nTrITSTPC16),values[AliDielectronVarManager::kZvPrim],5); + values[AliDielectronVarManager::kNaccItsPureEsd16Corr] = + AliDielectronHelper::GetNaccTrckltsCorrected(event,Double_t(nTrITSSA16),values[AliDielectronVarManager::kZvPrim],8); + } inline void AliDielectronVarManager::FillVarAODEvent(const AliAODEvent *event, Double_t * const values) @@ -886,11 +1321,49 @@ inline void AliDielectronVarManager::FillVarMCEvent(const AliMCEvent *event, Dou // // Fill common AliVEvent interface information - FillVarVEvent(event, values); - + // FillVarVEvent(event, values); + const AliVVertex* vtx = event->GetPrimaryVertex(); + values[AliDielectronVarManager::kXvPrim] = (vtx ? vtx->GetX() : 0.0); + values[AliDielectronVarManager::kYvPrim] = (vtx ? vtx->GetY() : 0.0); + values[AliDielectronVarManager::kZvPrim] = (vtx ? vtx->GetZ() : 0.0); // Fill AliMCEvent interface specific information - values[AliDielectronVarManager::kNch] = AliDielectronHelper::GetNch(event, 1.6); -} + values[AliDielectronVarManager::kNch] = AliDielectronHelper::GetNch(event, 1.6); + values[AliDielectronVarManager::kNch05] = AliDielectronHelper::GetNch(event, 0.5); + values[AliDielectronVarManager::kNch10] = AliDielectronHelper::GetNch(event, 1.0); + + values[AliDielectronVarManager::kNumberOfJPsis] = AliDielectronHelper::GetNMothers(event, 0.9, 443, 11); + values[AliDielectronVarManager::kNumberOfJPsisPrompt] = AliDielectronHelper::GetNMothers(event, 0.9, 443, 11, 1); + values[AliDielectronVarManager::kNumberOfJPsisNPrompt] = AliDielectronHelper::GetNMothers(event, 0.9, 443, 11, 0); +} + +inline void AliDielectronVarManager::FillVarTPCEventPlane(const AliEventplane *evplane, Double_t * const values) +{ + // + // Fill TPC event plane information after correction + // + if(!evplane) return; + TVector2 *qcorr = const_cast(evplane)->GetQVector(); // This is the "corrected" Q-Vector + TVector2 *qcsub1 = const_cast(evplane)->GetQsub1(); + TVector2 *qcsub2 = const_cast(evplane)->GetQsub2(); + if(!qcorr || !qcsub1 || !qcsub2) return; + + values[AliDielectronVarManager::kTPCxH2] = qcorr->X(); + values[AliDielectronVarManager::kTPCyH2] = qcorr->Y(); + values[AliDielectronVarManager::kTPCrpH2] = ((TMath::Abs(qcorr->X())>1.0e-10) ? TMath::ATan2(qcorr->Y(),qcorr->X())/2.0 : 0.0); + + values[AliDielectronVarManager::kTPCsub1xH2] = qcsub1->X(); + values[AliDielectronVarManager::kTPCsub1yH2] = qcsub1->Y(); + values[AliDielectronVarManager::kTPCsub1rpH2] = ((TMath::Abs(qcsub1->X())>1.0e-10) ? TMath::ATan2(qcsub1->Y(),qcsub1->X())/2.0 : 0.0); + + values[AliDielectronVarManager::kTPCsub2xH2] = qcsub2->X(); + values[AliDielectronVarManager::kTPCsub2yH2] = qcsub2->Y(); + values[AliDielectronVarManager::kTPCsub2rpH2] = ((TMath::Abs(qcsub2->X())>1.0e-10) ? TMath::ATan2(qcsub2->Y(),qcsub2->X())/2.0 : 0.0); + + values[AliDielectronVarManager::kTPCsub12DiffH2] = TMath::Cos( 2.*(values[AliDielectronVarManager::kTPCsub1rpH2] - + values[AliDielectronVarManager::kTPCsub2rpH2]) ); + values[AliDielectronVarManager::kTPCsub12DiffH2Sin] = TMath::Sin( 2.*(values[AliDielectronVarManager::kTPCsub1rpH2] - + values[AliDielectronVarManager::kTPCsub2rpH2]) ); +} inline void AliDielectronVarManager::InitESDpid(Int_t type) { @@ -958,6 +1431,124 @@ inline void AliDielectronVarManager::InitAODpidUtil(Int_t type) } +inline void AliDielectronVarManager::InitEstimatorAvg(const Char_t* filename) +{ + // + // initialize the profile histograms neccessary for the correction of the multiplicity estimators in pp collisions + // + + const Char_t* estimatorNames[9] = {"SPDmult05","SPDmult10","SPDmult16", + "ITSTPC05", "ITSTPC10", "ITSTPC16", + "ITSSA05", "ITSSA10", "ITSSA16"}; + const Char_t* periodNames[4] = {"LHC10b", "LHC10c", "LHC10d", "LHC10e"}; + TFile* file=TFile::Open(filename); + if(!file) return; + + for(Int_t ip=0; ip<4; ++ip) { + for(Int_t ie=0; ie<9; ++ie) { + fgMultEstimatorAvg[ip][ie] = (TProfile*)(file->Get(Form("%s_%s",estimatorNames[ie],periodNames[ip]))->Clone(Form("%s_%s_clone",estimatorNames[ie],periodNames[ip]))); + } + } +} + + +inline void AliDielectronVarManager::InitTRDpidEffHistograms(const Char_t* filename) +{ + // + // initialize the 3D histograms with the TRD pid efficiency histograms + // + + // reset the centrality ranges and the efficiency histograms + for(Int_t i=0; i<10; ++i) { // centrality ranges + for(Int_t j=0; j<4; ++j) fgTRDpidEffCentRanges[i][j] = -1.; + if(fgTRDpidEff[i][0]) { + delete fgTRDpidEff[i][0]; + fgTRDpidEff[i][0] = 0x0; + } + if(fgTRDpidEff[i][1]) { + delete fgTRDpidEff[i][1]; + fgTRDpidEff[i][1] = 0x0; + } + } + + TFile* file=TFile::Open(filename); + TList* keys=file->GetListOfKeys(); + Int_t idxp=0; Int_t idxn=0; + for(Int_t i=0; iGetEntries(); ++i) { + if(idxp>=10) continue; + if(idxn>=10) continue; + TString name=((TKey*)keys->At(i))->ReadObj()->GetName(); + // Name of histograms should be in the format: + // TRDeff_cent__ + // is either "BPLUS" or "BMINUS" + if(!(name.Contains("BPLUS") || name.Contains("BMINUS"))) continue; + TObjArray* arr = name.Tokenize("_"); + Bool_t isBplus = kTRUE; + if(name.Contains("BMINUS")) isBplus = kFALSE; + TString centMinStr = arr->At(2)->GetName(); TString centMaxStr = arr->At(3)->GetName(); + if(isBplus) { + fgTRDpidEffCentRanges[idxp][2] = centMinStr.Atof(); + fgTRDpidEffCentRanges[idxp][3] = centMaxStr.Atof(); + fgTRDpidEff[idxp][1] = (TH3D*)(file->Get(name.Data())->Clone(Form("%s_clone",name.Data()))); + ++idxp; + } + else { + fgTRDpidEffCentRanges[idxn][0] = centMinStr.Atof(); + fgTRDpidEffCentRanges[idxn][1] = centMaxStr.Atof(); + fgTRDpidEff[idxn][0] = (TH3D*)(file->Get(name.Data())->Clone(Form("%s_clone",name.Data()))); + ++idxn; + } + } +} + + +inline Double_t AliDielectronVarManager::GetTRDpidEfficiency(Int_t runNo, Double_t centrality, + Double_t eta, Double_t trdPhi, Double_t pout, + Double_t& effErr) { + // + // return the efficiency in the given phase space cell + // + // LHC10h data---------------------------------------------- + Bool_t isBplus = kTRUE; + if(runNo<=138275) isBplus = kFALSE; + // TODO: check magnetic polarity for runs in 2011 data + // --------------------------------------------------------- + Int_t centIdx = -1; + for(Int_t icent=0; icent<10; ++icent) { + if(isBplus) { + if(centrality>=fgTRDpidEffCentRanges[icent][2] && centrality=fgTRDpidEffCentRanges[icent][0] && centralityGetXaxis()->FindBin(eta); + if(etaGetXaxis()->GetXmin()) etaBin=1; + if(eta>effH->GetXaxis()->GetXmax()) etaBin=effH->GetXaxis()->GetNbins(); + Int_t phiBin = effH->GetYaxis()->FindBin(trdPhi); + if(trdPhiGetYaxis()->GetXmin()) phiBin=1; + if(trdPhi>effH->GetYaxis()->GetXmax()) phiBin=effH->GetYaxis()->GetNbins(); + Int_t poutBin = effH->GetZaxis()->FindBin(pout); + if(poutGetZaxis()->GetXmin()) poutBin=1; + if(pout>effH->GetZaxis()->GetXmax()) poutBin=effH->GetZaxis()->GetNbins(); + Double_t eff = effH->GetBinContent(etaBin, phiBin, poutBin); + effErr = effH->GetBinError(etaBin, phiBin, poutBin); + if(eff<-0.0001) { + effErr = 0.0; + eff = 1.0; + } + return eff; +} + + inline void AliDielectronVarManager::SetEvent(AliVEvent * const ev) { @@ -970,6 +1561,12 @@ inline void AliDielectronVarManager::SetEvent(AliVEvent * const ev) AliDielectronVarManager::Fill(fgEvent, fgData); } +inline void AliDielectronVarManager::SetEventData(const Double_t data[AliDielectronVarManager::kNMaxValues]) +{ + for (Int_t i=0; iGetVZEROData(); + + for(Int_t iChannel=0; iChannel<64; ++iChannel) { + if(iChannel<32 && sideOption==0) continue; + if(iChannel>=32 && sideOption==1) continue; + phi=iChannel%8; + //mult = vzero->GetMultiplicity(iChannel); + mult = event->GetVZEROEqMultiplicity(iChannel); + // 2nd harmonic + qvec[0] += mult*(2.0*TMath::Power(kX[phi],2.0)-1); + qvec[1] += mult*(2.0*kX[phi]*kY[phi]); + } // end loop over channels + if(TMath::Abs(qvec[0])>1.0e-10) + qvec[2] = TMath::ATan2(qvec[1],qvec[0])/2.0; +} + + + /* inline void AliDielectronVarManager::FillValues(const TParticle *particle, Double_t *values) {