Modified EveBase to compile with zmq<3.0. Corrected ThreadedSocket implementation.
[u/mrichter/AliRoot.git] / PWGGA / GammaConv / AliAnalysisTaskResolution.cxx
CommitLineData
72395bd9 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Authors: Friederike Bock *
5 * Version 1.0 *
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 * 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 **************************************************************************/
15
16////////////////////////////////////////////////
17//---------------------------------------------
18// QA Task for V0 Reader V1
19//---------------------------------------------
20////////////////////////////////////////////////
21
22#include "AliAnalysisTaskResolution.h"
23#include "TChain.h"
24#include "AliAnalysisManager.h"
25#include "TParticle.h"
26#include "TVectorF.h"
27#include "AliPIDResponse.h"
28#include "AliESDtrackCuts.h"
29#include "TFile.h"
30
31class iostream;
32
33using namespace std;
34
35ClassImp(AliAnalysisTaskResolution)
36
37//________________________________________________________________________
38AliAnalysisTaskResolution::AliAnalysisTaskResolution() : AliAnalysisTaskSE(),
63674d67 39 fV0Reader(NULL),
40 fConversionGammas(NULL),
41 fConversionCuts(NULL),
42 fTreeEvent(NULL),
43 fTreeResolution(NULL),
44 fPrimVtxZ(0.),
45 fNContrVtx(0),
46 fNESDtracksEta09(0),
47 fNESDtracksEta0914(0),
48 fNESDtracksEta14(0),
49 fGammaRecCoords(5),
50 fGammaMCCoords(5),
51 fChi2ndf(0),
52 fIsHeavyIon(0),
53 fIsMC(kFALSE),
54 fOutputList(NULL),
55 fEventList(NULL),
56 fResolutionList(NULL),
57 fESDEvent(NULL),
58 fMCEvent(NULL)
72395bd9 59{
60
61
62}
63
72395bd9 64//________________________________________________________________________
65AliAnalysisTaskResolution::AliAnalysisTaskResolution(const char *name) : AliAnalysisTaskSE(name),
63674d67 66 fV0Reader(NULL),
67 fConversionGammas(NULL),
68 fConversionCuts(NULL),
69 fTreeEvent(NULL),
70 fTreeResolution(NULL),
71 fPrimVtxZ(0.),
72 fNContrVtx(0),
73 fNESDtracksEta09(0),
74 fNESDtracksEta0914(0),
75 fNESDtracksEta14(0),
76 fGammaRecCoords(5),
77 fGammaMCCoords(5),
78 fChi2ndf(0),
79 fIsHeavyIon(0),
80 fIsMC(kFALSE),
81 fOutputList(NULL),
82 fEventList(NULL),
83 fResolutionList(NULL),
84 fESDEvent(NULL),
85 fMCEvent(NULL)
72395bd9 86{
63674d67 87 // Default constructor
72395bd9 88
63674d67 89 DefineInput(0, TChain::Class());
90 DefineOutput(1, TList::Class());
72395bd9 91}
92
93//________________________________________________________________________
94AliAnalysisTaskResolution::~AliAnalysisTaskResolution()
95{
63674d67 96 // default deconstructor
72395bd9 97
98}
99//________________________________________________________________________
100void AliAnalysisTaskResolution::UserCreateOutputObjects()
101{
63674d67 102 // Create User Output Objects
72395bd9 103
63674d67 104 if(fOutputList != NULL){
105 delete fOutputList;
106 fOutputList = NULL;
107 }
108 if(fOutputList == NULL){
109 fOutputList = new TList();
110 fOutputList->SetOwner(kTRUE);
111 }
112
113 fEventList = new TList();
114 fEventList->SetName("EventList");
115 fEventList->SetOwner(kTRUE);
116 fOutputList->Add(fEventList);
117
118 fTreeEvent = new TTree("Event","Event");
119 fTreeEvent->Branch("primVtxZ",&fPrimVtxZ,"fPrimVtxZ/F");
120 fTreeEvent->Branch("nContrVtx",&fNContrVtx,"fNContrVtx/I");
121 fTreeEvent->Branch("nGoodTracksEta09",&fNESDtracksEta09,"fNESDtracksEta09/I");
122 fTreeEvent->Branch("nGoodTracksEta0914",&fNESDtracksEta0914,"fNESDtracksEta0914/I");
123 fTreeEvent->Branch("nGoodTracksEta14",&fNESDtracksEta14,"fNESDtracksEta14/I");
124 fEventList->Add(fTreeEvent);
125
126 if (fIsMC){
127 fResolutionList = new TList();
128 fResolutionList->SetName("ResolutionList");
129 fResolutionList->SetOwner(kTRUE);
130 fOutputList->Add(fResolutionList);
131
132 fTreeResolution = new TTree("Resolution","Resolution");
133 fTreeResolution->Branch("RecCoords",&fGammaRecCoords);
134 fTreeResolution->Branch("MCCoords",&fGammaMCCoords);
135 fTreeResolution->Branch("chi2ndf",&fChi2ndf,"fChi2ndf/F");
136 fResolutionList->Add(fTreeResolution);
137 }
138 PostData(1, fOutputList);
72395bd9 139}
140
141//________________________________________________________________________
142void AliAnalysisTaskResolution::UserExec(Option_t *){
143
51ead2dc 144 fV0Reader=(AliV0ReaderV1*)AliAnalysisManager::GetAnalysisManager()->GetTask("V0ReaderV1");
72395bd9 145
51ead2dc 146 Int_t eventQuality = ((AliConversionCuts*)fV0Reader->GetConversionCuts())->GetEventQuality();
147 if(eventQuality != 0){// Event Not Accepted
148 return;
149 }
150 fESDEvent = (AliESDEvent*) InputEvent();
151 if (fESDEvent==NULL) return;
72395bd9 152 if(MCEvent()){
153 fMCEvent = MCEvent();
154 }
51ead2dc 155
63674d67 156 if(MCEvent()){
51ead2dc 157 // Process MC Particle
158 if(fConversionCuts->GetSignalRejection() != 0){
159// if(fESDEvent->IsA()==AliESDEvent::Class()){
160 fConversionCuts->GetNotRejectedParticles( fConversionCuts->GetSignalRejection(),
161 fConversionCuts->GetAcceptedHeader(),
162 fMCEvent);
163// }
164// else if(fInputEvent->IsA()==AliAODEvent::Class()){
165// fConversionCuts->GetNotRejectedParticles( fConversionCuts->GetSignalRejection(),
166// fConversionCuts->GetAcceptedHeader(),
167// fInputEvent);
168// }
169 }
170 }
171
172
63674d67 173 if(fIsHeavyIon > 0 && !fConversionCuts->IsCentralitySelected(fESDEvent)) return;
51ead2dc 174 fNESDtracksEta09 = CountTracks09(); // Estimate Event Multiplicity
175 fNESDtracksEta0914 = CountTracks0914(); // Estimate Event Multiplicity
176 fNESDtracksEta14 = fNESDtracksEta09 + fNESDtracksEta0914;
177 if(fESDEvent){
178 if(fESDEvent->GetPrimaryVertexTracks()->GetNContributors()>0) {
179 fNContrVtx = fESDEvent->GetPrimaryVertexTracks()->GetNContributors();
180 } else {
181 fNContrVtx = 0;
182 }
183// else if(fESDEvent->GetPrimaryVertexTracks()->GetNContributors()<1) {
184// if(fESDEvent->GetPrimaryVertexSPD()->GetNContributors()>0) {
185// fNContrVtx = fESDEvent->GetPrimaryVertexSPD()->GetNContributors();
186// } else if(fESDEvent->GetPrimaryVertexSPD()->GetNContributors()<1) {
187// fNContrVtx = 0;
188// }
189// }
190 }
191 fPrimVtxZ = fESDEvent->GetPrimaryVertex()->GetZ();
192
63674d67 193 if (fIsHeavyIon == 2){
51ead2dc 194 if (!(fNESDtracksEta09 > 20 && fNESDtracksEta09 < 80)) return;
195 }
196
197
198 if (fTreeEvent){
199 fTreeEvent->Fill();
200 }
201
202 fConversionGammas=fV0Reader->GetReconstructedGammas();
203 ProcessPhotons();
204 PostData(1, fOutputList);
72395bd9 205}
206
207
208///________________________________________________________________________
209void AliAnalysisTaskResolution::ProcessPhotons(){
210
211 // Fill Histograms for QA and MC
63674d67 212 for(Int_t firstGammaIndex=0;firstGammaIndex<fConversionGammas->GetEntriesFast();firstGammaIndex++){
213 AliAODConversionPhoton *gamma=dynamic_cast<AliAODConversionPhoton*>(fConversionGammas->At(firstGammaIndex));
214 if (gamma ==NULL) continue;
215 if(!fConversionCuts->PhotonIsSelected(gamma,fESDEvent)) continue;
216 fGammaRecCoords(0) = gamma->GetPhotonPt();
217 fGammaRecCoords(1) = gamma->GetPhotonPhi();
218 fGammaRecCoords(2) = gamma->GetPhotonEta();
219 fGammaRecCoords(3) = gamma->GetConversionRadius();
220 fGammaRecCoords(4) = gamma->GetConversionZ();
221 fChi2ndf = gamma->GetChi2perNDF();
222 if(MCEvent()){
72395bd9 223// cout << "generating MC stack"<< endl;
224 AliStack *fMCStack = fMCEvent->Stack();
225 if (!fMCStack) continue;
226 TParticle *posDaughter = gamma->GetPositiveMCDaughter(fMCStack);
227 TParticle *negDaughter = gamma->GetNegativeMCDaughter(fMCStack);
228// cout << "generate Daughters: "<<posDaughter << "\t" << negDaughter << endl;
51ead2dc 229 if(fMCStack && fConversionCuts->GetSignalRejection() != 0){
230 Int_t isPosFromMBHeader
231 = fConversionCuts->IsParticleFromBGEvent(gamma->GetMCLabelPositive(), fMCStack, fESDEvent);
232 Int_t isNegFromMBHeader
233 = fConversionCuts->IsParticleFromBGEvent(gamma->GetMCLabelNegative(), fMCStack, fESDEvent);
234 if( (isNegFromMBHeader < 1) || (isPosFromMBHeader < 1)) continue;
235 }
72395bd9 236
237 if(posDaughter == NULL || negDaughter == NULL){
238 continue;
239 } else if(posDaughter->GetMother(0) != negDaughter->GetMother(0) || (posDaughter->GetMother(0) == negDaughter->GetMother(0) && posDaughter->GetMother(0) ==-1)){
240 continue;
241 } else {
242// cout << "same mother" << endl;
243 Int_t pdgCodePos;
244 if (posDaughter->GetPdgCode()) pdgCodePos = posDaughter->GetPdgCode(); else continue;
245 Int_t pdgCodeNeg;
246 if (negDaughter->GetPdgCode()) pdgCodeNeg = negDaughter->GetPdgCode(); else continue;
247// cout << "PDG codes daughters: " << pdgCodePos << "\t" << pdgCodeNeg << endl;
248 Int_t pdgCode;
249 if (gamma->GetMCParticle(fMCStack)->GetPdgCode()) pdgCode = gamma->GetMCParticle(fMCStack)->GetPdgCode(); else continue;
250// cout << "PDG code: " << pdgCode << endl;
251 if(TMath::Abs(pdgCodePos)!=11 || TMath::Abs(pdgCodeNeg)!=11)
252 continue;
253 else if ( !(pdgCodeNeg==pdgCodePos)){
254 TParticle *truePhotonCanditate = gamma->GetMCParticle(fMCStack);
255 if(pdgCode == 111)
256 continue;
257 else if (pdgCode == 221)
258 continue;
259 else if (!(negDaughter->GetUniqueID() != 5 || posDaughter->GetUniqueID() !=5)){
260 if(pdgCode == 22){
63674d67 261 fGammaMCCoords(0) = truePhotonCanditate->Pt();
262 fGammaMCCoords(1) = gamma->GetNegativeMCDaughter(fMCStack)->Phi();
263 fGammaMCCoords(2) = gamma->GetNegativeMCDaughter(fMCStack)->Eta();
264 fGammaMCCoords(3) = gamma->GetNegativeMCDaughter(fMCStack)->R();
265 fGammaMCCoords(4) = gamma->GetNegativeMCDaughter(fMCStack)->Vz();
266
267 if (fTreeResolution){
268 fTreeResolution->Fill();
72395bd9 269 }
270 }
271 } else continue; //garbage
272 } else continue; //garbage
273 }
63674d67 274 }
72395bd9 275 }
276}
277
278//________________________________________________________________________
279Int_t AliAnalysisTaskResolution::CountTracks09(){
51ead2dc 280 Int_t fNumberOfESDTracks = 0;
281 if(fInputEvent->IsA()==AliESDEvent::Class()){
282 // Using standard function for setting Cuts
283
284 AliStack *fMCStack = NULL;
63674d67 285 if (MCEvent()){
51ead2dc 286 fMCStack= fMCEvent->Stack();
287 if (!fMCStack) return 0;
288 }
289
290 Bool_t selectPrimaries=kTRUE;
291 AliESDtrackCuts *EsdTrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2010(selectPrimaries);
292 EsdTrackCuts->SetMaxDCAToVertexZ(2);
293 EsdTrackCuts->SetEtaRange(-0.9, 0.9);
294 EsdTrackCuts->SetPtRange(0.15);
295
296 for(Int_t iTracks = 0; iTracks < fInputEvent->GetNumberOfTracks(); iTracks++){
297 AliESDtrack* curTrack = (AliESDtrack*) fInputEvent->GetTrack(iTracks);
298 if(!curTrack) continue;
299 if(EsdTrackCuts->AcceptTrack(curTrack) ){
300 if (fMCEvent){
301 if(fMCStack && fConversionCuts->GetSignalRejection() != 0){
302 Int_t isFromMBHeader
303 = fConversionCuts->IsParticleFromBGEvent(abs(curTrack->GetLabel()), fMCStack, fESDEvent);
304 if( (isFromMBHeader < 1) ) continue;
305 }
306 }
307 fNumberOfESDTracks++;
308 }
309 }
310 delete EsdTrackCuts;
311 EsdTrackCuts=0x0;
312 } else if(fInputEvent->IsA()==AliAODEvent::Class()){
313 for(Int_t iTracks = 0; iTracks<fInputEvent->GetNumberOfTracks(); iTracks++){
314 AliAODTrack* curTrack = (AliAODTrack*) fInputEvent->GetTrack(iTracks);
315 if(!curTrack->IsPrimaryCandidate()) continue;
316 if(abs(curTrack->Eta())>0.9) continue;
317 if(curTrack->Pt()<0.15) continue;
318 if(abs(curTrack->ZAtDCA())>2) continue;
319 fNumberOfESDTracks++;
320 }
321 }
72395bd9 322
51ead2dc 323 return fNumberOfESDTracks;
72395bd9 324}
325
51ead2dc 326//________________________________________________________________________
72395bd9 327Int_t AliAnalysisTaskResolution::CountTracks0914(){
51ead2dc 328 Int_t fNumberOfESDTracks = 0;
329 if(fInputEvent->IsA()==AliESDEvent::Class()){
330 // Using standard function for setting Cuts
331
332 AliStack *fMCStack = NULL;
63674d67 333 if (MCEvent()){
51ead2dc 334 fMCStack= fMCEvent->Stack();
335 if (!fMCStack) return 0;
336 }
72395bd9 337
51ead2dc 338 AliESDtrackCuts *EsdTrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
339 EsdTrackCuts->SetMaxDCAToVertexZ(5);
340 EsdTrackCuts->SetEtaRange(0.9, 1.4);
341 EsdTrackCuts->SetPtRange(0.15);
342
343 for(Int_t iTracks = 0; iTracks < fInputEvent->GetNumberOfTracks(); iTracks++){
344 AliESDtrack* curTrack = (AliESDtrack*) fInputEvent->GetTrack(iTracks);
345 if(!curTrack) continue;
346 if(EsdTrackCuts->AcceptTrack(curTrack) ){
347 if (fMCEvent){
348 if(fMCStack && fConversionCuts->GetSignalRejection() != 0){
349 Int_t isFromMBHeader
350 = fConversionCuts->IsParticleFromBGEvent(abs(curTrack->GetLabel()), fMCStack, fESDEvent);
351 if( (isFromMBHeader < 1) ) continue;
352 }
353 }
354 fNumberOfESDTracks++;
355 }
356 }
357 EsdTrackCuts->SetEtaRange(-1.4, -0.9);
358 for(Int_t iTracks = 0; iTracks < fESDEvent->GetNumberOfTracks(); iTracks++){
359 AliESDtrack* curTrack = fESDEvent->GetTrack(iTracks);
360 if(!curTrack) continue;
361 if(EsdTrackCuts->AcceptTrack(curTrack) ){
362 if (fMCEvent){
363 if(fMCStack && fConversionCuts->GetSignalRejection() != 0){
364 Int_t isFromMBHeader
365 = fConversionCuts->IsParticleFromBGEvent(abs(curTrack->GetLabel()), fMCStack, fESDEvent);
366 if( (isFromMBHeader < 1) ) continue;
367 }
368 }
369 fNumberOfESDTracks++;
370 }
371 }
372 delete EsdTrackCuts;
373 EsdTrackCuts=0x0;
374 } else if(fInputEvent->IsA()==AliAODEvent::Class()){
375 for(Int_t iTracks = 0; iTracks<fInputEvent->GetNumberOfTracks(); iTracks++){
376 AliAODTrack* curTrack = (AliAODTrack*) fInputEvent->GetTrack(iTracks);
377 if(abs(curTrack->Eta())<0.9 || abs(curTrack->Eta())>1.4 ) continue;
378 if(curTrack->Pt()<0.15) continue;
379 if(abs(curTrack->ZAtDCA())>5) continue;
380 fNumberOfESDTracks++;
381 }
382 }
383
384 return fNumberOfESDTracks;
72395bd9 385}
386
387//________________________________________________________________________
388void AliAnalysisTaskResolution::Terminate(Option_t *)
389{
63674d67 390
72395bd9 391}