X-Git-Url: http://git.uio.no/git/?p=u%2Fmrichter%2FAliRoot.git;a=blobdiff_plain;f=PWGHF%2Fhfe%2FAliAnalysisTaskHFECal.cxx;h=e7d2717092b901e67b6aa9ebdc2111ebd23768f8;hp=5dada5b6267e0266c308ee3ad4c278be1554b3d0;hb=32abf0b9c3fa3a93907a70b3b950f8a87a0bf667;hpb=f4bb86d044bd0fc2df2b843673458661018a439e diff --git a/PWGHF/hfe/AliAnalysisTaskHFECal.cxx b/PWGHF/hfe/AliAnalysisTaskHFECal.cxx index 5dada5b6267..e7d2717092b 100644 --- a/PWGHF/hfe/AliAnalysisTaskHFECal.cxx +++ b/PWGHF/hfe/AliAnalysisTaskHFECal.cxx @@ -29,6 +29,8 @@ #include "TFile.h" #include "TGraphErrors.h" +#include "TDatabasePDG.h" + #include "AliAnalysisTask.h" #include "AliAnalysisManager.h" @@ -93,7 +95,7 @@ AliAnalysisTaskHFECal::AliAnalysisTaskHFECal(const char *name) ,stack(0) ,fGeom(0) ,fOutputList(0) - ,fqahist(1) + ,fqahist(0) ,fTrackCuts(0) ,fCuts(0) ,fIdentifiedAsOutInz(kFALSE) @@ -105,8 +107,13 @@ AliAnalysisTaskHFECal::AliAnalysisTaskHFECal(const char *name) ,fPID(0) ,fPIDqa(0) ,fOpeningAngleCut(0.1) + ,fMimpTassCut(0.5) + ,fMimNsigassCut(-3) ,fInvmassCut(0) // no mass ,fSetMassConstraint(kTRUE) + ,fSetMassWidthCut(kFALSE) + ,fSetMassNonlinear(kFALSE) + ,fSetKFpart(kTRUE) ,fNoEvents(0) ,fEMCAccE(0) ,hEMCAccE(0) @@ -136,8 +143,11 @@ AliAnalysisTaskHFECal::AliAnalysisTaskHFECal(const char *name) ,fPhotoElecPt(0) ,fPhoElecPt(0) ,fPhoElecPtM20(0) + ,fPhoElecPtM20Mass(0) ,fSameElecPt(0) ,fSameElecPtM20(0) + ,fSameElecPtM20Mass(0) + ,fSemiElecPtM20(0) ,fTrackPtBefTrkCuts(0) ,fTrackPtAftTrkCuts(0) ,fTPCnsigma(0) @@ -169,8 +179,10 @@ AliAnalysisTaskHFECal::AliAnalysisTaskHFECal(const char *name) ,fIncpTMCM20pho(0) ,fPhoElecPtMC(0) ,fPhoElecPtMCM20(0) + ,fPhoElecPtMCM20Mass(0) ,fSameElecPtMC(0) ,fSameElecPtMCM20(0) + ,fSameElecPtMCM20Mass(0) ,fIncpTMCM20pho_pi0e(0) ,fPhoElecPtMCM20_pi0e(0) ,fSameElecPtMCM20_pi0e(0) @@ -180,8 +192,12 @@ AliAnalysisTaskHFECal::AliAnalysisTaskHFECal(const char *name) ,fIncpTMCpho_pi0e_TPC(0) ,fPhoElecPtMC_pi0e_TPC(0) ,fSameElecPtMC_pi0e_TPC(0) + ,fIncpTMCpho_eta_TPC(0) + ,fPhoElecPtMC_eta_TPC(0) + ,fSameElecPtMC_eta_TPC(0) ,CheckNclust(0) ,CheckNits(0) + ,CheckDCA(0) ,Hpi0pTcheck(0) ,HETApTcheck(0) ,HphopTcheck(0) @@ -194,9 +210,31 @@ AliAnalysisTaskHFECal::AliAnalysisTaskHFECal(const char *name) ,fFakeTrk0(0) ,fFakeTrk1(0) ,ftimingEle(0) + ,fIncMaxE(0) ,fIncReco(0) ,fPhoReco(0) ,fSamReco(0) + ,fIncRecoMaxE(0) + ,fPhoRecoMaxE(0) + ,fSamRecoMaxE(0) + ,fPhoVertexReco_TPC(0) + ,fPhoVertexReco_TPC_Invmass(0) + ,fPhoVertexReco_EMCal(0) + ,fPhoVertexReco_Invmass(0) + ,fPhoVertexReco_step0(0) + ,fPhoVertexReco_step1(0) + ,fMatchV0_0(0) + ,fMatchV0_1(0) + ,fMatchMC_0(0) + ,fMatchMC_1(0) + ,fpair(0) + ,fFakeRejection0(0) + ,fFakeRejection1(0) + ,fFakeRejection2(0) + ,EopFake(0) + ,EopTrue(0) + ,MatchFake(0) + ,MatchTrue(0) //,fnSigEtaCorr(NULL) { //Named constructor @@ -224,7 +262,7 @@ AliAnalysisTaskHFECal::AliAnalysisTaskHFECal() ,stack(0) ,fGeom(0) ,fOutputList(0) - ,fqahist(1) + ,fqahist(0) ,fTrackCuts(0) ,fCuts(0) ,fIdentifiedAsOutInz(kFALSE) @@ -236,8 +274,13 @@ AliAnalysisTaskHFECal::AliAnalysisTaskHFECal() ,fPID(0) ,fPIDqa(0) ,fOpeningAngleCut(0.1) + ,fMimpTassCut(0.5) + ,fMimNsigassCut(-3) ,fInvmassCut(0) // no mass ,fSetMassConstraint(kTRUE) + ,fSetMassWidthCut(kFALSE) + ,fSetMassNonlinear(kFALSE) + ,fSetKFpart(kTRUE) ,fNoEvents(0) ,fEMCAccE(0) ,hEMCAccE(0) @@ -267,8 +310,11 @@ AliAnalysisTaskHFECal::AliAnalysisTaskHFECal() ,fPhotoElecPt(0) ,fPhoElecPt(0) ,fPhoElecPtM20(0) + ,fPhoElecPtM20Mass(0) ,fSameElecPt(0) ,fSameElecPtM20(0) + ,fSameElecPtM20Mass(0) + ,fSemiElecPtM20(0) ,fTrackPtBefTrkCuts(0) ,fTrackPtAftTrkCuts(0) ,fTPCnsigma(0) @@ -300,8 +346,10 @@ AliAnalysisTaskHFECal::AliAnalysisTaskHFECal() ,fIncpTMCM20pho(0) ,fPhoElecPtMC(0) ,fPhoElecPtMCM20(0) + ,fPhoElecPtMCM20Mass(0) ,fSameElecPtMC(0) ,fSameElecPtMCM20(0) + ,fSameElecPtMCM20Mass(0) ,fIncpTMCM20pho_pi0e(0) ,fPhoElecPtMCM20_pi0e(0) ,fSameElecPtMCM20_pi0e(0) @@ -311,8 +359,12 @@ AliAnalysisTaskHFECal::AliAnalysisTaskHFECal() ,fIncpTMCpho_pi0e_TPC(0) ,fPhoElecPtMC_pi0e_TPC(0) ,fSameElecPtMC_pi0e_TPC(0) + ,fIncpTMCpho_eta_TPC(0) + ,fPhoElecPtMC_eta_TPC(0) + ,fSameElecPtMC_eta_TPC(0) ,CheckNclust(0) ,CheckNits(0) + ,CheckDCA(0) ,Hpi0pTcheck(0) ,HETApTcheck(0) ,HphopTcheck(0) @@ -325,9 +377,31 @@ AliAnalysisTaskHFECal::AliAnalysisTaskHFECal() ,fFakeTrk0(0) ,fFakeTrk1(0) ,ftimingEle(0) + ,fIncMaxE(0) ,fIncReco(0) ,fPhoReco(0) ,fSamReco(0) + ,fIncRecoMaxE(0) + ,fPhoRecoMaxE(0) + ,fSamRecoMaxE(0) + ,fPhoVertexReco_TPC(0) + ,fPhoVertexReco_TPC_Invmass(0) + ,fPhoVertexReco_EMCal(0) + ,fPhoVertexReco_Invmass(0) + ,fPhoVertexReco_step0(0) + ,fPhoVertexReco_step1(0) + ,fMatchV0_0(0) + ,fMatchV0_1(0) + ,fMatchMC_0(0) + ,fMatchMC_1(0) + ,fpair(0) + ,fFakeRejection0(0) + ,fFakeRejection1(0) + ,fFakeRejection2(0) + ,EopFake(0) + ,EopTrue(0) + ,MatchFake(0) + ,MatchTrue(0) //,fnSigEtaCorr(NULL) { //Default constructor @@ -417,8 +491,6 @@ void AliAnalysisTaskHFECal::UserExec(Option_t*) mcphoinfo[0] = cent; mcphoinfo[1] = pTMC; mcphoinfo[2] = iHijing; - //if(fPDG==111)Hpi0pTcheck->Fill(pTMC,iHijing); - //if(fPDG==221)HETApTcheck->Fill(pTMC,iHijing); if(fPDG==111)Hpi0pTcheck->Fill(mcphoinfo); if(fPDG==221)HETApTcheck->Fill(mcphoinfo); if(fabs(fPDG)==411 || fabs(fPDG)==413 || fabs(fPDG)==421 || fabs(fPDG)==423 || fabs(fPDG)==431)HDpTcheck->Fill(pTMC,iHijing); @@ -479,9 +551,6 @@ void AliAnalysisTaskHFECal::UserExec(Option_t*) fCFM->SetRecEventInfo(fESD); - //Float_t cent = -1.; - //AliCentrality *centrality = fESD->GetCentrality(); - //cent = centrality->GetCentralityPercentile("V0M"); fCent->Fill(cent); //if(cent>90.) return; @@ -491,6 +560,7 @@ void AliAnalysisTaskHFECal::UserExec(Option_t*) FindTriggerClusters(); // make EMCAL array + double maxE = 0.0; for(Int_t iCluster=0; iClusterGetNumberOfCaloClusters(); iCluster++) { AliESDCaloCluster *clust = fESD->GetCaloCluster(iCluster); @@ -505,9 +575,8 @@ void AliAnalysisTaskHFECal::UserExec(Option_t*) double calInfo[5]; calInfo[0] = emcphi; calInfo[1] = emceta; calInfo[2] = clustE; calInfo[3] = cent; calInfo[4] = clust->Chi2(); //fEMCAccE->Fill(calInfo); - //if(clustE>3.0)fEMCAccE->Fill(calInfo); - //if(fqahist==1 && clustE>1.5)fEMCAccE->Fill(calInfo); hEMCAccE->Fill(cent,clustE); + if(clustE>maxE)maxE = clustE; } } @@ -519,6 +588,8 @@ void AliAnalysisTaskHFECal::UserExec(Option_t*) continue; } + //--- Get MC informtion + int parentlabel = 99999; int parentPID = 99999; int grand_parentlabel = 99999; @@ -533,6 +604,7 @@ void AliAnalysisTaskHFECal::UserExec(Option_t*) double mcMompT = 0.0; //double mcGrandMompT = 0.0; double mcWeight = -10.0; + double conv_proR = -1.0; int iHijing = 1; int mcLabel = -1; @@ -540,7 +612,8 @@ void AliAnalysisTaskHFECal::UserExec(Option_t*) if(fmcData && fMC && stack) { Int_t label = TMath::Abs(track->GetLabel()); - mcLabel = track->GetLabel(); + //mcLabel = track->GetLabel(); + mcLabel = abs(track->GetLabel()); // check for conv. issue if(mcLabel>-1) { @@ -551,6 +624,7 @@ void AliAnalysisTaskHFECal::UserExec(Option_t*) TParticle* particle = stack->Particle(label); int mcpid = particle->GetPdgCode(); mcpT = particle->Pt(); + conv_proR = particle->R(); //printf("MCpid = %d",mcpid); if(particle->GetFirstMother()>-1) { @@ -656,15 +730,29 @@ void AliAnalysisTaskHFECal::UserExec(Option_t*) } // end of mcLabel>-1 - } // end of MC info. + } // <------ end of MC info. //cout << "Pi0 = " << mcOrgPi0 << " ; Eta = " << mcOrgEta << endl; //printf("weight = %f\n",mcWeight); if(TMath::Abs(track->Eta())>0.6) continue; - if(TMath::Abs(track->Pt()<2.0)) continue; + //if(TMath::Abs(track->Pt()<2.5)) continue; + if(TMath::Abs(track->Pt()<0.1)) continue; + int nITS = track->GetNcls(0); + fTrackPtBefTrkCuts->Fill(track->Pt()); + + /* + UChar_t itsPixel = track->GetITSClusterMap(); + cout << "nITS = " << nITS << endl; + if(itsPixel & BIT(0))cout << "1st layer hit" << endl; + if(itsPixel & BIT(1))cout << "2nd layer hit" << endl; + if(itsPixel & BIT(2))cout << "3rd layer hit" << endl; + if(itsPixel & BIT(3))cout << "4th layer hit" << endl; + if(itsPixel & BIT(4))cout << "5th layer hit" << endl; + */ + // RecKine: ITSTPC cuts if(!ProcessCutStep(AliHFEcuts::kStepRecKineITSTPC, track)) continue; @@ -682,47 +770,55 @@ void AliAnalysisTaskHFECal::UserExec(Option_t*) // HFE cuts: TPC PID cleanup if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsTPC, track)) continue; + if(mcPho && iHijing==0)fPhoVertexReco_step0->Fill(track->Pt(),conv_proR); // check MC vertex + if(mcPho && iHijing==1)fPhoVertexReco_step1->Fill(track->Pt(),conv_proR); // check MC vertex + int nTPCcl = track->GetTPCNcls(); //int nTPCclF = track->GetTPCNclsF(); // warnings - int nITS = track->GetNcls(0); - + //int nITS = track->GetNcls(0); + + fTrackPtAftTrkCuts->Fill(track->Pt()); Double_t mom = -999., eop=-999., pt = -999., dEdx=-999., fTPCnSigma=-10, phi=-999., eta=-999.; pt = track->Pt(); - if(pt<2.0)continue; - - // Track extrapolation + if(pt<2.5)continue; + if(pt<0.1)continue; - Int_t charge = track->Charge(); + //Int_t charge = track->Charge(); fTrkpt->Fill(pt); mom = track->P(); phi = track->Phi(); eta = track->Eta(); + float dca_xy; + float dca_z; + track->GetImpactParameters(dca_xy,dca_z); + dEdx = track->GetTPCsignal(); fTPCnSigma = fPID->GetPIDResponse() ? fPID->GetPIDResponse()->NumberOfSigmasTPC(track, AliPID::kElectron) : 1000; //cout << "nSigma correctoon-----" << endl; //cout << "org = " << fTPCnSigma << endl; + /* if(!fmcData) // nsigma eta correction { double nSigexpCorr = NsigmaCorrection(eta,cent); fTPCnSigma -= nSigexpCorr; } - + */ //cout << "correction = " << fTPCnSigma << endl; + //--- track cluster match + double ncells = -1.0; double m20 = -1.0; double m02 = -1.0; double disp = -1.0; double rmatch = -1.0; double nmatch = -1.0; - double oppstatus = 0.0; + //double oppstatus = 0.0; double emctof = 0.0; - - Bool_t fFlagPhotonicElec = kFALSE; - Bool_t fFlagConvinatElec = kFALSE; + Bool_t MaxEmatch = kFALSE; Int_t clsId = track->GetEMCALcluster(); if (clsId>0){ @@ -730,14 +826,8 @@ void AliAnalysisTaskHFECal::UserExec(Option_t*) if(clust && clust->IsEMCAL()){ double clustE = clust->E(); + if(clustE==maxE)MaxEmatch = kTRUE; eop = clustE/fabs(mom); - //cout << "eop org = "<< eop << endl; - if(mcLabel>-1.0) - { - double mceopcorr = MCEopMeanCorrection(pt,cent); - eop += mceopcorr; - } - //cout << "eop corr = " << eop << endl; //double clustT = clust->GetTOF(); ncells = clust->GetNCells(); @@ -750,28 +840,27 @@ void AliAnalysisTaskHFECal::UserExec(Option_t*) nmatch = clust->GetNTracksMatched(); emctof = clust->GetTOF(); //cout << "emctof = " << emctof << endl; - - if(fTPCnSigma>-1.5 && fTPCnSigma<3.0) - { - SelectPhotonicElectron(iTracks,cent,track,fFlagPhotonicElec,fFlagConvinatElec,fTPCnSigma,m20,eop,mcele,mcWeight,iHijing,mcOrgPi0,mcOrgEta); - } - if(fFlagPhotonicElec)oppstatus = 1.0; - if(fFlagConvinatElec)oppstatus = 2.0; - if(fFlagPhotonicElec && fFlagConvinatElec)oppstatus = 3.0; + cout << "eop org = "<< eop << endl; + double eoporg = eop; + if(fmcData) + { + double mceopcorr = MCEopMeanCorrection(pt,cent); + eop += mceopcorr; + } + cout << "eop corr = " << eop << endl; double valdedx[16]; valdedx[0] = pt; valdedx[1] = nITS; valdedx[2] = phi; valdedx[3] = eta; valdedx[4] = fTPCnSigma; - //valdedx[5] = eop; valdedx[6] = rmatch; valdedx[7] = ncells, valdedx[8] = nTPCclF; valdedx[9] = m20; valdedx[10] = mcpT; - valdedx[5] = eop; valdedx[6] = rmatch; valdedx[7] = ncells, valdedx[8] = nmatch; valdedx[9] = m20; valdedx[10] = mcpT; - valdedx[11] = cent; valdedx[12] = dEdx; valdedx[13] = oppstatus; valdedx[14] = nTPCcl; + valdedx[5] = eop; valdedx[6] = rmatch; valdedx[7] = dca_xy, valdedx[8] = dca_z; valdedx[9] = m20; valdedx[10] = mcpT; + valdedx[11] = cent; valdedx[12] = dEdx; valdedx[13] = eoporg; valdedx[14] = nTPCcl; valdedx[15] = mcele; - if(fqahist==1)fEleInfo->Fill(valdedx); - + fEleInfo->Fill(valdedx); } } //Get Cal info PID response + /* double eop2; double ss[4]; Double_t nSigmaEop = fPID->GetPIDResponse()->NumberOfSigmasEMCAL(track,AliPID::kElectron,eop2,ss); @@ -783,39 +872,26 @@ void AliAnalysisTaskHFECal::UserExec(Option_t*) valEop[2] = nSigmaEop; fElenSigma->Fill(valEop); } - - // ============ PID + */ + // --- tarck cut & e ID if(nITS<2.5)continue; if(nTPCcl<100)continue; - + + CheckNclust->Fill(nTPCcl); CheckNits->Fill(nITS); + CheckDCA->Fill(dca_xy,dca_z); + // check production vertex of photons + fdEdxBef->Fill(mom,fTPCnSigma); fTPCnsigma->Fill(mom,fTPCnSigma); if(fTPCnSigma >= -1.0 && fTPCnSigma <= 3)fTrkEovPBef->Fill(pt,eop); - Int_t pidpassed = 1; - - // check reco eff. with TPC - - double phoval[5]; - phoval[0] = cent; - phoval[1] = pt; - phoval[2] = fTPCnSigma; - phoval[3] = iHijing; - phoval[4] = mcMompT; - - if((fTPCnSigma >= -1.0 && fTPCnSigma <= 3) && mcele>-1 && mcPho && mcOrgPi0) - { - if(iHijing==1)mcWeight = 1.0; - fIncpTMCpho_pi0e_TPC->Fill(phoval,mcWeight); - if(fFlagPhotonicElec) fPhoElecPtMC_pi0e_TPC->Fill(phoval,mcWeight); - if(fFlagConvinatElec) fSameElecPtMC_pi0e_TPC->Fill(phoval,mcWeight); - } + //--- track accepted by HFE - //--- track accepted + Int_t pidpassed = 1; AliHFEpidObject hfetrack; hfetrack.SetAnalysisType(AliHFEpidObject::kESDanalysis); hfetrack.SetRecTrack(track); @@ -836,8 +912,91 @@ void AliAnalysisTaskHFECal::UserExec(Option_t*) hfetrack.SetPbPb(); if(!fPID->IsSelected(&hfetrack, NULL, "", fPIDqa)) pidpassed = 0; - if(pidpassed==0) continue; - + if(pidpassed==0) continue; // nSigma rejection + + //--- photonic ID + + Bool_t fFlagPhotonicTPC = kFALSE; + Bool_t fFlagConvinatTPC = kFALSE; + SelectPhotonicElectron(iTracks,cent,track,fFlagPhotonicTPC,fFlagConvinatTPC,fTPCnSigma,m20,eop,mcele,mcWeight,iHijing,mcOrgPi0,mcOrgEta,0); + + //--- check reco eff. with TPC + double phoval[5]; + phoval[0] = cent; + phoval[1] = pt; + phoval[2] = fTPCnSigma; + phoval[3] = iHijing; + phoval[4] = mcMompT; + + if((fTPCnSigma >= -5.0 && fTPCnSigma <= 5) && (mcOrgPi0 || mcOrgEta)) + { + if(iHijing==1)mcWeight = 1.0; + if(mcOrgPi0) + { + fIncpTMCpho_pi0e_TPC->Fill(phoval,mcWeight); + if(fFlagPhotonicTPC) fPhoElecPtMC_pi0e_TPC->Fill(phoval,mcWeight); + if(fFlagConvinatTPC) fSameElecPtMC_pi0e_TPC->Fill(phoval,mcWeight); + fPhoVertexReco_TPC->Fill(track->Pt(),conv_proR,mcWeight); // check MC vertex + if(fFlagPhotonicTPC)fPhoVertexReco_TPC_Invmass->Fill(track->Pt(),conv_proR,mcWeight); // check MC vertex + } + if(mcOrgEta) + { + fIncpTMCpho_eta_TPC->Fill(phoval,mcWeight); + if(fFlagPhotonicTPC) fPhoElecPtMC_eta_TPC->Fill(phoval,mcWeight); + if(fFlagConvinatTPC) fSameElecPtMC_eta_TPC->Fill(phoval,mcWeight); + } + } + + //--- matching check + double emcphimim = 1.396; + double emcphimax = 3.14; + if(phi>emcphimim && phiFill(pt); // data + if(mcele>2.1)fMatchMC_0->Fill(pt,mcWeight); // MC + + if(eop>=0.0) // have a match + { + if(fFlagPhotonicTPC)fMatchV0_1->Fill(pt); + if(mcele>2.1)fMatchMC_1->Fill(pt,mcWeight); + } + } + + // check fake rejection + if(mcOrgPi0 || mcOrgEta) + { + double phiacc0 = 80.0/180.0*acos(-1); + double phiacc1 = 180.0/180.0*acos(-1); + int TrStat = 0; + if(phi>phiacc0 && phiGetLabel()>0)TrStat = 1; + fFakeRejection0->Fill(TrStat,pt,mcWeight); + if(eop>-1.0)fFakeRejection1->Fill(TrStat,pt,mcWeight); // have match + if(eop>0.9 && eop<1.3)fFakeRejection2->Fill(TrStat,pt,mcWeight); // have PID + + if(TrStat==0) + { + EopFake->Fill(pt,eop,mcWeight); + MatchFake->Fill(pt,rmatch,mcWeight); + } + else + { + EopTrue->Fill(pt,eop,mcWeight); + MatchTrue->Fill(pt,rmatch,mcWeight); + } + + } + } + + //+++++++ E/p cut ++++++++++++++++ + + if(eop<0.9 || eop>1.3)continue; + + Bool_t fFlagPhotonicElec = kFALSE; + Bool_t fFlagConvinatElec = kFALSE; + SelectPhotonicElectron(iTracks,cent,track,fFlagPhotonicElec,fFlagConvinatElec,fTPCnSigma,m20,eop,mcele,mcWeight,iHijing,mcOrgPi0,mcOrgEta,1); + fTrkEovPAft->Fill(pt,eop); fdEdxAft->Fill(mom,fTPCnSigma); @@ -852,31 +1011,12 @@ void AliAnalysisTaskHFECal::UserExec(Option_t*) ftimingEle->Fill(pt,emctof); if(fFlagPhotonicElec) fPhoElecPtM20->Fill(cent,pt); if(fFlagConvinatElec) fSameElecPtM20->Fill(cent,pt); + if(!fFlagPhotonicElec) fSemiElecPtM20->Fill(cent,pt); } - //-------- - - double recopT = SumpT(iTracks,track); - - if(m20>0.0 && m20<0.3) - { - fIncpTM20->Fill(cent,pt); - ftimingEle->Fill(pt,emctof); - if(fFlagPhotonicElec) fPhoElecPtM20->Fill(cent,pt); - if(fFlagConvinatElec) fSameElecPtM20->Fill(cent,pt); - - if(pt>5.0) - { - fIncReco->Fill(cent,recopT); - if(fFlagPhotonicElec) fPhoReco->Fill(cent,recopT); - if(fFlagConvinatElec) fSamReco->Fill(cent,recopT); - } - } - // MC // check label for electron candidiates - int idlabel = 1; if(mcLabel==0)idlabel = 0; fLabelCheck->Fill(pt,idlabel); @@ -897,8 +1037,6 @@ void AliAnalysisTaskHFECal::UserExec(Option_t*) if(mcBtoE || mcDtoE) // select B->e & D->e { fIncpTMChfe->Fill(cent,pt); - //if(m20>0.0 && m20<0.3)fIncpTMCM20hfe->Fill(cent,pt); - //if(m20>0.0 && m20<0.3)fIncpTMCM20hfeCheck->Fill(cent,mcpT); if(m20>0.0 && m20<0.3) { //cout << "MC label = " << mcLabel << endl; @@ -929,20 +1067,22 @@ void AliAnalysisTaskHFECal::UserExec(Option_t*) fIncpTMCM20pho_pi0e->Fill(phoval,mcWeight); if(fFlagPhotonicElec) fPhoElecPtMCM20_pi0e->Fill(phoval,mcWeight); if(fFlagConvinatElec) fSameElecPtMCM20_pi0e->Fill(phoval,mcWeight); - //fIncpTMCM20pho_pi0e->Fill(phoval); // v5-04-02-AN & v5-04-06-AN - //if(fFlagPhotonicElec) fPhoElecPtMCM20_pi0e->Fill(phoval); - //if(fFlagConvinatElec) fSameElecPtMCM20_pi0e->Fill(phoval); - } + + // check production vertex + //fPhoVertexReco_EMCal->Fill(track->Pt(),conv_proR); + //if(fFlagPhotonicElec)fPhoVertexReco_Invmass->Fill(track->Pt(),conv_proR); + fPhoVertexReco_EMCal->Fill(track->Pt(),conv_proR,mcWeight); + if(fFlagPhotonicElec)fPhoVertexReco_Invmass->Fill(track->Pt(),conv_proR,mcWeight); + + } + // --- eta if(mcOrgEta) { fIncpTMCM20pho_eta->Fill(phoval,mcWeight); if(fFlagPhotonicElec) fPhoElecPtMCM20_eta->Fill(phoval,mcWeight); if(fFlagConvinatElec) fSameElecPtMCM20_eta->Fill(phoval,mcWeight); - //fIncpTMCM20pho_eta->Fill(phoval); - //if(fFlagPhotonicElec) fPhoElecPtMCM20_eta->Fill(phoval); - //if(fFlagConvinatElec) fSameElecPtMCM20_eta->Fill(phoval); } - // --- eta + // check production vertex } } } @@ -981,7 +1121,7 @@ void AliAnalysisTaskHFECal::UserCreateOutputObjects() if(!fPID->GetNumberOfPIDdetectors()) { fPID->AddDetector("TPC", 0); - fPID->AddDetector("EMCAL", 1); + //fPID->AddDetector("EMCAL", 1); <--- apply PID selection in task } Double_t params[4]; @@ -1010,12 +1150,15 @@ void AliAnalysisTaskHFECal::UserCreateOutputObjects() fCuts->SetTPCmodes(AliHFEextraCuts::kFound, AliHFEextraCuts::kFoundOverFindable); //fCuts->SetMinNClustersITS(3); fCuts->SetMinNClustersITS(2); + fCuts->SetProductionVertex(0,50,0,50); fCuts->SetCutITSpixel(AliHFEextraCuts::kAny); fCuts->SetCheckITSLayerStatus(kFALSE); fCuts->SetVertexRange(10.); fCuts->SetTOFPIDStep(kFALSE); - fCuts->SetPtRange(2, 50); - fCuts->SetMaxImpactParam(3.,3.); + //fCuts->SetPtRange(2, 50); + fCuts->SetPtRange(0.1, 50); + //fCuts->SetMaxImpactParam(3.,3.); + fCuts->SetMaxImpactParam(2.4,3.2); // standard in 2011 //--------Initialize correction Framework and Cuts fCFM = new AliCFManager; @@ -1043,7 +1186,7 @@ void AliAnalysisTaskHFECal::UserCreateOutputObjects() Double_t xminE[5] = {1.0, -1, 0.0, 0, -0.5}; Double_t xmaxE[5] = {3.5, 1, 100.0, 100, 9.5}; fEMCAccE = new THnSparseD("fEMCAccE","EMC acceptance & E;#eta;#phi;Energy;Centrality;trugCondition;",5,binsE,xminE,xmaxE); - if(fqahist==1)fOutputList->Add(fEMCAccE); + //if(fqahist==1)fOutputList->Add(fEMCAccE); hEMCAccE = new TH2F("hEMCAccE","Cluster Energy",200,0,100,100,0,20); fOutputList->Add(hEMCAccE); @@ -1094,43 +1237,43 @@ void AliAnalysisTaskHFECal::UserCreateOutputObjects() fInvmassULSmc = new THnSparseD("fInvmassULSmc", "Inv mass of ULS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2); nSigma; angle; m20cut; eop; MCele", 9, nBinspho,minpho, maxpho); if(fqahist==1)fOutputList->Add(fInvmassULSmc); - fInvmassLSreco = new TH2D("fInvmassLSreco", "Inv mass of LS (e,e) reco; cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,500,0,0.5 ); + fInvmassLSreco = new TH2D("fInvmassLSreco", "Inv mass of LS (e,e) reco; cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,600,-0.1,0.5 ); fInvmassLSreco->Sumw2(); fOutputList->Add(fInvmassLSreco); - fInvmassULSreco = new TH2D("fInvmassULSreco", "Inv mass of ULS (e,e) reco; cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,500,0,0.5 ); + fInvmassULSreco = new TH2D("fInvmassULSreco", "Inv mass of ULS (e,e) reco; cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,600,-0.1,0.5 ); fInvmassULSreco->Sumw2(); fOutputList->Add(fInvmassULSreco); - fInvmassLSmc0 = new TH2D("fInvmassLSmc0", "Inv mass of LS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,500,0,0.5 ); + fInvmassLSmc0 = new TH2D("fInvmassLSmc0", "Inv mass of LS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,600,-0.1,0.5 ); fInvmassLSmc0->Sumw2(); fOutputList->Add(fInvmassLSmc0); - fInvmassLSmc1 = new TH2D("fInvmassLSmc1", "Inv mass of LS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,500,0,0.5 ); + fInvmassLSmc1 = new TH2D("fInvmassLSmc1", "Inv mass of LS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,600,-0.1,0.5 ); fInvmassLSmc1->Sumw2(); fOutputList->Add(fInvmassLSmc1); - fInvmassLSmc2 = new TH2D("fInvmassLSmc2", "Inv mass of LS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,500,0,0.5 ); + fInvmassLSmc2 = new TH2D("fInvmassLSmc2", "Inv mass of LS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,600,-0.1,0.5 ); fInvmassLSmc2->Sumw2(); fOutputList->Add(fInvmassLSmc2); - fInvmassLSmc3 = new TH2D("fInvmassLSmc3", "Inv mass of LS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,500,0,0.5 ); + fInvmassLSmc3 = new TH2D("fInvmassLSmc3", "Inv mass of LS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,600,-0.1,0.5 ); fInvmassLSmc3->Sumw2(); fOutputList->Add(fInvmassLSmc3); - fInvmassULSmc0 = new TH2D("fInvmassULSmc0", "Inv mass of ULS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,500,0,0.5 ); + fInvmassULSmc0 = new TH2D("fInvmassULSmc0", "Inv mass of ULS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,600,-0.1,0.5 ); fInvmassULSmc0->Sumw2(); fOutputList->Add(fInvmassULSmc0); - fInvmassULSmc1 = new TH2D("fInvmassULSmc1", "Inv mass of ULS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,500,0,0.5 ); + fInvmassULSmc1 = new TH2D("fInvmassULSmc1", "Inv mass of ULS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,600,-0.1,0.5 ); fInvmassULSmc1->Sumw2(); fOutputList->Add(fInvmassULSmc1); - fInvmassULSmc2 = new TH2D("fInvmassULSmc2", "Inv mass of ULS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,500,0,0.5 ); + fInvmassULSmc2 = new TH2D("fInvmassULSmc2", "Inv mass of ULS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,600,-0.1,0.5 ); fInvmassULSmc2->Sumw2(); fOutputList->Add(fInvmassULSmc2); - fInvmassULSmc3 = new TH2D("fInvmassULSmc3", "Inv mass of ULS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,500,0,0.5 ); + fInvmassULSmc3 = new TH2D("fInvmassULSmc3", "Inv mass of ULS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,600,-0.1,0.5 ); fInvmassULSmc3->Sumw2(); fOutputList->Add(fInvmassULSmc3); @@ -1149,12 +1292,23 @@ void AliAnalysisTaskHFECal::UserCreateOutputObjects() fPhoElecPtM20 = new TH2F("fPhoElecPtM20", "Pho-inclusive electron pt with M20",200,0,100,100,0,50); fOutputList->Add(fPhoElecPtM20); + fPhoElecPtM20Mass = new TH2F("fPhoElecPtM20Mass", "Pho-inclusive electron pt with M20 Mass",200,0,100,100,0,50); + fPhoElecPtM20Mass->Sumw2(); + fOutputList->Add(fPhoElecPtM20Mass); + fSameElecPt = new TH2F("fSameElecPt", "Same-inclusive electron pt",200,0,100,100,0,50); fOutputList->Add(fSameElecPt); fSameElecPtM20 = new TH2F("fSameElecPtM20", "Same-inclusive electron pt with M20",200,0,100,100,0,50); fOutputList->Add(fSameElecPtM20); + fSameElecPtM20Mass = new TH2F("fSameElecPtM20Mass", "Same-inclusive electron pt with M20 Mass",200,0,100,100,0,50); + fSameElecPtM20Mass->Sumw2(); + fOutputList->Add(fSameElecPtM20Mass); + + fSemiElecPtM20 = new TH2F("fSemiElecPtM20", "Semi-inclusive electron pt with M20",200,0,100,100,0,50); + fOutputList->Add(fSemiElecPtM20); + fCent = new TH1F("fCent","Centrality",200,0,100) ; fOutputList->Add(fCent); @@ -1162,12 +1316,14 @@ void AliAnalysisTaskHFECal::UserCreateOutputObjects() const Double_t kMinP = 0.; const Double_t kMaxP = 20.; - // 1st histogram: TPC dEdx with/without EMCAL (p, pT, TPC Signal, phi, eta, Sig, e/p, ,match, cell, M02, M20, Disp, Centrality, select) - Int_t nBins[16] = { 100, 7, 60, 20, 90, 250, 25, 40, 10, 200, 100, 100, 500, 5, 100, 8}; - Double_t min[16] = {kMinP, -0.5, 1.0, -1.0, -5.0, 0, 0, 0, 0.0, 0.0, 0.0, 0, 0, -0.5, 80, -1.5}; - Double_t max[16] = {kMaxP, 6.5, 4.0, 1.0, 4.0, 2.5, 0.05, 40, 10, 2.0, 20.0, 100, 100, 4.5, 180, 6.5}; - fEleInfo = new THnSparseD("fEleInfo", "Electron Info; pT [GeV/c]; TPC signal;phi;eta;nSig; E/p;Rmatch;Ncell;clsF;M20;mcpT;Centrality;charge;opp;same;trigCond;MCele", 16, nBins, min, max); - if(fqahist==1)fOutputList->Add(fEleInfo); + //+++ 1st histogram: TPC dEdx with/without EMCAL (p, pT, TPC Signal, phi, eta, Sig, e/p, ,match, cell, M02, M20, Disp, Centrality, select) + // 1st histogram: TPC dEdx with/without EMCAL (p, pT, TPC Signal, phi, eta, Sig, e/p, ,match, dca_xy, dca_z, M20, Disp, Centrality, select) + Int_t nBins[16] = { 100, 7, 60, 20, 90, 100, 25, 60, 60, 100, 40, 10, 250, 100, 100, 8}; + Double_t min[16] = {kMinP, -0.5, 1.0, -1.0, -5.0, 0, 0, -3.0, -3.0, 0.0, 0.0, 0, 0, 0, 80, -1.5}; + Double_t max[16] = {kMaxP, 6.5, 4.0, 1.0, 4.0, 2.0, 0.05, 3.0, 3.0, 1.0, 20.0, 100, 100, 2.0, 180, 6.5}; + //fEleInfo = new THnSparseD("fEleInfo", "Electron Info; pT [GeV/c]; TPC signal;phi;eta;nSig; E/p;Rmatch;Ncell;clsF;M20;mcpT;Centrality;charge;opp;same;trigCond;MCele", 16, nBins, min, max); + fEleInfo = new THnSparseD("fEleInfo", "Electron Info; pT [GeV/c]; TPC signal;phi;eta;nSig; E/p;Rmatch;DCA_xy;DCA_z;M20;mcpT;Centrality;charge;opp;same;trigCond;MCele", 16, nBins, min, max); + //if(fqahist==1)fOutputList->Add(fEleInfo); // Make common binning Int_t nBinsEop[3] = { 10, 50, 100}; @@ -1254,12 +1410,18 @@ void AliAnalysisTaskHFECal::UserCreateOutputObjects() fPhoElecPtMCM20 = new THnSparseD("fPhoElecPtMCM20", "MC Pho-inclusive electron pt with M20",5,nBinspho2,minpho2,maxpho2); fOutputList->Add(fPhoElecPtMCM20); + fPhoElecPtMCM20Mass = new TH2D("fPhoElecPtMCM20Mass", "MC Pho-inclusive electron pt with M20 Mass",200,0,100,100,0,50); + fOutputList->Add(fPhoElecPtMCM20Mass); + fSameElecPtMC = new THnSparseD("fSameElecPtMC", "MC Same-inclusive electron pt",5,nBinspho2,minpho2,maxpho2); fOutputList->Add(fSameElecPtMC); fSameElecPtMCM20 = new THnSparseD("fSameElecPtMCM20", "MC Same-inclusive electron pt with M20",5,nBinspho2,minpho2,maxpho2); fOutputList->Add(fSameElecPtMCM20); + fSameElecPtMCM20Mass = new TH2D("fSameElecPtMCM20Mass", "MC Same-inclusive electron pt with M20 Mass",200,0,100,100,0,50); + fOutputList->Add(fSameElecPtMCM20Mass); + fIncpTMCM20pho_pi0e = new THnSparseD("fIncpTMCM20pho_pi0e","MC Pho pi0->e pid electro vs. centrality with M20",5,nBinspho2,minpho2,maxpho2); fIncpTMCM20pho_pi0e->Sumw2(); fOutputList->Add(fIncpTMCM20pho_pi0e); @@ -1295,6 +1457,20 @@ void AliAnalysisTaskHFECal::UserCreateOutputObjects() fSameElecPtMC_pi0e_TPC = new THnSparseD("fSameElecPtMC_pi0e_TPC", "MC Same-inclusive electron pt pi0->e",5,nBinspho2,minpho2,maxpho2); fSameElecPtMC_pi0e_TPC->Sumw2(); fOutputList->Add(fSameElecPtMC_pi0e_TPC); + + fIncpTMCpho_eta_TPC = new THnSparseD("fIncpTMCpho_eta_TPC","MC Pho pi0->e pid electro vs. centrality with M20",5,nBinspho2,minpho2,maxpho2); + fIncpTMCpho_eta_TPC->Sumw2(); + fOutputList->Add(fIncpTMCpho_eta_TPC); + + fPhoElecPtMC_eta_TPC = new THnSparseD("fPhoElecPtMC_eta_TPC", "MC Pho-inclusive electron pt with pi0->e",5,nBinspho2,minpho2,maxpho2); + fPhoElecPtMC_eta_TPC->Sumw2(); + fOutputList->Add(fPhoElecPtMC_eta_TPC); + + fSameElecPtMC_eta_TPC = new THnSparseD("fSameElecPtMC_eta_TPC", "MC Same-inclusive electron pt pi0->e",5,nBinspho2,minpho2,maxpho2); + fSameElecPtMC_eta_TPC->Sumw2(); + fOutputList->Add(fSameElecPtMC_eta_TPC); + + //------------- CheckNclust = new TH1D("CheckNclust","cluster check",200,0,200); @@ -1302,6 +1478,9 @@ void AliAnalysisTaskHFECal::UserCreateOutputObjects() CheckNits = new TH1D("CheckNits","ITS cluster check",8,-0.5,7.5); fOutputList->Add(CheckNits); + + CheckDCA = new TH2D("CheckDCA","DCA check",200,-5,5,200,-5,5); + fOutputList->Add(CheckDCA); /* Hpi0pTcheck = new TH2D("Hpi0pTcheck","Pi0 pT from Hijing",100,0,50,3,-0.5,2.5); fOutputList->Add(Hpi0pTcheck); @@ -1372,16 +1551,94 @@ void AliAnalysisTaskHFECal::UserCreateOutputObjects() fnSigEtaCorr[5] = new TGraphErrors(12,etaval,corr5); // 50-70 fnSigEtaCorr[6] = new TGraphErrors(12,etaval,corr6); // 70-90 + fIncMaxE = new TH2D("fIncMaxE","Inc",10,0,100,10,0,100); + fOutputList->Add(fIncMaxE); - fIncReco = new TH2D("fIncReco","Inc",10,0,100,100,0,100); + fIncReco = new TH2D("fIncReco","Inc",10,0,100,100,0,500); fOutputList->Add(fIncReco); - fPhoReco = new TH2D("fPhoReco","Pho",10,0,100,100,0,100); + fPhoReco = new TH2D("fPhoReco","Pho",10,0,100,100,0,500); fOutputList->Add(fPhoReco); - fSamReco = new TH2D("fSamReco","Same",10,0,100,100,0,100); + fSamReco = new TH2D("fSamReco","Same",10,0,100,100,0,500); fOutputList->Add(fSamReco); + fIncRecoMaxE = new TH2D("fIncRecoMaxE","Inc",10,0,100,100,0,500); + fOutputList->Add(fIncRecoMaxE); + + fPhoRecoMaxE = new TH2D("fPhoRecoMaxE","Pho",10,0,100,100,0,500); + fOutputList->Add(fPhoRecoMaxE); + + fSamRecoMaxE = new TH2D("fSamRecoMaxE","Same",10,0,100,100,0,500); + fOutputList->Add(fSamRecoMaxE); + + fPhoVertexReco_TPC = new TH2D("fPhoVertexReco_TPC","photon production Vertex mass selection TPC",40,0,20,250,0,50); + fPhoVertexReco_TPC->Sumw2(); + fOutputList->Add(fPhoVertexReco_TPC); + + fPhoVertexReco_TPC_Invmass = new TH2D("fPhoVertexReco_TPC_Invmass","photon production Vertex mass selection TPC Invmass",40,0,20,250,0,50); + fPhoVertexReco_TPC_Invmass->Sumw2(); + fOutputList->Add(fPhoVertexReco_TPC_Invmass); + + fPhoVertexReco_EMCal = new TH2D("fPhoVertexReco_EMCal","photon production Vertex mass selection",40,0,20,250,0,50); + fPhoVertexReco_EMCal->Sumw2(); + fOutputList->Add(fPhoVertexReco_EMCal); + + fPhoVertexReco_Invmass = new TH2D("fPhoVertexReco_Invmass","photon production Vertex mass selection",40,0,20,250,0,50); + fPhoVertexReco_Invmass->Sumw2(); + fOutputList->Add(fPhoVertexReco_Invmass); + + fPhoVertexReco_step0= new TH2D("fPhoVertexReco_step0","photon production Vertex mass selection",40,0,20,250,0,50); + fPhoVertexReco_step0->Sumw2(); + fOutputList->Add(fPhoVertexReco_step0); + + fPhoVertexReco_step1= new TH2D("fPhoVertexReco_step1","photon production Vertex mass selection",40,0,20,250,0,50); + fPhoVertexReco_step1->Sumw2(); + fOutputList->Add(fPhoVertexReco_step1); + + fMatchV0_0 = new TH1D("fMatchV0_0","V0 match",100,0,20); + fOutputList->Add(fMatchV0_0); + + fMatchV0_1 = new TH1D("fMatchV0_1","V0 match",100,0,20); + fOutputList->Add(fMatchV0_1); + + fMatchMC_0 = new TH1D("fMatchMC_0","MC match",100,0,20); + fOutputList->Add(fMatchMC_0); + + fMatchMC_1 = new TH1D("fMatchMC_1","MC match",100,0,20); + fOutputList->Add(fMatchMC_1); + + fpair = new TH2D("fpair","pair of associate",100,0,20,21,-10.5,10.5); + fOutputList->Add(fpair); + + fFakeRejection0 = new TH2D("fFakeRejection0","TPC PID",2,-0.5,1.5,100,0,20); + fFakeRejection0->Sumw2(); + fOutputList->Add(fFakeRejection0); + + fFakeRejection1 = new TH2D("fFakeRejection1","TPC PID + Tr match",2,-0.5,1.5,100,0,20); + fFakeRejection1->Sumw2(); + fOutputList->Add(fFakeRejection1); + + fFakeRejection2 = new TH2D("fFakeRejection2","TPC PID + Tr match + E/p",2,-0.5,1.5,100,0,20); + fFakeRejection2->Sumw2(); + fOutputList->Add(fFakeRejection2); + + EopFake = new TH2D("EopFake","negative track Eop",20,0,20,200,0,2); + EopFake->Sumw2(); + fOutputList->Add(EopFake); + + EopTrue = new TH2D("EopTrue","true track Eop",20,0,20,200,0,2); + EopTrue->Sumw2(); + fOutputList->Add(EopTrue); + + MatchFake = new TH2D("MatchFake","negative track Match",20,0,20,100,0,0.05); + MatchFake->Sumw2(); + fOutputList->Add(MatchFake); + + MatchTrue = new TH2D("MatchTrue","true track Match",20,0,20,100,0.05); + MatchTrue->Sumw2(); + fOutputList->Add(MatchTrue); + PostData(1,fOutputList); } @@ -1401,9 +1658,7 @@ Bool_t AliAnalysisTaskHFECal::ProcessCutStep(Int_t cutStep, AliVParticle *track) return kTRUE; } //_________________________________________ -//void AliAnalysisTaskHFECal::SelectPhotonicElectron(Int_t itrack, Double_t cent, AliESDtrack *track, Bool_t &fFlagPhotonicElec) -//void AliAnalysisTaskHFECal::SelectPhotonicElectron(Int_t itrack, Double_t cent, AliESDtrack *track, Bool_t &fFlagPhotonicElec, Bool_t &fFlagConvinatElec, Double_t nSig) -void AliAnalysisTaskHFECal::SelectPhotonicElectron(Int_t itrack, Double_t cent, AliESDtrack *track, Bool_t &fFlagPhotonicElec, Bool_t &fFlagConvinatElec, Double_t nSig, Double_t shower, Double_t ep, Double_t mce, Double_t w, Int_t ibgevent, Bool_t tagpi0, Bool_t tageta) +void AliAnalysisTaskHFECal::SelectPhotonicElectron(Int_t itrack, Double_t cent, AliESDtrack *track, Bool_t &fFlagPhotonic, Bool_t &fFlagConvinat, Double_t nSig, Double_t shower, Double_t ep, Double_t mce, Double_t w, Int_t ibgevent, Bool_t tagpi0, Bool_t tageta, Int_t iCal) { //Identify non-heavy flavour electrons using Invariant mass method @@ -1411,12 +1666,13 @@ void AliAnalysisTaskHFECal::SelectPhotonicElectron(Int_t itrack, Double_t cent, fTrackCuts->SetRequireTPCRefit(kTRUE); fTrackCuts->SetRequireITSRefit(kTRUE); fTrackCuts->SetEtaRange(-0.9,0.9); - //fTrackCuts->SetRequireSigmaToVertex(kTRUE); fTrackCuts->SetMaxChi2PerClusterTPC(3.5); fTrackCuts->SetMinNClustersTPC(90); - const AliESDVertex *pVtx = fESD->GetPrimaryVertex(); - + //const AliESDVertex *pVtx = fESD->GetPrimaryVertex(); + Double_t bfield = fESD->GetMagneticField(); + Double_t emass = 0.51*0.001; // (0.51 MeV) + double ptEle = track->Pt(); //add if(ibgevent==0 && w > 0.0) { @@ -1434,7 +1690,9 @@ void AliAnalysisTaskHFECal::SelectPhotonicElectron(Int_t itrack, Double_t cent, p1 = particle->GetFirstMother(); } - //for(Int_t jTracks = itrack+1; jTracksGetNumberOfTracks(); jTracks++){ + int numULS = 0; + int numLS = 0; + for(Int_t jTracks = 0; jTracksGetNumberOfTracks(); jTracks++){ AliESDtrack* trackAsso = fESD->GetTrack(jTracks); if (!trackAsso) { @@ -1468,9 +1726,13 @@ void AliAnalysisTaskHFECal::SelectPhotonicElectron(Int_t itrack, Double_t cent, Int_t charge = track->Charge(); - if(ptAsso <0.5) continue; + if(ptAsso AcceptTrack(trackAsso)) continue; - if(dEdxAsso <65 || dEdxAsso>100) continue; //11a pass1 + double fTPCnSigmaAss = fPID->GetPIDResponse() ? fPID->GetPIDResponse()->NumberOfSigmasTPC(trackAsso, AliPID::kElectron) : 1000; + //cout << "fTPCnSigmaAss = " << fTPCnSigmaAss << endl; + //cout << "fTPCnSigmaAss Cut = " << fMimNsigassCut << endl; + if(fTPCnSigmaAss 5) continue; + //cout << "fTPCnSigmaAss a.f. cut = " << fTPCnSigmaAss << endl; Int_t fPDGe1 = 11; Int_t fPDGe2 = 11; if(charge>0) fPDGe1 = -11; @@ -1483,17 +1745,29 @@ void AliAnalysisTaskHFECal::SelectPhotonicElectron(Int_t itrack, Double_t cent, //printf("fFlagLS = %d\n",fFlagLS); //printf("fFlagULS = %d\n",fFlagULS); - printf("\n"); + //printf("\n"); + AliKFParticle::SetField(bfield); AliKFParticle ge1(*track, fPDGe1); AliKFParticle ge2(*trackAsso, fPDGe2); + + if(fSetMassNonlinear) + { + ge1.SetNonlinearMassConstraint(emass); + ge2.SetNonlinearMassConstraint(emass); + } + AliKFParticle recg(ge1, ge2); + /* // vertex AliKFVertex primV(*pVtx); primV += recg; + primV -= ge1; + primV -= ge2; recg.SetProductionVertex(primV); - + */ + // check chi2 if(recg.GetNDF()<1) continue; Double_t chi2recg = recg.GetChi2()/recg.GetNDF(); @@ -1507,35 +1781,39 @@ void AliAnalysisTaskHFECal::SelectPhotonicElectron(Int_t itrack, Double_t cent, } recg.GetMass(mass,width); + if(fSetMassWidthCut && width>1e10)continue; + // angle openingAngle = ge1.GetAngle(ge2); if(fFlagLS) fOpeningAngleLS->Fill(openingAngle); if(fFlagULS) fOpeningAngleULS->Fill(openingAngle); - + double ishower = 0; if(shower>0.0 && shower<0.3)ishower = 1; - double phoinfo[9]; - phoinfo[0] = cent; - phoinfo[1] = ptPrim; - phoinfo[2] = mass; - //phoinfo[3] = nSig; - phoinfo[3] = dEdxAsso; - phoinfo[4] = openingAngle; - phoinfo[5] = ishower; - phoinfo[6] = ep; - phoinfo[7] = mce; - phoinfo[8] = ptAsso; - - if(fFlagLS) fInvmassLS->Fill(phoinfo); - if(fFlagULS) fInvmassULS->Fill(phoinfo); - if(fFlagLS && ibgevent==0 && jbgevent==0) fInvmassLSmc->Fill(phoinfo,w); - if(fFlagULS && ibgevent==0 && jbgevent==0) - { - fInvmassULSmc->Fill(phoinfo,w); - } - //printf("fInvmassCut %f\n",fInvmassCut); - //printf("openingAngle %f\n",fOpeningAngleCut); + if(iCal==1) + { + double phoinfo[9]; + phoinfo[0] = cent; + phoinfo[1] = ptPrim; + phoinfo[2] = mass; + phoinfo[3] = nSig; + phoinfo[4] = openingAngle; + phoinfo[5] = ishower; + phoinfo[6] = ep; + phoinfo[7] = mce; + phoinfo[8] = ptAsso; + + if(fFlagLS) fInvmassLS->Fill(phoinfo); + if(fFlagULS) fInvmassULS->Fill(phoinfo); + if(fFlagLS && ibgevent==0 && jbgevent==0) fInvmassLSmc->Fill(phoinfo,w); + if(fFlagULS && ibgevent==0 && jbgevent==0) + { + fInvmassULSmc->Fill(phoinfo,w); + } + //printf("fInvmassCut %f\n",fInvmassCut); + //printf("openingAngle %f\n",fOpeningAngleCut); + } // angle cut if(openingAngle > fOpeningAngleCut) continue; @@ -1546,50 +1824,76 @@ void AliAnalysisTaskHFECal::SelectPhotonicElectron(Int_t itrack, Double_t cent, if(fFlagLS ) fInvmassLSreco->Fill(ptPrim,mass); if(fFlagULS) fInvmassULSreco->Fill(ptPrim,mass); + // check double count + if(mass0.0 && shower<0.3 && iCal==1)fPhoElecPtM20Mass->Fill(cent,ptPrim); + if(mass0.0 && shower<0.3 && iCal==1)fSameElecPtM20Mass->Fill(cent,ptPrim); + + // flag + if(mass0.0) - { - //cout << "tagpi0 = " << tagpi0 << " ; tageta = " << tageta << endl; - if(fFlagLS && ibgevent==0 && jbgevent==0 && tagpi0) fInvmassLSmc0->Fill(ptPrim,mass); - if(fFlagULS && ibgevent==0 && jbgevent==0 && tagpi0) fInvmassULSmc0->Fill(ptPrim,mass); - if(fFlagLS && ibgevent==0 && jbgevent==0 && tageta) fInvmassLSmc1->Fill(ptPrim,mass); - if(fFlagULS && ibgevent==0 && jbgevent==0 && tageta) fInvmassULSmc1->Fill(ptPrim,mass); - if(fFlagLS && ibgevent==0 && jbgevent==0 && (p1==p2) && tagpi0) fInvmassLSmc2->Fill(ptPrim,mass); - if(fFlagULS && ibgevent==0 && jbgevent==0 && (p1==p2) && tagpi0) fInvmassULSmc2->Fill(ptPrim,mass); - if(fFlagLS && ibgevent==0 && jbgevent==0 && (p1==p2) && tageta) fInvmassLSmc3->Fill(ptPrim,mass); - if(fFlagULS && ibgevent==0 && jbgevent==0 && (p1==p2) && tageta) fInvmassULSmc3->Fill(ptPrim,mass); - } + { + // count + if(mass0.0 && shower<0.3 && mce>2.1 && iCal==1)fPhoElecPtMCM20Mass->Fill(cent,ptPrim,w); + if(mass0.0 && shower<0.3 && mce>2.1 && iCal==1)fSameElecPtMCM20Mass->Fill(cent,ptPrim,w); + + //printf("MC data\n"); + + if(w>0.0) + { + //cout << "tagpi0 = " << tagpi0 << " ; tageta = " << tageta << endl; + if(iCal==1) + { + if(fFlagLS && ibgevent==0 && jbgevent==0 && tagpi0) fInvmassLSmc0->Fill(ptPrim,mass); + if(fFlagULS && ibgevent==0 && jbgevent==0 && tagpi0) fInvmassULSmc0->Fill(ptPrim,mass); + if(fFlagLS && ibgevent==0 && jbgevent==0 && tageta) fInvmassLSmc1->Fill(ptPrim,mass); + if(fFlagULS && ibgevent==0 && jbgevent==0 && tageta) fInvmassULSmc1->Fill(ptPrim,mass); + if(fFlagLS && ibgevent==0 && jbgevent==0 && (p1==p2) && tagpi0) fInvmassLSmc2->Fill(ptPrim,mass); + if(fFlagULS && ibgevent==0 && jbgevent==0 && (p1==p2) && tagpi0) fInvmassULSmc2->Fill(ptPrim,mass); + if(fFlagLS && ibgevent==0 && jbgevent==0 && (p1==p2) && tageta) fInvmassLSmc3->Fill(ptPrim,mass); + if(fFlagULS && ibgevent==0 && jbgevent==0 && (p1==p2) && tageta) fInvmassULSmc3->Fill(ptPrim,mass); + } + } + + if(mass0 || numLS>0) + { + int numPair = numULS-numLS; + if(iCal==1)fpair->Fill(ptEle,numPair); + } + + fFlagPhotonic = flagPhotonicElec; + fFlagConvinat = flagConvinatElec; } + //------------------------------------------- void AliAnalysisTaskHFECal::FindMother(TParticle* part, int &label, int &pid) @@ -1933,6 +2237,9 @@ double AliAnalysisTaskHFECal::MCEopMeanCorrection(double pTmc, float central) shift = fcorr0->Eval(pTmc)-fcorr1->Eval(pTmc); + fcorr0->Delete(); + fcorr1->Delete(); + return shift; } @@ -1982,41 +2289,3 @@ double AliAnalysisTaskHFECal::NsigmaCorrection(double tmpeta, float central) } -double AliAnalysisTaskHFECal::SumpT(Int_t itrack, AliESDtrack* track) -{ - - fTrackCuts->SetAcceptKinkDaughters(kFALSE); - fTrackCuts->SetRequireTPCRefit(kTRUE); - fTrackCuts->SetRequireITSRefit(kTRUE); - fTrackCuts->SetEtaRange(-0.9,0.9); - //fTrackCuts->SetRequireSigmaToVertex(kTRUE); - fTrackCuts->SetMaxChi2PerClusterTPC(3.5); - fTrackCuts->SetMinNClustersTPC(90); - - double pTrecp = track->Pt(); - double phiorg = track->Phi(); - double etaorg = track->Eta(); - - for(Int_t jTracks = 0; jTracksGetNumberOfTracks(); jTracks++){ - AliESDtrack* trackAsso = fESD->GetTrack(jTracks); - if (!trackAsso) { - printf("ERROR: Could not receive track %d\n", jTracks); - continue; - } - if(itrack==jTracks)continue; - double pTAss = trackAsso->Pt(); - double etaAss = trackAsso->Eta(); - double phiAss = trackAsso->Phi(); - - double delphi = phiorg - phiAss; - double deleta = etaorg - etaAss; - - double R = sqrt(pow(deleta,2)+pow(delphi,2)); - if(pTAss<0.5)continue; - if(R<0.4)pTrecp+=pTAss; - - } - - return pTrecp; -} -