]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGLF/STRANGENESS/Hypernuclei/AliAnalysisTaskAntiHe4.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Updates from Debojit - PIDCorr
[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
20#include "AliAnalysisTaskAntiHe4.h"
21
22
23///////////////////////////////////////////////////////////
24// //
25// Analysis for the observation of anti-alpha particles. //
26// //
27///////////////////////////////////////////////////////////
28
29
30ClassImp(AliAnalysisTaskAntiHe4)
31
32//________________________________________________________________________
33AliAnalysisTaskAntiHe4::AliAnalysisTaskAntiHe4()
34 : AliAnalysisTaskSE(),
35 fEventHandler(0),
36 fESD(0),
37 fHistCentralityClass10(0),
38 fHistCentralityPercentile(0),
39 fHistTriggerStat(0),
40 fHistTriggerStatAfterEventSelection(0),
41 fHistDEDx(0),
42 fHistTOF3D(0),
43 fHistAlpha(0),
44 fHistAlphaSignal(0),
45 fGraphAlphaSignal(0),
46 fNCounter(0),
47 fHistDeDx(0),
48 fHistDeDxRegion(0),
49 fHistDeDxSharp(0),
50 fHistTOF2D(0),
51 fHistTOFnuclei(0x0),
52 fNTriggers(5),
53 fBBParametersLightParticles(),
54 fBBParametersNuclei(),
55 fMCtrue(0),
56 fTriggerFired(),
57 fESDtrackCuts(0),
58 fESDtrackCutsSharp(0),
59 fESDpid(0),
60 fAntiAlpha(0),
61 fTree(0),
62 fOutputContainer(0)
63{
64 // default Constructor
65
66
67 // Define input and output slots here
68}
69
70//________________________________________________________________________
71AliAnalysisTaskAntiHe4::AliAnalysisTaskAntiHe4(const char *name)
72 : AliAnalysisTaskSE(name),
73 fEventHandler(0),
74 fESD(0),
75 fHistCentralityClass10(0),
76 fHistCentralityPercentile(0),
77 fHistTriggerStat(0),
78 fHistTriggerStatAfterEventSelection(0),
79 fHistDEDx(0),
80 fHistTOF3D(0),
81 fHistAlpha(0),
82 fHistAlphaSignal(0),
83 fGraphAlphaSignal(0),
84 fNCounter(0),
85 fHistDeDx(0),
86 fHistDeDxRegion(0),
87 fHistDeDxSharp(0),
88 fHistTOF2D(0),
89 fHistTOFnuclei(0x0),
90 fNTriggers(5),
91 fBBParametersLightParticles(),
92 fBBParametersNuclei(),
93 fMCtrue(0),
94 fTriggerFired(),
95 fESDtrackCuts(0),
96 fESDtrackCutsSharp(0),
97 fESDpid(0),
98 fAntiAlpha(0),
99 fTree(0),
100 fOutputContainer(0)
101{
102 // Constructor
103
104 // Define input and output slots here
105 // Input slot #0 works with a TChain
106 DefineInput(0, TChain::Class());
107 // Output slot #0 writes into a TH1 container
108 DefineOutput(1, TObjArray::Class());
109 DefineOutput(2, TTree::Class());
110 //
111 // cuts for candidates
112 //
113 fESDtrackCuts = new AliESDtrackCuts("AliESDtrackCuts","AliESDtrackCuts");
114 //
115 fESDtrackCuts->SetAcceptKinkDaughters(kFALSE);
116 fESDtrackCuts->SetMinNClustersTPC(70);
117 fESDtrackCuts->SetMaxChi2PerClusterTPC(5);
118 fESDtrackCuts->SetMaxDCAToVertexXY(3);
119 fESDtrackCuts->SetMaxDCAToVertexZ(2);
120 fESDtrackCuts->SetRequireTPCRefit(kTRUE);
121 //fESDtrackCuts->SetRequireITSRefit(kTRUE);
122 fESDtrackCuts->SetMinNClustersITS(2);
123 fESDtrackCuts->SetEtaRange(-0.8,0.8);
124 //
125 // cuts for final plots
126 //
127 fESDtrackCutsSharp = new AliESDtrackCuts("AliESDtrackCuts","AliESDtrackCuts");
128 fESDtrackCutsSharp->SetAcceptKinkDaughters(kFALSE);
129 fESDtrackCutsSharp->SetMinNClustersTPC(80);
130 fESDtrackCutsSharp->SetMaxChi2PerClusterITS(10);// TO BE INVESTIGATED !!!!!!!!!!!!!!
131 fESDtrackCutsSharp->SetMaxChi2PerClusterTPC(5);
132 fESDtrackCutsSharp->SetRequireTPCRefit(kTRUE);
133 fESDtrackCutsSharp->SetRequireITSRefit(kTRUE);
134 fESDtrackCutsSharp->SetMinNClustersITS(2);
135 fESDtrackCutsSharp->SetMaxDCAToVertexXY(0.1);
136 fESDtrackCutsSharp->SetMaxDCAToVertexZ(0.5);
137 fESDtrackCutsSharp->SetEtaRange(-0.8,0.8);
138
139 //ESD Track cuts from TestFilterRawTask
140 /* fESDtrackCuts = new AliESDtrackCuts("AliESDtrackCuts","AliESDtrackCuts");
141 fESDtrackCuts->SetMinNClustersTPC(80);
142 fESDtrackCuts->SetMinNClustersITS(2);
143 fESDtrackCuts->SetAcceptKinkDaughters(kFALSE);
144 fESDtrackCuts->SetMaxChi2PerClusterTPC(5);
145 fESDtrackCuts->SetRequireTPCRefit(kTRUE);
146 fESDtrackCuts->SetRequireITSRefit(kTRUE);
147 fESDtrackCuts->SetEtaRange(-1,1);
148 fESDtrackCuts->SetMaxDCAToVertexXY(1);
149 fESDtrackCuts->SetMaxDCAToVertexZ(2);
150 //test strickter cuts
151 //fESDtrackCuts->SetMaxDCAToVertexXY(0.1);
152 //fESDtrackCuts->SetMaxDCAToVertexZ(0.5);
153 //fESDtrackCuts->SetEtaRange(-0.8,0.8);
154 */
155
156 fMCtrue = kTRUE;
157
158}
159
160//________________________________________________________________________
161void AliAnalysisTaskAntiHe4::UserCreateOutputObjects()
162{
163 // Create histograms
164 // Called once
165 //
166 //histogram to count number of events
167 //
168 fHistCentralityClass10 = new TH1F("fHistCentralityClass10", "centrality with class10", 11, -0.5, 10.5);
169 fHistCentralityClass10->GetXaxis()->SetTitle("Centrality");
170 fHistCentralityClass10->GetYaxis()->SetTitle("Entries");
171 //
172 fHistCentralityPercentile = new TH1F("fHistCentralityPercentile", "centrality with percentile", 101, -0.1, 100.1);
173 fHistCentralityPercentile->GetXaxis()->SetTitle("Centrality");
174 fHistCentralityPercentile->GetYaxis()->SetTitle("Entries");
175 //
176 //trigger statitics histogram
177 //
178 fHistTriggerStat = new TH1F("fHistTriggerStat","Trigger statistics", fNTriggers,-0.5,fNTriggers-0.5);
179 const Char_t* aTriggerNames[] = { "kMB", "kCentral", "kSemiCentral", "kEMCEJE", "kEMCEGA" };
180 for ( Int_t ii=0; ii < fNTriggers; ii++ )
181 fHistTriggerStat->GetXaxis()->SetBinLabel(ii+1, aTriggerNames[ii]);
182
183 fHistTriggerStatAfterEventSelection = new TH1F("fHistTriggerStatAfterEventSelection","Trigger statistics after event selection", fNTriggers,-0.5,fNTriggers-0.5);
184 for ( Int_t ii=0; ii < fNTriggers; ii++ )
185 fHistTriggerStatAfterEventSelection->GetXaxis()->SetBinLabel(ii+1, aTriggerNames[ii]);
186
187 //dE/dx performance
188 fHistDEDx= new TH3F("fHistDEDx","DEDx",1000,0.01,100,1000,1,2000,2,-2,2);
189 BinLogAxis(fHistDEDx, 0);
190 BinLogAxis(fHistDEDx, 1);
191 fHistDEDx->GetXaxis()->SetTitle("p_{tot}/sign");
192 fHistDEDx->GetYaxis()->SetTitle("TPC signal");
193
194 fHistDeDx = new TH2F("fHistDeDx", "dE/dx", 1000, 0.1, 6.0, 1000, 0.0, 1000);
195 fHistDeDx->GetYaxis()->SetTitle("TPC dE/dx signal (a.u.)");
196 fHistDeDx->GetXaxis()->SetTitle("#frac{p}{z} (GeV/c)");
197 BinLogAxis(fHistDeDx);
198
199 fHistDeDxRegion = new TH3F("fHistDeDxRegion", "dE/dx", 400, 0., 6.0, 300, 0., 3., 4, -0.5, 3.5);
200 fHistDeDxRegion->GetYaxis()->SetTitle("TPC dE/dx signal (a.u.)");
201 fHistDeDxRegion->GetXaxis()->SetTitle("#frac{p}{z} (GeV/c)");
202
203 fHistDeDxSharp = new TH2F("fHistDeDxSharp", "dE/dx", 1000, 0.1, 6.0, 1000, 0.0, 1000);
204 fHistDeDxSharp->GetYaxis()->SetTitle("TPC dE/dx signal (a.u.)");
205 fHistDeDxSharp->GetXaxis()->SetTitle("#frac{p}{z} (GeV/c)");
206 BinLogAxis(fHistDeDxSharp);
207
208 //TOF performance
209 fHistTOF3D = new TH3F("fHistTOF3D","mass determination from TOF",400,0.25,4.5,2,-0.5,1.5,2,-1,1);
210
211 fHistTOF2D = new TH2F("fHistTOF2D", "TOF2D", 500, 0.0, 10., 2250, 0.2, 1.1);
212 fHistTOF2D->GetYaxis()->SetTitle("#beta");
213 fHistTOF2D->GetXaxis()->SetTitle("p (GeV/c)");
214
215 fHistTOFnuclei = new TH2F("fHistTOFnuclei", "TOF", 500, 0.0, 10., 2250, 0.7, 1.);
216 fHistTOFnuclei->GetYaxis()->SetTitle("#beta");
217 fHistTOFnuclei->GetXaxis()->SetTitle("#frac{p}{z} (GeV/c)");
218
219 //alphas
220 fHistAlpha = new TH1F("fHistAlpha", "Anti-Alpha", 22, 1.12, 4.31);
221 fHistAlpha->GetYaxis()->SetTitle("Counts");
222 fHistAlpha->GetXaxis()->SetTitle("#frac{m^{2}}{z^{2}} (GeV^{2}/c^{4})");
223
224 fHistAlphaSignal = new TH1F("fHistAlphaSignal", "Anti-Alpha", 22, 1.12, 4.31);
225 fHistAlphaSignal->GetYaxis()->SetTitle("Counts");
226 fHistAlphaSignal->GetXaxis()->SetTitle("#frac{m^{2}}{z^{2}} (GeV^{2}/c^{4})");
227
228 fGraphAlphaSignal = new TGraph(20);
229 fNCounter = 0;
230 //
231 // (0.) dca, (1.) sign, (2.) particle Type, (3.) p_tot
232 //
233 TString axisNameMult[4]={"dca","sign","particleType","ptot"};
234 TString axisTitleMult[4]={"dca","sign","particleType","ptot"};
235 const Int_t kDcaBins = 76;
236 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};
237 //
238 // (0.) dca, (1.) sign, (2.) particle Type, (3.) p_tot
239 Int_t binsAntiAlpha[4] = {77, 2, 4, 100};
240 Double_t xminAntiAlpha[4] = {-3, -2, -0.5, 0};
241 Double_t xmaxAntiAlpha[4] = { 3, 2, 3.5, 6};
242 //
243 fAntiAlpha = new THnF("fAntiAlpha", "AntiAlpha; (0.) dca, (1.) sign, (2.) particle Type, (3.) p_tot", 4, binsAntiAlpha, xminAntiAlpha, xmaxAntiAlpha);
244 //
245 fAntiAlpha->GetAxis(0)->Set(kDcaBins, binsDca);
246 for (Int_t iaxis=0; iaxis<4;iaxis++){
247 fAntiAlpha->GetAxis(iaxis)->SetName(axisNameMult[iaxis]);
248 fAntiAlpha->GetAxis(iaxis)->SetTitle(axisTitleMult[iaxis]);
249 }
250 //
251 //
252 //
253 //
254 fOutputContainer = new TObjArray(1);
255 fOutputContainer->SetOwner(kTRUE);
256 fOutputContainer->SetName(GetName());
257 fOutputContainer->Add(fHistCentralityClass10);
258 fOutputContainer->Add(fHistCentralityPercentile);
259 fOutputContainer->Add(fHistTriggerStat);
260 fOutputContainer->Add(fHistTriggerStatAfterEventSelection);
261 fOutputContainer->Add(fHistDEDx);
262 fOutputContainer->Add(fHistDeDx);
263 fOutputContainer->Add(fHistDeDxRegion);
264 fOutputContainer->Add(fHistDeDxSharp);
265 fOutputContainer->Add(fAntiAlpha);
266 fOutputContainer->Add(fHistAlpha);
267 fOutputContainer->Add(fHistAlphaSignal);
268 fOutputContainer->Add(fGraphAlphaSignal);
269 fOutputContainer->Add(fHistTOF3D);
270 fOutputContainer->Add(fHistTOF2D);
271 fOutputContainer->Add(fHistTOFnuclei);
272
273
274 //------------ Tree and branch definitions ----------------//
275 fTree = new TTree("tree", " alpha tree");
276
277 //------------ Event variables ------------//
278 fTree->Branch("Name",Name,"Name/C");
279 fTree->Branch("Evnt",&evnt, "evnt/I");
280 fTree->Branch("itrk", &itrk, "itrk/I");
281
282 //-------------------------------------------//
283 //----------- Track variables --------------//
284
285 fTree->Branch("TrkPtot",TrkPtot,"TrkPtot[itrk]/F");
286 fTree->Branch("TPCPtot",TPCPtot,"TPCPtot[itrk]/F");
287 fTree->Branch("DeDx",DeDx,"DeDx[itrk]/F");
288 fTree->Branch("Sign",Sign,"Sign[itrk]/F");
289 fTree->Branch("DCAXY",DCAXY,"DCAXY[itrk]/F");
290 fTree->Branch("DCAZ",DCAZ,"DCAZ[itrk]/F");
291 fTree->Branch("ITSnCluster",ITSnCluster,"ITSnCluster[itrk]/F");
292 fTree->Branch("TPCNsignal",TPCNsignal,"TPCNsignal[itrk]/F");
293 fTree->Branch("Mass",Mass,"Mass[itrk]/F");
294 //
295 fTree->Branch("ITSRefit",ITSRefit,"ITSRefit[itrk]/F");
296 fTree->Branch("TOFtime",TOFtime,"TOFtime[itrk]/F");
297 fTree->Branch("TOFRefit",TOFRefit,"TOFRefit[itrk]/F");
298 fTree->Branch("TOFout",TOFout,"TOFout[itrk]/F");
299 //
300 fTree->Branch("ITSsignal",ITSsignal,"ITSsignal[itrk]/F");
301 fTree->Branch("SharedClusters",SharedClusters,"SharedClusters[itrk]/F");
302 fTree->Branch("fFileName",fFileName,"fFileName/C");
303 fTree->Branch("fEventNumber",fEventNumber,"fEventNumber/I");
304
305}
306
307//________________________________________________________________________
308void AliAnalysisTaskAntiHe4::UserExec(Option_t *)
309{
310 // Main loop
311 // Called for each event
312
313 //get Event-Handler for the trigger information
314 fEventHandler= dynamic_cast<AliInputEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
315 if (!fEventHandler) {
316 AliError("Could not get InputHandler");
317 //return -1;
318 return;
319 }
320
321 fESD=dynamic_cast<AliESDEvent*>(InputEvent());
322 if (!fESD) {
323 //Printf("ERROR: fESD not available");
324 return;
325 }
326
327 if (SetupEvent() < 0) {
328 PostData(1, fOutputContainer);
329 return;
330 }
331
332 const AliESDVertex *vertex = fESD->GetPrimaryVertexTracks();
333 if(vertex->GetNContributors()<1) {
334 // SPD vertex
335 vertex = fESD->GetPrimaryVertexSPD();
336 if(vertex->GetNContributors()<1) {
337 PostData(1, fOutputContainer);
338 return;
339 }
340
341 }
342 // check if event is selected by physics selection class
343 //
344 Bool_t isSelected = kFALSE;
345 isSelected = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
346 if (!isSelected || TMath::Abs(vertex->GetZv()) > 10) {
347 PostData(1, fOutputContainer);
348 return;
349 }
350
351 //
352 //centrality selection in PbPb
353 //
354 Int_t centralityClass10 = -1;
355 Int_t centralityPercentile = -1;
356 //
357 if (fESD->GetEventSpecie() == 4) { //species == 4 == PbPb
358 // PbPb
359 AliCentrality *esdCentrality = fESD->GetCentrality();
360 centralityClass10 = esdCentrality->GetCentralityClass10("V0M"); // centrality percentile determined with V0
361 centralityPercentile = esdCentrality->GetCentralityPercentile("V0M"); // centrality percentile determined with V0
362 if(!fMCtrue){
363 if (centralityPercentile < 0. || centralityPercentile > 80. ) return; //select only events with centralities between 0 and 80 %
364 }
365 }
366 //
367 fHistCentralityClass10->Fill(centralityClass10);
368 fHistCentralityPercentile->Fill(centralityPercentile);
369 //
370 if(fTriggerFired[0] == kTRUE) fHistTriggerStatAfterEventSelection->Fill(0);
371 if(fTriggerFired[1] == kTRUE) fHistTriggerStatAfterEventSelection->Fill(1);
372 if(fTriggerFired[2] == kTRUE) fHistTriggerStatAfterEventSelection->Fill(2);
373 if(fTriggerFired[3] == kTRUE) fHistTriggerStatAfterEventSelection->Fill(3);
374 if(fTriggerFired[4] == kTRUE) fHistTriggerStatAfterEventSelection->Fill(4);
375 //
376 //
377 if (!fESDpid) fESDpid = ((AliESDInputHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->GetESDpid();
378 if (!fESDpid) {
379 fESDpid = new AliESDpid(); // HACK FOR MC PBPB --> PLEASE REMOVE AS SOON AS POSSIBLE
380 fESDpid->GetTPCResponse().SetBetheBlochParameters(1.28778e+00/50., 3.13539e+01, TMath::Exp(-3.16327e+01), 1.87901e+00, 6.41583e+00);
381 }
382 //
383 Float_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z for th // for Anti-Alpha
384 evnt = fESD->GetEventNumberInFile();
385 sscanf(fInputHandler->GetTree()->GetCurrentFile()->GetName(),"%s", Name);
386 itrk = 0;
387 //
388 Int_t runNumber = 0;
389 runNumber = fESD->GetRunNumber();
390 //
391 Bool_t fillTree = kFALSE;
392 // Track loop to fill the spectram
393 for (Int_t iTracks = 0; iTracks < fESD->GetNumberOfTracks(); iTracks++) {
394 //
395 AliESDtrack* track = dynamic_cast<AliESDtrack*>(fESD->GetTrack(iTracks));
396 if (!fESDtrackCuts->AcceptTrack(track)) continue;
397 //
398 Double_t nClustersTPCPID=track->GetTPCsignalN();
399 if(nClustersTPCPID < 60) continue;
400 //
401 if (!track->GetInnerParam()) continue;
402 //
403 Double32_t signal[4] = {0,0,0,0};
404 Char_t ncl[3];
405 Char_t nrows[3];
406 if (track->GetTPCdEdxInfo()) {
407 track->GetTPCdEdxInfo()->GetTPCSignalRegionInfo(signal,ncl,nrows);
408 }
409 //
410 Double_t ptot = track->GetInnerParam()->GetP();
411 Double_t ptotInc = track->GetP(); // total momentum of the incoming particle
412 Double_t sign = track->GetSign();
413 //
414 Double_t eta = TMath::Abs(track->Eta());
415 TBits shared = track->GetTPCSharedMap();
416
417 track->GetImpactParameters(dca, cov);
418 Float_t dcaXY = TMath::Sqrt(dca[0]*dca[0]);
419 Float_t dcaXYsign = dca[0];
420 Float_t dcaZ = TMath::Sqrt(dca[1]*dca[1]);
421 //
422 //
423 Double_t tpcSignal = track->GetTPCsignal();
424 //
425 //PID via specific energy loss in the TPC
426 //define the arrays for the Bethe-Bloch-Parameters
427 //
428
429 SetBBParameters(runNumber);
430
431 //define expected signals for the various species
432 Double_t expSignalProton = AliExternalTrackParam::BetheBlochAleph(ptot/0.93827,fBBParametersLightParticles[0],fBBParametersLightParticles[1],fBBParametersLightParticles[2],fBBParametersLightParticles[3],fBBParametersLightParticles[4]);
433 Double_t expSignalDeuteron = AliExternalTrackParam::BetheBlochAleph(ptot/1.875612,fBBParametersLightParticles[0],fBBParametersLightParticles[1],fBBParametersLightParticles[2],fBBParametersLightParticles[3],fBBParametersLightParticles[4]);
434 Double_t expSignalTriton = AliExternalTrackParam::BetheBlochAleph(ptot/2.808921,fBBParametersLightParticles[0],fBBParametersLightParticles[1],fBBParametersLightParticles[2],fBBParametersLightParticles[3],fBBParametersLightParticles[4]);
435
436 Double_t expSignalHelium3 = 4*AliExternalTrackParam::BetheBlochAleph(2*ptot/2.80941,fBBParametersNuclei[0],fBBParametersNuclei[1],fBBParametersNuclei[2],fBBParametersNuclei[3],fBBParametersNuclei[4]);
437 Double_t expSignalHelium4 = 4*AliExternalTrackParam::BetheBlochAleph(2*ptot/3.728401,fBBParametersNuclei[0],fBBParametersNuclei[1],fBBParametersNuclei[2],fBBParametersNuclei[3],fBBParametersNuclei[4]);
438
439 //
440 Float_t mass = 0;
441 Float_t time = -1;
442 Float_t beta = 0;
443 Float_t gamma = -1;
444 Bool_t hasTOFout = kFALSE;
445 Bool_t hasTOFtime = kFALSE;
446 Float_t length = track->GetIntegratedLength();
447 UInt_t status = track->GetStatus();
448
449 if (length > 350) {
450 time = track->GetTOFsignal() - fESDpid->GetTOFResponse().GetStartTime(track->P());
451 if (time > 0) {
452 beta = length / (2.99792457999999984e-02 * time);
453 gamma = 1/TMath::Sqrt(1 - beta*beta);
454 mass = ptot/TMath::Sqrt(gamma*gamma - 1); // using inner TPC mom. as approx.
455 }
456 }
457 //
458 // fill tree and print candidates (for short list)
459 //
460 Float_t cut = 4*AliExternalTrackParam::BetheBlochAleph(2*ptot/(0.938*3),1.1931,
461 31.9806,
462 5.04114e-11,
463 2.13096,
464 2.38541);
465 if (eta < 0.8 && tpcSignal > 120 && tpcSignal > cut && tpcSignal < 1000 && track->GetTPCsignalN() > 60 && dcaZ < 15 && dcaXY < 15 && ptot > 1.0 && ptot < 20) {
466 //
467 cout << "AntiAlphaEvent" << " "
468 << AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()->GetTree()->GetCurrentFile()->GetName() << " "
469 << "event number in file: " << fESD->GetEventNumberInFile()
470 << " Index " << iTracks
471 << " ptot: " << ptot
472 << " sig: " << tpcSignal <<endl;
473 //
474 fillTree = kTRUE;
475 //
476 TrkPtot[itrk] = track->P();
477 TPCPtot[itrk] = ptot;
478 DeDx[itrk] = tpcSignal;
479 DCAZ[itrk] = dcaZ;
480 TPCNsignal[itrk] = track->GetTPCsignalN();
481 ITSnCluster[itrk] = track->GetNcls(0);
482 Sign[itrk] = sign;
483 DCAXY[itrk] = dcaXY;
484 Mass[itrk] = mass;
485 //
486 ITSsignal[itrk] = track->GetITSsignal();
487 SharedClusters[itrk] = shared.CountBits();
488 sscanf(fInputHandler->GetTree()->GetCurrentFile()->GetName(),"%s", fFileName);
489 //fFileName[itrk] = fInputHandler->GetTree()->GetCurrentFile()->GetName();
490 fEventNumber[itrk] = fESD->GetEventNumberInFile();
491 //
492 if(status&AliESDtrack::kITSrefit)
493 ITSRefit[itrk] = 1;
494 else ITSRefit[itrk] = 0;
495 //
496 if (time < 99998) {
497 TOFRefit[itrk] = 1;
498 } else {
499 TOFRefit[itrk] = 0;
500 }
501 //
502 hasTOFout = status&AliESDtrack::kTOFout;
503 hasTOFtime = status&AliESDtrack::kTIME;
504 //
505 TOFtime[itrk] = hasTOFtime;
506 TOFout[itrk] = hasTOFout;
507 itrk++;
508 }
509 //
510 // do pid fill histogram for raw ratios
511 //
512 // (0.) dca, (1.) sign, (2.) particle Type, (3.) p_tot
513 Int_t id = -1;
514 if (ptot > 0.2 && TMath::Abs(tpcSignal - expSignalProton)/expSignalProton < 0.2) id = 0;
515 if (ptot > 0.3 && TMath::Abs(tpcSignal - expSignalDeuteron)/expSignalDeuteron < 0.2) id = 1;
516 if (ptot > 0.7 && TMath::Abs(tpcSignal - expSignalTriton)/expSignalTriton < 0.2) id = 2;
517 if (ptot > 0.5 && (tpcSignal - expSignalHelium3)/expSignalHelium3 > -0.1 && (tpcSignal - expSignalHelium3)/expSignalHelium3 < 0.2) id = 3;
518 //
519 Double_t vecAntiAlpha[4] = {dcaXYsign, sign, id, ptotInc};
520 if (id != -1 && tpcSignal > 120) fAntiAlpha->Fill(vecAntiAlpha);
521 //
522 // fill final histograms
523 //
524 if (!fESDtrackCutsSharp->AcceptTrack(track) || shared.CountBits() > 1 || track->GetTPCsignalN() < 80 || track->GetKinkIndex(0) != 0 || track->GetTPCNclsIter1() < 80) continue;
525 //
526 fHistDEDx->Fill(ptot,track->GetTPCsignal(), sign);
527 if (TMath::Abs(tpcSignal - expSignalHelium4)/expSignalHelium4 < 0.2 && time < 99998) fHistTOF3D->Fill(mass,1,sign);
528 //
529 // remove overlapping tracks with special dE/dx cut
530 //
531 //if (signal[1] < 6) continue;
532 //if (signal[0]/signal[1] > 1.6 || signal[2]/signal[1] > 1.6 || signal[3]/signal[1] > 1.3) continue;
533 //
534 if(sign<0) {
535 fHistDeDx->Fill(ptot, track->GetTPCsignal());
536 if (track->GetTPCsignalN() > 100 &&
537 TMath::Abs(track->Eta()) < 0.8 &&
538 signal[3]/signal[1] > 0.6 &&
539 signal[0]/signal[1] > 0.5 &&
540 signal[3]/signal[1] < 1.2 &&
541 track->GetTPCNclsIter1() > 70 &&
542 track->GetTPCnclsS() < 10) {
543 fHistDeDxSharp->Fill(ptot, track->GetTPCsignal());
544 }
545 }
546 //
547 // dE/dx for specific regions
548 //
549 for(Int_t iSig = 0; iSig < 4; iSig++) {
550 if (signal[1] > 6) fHistDeDxRegion->Fill(ptot,signal[iSig]/signal[1],iSig);
551 }
552 //
553 // alpha TOF plot
554 //
555 if((track->GetTPCsignal()-expSignalHelium4)/expSignalHelium4 > -0.15 && (track->GetTPCsignal()-expSignalHelium4)/expSignalHelium4 < 0.15) {
556 //TOF
557 hasTOFout = status&AliESDtrack::kTOFout;
558 hasTOFtime = status&AliESDtrack::kTIME;
559 Bool_t hasTOF = kFALSE;
560
561 if (hasTOFout && hasTOFtime) hasTOF = kTRUE;
562 //
563 if (length < 350.) hasTOF = kFALSE;
564
565 Float_t time0 = fESDpid->GetTOFResponse().GetStartTime(track->P());//fESDpid->GetTOFResponse().GetTimeZero();
566 //
567 if (length > 350. && hasTOF == kTRUE && ptot < 4) {
568 time = track->GetTOFsignal() - time0;
569 if (time > 0) {
570 beta = length / (2.99792457999999984e-02 * time);
571 if (beta < 0.975) {
572 gamma = 1/TMath::Sqrt(1 - beta*beta);
573 mass = ptot/TMath::Sqrt(gamma*gamma - 1); // using inner TPC mom. as approx.
574 if (TMath::Sqrt(track->GetTOFsignalDz()*track->GetTOFsignalDz() + track->GetTOFsignalDx()*track->GetTOFsignalDx()) < 5.) {
575 fHistAlpha->Fill(mass*mass);
576 if (mass*mass > 3. && mass*mass < 4.) {
577 fHistAlphaSignal->Fill(mass*mass);
578 fGraphAlphaSignal->SetPoint(fNCounter, ptot, track->GetTPCsignal());
579 fNCounter++;
580 }
581 fHistTOFnuclei->Fill(ptot,beta);
582 }
583 }
584 }
585 }
586 }
587
588 }//end loop over tracks
589 if (fillTree) fTree->Fill();
590
591 // Post output data.
592 PostData(1, fOutputContainer);
593 PostData(2, fTree);
594}
595
596//________________________________________________________________________
597void AliAnalysisTaskAntiHe4::Terminate(const Option_t *)
598{
599 // Draw result to the screen
600 // Called once at the end of the query
601
602 //get output data and darw 'fHistPt'
603 if (!GetOutputData(0)) return;
604
605}
606
607
608//________________________________________________________________________
609void AliAnalysisTaskAntiHe4::BinLogAxis(const THn *h, Int_t axisNumber) {
610 //
611 // Method for the correct logarithmic binning of histograms
612 //
613 TAxis *axis = h->GetAxis(axisNumber);
614 int bins = axis->GetNbins();
615
616 Double_t from = axis->GetXmin();
617 Double_t to = axis->GetXmax();
618 Double_t *newBins = new Double_t[bins + 1];
619
620 newBins[0] = from;
621 Double_t factor = pow(to/from, 1./bins);
622
623 for (int i = 1; i <= bins; i++) {
624 newBins[i] = factor * newBins[i-1];
625 }
626 axis->Set(bins, newBins);
627 delete [] newBins;
628
629}
630
631//________________________________________________________________________
632void AliAnalysisTaskAntiHe4::BinLogAxis(const TH3 *h, Int_t axisNumber) {
633 //
634 // Method for the correct logarithmic binning of histograms
635 //
636 TAxis *axis = h->GetXaxis();
637 if (axisNumber == 1) axis = h->GetYaxis();
638 if (axisNumber == 2) axis = h->GetZaxis();
639 int bins = axis->GetNbins();
640
641 Double_t from = axis->GetXmin();
642 Double_t to = axis->GetXmax();
643 Double_t *newBins = new Double_t[bins + 1];
644
645 newBins[0] = from;
646 Double_t factor = pow(to/from, 1./bins);
647
648 for (int i = 1; i <= bins; i++) {
649 newBins[i] = factor * newBins[i-1];
650 }
651 axis->Set(bins, newBins);
652 delete [] newBins;
653
654}
655//________________________________________________________________________
656void AliAnalysisTaskAntiHe4::BinLogAxis(const TH1 *h) {
657 //
658 // Method for the correct logarithmic binning of histograms
659 //
660 TAxis *axis = h->GetXaxis();
661 int bins = axis->GetNbins();
662
663 Double_t from = axis->GetXmin();
664 Double_t to = axis->GetXmax();
665 Double_t *newBins = new Double_t[bins + 1];
666
667 newBins[0] = from;
668 Double_t factor = pow(to/from, 1./bins);
669
670 for (int i = 1; i <= bins; i++) {
671 newBins[i] = factor * newBins[i-1];
672 }
673 axis->Set(bins, newBins);
674 delete [] newBins;
675
676}
677//_____________________________________________________________________________
678Int_t AliAnalysisTaskAntiHe4::Initialize() {
679
680
681 // -- Reset Event
682 // ----------------
683 ResetEvent();
684
685 return 0;
686}
687//________________________________________________________________________
688Int_t AliAnalysisTaskAntiHe4::SetupEvent() {
689 // Setup Reading of event
690
691 ResetEvent();
692 // -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
693 // -- Get Trigger
694 // -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
695
696 //Bool_t isTriggered = IsTriggered();
697 IsTriggered();
698 return 0;
699}
700//_____________________________________________________________________________
701void AliAnalysisTaskAntiHe4::ResetEvent() {
702 // Reset event
703 // -- Reset QA
704 for (Int_t ii = 0; ii < fNTriggers; ++ii)
705 fTriggerFired[ii] = kFALSE;
706
707 return;
708}
709//________________________________________________________________________
710Bool_t AliAnalysisTaskAntiHe4::IsTriggered() {
711 // Check if Event is triggered and fill Trigger Histogram
712
713 if ((fEventHandler->IsEventSelected() & AliVEvent::kMB)) fTriggerFired[0] = kTRUE;
714 if ((fEventHandler->IsEventSelected() & AliVEvent::kCentral)) fTriggerFired[1] = kTRUE;
715 if ((fEventHandler->IsEventSelected() & AliVEvent::kSemiCentral)) fTriggerFired[2] = kTRUE;
716 if ((fEventHandler->IsEventSelected() & AliVEvent::kEMCEJE)) fTriggerFired[3] = kTRUE;
717 if ((fEventHandler->IsEventSelected() & AliVEvent::kEMCEGA)) fTriggerFired[4] = kTRUE;
718
719 Bool_t isTriggered = kFALSE;
720
721 for (Int_t ii=0; ii<fNTriggers; ++ii) {
722 if(fTriggerFired[ii]) {
723 isTriggered = kTRUE;
724 fHistTriggerStat->Fill(ii);
725 }
726 }
727
728 return isTriggered;
729}
730//________________________________________________________________________
731void AliAnalysisTaskAntiHe4::SetBBParameters(Int_t runNumber){
732
733 if(runNumber < 166500){ //LHC10h
734
735 fBBParametersLightParticles[0] = 1.45802;
736 fBBParametersLightParticles[1] = 27.4992;
737 fBBParametersLightParticles[2] = 4.00313e-15;
738 fBBParametersLightParticles[3] = 2.48485;
739 fBBParametersLightParticles[4] = 8.31768;
740
741 fBBParametersNuclei[0] = 1.69155;
742 fBBParametersNuclei[1] = 27.4992;
743 fBBParametersNuclei[2] = 4.00313e-15;
744 fBBParametersNuclei[3] = 2.48485;
745 fBBParametersNuclei[4] = 8.31768;
746
747 }
748
749 if(runNumber > 166500){ //LHC11h
750
751 fBBParametersLightParticles[0] = 4.69637;//1.11243;
752 fBBParametersLightParticles[1] = 7.51827;//26.1144;
753 fBBParametersLightParticles[2] = 0.0183746;//4.00313e-15;
754 fBBParametersLightParticles[3] = 2.60;//2.72969;
755 fBBParametersLightParticles[4] = 2.7;//9.15038;
756
757 fBBParametersNuclei[0] = 1.4906;
758 fBBParametersNuclei[1] = 27.9758;
759 fBBParametersNuclei[2] = 4.00313e-15;
760 fBBParametersNuclei[3] = 2.50804;
761 fBBParametersNuclei[4] = 8.31768;
762
763 }
764}