update: rejection of EMCAL triggers.
[u/mrichter/AliRoot.git] / PWGLF / STRANGENESS / Hypernuclei / AliAnalysisTaskAntiHe4.cxx
CommitLineData
110805a7 1#include "TChain.h"
2#include "TTree.h"
3#include "TF1.h"
4#include "TH1F.h"
5#include "TH2F.h"
6#include "TH3F.h"
7#include "THn.h"
8#include "TCanvas.h"
9#include "AliESDtrackCuts.h"
10#include "AliESDpid.h"
11#include "AliESDVertex.h"
12#include "TFile.h"
13
14#include "AliAnalysisTask.h"
15#include "AliAnalysisManager.h"
16#include "AliCentrality.h"
17#include "AliESDEvent.h"
18#include "AliESDInputHandler.h"
19
d1f5b077 20#include "AliMCEventHandler.h"
21#include "AliMCEvent.h"
22#include "AliStack.h"
23
110805a7 24#include "AliAnalysisTaskAntiHe4.h"
25
26
27///////////////////////////////////////////////////////////
28// //
29// Analysis for the observation of anti-alpha particles. //
30// //
31///////////////////////////////////////////////////////////
32
33
34ClassImp(AliAnalysisTaskAntiHe4)
35
36//________________________________________________________________________
37AliAnalysisTaskAntiHe4::AliAnalysisTaskAntiHe4()
38 : AliAnalysisTaskSE(),
39 fEventHandler(0),
40 fESD(0),
41 fHistCentralityClass10(0),
42 fHistCentralityPercentile(0),
43 fHistTriggerStat(0),
44 fHistTriggerStatAfterEventSelection(0),
45 fHistDEDx(0),
46 fHistTOF3D(0),
47 fHistAlpha(0),
48 fHistAlphaSignal(0),
49 fGraphAlphaSignal(0),
50 fNCounter(0),
51 fHistDeDx(0),
52 fHistDeDxRegion(0),
53 fHistDeDxSharp(0),
54 fHistTOF2D(0),
55 fHistTOFnuclei(0x0),
56 fNTriggers(5),
57 fBBParametersLightParticles(),
58 fBBParametersNuclei(),
59 fMCtrue(0),
60 fTriggerFired(),
61 fESDtrackCuts(0),
62 fESDtrackCutsSharp(0),
63 fESDpid(0),
64 fAntiAlpha(0),
d1f5b077 65 fHistHelium4PtGen(0),
66 fHistHelium4PtGenPrim(0),
67 fHistHelium4PtGenSec(0),
68 fHistHelium4PtGenEta(0),
69 fHistHelium4PtGenPrimEta(0),
70 fHistAntiHelium4PtGen(0),
71 fHistAntiHelium4PtGenPrim(0),
72 fHistAntiHelium4PtGenSec(0),
73 fHistAntiHelium4PtGenEta(0),
74 fHistHelium4PtAso(0),
75 fHistHelium4PtAsoPrim(0),
76 fHistHelium4PtAsoSec(0),
77 fHistAntiHelium4PtAso(0),
110805a7 78 fTree(0),
79 fOutputContainer(0)
80{
81 // default Constructor
82
83
84 // Define input and output slots here
85}
86
87//________________________________________________________________________
88AliAnalysisTaskAntiHe4::AliAnalysisTaskAntiHe4(const char *name)
89 : AliAnalysisTaskSE(name),
90 fEventHandler(0),
91 fESD(0),
92 fHistCentralityClass10(0),
93 fHistCentralityPercentile(0),
94 fHistTriggerStat(0),
95 fHistTriggerStatAfterEventSelection(0),
96 fHistDEDx(0),
97 fHistTOF3D(0),
98 fHistAlpha(0),
99 fHistAlphaSignal(0),
100 fGraphAlphaSignal(0),
101 fNCounter(0),
102 fHistDeDx(0),
103 fHistDeDxRegion(0),
104 fHistDeDxSharp(0),
105 fHistTOF2D(0),
106 fHistTOFnuclei(0x0),
107 fNTriggers(5),
108 fBBParametersLightParticles(),
109 fBBParametersNuclei(),
110 fMCtrue(0),
111 fTriggerFired(),
112 fESDtrackCuts(0),
113 fESDtrackCutsSharp(0),
114 fESDpid(0),
115 fAntiAlpha(0),
d1f5b077 116 fHistHelium4PtGen(0),
117 fHistHelium4PtGenPrim(0),
118 fHistHelium4PtGenSec(0),
119 fHistHelium4PtGenEta(0),
120 fHistHelium4PtGenPrimEta(0),
121 fHistAntiHelium4PtGen(0),
122 fHistAntiHelium4PtGenPrim(0),
123 fHistAntiHelium4PtGenSec(0),
124 fHistAntiHelium4PtGenEta(0),
125 fHistHelium4PtAso(0),
126 fHistHelium4PtAsoPrim(0),
127 fHistHelium4PtAsoSec(0),
128 fHistAntiHelium4PtAso(0),
110805a7 129 fTree(0),
130 fOutputContainer(0)
131{
132 // Constructor
133
134 // Define input and output slots here
135 // Input slot #0 works with a TChain
136 DefineInput(0, TChain::Class());
137 // Output slot #0 writes into a TH1 container
138 DefineOutput(1, TObjArray::Class());
139 DefineOutput(2, TTree::Class());
140 //
141 // cuts for candidates
142 //
143 fESDtrackCuts = new AliESDtrackCuts("AliESDtrackCuts","AliESDtrackCuts");
144 //
145 fESDtrackCuts->SetAcceptKinkDaughters(kFALSE);
146 fESDtrackCuts->SetMinNClustersTPC(70);
147 fESDtrackCuts->SetMaxChi2PerClusterTPC(5);
148 fESDtrackCuts->SetMaxDCAToVertexXY(3);
149 fESDtrackCuts->SetMaxDCAToVertexZ(2);
150 fESDtrackCuts->SetRequireTPCRefit(kTRUE);
151 //fESDtrackCuts->SetRequireITSRefit(kTRUE);
152 fESDtrackCuts->SetMinNClustersITS(2);
153 fESDtrackCuts->SetEtaRange(-0.8,0.8);
154 //
155 // cuts for final plots
156 //
157 fESDtrackCutsSharp = new AliESDtrackCuts("AliESDtrackCuts","AliESDtrackCuts");
158 fESDtrackCutsSharp->SetAcceptKinkDaughters(kFALSE);
159 fESDtrackCutsSharp->SetMinNClustersTPC(80);
160 fESDtrackCutsSharp->SetMaxChi2PerClusterITS(10);// TO BE INVESTIGATED !!!!!!!!!!!!!!
161 fESDtrackCutsSharp->SetMaxChi2PerClusterTPC(5);
162 fESDtrackCutsSharp->SetRequireTPCRefit(kTRUE);
163 fESDtrackCutsSharp->SetRequireITSRefit(kTRUE);
164 fESDtrackCutsSharp->SetMinNClustersITS(2);
165 fESDtrackCutsSharp->SetMaxDCAToVertexXY(0.1);
166 fESDtrackCutsSharp->SetMaxDCAToVertexZ(0.5);
167 fESDtrackCutsSharp->SetEtaRange(-0.8,0.8);
168
169 //ESD Track cuts from TestFilterRawTask
170 /* fESDtrackCuts = new AliESDtrackCuts("AliESDtrackCuts","AliESDtrackCuts");
171 fESDtrackCuts->SetMinNClustersTPC(80);
172 fESDtrackCuts->SetMinNClustersITS(2);
173 fESDtrackCuts->SetAcceptKinkDaughters(kFALSE);
174 fESDtrackCuts->SetMaxChi2PerClusterTPC(5);
175 fESDtrackCuts->SetRequireTPCRefit(kTRUE);
176 fESDtrackCuts->SetRequireITSRefit(kTRUE);
177 fESDtrackCuts->SetEtaRange(-1,1);
178 fESDtrackCuts->SetMaxDCAToVertexXY(1);
179 fESDtrackCuts->SetMaxDCAToVertexZ(2);
180 //test strickter cuts
181 //fESDtrackCuts->SetMaxDCAToVertexXY(0.1);
182 //fESDtrackCuts->SetMaxDCAToVertexZ(0.5);
183 //fESDtrackCuts->SetEtaRange(-0.8,0.8);
184 */
185
186 fMCtrue = kTRUE;
187
188}
189
190//________________________________________________________________________
191void AliAnalysisTaskAntiHe4::UserCreateOutputObjects()
192{
193 // Create histograms
194 // Called once
195 //
196 //histogram to count number of events
197 //
198 fHistCentralityClass10 = new TH1F("fHistCentralityClass10", "centrality with class10", 11, -0.5, 10.5);
199 fHistCentralityClass10->GetXaxis()->SetTitle("Centrality");
200 fHistCentralityClass10->GetYaxis()->SetTitle("Entries");
201 //
202 fHistCentralityPercentile = new TH1F("fHistCentralityPercentile", "centrality with percentile", 101, -0.1, 100.1);
203 fHistCentralityPercentile->GetXaxis()->SetTitle("Centrality");
204 fHistCentralityPercentile->GetYaxis()->SetTitle("Entries");
205 //
206 //trigger statitics histogram
207 //
208 fHistTriggerStat = new TH1F("fHistTriggerStat","Trigger statistics", fNTriggers,-0.5,fNTriggers-0.5);
209 const Char_t* aTriggerNames[] = { "kMB", "kCentral", "kSemiCentral", "kEMCEJE", "kEMCEGA" };
210 for ( Int_t ii=0; ii < fNTriggers; ii++ )
211 fHistTriggerStat->GetXaxis()->SetBinLabel(ii+1, aTriggerNames[ii]);
212
213 fHistTriggerStatAfterEventSelection = new TH1F("fHistTriggerStatAfterEventSelection","Trigger statistics after event selection", fNTriggers,-0.5,fNTriggers-0.5);
214 for ( Int_t ii=0; ii < fNTriggers; ii++ )
215 fHistTriggerStatAfterEventSelection->GetXaxis()->SetBinLabel(ii+1, aTriggerNames[ii]);
216
217 //dE/dx performance
218 fHistDEDx= new TH3F("fHistDEDx","DEDx",1000,0.01,100,1000,1,2000,2,-2,2);
219 BinLogAxis(fHistDEDx, 0);
220 BinLogAxis(fHistDEDx, 1);
221 fHistDEDx->GetXaxis()->SetTitle("p_{tot}/sign");
222 fHistDEDx->GetYaxis()->SetTitle("TPC signal");
223
224 fHistDeDx = new TH2F("fHistDeDx", "dE/dx", 1000, 0.1, 6.0, 1000, 0.0, 1000);
225 fHistDeDx->GetYaxis()->SetTitle("TPC dE/dx signal (a.u.)");
226 fHistDeDx->GetXaxis()->SetTitle("#frac{p}{z} (GeV/c)");
227 BinLogAxis(fHistDeDx);
228
229 fHistDeDxRegion = new TH3F("fHistDeDxRegion", "dE/dx", 400, 0., 6.0, 300, 0., 3., 4, -0.5, 3.5);
230 fHistDeDxRegion->GetYaxis()->SetTitle("TPC dE/dx signal (a.u.)");
231 fHistDeDxRegion->GetXaxis()->SetTitle("#frac{p}{z} (GeV/c)");
232
233 fHistDeDxSharp = new TH2F("fHistDeDxSharp", "dE/dx", 1000, 0.1, 6.0, 1000, 0.0, 1000);
234 fHistDeDxSharp->GetYaxis()->SetTitle("TPC dE/dx signal (a.u.)");
235 fHistDeDxSharp->GetXaxis()->SetTitle("#frac{p}{z} (GeV/c)");
236 BinLogAxis(fHistDeDxSharp);
237
238 //TOF performance
239 fHistTOF3D = new TH3F("fHistTOF3D","mass determination from TOF",400,0.25,4.5,2,-0.5,1.5,2,-1,1);
240
241 fHistTOF2D = new TH2F("fHistTOF2D", "TOF2D", 500, 0.0, 10., 2250, 0.2, 1.1);
242 fHistTOF2D->GetYaxis()->SetTitle("#beta");
243 fHistTOF2D->GetXaxis()->SetTitle("p (GeV/c)");
244
245 fHistTOFnuclei = new TH2F("fHistTOFnuclei", "TOF", 500, 0.0, 10., 2250, 0.7, 1.);
246 fHistTOFnuclei->GetYaxis()->SetTitle("#beta");
247 fHistTOFnuclei->GetXaxis()->SetTitle("#frac{p}{z} (GeV/c)");
248
249 //alphas
250 fHistAlpha = new TH1F("fHistAlpha", "Anti-Alpha", 22, 1.12, 4.31);
251 fHistAlpha->GetYaxis()->SetTitle("Counts");
252 fHistAlpha->GetXaxis()->SetTitle("#frac{m^{2}}{z^{2}} (GeV^{2}/c^{4})");
253
254 fHistAlphaSignal = new TH1F("fHistAlphaSignal", "Anti-Alpha", 22, 1.12, 4.31);
255 fHistAlphaSignal->GetYaxis()->SetTitle("Counts");
256 fHistAlphaSignal->GetXaxis()->SetTitle("#frac{m^{2}}{z^{2}} (GeV^{2}/c^{4})");
257
258 fGraphAlphaSignal = new TGraph(20);
259 fNCounter = 0;
260 //
261 // (0.) dca, (1.) sign, (2.) particle Type, (3.) p_tot
262 //
263 TString axisNameMult[4]={"dca","sign","particleType","ptot"};
264 TString axisTitleMult[4]={"dca","sign","particleType","ptot"};
265 const Int_t kDcaBins = 76;
266 Double_t binsDca[kDcaBins+1] = {-3,-2.85,-2.7,-2.55,-2.4,-2.25,-2.1,-1.95,-1.8,-1.65,-1.5,-1.35,-1.2,-1.05,-0.9,-0.75,-0.6,-0.45,-0.3,-0.285,-0.27,-0.255,-0.24,-0.225,-0.21,-0.195,-0.18,-0.165,-0.15,-0.135,-0.12,-0.105,-0.09,-0.075,-0.06,-0.045,-0.03,-0.015,0,0.015,0.03,0.045,0.06,0.075,0.09,0.105,0.12,0.135,0.15,0.165,0.18,0.195,0.21,0.225,0.24,0.255,0.27,0.285,0.3,0.45,0.6,0.75,0.9,1.05,1.2,1.35,1.5,1.65,1.8,1.95,2.1,2.25,2.4,2.55,2.7,2.85,3};
267 //
268 // (0.) dca, (1.) sign, (2.) particle Type, (3.) p_tot
269 Int_t binsAntiAlpha[4] = {77, 2, 4, 100};
270 Double_t xminAntiAlpha[4] = {-3, -2, -0.5, 0};
271 Double_t xmaxAntiAlpha[4] = { 3, 2, 3.5, 6};
272 //
273 fAntiAlpha = new THnF("fAntiAlpha", "AntiAlpha; (0.) dca, (1.) sign, (2.) particle Type, (3.) p_tot", 4, binsAntiAlpha, xminAntiAlpha, xmaxAntiAlpha);
274 //
275 fAntiAlpha->GetAxis(0)->Set(kDcaBins, binsDca);
276 for (Int_t iaxis=0; iaxis<4;iaxis++){
277 fAntiAlpha->GetAxis(iaxis)->SetName(axisNameMult[iaxis]);
278 fAntiAlpha->GetAxis(iaxis)->SetTitle(axisTitleMult[iaxis]);
279 }
280 //
d1f5b077 281 //
282 //Create histograms for MC
283 //generated histogramms
284 fHistHelium4PtGen = new TH1F("fHistHelium4PtGen", "Generated ^{4}He", 200, 0.0, 10.0);
285 fHistHelium4PtGen->GetYaxis()->SetTitle("Counts");
286 fHistHelium4PtGen->GetXaxis()->SetTitle("p_{T} (GeV/c)");
287
288 fHistHelium4PtGenPrim = new TH1F("fHistHelium4PtGenPrim", "Primary generated ^{4}He", 200, 0.0, 10.0);
289 fHistHelium4PtGenPrim->GetYaxis()->SetTitle("Counts");
290 fHistHelium4PtGenPrim->GetXaxis()->SetTitle("p_{T} (GeV/c)");
291
292 fHistHelium4PtGenSec = new TH1F("fHistHelium4PtGenSec", "Secondary generated ^{4}He", 200, 0.0, 10.0);
293 fHistHelium4PtGenSec->GetYaxis()->SetTitle("Counts");
294 fHistHelium4PtGenSec->GetXaxis()->SetTitle("p_{T} (GeV/c)");
295
296 fHistHelium4PtGenEta = new TH1F("fHistHelium4PtGenEta", "Generated ^{4}He in #eta < 0.8", 200, 0.0, 10.0);
297 fHistHelium4PtGenEta->GetYaxis()->SetTitle("Counts");
298 fHistHelium4PtGenEta->GetXaxis()->SetTitle("p_{T} (GeV/c)");
299
300 fHistHelium4PtGenPrimEta = new TH1F("fHistHelium4PtGenPrimEta", "Primary generated ^{4}He in #eta < 0.8", 200, 0.0, 10.0);
301 fHistHelium4PtGenPrimEta->GetYaxis()->SetTitle("Counts");
302 fHistHelium4PtGenPrimEta->GetXaxis()->SetTitle("p_{T} (GeV/c)");
303
304 fHistAntiHelium4PtGen = new TH1F("fHistAntiHelium4PtGen", "Generated #bar{^{4}He}", 200, 0.0, 10.0);
305 fHistAntiHelium4PtGen->GetYaxis()->SetTitle("Counts");
306 fHistAntiHelium4PtGen->GetXaxis()->SetTitle("p_{T} (GeV/c)");
307
308 fHistAntiHelium4PtGenPrim = new TH1F("fHistAntiHelium4PtGenPrim", "Primary generated #bar{^{4}He}", 200, 0.0, 10.0);
309 fHistAntiHelium4PtGenPrim->GetYaxis()->SetTitle("Counts");
310 fHistAntiHelium4PtGenPrim->GetXaxis()->SetTitle("p_{T} (GeV/c)");
311
312 fHistAntiHelium4PtGenSec = new TH1F("fHistAntiHelium4PtGenSec", "Secondary generated #bar{^{4}He}", 200, 0.0, 10.0);
313 fHistAntiHelium4PtGenSec->GetYaxis()->SetTitle("Counts");
314 fHistAntiHelium4PtGenSec->GetXaxis()->SetTitle("p_{T} (GeV/c)");
315
316 fHistAntiHelium4PtGenEta = new TH1F("fHistAntiHelium4PtGenEta", "Generated #bar{^{4}He} in #eta < 0.8", 200, 0.0, 10.0);
317 fHistAntiHelium4PtGenEta->GetYaxis()->SetTitle("Counts");
318 fHistAntiHelium4PtGenEta->GetXaxis()->SetTitle("p_{T} (GeV/c)");
319
320 //associated histograms
321 fHistHelium4PtAso = new TH1F("fHistHelium4PtAso", "Associated ^{4}He", 200, 0.0, 10.0);
322 fHistHelium4PtAso->GetYaxis()->SetTitle("Counts");
323 fHistHelium4PtAso->GetXaxis()->SetTitle("p_{T} (GeV/c)");
324
325 fHistHelium4PtAsoPrim = new TH1F("fHistHelium4PtAsoPrim", "Associated prim ^{4}He", 200, 0.0, 10.0);
326 fHistHelium4PtAsoPrim->GetYaxis()->SetTitle("Counts");
327 fHistHelium4PtAsoPrim->GetXaxis()->SetTitle("p_{T} (GeV/c)");
328
329 fHistHelium4PtAsoSec = new TH1F("fHistHelium4PtAsoSec", "Associated sec ^{4}He", 200, 0.0, 10.0);
330 fHistHelium4PtAsoSec->GetYaxis()->SetTitle("Counts");
331 fHistHelium4PtAsoSec->GetXaxis()->SetTitle("p_{T} (GeV/c)");
332
333 fHistAntiHelium4PtAso = new TH1F("fHistAntiHelium4PtAso", "Associated #bar{^{4}He}", 200, 0.0, 10.0);
334 fHistAntiHelium4PtAso->GetYaxis()->SetTitle("Counts");
335 fHistAntiHelium4PtAso->GetXaxis()->SetTitle("p_{T} (GeV/c)");
336
110805a7 337 //
338 //
339 fOutputContainer = new TObjArray(1);
340 fOutputContainer->SetOwner(kTRUE);
341 fOutputContainer->SetName(GetName());
342 fOutputContainer->Add(fHistCentralityClass10);
343 fOutputContainer->Add(fHistCentralityPercentile);
344 fOutputContainer->Add(fHistTriggerStat);
345 fOutputContainer->Add(fHistTriggerStatAfterEventSelection);
346 fOutputContainer->Add(fHistDEDx);
347 fOutputContainer->Add(fHistDeDx);
348 fOutputContainer->Add(fHistDeDxRegion);
349 fOutputContainer->Add(fHistDeDxSharp);
350 fOutputContainer->Add(fAntiAlpha);
351 fOutputContainer->Add(fHistAlpha);
352 fOutputContainer->Add(fHistAlphaSignal);
353 fOutputContainer->Add(fGraphAlphaSignal);
354 fOutputContainer->Add(fHistTOF3D);
355 fOutputContainer->Add(fHistTOF2D);
356 fOutputContainer->Add(fHistTOFnuclei);
d1f5b077 357 fOutputContainer->Add(fHistHelium4PtGen);
358 fOutputContainer->Add(fHistHelium4PtGenPrim);
359 fOutputContainer->Add(fHistHelium4PtGenSec);
360 fOutputContainer->Add(fHistHelium4PtGenEta);
361 fOutputContainer->Add(fHistHelium4PtGenPrimEta);
362 fOutputContainer->Add(fHistAntiHelium4PtGen);
363 fOutputContainer->Add(fHistAntiHelium4PtGenPrim);
364 fOutputContainer->Add(fHistAntiHelium4PtGenSec);
365 fOutputContainer->Add(fHistAntiHelium4PtGenEta);
366 fOutputContainer->Add(fHistHelium4PtAso);
367 fOutputContainer->Add(fHistHelium4PtAsoPrim);
368 fOutputContainer->Add(fHistHelium4PtAsoSec);
369 fOutputContainer->Add(fHistAntiHelium4PtAso);
110805a7 370
371
372 //------------ Tree and branch definitions ----------------//
373 fTree = new TTree("tree", " alpha tree");
374
375 //------------ Event variables ------------//
376 fTree->Branch("Name",Name,"Name/C");
377 fTree->Branch("Evnt",&evnt, "evnt/I");
378 fTree->Branch("itrk", &itrk, "itrk/I");
379
380 //-------------------------------------------//
381 //----------- Track variables --------------//
382
383 fTree->Branch("TrkPtot",TrkPtot,"TrkPtot[itrk]/F");
384 fTree->Branch("TPCPtot",TPCPtot,"TPCPtot[itrk]/F");
385 fTree->Branch("DeDx",DeDx,"DeDx[itrk]/F");
386 fTree->Branch("Sign",Sign,"Sign[itrk]/F");
387 fTree->Branch("DCAXY",DCAXY,"DCAXY[itrk]/F");
388 fTree->Branch("DCAZ",DCAZ,"DCAZ[itrk]/F");
389 fTree->Branch("ITSnCluster",ITSnCluster,"ITSnCluster[itrk]/F");
390 fTree->Branch("TPCNsignal",TPCNsignal,"TPCNsignal[itrk]/F");
391 fTree->Branch("Mass",Mass,"Mass[itrk]/F");
392 //
393 fTree->Branch("ITSRefit",ITSRefit,"ITSRefit[itrk]/F");
394 fTree->Branch("TOFtime",TOFtime,"TOFtime[itrk]/F");
395 fTree->Branch("TOFRefit",TOFRefit,"TOFRefit[itrk]/F");
396 fTree->Branch("TOFout",TOFout,"TOFout[itrk]/F");
397 //
398 fTree->Branch("ITSsignal",ITSsignal,"ITSsignal[itrk]/F");
399 fTree->Branch("SharedClusters",SharedClusters,"SharedClusters[itrk]/F");
400 fTree->Branch("fFileName",fFileName,"fFileName/C");
401 fTree->Branch("fEventNumber",fEventNumber,"fEventNumber/I");
d1f5b077 402 //
403 fTree->Branch("fAssociated",fAssociated,"fAssociated[itrk]/O");
404 fTree->Branch("fTrackPt",fTrackPt,"fTrackPt[itrk]/F");
110805a7 405
406}
407
408//________________________________________________________________________
409void AliAnalysisTaskAntiHe4::UserExec(Option_t *)
410{
411 // Main loop
412 // Called for each event
413
414 //get Event-Handler for the trigger information
415 fEventHandler= dynamic_cast<AliInputEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
416 if (!fEventHandler) {
417 AliError("Could not get InputHandler");
418 //return -1;
419 return;
420 }
d1f5b077 421
422 // Monte Carlo
423 AliMCEventHandler* eventHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
424 if (!eventHandler){
425 //Printf("ERROR: Could not retrieve MC event handler");
426 fMCtrue = kFALSE;
427 }
428
429 AliMCEvent* mcEvent = 0x0;
430 AliStack* stack = 0x0;
431 if (eventHandler) mcEvent = eventHandler->MCEvent();
432 if (!mcEvent){
433 //Printf("ERROR: Could not retrieve MC event");
434 if (fMCtrue) return;
435 }
436
437 if (fMCtrue){
438 stack = mcEvent->Stack();
439 if (!stack) return;
440 }
441
442 //look for the generated particles
443 MCGenerated(stack);
110805a7 444
445 fESD=dynamic_cast<AliESDEvent*>(InputEvent());
446 if (!fESD) {
447 //Printf("ERROR: fESD not available");
448 return;
449 }
450
451 if (SetupEvent() < 0) {
452 PostData(1, fOutputContainer);
453 return;
454 }
455
456 const AliESDVertex *vertex = fESD->GetPrimaryVertexTracks();
457 if(vertex->GetNContributors()<1) {
458 // SPD vertex
459 vertex = fESD->GetPrimaryVertexSPD();
460 if(vertex->GetNContributors()<1) {
461 PostData(1, fOutputContainer);
462 return;
463 }
464
465 }
466 // check if event is selected by physics selection class
467 //
468 Bool_t isSelected = kFALSE;
469 isSelected = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
470 if (!isSelected || TMath::Abs(vertex->GetZv()) > 10) {
471 PostData(1, fOutputContainer);
472 return;
473 }
474
475 //
476 //centrality selection in PbPb
477 //
478 Int_t centralityClass10 = -1;
479 Int_t centralityPercentile = -1;
480 //
481 if (fESD->GetEventSpecie() == 4) { //species == 4 == PbPb
482 // PbPb
483 AliCentrality *esdCentrality = fESD->GetCentrality();
484 centralityClass10 = esdCentrality->GetCentralityClass10("V0M"); // centrality percentile determined with V0
485 centralityPercentile = esdCentrality->GetCentralityPercentile("V0M"); // centrality percentile determined with V0
486 if(!fMCtrue){
487 if (centralityPercentile < 0. || centralityPercentile > 80. ) return; //select only events with centralities between 0 and 80 %
488 }
489 }
490 //
d1f5b077 491 if (!fTriggerFired[0] && !fTriggerFired[1] && !fTriggerFired[2]) return; // select only events which pass kMB, kCentral, kSemiCentral
492 //
110805a7 493 fHistCentralityClass10->Fill(centralityClass10);
494 fHistCentralityPercentile->Fill(centralityPercentile);
495 //
496 if(fTriggerFired[0] == kTRUE) fHistTriggerStatAfterEventSelection->Fill(0);
497 if(fTriggerFired[1] == kTRUE) fHistTriggerStatAfterEventSelection->Fill(1);
498 if(fTriggerFired[2] == kTRUE) fHistTriggerStatAfterEventSelection->Fill(2);
499 if(fTriggerFired[3] == kTRUE) fHistTriggerStatAfterEventSelection->Fill(3);
500 if(fTriggerFired[4] == kTRUE) fHistTriggerStatAfterEventSelection->Fill(4);
501 //
502 //
503 if (!fESDpid) fESDpid = ((AliESDInputHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->GetESDpid();
504 if (!fESDpid) {
505 fESDpid = new AliESDpid(); // HACK FOR MC PBPB --> PLEASE REMOVE AS SOON AS POSSIBLE
506 fESDpid->GetTPCResponse().SetBetheBlochParameters(1.28778e+00/50., 3.13539e+01, TMath::Exp(-3.16327e+01), 1.87901e+00, 6.41583e+00);
507 }
508 //
509 Float_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z for th // for Anti-Alpha
510 evnt = fESD->GetEventNumberInFile();
511 sscanf(fInputHandler->GetTree()->GetCurrentFile()->GetName(),"%s", Name);
512 itrk = 0;
513 //
514 Int_t runNumber = 0;
515 runNumber = fESD->GetRunNumber();
516 //
517 Bool_t fillTree = kFALSE;
518 // Track loop to fill the spectram
519 for (Int_t iTracks = 0; iTracks < fESD->GetNumberOfTracks(); iTracks++) {
520 //
521 AliESDtrack* track = dynamic_cast<AliESDtrack*>(fESD->GetTrack(iTracks));
522 if (!fESDtrackCuts->AcceptTrack(track)) continue;
523 //
524 Double_t nClustersTPCPID=track->GetTPCsignalN();
525 if(nClustersTPCPID < 60) continue;
526 //
527 if (!track->GetInnerParam()) continue;
528 //
529 Double32_t signal[4] = {0,0,0,0};
530 Char_t ncl[3];
531 Char_t nrows[3];
532 if (track->GetTPCdEdxInfo()) {
533 track->GetTPCdEdxInfo()->GetTPCSignalRegionInfo(signal,ncl,nrows);
534 }
535 //
536 Double_t ptot = track->GetInnerParam()->GetP();
537 Double_t ptotInc = track->GetP(); // total momentum of the incoming particle
538 Double_t sign = track->GetSign();
539 //
540 Double_t eta = TMath::Abs(track->Eta());
541 TBits shared = track->GetTPCSharedMap();
542
543 track->GetImpactParameters(dca, cov);
544 Float_t dcaXY = TMath::Sqrt(dca[0]*dca[0]);
545 Float_t dcaXYsign = dca[0];
546 Float_t dcaZ = TMath::Sqrt(dca[1]*dca[1]);
547 //
548 //
549 Double_t tpcSignal = track->GetTPCsignal();
550 //
551 //PID via specific energy loss in the TPC
552 //define the arrays for the Bethe-Bloch-Parameters
553 //
554
555 SetBBParameters(runNumber);
556
557 //define expected signals for the various species
558 Double_t expSignalProton = AliExternalTrackParam::BetheBlochAleph(ptot/0.93827,fBBParametersLightParticles[0],fBBParametersLightParticles[1],fBBParametersLightParticles[2],fBBParametersLightParticles[3],fBBParametersLightParticles[4]);
559 Double_t expSignalDeuteron = AliExternalTrackParam::BetheBlochAleph(ptot/1.875612,fBBParametersLightParticles[0],fBBParametersLightParticles[1],fBBParametersLightParticles[2],fBBParametersLightParticles[3],fBBParametersLightParticles[4]);
560 Double_t expSignalTriton = AliExternalTrackParam::BetheBlochAleph(ptot/2.808921,fBBParametersLightParticles[0],fBBParametersLightParticles[1],fBBParametersLightParticles[2],fBBParametersLightParticles[3],fBBParametersLightParticles[4]);
561
562 Double_t expSignalHelium3 = 4*AliExternalTrackParam::BetheBlochAleph(2*ptot/2.80941,fBBParametersNuclei[0],fBBParametersNuclei[1],fBBParametersNuclei[2],fBBParametersNuclei[3],fBBParametersNuclei[4]);
563 Double_t expSignalHelium4 = 4*AliExternalTrackParam::BetheBlochAleph(2*ptot/3.728401,fBBParametersNuclei[0],fBBParametersNuclei[1],fBBParametersNuclei[2],fBBParametersNuclei[3],fBBParametersNuclei[4]);
564
565 //
566 Float_t mass = 0;
567 Float_t time = -1;
568 Float_t beta = 0;
569 Float_t gamma = -1;
570 Bool_t hasTOFout = kFALSE;
571 Bool_t hasTOFtime = kFALSE;
572 Float_t length = track->GetIntegratedLength();
573 UInt_t status = track->GetStatus();
d1f5b077 574 Bool_t isAssociated = kFALSE;
110805a7 575
576 if (length > 350) {
577 time = track->GetTOFsignal() - fESDpid->GetTOFResponse().GetStartTime(track->P());
578 if (time > 0) {
579 beta = length / (2.99792457999999984e-02 * time);
580 gamma = 1/TMath::Sqrt(1 - beta*beta);
581 mass = ptot/TMath::Sqrt(gamma*gamma - 1); // using inner TPC mom. as approx.
582 }
583 }
584 //
585 // fill tree and print candidates (for short list)
586 //
587 Float_t cut = 4*AliExternalTrackParam::BetheBlochAleph(2*ptot/(0.938*3),1.1931,
588 31.9806,
589 5.04114e-11,
590 2.13096,
591 2.38541);
d1f5b077 592 if (fMCtrue) cut = 4*AliExternalTrackParam::BetheBlochAleph(2*ptot/(0.938*3),0.9931,
593 31.9806,
594 5.04114e-11,
595 2.13096,
596 2.38541);
110805a7 597 if (eta < 0.8 && tpcSignal > 120 && tpcSignal > cut && tpcSignal < 1000 && track->GetTPCsignalN() > 60 && dcaZ < 15 && dcaXY < 15 && ptot > 1.0 && ptot < 20) {
598 //
d1f5b077 599 cout << "AntiAlphaEvent" << " "
110805a7 600 << AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()->GetTree()->GetCurrentFile()->GetName() << " "
601 << "event number in file: " << fESD->GetEventNumberInFile()
602 << " Index " << iTracks
603 << " ptot: " << ptot
604 << " sig: " << tpcSignal <<endl;
605 //
606 fillTree = kTRUE;
607 //
608 TrkPtot[itrk] = track->P();
609 TPCPtot[itrk] = ptot;
610 DeDx[itrk] = tpcSignal;
611 DCAZ[itrk] = dcaZ;
612 TPCNsignal[itrk] = track->GetTPCsignalN();
613 ITSnCluster[itrk] = track->GetNcls(0);
614 Sign[itrk] = sign;
615 DCAXY[itrk] = dcaXY;
616 Mass[itrk] = mass;
617 //
618 ITSsignal[itrk] = track->GetITSsignal();
619 SharedClusters[itrk] = shared.CountBits();
620 sscanf(fInputHandler->GetTree()->GetCurrentFile()->GetName(),"%s", fFileName);
621 //fFileName[itrk] = fInputHandler->GetTree()->GetCurrentFile()->GetName();
622 fEventNumber[itrk] = fESD->GetEventNumberInFile();
623 //
624 if(status&AliESDtrack::kITSrefit)
625 ITSRefit[itrk] = 1;
626 else ITSRefit[itrk] = 0;
627 //
628 if (time < 99998) {
629 TOFRefit[itrk] = 1;
630 } else {
631 TOFRefit[itrk] = 0;
632 }
d1f5b077 633 //itrk++;
110805a7 634 //
635 hasTOFout = status&AliESDtrack::kTOFout;
636 hasTOFtime = status&AliESDtrack::kTIME;
637 //
638 TOFtime[itrk] = hasTOFtime;
639 TOFout[itrk] = hasTOFout;
d1f5b077 640
641 if (fMCtrue){ //associated
642
643 Int_t label = track->GetLabel();
644 TParticle *tparticle = stack->Particle(TMath::Abs(label));
645
646 Bool_t isPrimary = stack->IsPhysicalPrimary(TMath::Abs(label));
647 Bool_t isSecondary = stack->IsSecondaryFromMaterial(TMath::Abs(label));
648
649 Long_t pdgCode = tparticle->GetPdgCode();
650 Double_t pT =(track->Pt())*2;
651
652 if(pdgCode == 1000020040)
653 {
654 fHistHelium4PtAso->Fill(pT);
655 if(isPrimary) fHistHelium4PtAsoPrim->Fill(pT);
656 if(isSecondary) fHistHelium4PtAsoSec->Fill(pT);
657 isAssociated = kTRUE;
658 }
659
660 if(pdgCode == -1000020040)
661 {
662 fHistAntiHelium4PtAso->Fill(pT);
663 isAssociated = kTRUE;
664 }
665
666 }
667
668 fAssociated[itrk] = isAssociated;
669 fTrackPt[itrk] = track->Pt();
670
110805a7 671 itrk++;
672 }
673 //
674 // do pid fill histogram for raw ratios
675 //
676 // (0.) dca, (1.) sign, (2.) particle Type, (3.) p_tot
677 Int_t id = -1;
678 if (ptot > 0.2 && TMath::Abs(tpcSignal - expSignalProton)/expSignalProton < 0.2) id = 0;
679 if (ptot > 0.3 && TMath::Abs(tpcSignal - expSignalDeuteron)/expSignalDeuteron < 0.2) id = 1;
680 if (ptot > 0.7 && TMath::Abs(tpcSignal - expSignalTriton)/expSignalTriton < 0.2) id = 2;
681 if (ptot > 0.5 && (tpcSignal - expSignalHelium3)/expSignalHelium3 > -0.1 && (tpcSignal - expSignalHelium3)/expSignalHelium3 < 0.2) id = 3;
682 //
683 Double_t vecAntiAlpha[4] = {dcaXYsign, sign, id, ptotInc};
684 if (id != -1 && tpcSignal > 120) fAntiAlpha->Fill(vecAntiAlpha);
685 //
686 // fill final histograms
687 //
688 if (!fESDtrackCutsSharp->AcceptTrack(track) || shared.CountBits() > 1 || track->GetTPCsignalN() < 80 || track->GetKinkIndex(0) != 0 || track->GetTPCNclsIter1() < 80) continue;
689 //
690 fHistDEDx->Fill(ptot,track->GetTPCsignal(), sign);
691 if (TMath::Abs(tpcSignal - expSignalHelium4)/expSignalHelium4 < 0.2 && time < 99998) fHistTOF3D->Fill(mass,1,sign);
692 //
693 // remove overlapping tracks with special dE/dx cut
694 //
695 //if (signal[1] < 6) continue;
696 //if (signal[0]/signal[1] > 1.6 || signal[2]/signal[1] > 1.6 || signal[3]/signal[1] > 1.3) continue;
697 //
698 if(sign<0) {
699 fHistDeDx->Fill(ptot, track->GetTPCsignal());
700 if (track->GetTPCsignalN() > 100 &&
701 TMath::Abs(track->Eta()) < 0.8 &&
702 signal[3]/signal[1] > 0.6 &&
703 signal[0]/signal[1] > 0.5 &&
704 signal[3]/signal[1] < 1.2 &&
705 track->GetTPCNclsIter1() > 70 &&
706 track->GetTPCnclsS() < 10) {
707 fHistDeDxSharp->Fill(ptot, track->GetTPCsignal());
708 }
709 }
710 //
711 // dE/dx for specific regions
712 //
713 for(Int_t iSig = 0; iSig < 4; iSig++) {
714 if (signal[1] > 6) fHistDeDxRegion->Fill(ptot,signal[iSig]/signal[1],iSig);
715 }
716 //
717 // alpha TOF plot
718 //
719 if((track->GetTPCsignal()-expSignalHelium4)/expSignalHelium4 > -0.15 && (track->GetTPCsignal()-expSignalHelium4)/expSignalHelium4 < 0.15) {
720 //TOF
721 hasTOFout = status&AliESDtrack::kTOFout;
722 hasTOFtime = status&AliESDtrack::kTIME;
723 Bool_t hasTOF = kFALSE;
724
725 if (hasTOFout && hasTOFtime) hasTOF = kTRUE;
726 //
727 if (length < 350.) hasTOF = kFALSE;
728
729 Float_t time0 = fESDpid->GetTOFResponse().GetStartTime(track->P());//fESDpid->GetTOFResponse().GetTimeZero();
730 //
731 if (length > 350. && hasTOF == kTRUE && ptot < 4) {
732 time = track->GetTOFsignal() - time0;
733 if (time > 0) {
734 beta = length / (2.99792457999999984e-02 * time);
735 if (beta < 0.975) {
736 gamma = 1/TMath::Sqrt(1 - beta*beta);
737 mass = ptot/TMath::Sqrt(gamma*gamma - 1); // using inner TPC mom. as approx.
738 if (TMath::Sqrt(track->GetTOFsignalDz()*track->GetTOFsignalDz() + track->GetTOFsignalDx()*track->GetTOFsignalDx()) < 5.) {
739 fHistAlpha->Fill(mass*mass);
740 if (mass*mass > 3. && mass*mass < 4.) {
741 fHistAlphaSignal->Fill(mass*mass);
742 fGraphAlphaSignal->SetPoint(fNCounter, ptot, track->GetTPCsignal());
743 fNCounter++;
744 }
745 fHistTOFnuclei->Fill(ptot,beta);
746 }
747 }
748 }
749 }
750 }
751
752 }//end loop over tracks
d1f5b077 753
110805a7 754 if (fillTree) fTree->Fill();
755
756 // Post output data.
757 PostData(1, fOutputContainer);
758 PostData(2, fTree);
759}
760
761//________________________________________________________________________
762void AliAnalysisTaskAntiHe4::Terminate(const Option_t *)
763{
764 // Draw result to the screen
765 // Called once at the end of the query
766
767 //get output data and darw 'fHistPt'
768 if (!GetOutputData(0)) return;
769
770}
771
772
773//________________________________________________________________________
774void AliAnalysisTaskAntiHe4::BinLogAxis(const THn *h, Int_t axisNumber) {
775 //
776 // Method for the correct logarithmic binning of histograms
777 //
778 TAxis *axis = h->GetAxis(axisNumber);
779 int bins = axis->GetNbins();
780
781 Double_t from = axis->GetXmin();
782 Double_t to = axis->GetXmax();
783 Double_t *newBins = new Double_t[bins + 1];
784
785 newBins[0] = from;
786 Double_t factor = pow(to/from, 1./bins);
787
788 for (int i = 1; i <= bins; i++) {
789 newBins[i] = factor * newBins[i-1];
790 }
791 axis->Set(bins, newBins);
792 delete [] newBins;
793
794}
795
796//________________________________________________________________________
797void AliAnalysisTaskAntiHe4::BinLogAxis(const TH3 *h, Int_t axisNumber) {
798 //
799 // Method for the correct logarithmic binning of histograms
800 //
801 TAxis *axis = h->GetXaxis();
802 if (axisNumber == 1) axis = h->GetYaxis();
803 if (axisNumber == 2) axis = h->GetZaxis();
804 int bins = axis->GetNbins();
805
806 Double_t from = axis->GetXmin();
807 Double_t to = axis->GetXmax();
808 Double_t *newBins = new Double_t[bins + 1];
809
810 newBins[0] = from;
811 Double_t factor = pow(to/from, 1./bins);
812
813 for (int i = 1; i <= bins; i++) {
814 newBins[i] = factor * newBins[i-1];
815 }
816 axis->Set(bins, newBins);
817 delete [] newBins;
818
819}
820//________________________________________________________________________
821void AliAnalysisTaskAntiHe4::BinLogAxis(const TH1 *h) {
822 //
823 // Method for the correct logarithmic binning of histograms
824 //
825 TAxis *axis = h->GetXaxis();
826 int bins = axis->GetNbins();
827
828 Double_t from = axis->GetXmin();
829 Double_t to = axis->GetXmax();
830 Double_t *newBins = new Double_t[bins + 1];
831
832 newBins[0] = from;
833 Double_t factor = pow(to/from, 1./bins);
834
835 for (int i = 1; i <= bins; i++) {
836 newBins[i] = factor * newBins[i-1];
837 }
838 axis->Set(bins, newBins);
839 delete [] newBins;
840
841}
842//_____________________________________________________________________________
843Int_t AliAnalysisTaskAntiHe4::Initialize() {
844
845
846 // -- Reset Event
847 // ----------------
848 ResetEvent();
849
850 return 0;
851}
852//________________________________________________________________________
853Int_t AliAnalysisTaskAntiHe4::SetupEvent() {
854 // Setup Reading of event
855
856 ResetEvent();
857 // -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
858 // -- Get Trigger
859 // -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
860
861 //Bool_t isTriggered = IsTriggered();
862 IsTriggered();
863 return 0;
864}
865//_____________________________________________________________________________
866void AliAnalysisTaskAntiHe4::ResetEvent() {
867 // Reset event
868 // -- Reset QA
869 for (Int_t ii = 0; ii < fNTriggers; ++ii)
870 fTriggerFired[ii] = kFALSE;
871
872 return;
873}
874//________________________________________________________________________
875Bool_t AliAnalysisTaskAntiHe4::IsTriggered() {
876 // Check if Event is triggered and fill Trigger Histogram
877
878 if ((fEventHandler->IsEventSelected() & AliVEvent::kMB)) fTriggerFired[0] = kTRUE;
879 if ((fEventHandler->IsEventSelected() & AliVEvent::kCentral)) fTriggerFired[1] = kTRUE;
880 if ((fEventHandler->IsEventSelected() & AliVEvent::kSemiCentral)) fTriggerFired[2] = kTRUE;
881 if ((fEventHandler->IsEventSelected() & AliVEvent::kEMCEJE)) fTriggerFired[3] = kTRUE;
882 if ((fEventHandler->IsEventSelected() & AliVEvent::kEMCEGA)) fTriggerFired[4] = kTRUE;
883
884 Bool_t isTriggered = kFALSE;
885
886 for (Int_t ii=0; ii<fNTriggers; ++ii) {
887 if(fTriggerFired[ii]) {
888 isTriggered = kTRUE;
889 fHistTriggerStat->Fill(ii);
890 }
891 }
892
893 return isTriggered;
894}
895//________________________________________________________________________
896void AliAnalysisTaskAntiHe4::SetBBParameters(Int_t runNumber){
897
898 if(runNumber < 166500){ //LHC10h
899
900 fBBParametersLightParticles[0] = 1.45802;
901 fBBParametersLightParticles[1] = 27.4992;
902 fBBParametersLightParticles[2] = 4.00313e-15;
903 fBBParametersLightParticles[3] = 2.48485;
904 fBBParametersLightParticles[4] = 8.31768;
905
906 fBBParametersNuclei[0] = 1.69155;
907 fBBParametersNuclei[1] = 27.4992;
908 fBBParametersNuclei[2] = 4.00313e-15;
909 fBBParametersNuclei[3] = 2.48485;
910 fBBParametersNuclei[4] = 8.31768;
911
912 }
913
914 if(runNumber > 166500){ //LHC11h
915
916 fBBParametersLightParticles[0] = 4.69637;//1.11243;
917 fBBParametersLightParticles[1] = 7.51827;//26.1144;
918 fBBParametersLightParticles[2] = 0.0183746;//4.00313e-15;
919 fBBParametersLightParticles[3] = 2.60;//2.72969;
920 fBBParametersLightParticles[4] = 2.7;//9.15038;
921
922 fBBParametersNuclei[0] = 1.4906;
923 fBBParametersNuclei[1] = 27.9758;
924 fBBParametersNuclei[2] = 4.00313e-15;
925 fBBParametersNuclei[3] = 2.50804;
926 fBBParametersNuclei[4] = 8.31768;
927
928 }
929}
d1f5b077 930//______________________________________________________________________________
931void AliAnalysisTaskAntiHe4::MCGenerated(AliStack* stack)
932{
933
934 // Monte Carlo for genenerated particles
935 if (fMCtrue) //MC loop
936 {
937
938 Int_t stackN = 0;
939
940 for(stackN = 0; stackN < stack->GetNtrack(); stackN++) //loop over stack
941 {
942
943 Bool_t isPrimary = stack->IsPhysicalPrimary(stackN);
944 Bool_t isSecondary = stack->IsSecondaryFromMaterial(stackN);
945
946 const TParticle *tparticle = stack->Particle(stackN);
947 Long_t pdgCode = tparticle->GetPdgCode();
948
949 Double_t pTGen = tparticle->Pt();
950 Double_t eta = tparticle->Eta();
951 //check which particle it is
952
953 //Alpha
954 if(pdgCode == 1000020040)
955 {
956 fHistHelium4PtGen->Fill(pTGen);
957 if(isPrimary) fHistHelium4PtGenPrim->Fill(pTGen);
958 if(isSecondary) fHistHelium4PtGenSec->Fill(pTGen);
959 if(eta < 0.8)fHistHelium4PtGenEta->Fill(pTGen);
960 if(isPrimary && eta < 0.8)fHistHelium4PtGenPrimEta->Fill(pTGen);
961 }
962
963 //Anti-Alpha
964 if(pdgCode == -1000020040)
965 {
966 fHistAntiHelium4PtGen->Fill(pTGen);
967 if(isPrimary) fHistAntiHelium4PtGenPrim->Fill(pTGen);
968 if(isSecondary) fHistAntiHelium4PtGenSec->Fill(pTGen);
969 if(eta < 0.8)fHistAntiHelium4PtGenEta->Fill(pTGen);
970 }
971
972
973 }//end loop over stack
974
975
976 }//end MC
977}