]> git.uio.no Git - u/mrichter/AliRoot.git/blame_incremental - PWGDQ/dielectron/AliAnalysisTaskReducedTree.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame_incremental
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (non-incremental) | history | HEAD
MFT track shit tool added
[u/mrichter/AliRoot.git] / PWGDQ / dielectron / AliAnalysisTaskReducedTree.cxx
... / ...
CommitLineData
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* 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// Analysis task for creating a reduced data tree //
19// //
20///////////////////////////////////////////////////////////////////////////
21
22
23#include <TChain.h>
24#include <TH1D.h>
25#include <TFile.h>
26
27#include <AliCFContainer.h>
28#include <AliInputEventHandler.h>
29#include <AliESDInputHandler.h>
30#include <AliAODInputHandler.h>
31#include <AliAnalysisManager.h>
32#include <AliVEvent.h>
33#include <AliESDEvent.h>
34#include <AliAODEvent.h>
35#include <AliESDHeader.h>
36#include <AliAODHeader.h>
37#include <AliAODTrack.h>
38//#include <AliAODForwardMult.h>
39//#include <AliForwardUtil.h>
40#include <AliTriggerAnalysis.h>
41#include <AliESDtrackCuts.h>
42#include <AliVZDC.h>
43#include <AliESDv0.h>
44#include <AliESDv0Cuts.h>
45#include <AliESDv0KineCuts.h>
46#include <AliESDFMD.h>
47#include <AliVCluster.h>
48#include <AliAODTracklets.h>
49#include <AliMultiplicity.h>
50#include <AliPIDResponse.h>
51#include "AliDielectron.h"
52#include "AliDielectronHistos.h"
53#include "AliDielectronMC.h"
54#include "AliDielectronVarManager.h"
55#include "AliFlowTrackCuts.h"
56#include "AliFlowBayesianPID.h"
57
58#include "AliReducedEvent.h"
59#include "AliAnalysisTaskReducedTree.h"
60
61#include <iostream>
62using std::cout;
63using std::endl;
64
65ClassImp(AliAnalysisTaskReducedTree)
66
67
68//_________________________________________________________________________________
69AliAnalysisTaskReducedTree::AliAnalysisTaskReducedTree() :
70 AliAnalysisTaskSE(),
71 fListDielectron(),
72 fListHistos(),
73 fSelectPhysics(kFALSE),
74 fTriggerMask(AliVEvent::kAny),
75 fRejectPileup(kFALSE),
76 fFillTrackInfo(kTRUE),
77 fFillDielectronInfo(kTRUE),
78 fFillV0Info(kTRUE),
79 fFillGammaConversions(kTRUE),
80 fFillK0s(kTRUE),
81 fFillLambda(kTRUE),
82 fFillALambda(kTRUE),
83 fFillCaloClusterInfo(kTRUE),
84 fFillFMDSectorInfo(kFALSE),
85 fFillFMDChannelInfo(kFALSE),
86 //fFillCorrectedFMDInfo(kTRUE),
87 fFillFriendInfo(kTRUE),
88 fEventFilter(0x0),
89 fTrackFilter(0x0),
90 fFlowTrackFilter(0x0),
91 fK0sCuts(0x0),
92 fLambdaCuts(0x0),
93 fGammaConvCuts(0x0),
94 fK0sPionCuts(0x0),
95 fLambdaProtonCuts(0x0),
96 fLambdaPionCuts(0x0),
97 fGammaElectronCuts(0x0),
98 fV0OpenCuts(0x0),
99 fV0StrongCuts(0x0),
100 fK0sMassRange(),
101 fLambdaMassRange(),
102 fGammaMassRange(),
103 fV0Histos(0x0),
104 fTreeFile(0x0),
105 fTree(0x0),
106 fFriendTreeFile(0x0),
107 fFriendTree(0x0),
108 fReducedEvent(0x0),
109 fReducedEventFriend(0x0)
110{
111 //
112 // Constructor
113 //
114}
115
116//_________________________________________________________________________________
117AliAnalysisTaskReducedTree::AliAnalysisTaskReducedTree(const char *name) :
118 AliAnalysisTaskSE(name),
119 fListDielectron(),
120 fListHistos(),
121 fSelectPhysics(kFALSE),
122 fTriggerMask(AliVEvent::kAny),
123 fRejectPileup(kFALSE),
124 fFillTrackInfo(kTRUE),
125 fFillDielectronInfo(kTRUE),
126 fFillV0Info(kTRUE),
127 fFillGammaConversions(kTRUE),
128 fFillK0s(kTRUE),
129 fFillLambda(kTRUE),
130 fFillALambda(kTRUE),
131 fFillCaloClusterInfo(kTRUE),
132 fFillFMDSectorInfo(kFALSE),
133 fFillFMDChannelInfo(kFALSE),
134 //fFillCorrectedFMDInfo(kTRUE),
135 fFillFriendInfo(kTRUE),
136 fEventFilter(0x0),
137 fTrackFilter(0x0),
138 fFlowTrackFilter(0x0),
139 fK0sCuts(0x0),
140 fLambdaCuts(0x0),
141 fGammaConvCuts(0x0),
142 fK0sPionCuts(0x0),
143 fLambdaProtonCuts(0x0),
144 fLambdaPionCuts(0x0),
145 fGammaElectronCuts(0x0),
146 fV0OpenCuts(0x0),
147 fV0StrongCuts(0x0),
148 fK0sMassRange(),
149 fLambdaMassRange(),
150 fGammaMassRange(),
151 fV0Histos(0x0),
152 fTreeFile(0x0),
153 fTree(0x0),
154 fFriendTreeFile(0x0),
155 fFriendTree(0x0),
156 fReducedEvent(0x0),
157 fReducedEventFriend(0x0)
158{
159 //
160 // Constructor
161 //
162 fK0sMassRange[0] = 0.4; fK0sMassRange[1] = 0.6;
163 fLambdaMassRange[0] = 1.08; fLambdaMassRange[1] = 1.15;
164 fGammaMassRange[0] = 0.0; fGammaMassRange[1] = 0.1;
165
166 DefineInput(0,TChain::Class());
167 //DefineInput(2,AliAODForwardMult::Class());
168 DefineOutput(1, TList::Class()); // QA histograms
169 DefineOutput(2, TTree::Class()); // reduced information tree
170 //if(fFillFriendInfo) DefineOutput(3, TTree::Class()); // reduced information tree with friends
171 //DefineOutput(2, TTree::Class()); // reduced information tree with friends
172 //DefineOutput(2, TTree::Class()); // reduced information tree
173
174 fListHistos.SetName("QAhistograms");
175 fListDielectron.SetOwner();
176 fListHistos.SetOwner(kFALSE);
177}
178
179
180//_________________________________________________________________________________
181void AliAnalysisTaskReducedTree::UserCreateOutputObjects()
182{
183 //
184 // Add all histogram manager histogram lists to the output TList
185 //
186
187 if (!fListHistos.IsEmpty() || fTree || fFriendTree) return; //already initialised
188
189 TIter nextDie(&fListDielectron);
190 AliDielectron *die=0;
191 while ( (die=static_cast<AliDielectron*>(nextDie())) ){
192 die->Init();
193 if (die->GetHistogramList()) fListHistos.Add(const_cast<THashList*>(die->GetHistogramList()));
194 }
195 if(fV0Histos) fListHistos.Add(const_cast<THashList*>(fV0Histos->GetHistogramList()));
196
197 if(fFillFriendInfo) {
198 fFriendTree = new TTree("DstFriendTree","Reduced ESD information");
199 fReducedEventFriend = new AliReducedEventFriend();
200 fFriendTree->Branch("Event",&fReducedEventFriend,16000,99);
201 }
202
203 //fTreeFile = new TFile("dstTree.root", "RECREATE");
204 OpenFile(2);
205 fTree = new TTree("DstTree","Reduced ESD information");
206 fReducedEvent = new AliReducedEvent("DstEvent");
207 fTree->Branch("Event",&fReducedEvent,16000,99);
208
209 PostData(1, &fListHistos);
210 PostData(2, fTree);
211 //if(fFillFriendInfo) PostData(3, fFriendTree);
212 //PostData(2, fFriendTree);
213 //PostData(1, fTree);
214}
215
216//_________________________________________________________________________________
217void AliAnalysisTaskReducedTree::UserExec(Option_t *)
218{
219 //
220 // Main loop. Called for every event
221 //
222 AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
223 Bool_t isESD=man->GetInputEventHandler()->IsA()==AliESDInputHandler::Class();
224 Bool_t isAOD=man->GetInputEventHandler()->IsA()==AliAODInputHandler::Class();
225
226 AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
227 if (!inputHandler) return;
228
229 if ( inputHandler->GetPIDResponse() ){
230 AliDielectronVarManager::SetPIDResponse( inputHandler->GetPIDResponse() );
231 } else {
232 AliFatal("This task needs the PID response attached to the input event handler!");
233 }
234
235 // Was event selected ?
236 UInt_t isSelected = AliVEvent::kAny;
237 if(fSelectPhysics && inputHandler){
238 if((isESD && inputHandler->GetEventSelection()) || isAOD){
239 isSelected = inputHandler->IsEventSelected();
240 isSelected&=fTriggerMask;
241 }
242 }
243
244 if(isSelected==0) {
245 return;
246 }
247/*
248 cout<<"get AOD"<<endl;
249
250 AliAODEvent* aodEvent = AliForwardUtil::GetAODEvent(this);
251 if (!aodEvent) return;
252 cout<<"got AOD"<<endl;
253
254 TObject* obj = aodEvent->FindListObject("Forward");
255 if (!obj) return;
256 cout<<"got AOD forward"<<endl;
257
258 AliAODForwardMult* aodForward = static_cast<AliAODForwardMult*>(obj);
259
260 //if (!aodForward->CheckEvent(mask,ipZmin,ipZmax,cMin,cMax)) return 0;
261
262 Double_t ret = 0;
263 const TH2D& d2Ndetadphi = aodForward->GetHistogram();
264
265
266 cout<<d2Ndetadphi.GetXaxis()->GetNbins()<<endl;
267*/
268
269 // FMD corrections
270// AliAODEvent* aodEvent;
271// TObject* obj;
272// AliAODForwardMult* aodForward=NULL;
273//
274// if(fFillCorrectedFMDInfo){
275// aodEvent = AliForwardUtil::GetAODEvent(this);
276// if (!aodEvent) return;
277//
278// obj = aodEvent->FindListObject("Forward");
279// if (!obj) return;}
280//
281// aodForward = static_cast<AliAODForwardMult*>(obj);
282//
283// //if (!aodForward->CheckE vent(fTriggerMask,ipZmin,ipZmax,cMin,cMax)) return 0;
284//
285// const TH2D& d2Ndetadphi = aodForward->GetHistogram();
286
287
288 // fill event histograms before event filter
289 Double_t values[AliDielectronVarManager::kNMaxValues]={0};
290 AliDielectronVarManager::Fill(InputEvent(),values);
291
292 TIter nextDie(&fListDielectron);
293 AliDielectron *die=0;
294 while ( (die=static_cast<AliDielectron*>(nextDie())) ){
295 AliDielectronHistos *h=die->GetHistoManager();
296 if (h){
297 if (h->GetHistogramList()->FindObject("Event_noCuts"))
298 h->FillClass("Event_noCuts",AliDielectronVarManager::kNMaxValues,values);
299 }
300 }
301 nextDie.Reset();
302
303 //event filter
304 if (fEventFilter) {
305 if (!fEventFilter->IsSelected(InputEvent())) return;
306 }
307
308 //pileup
309 if (fRejectPileup){
310 if (InputEvent()->IsPileupFromSPD(3,0.8,3.,2.,5.)) return;
311 }
312
313 //bz for AliKF
314 Double_t bz = InputEvent()->GetMagneticField();
315 AliKFParticle::SetField( bz );
316
317 //Process event in all AliDielectron instances
318 fReducedEvent->ClearEvent();
319 if(fFillFriendInfo) fReducedEventFriend->ClearEvent();
320 FillEventInfo();
321 //if(fFillCorrectedFMDInfo) FillCorrectedFMDInfo(d2Ndetadphi); //Fill corrected FMD info
322 if(fFillV0Info) FillV0PairInfo();
323
324 Short_t idie=0;
325 if(fFillDielectronInfo) {
326 while((die=static_cast<AliDielectron*>(nextDie()))){
327 die->Process(InputEvent());
328 FillDielectronPairInfo(die, idie);
329 ++idie;
330 }
331 }
332 nextDie.Reset();
333
334 if(fFillTrackInfo) FillTrackInfo();
335 if(fFillFriendInfo) FillFriendEventInfo(); // Q-vector calculation
336
337 fTree->Fill();
338 if(fFillFriendInfo) fFriendTree->Fill();
339
340 // if there are candidate pairs, add the information to the reduced tree
341 PostData(1, &fListHistos);
342 PostData(2, fTree);
343 //if(fFillFriendInfo) PostData(3, fFriendTree);
344 //PostData(2, fFriendTree);
345 //PostData(2, fTree);
346}
347
348
349//_________________________________________________________________________________
350void AliAnalysisTaskReducedTree::FillEventInfo()
351{
352 //
353 // fill reduced event information
354 //
355 AliVEvent* event = InputEvent();
356 // Was event selected ?
357 AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
358 Bool_t isESD = (event->IsA()==AliESDEvent::Class());
359 Bool_t isAOD = (event->IsA()==AliAODEvent::Class());
360
361 AliESDEvent* esdEvent = 0x0;
362 if(isESD) esdEvent = static_cast<AliESDEvent*>(event);
363 AliAODEvent* aodEvent = 0x0;
364 if(isAOD) aodEvent = static_cast<AliAODEvent*>(event);
365
366 AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
367 UInt_t isSelected = AliVEvent::kAny;
368 if(inputHandler){
369 if((isESD && inputHandler->GetEventSelection()) || isAOD){
370 isSelected = inputHandler->IsEventSelected();
371 isSelected&=fTriggerMask;
372 }
373 }
374 Double_t values[AliDielectronVarManager::kNMaxValues];
375 AliDielectronVarManager::Fill(event, values);
376
377 fReducedEvent->fRunNo = event->GetRunNumber();
378 fReducedEvent->fBC = event->GetBunchCrossNumber();
379 fReducedEvent->fEventType = event->GetEventType();
380 fReducedEvent->fTriggerMask = event->GetTriggerMask();
381 fReducedEvent->fIsPhysicsSelection = (isSelected!=0 ? kTRUE : kFALSE);
382 fReducedEvent->fIsSPDPileup = event->IsPileupFromSPD(3,0.8,3.,2.,5.);
383 fReducedEvent->fIsSPDPileupMultBins = event->IsPileupFromSPDInMultBins();
384 AliVVertex* eventVtx = 0x0;
385 if(isESD) eventVtx = const_cast<AliESDVertex*>(esdEvent->GetPrimaryVertexTracks());
386 if(isAOD) eventVtx = const_cast<AliAODVertex*>(aodEvent->GetPrimaryVertex());
387 if(eventVtx) {
388 fReducedEvent->fVtx[0] = (isESD ? ((AliESDVertex*)eventVtx)->GetX() : ((AliAODVertex*)eventVtx)->GetX());
389 fReducedEvent->fVtx[1] = (isESD ? ((AliESDVertex*)eventVtx)->GetY() : ((AliAODVertex*)eventVtx)->GetY());
390 fReducedEvent->fVtx[2] = (isESD ? ((AliESDVertex*)eventVtx)->GetZ() : ((AliAODVertex*)eventVtx)->GetZ());
391 fReducedEvent->fNVtxContributors = eventVtx->GetNContributors();
392 }
393 if(isESD) {
394 eventVtx = const_cast<AliESDVertex*>(esdEvent->GetPrimaryVertexTPC());
395 fReducedEvent->fEventNumberInFile = esdEvent->GetEventNumberInFile();
396 fReducedEvent->fL0TriggerInputs = esdEvent->GetHeader()->GetL0TriggerInputs();
397 fReducedEvent->fL1TriggerInputs = esdEvent->GetHeader()->GetL1TriggerInputs();
398 fReducedEvent->fL2TriggerInputs = esdEvent->GetHeader()->GetL2TriggerInputs();
399 fReducedEvent->fIRIntClosestIntMap[0] = esdEvent->GetHeader()->GetIRInt1ClosestInteractionMap();
400 fReducedEvent->fIRIntClosestIntMap[1] = esdEvent->GetHeader()->GetIRInt2ClosestInteractionMap();
401 if(eventVtx) {
402 fReducedEvent->fVtxTPC[0] = ((AliESDVertex*)eventVtx)->GetX();
403 fReducedEvent->fVtxTPC[1] = ((AliESDVertex*)eventVtx)->GetY();
404 fReducedEvent->fVtxTPC[2] = ((AliESDVertex*)eventVtx)->GetZ();
405 fReducedEvent->fNVtxTPCContributors = eventVtx->GetNContributors();
406 }
407 fReducedEvent->fTimeStamp = esdEvent->GetTimeStamp();
408 fReducedEvent->fNpileupSPD = esdEvent->GetNumberOfPileupVerticesSPD();
409 fReducedEvent->fNpileupTracks = esdEvent->GetNumberOfPileupVerticesTracks();
410 fReducedEvent->fNPMDtracks = esdEvent->GetNumberOfPmdTracks();
411 fReducedEvent->fNTRDtracks = esdEvent->GetNumberOfTrdTracks();
412 fReducedEvent->fNTRDtracklets = esdEvent->GetNumberOfTrdTracklets();
413
414 AliESDZDC* zdc = esdEvent->GetESDZDC();
415 if(zdc) {
416 for(Int_t i=0; i<5; ++i) fReducedEvent->fZDCnEnergy[i] = zdc->GetZN1TowerEnergy()[i];
417 for(Int_t i=5; i<10; ++i) fReducedEvent->fZDCnEnergy[i] = zdc->GetZN2TowerEnergy()[i-5];
418 for(Int_t i=0; i<5; ++i) fReducedEvent->fZDCpEnergy[i] = zdc->GetZP1TowerEnergy()[i];
419 for(Int_t i=5; i<10; ++i) fReducedEvent->fZDCpEnergy[i] = zdc->GetZP2TowerEnergy()[i-5];
420 }
421 }
422 if(isAOD) {
423 AliAODHeader * header = dynamic_cast<AliAODHeader*>(aodEvent->GetHeader());
424 if(!header) AliFatal("Not a standard AOD");
425
426
427 fReducedEvent->fIRIntClosestIntMap[0] = header->GetIRInt1ClosestInteractionMap();
428 fReducedEvent->fIRIntClosestIntMap[1] = header->GetIRInt2ClosestInteractionMap();
429 fReducedEvent->fEventNumberInFile = header->GetEventNumberESDFile();
430 fReducedEvent->fL0TriggerInputs = header->GetL0TriggerInputs();
431 fReducedEvent->fL1TriggerInputs = header->GetL1TriggerInputs();
432 fReducedEvent->fL2TriggerInputs = header->GetL2TriggerInputs();
433 fReducedEvent->fTimeStamp = 0;
434 fReducedEvent->fNpileupSPD = aodEvent->GetNumberOfPileupVerticesSPD();
435 fReducedEvent->fNpileupTracks = aodEvent->GetNumberOfPileupVerticesTracks();
436 fReducedEvent->fNPMDtracks = aodEvent->GetNPmdClusters();
437 fReducedEvent->fNTRDtracks = 0;
438 fReducedEvent->fNTRDtracklets = 0;
439
440 AliAODZDC* zdc = aodEvent->GetZDCData();
441 if(zdc) {
442 for(Int_t i=0; i<5; ++i) fReducedEvent->fZDCnEnergy[i] = zdc->GetZNATowerEnergy()[i];
443 for(Int_t i=5; i<10; ++i) fReducedEvent->fZDCnEnergy[i] = zdc->GetZNCTowerEnergy()[i-5];
444 for(Int_t i=0; i<5; ++i) fReducedEvent->fZDCpEnergy[i] = zdc->GetZPATowerEnergy()[i];
445 for(Int_t i=5; i<10; ++i) fReducedEvent->fZDCpEnergy[i] = zdc->GetZPCTowerEnergy()[i-5];
446 }
447 }
448
449 // Fill TZERO information
450 if(isESD) {
451 const AliESDTZERO* tzero = esdEvent->GetESDTZERO();
452 if(tzero) {
453 fReducedEvent->fT0start = tzero->GetT0();
454 fReducedEvent->fT0zVertex = tzero->GetT0zVertex();
455 for(Int_t i = 0;i<24;i++)
456 fReducedEvent->fT0amplitude[i] = tzero->GetT0amplitude()[i];
457 for(Int_t i = 0;i<3;i++)
458 fReducedEvent->fT0TOF[i] = tzero->GetT0TOF()[i];
459 for(Int_t i = 0;i<3;i++)
460 fReducedEvent->fT0TOFbest[i] = tzero->GetT0TOFbest()[i];
461 fReducedEvent->fT0pileup = tzero->GetPileupFlag();
462 fReducedEvent->fT0sattelite = tzero->GetSatellite();
463 }
464 }
465 if(isAOD) {
466 AliAODTZERO* tzero = aodEvent->GetTZEROData();
467 if(tzero) {
468 fReducedEvent->fT0start = -999.; // not available
469 fReducedEvent->fT0zVertex = tzero->GetT0zVertex();
470 for(Int_t i = 0;i<26;i++)
471 fReducedEvent->fT0amplitude[i] = tzero->GetAmp(i);
472 for(Int_t i = 0;i<3;i++)
473 fReducedEvent->fT0TOF[i] = tzero->GetT0TOF()[i];
474 for(Int_t i = 0;i<3;i++)
475 fReducedEvent->fT0TOFbest[i] = tzero->GetT0TOFbest()[i];
476 fReducedEvent->fT0pileup = tzero->GetPileupFlag();
477 fReducedEvent->fT0sattelite = tzero->GetSatellite();
478 }
479 }
480
481 if(fFillFMDChannelInfo&&isESD) fReducedEvent->fIsFMDReduced = kFALSE;
482 if((fFillFMDSectorInfo||fFillFMDChannelInfo)&&isESD) FillFMDInfo();
483
484 AliCentrality *centrality = event->GetCentrality();
485 if(centrality) {
486 fReducedEvent->fCentrality[0] = centrality->GetCentralityPercentile("V0M");
487 fReducedEvent->fCentrality[1] = centrality->GetCentralityPercentile("CL1");
488 fReducedEvent->fCentrality[2] = centrality->GetCentralityPercentile("TRK");
489 fReducedEvent->fCentrality[3] = centrality->GetCentralityPercentile("ZEMvsZDC");
490 fReducedEvent->fCentQuality = centrality->GetQuality();
491 }
492
493 //cout << "event vtxZ/cent: " << fReducedEvent->fVtx[2] << "/" << fReducedEvent->fCentrality[0] << endl;
494
495 fReducedEvent->fNtracks[0] = event->GetNumberOfTracks();
496 fReducedEvent->fSPDntracklets = GetSPDTrackletMultiplicity(event, -1.0, 1.0);
497 for(Int_t ieta=0; ieta<32; ++ieta)
498 fReducedEvent->fSPDntrackletsEta[ieta] = GetSPDTrackletMultiplicity(event, -1.6+0.1*ieta, -1.6+0.1*(ieta+1));
499
500 AliVVZERO* vzero = event->GetVZEROData();
501 for(Int_t i=0;i<64;++i)
502 fReducedEvent->fVZEROMult[i] = vzero->GetMultiplicity(i);
503
504 // EMCAL/PHOS clusters
505 if(fFillCaloClusterInfo) FillCaloClusters();
506
507 // TODO FMD multiplicities
508
509}
510
511
512//_________________________________________________________________________________
513void AliAnalysisTaskReducedTree::FillCaloClusters() {
514 //
515 // Fill info about the calorimeter clusters
516 //
517 AliVEvent* event = InputEvent();
518 Int_t nclusters = event->GetNumberOfCaloClusters();
519
520 fReducedEvent->fNCaloClusters = 0;
521 for(Int_t iclus=0; iclus<nclusters; ++iclus) {
522 AliVCluster* cluster = event->GetCaloCluster(iclus);
523
524 TClonesArray& clusters = *(fReducedEvent->fCaloClusters);
525 AliReducedCaloCluster *reducedCluster=new(clusters[fReducedEvent->fNCaloClusters]) AliReducedCaloCluster();
526
527 reducedCluster->fType = (cluster->IsEMCAL() ? AliReducedCaloCluster::kEMCAL : AliReducedCaloCluster::kPHOS);
528 reducedCluster->fEnergy = cluster->E();
529 reducedCluster->fTrackDx = cluster->GetTrackDx();
530 reducedCluster->fTrackDz = cluster->GetTrackDz();
531 reducedCluster->fM20 = cluster->GetM20();
532 reducedCluster->fM02 = cluster->GetM02();
533 reducedCluster->fDispersion = cluster->GetDispersion();
534 fReducedEvent->fNCaloClusters += 1;
535 } // end loop over clusters
536}
537
538
539//_________________________________________________________________________________
540void AliAnalysisTaskReducedTree::FillFriendEventInfo() {
541 //
542 // Fill event info into the friend tree
543 //
544 // Add here calculated Q-vector components from all detectors
545 for(Int_t idet=0; idet<AliReducedEventFriend::kNdetectors; ++idet) {
546 fReducedEvent->GetQvector(fReducedEventFriend->fQvector[idet], idet);
547 for(Int_t ih=0; ih<fgkNMaxHarmonics; ++ih)
548 fReducedEventFriend->fEventPlaneStatus[idet][ih] = AliReducedEventFriend::kRaw;
549 }
550}
551
552
553//_________________________________________________________________________________
554void AliAnalysisTaskReducedTree::FillFMDInfo()
555{
556 //
557 // fill reduced FMD information
558 //
559 AliVEvent* event = InputEvent();
560 Bool_t isESD = (event->IsA()==AliESDEvent::Class());
561
562 if(!isESD) return;
563
564 AliESDEvent* esdEvent = 0x0;
565 if(isESD) esdEvent = static_cast<AliESDEvent*>(event);
566
567 AliESDFMD* esdFmd = esdEvent->GetFMDData();
568 //const AliFMDFloatMap multMap = esdFmd->MultiplicityMap();
569 //const AliFMDFloatMap etaMap = esdFmd->EtaMap();
570
571 //esdFmd->Print();
572 Int_t nFMD=0;
573 Int_t id=-1;
574 Int_t maxDet=3;
575 Int_t maxRing=2;
576 Int_t maxSector;
577 Int_t maxStrip;
578 Float_t m=0.0;
579 Double_t phi;
580 Char_t ring;
581 Float_t fmdMult;
582 Float_t msum=0;
583 UShort_t fmdDet=0;
584 Int_t phiBin=0;
585
586 for(UShort_t det = 1; det <= maxDet; ++det) {
587 (det == 1 ? maxRing=1 : maxRing=2);
588 for(UShort_t ir = 0; ir < maxRing; ++ir) {
589 ring = (ir == 0 ? 'I' : 'O');
590 (ir == 0 ? maxSector=20 : maxSector=40);
591 (ir == 0 ? maxStrip=512 : maxStrip=256);
592 nFMD=-1;
593 for(UShort_t sec = 0; sec < maxSector; ++sec) {
594 phi = esdFmd->Phi(det, ring, sec, 0)/180.*TMath::Pi();
595 phiBin = Int_t (phi/2/TMath::Pi()*maxSector);
596 fmdMult = 0;
597 for(UShort_t str = 0; str < maxStrip; ++str) {
598 ++id;
599 m = esdFmd->Multiplicity(det, ring, sec, str);
600 //cout << "det/ir/sec/str/m :: " << det << "/" << ir << "/" << sec << "/" << str << "/" << m << endl;
601 if(fFillFMDChannelInfo)
602 if(m<1.e-6) continue;
603 if(m == AliESDFMD::kInvalidMult) m=0;
604 fmdMult += m;
605 msum+=m;
606 ++nFMD;
607
608 if(fFillFMDChannelInfo){
609 if(m>15.) m=15.;
610 m = UShort_t(m*4369+0.5);
611 TClonesArray& fmd = *(fReducedEvent->GetFMD(fmdDet));
612 AliReducedFMD *reducedFMD=new(fmd[nFMD]) AliReducedFMD();
613 fReducedEvent->fNFMDchannels[fmdDet] += 1;
614 reducedFMD->fMultiplicity = m;
615 //reducedFMD->fEta = esdFmd->Eta(det, ring, 0, str);
616 reducedFMD->fId = id;
617 }
618 } // end loop over strips
619
620 if(fFillFMDSectorInfo) {
621 TClonesArray& fmd = *(fReducedEvent->GetFMD(fmdDet));
622 AliReducedFMD *reducedFMD=new(fmd[phiBin]) AliReducedFMD();
623 reducedFMD->fMultiplicity = fmdMult;
624 fReducedEvent->fNFMDchannels[fmdDet] += 1;
625 //cout<<sec<<" "<<fmdMult<<endl;
626 fmdMult=0.0;
627 }
628 } // end loop over sectors
629
630 fReducedEvent->fFMDtotalMult[fmdDet] = msum;
631 msum=0.0;
632 ++fmdDet;
633 id=-1;
634 } // end loop over rings
635 } // end loop over detectors
636}
637
638////_________________________________________________________________________________
639//void AliAnalysisTaskReducedTree::FillCorrectedFMDInfo(const TH2D& fmdhist)
640//{
641//
642//
643//Int_t nEta = fmdhist.GetXaxis()->GetNbins();
644//Int_t nPhi = fmdhist.GetYaxis()->GetNbins();
645////Int_t nBins= fmdhist.GetNbins();
646//Float_t eta=0.0;
647//Float_t phi=0.0;
648//Float_t mult=0.0;
649//
650//
651//Int_t nFMD=0;
652////fReducedEvent->fNCorFmdChannels = 0;
653//for (Int_t e = 1; e <= nEta; e++) {
654// eta = fmdhist.GetXaxis()->GetBinCenter(e);
655// for (Int_t p = 1; p <= nPhi; p++) {
656// phi = fmdhist.GetYaxis()->GetBinCenter(p);
657// mult = fmdhist.GetBinContent(e, p);
658// //TClonesArray& Corfmd = *(fReducedEvent->fCorFMD);
659// //AliReducedFMD *reducedCorFMD=new(Corfmd[nFMD]) AliReducedCorFMD();
660// std::cout<<mult<<" "<<eta<<" "<<phi<<std::endl;
661// }
662//
663// //fReducedEvent->fNCorFmdChannels += 1;
664// nFMD += 1;
665////reducedFMD->fCorMultiplicity = mult;
666////reducedFMD->fCorEta = eta;
667////reducedFMD->fCorPhi = phi;
668//
669//
670//}
671//
672//
673//}
674
675
676//_________________________________________________________________________________
677void AliAnalysisTaskReducedTree::FillTrackInfo()
678{
679 //
680 // fill reduced track information
681 //
682 AliVEvent* event = InputEvent();
683 Bool_t isESD = (event->IsA()==AliESDEvent::Class());
684 Bool_t isAOD = (event->IsA()==AliAODEvent::Class());
685
686 AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
687 AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
688 AliPIDResponse* pidResponse = inputHandler->GetPIDResponse();
689
690 // find all the tracks which belong to a V0 stored in the reduced event
691 UShort_t trackIdsV0[4][20000]={{0}};
692 UShort_t trackIdsPureV0[4][20000]={{0}};
693 Int_t nV0LegsTagged[4] = {0}; Int_t nPureV0LegsTagged[4] = {0};
694 Bool_t leg1Found[4]; Bool_t leg2Found[4];
695 for(Int_t iv0=0;iv0<fReducedEvent->fNV0candidates[1];++iv0) {
696 AliReducedPair* pair = fReducedEvent->GetV0Pair(iv0);
697 if(!pair) continue;
698 Int_t pairId = 0; Bool_t isPureV0 = kFALSE;
699 if(pair->fCandidateId==AliReducedPair::kGammaConv) {
700 pairId=0;
701 if(pair->IsPureV0Gamma()) isPureV0 = kTRUE;
702 }
703 if(pair->fCandidateId==AliReducedPair::kK0sToPiPi) {
704 pairId=1;
705 if(pair->IsPureV0K0s()) isPureV0 = kTRUE;
706 }
707 if(pair->fCandidateId==AliReducedPair::kLambda0ToPPi) {
708 pairId=2;
709 if(pair->IsPureV0Lambda()) isPureV0 = kTRUE;
710 }
711 if(pair->fCandidateId==AliReducedPair::kALambda0ToPPi) {
712 pairId=3;
713 if(pair->IsPureV0ALambda()) isPureV0 = kTRUE;
714 }
715
716 leg1Found[pairId] = kFALSE; leg2Found[pairId] = kFALSE;
717 for(Int_t it=0;it<nV0LegsTagged[pairId];++it) {
718 if(trackIdsV0[pairId][it]==pair->fLegIds[0]) leg1Found[pairId]=kTRUE;
719 if(trackIdsV0[pairId][it]==pair->fLegIds[1]) leg2Found[pairId]=kTRUE;
720 }
721 // if the legs of this V0 were not already stored then add them now to the list
722 if(!leg1Found[pairId]) {trackIdsV0[pairId][nV0LegsTagged[pairId]] = pair->fLegIds[0]; ++nV0LegsTagged[pairId];}
723 if(!leg2Found[pairId]) {trackIdsV0[pairId][nV0LegsTagged[pairId]] = pair->fLegIds[1]; ++nV0LegsTagged[pairId];}
724
725 if(isPureV0) {
726 leg1Found[pairId] = kFALSE; leg2Found[pairId] = kFALSE;
727 for(Int_t it=0;it<nPureV0LegsTagged[pairId];++it) {
728 if(trackIdsPureV0[pairId][it]==pair->fLegIds[0]) leg1Found[pairId]=kTRUE;
729 if(trackIdsPureV0[pairId][it]==pair->fLegIds[1]) leg2Found[pairId]=kTRUE;
730 }
731 // if the legs of this pure V0 were not already stored then add them now to the list
732 if(!leg1Found[pairId]) {trackIdsPureV0[pairId][nPureV0LegsTagged[pairId]] = pair->fLegIds[0]; ++nPureV0LegsTagged[pairId];}
733 if(!leg2Found[pairId]) {trackIdsPureV0[pairId][nPureV0LegsTagged[pairId]] = pair->fLegIds[1]; ++nPureV0LegsTagged[pairId];}
734 }
735 }
736
737 // find all the tracks which belong to a stored dielectron pair
738 UShort_t trackIdsDiele[20000]={0};
739 Int_t nDieleLegsTagged = 0;
740 for(Int_t idie=0;idie<fReducedEvent->NDielectrons();++idie) {
741 AliReducedPair* pair = fReducedEvent->GetDielectronPair(idie);
742 leg1Found[0]=kFALSE; leg2Found[0]=kFALSE;
743 for(Int_t it=0; it<nDieleLegsTagged; ++it) {
744 if(trackIdsDiele[it]==pair->fLegIds[0]) leg1Found[0]=kTRUE;
745 if(trackIdsDiele[it]==pair->fLegIds[1]) leg2Found[0]=kTRUE;
746 }
747 // if the legs of this dielectron were not already stored then add them now to the list
748 if(!leg1Found[0]) {trackIdsDiele[nDieleLegsTagged] = pair->fLegIds[0]; ++nDieleLegsTagged;}
749 if(!leg2Found[0]) {trackIdsDiele[nDieleLegsTagged] = pair->fLegIds[1]; ++nDieleLegsTagged;}
750 }
751
752 AliESDtrack* esdTrack=0;
753 AliAODTrack* aodTrack=0;
754 Int_t ntracks=event->GetNumberOfTracks();
755 Int_t trackId = 0;
756 Bool_t usedForV0[4] = {kFALSE};
757 Bool_t usedForPureV0[4] = {kFALSE};
758 Bool_t usedForV0Or = kFALSE;
759 Bool_t usedForDielectron = kFALSE;
760 for(Int_t itrack=0; itrack<ntracks; ++itrack){
761 AliVParticle *particle=event->GetTrack(itrack);
762 if(isESD) {
763 esdTrack=static_cast<AliESDtrack*>(particle);
764 trackId = esdTrack->GetID();
765 }
766 if(isAOD) {
767 aodTrack=static_cast<AliAODTrack*>(particle);
768 trackId = aodTrack->GetID();
769 }
770 // check whether this track belongs to a V0 stored in the reduced event
771 usedForV0Or = kFALSE;
772 for(Int_t i=0; i<4; ++i) {
773 usedForV0[i] = kFALSE;
774 for(Int_t ii=0; ii<nV0LegsTagged[i]; ++ii) {
775 if(UShort_t(trackId)==trackIdsV0[i][ii]) {
776 usedForV0[i] = kTRUE;
777 //cout << "track " << trackId << " used for V0 type " << i << endl;
778 break;
779 }
780 }
781 usedForV0Or = usedForV0Or || usedForV0[i];
782 usedForPureV0[i] = kFALSE;
783 for(Int_t ii=0; ii<nPureV0LegsTagged[i]; ++ii) {
784 if(UShort_t(trackId)==trackIdsPureV0[i][ii]) {
785 usedForPureV0[i] = kTRUE;
786 //cout << "track " << trackId << " used for pure V0 type " << i << endl;
787 break;
788 }
789 }
790 }
791 // check whether this track belongs to a dielectron stored in the reduced event
792 usedForDielectron = kFALSE;
793 for(Int_t ii=0; ii<nDieleLegsTagged; ++ii) {
794 if(UShort_t(trackId)==trackIdsDiele[ii]) {
795 usedForDielectron = kTRUE;
796 break;
797 }
798 }
799
800 ULong_t status = (isESD ? esdTrack->GetStatus() : aodTrack->GetStatus());
801 //cout << "TRACK" << endl;
802 for(Int_t ibit=0; ibit<32; ++ibit) {
803 if(status & (ULong_t(1)<<ibit)) {
804 //cout << "bit " << ibit << endl;
805 fReducedEvent->fNtracksPerTrackingFlag[ibit] += 1;
806 }
807 }
808
809 //apply track cuts
810 if(!usedForV0Or && !usedForDielectron && fTrackFilter && !fTrackFilter->IsSelected(particle)) continue;
811 //cout << "storing track " << trackId << endl;
812
813 TClonesArray& tracks = *(fReducedEvent->fTracks);
814 AliReducedTrack *reducedParticle=new(tracks[fReducedEvent->fNtracks[1]]) AliReducedTrack();
815
816 Double_t values[AliDielectronVarManager::kNMaxValues];
817 AliDielectronVarManager::Fill(particle, values);
818 reducedParticle->fStatus = status;//(ULong_t)values[AliDielectronVarManager::kTrackStatus];
819 reducedParticle->fGlobalPhi = values[AliDielectronVarManager::kPhi];
820 reducedParticle->fGlobalPt = values[AliDielectronVarManager::kPt]*values[AliDielectronVarManager::kCharge];
821 reducedParticle->fGlobalEta = values[AliDielectronVarManager::kEta];
822 reducedParticle->fMomentumInner = values[AliDielectronVarManager::kPIn];
823 reducedParticle->fDCA[0] = values[AliDielectronVarManager::kImpactParXY];
824 reducedParticle->fDCA[1] = values[AliDielectronVarManager::kImpactParZ];
825 reducedParticle->fTrackLength = values[AliDielectronVarManager::kTrackLength];
826
827 reducedParticle->fITSclusterMap = (UChar_t)values[AliDielectronVarManager::kITSclusterMap];
828 reducedParticle->fITSsignal = values[AliDielectronVarManager::kITSsignal];
829 reducedParticle->fITSnSig[0] = values[AliDielectronVarManager::kITSnSigmaEle];
830 reducedParticle->fITSnSig[1] = values[AliDielectronVarManager::kITSnSigmaPio];
831 reducedParticle->fITSnSig[2] = values[AliDielectronVarManager::kITSnSigmaKao];
832 reducedParticle->fITSnSig[3] = values[AliDielectronVarManager::kITSnSigmaPro];
833 reducedParticle->fITSchi2 = values[AliDielectronVarManager::kITSchi2Cl];
834
835 reducedParticle->fTPCNcls = (UChar_t)values[AliDielectronVarManager::kNclsTPC];
836 reducedParticle->fTPCNclsF = (UChar_t)values[AliDielectronVarManager::kNFclsTPC];
837 reducedParticle->fTPCNclsIter1 = (UChar_t)values[AliDielectronVarManager::kNclsTPCiter1];
838 reducedParticle->fTPCsignal = values[AliDielectronVarManager::kTPCsignal];
839 reducedParticle->fTPCsignalN = values[AliDielectronVarManager::kTPCsignalN];
840 reducedParticle->fTPCnSig[0] = values[AliDielectronVarManager::kTPCnSigmaEle];
841 reducedParticle->fTPCnSig[1] = values[AliDielectronVarManager::kTPCnSigmaPio];
842 reducedParticle->fTPCnSig[2] = values[AliDielectronVarManager::kTPCnSigmaKao];
843 reducedParticle->fTPCnSig[3] = values[AliDielectronVarManager::kTPCnSigmaPro];
844 reducedParticle->fTPCClusterMap = EncodeTPCClusterMap(particle, isAOD);
845 reducedParticle->fTPCchi2 = values[AliDielectronVarManager::kTPCchi2Cl];
846
847 reducedParticle->fTOFbeta = values[AliDielectronVarManager::kTOFbeta];
848 reducedParticle->fTOFnSig[0] = values[AliDielectronVarManager::kTOFnSigmaEle];
849 reducedParticle->fTOFnSig[1] = values[AliDielectronVarManager::kTOFnSigmaPio];
850 reducedParticle->fTOFnSig[2] = values[AliDielectronVarManager::kTOFnSigmaKao];
851 reducedParticle->fTOFnSig[3] = values[AliDielectronVarManager::kTOFnSigmaPro];
852
853 Double_t trdProbab[AliPID::kSPECIES]={0.0};
854
855 if(fFlowTrackFilter) {
856 // switch on the first bit if this particle should be used for the event plane
857 if(fFlowTrackFilter->IsSelected(particle)) reducedParticle->fFlags |= (UShort_t(1)<<0);
858 }
859 for(Int_t iV0type=0;iV0type<4;++iV0type) {
860 if(usedForV0[iV0type]) reducedParticle->fFlags |= (UShort_t(1)<<(iV0type+1));
861 if(usedForPureV0[iV0type]) reducedParticle->fFlags |= (UShort_t(1)<<(iV0type+8));
862 }
863
864 if(isESD){
865 //AliESDtrack *track=static_cast<AliESDtrack*>(particle);
866 reducedParticle->fTrackId = (UShort_t)esdTrack->GetID();
867 reducedParticle->fTPCCrossedRows = (UChar_t)esdTrack->GetTPCCrossedRows();
868 //reducedParticle->fTPCClusterMap = EncodeTPCClusterMap(esdTrack);
869 const AliExternalTrackParam* tpcInner = esdTrack->GetInnerParam();
870 reducedParticle->fTPCPhi = (tpcInner ? tpcInner->Phi() : 0.0);
871 reducedParticle->fTPCPt = (tpcInner ? tpcInner->Pt() : 0.0);
872 reducedParticle->fTPCEta = (tpcInner ? tpcInner->Eta() : 0.0);
873
874 reducedParticle->fTOFdeltaBC = esdTrack->GetTOFDeltaBC();
875
876 reducedParticle->fTRDntracklets[0] = esdTrack->GetTRDntracklets();
877 reducedParticle->fTRDntracklets[1] = esdTrack->GetTRDntrackletsPID();
878 pidResponse->ComputeTRDProbability(esdTrack,AliPID::kSPECIES,trdProbab,AliTRDPIDResponse::kLQ1D);
879 reducedParticle->fTRDpid[0] = trdProbab[AliPID::kElectron];
880 reducedParticle->fTRDpid[1] = trdProbab[AliPID::kPion];
881 pidResponse->ComputeTRDProbability(esdTrack,AliPID::kSPECIES,trdProbab,AliTRDPIDResponse::kLQ2D);
882 reducedParticle->fTRDpidLQ2D[0] = trdProbab[AliPID::kElectron];
883 reducedParticle->fTRDpidLQ2D[1] = trdProbab[AliPID::kPion];
884
885 for(Int_t idx=0; idx<3; ++idx) if(esdTrack->GetKinkIndex(idx)>0) reducedParticle->fFlags |= (1<<(5+idx));
886 if(esdTrack->IsEMCAL()) reducedParticle->fCaloClusterId = esdTrack->GetEMCALcluster();
887 if(esdTrack->IsPHOS()) reducedParticle->fCaloClusterId = esdTrack->GetPHOScluster();
888 }
889 if(isAOD) {
890 //AliAODTrack *track=static_cast<AliAODTrack*>(particle);
891 const AliExternalTrackParam* tpcInner = aodTrack->GetInnerParam();
892 reducedParticle->fTPCPhi = (tpcInner ? tpcInner->Phi() : 0.0);
893 reducedParticle->fTPCPt = (tpcInner ? tpcInner->Pt() : 0.0);
894 reducedParticle->fTPCEta = (tpcInner ? tpcInner->Eta() : 0.0);
895
896 reducedParticle->fTrackId = aodTrack->GetID();
897 reducedParticle->fITSsignal = aodTrack->GetITSsignal();
898 if(pidResponse) {
899 reducedParticle->fITSnSig[0] = pidResponse->NumberOfSigmasITS(aodTrack,AliPID::kElectron);
900 reducedParticle->fITSnSig[1] = pidResponse->NumberOfSigmasITS(aodTrack,AliPID::kPion);
901 reducedParticle->fITSnSig[2] = pidResponse->NumberOfSigmasITS(aodTrack,AliPID::kKaon);
902 reducedParticle->fITSnSig[3] = pidResponse->NumberOfSigmasITS(aodTrack,AliPID::kProton);
903 }
904 reducedParticle->fTRDntracklets[0] = aodTrack->GetTRDntrackletsPID();
905 reducedParticle->fTRDntracklets[1] = aodTrack->GetTRDntrackletsPID();
906 pidResponse->ComputeTRDProbability(aodTrack,AliPID::kSPECIES,trdProbab,AliTRDPIDResponse::kLQ1D);
907 reducedParticle->fTRDpid[0] = trdProbab[AliPID::kElectron];
908 reducedParticle->fTRDpid[1] = trdProbab[AliPID::kPion];
909 pidResponse->ComputeTRDProbability(aodTrack,AliPID::kSPECIES,trdProbab,AliTRDPIDResponse::kLQ2D);
910 reducedParticle->fTRDpidLQ2D[0] = trdProbab[AliPID::kElectron];
911 reducedParticle->fTRDpidLQ2D[1] = trdProbab[AliPID::kPion];
912
913 if(aodTrack->IsEMCAL()) reducedParticle->fCaloClusterId = aodTrack->GetEMCALcluster();
914 if(aodTrack->IsPHOS()) reducedParticle->fCaloClusterId = aodTrack->GetPHOScluster();
915 if(values[AliDielectronVarManager::kKinkIndex0]>0.0) reducedParticle->fFlags |= (1<<5);
916 }
917
918 fReducedEvent->fNtracks[1] += 1;
919 }
920}
921
922
923//_________________________________________________________________________________
924void AliAnalysisTaskReducedTree::FillDielectronPairInfo(AliDielectron* die, Short_t iDie)
925{
926 //
927 // fill reduced pair information
928 //
929 Bool_t hasMC=AliDielectronMC::Instance()->HasMC();
930
931 for(Int_t iType=0; iType<3; ++iType) {
932
933 const TObjArray* array = die->GetPairArray(iType);
934 if(!array || array->GetEntriesFast()==0) continue;
935
936 for(Int_t iCandidate=0; iCandidate<array->GetEntriesFast(); ++iCandidate) {
937 AliDielectronPair* pair = (AliDielectronPair*)array->At(iCandidate);
938 Double_t values[AliDielectronVarManager::kNMaxValues];
939 AliDielectronVarManager::Fill(pair, values);
940
941 TClonesArray& tracks = *(fReducedEvent->fCandidates);
942 AliReducedPair *reducedParticle=
943 new (tracks[fReducedEvent->fNV0candidates[1]+fReducedEvent->fNDielectronCandidates]) AliReducedPair();
944 // !!! hardcoded flag for dielectron id
945 reducedParticle->fCandidateId = (iDie==0 ? AliReducedPair::kJpsiToEE : AliReducedPair::kPhiToKK);
946 reducedParticle->fPairType = (Char_t)values[AliDielectronVarManager::kPairType];
947 reducedParticle->fLegIds[0] = (UShort_t)(static_cast<AliVTrack*>(pair->GetFirstDaughterP()))->GetID();
948 reducedParticle->fLegIds[1] = (UShort_t)(static_cast<AliVTrack*>(pair->GetSecondDaughterP()))->GetID();
949 reducedParticle->fMass[0] = values[AliDielectronVarManager::kM];
950 reducedParticle->fMass[1] = -999.;
951 reducedParticle->fMass[2] = -999.;
952 reducedParticle->fMass[3] = -999.;
953 reducedParticle->fPhi = values[AliDielectronVarManager::kPhi]; // in the [-pi,pi] interval
954 if(reducedParticle->fPhi<0.0) reducedParticle->fPhi = 2.0*TMath::Pi() + reducedParticle->fPhi; // converted to [0,2pi]
955 reducedParticle->fPt = values[AliDielectronVarManager::kPt];
956 reducedParticle->fEta = values[AliDielectronVarManager::kEta];
957 reducedParticle->fLxy = values[AliDielectronVarManager::kPseudoProperTime];
958 reducedParticle->fLxyErr = values[AliDielectronVarManager::kPseudoProperTimeErr];
959 reducedParticle->fPointingAngle = values[AliDielectronVarManager::kCosPointingAngle];
960
961 reducedParticle->fMCid = 0;
962 if(hasMC) {
963 AliDielectronMC::Instance()->ConnectMCEvent();
964 const TObjArray* mcSignals = die->GetMCSignals();
965 for(Int_t iSig=0; iSig<mcSignals->GetEntries(); ++iSig) {
966 if(iSig>31) break;
967 AliDielectronMC *mc=AliDielectronMC::Instance();
968 if(mc->IsMCTruth(pair, (AliDielectronSignalMC*)mcSignals->At(iSig))) {
969 reducedParticle->fMCid = reducedParticle->fMCid | (1<<iSig);
970 }
971 }
972 } // end if has MC
973 fReducedEvent->fNDielectronCandidates += 1;
974 } // end loop over candidates
975 } // end loop over pair type
976}
977
978
979//_________________________________________________________________________________
980void AliAnalysisTaskReducedTree::FillV0PairInfo()
981{
982 //
983 // fill reduced pair information
984 //
985 AliESDEvent* esd = (AliESDEvent*)InputEvent();
986 const AliESDVertex *primaryVertex = esd->GetPrimaryVertex();
987 AliKFVertex primaryVertexKF(*primaryVertex);
988
989 fReducedEvent->fNV0candidates[0] = InputEvent()->GetNumberOfV0s();
990
991 if(!(fFillK0s || fFillLambda || fFillALambda || fFillGammaConversions)) return;
992
993 Double_t valuesPos[AliDielectronVarManager::kNMaxValues];
994 Double_t valuesNeg[AliDielectronVarManager::kNMaxValues];
995
996 if(fV0OpenCuts) {
997 fV0OpenCuts->SetEvent(esd);
998 fV0OpenCuts->SetPrimaryVertex(&primaryVertexKF);
999 }
1000 if(fV0StrongCuts) {
1001 fV0StrongCuts->SetEvent(esd);
1002 fV0StrongCuts->SetPrimaryVertex(&primaryVertexKF);
1003 }
1004
1005 Int_t pdgV0=0; Int_t pdgP=0; Int_t pdgN=0;
1006 for(Int_t iV0=0; iV0<InputEvent()->GetNumberOfV0s(); ++iV0) { // loop over V0s
1007 AliESDv0 *v0 = esd->GetV0(iV0);
1008
1009 AliESDtrack* legPos = esd->GetTrack(v0->GetPindex());
1010 AliESDtrack* legNeg = esd->GetTrack(v0->GetNindex());
1011
1012 if(legPos->GetSign() == legNeg->GetSign()) {
1013 continue;
1014 }
1015
1016 Bool_t v0ChargesAreCorrect = (legPos->GetSign()==+1 ? kTRUE : kFALSE);
1017 legPos = (!v0ChargesAreCorrect ? esd->GetTrack(v0->GetNindex()) : legPos);
1018 legNeg = (!v0ChargesAreCorrect ? esd->GetTrack(v0->GetPindex()) : legNeg);
1019
1020 pdgV0=0; pdgP=0; pdgN=0;
1021 Bool_t goodK0s = kTRUE; Bool_t goodLambda = kTRUE; Bool_t goodALambda = kTRUE; Bool_t goodGamma = kTRUE;
1022 if(fV0OpenCuts) {
1023 goodK0s = kFALSE; goodLambda = kFALSE; goodALambda = kFALSE; goodGamma = kFALSE;
1024 Bool_t processV0 = fV0OpenCuts->ProcessV0(v0, pdgV0, pdgP, pdgN);
1025 if(processV0 && TMath::Abs(pdgV0)==310 && TMath::Abs(pdgP)==211 && TMath::Abs(pdgN)==211) {
1026 goodK0s = kTRUE;
1027 if(fK0sPionCuts && (!fK0sPionCuts->IsSelected(legPos) || !fK0sPionCuts->IsSelected(legNeg))) goodK0s = kFALSE;
1028 }
1029 if(processV0 && pdgV0==3122 && (TMath::Abs(pdgP)==211 || TMath::Abs(pdgP)==2212) && (TMath::Abs(pdgN)==211 || TMath::Abs(pdgN)==2212)) {
1030 goodLambda = kTRUE;
1031 if(fLambdaProtonCuts && !fLambdaProtonCuts->IsSelected(legPos)) goodLambda = kFALSE;
1032 if(fLambdaPionCuts && !fLambdaPionCuts->IsSelected(legNeg)) goodLambda = kFALSE;
1033 }
1034 if(processV0 && pdgV0==-3122 && (TMath::Abs(pdgP)==211 || TMath::Abs(pdgP)==2212) && (TMath::Abs(pdgN)==211 || TMath::Abs(pdgN)==2212)) {
1035 goodALambda = kTRUE;
1036 if(fLambdaProtonCuts && !fLambdaProtonCuts->IsSelected(legNeg)) goodALambda = kFALSE;
1037 if(fLambdaPionCuts && !fLambdaPionCuts->IsSelected(legPos)) goodALambda = kFALSE;
1038 }
1039 if(processV0 && TMath::Abs(pdgV0)==22 && TMath::Abs(pdgP)==11 && TMath::Abs(pdgN)==11) {
1040 goodGamma = kTRUE;
1041 if(fGammaElectronCuts && (!fGammaElectronCuts->IsSelected(legPos) || !fGammaElectronCuts->IsSelected(legNeg))) goodGamma = kFALSE;
1042 }
1043 //cout << "open cuts pdgV0/pdgP/pdgN/processV0 : " << pdgV0 << "/" << pdgP << "/" << pdgN << "/" << processV0 << endl;
1044 //cout << "good K0s/Lambda/ALambda/Gamma : " << goodK0s << "/" << goodLambda << "/" << goodALambda << "/" << goodGamma << endl;
1045 }
1046
1047 Bool_t veryGoodK0s = kFALSE; Bool_t veryGoodLambda = kFALSE; Bool_t veryGoodALambda = kFALSE; Bool_t veryGoodGamma = kFALSE;
1048 if(fV0StrongCuts && (goodK0s || goodLambda || goodALambda || goodGamma)) {
1049 pdgV0=0; pdgP=0; pdgN=0;
1050 Bool_t processV0 = fV0StrongCuts->ProcessV0(v0, pdgV0, pdgP, pdgN);
1051 if(processV0 && goodK0s && TMath::Abs(pdgV0)==310 && TMath::Abs(pdgP)==211 && TMath::Abs(pdgN)==211)
1052 veryGoodK0s = kTRUE;
1053 if(processV0 && goodLambda && pdgV0==3122 && (TMath::Abs(pdgP)==211 || TMath::Abs(pdgP)==2212) && (TMath::Abs(pdgN)==211 || TMath::Abs(pdgN)==2212))
1054 veryGoodLambda = kTRUE;
1055 if(processV0 && goodALambda && pdgV0==-3122 && (TMath::Abs(pdgP)==211 || TMath::Abs(pdgP)==2212) && (TMath::Abs(pdgN)==211 || TMath::Abs(pdgN)==2212))
1056 veryGoodALambda = kTRUE;
1057 if(processV0 && goodGamma && TMath::Abs(pdgV0)==22 && TMath::Abs(pdgP)==11 && TMath::Abs(pdgN)==11)
1058 veryGoodGamma = kTRUE;
1059 //cout << "strong cuts pdgV0/pdgP/pdgN/processV0 : " << pdgV0 << "/" << pdgP << "/" << pdgN << "/" << processV0 << endl;
1060 //cout << "very good K0s/Lambda/ALambda/Gamma : " << veryGoodK0s << "/" << veryGoodLambda << "/" << veryGoodALambda << "/" << veryGoodGamma << endl;
1061 }
1062
1063 if(!((goodK0s && fFillK0s) ||
1064 (goodLambda && fFillLambda) ||
1065 (goodALambda && fFillALambda) ||
1066 (goodGamma && fFillGammaConversions))) continue;
1067
1068 // Fill the V0 information into the tree for 4 hypothesis: K0s, Lambda, Anti-Lambda and gamma conversion
1069 AliReducedPair* k0sReducedPair = FillV0PairInfo(v0, AliReducedPair::kK0sToPiPi, legPos, legNeg, &primaryVertexKF, v0ChargesAreCorrect);
1070 AliReducedPair* lambdaReducedPair = FillV0PairInfo(v0, AliReducedPair::kLambda0ToPPi, legPos, legNeg, &primaryVertexKF, v0ChargesAreCorrect);
1071 AliReducedPair* alambdaReducedPair = FillV0PairInfo(v0, AliReducedPair::kALambda0ToPPi, legPos, legNeg, &primaryVertexKF, v0ChargesAreCorrect);
1072 AliReducedPair* gammaReducedPair = FillV0PairInfo(v0, AliReducedPair::kGammaConv, legPos, legNeg, &primaryVertexKF, v0ChargesAreCorrect);
1073
1074 if(fFillK0s && goodK0s && k0sReducedPair->fMass[0]>fK0sMassRange[0] && k0sReducedPair->fMass[0]<fK0sMassRange[1]) {
1075 TClonesArray& tracks = *(fReducedEvent->fCandidates);
1076 AliReducedPair *goodK0sPair = new (tracks[fReducedEvent->fNV0candidates[1]]) AliReducedPair(*k0sReducedPair);
1077 goodK0sPair->fMass[0] = k0sReducedPair->fMass[0];
1078 goodK0sPair->fMass[1] = lambdaReducedPair->fMass[0];
1079 goodK0sPair->fMass[2] = alambdaReducedPair->fMass[0];
1080 goodK0sPair->fMass[3] = gammaReducedPair->fMass[0];
1081 if(veryGoodK0s) goodK0sPair->fMCid |= (UInt_t(1)<<1);
1082 fReducedEvent->fNV0candidates[1] += 1;
1083 } else {goodK0s=kFALSE;}
1084 if(fFillLambda && goodLambda && lambdaReducedPair->fMass[0]>fLambdaMassRange[0] && lambdaReducedPair->fMass[0]<fLambdaMassRange[1]) {
1085 TClonesArray& tracks = *(fReducedEvent->fCandidates);
1086 AliReducedPair *goodLambdaPair = new (tracks[fReducedEvent->fNV0candidates[1]]) AliReducedPair(*lambdaReducedPair);
1087 goodLambdaPair->fMass[0] = k0sReducedPair->fMass[0];
1088 goodLambdaPair->fMass[1] = lambdaReducedPair->fMass[0];
1089 goodLambdaPair->fMass[2] = alambdaReducedPair->fMass[0];
1090 goodLambdaPair->fMass[3] = gammaReducedPair->fMass[0];
1091 if(veryGoodLambda) goodLambdaPair->fMCid |= (UInt_t(1)<<2);
1092 fReducedEvent->fNV0candidates[1] += 1;
1093 } else {goodLambda=kFALSE;}
1094 if(fFillALambda && goodALambda && alambdaReducedPair->fMass[0]>fLambdaMassRange[0] && alambdaReducedPair->fMass[0]<fLambdaMassRange[1]) {
1095 TClonesArray& tracks = *(fReducedEvent->fCandidates);
1096 AliReducedPair *goodALambdaPair = new (tracks[fReducedEvent->fNV0candidates[1]]) AliReducedPair(*alambdaReducedPair);
1097 goodALambdaPair->fMass[0] = k0sReducedPair->fMass[0];
1098 goodALambdaPair->fMass[1] = lambdaReducedPair->fMass[0];
1099 goodALambdaPair->fMass[2] = alambdaReducedPair->fMass[0];
1100 goodALambdaPair->fMass[3] = gammaReducedPair->fMass[0];
1101 if(veryGoodALambda) goodALambdaPair->fMCid |= (UInt_t(1)<<3);
1102 fReducedEvent->fNV0candidates[1] += 1;
1103 } else {goodALambda = kFALSE;}
1104 //cout << "gamma mass: " << gammaReducedPair->fMass[0] << endl;
1105 if(fFillGammaConversions && goodGamma && gammaReducedPair->fMass[0]>fGammaMassRange[0] && gammaReducedPair->fMass[0]<fGammaMassRange[1]) {
1106 TClonesArray& tracks = *(fReducedEvent->fCandidates);
1107 AliReducedPair *goodGammaPair = new (tracks[fReducedEvent->fNV0candidates[1]]) AliReducedPair(*gammaReducedPair);
1108 goodGammaPair->fMass[0] = k0sReducedPair->fMass[0];
1109 goodGammaPair->fMass[1] = lambdaReducedPair->fMass[0];
1110 goodGammaPair->fMass[2] = alambdaReducedPair->fMass[0];
1111 goodGammaPair->fMass[3] = gammaReducedPair->fMass[0];
1112 if(veryGoodGamma) goodGammaPair->fMCid |= (UInt_t(1)<<4);
1113 fReducedEvent->fNV0candidates[1] += 1;
1114 } else {goodGamma=kFALSE;}
1115 delete k0sReducedPair;
1116 delete lambdaReducedPair;
1117 delete alambdaReducedPair;
1118 delete gammaReducedPair;
1119
1120 if(!(goodK0s || goodLambda || goodALambda || goodGamma)) continue;
1121
1122 // Fill histograms and the CF container
1123 AliDielectronVarManager::Fill(legPos, valuesPos);
1124 AliDielectronVarManager::Fill(legNeg, valuesNeg);
1125
1126 if(fV0Histos && fV0Histos->GetHistogramList()->FindObject("V0Track_Pos"))
1127 fV0Histos->FillClass("V0Track_Pos", AliDielectronVarManager::kNMaxValues, valuesPos);
1128 if(fV0Histos && fV0Histos->GetHistogramList()->FindObject("V0Track_Neg"))
1129 fV0Histos->FillClass("V0Track_Neg", AliDielectronVarManager::kNMaxValues, valuesNeg);
1130 } // end loop over V0s
1131}
1132
1133
1134//_________________________________________________________________________________
1135AliReducedPair* AliAnalysisTaskReducedTree::FillV0PairInfo(AliESDv0* v0, Int_t id,
1136 AliESDtrack* legPos, AliESDtrack* legNeg,
1137 AliKFVertex* vtxKF, Bool_t chargesAreCorrect) {
1138 //
1139 // Create a reduced V0 object and fill it
1140 //
1141 AliReducedPair* reducedPair=new AliReducedPair();
1142 reducedPair->fCandidateId = id;
1143 reducedPair->fPairType = v0->GetOnFlyStatus(); // on the fly status
1144 reducedPair->fLegIds[0] = legPos->GetID();
1145 reducedPair->fLegIds[1] = legNeg->GetID();
1146 if(!reducedPair->fPairType) { // offline
1147 UInt_t pidPos = AliPID::kPion;
1148 if(id==AliReducedPair::kLambda0ToPPi) pidPos = AliPID::kProton;
1149 if(id==AliReducedPair::kGammaConv) pidPos = AliPID::kElectron;
1150 UInt_t pidNeg = AliPID::kPion;
1151 if(id==AliReducedPair::kALambda0ToPPi) pidNeg = AliPID::kProton;
1152 if(id==AliReducedPair::kGammaConv) pidNeg = AliPID::kElectron;
1153 reducedPair->fMass[0] = v0->GetEffMass(pidPos, pidNeg);
1154 reducedPair->fPhi = v0->Phi();
1155 if(reducedPair->fPhi<0.0) reducedPair->fPhi = 2.0*TMath::Pi() + reducedPair->fPhi; // converted to [0,2pi]
1156 reducedPair->fPt = v0->Pt();
1157 reducedPair->fEta = v0->Eta();
1158 reducedPair->fLxy = v0->GetRr();
1159 reducedPair->fPointingAngle = v0->GetV0CosineOfPointingAngle(vtxKF->GetX(), vtxKF->GetY(), vtxKF->GetZ());
1160 reducedPair->fChisquare = v0->GetChi2V0();
1161 }
1162 else {
1163 const AliExternalTrackParam *negHelix=v0->GetParamN();
1164 const AliExternalTrackParam *posHelix=v0->GetParamP();
1165 if(!chargesAreCorrect) {
1166 negHelix = v0->GetParamP();
1167 posHelix = v0->GetParamN();
1168 }
1169 Int_t pdgPos = 211;
1170 if(id==AliReducedPair::kLambda0ToPPi) pdgPos = 2212;
1171 if(id==AliReducedPair::kGammaConv) pdgPos = -11;
1172 Int_t pdgNeg = -211;
1173 if(id==AliReducedPair::kALambda0ToPPi) pdgNeg = -2212;
1174 if(id==AliReducedPair::kGammaConv) pdgNeg = 11;
1175 AliKFParticle negKF(*(negHelix), pdgPos);
1176 AliKFParticle posKF(*(posHelix), pdgNeg);
1177 AliKFParticle v0Refit;
1178 v0Refit += negKF;
1179 v0Refit += posKF;
1180 Double_t massFit=0.0, massErrFit=0.0;
1181 v0Refit.GetMass(massFit,massErrFit);
1182 reducedPair->fMass[0] = massFit;
1183 reducedPair->fPhi = v0Refit.GetPhi();
1184 if(reducedPair->fPhi<0.0) reducedPair->fPhi = 2.0*TMath::Pi() + reducedPair->fPhi; // converted to [0,2pi]
1185 reducedPair->fPt = v0Refit.GetPt();
1186 reducedPair->fEta = v0Refit.GetEta();
1187 reducedPair->fLxy = v0Refit.GetPseudoProperDecayTime(*vtxKF, massFit);
1188 Double_t deltaPos[3];
1189 deltaPos[0] = v0Refit.GetX() - vtxKF->GetX(); deltaPos[1] = v0Refit.GetY() - vtxKF->GetY(); deltaPos[2] = v0Refit.GetZ() - vtxKF->GetZ();
1190 Double_t momV02 = v0Refit.GetPx()*v0Refit.GetPx() + v0Refit.GetPy()*v0Refit.GetPy() + v0Refit.GetPz()*v0Refit.GetPz();
1191 Double_t deltaPos2 = deltaPos[0]*deltaPos[0] + deltaPos[1]*deltaPos[1] + deltaPos[2]*deltaPos[2];
1192 reducedPair->fPointingAngle = (deltaPos[0]*v0Refit.GetPx() + deltaPos[1]*v0Refit.GetPy() + deltaPos[2]*v0Refit.GetPz()) /
1193 TMath::Sqrt(momV02*deltaPos2);
1194 reducedPair->fChisquare = v0Refit.GetChi2();
1195 }
1196 return reducedPair;
1197}
1198
1199
1200//_________________________________________________________________________________
1201UChar_t AliAnalysisTaskReducedTree::EncodeTPCClusterMap(AliVParticle* track, Bool_t isAOD) {
1202 //
1203 // Encode the TPC cluster map into an UChar_t
1204 // Divide the 159 bits from the bit map into 8 groups of adiacent clusters
1205 // For each group enable its corresponding bit if in that group there are more clusters compared to
1206 // a threshold.
1207 //
1208 AliESDtrack* esdTrack=0x0;
1209 AliAODTrack* aodTrack=0x0;
1210 if(isAOD)
1211 aodTrack=static_cast<AliAODTrack*>(track);
1212 else
1213 esdTrack=static_cast<AliESDtrack*>(track);
1214
1215 const UChar_t threshold=5;
1216 TBits tpcClusterMap = (isAOD ? aodTrack->GetTPCClusterMap() : esdTrack->GetTPCClusterMap());
1217 UChar_t map=0;
1218 UChar_t n=0;
1219 UChar_t j=0;
1220 for(UChar_t i=0; i<8; ++i) {
1221 n=0;
1222 for(j=i*20; j<(i+1)*20 && j<159; ++j) n+=tpcClusterMap.TestBitNumber(j);
1223 if(n>=threshold) map |= (1<<i);
1224 }
1225 return map;
1226}
1227
1228
1229//_________________________________________________________________________________
1230Int_t AliAnalysisTaskReducedTree::GetSPDTrackletMultiplicity(AliVEvent* event, Float_t lowEta, Float_t highEta) {
1231 //
1232 // Count the number of SPD tracklets in a given eta range
1233 //
1234 if (!event) return -1;
1235
1236 Int_t nTracklets = 0;
1237 Int_t nAcc = 0;
1238
1239 if(event->IsA() == AliAODEvent::Class()) {
1240 AliAODTracklets *tracklets = ((AliAODEvent*)event)->GetTracklets();
1241 nTracklets = tracklets->GetNumberOfTracklets();
1242 for(Int_t nn=0; nn<nTracklets; ++nn) {
1243 Double_t theta = tracklets->GetTheta(nn);
1244 Double_t eta = -TMath::Log(TMath::Tan(theta/2.0));
1245 if(eta < lowEta) continue;
1246 if(eta > highEta) continue;
1247 ++nAcc;
1248 }
1249 } else if(event->IsA() == AliESDEvent::Class()) {
1250 nTracklets = ((AliESDEvent*)event)->GetMultiplicity()->GetNumberOfTracklets();
1251 for(Int_t nn=0; nn<nTracklets; ++nn) {
1252 Double_t eta = ((AliESDEvent*)event)->GetMultiplicity()->GetEta(nn);
1253 if(eta < lowEta) continue;
1254 if(eta > highEta) continue;
1255 ++nAcc;
1256 }
1257 } else return -1;
1258
1259 return nAcc;
1260}
1261
1262
1263//_________________________________________________________________________________
1264void AliAnalysisTaskReducedTree::FinishTaskOutput()
1265{
1266 //
1267 // Finish Task
1268 //
1269
1270 //fTreeFile->Write();
1271 //fTreeFile->Close();
1272}