1 /**************************************************************************
2 * Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
18 //*************************************************************************
19 // Class AliAnalysisTaskSEDvsMultiplicity
20 // AliAnalysisTaskSE for the D meson vs. multiplcity analysis
21 // Authors: Renu Bala, Zaida Conesa del Valle, Francesco Prino
22 /////////////////////////////////////////////////////////////
24 #include <TClonesArray.h>
28 #include <TDatabasePDG.h>
32 #include <THnSparse.h>
34 #include "AliAnalysisManager.h"
35 #include "AliRDHFCuts.h"
36 #include "AliRDHFCutsDplustoKpipi.h"
37 #include "AliRDHFCutsDStartoKpipi.h"
38 #include "AliRDHFCutsD0toKpi.h"
39 #include "AliAODHandler.h"
40 #include "AliAODEvent.h"
41 #include "AliAODVertex.h"
42 #include "AliAODTrack.h"
43 #include "AliAODRecoDecayHF.h"
44 #include "AliAODRecoCascadeHF.h"
45 #include "AliAnalysisVertexingHF.h"
46 #include "AliAnalysisTaskSE.h"
47 #include "AliAnalysisTaskSEDvsMultiplicity.h"
48 #include "AliNormalizationCounter.h"
49 #include "AliVertexingHFUtils.h"
50 ClassImp(AliAnalysisTaskSEDvsMultiplicity)
53 //________________________________________________________________________
54 AliAnalysisTaskSEDvsMultiplicity::AliAnalysisTaskSEDvsMultiplicity():
60 fHistNtrEta16vsNtrEta1(0),
61 fHistNtrCorrEta1vsNtrRawEta1(0),
63 fHistNtrCorrVsZvtx(0),
65 fHistNtrCorrEvWithCand(0),
66 fHistNtrCorrEvWithD(0),
68 fPtVsMassVsMultNoPid(0),
69 fPtVsMassVsMultUncorr(0),
70 fPtVsMassVsMultPart(0),
71 fPtVsMassVsMultAntiPart(0),
88 // Default constructor
89 for(Int_t i=0; i<5; i++) fHistMassPtImpPar[i]=0;
90 for(Int_t i=0; i<4; i++) fMultEstimatorAvg[i]=0;
93 //________________________________________________________________________
94 AliAnalysisTaskSEDvsMultiplicity::AliAnalysisTaskSEDvsMultiplicity(const char *name, Int_t pdgMeson,AliRDHFCuts *cuts):
95 AliAnalysisTaskSE(name),
100 fHistNtrEta16vsNtrEta1(0),
101 fHistNtrCorrEta1vsNtrRawEta1(0),
103 fHistNtrCorrVsZvtx(0),
104 fHistNtrCorrEvSel(0),
105 fHistNtrCorrEvWithCand(0),
106 fHistNtrCorrEvWithD(0),
108 fPtVsMassVsMultNoPid(0),
109 fPtVsMassVsMultUncorr(0),
110 fPtVsMassVsMultPart(0),
111 fPtVsMassVsMultAntiPart(0),
113 fLowmasslimit(1.765),
115 fRDCutsAnalysis(cuts),
120 fLowerImpPar(-2000.),
121 fHigherImpPar(2000.),
129 // Standard constructor
131 for(Int_t i=0; i<5; i++) fHistMassPtImpPar[i]=0;
132 for(Int_t i=0; i<4; i++) fMultEstimatorAvg[i]=0;
134 fNMassBins=200; // FIXME
135 SetMassLimits(0.,0.2); // FIXME
138 SetMassLimits(fPdgMeson,0.1);
140 // Default constructor
141 // Output slot #1 writes into a TList container
142 DefineOutput(1,TList::Class()); //My private output
143 // Output slot #2 writes cut to private output
144 DefineOutput(2,TList::Class());
145 // Output slot #3 writes cut to private output
146 DefineOutput(3,TList::Class());
148 //________________________________________________________________________
149 AliAnalysisTaskSEDvsMultiplicity::~AliAnalysisTaskSEDvsMultiplicity()
157 delete fRDCutsAnalysis;
160 for(Int_t i=0; i<4; i++) delete fMultEstimatorAvg[i];
161 for(Int_t i=0; i<5; i++){
162 delete fHistMassPtImpPar[i];
166 //_________________________________________________________________
167 void AliAnalysisTaskSEDvsMultiplicity::SetMassLimits(Double_t lowlimit, Double_t uplimit){
168 // set invariant mass limits
169 if(uplimit>lowlimit){
170 fUpmasslimit = lowlimit;
171 fLowmasslimit = uplimit;
173 AliError("Wrong mass limits: upper value should be larger than lower one");
176 //_________________________________________________________________
177 void AliAnalysisTaskSEDvsMultiplicity::SetMassLimits(Int_t pdg, Double_t range){
178 // set invariant mass limits
179 Double_t mass=TDatabasePDG::Instance()->GetParticle(TMath::Abs(pdg))->Mass();
180 SetMassLimits(mass-range,mass+range);
182 //________________________________________________________________________
183 void AliAnalysisTaskSEDvsMultiplicity::Init(){
187 printf("AnalysisTaskSEDvsMultiplicity::Init() \n");
189 fListCuts=new TList();
192 AliRDHFCutsDplustoKpipi* copycut=new AliRDHFCutsDplustoKpipi(*(static_cast<AliRDHFCutsDplustoKpipi*>(fRDCutsAnalysis)));
193 copycut->SetName("AnalysisCutsDplus");
194 fListCuts->Add(copycut);
195 }else if(fPdgMeson==421){
196 AliRDHFCutsD0toKpi* copycut=new AliRDHFCutsD0toKpi(*(static_cast<AliRDHFCutsD0toKpi*>(fRDCutsAnalysis)));
197 copycut->SetName("AnalysisCutsDzero");
198 fListCuts->Add(copycut);
199 }else if(fPdgMeson==413){
200 AliRDHFCutsDStartoKpipi* copycut=new AliRDHFCutsDStartoKpipi(*(static_cast<AliRDHFCutsDStartoKpipi*>(fRDCutsAnalysis)));
201 copycut->SetName("AnalysisCutsDStar");
202 fListCuts->Add(copycut);
204 PostData(2,fListCuts);
209 //________________________________________________________________________
210 void AliAnalysisTaskSEDvsMultiplicity::UserCreateOutputObjects()
212 // Create the output container
214 if(fDebug > 1) printf("AnalysisTaskSEDvsMultiplicity::UserCreateOutputObjects() \n");
216 // Several histograms are more conveniently managed in a TList
217 fOutput = new TList();
219 fOutput->SetName("OutputHistos");
221 fHistNtrCorrEvSel = new TH1F("hNtrCorrEvSel","Tracklets multiplicity for selected events; Tracklets ; Entries",200,0.,200.);
222 fHistNtrCorrEvWithCand = new TH1F("hNtrCorrEvWithCand", "Tracklets multiplicity for events with D candidates; Tracklets ; Entries",200,0.,200.);// Total multiplicity
223 fHistNtrCorrEvWithD = new TH1F("hNtrCorrEvWithD", "Tracklets multiplicity for events with D in mass region ; Tracklets ; Entries",200,0.,200.); //
224 fHistNtrEta16vsNtrEta1 = new TH2F("hNtrEta16vsNtrEta1","Uncorrected Eta1.6 vs Eta1.0; Ntracklets #eta<1.0; Ntracklets #eta<1.6",200,-0.5,199.5,200,-0.5,199.5); //eta 1.6 vs eta 1.0 histogram
225 fHistNtrCorrEta1vsNtrRawEta1 = new TH2F("hNtrCorrEta1vsNtrRawEta1","Corrected Eta1 vs Eta1.0; Ntracklets #eta<1.0 corrected; Ntracklets #eta<1",200,-0.,200.,200,-0.5,199.5); //eta 1.6 vs eta 1.0 histogram
226 fHistNtrVsZvtx = new TH2F("hNtrVsZvtx","Ntracklet vs VtxZ; VtxZ;N_{tracklet};",300,-15,15,200,0,200.); //
227 fHistNtrCorrVsZvtx = new TH2F("hNtrCorrVsZvtx","Ntracklet vs VtxZ; VtxZ;N_{tracklet};",300,-15,15,200,0,200.); //
230 fHistNtrCorrEvSel->Sumw2();
231 fHistNtrCorrEvWithCand->Sumw2();
232 fHistNtrCorrEvWithD->Sumw2();
234 fOutput->Add(fHistNtrCorrEvSel);
235 fOutput->Add(fHistNtrCorrEvWithCand);
236 fOutput->Add(fHistNtrCorrEvWithD);
237 fOutput->Add(fHistNtrEta16vsNtrEta1);
238 fOutput->Add(fHistNtrCorrEta1vsNtrRawEta1);
239 fOutput->Add(fHistNtrVsZvtx);
240 fOutput->Add(fHistNtrCorrVsZvtx);
243 fHistNEvents = new TH1F("fHistNEvents", "number of events ",11,-0.5,10.5);
244 fHistNEvents->GetXaxis()->SetBinLabel(1,"nEvents total");
245 fHistNEvents->GetXaxis()->SetBinLabel(2,"nEvents with Z vertex");
246 fHistNEvents->GetXaxis()->SetBinLabel(3,"nEvents selected");
247 fHistNEvents->GetXaxis()->SetBinLabel(4,"Rejected due to trigger");
248 fHistNEvents->GetXaxis()->SetBinLabel(5,"Rejected due to phys sel");
249 fHistNEvents->GetXaxis()->SetBinLabel(6,"Rejected due to vertex cuts");
250 fHistNEvents->GetXaxis()->SetBinLabel(7,"Rejected due to pileup");
251 fHistNEvents->GetXaxis()->SetBinLabel(8,"Total no. of candidate");
252 fHistNEvents->GetXaxis()->SetBinLabel(9,"no. of cand wo bitmask");
253 fHistNEvents->GetXaxis()->SetBinLabel(10,"D after cuts (No PID)");
254 fHistNEvents->GetXaxis()->SetBinLabel(11,"D after cuts + PID)");
255 fHistNEvents->GetXaxis()->SetNdivisions(1,kFALSE);
256 fHistNEvents->Sumw2();
257 fHistNEvents->SetMinimum(0);
258 fOutput->Add(fHistNEvents);
260 fPtVsMassVsMult=new TH3F("hPtVsMassvsMult", "D candidates: p_{t} vs mass vs tracklets multiplicity; Tracklets; Mass M [GeV/c^{2}]; p_{t} [GeV/c]",200,0.,200.,fNMassBins,fLowmasslimit,fUpmasslimit,48,0.,24.);
262 fPtVsMassVsMultNoPid=new TH3F("hPtVsMassvsMultNoPid", "D candidates: p_{t} vs mass vs tracklets multiplicity; Tracklets; Mass M [GeV/c^{2}]; p_{t} [GeV/c]",200,0.,200.,fNMassBins,fLowmasslimit,fUpmasslimit,48,0.,24.);
264 fPtVsMassVsMultUncorr=new TH3F("hPtVsMassvsMultUncorr", "D candidates: p_{t} vs mass vs tracklets multiplicity; Tracklets; Mass M [GeV/c^{2}]; p_{t} [GeV/c]",200,-0.5,199.5,fNMassBins,fLowmasslimit,fUpmasslimit,48,0.,24.);
266 fPtVsMassVsMultPart=new TH3F("hPtVsMassvsMultPart", "D candidates: p_{t} vs mass vs tracklets multiplicity; Tracklets; Mass M [GeV/c^{2}]; p_{t} [GeV/c]",200,0.,200.,fNMassBins,fLowmasslimit,fUpmasslimit,48,0.,24.);
268 fPtVsMassVsMultAntiPart=new TH3F("hPtVsMassvsMultAntiPart", "D candidates: p_{t} vs mass vs tracklets multiplicity; Tracklets; Mass M [GeV/c^{2}]; p_{t} [GeV/c]",200,0.,200.,fNMassBins,fLowmasslimit,fUpmasslimit,48,0.,24.);
270 fOutput->Add(fPtVsMassVsMult);
271 fOutput->Add(fPtVsMassVsMultUncorr);
272 fOutput->Add(fPtVsMassVsMultNoPid);
273 fOutput->Add(fPtVsMassVsMultPart);
274 fOutput->Add(fPtVsMassVsMultAntiPart);
276 if(fDoImpPar) CreateImpactParameterHistos();
278 fCounter = new AliNormalizationCounter("NormCounterCorrMult");
279 fCounter->SetStudyMultiplicity(kTRUE,1.);
282 fCounterU = new AliNormalizationCounter("NormCounterUnCorrMult");
283 fCounterU->SetStudyMultiplicity(kTRUE,1.);
286 fOutputCounters = new TList();
287 fOutputCounters->SetOwner();
288 fOutputCounters->SetName("OutputCounters");
289 fOutputCounters->Add(fCounter);
290 fOutputCounters->Add(fCounterU);
293 PostData(2,fListCuts);
294 PostData(3,fOutputCounters);
301 //________________________________________________________________________
302 void AliAnalysisTaskSEDvsMultiplicity::UserExec(Option_t */*option*/)
304 // Execute analysis for current event:
305 // heavy flavor candidates association to MC truth
307 AliAODEvent *aod = dynamic_cast<AliAODEvent*> (InputEvent());
309 // AliAODTracklets* tracklets = aod->GetTracklets();
310 //Int_t ntracklets = tracklets->GetNumberOfTracklets();
313 TClonesArray *arrayCand = 0;
314 TString arrayName="";
319 arrayName="Charm3Prong";
320 pdgDau[0]=211; pdgDau[1]=321; pdgDau[2]=211;
322 selbit=AliRDHFCuts::kDplusCuts;
323 }else if(fPdgMeson==421){
325 pdgDau[0]=211; pdgDau[1]=321; pdgDau[2]=0;
327 selbit=AliRDHFCuts::kD0toKpiCuts;
328 }else if(fPdgMeson==413){
330 pdgDau[0]=321; pdgDau[1]=211; pdgDau[2]=211;
332 selbit=AliRDHFCuts::kDstarCuts;
335 if(!aod && AODEvent() && IsStandardAOD()) {
336 // In case there is an AOD handler writing a standard AOD, use the AOD
337 // event in memory rather than the input (ESD) event.
338 aod = dynamic_cast<AliAODEvent*> (AODEvent());
339 // in this case the braches in the deltaAOD (AliAOD.VertexingHF.root)
340 // have to taken from the AOD event hold by the AliAODExtension
341 AliAODHandler* aodHandler = (AliAODHandler*)
342 ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());
343 if(aodHandler->GetExtensions()) {
344 AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");
345 AliAODEvent *aodFromExt = ext->GetAOD();
346 arrayCand=(TClonesArray*)aodFromExt->GetList()->FindObject(arrayName.Data());
349 arrayCand=(TClonesArray*)aod->GetList()->FindObject(arrayName.Data());
352 if(!aod || !arrayCand) {
353 printf("AliAnalysisTaskSEDvsMultiplicity::UserExec: Charm3Prong branch not found!\n");
359 // fix for temporary bug in ESDfilter
360 // the AODs with null vertex pointer didn't pass the PhysSel
361 if(!aod->GetPrimaryVertex()||TMath::Abs(aod->GetMagneticField())<0.001) return;
363 Int_t countTreta1=AliVertexingHFUtils::GetNumberOfTrackletsInEtaRange(aod,-1.,1.);
364 Int_t countTr=AliVertexingHFUtils::GetNumberOfTrackletsInEtaRange(aod,-1.6,1.6);
366 fCounterU->StoreEvent(aod,fRDCutsAnalysis,fReadMC,countTreta1);
367 fHistNEvents->Fill(0); // count event
369 Double_t countTreta1corr=countTreta1;
370 AliAODVertex *vtx1 = (AliAODVertex*)aod->GetPrimaryVertex();
372 if(vtx1->GetNContributors()>0){
373 fHistNEvents->Fill(1);
374 TProfile* estimatorAvg = GetEstimatorHistogram(aod);
376 countTreta1corr=AliVertexingHFUtils::GetCorrectedNtracklets(estimatorAvg,countTreta1,vtx1->GetZ(),fRefMult);
381 fCounter->StoreEvent(aod,fRDCutsAnalysis,fReadMC,(Int_t)countTreta1corr);
383 Bool_t isEvSel=fRDCutsAnalysis->IsEventSelected(aod);
385 if(fRDCutsAnalysis->GetWhyRejection()==5) fHistNEvents->Fill(3);
386 if(fRDCutsAnalysis->GetWhyRejection()==7) fHistNEvents->Fill(4);
387 if(fRDCutsAnalysis->GetWhyRejection()==6) fHistNEvents->Fill(5);
388 if(fRDCutsAnalysis->GetWhyRejection()==1) fHistNEvents->Fill(6);
392 fHistNtrEta16vsNtrEta1->Fill(countTreta1,countTr);
393 fHistNtrCorrEta1vsNtrRawEta1->Fill(countTreta1,countTreta1corr);
395 fHistNtrVsZvtx->Fill(vtx1->GetZ(),countTreta1);
396 fHistNtrCorrVsZvtx->Fill(vtx1->GetZ(),countTreta1corr);
399 TClonesArray *arrayMC=0;
400 AliAODMCHeader *mcHeader=0;
405 arrayMC = (TClonesArray*)aod->GetList()->FindObject(AliAODMCParticle::StdBranchName());
407 printf("AliAnalysisTaskSEDvsMultiplicity::UserExec: MC particles branch not found!\n");
411 mcHeader = (AliAODMCHeader*)aod->GetList()->FindObject(AliAODMCHeader::StdBranchName());
413 printf("AliAnalysisTaskSEDvsMultiplicity::UserExec: MC header branch not found!\n");
418 Int_t nCand = arrayCand->GetEntriesFast();
419 Int_t nSelectedNoPID=0,nSelectedPID=0,nSelectedInMassPeak=0;
420 Double_t mD0PDG = TDatabasePDG::Instance()->GetParticle(421)->Mass();
421 Double_t mDplusPDG = TDatabasePDG::Instance()->GetParticle(411)->Mass();
422 Double_t mDstarPDG = TDatabasePDG::Instance()->GetParticle(413)->Mass();
424 for (Int_t iCand = 0; iCand < nCand; iCand++) {
425 AliAODRecoDecayHF *d = (AliAODRecoDecayHF*)arrayCand->UncheckedAt(iCand);
426 fHistNEvents->Fill(7);
427 if(fUseBit && !d->HasSelectionBit(selbit)){
428 fHistNEvents->Fill(8);
432 Double_t ptCand = d->Pt();
433 Double_t rapid=d->Y(fPdgMeson);
434 Bool_t isFidAcc=fRDCutsAnalysis->IsInFiducialAcceptance(ptCand,rapid);
435 Int_t passAllCuts=fRDCutsAnalysis->IsSelected(d,AliRDHFCuts::kCandidate,aod);
436 Int_t passTopolCuts=fRDCutsAnalysis->GetIsSelectedCuts();
437 if(passTopolCuts==0) continue;
439 fHistNEvents->Fill(9);
442 fHistNEvents->Fill(10);
444 Bool_t isPrimary=kTRUE;
446 Double_t trueImpParXY=9999.;
447 Double_t impparXY=d->ImpParXY()*10000.;
448 Double_t dlen=0.1; //FIXME
451 mass[0]=d->InvMass(nDau,pdgDau);
453 if(TMath::Abs(mass[0]-mDplusPDG)<0.02) nSelectedInMassPeak++; //20 MeV for now... FIXME
454 }else if(fPdgMeson==421){
455 UInt_t pdgdaughtersD0[2]={211,321};//pi,K
456 UInt_t pdgdaughtersD0bar[2]={321,211};//K,pi
457 mass[0]=d->InvMass(2,pdgdaughtersD0);
458 mass[1]=d->InvMass(2,pdgdaughtersD0bar);
459 if(TMath::Abs(mass[0]-mD0PDG)<0.02 || TMath::Abs(mass[1]-mD0PDG)<0.02 ) nSelectedInMassPeak++; //20 MeV for now... FIXME
460 }else if(fPdgMeson==413){
462 AliAODRecoCascadeHF* temp = (AliAODRecoCascadeHF*)d;
463 mass[0]=temp->DeltaInvMass();
465 if(TMath::Abs(mass[0]-(mDstarPDG-mD0PDG))<0.001) nSelectedInMassPeak++; //1 MeV for now... FIXME
467 for(Int_t iHyp=0; iHyp<2; iHyp++){
468 if(mass[iHyp]<0.) continue; // for D+ and D* we have 1 mass hypothesis
469 Double_t invMass=mass[iHyp];
470 Double_t arrayForSparse[5]={invMass,ptCand,impparXY,dlen,countTreta1corr};
474 labD = d->MatchToMC(fPdgMeson,arrayMC,nDau,(Int_t*)pdgDau);
475 Bool_t fillHisto=fDoImpPar;
477 AliAODMCParticle *partD = (AliAODMCParticle*)arrayMC->At(labD);
478 Int_t code=partD->GetPdgCode();
479 if(CheckOrigin(arrayMC,partD)==5) isPrimary=kFALSE;
480 if(code<0 && iHyp==0) fillHisto=kFALSE;
481 if(code>0 && iHyp==1) fillHisto=kFALSE;
484 trueImpParXY=AliVertexingHFUtils::GetTrueImpactParameterDplus(mcHeader,arrayMC,partD)*10000.;
485 }else if(fPdgMeson==421){
486 trueImpParXY=AliVertexingHFUtils::GetTrueImpactParameterDzero(mcHeader,arrayMC,partD)*10000.;
487 }else if(fPdgMeson==413){
488 trueImpParXY=0.; /// FIXME
490 Double_t arrayForSparseTrue[5]={invMass,ptCand,trueImpParXY,dlen,countTreta1corr};
491 if(fillHisto && isFidAcc && passAllCuts){
492 fHistMassPtImpPar[2]->Fill(arrayForSparse);
493 fHistMassPtImpPar[3]->Fill(arrayForSparseTrue);
496 if(fillHisto && isFidAcc && passAllCuts) fHistMassPtImpPar[1]->Fill(arrayForSparse);
499 if(fillHisto && isFidAcc && passAllCuts)fHistMassPtImpPar[4]->Fill(arrayForSparse);
502 if(TMath::Abs(labD)==fPdgMeson && fMCOption==2) continue;
503 if(TMath::Abs(labD)!=fPdgMeson && fMCOption==1) continue;
507 if(iHyp==0 && !(passTopolCuts&1)) continue; // candidate not passing as D0
508 if(iHyp==1 && !(passTopolCuts&2)) continue; // candidate not passing as D0bar
511 fPtVsMassVsMultNoPid->Fill(countTreta1corr,invMass,ptCand);
513 fPtVsMassVsMult->Fill(countTreta1corr,invMass,ptCand);
514 fPtVsMassVsMultUncorr->Fill(countTreta1,invMass,ptCand);
515 // Add separation between part antipart
517 if(d->GetCharge()>0) fPtVsMassVsMultPart->Fill(countTreta1corr,invMass,ptCand);
518 else fPtVsMassVsMultAntiPart->Fill(countTreta1corr,invMass,ptCand);
519 }else if(fPdgMeson==421){
520 if(passTopolCuts&1) fPtVsMassVsMultPart->Fill(countTreta1corr,invMass,ptCand);
521 if(passTopolCuts&2) fPtVsMassVsMultAntiPart->Fill(countTreta1corr,invMass,ptCand);
522 }else if(fPdgMeson==413){
523 // FIXME ADD Dstar!!!!!!!!
528 fHistMassPtImpPar[0]->Fill(arrayForSparse);
536 fCounter->StoreCandidates(aod,nSelectedNoPID,kTRUE);
537 fCounter->StoreCandidates(aod,nSelectedPID,kFALSE);
538 fHistNtrCorrEvSel->Fill(countTreta1corr);
539 if(nSelectedPID>0) fHistNtrCorrEvWithCand->Fill(countTreta1corr);
540 if(nSelectedInMassPeak>0) fHistNtrCorrEvWithD->Fill(countTreta1corr);
543 PostData(2,fListCuts);
549 //________________________________________________________________________
550 void AliAnalysisTaskSEDvsMultiplicity::CreateImpactParameterHistos(){
551 // Histos for impact paramter study
552 // mass . pt , impact parameter , decay length , multiplicity
554 Int_t nbins[5]={fNMassBins,200,fNImpParBins,50,100};
555 Double_t xmin[5]={fLowmasslimit,0.,fLowerImpPar,0.,0.};
556 Double_t xmax[5]={fUpmasslimit,20.,fHigherImpPar,1.,100.};
558 fHistMassPtImpPar[0]=new THnSparseF("hMassPtImpParAll",
559 "Mass vs. pt vs.imppar - All",
561 fHistMassPtImpPar[1]=new THnSparseF("hMassPtImpParPrompt",
562 "Mass vs. pt vs.imppar - promptD",
564 fHistMassPtImpPar[2]=new THnSparseF("hMassPtImpParBfeed",
565 "Mass vs. pt vs.imppar - DfromB",
567 fHistMassPtImpPar[3]=new THnSparseF("hMassPtImpParTrueBfeed",
568 "Mass vs. pt vs.true imppar -DfromB",
570 fHistMassPtImpPar[4]=new THnSparseF("hMassPtImpParBkg",
571 "Mass vs. pt vs.imppar - backgr.",
573 for(Int_t i=0; i<5;i++){
574 fOutput->Add(fHistMassPtImpPar[i]);
578 //________________________________________________________________________
579 void AliAnalysisTaskSEDvsMultiplicity::Terminate(Option_t */*option*/)
581 // Terminate analysis
583 if(fDebug > 1) printf("AnalysisTaskSEDvsMultiplicity: Terminate() \n");
585 fOutput = dynamic_cast<TList*> (GetOutputData(1));
587 printf("ERROR: fOutput not available\n");
590 fHistNEvents = dynamic_cast<TH1F*>(fOutput->FindObject("fHistNEvents"));
591 printf("Number of Analyzed Events = %d\n",(Int_t)fHistNEvents->GetBinContent(3));
596 //_________________________________________________________________________________________________
597 Int_t AliAnalysisTaskSEDvsMultiplicity::CheckOrigin(TClonesArray* arrayMC, AliAODMCParticle *mcPartCandidate) const {
599 // checking whether the mother of the particles come from a charm or a bottom quark
604 mother = mcPartCandidate->GetMother();
606 Int_t abspdgGranma =0;
607 Bool_t isFromB=kFALSE;
608 Bool_t isQuarkFound=kFALSE;
611 AliAODMCParticle* mcGranma = dynamic_cast<AliAODMCParticle*>(arrayMC->At(mother));
613 pdgGranma = mcGranma->GetPdgCode();
614 abspdgGranma = TMath::Abs(pdgGranma);
615 if ((abspdgGranma > 500 && abspdgGranma < 600) || (abspdgGranma > 5000 && abspdgGranma < 6000)){
618 if(abspdgGranma==4 || abspdgGranma==5) isQuarkFound=kTRUE;
619 mother = mcGranma->GetMother();
621 AliError("Failed casting the mother particle!");
626 if(isFromB) return 5;
632 //____________________________________________________________________________
633 TProfile* AliAnalysisTaskSEDvsMultiplicity::GetEstimatorHistogram(const AliVEvent* event){
634 // Get Estimator Histogram from period event->GetRunNumber();
636 // If you select SPD tracklets in |eta|<1 you should use type == 1
639 Int_t runNo = event->GetRunNumber();
640 Int_t period = -1; // 0-LHC10b, 1-LHC10c, 2-LHC10d, 3-LHC10e
641 if(runNo>114930 && runNo<117223) period = 0;
642 if(runNo>119158 && runNo<120830) period = 1;
643 if(runNo>122373 && runNo<126438) period = 2;
644 if(runNo>127711 && runNo<130841) period = 3;
645 if(period<0 || period>3) return 0;
647 return fMultEstimatorAvg[period];