fUseDmesonEfficiencyCorrection(kFALSE),
fUseCentrality(kFALSE),
fUseHadronicChannelAtKineLevel(kFALSE),
+fRemoveMoreThanOneDmesonCandidate(kFALSE),
fPhiBins(32),
fEvents(0),
fDebugLevel(0),
fUseDmesonEfficiencyCorrection(kFALSE),
fUseCentrality(kFALSE),
fUseHadronicChannelAtKineLevel(kFALSE),
+fRemoveMoreThanOneDmesonCandidate(kFALSE),
fPhiBins(32),
fEvents(0),
fDebugLevel(0),
Info("AliAnalysisTaskDStarCorrelations","Calling Constructor");
SetDim();
- if(fSystem == pp) fUseCentrality = kFALSE; else fUseCentrality = kTRUE;
+ if(fSystem == AA) fUseCentrality = kTRUE; else fUseCentrality = kFALSE;
// if(fSystem == AA) fBkgMethod = kDStarSB; else fBkgMethod = kDZeroSB;
// loop over D meson candidates
for (Int_t iDStartoD0pi = 0; iDStartoD0pi<looponDCands; iDStartoD0pi++) {
- //if(ncandidates) continue; // if there is more than one D candidate, skip it
+ if(fRemoveMoreThanOneDmesonCandidate && ncandidates) continue; // if there is more than one D candidate, skip it
// initialize variables
isInPeak = kFALSE;
// get D candidate variables
ptDStar = dstarD0pi->Pt();
+ if(ptDStar<fCuts->GetMinPtCandidate()||ptDStar>fCuts->GetMaxPtCandidate()) continue;
+
phiDStar = dstarD0pi->Phi();
etaDStar = dstarD0pi->Eta();
if(TMath::Abs(etaDStar) > fMaxEtaDStar) continue;
// cout << "crash here 1" << endl;
// plot D0 vs DStar mass
if(!fmixing){
- cout << Form("histDZerovsDStarMass_%d",ptbin) <<endl;
+ // cout << Form("histDZerovsDStarMass_%d",ptbin) <<endl;
((TH2D*)fDmesonOutput->FindObject(Form("histDZerovsDStarMass_%d",ptbin)))->Fill(invMassDZero,deltainvMDStar);
if(fUseDmesonEfficiencyCorrection) ((TH2D*)fDmesonOutput->FindObject(Form("histDZerovsDStarMassWeight_%d",ptbin)))->Fill(invMassDZero,deltainvMDStar,DmesonWeight);
}
} // end if !mixing
+ // Check for D* and D0 invariant candidates for phi and eta
+ if(fFullmode){
+ Double_t arrayDStar[5];
+ arrayDStar[0] = deltainvMDStar;
+ arrayDStar[1] = invMassDZero;
+ arrayDStar[2] = phiDStar;
+ arrayDStar[3] = etaDStar;
+ arrayDStar[4] = ptDStar;
+ if(!fmixing)((THnSparseF*)fDmesonOutput->FindObject("DmesonMassCheck"))->Fill(arrayDStar);
+ }
// good D0 canidates
phiDStar = DStarMC->Phi();
etaDStar = DStarMC->Eta();
+ phiDStar = fCorrelator->SetCorrectPhiRange(phiDStar);
+
if(TMath::Abs(etaDStar) > fMaxEtaDStar) continue;
// cout << "Dstars are selected" << endl;
}// end if pure MC information
-
-
+ // check kinematics for tagged D*
+ if(fmontecarlo && !fmixing) ((TH2F*)fOutputMC->FindObject("MCTagEtaInclusiveDStar"))->Fill(etaDStar,ptDStar); // fill phi, eta
+ if(fmontecarlo && !fmixing) ((TH2F*)fOutputMC->FindObject("MCTagPhiInclusiveDStar"))->Fill(phiDStar,ptDStar); // fill phi, eta
// getting the number of triggers in the MCtag D* case
- if(fmontecarlo && isDStarMCtag) ((TH1F*)fOutputMC->FindObject("MCtagPtDStar"))->Fill(ptDStar);
- if(fmontecarlo && isDStarMCtag && !isDfromB) ((TH1D*)fOutputMC->FindObject("MCtagPtDStarfromCharm"))->Fill(ptDStar);
- if(fmontecarlo && isDStarMCtag && isDfromB) ((TH1D*)fOutputMC->FindObject("MCtagPtDStarfromBeauty"))->Fill(ptDStar);
+ if(fmontecarlo && isDStarMCtag) ((TH1F*)fOutputMC->FindObject("MCtagPtDStar"))->Fill(ptDStar,DmesonWeight);
+ if(fmontecarlo && isDStarMCtag && !isDfromB) ((TH1D*)fOutputMC->FindObject("MCtagPtDStarfromCharm"))->Fill(ptDStar,DmesonWeight);
+ if(fmontecarlo && isDStarMCtag && isDfromB) ((TH1D*)fOutputMC->FindObject("MCtagPtDStarfromBeauty"))->Fill(ptDStar,DmesonWeight);
fCorrelator->SetTriggerParticleProperties(ptDStar,phiDStar,etaDStar); // pass to the object the necessary trigger part parameters
Int_t NofTracks = fCorrelator->GetNofTracks(); // number of tracks in event
+
+ if(!fmixing){
+ if(EventHasDStarCandidate) ((TH1D*)fTracksOutput->FindObject("TracksPerDcandidate"))->Fill(NofTracks);
+ if(fBkgMethod == kDStarSB && EventHasDStarSideBandCandidate) ((TH1D*)fTracksOutput->FindObject("TracksPerSBcandidate"))->Fill(NofTracks);
+ if(fBkgMethod == kDZeroSB && EventHasDZeroSideBandCandidate) ((TH1D*)fTracksOutput->FindObject("TracksPerSBcandidate"))->Fill(NofTracks);
+
+ if(isDStarMCtag && fmontecarlo) ((TH1D*)fTracksOutput->FindObject("TracksPerDMC"))->Fill(NofTracks);
+ }
+
//************************************************** LOOP ON TRACKS *****************************************************
// cout << "crash here 6" << endl;
for(Int_t iTrack = 0; iTrack<NofTracks; iTrack++){ // looping on track candidates
Int_t label = hadron->GetLabel();
- if(!fmixing && !fReco){ // for kinematic MC
+
+
+ if(fmontecarlo){
AliAODMCParticle *part = (AliAODMCParticle*)fmcArray->At(label);
if(!part) continue;
- if(IsDDaughter(DStarMC, part)) continue;
- cout << "Skipping DStar daugheter " << endl;
+ if (!part->IsPhysicalPrimary()) continue; // removing secondary tracks
+ if(!fmixing && !fReco){
+ if(IsDDaughter(DStarMC, part)) continue; // skipping D* daughter
+ }
}
+
+
// if it is ok, then do the rest
Double_t ptHad = hadron->Pt();
- if(fmontecarlo) trackOrigin = fAssocCuts->IsMCpartFromHF(label,fmcArray);
+ if(fmontecarlo) {trackOrigin = fAssocCuts->IsMCpartFromHF(label,fmcArray);
+ if(trackOrigin[4]) {cout << "Something is wrong with hadron in MC - skipping" << endl; continue;}
+ }
// cout << "crash correlation 3" << endl;
if(!isTrackForPeakFilled && !fmixing && EventHasDStarCandidate){
MCarraytofill[3] = ptHad;
MCarraytofill[4] = poolbin;
- if(trackOrigin[0]) MCarraytofill[5] = 1 ;
- else if(trackOrigin[1]) MCarraytofill[5] = 2 ;
- else MCarraytofill[5] = 0;
+ if(trackOrigin[0]) MCarraytofill[5] = 1 ; // is from charm
+ else if(trackOrigin[1]) MCarraytofill[5] = 2 ; // is from beauty
+ else MCarraytofill[5] = 0; // non HF track
}
// filling signal
// if(isInPeak){
- if(EventHasDStarCandidate){
+ if(isInPeak){
// cout << "Filling signal " << endl;
// if(!fReco && TMath::Abs(etaHad)>0.8) continue;
//cout ("CorrelationsDStarHadron_%d",ptbin)
while (mother > 0){
AliAODMCParticle* mcMoth = dynamic_cast<AliAODMCParticle*>(fmcArray->At(mother)); //it's the mother of the track!
if (mcMoth){
- if (mcMoth->GetLabel() == labelD0) isDaughter = kTRUE;
+ if (mcMoth->GetLabel() == labelD0) {isDaughter = kTRUE; return isDaughter;}
mother = mcMoth->GetMother(); //goes back by one
} else{
AliError("Failed casting the mother particle!");
} // end loop on bins
+ // make Dstar mass ThnSparse
+
+ Int_t NDstarSparse[5]= {30 , 15 , 32, 18,14};
+ Double_t binLowDstarSparse[5]={0.14,1.7 ,-0.5*TMath::Pi(), -0.9,2};
+ Double_t binUpDstarSparse[5]= {0.16,2.0 , 1.5*TMath::Pi(), 0.9,16};
+ THnSparseF * DmesonMassCheck = new THnSparseF("DmesonMassCheck","Check D meson mass; #Delta Inv Mass GeV/c^{2}; M(K#pi) GeV/c^{2}; #varphi; #eta; p_{T} GeV/c",5,NDstarSparse,binLowDstarSparse,binUpDstarSparse);
+ if(!fmixing) fDmesonOutput->Add(DmesonMassCheck);
}
fOutput->Add(EventPropsCheck);
//event properties
- TH1D * SPDMultiplicty = new TH1D("MultiplicityOnlyCheck","Properties of the event; SPD Multiplicity",1000,0,1000);
+ TH1D * SPDMultiplicty = new TH1D("MultiplicityOnlyCheck","Properties of the event; SPD Multiplicity",20000,0,20000);
fOutput->Add(SPDMultiplicty);
if(!fmixing) fDmesonOutput->Add(NumberOfSBCandidates);
// phi distribution
- TH2F * PhiInclusiveDStar = new TH2F("PhiInclusiveDStar","Azimuthal distributions of Inclusive Dmesons; #phi; pT;Entries",nbinscorr,lowcorrbin,upcorrbin,30,0,30);
- TH2F * PhiSidebandDStar = new TH2F("PhiSidebandDStar","Azimuthal distributions of Sideband Dmesons; #phi; pT;Entries",nbinscorr,lowcorrbin,upcorrbin,30,0,30);
+ TH2F * PhiInclusiveDStar = new TH2F("PhiInclusiveDStar","Azimuthal distributions of Inclusive Dmesons; #phi; pT;Entries",nbinscorr,lowcorrbin,upcorrbin,50,0,50);
+ TH2F * PhiSidebandDStar = new TH2F("PhiSidebandDStar","Azimuthal distributions of Sideband Dmesons; #phi; pT;Entries",nbinscorr,lowcorrbin,upcorrbin,50,0,50);
// eta distribution
- TH2F * EtaInclusiveDStar = new TH2F("EtaInclusiveDStar","Azimuthal distributions of Inclusive Dmesons; #eta; pT;Entries",20,-1,1,30,0,30);
- TH2F * EtaSidebandDStar = new TH2F("EtaSidebandDStar","Azimuthal distributions of Sideband Dmesons; #eta; pT;Entries",20,-1,1,30,0,30);
+ TH2F * EtaInclusiveDStar = new TH2F("EtaInclusiveDStar","Azimuthal distributions of Inclusive Dmesons; #eta; pT;Entries",20,-1,1,50,0,50);
+ TH2F * EtaSidebandDStar = new TH2F("EtaSidebandDStar","Azimuthal distributions of Sideband Dmesons; #eta; pT;Entries",20,-1,1,50,0,50);
if(!fmixing) fDmesonOutput->Add(PhiInclusiveDStar);
if(!fmixing) fDmesonOutput->Add(PhiSidebandDStar);
TH2F * EtaInclusiveTracks = new TH2F("EtaInclusiveTracks","Azimuthal distributions of tracks if Inclusive Dmesons; #eta; pT;Entries",20,-1,1,100,0,10);
TH2F * EtaSidebandTracks = new TH2F("EtaSidebandTracks","Azimuthal distributions of tracks if Sideband Dmesons; #eta; pT;Entries",20,-1,1,100,0,10);
- TH1D * TracksPerDcandidate = new TH1D("TracksPerDcandidate","Distribution of number of tracks per D meson candidate; N tracks; counts",200,-0.5,199.5);
- TH1D * TracksPerSBcandidate = new TH1D("TracksPerSBcandidate","Distribution of number of tracks per sideband candidate; N tracks; counts",200,-0.5,199.5);
+ TH1D * TracksPerDcandidate = new TH1D("TracksPerDcandidate","Distribution of number of tracks per D meson candidate; N tracks; counts",20000,-0.5,19999.5);
+ TH1D * TracksPerSBcandidate = new TH1D("TracksPerSBcandidate","Distribution of number of tracks per sideband candidate; N tracks; counts",20000,-0.5,19999.5);
+ TH1D * TracksPerDMC = new TH1D("TracksPerDMC","Distribution of number of tracks per tagged D ; N tracks; counts",20000,-0.5,19999.5);
if(!fmixing) fTracksOutput->Add(PhiInclusiveTracks);
if(!fmixing) fTracksOutput->Add(PhiSidebandTracks);
if(!fmixing) fTracksOutput->Add(EtaInclusiveTracks);
if(!fmixing) fTracksOutput->Add(EtaSidebandTracks);
+
if(!fmixing) fTracksOutput->Add(TracksPerDcandidate);
if(!fmixing) fTracksOutput->Add(TracksPerSBcandidate);
+ if(!fmixing && fmontecarlo) fTracksOutput->Add(TracksPerDMC);
// Montecarlo for D*
TH1F *MCtagPtDStar = new TH1F("MCtagPtDStar","RECO pt of MCtagged DStars side bands",50,0,50);
if(fmontecarlo) fOutputMC->Add(MCtagPtDStar);
+ TH2F * MCTagEtaInclusiveDStar = new TH2F("MCTagEtaInclusiveDStar","MC Tag eta distributions of Inclusive Dmesons; #eta; pT;Entries",20,-1,1,50,0,50);
+ if(fmontecarlo && !fmixing) fOutputMC->Add(MCTagEtaInclusiveDStar);
+
+ TH2F * MCTagPhiInclusiveDStar = new TH2F("MCTagPhiInclusiveDStar","MC Tag Azimuthal distributions of Inclusive Dmesons; #eta; pT;Entries",20,-0.5*TMath::Pi(),1.5*TMath::Pi(),50,0,50);
+ if(fmontecarlo && !fmixing) fOutputMC->Add(MCTagPhiInclusiveDStar);
+
+
+
// event mixing histograms
- if(pt < fAODTrackCuts[0]) return kFALSE;
- if(pt > fAODTrackCuts[1]) return kFALSE;
- if(d0 < fAODTrackCuts[2]) return kFALSE;
- if(d0 > fAODTrackCuts[3]) return kFALSE;
+ if(pt < fAODTrackCuts[0]) {return kFALSE; cout << "reject min pt" << endl;}
+ if(pt > fAODTrackCuts[1]) {return kFALSE; cout << "reject max pt" << endl;}
+ if(d0 < fAODTrackCuts[2]) {return kFALSE; cout << "reject min d0" << endl;}
+ if(d0 > fAODTrackCuts[3]) {return kFALSE; cout << "reject max d0" << endl;}
return kTRUE;
//--------------------------------------------------------------------------
Bool_t *AliHFAssociatedTrackCuts::IsMCpartFromHF(Int_t label, TClonesArray*mcArray){
// Check origin in MC
-
- AliAODMCParticle* mcParticle;
- Int_t pdgCode = -1;
-
+
Bool_t isCharmy = kFALSE;
Bool_t isBeauty = kFALSE;
Bool_t isD = kFALSE;
Bool_t isB = kFALSE;
- Bool_t *originvect = new Bool_t[4];
+ Bool_t *originvect = new Bool_t[5];
- originvect[0] = kFALSE;
- originvect[1] = kFALSE;
- originvect[2] = kFALSE;
- originvect[3] = kFALSE;
-
- if (label<0) return originvect;
-
- while(pdgCode!=2212){ // loops back to the collision to check the particle source
-
- mcParticle = dynamic_cast<AliAODMCParticle*>(mcArray->At(label));
- if(!mcParticle) {AliError("NO MC PARTICLE"); break;}
- pdgCode = TMath::Abs(mcParticle->GetPdgCode());
-
- label = mcParticle->GetMother();
-
-
- if((pdgCode>=400 && pdgCode <500) || (pdgCode>=4000 && pdgCode<5000 )) isD = kTRUE;
- if((pdgCode>=500 && pdgCode <600) || (pdgCode>=5000 && pdgCode<6000 )) {isD = kFALSE; isB = kTRUE;}
-
-
- if(pdgCode == 4) isCharmy = kTRUE;
- if(pdgCode == 5) {isBeauty = kTRUE; isCharmy = kFALSE;}
- if(label<0) break;
-
- }
+ originvect[0] = kFALSE; // is from charm
+ originvect[1] = kFALSE; // is from beauty
+ originvect[2] = kFALSE; // is from D
+ originvect[3] = kFALSE; // is from B
+ originvect[4] = kFALSE; // did something go wrong? (kTRUE yes, kFALSE no))
+
+ //__________________________________
-
+ if (label<0) {originvect[4] = kTRUE; return originvect;}
+
+ AliAODMCParticle* mcParticle = dynamic_cast<AliAODMCParticle*>(mcArray->At(label));
+ if(!mcParticle) {originvect[4] = kTRUE; return originvect;}
+ Int_t pdgCode = -1;
+ Int_t mother = mcParticle->GetMother();
+
+ if(mother <0) {originvect[4] = kTRUE; return originvect;}
+
+ while (mother >= 0){
+ mcParticle = dynamic_cast<AliAODMCParticle*>(mcArray->At(mother));
+ if(!mcParticle) {AliError("NO MC PARTICLE"); break;}
+ pdgCode = TMath::Abs(mcParticle->GetPdgCode());
+
+ mother = mcParticle->GetMother();
+
+ if((pdgCode>=400 && pdgCode <500) || (pdgCode>=4000 && pdgCode<5000 )) isD = kTRUE;
+ if((pdgCode>=500 && pdgCode <600) || (pdgCode>=5000 && pdgCode<6000 )) {isD = kFALSE; isB = kTRUE;}
+ if(pdgCode == 4) isCharmy = kTRUE;
+ if(pdgCode == 5) {isBeauty = kTRUE; isCharmy = kFALSE;}
+ if(mother<0) break;
+ }
+
originvect[0] = isCharmy;
originvect[1] = isBeauty;
originvect[2] = isD;
/* $Id$ */
-AliAnalysisTaskDStarCorrelations *AddTaskDStarCorrelations(AliAnalysisTaskDStarCorrelations::CollSyst syst,
- Bool_t theMCon, Bool_t mixing, Bool_t UseReco=kTRUE,Bool_t UseHadChannelinMC,Bool_t fullmode = kFALSE,Bool_t UseEffic=kFALSE,Bool_t UseDEffic = kFALSE, Bool_t useDStarSidebands = kFALSE,
- AliAnalysisTaskDStarCorrelations::DEffVariable var,
- Int_t trackselect =1, Int_t usedispl =0, Int_t nbins, Float_t DStarSigma, Float_t D0Sigma, Float_t D0SBSigmaLow, Float_t D0SBSigmaHigh, Float_t eta,
- TString DStarCutsFile, TString TrackCutsFile,
- TString suffix = "")
+AliAnalysisTaskDStarCorrelations *AddTaskDStarCorrelations(AliAnalysisTaskDStarCorrelations::CollSyst syst, // set collisional system (pp, pA, AA)
+ Bool_t theMCon, // flag for Data (kFALSE) or MC (kTRUE) analysis
+ Bool_t mixing, // flag for Single Event (kFALSE) or Mixed Event (kTRUE) analysis
+ Bool_t UseReco=kTRUE, // flag for Kine/pure MC (kFALSE) or Reconstruction (kTRUE) analysis - in data, kTRUE by default
+ Bool_t UseHadChannelinMC, // flag to use D*->Kpipi (kTRUE) or any decay (kFALSE) in MC kine
+ Bool_t fullmode = kFALSE, // flag to run in fast mode (kFALSE) or slow and detailed (kTRUE)
+ Bool_t UseEffic=kFALSE, // flag to use associated track eff (kTRUE = YES, kFALSE = NO)
+ Bool_t UseDEffic = kFALSE, // flag to use Dmeson track eff (kTRUE = YES, kFALSE = NO)
+ Bool_t useDStarSidebands = kFALSE, // flag to use sidebands from D0 (kFALSE) or sidebands from Dstar (kTRUE)
+ Bool_t useOnlyOneDStarPerEvent, // use only one D* per event (kTRUE) or all of them (kFALSE)
+ AliAnalysisTaskDStarCorrelations::DEffVariable var, // variable to use in D mesn efficiency correction besides pt
+ Int_t trackselect =1, // correlate with hadrons (1), kaons (2), kzeros (3)
+ Int_t usedispl =0, // don't use displacement (0), use absolute displacement (1) or relative (normalized by impac.par resol.) displacement (2)
+ Int_t nbins, //number of bins in correlations histo
+ Float_t DStarSigma, // number of sigmas in dstar selection (for pt shape, eta, phi distr. studies)
+ Float_t D0Sigma, // number of sigmas in dzero selection (for pt shape, eta, phi distr. studies)
+ Float_t D0SBSigmaLow, // number of sigmas in dzero sb selection (for pt shape, eta, phi distr. studies)
+ Float_t D0SBSigmaHigh, // number of sigmas in dzero sb selection (for pt shape, eta, phi distr. studies)
+ Float_t eta, // maximum D* eta
+ Float_t minDStarPt, // set minimum pt for Dstar
+ TString DStarCutsFile, // path of Dmeson cut object
+ TString TrackCutsFile, // path of associated cut object
+ TString suffix = "" // suffix for output
+ )
{
// RDHFDStartoKpipi->SetTriggerClass("");
// RDHFDStartoKpipi->SetTriggerMask(AliVEvent::kCentral);
-
+ if(minDStarPt>RDHFDStartoKpipi->GetMinPtCandidate()){
+ cout << "Changing minimum pT of DStar from " << RDHFDStartoKpipi->GetMinPtCandidate() << " to " << minDStarPt << endl;
+ RDHFDStartoKpipi->SetMinPtCandidate(minDStarPt);
+ }
// ******************************** OPENING THE ASSOCIATED TRACK CUTS ************************************
cout << "Getting associated track cut object from file \n" << TrackCutsFile.Data() << "\n " << endl;
return;
}
+ if(UseEffic && !corrCuts->IsTrackEffMap()){
+ cout << "You are trying to use the single track efficiency, but there is no map loaded in the cut object " << endl;
+ return;
+ }
+
+
+ if(UseDEffic && var == AliAnalysisTaskDStarCorrelations::kNone && (!corrCuts->IsTrigEffMap1D()&&(!corrCuts->IsTrigEffMap1DB()))){
+ cout << "You are trying to use the DStar efficiency vs pt only, but there is no map loaded in the cut object " << endl;
+ return;
+ }
+
+ if(UseDEffic && (var != AliAnalysisTaskDStarCorrelations::kNone) && (!corrCuts->IsTrigEffMap2D()&&(!corrCuts->IsTrigEffMap2DB()))){
+ cout << "You are trying to use the DStar efficiency vs pt only, but there is no map loaded in the cut object " << endl;
+ return;
+ }
+
// ******************************** SELECTING THE MC PROCESS ************************************
TString selectMCproc = "";
//task->SetDMesonSigmas(sigmas);
task->SetUseEfficiencyCorrection(UseEffic);
task->SetUseDmesonEfficiencyCorrection(UseDEffic);
+ task->SetUseRemoveMoreThanOneCDmesonCandidate(useOnlyOneDStarPerEvent);
task->SetEfficiencyVariable(var);
task->SetMaxDStarEta(eta);