AliPhysics  e34b7ac (e34b7ac)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliAnalysisTaskSEDStarSpectra.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
3  * *
4  * Author: The ALICE Off-line Project. *
5  * Contributors are mentioned in the code where appropriate. *
6  * *
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  * appeuear 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  **************************************************************************/
15 
16 /* $Id$ */
17 
18 //
19 //
20 // Base class for DStar Analysis
21 //
22 //
23 // The D* spectra study is done in pt bins:
24 // [0,0.5] [0.5,1] [1,2] [2,3] [3,4] [4,5] [5,6] [6,7] [7,8],
25 // [8,10],[10,12], [12,16], [16,20] and [20,24]
26 //
27 // Cuts arew centralized in AliRDHFCutsDStartoKpipi
28 // Side Band and like sign background are implemented in the macro
29 //
30 //-----------------------------------------------------------------------
31 //
32 // Author A.Grelli
33 // ERC-QGP Utrecht University - a.grelli@uu.nl,
34 // Author Y.Wang
35 // University of Heidelberg - yifei@physi.uni-heidelberg.de
36 // Author C.Ivan
37 // ERC-QGP Utrecht University - c.ivan@uu.nl,
38 //
39 //-----------------------------------------------------------------------
40 
41 #include <TSystem.h>
42 #include <TParticle.h>
43 #include <TH1I.h>
44 #include "TROOT.h"
45 #include <TDatabasePDG.h>
46 #include <AliAnalysisDataSlot.h>
47 #include <AliAnalysisDataContainer.h>
49 #include "AliStack.h"
50 #include "AliMCEvent.h"
51 #include "AliAnalysisManager.h"
52 #include "AliAODMCHeader.h"
53 #include "AliAODHandler.h"
54 #include "AliLog.h"
55 #include "AliAODVertex.h"
56 #include "AliAODRecoDecay.h"
57 #include "AliAODRecoDecayHF.h"
58 #include "AliAODRecoCascadeHF.h"
60 #include "AliAnalysisVertexingHF.h"
61 #include "AliESDtrack.h"
62 #include "AliAODMCParticle.h"
64 #include "AliAODEvent.h"
66 
70 
71 //__________________________________________________________________________
74  fEvents(0),
75  fAnalysis(0),
76  fD0Window(0),
77  fPeakWindow(0),
78  fUseMCInfo(kFALSE),
79  fDoSearch(kFALSE),
80  fOutput(0),
81  fOutputAll(0),
82  fOutputPID(0),
83  fNSigma(3),
84  fCuts(0),
85  fCEvents(0),
86  fTrueDiff2(0),
87  fDeltaMassD1(0),
88  fCounter(0),
89  fDoImpParDstar(kFALSE),
90  fNImpParBins(400),
91  fLowerImpPar(-2000.),
92  fHigherImpPar(2000.),
93  fDoDStarVsY(kFALSE)
94 {
95  //
97  //
98  for(Int_t i=0;i<5;i++) fHistMassPtImpParTCDs[i]=0;
99 }
100 //___________________________________________________________________________
102  AliAnalysisTaskSE(name),
103  fEvents(0),
104  fAnalysis(0),
105  fD0Window(0),
106  fPeakWindow(0),
107  fUseMCInfo(kFALSE),
108  fDoSearch(kFALSE),
109  fOutput(0),
110  fOutputAll(0),
111  fOutputPID(0),
112  fNSigma(3),
113  fCuts(0),
114  fCEvents(0),
115  fTrueDiff2(0),
116  fDeltaMassD1(0),
117  fCounter(0),
118  fDoImpParDstar(kFALSE),
119  fNImpParBins(400),
120  fLowerImpPar(-2000.),
121  fHigherImpPar(2000.),
122  fDoDStarVsY(kFALSE)
123 {
124  //
126  //
127  Info("AliAnalysisTaskSEDStarSpectra","Calling Constructor");
128 
129  fCuts=cuts;
130  for(Int_t i=0;i<5;i++) fHistMassPtImpParTCDs[i]=0;
131 
132  DefineOutput(1,TList::Class()); //conters
133  DefineOutput(2,TList::Class()); //All Entries output
134  DefineOutput(3,TList::Class()); //3sigma PID output
135  DefineOutput(4,AliRDHFCutsDStartoKpipi::Class()); //My private output
136  DefineOutput(5,AliNormalizationCounter::Class()); // normalization
137 }
138 
139 //___________________________________________________________________________
141  //
143  //
144  Info("~AliAnalysisTaskSEDStarSpectra","Calling Destructor");
145 
146  delete fOutput;
147  delete fOutputAll;
148  delete fOutputPID;
149  delete fCuts;
150  delete fCEvents;
151  delete fDeltaMassD1;
152  for(Int_t i=0; i<5; i++){
153  delete fHistMassPtImpParTCDs[i];
154  }
155 }
156 //_________________________________________________
158  //
160  //
161 
162  if(fDebug > 1) printf("AnalysisTaskSEDStarSpectra::Init() \n");
164  // Post the data
165  PostData(4,copyfCuts);
166 
167  return;
168 }
169 
170 //_________________________________________________
172 {
174  if (!fInputEvent) {
175  Error("UserExec","NO EVENT FOUND!");
176  return;
177  }
178 
179  fEvents++;
180 
181  AliAODEvent* aodEvent = dynamic_cast<AliAODEvent*>(fInputEvent);
182  TClonesArray *arrayDStartoD0pi=0;
183 
184  fCEvents->Fill(1);
185 
186  if(!aodEvent && AODEvent() && IsStandardAOD()) {
187  // In case there is an AOD handler writing a standard AOD, use the AOD
188  // event in memory rather than the input (ESD) event.
189  aodEvent = dynamic_cast<AliAODEvent*> (AODEvent());
190  // in this case the braches in the deltaAOD (AliAOD.VertexingHF.root)
191  // have to taken from the AOD event hold by the AliAODExtension
192  AliAODHandler* aodHandler = (AliAODHandler*)
193  ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());
194  if(aodHandler->GetExtensions()) {
195  AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");
196  AliAODEvent *aodFromExt = ext->GetAOD();
197  arrayDStartoD0pi=(TClonesArray*)aodFromExt->GetList()->FindObject("Dstar");
198  }
199  } else {
200  arrayDStartoD0pi=(TClonesArray*)aodEvent->GetList()->FindObject("Dstar");
201  }
202 
203  // fix for temporary bug in ESDfilter
204  // the AODs with null vertex pointer didn't pass the PhysSel
205  if(!aodEvent->GetPrimaryVertex() || TMath::Abs(aodEvent->GetMagneticField())<0.001) return;
206  fCEvents->Fill(2);
207 
208  fCounter->StoreEvent(aodEvent,fCuts,fUseMCInfo);
209 
210  // trigger class for PbPb C0SMH-B-NOPF-ALLNOTRD
211  TString trigclass=aodEvent->GetFiredTriggerClasses();
212  if(trigclass.Contains("C0SMH-B-NOPF-ALLNOTRD")||trigclass.Contains("C0SMH-B-NOPF-ALL")) fCEvents->Fill(5);
213 
214  if(!fCuts->IsEventSelected(aodEvent)) {
215  if(fCuts->GetWhyRejection()==6) // rejected for Z vertex
216  fCEvents->Fill(6);
217  return;
218  }
219 
220  Bool_t isEvSel=fCuts->IsEventSelected(aodEvent);
221  fCEvents->Fill(3);
222  if(!isEvSel) return;
223 
224  // Load the event
225  // AliInfo(Form("Event %d",fEvents));
226  //if (fEvents%10000 ==0) AliInfo(Form("Event %d",fEvents));
227 
228  // counters for efficiencies
229  Int_t icountReco = 0;
230 
231  //D* and D0 prongs needed to MatchToMC method
232  Int_t pdgDgDStartoD0pi[2]={421,211};
233  Int_t pdgDgD0toKpi[2]={321,211};
234 
235  // AOD primary vertex
236  AliAODVertex *vtx1 = (AliAODVertex*)aodEvent->GetPrimaryVertex();
237  if(!vtx1) return;
238  if(vtx1->GetNContributors()<1) return;
239  fCEvents->Fill(4);
240 
241  if (!arrayDStartoD0pi){
242  AliInfo("Could not find array of HF vertices, skipping the event");
243  return;
244  }else AliDebug(2, Form("Found %d vertices",arrayDStartoD0pi->GetEntriesFast()));
245 
246  Int_t nSelectedAna =0;
247  Int_t nSelectedProd =0;
248 
249  // vHF object is needed to call the method that refills the missing info of the candidates
250  // if they have been deleted in dAOD reconstruction phase
251  // in order to reduce the size of the file
253 
254  // loop over the tracks to search for candidates soft pion
255  for (Int_t iDStartoD0pi = 0; iDStartoD0pi<arrayDStartoD0pi->GetEntriesFast(); iDStartoD0pi++) {
256 
257  // D* candidates and D0 from D*
258  AliAODRecoCascadeHF* dstarD0pi = (AliAODRecoCascadeHF*)arrayDStartoD0pi->At(iDStartoD0pi);
259  Bool_t isDStarCand =kTRUE;
260  if(!(vHF->FillRecoCasc(aodEvent,dstarD0pi,isDStarCand))) {//Fill the data members of the candidate only if they are empty.
261  fCEvents->Fill(12); //monitor how often this fails
262  continue;
263  }
264  if(!dstarD0pi->GetSecondaryVtx()) continue;
265  AliAODRecoDecayHF2Prong* theD0particle = (AliAODRecoDecayHF2Prong*)dstarD0pi->Get2Prong();
266  if (!theD0particle) continue;
267 
268  Int_t isDStar = 0;
269  TClonesArray *mcArray = 0;
270  AliAODMCHeader *mcHeader=0;
271 
272  Bool_t isPrimary=kTRUE;
273  Float_t pdgCode=-2;
274  Float_t trueImpParXY=0.;
275 
276  // mc analysis
277  if(fUseMCInfo){
278  //MC array need for maching
279  mcArray = dynamic_cast<TClonesArray*>(aodEvent->FindListObject(AliAODMCParticle::StdBranchName()));
280  if (!mcArray) {
281  AliError("Could not find Monte-Carlo in AOD");
282  return;
283  }
284  // load MC header
285  mcHeader = (AliAODMCHeader*)aodEvent->GetList()->FindObject(AliAODMCHeader::StdBranchName());
286  if(!mcHeader) {
287  printf("AliAnalysisTaskSEDplus::UserExec: MC header branch not found!\n");
288  return;
289  }
290  // find associated MC particle for D* ->D0toKpi
291  Int_t mcLabel = dstarD0pi->MatchToMC(413,421,pdgDgDStartoD0pi,pdgDgD0toKpi,mcArray);
292  if(mcLabel>=0){
293 
294  AliAODMCParticle *partDSt = (AliAODMCParticle*)mcArray->At(mcLabel);
295  Int_t checkOrigin = CheckOrigin(mcArray,partDSt);
296  if(checkOrigin==5) isPrimary=kFALSE;
297  AliAODMCParticle *dg0 = (AliAODMCParticle*)mcArray->At(partDSt->GetDaughter(0));
298  // AliAODMCParticle *dg01 = (AliAODMCParticle*)mcArray->At(dg0->GetDaughter(0));
299  pdgCode=TMath::Abs(partDSt->GetPdgCode());
300  if(!isPrimary){
301  trueImpParXY=GetTrueImpactParameterD0(mcHeader,mcArray,dg0)*1000.;
302  }
303  isDStar = 1;
304  }else{
305  pdgCode=-1;
306  }
307  }
308 
309  if(pdgCode==-1) AliDebug(2,"No particle assigned! check\n");
310 
311  Int_t ptbin=fCuts->PtBin(dstarD0pi->Pt());
312 
313  // quality selction on tracks and region of interest
314  Int_t isTkSelected = fCuts->IsSelected(dstarD0pi,AliRDHFCuts::kTracks); // quality cuts on tracks
315  if(!isTkSelected) continue;
316 
317  if(!fCuts->IsInFiducialAcceptance(dstarD0pi->Pt(),dstarD0pi->YDstar())) continue;
318 
319 
320  //histos for impact par studies - D0!!!
321  Double_t ptCand = dstarD0pi->Get2Prong()->Pt();
322  Double_t invMass=dstarD0pi->InvMassD0();
323  Double_t impparXY=dstarD0pi->Get2Prong()->ImpParXY()*10000.;
324 
325  Double_t arrayForSparse[3]={invMass,ptCand,impparXY};
326  Double_t arrayForSparseTrue[3]={invMass,ptCand,trueImpParXY};
327 
328  // set the D0 search window bin by bin - useful to calculate side band bkg
329  if (ptbin==0){
330  if(fAnalysis==1){
331  fD0Window=0.035;
332  fPeakWindow=0.03;
333  }else{
334  fD0Window=0.020;
335  fPeakWindow=0.0018;
336  }
337  }
338  if (ptbin==1){
339  if(fAnalysis==1){
340  fD0Window=0.035;
341  fPeakWindow=0.03;
342  }else{
343  fD0Window=0.020;
344  fPeakWindow=0.0018;
345  }
346  }
347  if (ptbin==2){
348  if(fAnalysis==1){
349  fD0Window=0.035;
350  fPeakWindow=0.03;
351  }else{
352  fD0Window=0.020;
353  fPeakWindow=0.0018;
354  }
355  }
356  if (ptbin==3){
357  if(fAnalysis==1){
358  fD0Window=0.035;
359  fPeakWindow=0.03;
360  }else{
361  fD0Window=0.022;
362  fPeakWindow=0.0016;
363  }
364  }
365  if (ptbin==4){
366  if(fAnalysis==1){
367  fD0Window=0.035;
368  fPeakWindow=0.03;
369  }else{
370  fD0Window=0.026;
371  fPeakWindow=0.0014;
372  }
373  }
374  if (ptbin==5){
375  if(fAnalysis==1){
376  fD0Window=0.045;
377  fPeakWindow=0.03;
378  }else{
379  fD0Window=0.026;
380  fPeakWindow=0.0014;
381  }
382  }
383  if (ptbin==6){
384  if(fAnalysis==1){
385  fD0Window=0.045;
386  fPeakWindow=0.03;
387  }else{
388  fD0Window=0.026;
389  fPeakWindow=0.006;
390  }
391  }
392  if (ptbin==7){
393  if(fAnalysis==1){
394  fD0Window=0.055;
395  fPeakWindow=0.03;
396  }else{
397  fD0Window=0.026;
398  fPeakWindow=0.006;
399  }
400  }
401  if (ptbin>7){
402  if(fAnalysis==1){
403  fD0Window=0.074;
404  fPeakWindow=0.03;
405  }else{
406  fD0Window=0.026;
407  fPeakWindow=0.006;
408  }
409  }
410 
411  nSelectedProd++;
412  nSelectedAna++;
413 
414  // check that we are close to signal in the DeltaM - here to save time for PbPb
415  Double_t mPDGD0=TDatabasePDG::Instance()->GetParticle(421)->Mass();
416  Double_t mPDGDstar=TDatabasePDG::Instance()->GetParticle(413)->Mass();
417  Double_t invmassDelta = dstarD0pi->DeltaInvMass();
418 
419  if (TMath::Abs(invmassDelta-(mPDGDstar-mPDGD0))>fPeakWindow) continue;
420  Int_t isSelected=fCuts->IsSelected(dstarD0pi,AliRDHFCuts::kCandidate,aodEvent); //selected
421 
422  // after cuts
423  if(fDoImpParDstar && isSelected){
424  fHistMassPtImpParTCDs[0]->Fill(arrayForSparse);
425  if(isPrimary) fHistMassPtImpParTCDs[1]->Fill(arrayForSparse);
426  else{
427  fHistMassPtImpParTCDs[2]->Fill(arrayForSparse);
428  fHistMassPtImpParTCDs[3]->Fill(arrayForSparseTrue);
429  }
430  }
431 
432  if (fDoDStarVsY && isSelected){
433  ((TH3F*) (fOutputPID->FindObject("deltamassVsyVsPt")))->Fill(dstarD0pi->DeltaInvMass(),dstarD0pi->YDstar(),dstarD0pi->Pt() );
434  }
435 
436 
437  // fill PID
438  FillSpectrum(dstarD0pi,isDStar,fCuts,isSelected,fOutputPID);
439  SideBandBackground(dstarD0pi,fCuts,isSelected, fOutputPID);
440  //WrongSignForDStar(dstarD0pi,fCuts,fOutputPID);
441 
442  //swich off the PID selection
443  fCuts->SetUsePID(kFALSE);
444  Int_t isSelectedNoPID=fCuts->IsSelected(dstarD0pi,AliRDHFCuts::kCandidate, aodEvent); //selected
445  fCuts->SetUsePID(kTRUE);
446 
447  FillSpectrum(dstarD0pi,isDStar,fCuts,isSelectedNoPID,fOutputAll);
448  // SideBandBackground(dstarD0pi,fCuts,isSelectedNoPID, fOutputAll);
449 
450  // rare D search ------
451  if(fDoSearch){
452  TLorentzVector lorentzTrack1(0,0,0,0); // lorentz 4 vector
453  TLorentzVector lorentzTrack2(0,0,0,0); // lorentz 4 vector
454 
455  for (Int_t i=0; i<aodEvent->GetNumberOfTracks(); i++){
456 
457  AliAODTrack* aodTrack = dynamic_cast<AliAODTrack*>(aodEvent->GetTrack(i));
458  if(!aodTrack) AliFatal("Not a standard AOD");
459 
460  if(dstarD0pi->Charge() == aodTrack->Charge()) continue;
461  if((!(aodTrack->GetStatus()&AliESDtrack::kITSrefit)|| (!(aodTrack->GetStatus()&AliESDtrack::kTPCrefit)))) continue;
462  if (TMath::Abs(invmassDelta-(mPDGDstar-mPDGD0))>0.02) continue;
463 
464  //build the D1 mass
465  Double_t mass = TDatabasePDG::Instance()->GetParticle(211)->Mass();
466 
467  lorentzTrack1.SetPxPyPzE( dstarD0pi->Px(),dstarD0pi->Py(), dstarD0pi->Pz(), dstarD0pi->E(413) );
468  lorentzTrack2.SetPxPyPzE( aodTrack->Px(),aodTrack->Py(), aodTrack->Pz(),aodTrack->E(mass) );
469 
470  //D1 mass
471  Double_t d1mass = ((lorentzTrack1+lorentzTrack2).M());
472  //mass difference - at 0.4117 and 0.4566
473  fDeltaMassD1->Fill(d1mass-dstarD0pi->InvMassDstarKpipi());
474  }
475  }
476 
477  if(isDStar == 1) {
478  fTrueDiff2->Fill(dstarD0pi->Pt(),dstarD0pi->DeltaInvMass());
479  }
480 
481  }
482 
483  fCounter->StoreCandidates(aodEvent,nSelectedProd,kTRUE);
484  fCounter->StoreCandidates(aodEvent,nSelectedAna,kFALSE);
485 
486  delete vHF;
487 
488  AliDebug(2, Form("Found %i Reco particles that are D*!!",icountReco));
489 
490  PostData(1,fOutput);
491  PostData(2,fOutputAll);
492  PostData(3,fOutputPID);
493  PostData(5,fCounter);
494 
495 }
496 //________________________________________ terminate ___________________________
498 {
502 
503  //Info("Terminate","");
504  AliAnalysisTaskSE::Terminate();
505 
506  fOutput = dynamic_cast<TList*> (GetOutputData(1));
507  if (!fOutput) {
508  printf("ERROR: fOutput not available\n");
509  return;
510  }
511 
512  fCEvents = dynamic_cast<TH1F*>(fOutput->FindObject("fCEvents"));
513  fDeltaMassD1 = dynamic_cast<TH1F*>(fOutput->FindObject("fDeltaMassD1"));
514  fTrueDiff2 = dynamic_cast<TH2F*>(fOutput->FindObject("fTrueDiff2"));
515 
516  fOutputAll = dynamic_cast<TList*> (GetOutputData(1));
517  if (!fOutputAll) {
518  printf("ERROR: fOutputAll not available\n");
519  return;
520  }
521  fOutputPID = dynamic_cast<TList*> (GetOutputData(2));
522  if (!fOutputPID) {
523  printf("ERROR: fOutputPID not available\n");
524  return;
525  }
526  return;
527 }
528 //___________________________________________________________________________
531  Info("UserCreateOutputObjects","CreateOutputObjects of task %s\n", GetName());
532 
533  //slot #1
534  //OpenFile(1);
535  fOutput = new TList();
536  fOutput->SetOwner();
537  fOutput->SetName("chist0");
538 
539  fOutputAll = new TList();
540  fOutputAll->SetOwner();
541  fOutputAll->SetName("listAll");
542 
543  fOutputPID = new TList();
544  fOutputPID->SetOwner();
545  fOutputPID->SetName("listPID");
546 
547  // define histograms
549 
550  //Counter for Normalization
551  fCounter = new AliNormalizationCounter(Form("%s",GetOutputSlot(5)->GetContainer()->GetName()));
552  fCounter->Init();
553 
555 
556  PostData(1,fOutput);
557  PostData(2,fOutputAll);
558  PostData(3,fOutputPID);
559 
560  return;
561 }
562 //___________________________________ hiostograms _______________________________________
565 
566  fCEvents = new TH1F("fCEvents","conter",13,0,13);
567  fCEvents->SetStats(kTRUE);
568  fCEvents->GetXaxis()->SetTitle("1");
569  fCEvents->GetYaxis()->SetTitle("counts");
570  fCEvents->GetXaxis()->SetBinLabel(2,"no. of events");
571  fCEvents->GetXaxis()->SetBinLabel(3,"good prim vtx and B field");
572  fCEvents->GetXaxis()->SetBinLabel(4,"no event selected");
573  fCEvents->GetXaxis()->SetBinLabel(5,"no vtx contributors");
574  fCEvents->GetXaxis()->SetBinLabel(6,"trigger for PbPb");
575  fCEvents->GetXaxis()->SetBinLabel(7,"no z vtx");
576  fCEvents->GetXaxis()->SetBinLabel(12,"no. of D0 fail to be rec");
577  fOutput->Add(fCEvents);
578 
579  fTrueDiff2 = new TH2F("DiffDstar_pt","True Reco diff vs pt",200,0,15,900,0,0.3);
580  fOutput->Add(fTrueDiff2);
581 
582  fDeltaMassD1 = new TH1F("DeltaMassD1","delta mass d1",600,0,0.8);
583  fOutput->Add(fDeltaMassD1);
584 
585  const Int_t nhist=14;
586  TString nameMass=" ", nameSgn=" ", nameBkg=" ";
587 
588  for(Int_t i=-2;i<nhist;i++){
589  nameMass="histDeltaMass_";
590  nameMass+=i+1;
591  nameSgn="histDeltaSgn_";
592  nameSgn+=i+1;
593  nameBkg="histDeltaBkg_";
594  nameBkg+=i+1;
595 
596  if (i==-2) {
597  nameMass="histDeltaMass";
598  nameSgn="histDeltaSgn";
599  nameBkg="histDeltaBkg";
600  }
601 
602  TH1F* spectrumMass = new TH1F(nameMass.Data(),"D^{*}-D^{0} invariant mass; #DeltaM [GeV/c^{2}]; Entries",700,0.13,0.2);
603  TH1F* spectrumSgn = new TH1F(nameSgn.Data(), "D^{*}-D^{0} Signal invariant mass - MC; #DeltaM [GeV/c^{2}]; Entries",700,0.13,0.2);
604  TH1F* spectrumBkg = new TH1F(nameBkg.Data(), "D^{*}-D^{0} Background invariant mass - MC; #DeltaM [GeV/c^{2}]; Entries",700,0.13,0.2);
605 
606  nameMass="histD0Mass_";
607  nameMass+=i+1;
608  nameSgn="histD0Sgn_";
609  nameSgn+=i+1;
610  nameBkg="histD0Bkg_";
611  nameBkg+=i+1;
612 
613  if (i==-2) {
614  nameMass="histD0Mass";
615  nameSgn="histD0Sgn";
616  nameBkg="histD0Bkg";
617  }
618 
619  TH1F* spectrumD0Mass = new TH1F(nameMass.Data(),"D^{0} invariant mass; M(D^{0}) [GeV/c^{2}]; Entries",200,1.75,1.95);
620  TH1F* spectrumD0Sgn = new TH1F(nameSgn.Data(), "D^{0} Signal invariant mass - MC; M(D^{0}) [GeV/c^{2}]; Entries",200,1.75,1.95);
621  TH1F* spectrumD0Bkg = new TH1F(nameBkg.Data(), "D^{0} Background invariant mass - MC; M(D^{0}) [GeV/c^{2}]; Entries",200,1.75,1.95);
622 
623  nameMass="histDstarMass_";
624  nameMass+=i+1;
625  nameSgn="histDstarSgn_";
626  nameSgn+=i+1;
627  nameBkg="histDstarBkg_";
628  nameBkg+=i+1;
629 
630  if (i==-2) {
631  nameMass="histDstarMass";
632  nameSgn="histDstarSgn";
633  nameBkg="histDstarBkg";
634  }
635 
636  TH1F* spectrumDstarMass = new TH1F(nameMass.Data(),"D^{*} invariant mass; M(D^{*}) [GeV/c^{2}]; Entries",200,1.9,2.1);
637  TH1F* spectrumDstarSgn = new TH1F(nameSgn.Data(), "D^{*} Signal invariant mass - MC; M(D^{*}) [GeV/c^{2}]; Entries",200,1.9,2.1);
638  TH1F* spectrumDstarBkg = new TH1F(nameBkg.Data(), "D^{*} Background invariant mass - MC; M(D^{*}) [GeV/c^{2}]; Entries",200,1.9,2.1);
639 
640  nameMass="histSideBandMass_";
641  nameMass+=i+1;
642  if (i==-2) {
643  nameMass="histSideBandMass";
644  }
645 
646  TH1F* spectrumSideBandMass = new TH1F(nameMass.Data(),"D^{*}-D^{0} sideband mass; M(D^{*}) [GeV/c^{2}]; Entries",200,0.1,0.2);
647 
648  nameMass="histWrongSignMass_";
649  nameMass+=i+1;
650  if (i==-2) {
651  nameMass="histWrongSignMass";
652  }
653 
654  TH1F* spectrumWrongSignMass = new TH1F(nameMass.Data(),"D^{*}-D^{0} wrongsign mass; M(D^{*}) [GeV/c^{2}]; Entries",200,0.1,0.2);
655 
656 
657  spectrumMass->Sumw2();
658  spectrumSgn->Sumw2();
659  spectrumBkg->Sumw2();
660 
661  spectrumMass->SetLineColor(6);
662  spectrumSgn->SetLineColor(2);
663  spectrumBkg->SetLineColor(4);
664 
665  spectrumMass->SetMarkerStyle(20);
666  spectrumSgn->SetMarkerStyle(20);
667  spectrumBkg->SetMarkerStyle(20);
668  spectrumMass->SetMarkerSize(0.6);
669  spectrumSgn->SetMarkerSize(0.6);
670  spectrumBkg->SetMarkerSize(0.6);
671  spectrumMass->SetMarkerColor(6);
672  spectrumSgn->SetMarkerColor(2);
673  spectrumBkg->SetMarkerColor(4);
674 
675  spectrumD0Mass->Sumw2();
676  spectrumD0Sgn->Sumw2();
677  spectrumD0Bkg->Sumw2();
678 
679  spectrumD0Mass->SetLineColor(6);
680  spectrumD0Sgn->SetLineColor(2);
681  spectrumD0Bkg->SetLineColor(4);
682 
683  spectrumD0Mass->SetMarkerStyle(20);
684  spectrumD0Sgn->SetMarkerStyle(20);
685  spectrumD0Bkg->SetMarkerStyle(20);
686  spectrumD0Mass->SetMarkerSize(0.6);
687  spectrumD0Sgn->SetMarkerSize(0.6);
688  spectrumD0Bkg->SetMarkerSize(0.6);
689  spectrumD0Mass->SetMarkerColor(6);
690  spectrumD0Sgn->SetMarkerColor(2);
691  spectrumD0Bkg->SetMarkerColor(4);
692 
693  spectrumDstarMass->Sumw2();
694  spectrumDstarSgn->Sumw2();
695  spectrumDstarBkg->Sumw2();
696 
697  spectrumDstarMass->SetLineColor(6);
698  spectrumDstarSgn->SetLineColor(2);
699  spectrumDstarBkg->SetLineColor(4);
700 
701  spectrumDstarMass->SetMarkerStyle(20);
702  spectrumDstarSgn->SetMarkerStyle(20);
703  spectrumDstarBkg->SetMarkerStyle(20);
704  spectrumDstarMass->SetMarkerSize(0.6);
705  spectrumDstarSgn->SetMarkerSize(0.6);
706  spectrumDstarBkg->SetMarkerSize(0.6);
707  spectrumDstarMass->SetMarkerColor(6);
708  spectrumDstarSgn->SetMarkerColor(2);
709  spectrumDstarBkg->SetMarkerColor(4);
710 
711  spectrumSideBandMass->Sumw2();
712  spectrumSideBandMass->SetLineColor(4);
713  spectrumSideBandMass->SetMarkerStyle(20);
714  spectrumSideBandMass->SetMarkerSize(0.6);
715  spectrumSideBandMass->SetMarkerColor(4);
716 
717  spectrumWrongSignMass->Sumw2();
718  spectrumWrongSignMass->SetLineColor(4);
719  spectrumWrongSignMass->SetMarkerStyle(20);
720  spectrumWrongSignMass->SetMarkerSize(0.6);
721  spectrumWrongSignMass->SetMarkerColor(4);
722 
723  TH1F* allMass = (TH1F*)spectrumMass->Clone();
724  TH1F* allSgn = (TH1F*)spectrumSgn->Clone();
725  TH1F* allBkg = (TH1F*)spectrumBkg->Clone();
726 
727  TH1F* pidMass = (TH1F*)spectrumMass->Clone();
728  TH1F* pidSgn = (TH1F*)spectrumSgn->Clone();
729  TH1F* pidBkg = (TH1F*)spectrumBkg->Clone();
730 
731  fOutputAll->Add(allMass);
732  fOutputAll->Add(allSgn);
733  fOutputAll->Add(allBkg);
734 
735  fOutputPID->Add(pidMass);
736  fOutputPID->Add(pidSgn);
737  fOutputPID->Add(pidBkg);
738 
739  TH1F* allD0Mass = (TH1F*)spectrumD0Mass->Clone();
740  TH1F* allD0Sgn = (TH1F*)spectrumD0Sgn->Clone();
741  TH1F* allD0Bkg = (TH1F*)spectrumD0Bkg->Clone();
742 
743  TH1F* pidD0Mass = (TH1F*)spectrumD0Mass->Clone();
744  TH1F* pidD0Sgn = (TH1F*)spectrumD0Sgn->Clone();
745  TH1F* pidD0Bkg = (TH1F*)spectrumD0Bkg->Clone();
746 
747  fOutputAll->Add(allD0Mass);
748  fOutputAll->Add(allD0Sgn);
749  fOutputAll->Add(allD0Bkg);
750 
751  fOutputPID->Add(pidD0Mass);
752  fOutputPID->Add(pidD0Sgn);
753  fOutputPID->Add(pidD0Bkg);
754 
755  TH1F* allDstarMass = (TH1F*)spectrumDstarMass->Clone();
756  TH1F* allDstarSgn = (TH1F*)spectrumDstarSgn->Clone();
757  TH1F* allDstarBkg = (TH1F*)spectrumDstarBkg->Clone();
758 
759  TH1F* pidDstarMass = (TH1F*)spectrumDstarMass->Clone();
760  TH1F* pidDstarSgn = (TH1F*)spectrumDstarSgn->Clone();
761  TH1F* pidDstarBkg = (TH1F*)spectrumDstarBkg->Clone();
762 
763  fOutputAll->Add(allDstarMass);
764  fOutputAll->Add(allDstarSgn);
765  fOutputAll->Add(allDstarBkg);
766 
767  fOutputPID->Add(pidDstarMass);
768  fOutputPID->Add(pidDstarSgn);
769  fOutputPID->Add(pidDstarBkg);
770 
771  TH1F* allSideBandMass = (TH1F*)spectrumSideBandMass->Clone();
772  TH1F* pidSideBandMass = (TH1F*)spectrumSideBandMass->Clone();
773 
774  fOutputAll->Add(allSideBandMass);
775  fOutputPID->Add(pidSideBandMass);
776 
777  TH1F* allWrongSignMass = (TH1F*)spectrumWrongSignMass->Clone();
778  TH1F* pidWrongSignMass = (TH1F*)spectrumWrongSignMass->Clone();
779 
780  fOutputAll->Add(allWrongSignMass);
781  fOutputPID->Add(pidWrongSignMass);
782 
783  }
784 
785  // pt spectra
786  nameMass="ptMass";
787  nameSgn="ptSgn";
788  nameBkg="ptBkg";
789 
790  TH1F* ptspectrumMass = new TH1F(nameMass.Data(),"D^{*} p_{T}; p_{T} [GeV]; Entries",200,0,10);
791  TH1F* ptspectrumSgn = new TH1F(nameSgn.Data(), "D^{*} Signal p_{T} - MC; p_{T} [GeV]; Entries",200,0,10);
792  TH1F* ptspectrumBkg = new TH1F(nameBkg.Data(), "D^{*} Background p_{T} - MC; p_{T} [GeV]; Entries",200,0,10);
793 
794  ptspectrumMass->Sumw2();
795  ptspectrumSgn->Sumw2();
796  ptspectrumBkg->Sumw2();
797 
798  ptspectrumMass->SetLineColor(6);
799  ptspectrumSgn->SetLineColor(2);
800  ptspectrumBkg->SetLineColor(4);
801 
802  ptspectrumMass->SetMarkerStyle(20);
803  ptspectrumSgn->SetMarkerStyle(20);
804  ptspectrumBkg->SetMarkerStyle(20);
805  ptspectrumMass->SetMarkerSize(0.6);
806  ptspectrumSgn->SetMarkerSize(0.6);
807  ptspectrumBkg->SetMarkerSize(0.6);
808  ptspectrumMass->SetMarkerColor(6);
809  ptspectrumSgn->SetMarkerColor(2);
810  ptspectrumBkg->SetMarkerColor(4);
811 
812  TH1F* ptallMass = (TH1F*)ptspectrumMass->Clone();
813  TH1F* ptallSgn = (TH1F*)ptspectrumSgn->Clone();
814  TH1F* ptallBkg = (TH1F*)ptspectrumBkg->Clone();
815 
816  TH1F* ptpidMass = (TH1F*)ptspectrumMass->Clone();
817  TH1F* ptpidSgn = (TH1F*)ptspectrumSgn->Clone();
818  TH1F* ptpidBkg = (TH1F*)ptspectrumBkg->Clone();
819 
820  fOutputAll->Add(ptallMass);
821  fOutputAll->Add(ptallSgn);
822  fOutputAll->Add(ptallBkg);
823 
824  fOutputPID->Add(ptpidMass);
825  fOutputPID->Add(ptpidSgn);
826  fOutputPID->Add(ptpidBkg);
827 
828  // eta spectra
829  nameMass="etaMass";
830  nameSgn="etaSgn";
831  nameBkg="etaBkg";
832 
833  TH1F* etaspectrumMass = new TH1F(nameMass.Data(),"D^{*} #eta; #eta; Entries",200,-1,1);
834  TH1F* etaspectrumSgn = new TH1F(nameSgn.Data(), "D^{*} Signal #eta - MC; #eta; Entries",200,-1,1);
835  TH1F* etaspectrumBkg = new TH1F(nameBkg.Data(), "D^{*} Background #eta - MC; #eta; Entries",200,-1,1);
836 
837  etaspectrumMass->Sumw2();
838  etaspectrumSgn->Sumw2();
839  etaspectrumBkg->Sumw2();
840 
841  etaspectrumMass->SetLineColor(6);
842  etaspectrumSgn->SetLineColor(2);
843  etaspectrumBkg->SetLineColor(4);
844 
845  etaspectrumMass->SetMarkerStyle(20);
846  etaspectrumSgn->SetMarkerStyle(20);
847  etaspectrumBkg->SetMarkerStyle(20);
848  etaspectrumMass->SetMarkerSize(0.6);
849  etaspectrumSgn->SetMarkerSize(0.6);
850  etaspectrumBkg->SetMarkerSize(0.6);
851  etaspectrumMass->SetMarkerColor(6);
852  etaspectrumSgn->SetMarkerColor(2);
853  etaspectrumBkg->SetMarkerColor(4);
854 
855  TH1F* etaallMass = (TH1F*)etaspectrumMass->Clone();
856  TH1F* etaallSgn = (TH1F*)etaspectrumSgn->Clone();
857  TH1F* etaallBkg = (TH1F*)etaspectrumBkg->Clone();
858 
859  TH1F* etapidMass = (TH1F*)etaspectrumMass->Clone();
860  TH1F* etapidSgn = (TH1F*)etaspectrumSgn->Clone();
861  TH1F* etapidBkg = (TH1F*)etaspectrumBkg->Clone();
862 
863  fOutputAll->Add(etaallMass);
864  fOutputAll->Add(etaallSgn);
865  fOutputAll->Add(etaallBkg);
866 
867  fOutputPID->Add(etapidMass);
868  fOutputPID->Add(etapidSgn);
869  fOutputPID->Add(etapidBkg);
870 
871  if (fDoDStarVsY){
872  TH3F* deltamassVsyVsPtPID = new TH3F("deltamassVsyVsPt", "delta mass Vs y Vs pT; #DeltaM [GeV/c^{2}]; y; p_{T} [GeV/c]", 700,0.13,0.2, 40, -1, 1, 36, 0., 36.);
873  fOutputPID->Add(deltamassVsyVsPtPID);
874  }
875  return;
876 }
877 //________________________________________________________________________
879  //
881  //
882 
883  if(!isSel) return;
884 
885  // D0 window
886  Double_t mPDGD0=TDatabasePDG::Instance()->GetParticle(421)->Mass();
887  Double_t invmassD0 = part->InvMassD0();
888  if (TMath::Abs(invmassD0-mPDGD0)>fD0Window) return;
889 
890 
891  Int_t ptbin=cuts->PtBin(part->Pt());
892  Double_t pt = part->Pt();
893  Double_t eta = part->Eta();
894 
895  Double_t invmassDelta = part->DeltaInvMass();
896  Double_t invmassDstar = part->InvMassDstarKpipi();
897 
898  TString fillthis="";
899  Bool_t massInRange=kFALSE;
900 
901  Double_t mPDGDstar=TDatabasePDG::Instance()->GetParticle(413)->Mass();
902 
903  // delta M(Kpipi)-M(Kpi)
904  if (TMath::Abs(invmassDelta-(mPDGDstar-mPDGD0))<fPeakWindow) massInRange=kTRUE;
905 
906  if(fUseMCInfo) {
907  if(isDStar==1) {
908  fillthis="histD0Sgn_";
909  fillthis+=ptbin;
910  ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0);
911  fillthis="histD0Sgn";
912  ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0);
913  fillthis="histDstarSgn_";
914  fillthis+=ptbin;
915  ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassDstar);
916  fillthis="histDstarSgn";
917  ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassDstar);
918  fillthis="histDeltaSgn_";
919  fillthis+=ptbin;
920  ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassDelta);
921  fillthis="histDeltaSgn";
922  ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassDelta);
923  if (massInRange) {
924  fillthis="ptSgn";
925  ((TH1F*)(listout->FindObject(fillthis)))->Fill(pt);
926  fillthis="etaSgn";
927  ((TH1F*)(listout->FindObject(fillthis)))->Fill(eta);
928  }
929  }
930  else {//background
931  fillthis="histD0Bkg_";
932  fillthis+=ptbin;
933  ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0);
934  fillthis="histD0Bkg";
935  ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0);
936  fillthis="histDstarBkg_";
937  fillthis+=ptbin;
938  ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassDstar);
939  fillthis="histDstarBkg";
940  ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassDstar);
941  fillthis="histDeltaBkg_";
942  fillthis+=ptbin;
943  ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassDelta);
944  fillthis="histDeltaBkg";
945  ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassDelta);
946  if (massInRange) {
947  fillthis="ptBkg";
948  ((TH1F*)(listout->FindObject(fillthis)))->Fill(pt);
949  fillthis="etaBkg";
950  ((TH1F*)(listout->FindObject(fillthis)))->Fill(eta);
951  }
952  }
953  }
954  //no MC info, just cut selection
955  fillthis="histD0Mass_";
956  fillthis+=ptbin;
957  ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0);
958  fillthis="histD0Mass";
959  ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0);
960  fillthis="histDstarMass_";
961  fillthis+=ptbin;
962  ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassDstar);
963  fillthis="histDstarMass";
964  ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassDstar);
965  fillthis="histDeltaMass_";
966  fillthis+=ptbin;
967  ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassDelta);
968  fillthis="histDeltaMass";
969  ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassDelta);
970 
971  if (massInRange) {
972  fillthis="ptMass";
973  ((TH1F*)(listout->FindObject(fillthis)))->Fill(pt);
974  fillthis="etaMass";
975  ((TH1F*)(listout->FindObject(fillthis)))->Fill(eta);
976  }
977 
978  return;
979 }
980 //______________________________ side band background for D*___________________________________
982 
986 
987  if(!isSel) return;
988 
989  Int_t ptbin=cuts->PtBin(part->Pt());
990 
991  // select the side bands intervall
992  Double_t invmassD0 = part->InvMassD0();
993  if(TMath::Abs(invmassD0-1.865)>4*fD0Window && TMath::Abs(invmassD0-1.865)<8*fD0Window){
994 
995  // for pt and eta
996  Double_t invmassDelta = part->DeltaInvMass();
997 
998  TString fillthis="";
999  fillthis="histSideBandMass_";
1000  fillthis+=ptbin;
1001  ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassDelta);
1002  fillthis="histSideBandMass";
1003  ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassDelta);
1004 
1005  }
1006 }
1007 //________________________________________________________________________________________________________________
1009  //
1011  //
1012  Int_t ptbin=cuts->PtBin(part->Pt());
1013 
1014  Double_t mPDGD0=TDatabasePDG::Instance()->GetParticle(421)->Mass();
1015  Double_t invmassD0 = part->InvMassD0();
1016  if (TMath::Abs(invmassD0-mPDGD0)>fD0Window) return;
1017 
1018  AliAODRecoDecayHF2Prong* theD0particle = (AliAODRecoDecayHF2Prong*)part->Get2Prong();
1019 
1020  Int_t okD0WrongSign;
1021  Double_t wrongMassD0=0.;
1022 
1023  Int_t isSelected=cuts->IsSelected(part,AliRDHFCuts::kCandidate); //selected
1024  if (!isSelected){
1025  return;
1026  }
1027 
1028  okD0WrongSign = 1;
1029 
1030  //if is D*+ than assume D0bar
1031  if(part->Charge()>0 && (isSelected ==1)) {
1032  okD0WrongSign = 0;
1033  }
1034 
1035  // assign the wrong mass in case the cuts return both D0 and D0bar
1036  if(part->Charge()>0 && (isSelected ==3)) {
1037  okD0WrongSign = 0;
1038  }
1039 
1040  //wrong D0 inv mass
1041  if(okD0WrongSign!=0){
1042  wrongMassD0 = theD0particle->InvMassD0();
1043  }else if(okD0WrongSign==0){
1044  wrongMassD0 = theD0particle->InvMassD0bar();
1045  }
1046 
1047  if(TMath::Abs(wrongMassD0-1.865)<fD0Window){
1048 
1049  // wrong D* inv mass
1050  Double_t e[3];
1051  if (part->Charge()>0){
1052  e[0]=theD0particle->EProng(0,321);
1053  e[1]=theD0particle->EProng(1,211);
1054  }else{
1055  e[0]=theD0particle->EProng(0,211);
1056  e[1]=theD0particle->EProng(1,321);
1057  }
1058  e[2]=part->EProng(0,211);
1059 
1060  Double_t esum = e[0]+e[1]+e[2];
1061  Double_t pds = part->P();
1062 
1063  Double_t wrongMassDstar = TMath::Sqrt(esum*esum-pds*pds);
1064 
1065  TString fillthis="";
1066  fillthis="histWrongSignMass_";
1067  fillthis+=ptbin;
1068  ((TH1F*)(listout->FindObject(fillthis)))->Fill(wrongMassDstar-wrongMassD0);
1069  fillthis="histWrongSignMass";
1070  ((TH1F*)(listout->FindObject(fillthis)))->Fill(wrongMassDstar-wrongMassD0);
1071 
1072  }
1073 }
1074 
1075 //-------------------------------------------------------------------------------
1076 Int_t AliAnalysisTaskSEDStarSpectra::CheckOrigin(TClonesArray* arrayMC, const AliAODMCParticle *mcPartCandidate) const {
1077  //
1078  // checking whether the mother of the particles come from a charm or a bottom quark
1079  //
1080 
1081  Int_t pdgGranma = 0;
1082  Int_t mother = 0;
1083  mother = mcPartCandidate->GetMother();
1084  Int_t istep = 0;
1085  Int_t abspdgGranma =0;
1086  Bool_t isFromB=kFALSE;
1087  while (mother >0 ){
1088  istep++;
1089  AliAODMCParticle* mcGranma = dynamic_cast<AliAODMCParticle*>(arrayMC->At(mother));
1090  if (mcGranma){
1091  pdgGranma = mcGranma->GetPdgCode();
1092  abspdgGranma = TMath::Abs(pdgGranma);
1093  if ((abspdgGranma > 500 && abspdgGranma < 600) || (abspdgGranma > 5000 && abspdgGranma < 6000)){
1094  isFromB=kTRUE;
1095  }
1096  mother = mcGranma->GetMother();
1097  }else{
1098  AliError("Failed casting the mother particle!");
1099  break;
1100  }
1101  }
1102 
1103  if(isFromB) return 5;
1104  else return 4;
1105 }
1106 //-------------------------------------------------------------------------------------
1107 Float_t AliAnalysisTaskSEDStarSpectra::GetTrueImpactParameterD0(const AliAODMCHeader *mcHeader, TClonesArray* arrayMC, const AliAODMCParticle *partDp) const {
1109 
1110  Double_t vtxTrue[3];
1111  mcHeader->GetVertex(vtxTrue);
1112  Double_t origD[3];
1113  partDp->XvYvZv(origD);
1114  Short_t charge=partDp->Charge();
1115  Double_t pXdauTrue[3],pYdauTrue[3],pZdauTrue[3];
1116  Int_t labelFirstDau = partDp->GetDaughter(0);
1117 
1118  Int_t nDau=partDp->GetNDaughters();
1119 
1120  Int_t theDau=0;
1121  if(nDau==2){
1122  for(Int_t iDau=0; iDau<2; iDau++){
1123  Int_t ind = labelFirstDau+iDau;
1124  AliAODMCParticle* part = dynamic_cast<AliAODMCParticle*>(arrayMC->At(ind));
1125  if(!part){
1126  AliError("Daughter particle not found in MC array");
1127  return 99999.;
1128  }
1129  Int_t pdgCode=TMath::Abs(part->GetPdgCode());
1130  if(pdgCode==211 || pdgCode==321){
1131  pXdauTrue[theDau]=part->Px();
1132  pYdauTrue[theDau]=part->Py();
1133  pZdauTrue[theDau]=part->Pz();
1134  ++theDau;
1135  }
1136  }
1137  }
1138  if(theDau!=2){
1139  AliError("Wrong number of decay prongs");
1140  return 99999.;
1141  }
1142 
1143  Double_t d0dummy[3]={0.,0.,0.};
1144  AliAODRecoDecayHF aodD0MC(vtxTrue,origD,3,charge,pXdauTrue,pYdauTrue,pZdauTrue,d0dummy);
1145  return aodD0MC.ImpParXY();
1146 
1147 }
1148 //______________________________________________________-
1151 
1152  Int_t nbins[3]={400,200,fNImpParBins};
1153  Double_t xmin[3]={1.75,0.,fLowerImpPar};
1154  Double_t xmax[3]={1.98,20.,fHigherImpPar};
1155 
1156  fHistMassPtImpParTCDs[0]=new THnSparseF("hMassPtImpParAll",
1157  "Mass vs. pt vs.imppar - All",
1158  3,nbins,xmin,xmax);
1159  fHistMassPtImpParTCDs[1]=new THnSparseF("hMassPtImpParPrompt",
1160  "Mass vs. pt vs.imppar - promptD",
1161  3,nbins,xmin,xmax);
1162  fHistMassPtImpParTCDs[2]=new THnSparseF("hMassPtImpParBfeed",
1163  "Mass vs. pt vs.imppar - DfromB",
1164  3,nbins,xmin,xmax);
1165  fHistMassPtImpParTCDs[3]=new THnSparseF("hMassPtImpParTrueBfeed",
1166  "Mass vs. pt vs.true imppar -DfromB",
1167  3,nbins,xmin,xmax);
1168  fHistMassPtImpParTCDs[4]=new THnSparseF("hMassPtImpParBkg",
1169  "Mass vs. pt vs.imppar - backgr.",
1170  3,nbins,xmin,xmax);
1171 
1172  for(Int_t i=0; i<5;i++){
1173  fOutput->Add(fHistMassPtImpParTCDs[i]);
1174  }
1175 }
1176 
Int_t charge
TH1F * fCEvents
Cuts - sent to output slot 3.
double Double_t
Definition: External.C:58
Definition: External.C:260
void StoreCandidates(AliVEvent *, Int_t nCand=0, Bool_t flagFilter=kTRUE)
Definition: External.C:236
Double_t DeltaInvMass() const
Float_t fHigherImpPar
lower limit in impact parameter (um)
Double_t YDstar() const
Int_t CheckOrigin(TClonesArray *arrayMC, const AliAODMCParticle *mcPartCandidate) const
Int_t MatchToMC(Int_t pdgabs, Int_t pdgabs2prong, Int_t *pdgDg, Int_t *pdgDg2prong, TClonesArray *mcArray, Bool_t isV0=kFALSE) const
Float_t fLowerImpPar
nunber of bins in impact parameter histos
void FillSpectrum(AliAODRecoCascadeHF *part, Int_t isDStar, AliRDHFCutsDStartoKpipi *cuts, Int_t isSel, TList *listout)
histos
AliRDHFCutsDStartoKpipi * fCuts
n sigma for kaon PID
Double_t mass
Double_t ImpParXY() const
char Char_t
Definition: External.C:18
Int_t GetWhyRejection() const
Definition: AliRDHFCuts.h:297
Bool_t fDoDStarVsY
higher limit in impact parameter (um)
Float_t GetTrueImpactParameterD0(const AliAODMCHeader *mcHeader, TClonesArray *arrayMC, const AliAODMCParticle *partDp) const
AliNormalizationCounter * fCounter
!Counter for normalization slot 4
Bool_t FillRecoCasc(AliVEvent *event, AliAODRecoCascadeHF *rc, Bool_t isDStar, Bool_t recoSecVtx=kFALSE)
virtual void UserCreateOutputObjects()
Implementation of interface methods.
int Int_t
Definition: External.C:63
float Float_t
Definition: External.C:68
virtual void UserExec(Option_t *option)
virtual Int_t IsSelected(TObject *obj, Int_t selectionLevel, AliAODEvent *aod)
virtual void Terminate(Option_t *option)
void SideBandBackground(AliAODRecoCascadeHF *part, AliRDHFCutsDStartoKpipi *cuts, Int_t isSel, TList *listout)
Background simulation.
Double_t fPeakWindow
select width on D0Mass
short Short_t
Definition: External.C:23
void WrongSignForDStar(AliAODRecoCascadeHF *part, AliRDHFCutsDStartoKpipi *cuts, TList *listout)
Bool_t IsEventSelected(AliVEvent *event)
void StoreEvent(AliVEvent *, AliRDHFCuts *, Bool_t mc=kFALSE, Int_t multiplicity=-9999, Double_t spherocity=-99.)
void SetUsePID(Bool_t flag=kTRUE)
Definition: AliRDHFCuts.h:204
virtual Bool_t IsInFiducialAcceptance(Double_t pt, Double_t y) const
ClassImp(AliAnalysisTaskCRC) AliAnalysisTaskCRC
THnSparseF * fHistMassPtImpParTCDs[5]
flag to enable D* vs y
Double_t InvMassD0() const
const char Option_t
Definition: External.C:48
Double_t InvMassDstarKpipi() const
const Int_t nbins
bool Bool_t
Definition: External.C:53
AliAODRecoDecayHF2Prong * Get2Prong() const
Int_t PtBin(Double_t pt) const
Bool_t fUseMCInfo
select width on DstarMass