1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Authors: Friederike Bock *
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 **************************************************************************/
16 ////////////////////////////////////////////////
17 //---------------------------------------------
18 // QA Task for V0 Reader V1
19 //---------------------------------------------
20 ////////////////////////////////////////////////
22 #include "AliAnalysisTaskMaterial.h"
24 #include "AliAnalysisManager.h"
25 #include "TParticle.h"
27 #include "AliPIDResponse.h"
28 #include "AliESDtrackCuts.h"
35 ClassImp(AliAnalysisTaskMaterial)
37 AliAnalysisTaskMaterial::AliAnalysisTaskMaterial() : AliAnalysisTaskSE(),
39 fConversionGammas(NULL),
40 fConversionCuts(NULL),
45 fAllMCGammaList(NULL),
46 fAllMCConvGammaList(NULL),
48 fTreeMaterialRec(NULL),
49 fTreeMaterialAllGamma(NULL),
50 fTreeMaterialConvGamma(NULL),
54 fNESDtracksEta0914(0),
59 fGammaMCConvTheta(0.),
61 fMCConvDaughterProp(4),
77 //________________________________________________________________________
78 AliAnalysisTaskMaterial::AliAnalysisTaskMaterial(const char *name) : AliAnalysisTaskSE(name),
80 fConversionGammas(NULL),
81 fConversionCuts(NULL),
86 fAllMCGammaList(NULL),
87 fAllMCConvGammaList(NULL),
89 fTreeMaterialRec(NULL),
90 fTreeMaterialAllGamma(NULL),
91 fTreeMaterialConvGamma(NULL),
95 fNESDtracksEta0914(0),
100 fGammaMCConvTheta(0.),
102 fMCConvDaughterProp(4),
114 // Default constructor
116 DefineInput(0, TChain::Class());
117 DefineOutput(1, TList::Class());
120 //________________________________________________________________________
121 AliAnalysisTaskMaterial::~AliAnalysisTaskMaterial()
123 // default deconstructor
125 //________________________________________________________________________
126 void AliAnalysisTaskMaterial::UserCreateOutputObjects()
128 // Create User Output Objects
130 if(fOutputList != NULL){
134 if(fOutputList == NULL){
135 fOutputList = new TList();
136 fOutputList->SetOwner(kTRUE);
139 fEventList = new TList();
140 fEventList->SetName("EventList");
141 fEventList->SetOwner(kTRUE);
142 fOutputList->Add(fEventList);
144 fTreeEvent = new TTree("Event","Event");
145 fTreeEvent->Branch("primVtxZ",&fPrimVtxZ,"fPrimVtxZ/F");
146 fTreeEvent->Branch("nContrVtx",&fNContrVtx,"fNContrVtx/I");
147 fTreeEvent->Branch("nGoodTracksEta09",&fNESDtracksEta09,"fNESDtracksEta09/I");
148 fTreeEvent->Branch("nGoodTracksEta14",&fNESDtracksEta14,"fNESDtracksEta14/I");
149 fEventList->Add(fTreeEvent);
151 fRecGammaList= new TList();
152 fRecGammaList->SetName("RecGammaList");
153 fRecGammaList->SetOwner(kTRUE);
154 fOutputList->Add(fRecGammaList);
156 fTreeMaterialRec = new TTree("ConvPointRec","ConvPointRec");
157 fTreeMaterialRec->Branch("nGoodTracksEta09",&fNESDtracksEta09,"fNESDtracksEta09/I");
158 fTreeMaterialRec->Branch("recCords",&fRecCords);
159 fTreeMaterialRec->Branch("daughterProp",&fDaughterProp);
160 fTreeMaterialRec->Branch("pt",&fGammaPt,"fGammaPt/F");
161 fTreeMaterialRec->Branch("theta",&fGammaTheta,"fGammaTheta/F");
162 fTreeMaterialRec->Branch("chi2ndf",&fGammaChi2NDF,"fGammaChi2NDF/F");
164 fTreeMaterialRec->Branch("kind",&fKind,"fKind/b");
166 fRecGammaList->Add(fTreeMaterialRec);
169 // fAllMCGammaList = new TList();
170 // fAllMCGammaList->SetName("AllMCGammaList");
171 // fAllMCGammaList->SetOwner(kTRUE);
172 // fOutputList->Add(fAllMCGammaList);
174 // fTreeMaterialAllGamma = new TTree("AllGamma","AllGamma");
175 // fTreeMaterialAllGamma->Branch("pt",&fGammaMCPt,"fGammaMCPt/F");
176 // fTreeMaterialAllGamma->Branch("theta",&fGammaMCTheta,"fGammaMCTheta/F");
177 // fAllMCGammaList->Add(fTreeMaterialAllGamma);
179 // fAllMCConvGammaList = new TList();
180 // fAllMCConvGammaList->SetName("AllMCGammaConvList");
181 // fAllMCConvGammaList->SetOwner(kTRUE);
182 // fOutputList->Add(fAllMCConvGammaList);
184 // // fMCConvCords = new Float_t[5];
185 // // fMCConvDaughterProp = new Float_t[4];
188 // fTreeMaterialConvGamma = new TTree("ConvGammaMC","ConvGammaMC");
189 // fTreeMaterialConvGamma->Branch("Cords",&fMCConvCords);
190 // fTreeMaterialConvGamma->Branch("daughterProp",&fMCConvDaughterProp);
191 // fTreeMaterialConvGamma->Branch("Pt",&fGammaMCConvPt,"fGammaMCConvPt/F");
192 // fTreeMaterialConvGamma->Branch("Theta",&fGammaMCConvTheta,"fGammaMCConvTheta/F");
193 // fAllMCConvGammaList->Add(fTreeMaterialConvGamma);
196 PostData(1, fOutputList);
200 //________________________________________________________________________
201 void AliAnalysisTaskMaterial::UserExec(Option_t *){
203 fV0Reader=(AliV0ReaderV1*)AliAnalysisManager::GetAnalysisManager()->GetTask("V0ReaderV1");
205 Int_t eventQuality = ((AliConvEventCuts*)fV0Reader->GetEventCuts())->GetEventQuality();
206 if(eventQuality != 0){// Event Not Accepted
209 fESDEvent = (AliESDEvent*) InputEvent();
210 if (fESDEvent==NULL) return;
212 fMCEvent = MCEvent();
216 // Process MC Particle
217 if(fEventCuts->GetSignalRejection() != 0){
218 // if(fESDEvent->IsA()==AliESDEvent::Class()){
219 fEventCuts->GetNotRejectedParticles( fEventCuts->GetSignalRejection(),
220 fEventCuts->GetAcceptedHeader(),
223 // else if(fInputEvent->IsA()==AliAODEvent::Class()){
224 // fEventCuts->GetNotRejectedParticles( fEventCuts->GetSignalRejection(),
225 // fEventCuts->GetAcceptedHeader(),
231 if(fIsHeavyIon > 0 && !fEventCuts->IsCentralitySelected(fESDEvent)) return;
232 fNESDtracksEta09 = CountTracks09(); // Estimate Event Multiplicity
233 fNESDtracksEta0914 = CountTracks0914(); // Estimate Event Multiplicity
234 fNESDtracksEta14 = fNESDtracksEta09 + fNESDtracksEta0914;
236 if(fESDEvent->GetPrimaryVertexTracks()->GetNContributors()>0) {
237 fNContrVtx = fESDEvent->GetPrimaryVertexTracks()->GetNContributors();
241 // else if(fESDEvent->GetPrimaryVertexTracks()->GetNContributors()<1) {
242 // if(fESDEvent->GetPrimaryVertexSPD()->GetNContributors()>0) {
243 // fNContrVtx = fESDEvent->GetPrimaryVertexSPD()->GetNContributors();
244 // } else if(fESDEvent->GetPrimaryVertexSPD()->GetNContributors()<1) {
249 fPrimVtxZ = fESDEvent->GetPrimaryVertex()->GetZ();
251 // if (fIsHeavyIon == 2){
252 // if (!(fNESDtracksEta09 > 20 && fNESDtracksEta09 < 80)) return;
260 fConversionGammas=fV0Reader->GetReconstructedGammas();
263 // ProcessMCPhotons();
265 PostData(1, fOutputList);
268 // ///________________________________________________________________________
269 // void AliAnalysisTaskMaterial::FillMCTree(Int_t stackPos){
270 // AliStack *MCStack = fMCEvent->Stack();
271 // TParticle* candidate = (TParticle *)MCStack->Particle(stackPos);
273 // if(fConversionCuts->PhotonIsSelectedMC(candidate,MCStack,kFALSE)){
274 // fGammaMCPt = candidate->Pt();
275 // fGammaMCTheta = candidate->Theta();
276 // if (fTreeMaterialAllGamma){
277 // fTreeMaterialAllGamma->Fill();
280 // if(fConversionCuts->PhotonIsSelectedMC(candidate,MCStack,kTRUE)){
281 // fGammaMCConvPt = candidate->Pt();
282 // fGammaMCConvTheta = candidate->Theta();
283 // TParticle* daughter1 = (TParticle *)MCStack->Particle(candidate->GetFirstDaughter());
284 // TParticle* daughter2 = (TParticle *)MCStack->Particle(candidate->GetLastDaughter());
285 // fMCConvCords(0) = (Float_t)daughter1->Vx();
286 // fMCConvCords(1) = (Float_t)daughter1->Vy();
287 // fMCConvCords(2) = (Float_t)daughter1->Vz();
288 // fMCConvCords(3) = (Float_t)daughter1->R();
289 // fMCConvCords(4) = (Float_t)daughter1->Phi();
291 // fMCConvDaughterProp(0) = (Float_t)daughter1->Pt();
292 // fMCConvDaughterProp(1) = (Float_t)daughter1->Theta();
293 // fMCConvDaughterProp(2) = (Float_t)daughter2->Pt();
294 // fMCConvDaughterProp(3) = (Float_t)daughter2->Theta();
296 // if (fTreeMaterialConvGamma){
297 // fTreeMaterialConvGamma->Fill();
299 // } // Converted MC Gamma
302 // ///________________________________________________________________________
303 // void AliAnalysisTaskMaterial::ProcessMCPhotons(){
304 // // Loop over all primary MC particle
305 // AliStack *ffMCStack = fMCEvent->Stack();
306 // for(Int_t i = 0; i < ffMCStack->GetNprimary(); i++) {
307 // TParticle* particle = (TParticle *)ffMCStack->Particle(i);
308 // if (!particle) continue;
309 // Int_t isMCFromMBHeader = -1;
310 // if(fConversionCuts->GetSignalRejection() != 0){
312 // = fConversionCuts->IsParticleFromBGEvent(i, ffMCStack, fESDEvent);
313 // if(isMCFromMBHeader == 0) continue;
315 // if (particle->GetPdgCode() == 111 && particle->GetFirstDaughter() >= ffMCStack->GetNprimary()){
316 // // cout << "Undecayed pi0 found with mother: " << particle->GetMother(0) << endl;
317 // for (Int_t j = 0; j < 2 ; j++){
318 // FillMCTree(particle->GetDaughter(j));
326 ///________________________________________________________________________
327 void AliAnalysisTaskMaterial::ProcessPhotons(){
329 // Fill Histograms for QA and MC
330 for(Int_t firstGammaIndex=0;firstGammaIndex<fConversionGammas->GetEntriesFast();firstGammaIndex++){
331 AliAODConversionPhoton *gamma=dynamic_cast<AliAODConversionPhoton*>(fConversionGammas->At(firstGammaIndex));
332 if (gamma ==NULL) continue;
333 if(!fConversionCuts->PhotonIsSelected(gamma,fESDEvent)) continue;
335 fNESDtracksEta09 = CountTracks09();
336 fGammaPt = gamma->GetPhotonPt();
337 fGammaTheta = gamma->GetPhotonTheta();
338 fGammaChi2NDF = gamma->GetChi2perNDF();
339 fRecCords(0) = (Float_t)gamma->GetConversionX();
340 fRecCords(1) = (Float_t)gamma->GetConversionY();
341 fRecCords(2) = (Float_t)gamma->GetConversionZ();
342 fRecCords(3) = (Float_t)gamma->GetConversionRadius();
343 fRecCords(4) = (Float_t)gamma->GetPhotonPhi();
345 AliESDtrack * negTrack = fConversionCuts->GetESDTrack(fESDEvent, gamma->GetTrackLabelNegative());
346 AliESDtrack * posTrack = fConversionCuts->GetESDTrack(fESDEvent, gamma->GetTrackLabelPositive());
347 fDaughterProp(0) = (Float_t)posTrack->Pt();
348 fDaughterProp(1) = (Float_t)posTrack->Theta();
349 fDaughterProp(2) = (Float_t)negTrack->Pt();
350 fDaughterProp(3) = (Float_t)negTrack->Theta();
355 // cout << "generating MC stack"<< endl;
356 AliStack *fMCStack = fMCEvent->Stack();
357 if (!fMCStack) continue;
358 TParticle *posDaughter = gamma->GetPositiveMCDaughter(fMCStack);
359 TParticle *negDaughter = gamma->GetNegativeMCDaughter(fMCStack);
360 // cout << "generate Daughters: "<<posDaughter << "\t" << negDaughter << endl;
362 if(fMCStack && fEventCuts->GetSignalRejection() != 0){
363 Int_t isPosFromMBHeader
364 = fEventCuts->IsParticleFromBGEvent(gamma->GetMCLabelPositive(), fMCStack, fESDEvent);
365 Int_t isNegFromMBHeader
366 = fEventCuts->IsParticleFromBGEvent(gamma->GetMCLabelNegative(), fMCStack, fESDEvent);
367 if( (isNegFromMBHeader < 1) || (isPosFromMBHeader < 1)) continue;
370 if(posDaughter == NULL || negDaughter == NULL){
371 fKind = 9; // garbage
372 // cout << "one of the daughters not available" << endl;
373 } else if(posDaughter->GetMother(0) != negDaughter->GetMother(0) || (posDaughter->GetMother(0) == negDaughter->GetMother(0) && posDaughter->GetMother(0) ==-1)){
374 // Not Same Mother == Combinatorial Bck
376 // cout << "not the same mother" << endl;
378 if (posDaughter->GetPdgCode()) pdgCodePos = posDaughter->GetPdgCode(); else continue;
380 if (negDaughter->GetPdgCode()) pdgCodeNeg = negDaughter->GetPdgCode(); else continue;
381 // cout << "PDG codes daughters: " << pdgCodePos << "\t" << pdgCodeNeg << endl;
382 if(TMath::Abs(pdgCodePos)==11 && TMath::Abs(pdgCodeNeg)==11)
383 fKind = 10; //Electron Combinatorial
384 if(TMath::Abs(pdgCodePos)==11 && TMath::Abs(pdgCodeNeg)==11 && (posDaughter->GetMother(0) == negDaughter->GetMother(0) && posDaughter->GetMother(0) ==-1))
385 fKind = 15; //direct Electron Combinatorial
387 if(TMath::Abs(pdgCodePos)==211 && TMath::Abs(pdgCodeNeg)==211)
388 fKind = 11; //Pion Combinatorial
389 if((TMath::Abs(pdgCodePos)==211 && TMath::Abs(pdgCodeNeg)==2212) ||
390 (TMath::Abs(pdgCodePos)==2212 && TMath::Abs(pdgCodeNeg)==211))
391 fKind = 12; //Pion, Proton Combinatorics
392 if((TMath::Abs(pdgCodePos)==211 && TMath::Abs(pdgCodeNeg)==11) ||
393 (TMath::Abs(pdgCodePos)==11 && TMath::Abs(pdgCodeNeg)==211))
394 fKind = 13; //Pion, Electron Combinatorics
395 if (TMath::Abs(pdgCodePos)==321 || TMath::Abs(pdgCodeNeg)==321)
396 fKind = 14; //Kaon combinatorics
399 // cout << "same mother" << endl;
401 if (posDaughter->GetPdgCode()) pdgCodePos = posDaughter->GetPdgCode(); else continue;
403 if (negDaughter->GetPdgCode()) pdgCodeNeg = negDaughter->GetPdgCode(); else continue;
404 // cout << "PDG codes daughters: " << pdgCodePos << "\t" << pdgCodeNeg << endl;
406 if (gamma->GetMCParticle(fMCStack)->GetPdgCode()) pdgCode = gamma->GetMCParticle(fMCStack)->GetPdgCode(); else continue;
407 // cout << "PDG code: " << pdgCode << endl;
408 if(TMath::Abs(pdgCodePos)!=11 || TMath::Abs(pdgCodeNeg)!=11)
409 fKind = 2; // combinatorics from hadronic decays
410 else if ( !(pdgCodeNeg==pdgCodePos)){
411 TParticle *truePhotonCanditate = gamma->GetMCParticle(fMCStack);
412 Int_t motherLabelPhoton = truePhotonCanditate->GetMother(0);
414 fKind = 3; // pi0 Dalitz
415 else if (pdgCode == 221)
416 fKind = 4; // eta Dalitz
417 else if (!(negDaughter->GetUniqueID() != 5 || posDaughter->GetUniqueID() !=5)){
418 if(pdgCode == 22 && motherLabelPhoton < fMCStack->GetNprimary()){
419 fKind = 0; // primary photons
420 } else if (pdgCode == 22){
421 fKind = 5; //secondary photons
423 } else fKind = 9; //garbage
424 } else fKind = 9; //garbage
428 if (fTreeMaterialRec){
429 fTreeMaterialRec->Fill();
434 //________________________________________________________________________
435 Int_t AliAnalysisTaskMaterial::CountTracks09(){
436 Int_t fNumberOfESDTracks = 0;
437 if(fInputEvent->IsA()==AliESDEvent::Class()){
438 // Using standard function for setting Cuts
440 AliStack *fMCStack = NULL;
442 fMCStack= fMCEvent->Stack();
443 if (!fMCStack) return 0;
446 Bool_t selectPrimaries=kTRUE;
447 AliESDtrackCuts *EsdTrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2010(selectPrimaries);
448 EsdTrackCuts->SetMaxDCAToVertexZ(2);
449 EsdTrackCuts->SetEtaRange(-0.9, 0.9);
450 EsdTrackCuts->SetPtRange(0.15);
452 for(Int_t iTracks = 0; iTracks < fInputEvent->GetNumberOfTracks(); iTracks++){
453 AliESDtrack* curTrack = (AliESDtrack*) fInputEvent->GetTrack(iTracks);
454 if(!curTrack) continue;
455 if(EsdTrackCuts->AcceptTrack(curTrack) ){
457 if(fMCStack && fEventCuts->GetSignalRejection() != 0){
458 Int_t isFromMBHeader = fEventCuts->IsParticleFromBGEvent(abs(curTrack->GetLabel()), fMCStack, fESDEvent);
459 if( (isFromMBHeader < 1) ) continue;
462 fNumberOfESDTracks++;
468 return fNumberOfESDTracks;
471 //________________________________________________________________________
472 Int_t AliAnalysisTaskMaterial::CountTracks0914(){
473 Int_t fNumberOfESDTracks = 0;
474 if(fInputEvent->IsA()==AliESDEvent::Class()){
475 // Using standard function for setting Cuts
477 AliStack *fMCStack = NULL;
479 fMCStack= fMCEvent->Stack();
480 if (!fMCStack) return 0;
483 AliESDtrackCuts *EsdTrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
484 EsdTrackCuts->SetMaxDCAToVertexZ(5);
485 EsdTrackCuts->SetEtaRange(0.9, 1.4);
486 EsdTrackCuts->SetPtRange(0.15);
488 for(Int_t iTracks = 0; iTracks < fInputEvent->GetNumberOfTracks(); iTracks++){
489 AliESDtrack* curTrack = (AliESDtrack*) fInputEvent->GetTrack(iTracks);
490 if(!curTrack) continue;
491 if(EsdTrackCuts->AcceptTrack(curTrack) ){
493 if(fMCStack && fEventCuts->GetSignalRejection() != 0){
494 Int_t isFromMBHeader = fEventCuts->IsParticleFromBGEvent(abs(curTrack->GetLabel()), fMCStack, fESDEvent);
495 if( (isFromMBHeader < 1) ) continue;
498 fNumberOfESDTracks++;
501 EsdTrackCuts->SetEtaRange(-1.4, -0.9);
502 for(Int_t iTracks = 0; iTracks < fESDEvent->GetNumberOfTracks(); iTracks++){
503 AliESDtrack* curTrack = fESDEvent->GetTrack(iTracks);
504 if(!curTrack) continue;
505 if(EsdTrackCuts->AcceptTrack(curTrack) ){
507 if(fMCStack && fEventCuts->GetSignalRejection() != 0){
508 Int_t isFromMBHeader = fEventCuts->IsParticleFromBGEvent(abs(curTrack->GetLabel()), fMCStack, fESDEvent);
509 if( (isFromMBHeader < 1) ) continue;
512 fNumberOfESDTracks++;
519 return fNumberOfESDTracks;
522 //________________________________________________________________________
523 void AliAnalysisTaskMaterial::Terminate(Option_t *)
525 // if (fStreamMaterial){
526 // fStreamMaterial->GetFile()->Write();