]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGHF/vertexingHF/AliAnalysisTaskCheckHFMCProd.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
updated for photonic electron efficiency
[u/mrichter/AliRoot.git] / PWGHF / vertexingHF / AliAnalysisTaskCheckHFMCProd.cxx
CommitLineData
7b6a4dcd 1#include "AliAnalysisTaskSE.h"
2#include "AliAnalysisManager.h"
3#include "AliAnalysisDataContainer.h"
4#include "AliESDEvent.h"
5#include "AliStack.h"
6#include "AliCentrality.h"
7#include "AliMCEventHandler.h"
8#include "AliMCEvent.h"
9#include "AliMultiplicity.h"
10#include <TParticle.h>
11#include <TSystem.h>
12#include <TTree.h>
13#include <TNtuple.h>
14#include <TH1F.h>
15#include <TH2F.h>
16#include <TChain.h>
17#include "AliESDInputHandlerRP.h"
18#include "AliAnalysisTaskCheckHFMCProd.h"
19
20/**************************************************************************
21 * Copyright(c) 1998-2012, ALICE Experiment at CERN, All rights reserved. *
22 * *
23 * Author: The ALICE Off-line Project. *
24 * Contributors are mentioned in the code where appropriate. *
25 * *
26 * Permission to use, copy, modify and distribute this software and its *
27 * documentation strictly for non-commercial purposes is hereby granted *
28 * without fee, provided that the above copyright notice appears in all *
29 * copies and that both the copyright notice and this permission notice *
30 * appear in the supporting documentation. The authors make no claims *
31 * about the suitability of this software for any purpose. It is *
32 * provided "as is" without express or implied warranty. *
33 **************************************************************************/
34
35/* $Id$ */
36
37//*************************************************************************
38// Implementation of class AliAnalysisTaskCheckHFMCProd
39// AliAnalysisTask to check MC production at ESD+Kine level
40//
41//
42// Authors: F. Prino, prino@to.infn.it
43//
44//*************************************************************************
45
46ClassImp(AliAnalysisTaskCheckHFMCProd)
47//______________________________________________________________________________
48AliAnalysisTaskCheckHFMCProd::AliAnalysisTaskCheckHFMCProd() : AliAnalysisTaskSE("HFMCChecks"),
49 fOutput(0),
50 fHistoNEvents(0),
51 fHistoTracks(0),
52 fHistoSelTracks(0),
53 fHistoTracklets(0),
54 fHistoSPD3DVtxX(0),
55 fHistoSPD3DVtxY(0),
56 fHistoSPD3DVtxZ(0),
57 fHistoSPDZVtxX(0),
58 fHistoSPDZVtxY(0),
59 fHistoSPDZVtxZ(0),
60 fHistoTRKVtxX(0),
61 fHistoTRKVtxY(0),
62 fHistoTRKVtxZ(0),
63 fHistoNcharmed(0),
64 fHistoNbVsNc(0),
65 fPbPb(kFALSE),
66 fReadMC(kTRUE)
67{
68 //
69 DefineInput(0, TChain::Class());
70 DefineOutput(1, TList::Class());
71}
72
73
74//___________________________________________________________________________
75AliAnalysisTaskCheckHFMCProd::~AliAnalysisTaskCheckHFMCProd(){
76 //
77 if (AliAnalysisManager::GetAnalysisManager()->IsProofMode()) return;
78 if (fOutput) {
79 delete fOutput;
80 fOutput = 0;
81 }
82}
83
84//___________________________________________________________________________
85void AliAnalysisTaskCheckHFMCProd::UserCreateOutputObjects() {
86 // create output histos
87
88 fOutput = new TList();
89 fOutput->SetOwner();
90 fOutput->SetName("OutputHistos");
91
92 fHistoNEvents = new TH1F("hNEvents", "Number of processed events",3,-0.5,2.5);
93 fHistoNEvents->Sumw2();
94 fHistoNEvents->SetMinimum(0);
95 fOutput->Add(fHistoNEvents);
96
97 Double_t maxMult=100.;
98 if(fPbPb) maxMult=10000.;
99 fHistoTracks = new TH1F("hTracks","",100,0.,maxMult*2);
100 fHistoTracks->Sumw2();
101 fOutput->Add(fHistoTracks);
102 fHistoSelTracks = new TH1F("hSelTracks","",100,0.,maxMult);
103 fHistoSelTracks->Sumw2();
104 fOutput->Add(fHistoSelTracks);
105 fHistoTracklets = new TH1F("hTracklets","",100,0.,maxMult);
106 fHistoTracklets->Sumw2();
107 fOutput->Add(fHistoTracklets);
108
109 fHistoSPD3DVtxX = new TH1F("hSPD3DvX","",100,-1.,1.);
110 fHistoSPD3DVtxX->Sumw2();
111 fOutput->Add(fHistoSPD3DVtxX);
112 fHistoSPD3DVtxY = new TH1F("hSPD3DvY","",100,-1.,1.);
113 fHistoSPD3DVtxY->Sumw2();
114 fOutput->Add(fHistoSPD3DVtxY);
115 fHistoSPD3DVtxZ = new TH1F("hSPD3DvZ","",100,-15.,15.);
116 fHistoSPD3DVtxZ->Sumw2();
117 fOutput->Add(fHistoSPD3DVtxZ);
118
119 fHistoSPDZVtxX = new TH1F("hSPDZvX","",100,-1.,1.);
120 fHistoSPDZVtxX->Sumw2();
121 fOutput->Add(fHistoSPDZVtxX);
122 fHistoSPDZVtxY = new TH1F("hSPDZvY","",100,-1.,1.);
123 fHistoSPDZVtxY->Sumw2();
124 fOutput->Add(fHistoSPDZVtxY);
125 fHistoSPDZVtxZ = new TH1F("hSPDZvZ","",100,-15.,15.);
126 fHistoSPDZVtxZ->Sumw2();
127 fOutput->Add(fHistoSPDZVtxZ);
128
129
130 fHistoTRKVtxX = new TH1F("hTRKvX","",100,-1.,1.);
131 fHistoTRKVtxX->Sumw2();
132 fOutput->Add(fHistoTRKVtxX);
133 fHistoTRKVtxY = new TH1F("hTRKvY","",100,-1.,1.);
134 fHistoTRKVtxY->Sumw2();
135 fOutput->Add(fHistoTRKVtxY);
136 fHistoTRKVtxZ = new TH1F("hTRKvZ","",100,-15.,15.);
137 fHistoTRKVtxZ->Sumw2();
138 fOutput->Add(fHistoTRKVtxZ);
139
140 Int_t nBinscb=11;
141 if(fPbPb) nBinscb=200;
142 Double_t maxncn=nBinscb-0.5;
143 fHistoNcharmed = new TH2F("hncharmed","",100,0.,maxMult,nBinscb,-0.5,maxncn);
144 fHistoNcharmed->Sumw2();
145 fOutput->Add(fHistoNcharmed);
146 fHistoNbVsNc = new TH2F("hnbvsnc","",nBinscb,-0.5,maxncn,nBinscb,-0.5,maxncn);
147 fHistoNbVsNc->Sumw2();
148 fOutput->Add(fHistoNbVsNc);
149
150 fHistYPtPrompt[0] = new TH2F("hyptd0prompt","D0 - Prompt",20,0.,20.,20,-2.,2.);
151 fHistYPtPrompt[1] = new TH2F("hyptdplusprompt","Dplus - Prompt",20,0.,20.,20,-2.,2.);
152 fHistYPtPrompt[2] = new TH2F("hyptdstarprompt","Dstar - Prompt",20,0.,20.,20,-2.,2.);
153 fHistYPtPrompt[3] = new TH2F("hyptdsprompt","Ds - Prompt",20,0.,20.,20,-2.,2.);
154 fHistYPtPrompt[4] = new TH2F("hyptlcprompt","Lc - Prompt",20,0.,20.,20,-2.,2.);
155
5472ae5b 156 fHistYPtAllDecay[0] = new TH2F("hyptd0AllDecay","D0 - All",20,0.,20.,40,-2.,2.);
157 fHistYPtAllDecay[1] = new TH2F("hyptdplusAllDecay","Dplus - All",20,0.,20.,40,-2.,2.);
158 fHistYPtAllDecay[2] = new TH2F("hyptdstarAllDecay","Dstar - All",20,0.,20.,40,-2.,2.);
159 fHistYPtAllDecay[3] = new TH2F("hyptdsAllDecay","Ds - All",20,0.,20.,40,-2.,2.);
160 fHistYPtAllDecay[4] = new TH2F("hyptlcAllDecay","Lc - All",20,0.,20.,40,-2.,2.);
161
162 fHistYPtPromptAllDecay[0] = new TH2F("hyptd0promptAllDecay","D0 - Prompt",20,0.,20.,40,-2.,2.);
163 fHistYPtPromptAllDecay[1] = new TH2F("hyptdpluspromptAllDecay","Dplus - Prompt",20,0.,20.,40,-2.,2.);
164 fHistYPtPromptAllDecay[2] = new TH2F("hyptdstarpromptAllDecay","Dstar - Prompt",20,0.,20.,40,-2.,2.);
165 fHistYPtPromptAllDecay[3] = new TH2F("hyptdspromptAllDecay","Ds - Prompt",20,0.,20.,40,-2.,2.);
166 fHistYPtPromptAllDecay[4] = new TH2F("hyptlcpromptAllDecay","Lc - Prompt",20,0.,20.,40,-2.,2.);
167
168 fHistYPtFeeddownAllDecay[0] = new TH2F("hyptd0feeddownAllDecay","D0 - FromB",20,0.,20.,40,-2.,2.);
169 fHistYPtFeeddownAllDecay[1] = new TH2F("hyptdplusfeeddownAllDecay","Dplus - FromB",20,0.,20.,40,-2.,2.);
170 fHistYPtFeeddownAllDecay[2] = new TH2F("hyptdstarfeeddownAllDecay","Dstar - FromB",20,0.,20.,40,-2.,2.);
171 fHistYPtFeeddownAllDecay[3] = new TH2F("hyptdsfeeddownAllDecay","Ds - FromB",20,0.,20.,40,-2.,2.);
172 fHistYPtFeeddownAllDecay[4] = new TH2F("hyptlcfeeddownAllDecay","Lc - FromB",20,0.,20.,40,-2.,2.);
173
174
175 fHistYPtFeeddown[0] = new TH2F("hyptd0feeddown","D0 - Feeddown",20,0.,20.,20,-2.,2.);
7b6a4dcd 176 fHistYPtFeeddown[1] = new TH2F("hyptdplusfeeddown","Dplus - Feeddown",20,0.,20.,20,-2.,2.);
177 fHistYPtFeeddown[2] = new TH2F("hyptdstarfeedown","Dstar - Feeddown",20,0.,20.,20,-2.,2.);
178 fHistYPtFeeddown[3] = new TH2F("hyptdsfeedown","Ds - Feeddown",20,0.,20.,20,-2.,2.);
179 fHistYPtFeeddown[4] = new TH2F("hyptlcfeedown","Lc - Feeddown",20,0.,20.,20,-2.,2.);
180
181 for(Int_t ih=0; ih<5; ih++){
5472ae5b 182 fHistYPtAllDecay[ih]->Sumw2();
183 fHistYPtAllDecay[ih]->SetMinimum(0);
184 fOutput->Add(fHistYPtAllDecay[ih]);
185 fHistYPtPromptAllDecay[ih]->Sumw2();
186 fHistYPtPromptAllDecay[ih]->SetMinimum(0);
187 fOutput->Add(fHistYPtPromptAllDecay[ih]);
188 fHistYPtFeeddownAllDecay[ih]->Sumw2();
189 fHistYPtFeeddownAllDecay[ih]->SetMinimum(0);
190 fOutput->Add(fHistYPtFeeddownAllDecay[ih]);
7b6a4dcd 191 fHistYPtPrompt[ih]->Sumw2();
192 fHistYPtPrompt[ih]->SetMinimum(0);
193 fOutput->Add(fHistYPtPrompt[ih]);
194 fHistYPtFeeddown[ih]->Sumw2();
195 fHistYPtFeeddown[ih]->SetMinimum(0);
196 fOutput->Add(fHistYPtFeeddown[ih]);
197 }
198
199 fHistYPtD0byDecChannel[0] = new TH2F("hyptD02","D0 - 2prong",20,0.,20.,20,-2.,2.);
200 fHistYPtD0byDecChannel[1] = new TH2F("hyptD04","D0 - 4prong",20,0.,20.,20,-2.,2.);
201 fHistYPtDplusbyDecChannel[0] = new TH2F("hyptDplusnonreson","Dplus - non reson",20,0.,20.,20,-2.,2.);
202 fHistYPtDplusbyDecChannel[1] = new TH2F("hyptDplusreson","Dplus - reson via K0*",20,0.,20.,20,-2.,2.);
203 fHistYPtDsbyDecChannel[0] = new TH2F("hyptdsphi","Ds - vis Phi",20,0.,20.,20,-2.,2.);
204 fHistYPtDsbyDecChannel[1] = new TH2F("hyptdsk0st","Ds - via k0*",20,0.,20.,20,-2.,2.);
205
206 for(Int_t ih=0; ih<2; ih++){
207
208 fHistYPtD0byDecChannel[ih]->Sumw2();
209 fHistYPtD0byDecChannel[ih]->SetMinimum(0);
210 fOutput->Add(fHistYPtD0byDecChannel[ih]);
211 fHistYPtDplusbyDecChannel[ih]->Sumw2();
212 fHistYPtDplusbyDecChannel[ih]->SetMinimum(0);
213 fOutput->Add(fHistYPtDplusbyDecChannel[ih]);
214 fHistYPtDsbyDecChannel[ih]->Sumw2();
215 fHistYPtDsbyDecChannel[ih]->SetMinimum(0);
216 fOutput->Add(fHistYPtDsbyDecChannel[ih]);
217 }
218
219
220 PostData(1,fOutput);
221
222}
223//______________________________________________________________________________
224void AliAnalysisTaskCheckHFMCProd::UserExec(Option_t *)
225{
226 //
227
228 AliESDEvent *esd = (AliESDEvent*) (InputEvent());
229
230
231 if(!esd) {
232 printf("AliAnalysisTaskSDDRP::Exec(): bad ESD\n");
233 return;
234 }
235
236 fHistoNEvents->Fill(0);
237
238 Int_t nTracks=esd->GetNumberOfTracks();
239 fHistoTracks->Fill(nTracks);
240 Int_t nSelTracks=0;
241 for(Int_t it=0; it<nTracks; it++){
242 AliESDtrack* tr=esd->GetTrack(it);
243 UInt_t status=tr->GetStatus();
244 if(!(status&AliESDtrack::kITSrefit)) continue;
245 if(!(status&AliESDtrack::kTPCin)) continue;
246 nSelTracks++;
247 }
248 fHistoSelTracks->Fill(nSelTracks);
249
250 const AliMultiplicity* mult=esd->GetMultiplicity();
251 Int_t nTracklets=mult->GetNumberOfTracklets();
252 fHistoTracklets->Fill(nTracklets);
253
254 const AliESDVertex *spdv=esd->GetVertex();
255 if(spdv && spdv->IsFromVertexer3D()){
256 fHistoSPD3DVtxX->Fill(spdv->GetXv());
257 fHistoSPD3DVtxY->Fill(spdv->GetYv());
258 fHistoSPD3DVtxZ->Fill(spdv->GetZv());
259 }
260 if(spdv && spdv->IsFromVertexerZ()){
261 fHistoSPDZVtxX->Fill(spdv->GetXv());
262 fHistoSPDZVtxY->Fill(spdv->GetYv());
263 fHistoSPDZVtxZ->Fill(spdv->GetZv());
264 }
265 const AliESDVertex *trkv=esd->GetPrimaryVertex();
266 if(trkv && trkv->GetNContributors()>1){
267 fHistoTRKVtxX->Fill(trkv->GetXv());
268 fHistoTRKVtxY->Fill(trkv->GetYv());
269 fHistoTRKVtxZ->Fill(trkv->GetZv());
270 }
271
272 AliStack* stack=0;
273
274 if(fReadMC){
275 AliMCEventHandler* eventHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
276 if (!eventHandler) {
277 Printf("ERROR: Could not retrieve MC event handler");
278 return;
279 }
280 AliMCEvent* mcEvent = eventHandler->MCEvent();
281 if (!mcEvent) {
282 Printf("ERROR: Could not retrieve MC event");
283 return;
284 }
285 stack = mcEvent->Stack();
286 if (!stack) {
287 Printf("ERROR: stack not available");
288 return;
289 }
290
291
292 Int_t nParticles=stack->GetNtrack();
293 Double_t dNchdy = 0.;
294 Int_t nb = 0, nc=0;
a0dabf3d 295 Int_t nCharmed=0;
7b6a4dcd 296 for (Int_t i=0;i<nParticles;i++){
297 TParticle* part = (TParticle*)stack->Particle(i);
298 Int_t absPdg=TMath::Abs(part->GetPdgCode());
299 if(absPdg==4) nc++;
300 if(absPdg==5) nb++;
301 if(stack->IsPhysicalPrimary(i)){
302 Double_t eta=part->Eta();
303 if(TMath::Abs(eta)<0.5) dNchdy+=0.6666; // 2/3 for the ratio charged/all
304 }
5f7f126e 305 Float_t rapid=-999.;
306 if (part->Energy() != TMath::Abs(part->Pz())){
307 rapid=0.5*TMath::Log((part->Energy()+part->Pz())/(part->Energy()-part->Pz()));
308 }
5472ae5b 309
310
7b6a4dcd 311 Int_t iPart=-1;
312 Int_t iType=0;
5472ae5b 313 Int_t iSpecies=-1;
314 if(absPdg==421){
315 iSpecies=0;
7b6a4dcd 316 iType=CheckD0Decay(i,stack);
317 if(iType>=0) iPart=0;
318 }
319 else if(absPdg==411){
5472ae5b 320 iSpecies=1;
7b6a4dcd 321 iType=CheckDplusDecay(i,stack);
322 if(iType>=0) iPart=1;
323 }
324 else if(absPdg==413){
5472ae5b 325 iSpecies=2;
7b6a4dcd 326 iType=CheckDstarDecay(i,stack);
327 if(iType>=0) iPart=2;
328 }
329 else if(absPdg==431){
5472ae5b 330 iSpecies=3;
7b6a4dcd 331 iType=CheckDsDecay(i,stack);
332 if(iType==0 || iType==1) iPart=3;
333 }
334 else if(absPdg==4122){
5472ae5b 335 iSpecies=4;
7b6a4dcd 336 iType=CheckLcDecay(i,stack);
337 if(iType>=0) iPart=4;
338 }
5472ae5b 339 if(iSpecies<0) continue;
340 fHistYPtAllDecay[iSpecies]->Fill(part->Pt(),rapid);
341
7b6a4dcd 342 TParticle* runningpart=part;
343 Int_t iFromB=-1;
5472ae5b 344 Int_t pdgmoth=-1;
7b6a4dcd 345 while(1){
346 Int_t labmoth=runningpart->GetFirstMother();
347 if(labmoth==-1) break;
348 TParticle *mot=(TParticle*)stack->Particle(labmoth);
5472ae5b 349 pdgmoth=TMath::Abs(mot->GetPdgCode());
7b6a4dcd 350 if(pdgmoth==5){
351 iFromB=1;
352 break;
353 }else if(pdgmoth==4){
354 iFromB=0;
5472ae5b 355 break;
7b6a4dcd 356 }
357 runningpart=mot;
358 }
359 if(iFromB<0) continue;
5472ae5b 360 else if(iFromB==0) fHistYPtPromptAllDecay[iSpecies]->Fill(part->Pt(),rapid);
361 else if(iFromB==1) fHistYPtFeeddownAllDecay[iSpecies]->Fill(part->Pt(),rapid);
362
363 if(iPart<0) continue;
364 if(iType<0) continue;
365 nCharmed++;
366 if(iPart==0 && iType<=1){
367 fHistYPtD0byDecChannel[iType]->Fill(part->Pt(),rapid);
368 }else if(iPart==1 && iType<=1){
369 fHistYPtDplusbyDecChannel[iType]->Fill(part->Pt(),rapid);
370 }else if(iPart==3 && iType<=1){
371 fHistYPtDsbyDecChannel[iType]->Fill(part->Pt(),rapid);
372 }
373
7b6a4dcd 374 if(iFromB==0 && iPart>=0 && iPart<5) fHistYPtPrompt[iPart]->Fill(part->Pt(),rapid);
375 else if(iFromB==1 && iPart>=0 && iPart<5) fHistYPtFeeddown[iPart]->Fill(part->Pt(),rapid);
376
377 }
7b6a4dcd 378 fHistoNcharmed->Fill(dNchdy,nCharmed);
379 fHistoNbVsNc->Fill(nc,nb);
380 }
381
382 PostData(1,fOutput);
383
384}
385//______________________________________________________________________________
386void AliAnalysisTaskCheckHFMCProd::Terminate(Option_t */*option*/)
387{
388 // Terminate analysis
389 fOutput = dynamic_cast<TList*> (GetOutputData(1));
390 if (!fOutput) {
391 printf("ERROR: fOutput not available\n");
392 return;
393 }
394
395 return;
396}
397
398
399
400
401//______________________________________________________________________________
402Int_t AliAnalysisTaskCheckHFMCProd::CheckD0Decay(Int_t labD0, AliStack* stack) const{
403
404 if(labD0<0) return -1;
405 TParticle* dp = (TParticle*)stack->Particle(labD0);
406 Int_t pdgdp=dp->GetPdgCode();
407 Int_t nDau=dp->GetNDaughters();
7b6a4dcd 408
409 if(nDau==2){
410 Int_t nKaons=0;
411 Int_t nPions=0;
412 for(Int_t iDau=dp->GetFirstDaughter(); iDau<=dp->GetLastDaughter(); iDau++){
413 if(iDau<0) return -1;
414 TParticle* dau=(TParticle*)stack->Particle(iDau);
415 Int_t pdgdau=dau->GetPdgCode();
416 if(TMath::Abs(pdgdau)==321){
417 if(pdgdp>0 && pdgdau>0) return -1;
418 if(pdgdp<0 && pdgdau<0) return -1;
419 nKaons++;
420 }else if(TMath::Abs(pdgdau)==211){
421 if(pdgdp<0 && pdgdau>0) return -1;
422 if(pdgdp>0 && pdgdau<0) return -1;
423 nPions++;
424 }else{
425 return -1;
426 }
427 }
428 if(nPions!=1) return -1;
429 if(nKaons!=1) return -1;
430 for(Int_t iDau=dp->GetFirstDaughter(); iDau<=dp->GetLastDaughter(); iDau++){
431 if(iDau<0) return -1;
432 }
433 return 0;
434 }
435
436 if(nDau==3 || nDau==4){
437 Int_t nKaons=0;
438 Int_t nPions=0;
439 for(Int_t iDau=dp->GetFirstDaughter(); iDau<=dp->GetLastDaughter(); iDau++){
440 if(iDau<0) return -1;
441 TParticle* dau=(TParticle*)stack->Particle(iDau);
442 Int_t pdgdau=dau->GetPdgCode();
443 if(TMath::Abs(pdgdau)==321){
444 if(pdgdp>0 && pdgdau>0) return -1;
445 if(pdgdp<0 && pdgdau<0) return -1;
446 nKaons++;
447 }else if(TMath::Abs(pdgdau)==113 || TMath::Abs(pdgdau)==313){
448 for(Int_t resDau=dau->GetFirstDaughter(); resDau<=dau->GetLastDaughter(); resDau++){
449 if(resDau<0) return -1;
450 TParticle* resdau=(TParticle*)stack->Particle(resDau);
451 Int_t pdgresdau=resdau->GetPdgCode();
452 if(TMath::Abs(pdgresdau)==321){
453 if(pdgdp>0 && pdgresdau>0) return -1;
454 if(pdgdp<0 && pdgresdau<0) return -1;
455 nKaons++;
456 }
457 if(TMath::Abs(pdgresdau)==211){
458 nPions++;
459 }
460 }
461 }else if(TMath::Abs(pdgdau)==211){
462 nPions++;
463 }else{
464 return -1;
465 }
466 }
467 if(nPions!=3) return -1;
468 if(nKaons!=1) return -1;
469 for(Int_t iDau=dp->GetFirstDaughter(); iDau<=dp->GetLastDaughter(); iDau++){
470 if(iDau<0) return -1;
471 }
472 return 1;
473 }
474
475 return -1;
476}
477
478
479//______________________________________________________________________________
480Int_t AliAnalysisTaskCheckHFMCProd::CheckDplusDecay(Int_t labDplus, AliStack* stack) const{
481
482 if(labDplus<0) return -1;
483 TParticle* dp = (TParticle*)stack->Particle(labDplus);
484 Int_t pdgdp=dp->GetPdgCode();
485 Int_t nDau=dp->GetNDaughters();
7b6a4dcd 486
487 if(nDau==3){
488 Int_t nKaons=0;
489 Int_t nPions=0;
490 for(Int_t iDau=dp->GetFirstDaughter(); iDau<=dp->GetLastDaughter(); iDau++){
491 if(iDau<0) return -1;
492 TParticle* dau=(TParticle*)stack->Particle(iDau);
493 Int_t pdgdau=dau->GetPdgCode();
494 if(TMath::Abs(pdgdau)==321){
495 if(pdgdp>0 && pdgdau>0) return -1;
496 if(pdgdp<0 && pdgdau<0) return -1;
497 nKaons++;
498 }else if(TMath::Abs(pdgdau)==211){
499 if(pdgdp<0 && pdgdau>0) return -1;
500 if(pdgdp>0 && pdgdau<0) return -1;
501 nPions++;
502 }else{
503 return -1;
504 }
505 }
506 if(nPions!=2) return -1;
507 if(nKaons!=1) return -1;
508 for(Int_t iDau=dp->GetFirstDaughter(); iDau<=dp->GetLastDaughter(); iDau++){
509 if(iDau<0) return -1;
510 }
511 return 0;
512 }
513
514 if(nDau==2){
515 Int_t nKaons=0;
516 Int_t nPions=0;
517 for(Int_t iDau=dp->GetFirstDaughter(); iDau<=dp->GetLastDaughter(); iDau++){
518 if(iDau<0) return -1;
519 TParticle* dau=(TParticle*)stack->Particle(iDau);
520 Int_t pdgdau=dau->GetPdgCode();
521 if(TMath::Abs(pdgdau)==313){
522 for(Int_t resDau=dau->GetFirstDaughter(); resDau<=dau->GetLastDaughter(); resDau++){
523 if(resDau<0) return -1;
524 TParticle* resdau=(TParticle*)stack->Particle(resDau);
525 Int_t pdgresdau=resdau->GetPdgCode();
526 if(TMath::Abs(pdgresdau)==321){
527 if(pdgdp>0 && pdgresdau>0) return -1;
528 if(pdgdp<0 && pdgresdau<0) return -1;
529 nKaons++;
530 }
531 if(TMath::Abs(pdgresdau)==211){
532 if(pdgdp<0 && pdgresdau>0) return -1;
533 if(pdgdp>0 && pdgresdau<0) return -1;
534 nPions++;
535 }
536 }
537 }else if(TMath::Abs(pdgdau)==211){
538 if(pdgdp<0 && pdgdau>0) return -1;
539 if(pdgdp>0 && pdgdau<0) return -1;
540 nPions++;
541 }else{
542 return -1;
543 }
544 }
545 if(nPions!=2) return -1;
546 if(nKaons!=1) return -1;
547 for(Int_t iDau=dp->GetFirstDaughter(); iDau<=dp->GetLastDaughter(); iDau++){
548 if(iDau<0) return -1;
549 }
550 return 1;
551 }
552 return -1;
553}
554
555//______________________________________________________________________________
556Int_t AliAnalysisTaskCheckHFMCProd::CheckDsDecay(Int_t labDs, AliStack* stack) const{
557 // Ds decay
558 if(labDs<0) return -1;
559 TParticle* dp = (TParticle*)stack->Particle(labDs);
560 Int_t pdgdp=dp->GetPdgCode();
561 Int_t nDau=dp->GetNDaughters();
7b6a4dcd 562
563 if(nDau==3){
564 Int_t nKaons=0;
565 Int_t nPions=0;
566 for(Int_t iDau=dp->GetFirstDaughter(); iDau<=dp->GetLastDaughter(); iDau++){
567 if(iDau<0) return -1;
568 TParticle* dau=(TParticle*)stack->Particle(iDau);
569 Int_t pdgdau=dau->GetPdgCode();
570 if(TMath::Abs(pdgdau)==321){
571 nKaons++;
572 }else if(TMath::Abs(pdgdau)==211){
573 if(pdgdp<0 && pdgdau>0) return -1;
574 if(pdgdp>0 && pdgdau<0) return -1;
575 nPions++;
576 }else{
577 return -1;
578 }
579 }
580 if(nPions!=1) return -1;
581 if(nKaons!=2) return -1;
582 for(Int_t iDau=dp->GetFirstDaughter(); iDau<=dp->GetLastDaughter(); iDau++){
583 if(iDau<0) return -1;
584 }
585 return 2;
586 }
587
588 if(nDau==2){
589 Int_t nKaons=0;
590 Int_t nPions=0;
591 Bool_t isPhi=kFALSE;
592 Bool_t isk0st=kFALSE;
593 for(Int_t iDau=dp->GetFirstDaughter(); iDau<=dp->GetLastDaughter(); iDau++){
594 if(iDau<0) return -1;
595 TParticle* dau=(TParticle*)stack->Particle(iDau);
596 Int_t pdgdau=dau->GetPdgCode();
597 if(TMath::Abs(pdgdau)==313){
598 isk0st=kTRUE;
599 for(Int_t resDau=dau->GetFirstDaughter(); resDau<=dau->GetLastDaughter(); resDau++){
600 if(resDau<0) return -1;
601 TParticle* resdau=(TParticle*)stack->Particle(resDau);
602 Int_t pdgresdau=resdau->GetPdgCode();
603 if(TMath::Abs(pdgresdau)==321){
604 nKaons++;
605 }
606 if(TMath::Abs(pdgresdau)==211){
607 if(pdgdp<0 && pdgresdau>0) return -1;
608 if(pdgdp>0 && pdgresdau<0) return -1;
609 nPions++;
610 }
611 }
612 }else if(TMath::Abs(pdgdau)==333){
613 isPhi=kTRUE;
614 for(Int_t resDau=dau->GetFirstDaughter(); resDau<=dau->GetLastDaughter(); resDau++){
615 if(resDau<0) return -1;
616 TParticle* resdau=(TParticle*)stack->Particle(resDau);
48709ec3 617 if(!resdau) return -1;
7b6a4dcd 618 Int_t pdgresdau=resdau->GetPdgCode();
619 if(TMath::Abs(pdgresdau)==321){
620 nKaons++;
621 }else{
622 return -1;
623 }
624 }
625 }else if(TMath::Abs(pdgdau)==211){
626 if(pdgdp<0 && pdgdau>0) return -1;
627 if(pdgdp>0 && pdgdau<0) return -1;
628 nPions++;
629 }else if(TMath::Abs(pdgdau)==321){
630 nKaons++;
631 }else{
632 return -1;
633 }
634 }
635 if(nPions!=1) return -1;
636 if(nKaons!=2) return -1;
637 for(Int_t iDau=dp->GetFirstDaughter(); iDau<=dp->GetLastDaughter(); iDau++){
638 if(iDau<0) return -1;
639 }
640 if(isk0st) return 1;
641 else if(isPhi) return 0;
642 else return 3;
643 }
644 return -1;
645}
646
647//______________________________________________________________________________
648Int_t AliAnalysisTaskCheckHFMCProd::CheckDstarDecay(Int_t labDstar, AliStack* stack) const{
649
650 if(labDstar<0) return -1;
651 TParticle* dp = (TParticle*)stack->Particle(labDstar);
652 Int_t pdgdp=dp->GetPdgCode();
653 Int_t nDau=dp->GetNDaughters();
654 if(nDau!=2) return -1;
655
656 Int_t nKaons=0;
657 Int_t nPions=0;
658 for(Int_t iDau=dp->GetFirstDaughter(); iDau<=dp->GetLastDaughter(); iDau++){
659 if(iDau<0) return -1;
660 TParticle* dau=(TParticle*)stack->Particle(iDau);
661 Int_t pdgdau=dau->GetPdgCode();
662 if(TMath::Abs(pdgdau)==421){
663 for(Int_t resDau=dau->GetFirstDaughter(); resDau<=dau->GetLastDaughter(); resDau++){
664 if(resDau<0) return -1;
665 TParticle* resdau=(TParticle*)stack->Particle(resDau);
666 Int_t pdgresdau=resdau->GetPdgCode();
667 if(TMath::Abs(pdgresdau)==321){
668 if(pdgdp>0 && pdgresdau>0) return -1;
669 if(pdgdp<0 && pdgresdau<0) return -1;
670 nKaons++;
671 }
672 if(TMath::Abs(pdgresdau)==211){
673 if(pdgdp<0 && pdgresdau>0) return -1;
674 if(pdgdp>0 && pdgresdau<0) return -1;
675 nPions++;
676 }
677 }
678 }else if(TMath::Abs(pdgdau)==211){
679 if(pdgdp<0 && pdgdau>0) return -1;
680 if(pdgdp>0 && pdgdau<0) return -1;
681 nPions++;
682 }else{
683 return -1;
684 }
685 }
686 if(nPions!=2) return -1;
687 if(nKaons!=1) return -1;
688 for(Int_t iDau=dp->GetFirstDaughter(); iDau<=dp->GetLastDaughter(); iDau++){
689 if(iDau<0) return -1;
690 }
691 return 0;
692
693}
694
695//______________________________________________________________________________
696Int_t AliAnalysisTaskCheckHFMCProd::CheckLcDecay(Int_t labLc, AliStack* stack) const{
697 if(labLc<0) return -1;
698 TParticle* dp = (TParticle*)stack->Particle(labLc);
699 Int_t pdgdp=dp->GetPdgCode();
700 Int_t nDau=dp->GetNDaughters();
7b6a4dcd 701
702 if(nDau==3){
703 Int_t nKaons=0;
704 Int_t nPions=0;
705 Int_t nProtons=0;
706 for(Int_t iDau=dp->GetFirstDaughter(); iDau<=dp->GetLastDaughter(); iDau++){
707 if(iDau<0) return -1;
708 TParticle* dau=(TParticle*)stack->Particle(iDau);
709 Int_t pdgdau=dau->GetPdgCode();
710 if(TMath::Abs(pdgdau)==321){
711 if(pdgdp>0 && pdgdau>0) return -1;
712 if(pdgdp<0 && pdgdau<0) return -1;
713 nKaons++;
714 }else if(TMath::Abs(pdgdau)==211){
715 if(pdgdp<0 && pdgdau>0) return -1;
716 if(pdgdp>0 && pdgdau<0) return -1;
717 nPions++;
718 }else if(TMath::Abs(pdgdau)==2212){
719 if(pdgdp<0 && pdgdau>0) return -1;
720 if(pdgdp>0 && pdgdau<0) return -1;
721 nProtons++;
722 }else{
723 return -1;
724 }
725 }
726 if(nPions!=1) return -1;
727 if(nKaons!=1) return -1;
728 if(nProtons!=1) return -1;
729 for(Int_t iDau=dp->GetFirstDaughter(); iDau<=dp->GetLastDaughter(); iDau++){
730 if(iDau<0) return -1;
731 }
732 return 0;
733 }
734
735 if(nDau==2){
736 Int_t nKaons=0;
737 Int_t nPions=0;
738 Int_t nProtons=0;
739 for(Int_t iDau=dp->GetFirstDaughter(); iDau<=dp->GetLastDaughter(); iDau++){
740 if(iDau<0) return -1;
741 TParticle* dau=(TParticle*)stack->Particle(iDau);
742 Int_t pdgdau=dau->GetPdgCode();
743 if(TMath::Abs(pdgdau)==313 || TMath::Abs(pdgdau)==3124 || TMath::Abs(pdgdau)==2224
744 || TMath::Abs(pdgdau)==3122 || TMath::Abs(pdgdau)==311){
745 Int_t nDauRes=dau->GetNDaughters();
746 if(nDauRes!=2) return -1;
747 for(Int_t resDau=dau->GetFirstDaughter(); resDau<=dau->GetLastDaughter(); resDau++){
748 if(resDau<0) return -1;
749 TParticle* resdau=(TParticle*)stack->Particle(resDau);
750 Int_t pdgresdau=resdau->GetPdgCode();
751 if(TMath::Abs(pdgresdau)==321){
752 if(pdgdp>0 && pdgresdau>0) return -1;
753 if(pdgdp<0 && pdgresdau<0) return -1;
754 nKaons++;
755 }
756 if(TMath::Abs(pdgresdau)==211){
757 if(pdgdp<0 && pdgresdau>0) return -1;
758 if(pdgdp>0 && pdgresdau<0) return -1;
759 nPions++;
760 }
761 if(TMath::Abs(pdgresdau)==2212){
762 if(pdgdp<0 && pdgresdau>0) return -1;
763 if(pdgdp>0 && pdgresdau<0) return -1;
764 nProtons++;
765 }
766 }
767 }else if(TMath::Abs(pdgdau)==321){
768 if(pdgdp>0 && pdgdau>0) return -1;
769 if(pdgdp<0 && pdgdau<0) return -1;
770 nKaons++;
771 }else if(TMath::Abs(pdgdau)==211){
772 if(pdgdp<0 && pdgdau>0) return -1;
773 if(pdgdp>0 && pdgdau<0) return -1;
774 nPions++;
775 }else if(TMath::Abs(pdgdau)==2212){
776 if(pdgdp<0 && pdgdau>0) return -1;
777 if(pdgdp>0 && pdgdau<0) return -1;
778 nProtons++;
779 }else{
780 return -1;
781 }
782 }
783 if(nPions!=1) return -1;
784 if(nKaons!=1) return -1;
785 if(nProtons!=1) return -1;
786 for(Int_t iDau=dp->GetFirstDaughter(); iDau<=dp->GetLastDaughter(); iDau++){
787 if(iDau<0) return -1;
788 }
789 return 1;
790 }
791 return -1;
792}
793