]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG2/SPECTRA/AliAnalysisTaskSEITSsaSpectra.cxx
Size of fTrackLabel array increased to 4
[u/mrichter/AliRoot.git] / PWG2 / SPECTRA / AliAnalysisTaskSEITSsaSpectra.cxx
CommitLineData
36be14b3 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// AliAnalysisTaskSE for the extraction of the various histograms to
18// study the pt spectra of identified hadrons:
19// - log(dEdx)-log(dEdxBB) distributions for pions, kaons and protons in pt bins
20// - Pt distributions of pions, kaons and protons with nSigma PID
21// Authors:
22// E. Biolcati, biolcati@to.infn.it
23// L. Milano, milano@to.infn.it
24// F. Prino, prino@to.infn.it
25///////////////////////////////////////////////////////////////////////////
26
36be14b3 27#include <TH1F.h>
28#include <TRandom3.h>
29#include <TH2F.h>
30#include <TChain.h>
31#include <TNtuple.h>
32#include <TParticle.h>
33#include <Rtypes.h>
34#include "AliAnalysisTaskSE.h"
35#include "AliAnalysisManager.h"
36#include "AliAnalysisDataContainer.h"
37#include "AliESDEvent.h"
38#include "AliESDInputHandler.h"
39#include "AliESDtrack.h"
40#include "AliStack.h"
41#include "AliMCEventHandler.h"
42#include "AliMCEvent.h"
43#include "AliPhysicsSelection.h"
44#include "AliAnalysisTaskSEITSsaSpectra.h"
c841b89f 45#include "AliESDtrackCuts.h"
1f75ac7d 46#include "AliCentrality.h"
47#include "AliMultiplicity.h"
48#include "AliESDUtils.h"
36be14b3 49
50ClassImp(AliAnalysisTaskSEITSsaSpectra)
6b77d2c0 51/* $Id$ */
36be14b3 52
53//________________________________________________________________________
54AliAnalysisTaskSEITSsaSpectra::AliAnalysisTaskSEITSsaSpectra():
1f75ac7d 55AliAnalysisTaskSE("Task CFit"),
36be14b3 56 fESD(0),
7f8f788e 57 fesdTrackCutsMult(0),
36be14b3 58 fOutput(0),
59 fHistNEvents(0),
c841b89f 60 fHistMult(0),
1f75ac7d 61 fHistCen(0),
36be14b3 62 fHistNTracks(0),
c841b89f 63 fHistNTracksPos(0),
64 fHistNTracksNeg(0),
36be14b3 65 fHistDEDX(0),
66 fHistDEDXdouble(0),
67 fHistBeforeEvSel(0),
68 fHistAfterEvSel(0),
69 fMinSPDPts(1),
70 fMinNdEdxSamples(3),
71 fMindEdx(0.),
c841b89f 72 fMinNSigma(1.5),
36be14b3 73 fMaxY(0.5),
c841b89f 74 fMaxChi2Clu(2.5),
36be14b3 75 fNSigmaDCAxy(7.),
76 fNSigmaDCAz(7.),
dac69cb0 77 fEtaRange(0.8),
1f75ac7d 78 fLowMult(-1),
79 fUpMult(-1),
80 fLowCentrality(-1.0),
81 fUpCentrality(-1.0),
82 fSPD(0),
83 fHImode(0),
6b77d2c0 84 fYear(2010),
36be14b3 85 fMC(kFALSE),
86 fSmearMC(kFALSE),
87 fSmearP(0.),
88 fSmeardEdx(0.),
89 fRandGener(0),
90 fFillNtuple(kFALSE),
91 fNtupleNSigma(0),
92 fNtupleMC(0)
93{
94 // Constructor
95 Double_t xbins[kNbins+1]={0.08,0.10,0.12,0.14,0.16,0.18,0.20,0.25,0.30,0.35,0.40,0.45,0.50,0.55,0.60,0.65,0.70,0.75,0.80,0.85,0.90,0.95,1.0};
96 for(Int_t iBin=0; iBin<kNbins+1; iBin++) fPtBinLimits[iBin]=xbins[iBin];
97 fRandGener=new TRandom3(0);
7f8f788e 98 fesdTrackCutsMult = new AliESDtrackCuts;
99 // TPC
100 fesdTrackCutsMult->SetMinNClustersTPC(70);
101 fesdTrackCutsMult->SetMaxChi2PerClusterTPC(4);
102 fesdTrackCutsMult->SetAcceptKinkDaughters(kFALSE);
103 fesdTrackCutsMult->SetRequireTPCRefit(kTRUE);
104 // ITS
105 fesdTrackCutsMult->SetRequireITSRefit(kTRUE);
106 fesdTrackCutsMult->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
107 AliESDtrackCuts::kAny);
108 fesdTrackCutsMult->SetDCAToVertex2D(kFALSE);
109 fesdTrackCutsMult->SetRequireSigmaToVertex(kFALSE);
110 fesdTrackCutsMult->SetEtaRange(-0.8,+0.8);
111 fesdTrackCutsMult->SetPtRange(0.15, 1e10);
112 SetYear(2010);
36be14b3 113
114 DefineInput(0, TChain::Class());
115 DefineOutput(1, TList::Class());
116 Printf("end of AliAnalysisTaskSEITSsaSpectra");
117}
118
119//___________________________________________________________________________
120AliAnalysisTaskSEITSsaSpectra::~AliAnalysisTaskSEITSsaSpectra(){
121 // Destructor
122 if (fOutput) {
123 delete fOutput;
124 fOutput = 0;
125 }
126 if(fRandGener) delete fRandGener;
127}
128
129//________________________________________________________________________
4d668cf7 130Double_t AliAnalysisTaskSEITSsaSpectra::CookdEdx(Double_t *s) const {
36be14b3 131 // truncated mean for the dEdx
132 Int_t nc=0;
133 Double_t dedx[4]={0.,0.,0.,0.};
134 for (Int_t il=0; il<4; il++) { // count good (>0) dE/dx values
135 if(s[il]>fMindEdx){
136 dedx[nc]= s[il];
137 nc++;
138 }
139 }
140 if(nc<fMinNdEdxSamples) return -1.;
141
142 Double_t tmp;
143 Int_t swap; // sort in ascending order
144 do {
145 swap=0;
146 for (Int_t i=0; i<nc-1; i++) {
147 if (dedx[i]<=dedx[i+1]) continue;
148 tmp=dedx[i];
149 dedx[i]=dedx[i+1];
150 dedx[i+1]=tmp;
151 swap++;
152 }
153 } while (swap);
154
155 Double_t sumamp=0,sumweight=0;
156 Double_t weight[4]={1.,1.,0.,0.};
157 if(nc==3) weight[1]=0.5;
158 else if(nc<3) weight[1]=0.;
159 for (Int_t i=0; i<nc; i++) {
160 sumamp+= dedx[i]*weight[i];
161 sumweight+=weight[i];
162 }
163 return sumamp/sumweight;
164}
165
166
7f8f788e 167//________________________________________________________________________
168void AliAnalysisTaskSEITSsaSpectra::SetYear(Int_t year){
169 // Set year dependent quantities
170 fYear=year;
171 if(fYear==2009){
172 fesdTrackCutsMult->SetMaxDCAToVertexXYPtDep("0.0350+0.0420/TMath::Power(pt,0.9)"); //2009 standard cut
173 fesdTrackCutsMult->SetMaxDCAToVertexZ(20); //2009 standard cut
174 }else{
175 fesdTrackCutsMult->SetMaxDCAToVertexXYPtDep("0.0182+0.0350/pt^1.01"); //2010 standard cut
176 fesdTrackCutsMult->SetMaxDCAToVertexZ(2); //2010 standard cut
177 }
178}
179
36be14b3 180//________________________________________________________________________
4d668cf7 181Bool_t AliAnalysisTaskSEITSsaSpectra::DCAcut(Double_t impactXY, Double_t impactZ, Double_t pt, Bool_t optMC) const {
36be14b3 182 // cut on transverse impact parameter updaated on 20-5-2010
183 // from the study of L. Milano, F. Prino on the ITS standalone tracks
184 // using the common binning of the TPC tracks
185 Double_t xyP[3];
186 Double_t zP[3];
187 if(optMC){
6b77d2c0 188 if(fYear==2009){
189 xyP[0]=88.63; //MC LHC10a12
190 xyP[1]=19.57;
191 xyP[2]=1.65;
192 zP[0]=140.98;
193 zP[1]=62.33;
194 zP[2]=1.15;
195 }else{
196 xyP[0]=36.; //MC LHC10d1
197 xyP[1]=43.9;
198 xyP[2]=1.3;
199 zP[0]=111.9;
200 zP[1]=59.8;
201 zP[2]=1.2;
202 }
36be14b3 203 }
204 else{
6b77d2c0 205 if(fYear==2009){
206 xyP[0]=85.28;//DATA 900 GeV pass6
207 xyP[1]=25.78;
208 xyP[2]=1.55;
209 zP[0]=146.80;
210 zP[1]=70.07;
211 zP[2]=1.11;
212 }else{
213 xyP[0]=32.7;//DATA 7 TeV pass2
214 xyP[1]=44.8;
215 xyP[2]=1.3;
216 zP[0]=117.3;
217 zP[1]=66.8;
218 zP[2]=1.2;
219 }
36be14b3 220 }
221 Double_t xySigma = xyP[0] + xyP[1]/TMath::Power(TMath::Abs(pt),xyP[2]);
222 Double_t xyMax = fNSigmaDCAxy*xySigma; //in micron
223 if((TMath::Abs(impactXY)*10000)>xyMax) return kFALSE;
224
225 Double_t zSigma = zP[0] + zP[1]/TMath::Power(TMath::Abs(pt),zP[2]);
226 Double_t zMax = fNSigmaDCAz*zSigma; //in micron
227 if((TMath::Abs(impactZ)*10000)>zMax) return kFALSE;
228
229 return kTRUE;
230}
231
232//________________________________________________________________________
4d668cf7 233Double_t AliAnalysisTaskSEITSsaSpectra::Eta2y(Double_t pt, Double_t m, Double_t eta) const {
234 // convert eta to y
36be14b3 235 Double_t mt = TMath::Sqrt(m*m + pt*pt);
236 return TMath::ASinH(pt/mt*TMath::SinH(eta));
237}
238
239
240//________________________________________________________________________
4d668cf7 241Double_t AliAnalysisTaskSEITSsaSpectra::BetheBloch(Double_t bg,Bool_t optMC) const{
36be14b3 242 // BB PHOBOS parametrization tuned by G. Ortona on 900 GeV pp data of 2009
243 Double_t par[5];
244 if(optMC){//Double_t par[5]={139.1,23.36,0.06052,0.2043,-0.0004999};
10ba520e 245
246 par[0]=-2.48; //new parameter from LHC10d2
247 par[1]=23.13;
248 par[2]=1.161;
249 par[3]=0.93;
250 par[4]=-1.2973;
251
36be14b3 252 }else {
253 //Double_t par[5]={5.33458e+04,1.65303e+01,2.60065e-03,3.59533e-04,7.51168e-05};}
254 par[0]=5.33458e+04;
255 par[1]=1.65303e+01;
256 par[2]=2.60065e-03;
257 par[3]=3.59533e-04;
258 par[4]=7.51168e-05;
259 }
260 Double_t beta = bg/TMath::Sqrt(1.+ bg*bg);
261 Double_t gamma=bg/beta;
262 Double_t eff=1.0;
263 if(bg<par[2]) eff=(bg-par[3])*(bg-par[3])+par[4];
264 else eff=(par[2]-par[3])*(par[2]-par[3])+par[4];
265 return (par[1]+2.0*TMath::Log(gamma)-beta*beta)*(par[0]/(beta*beta))*eff;
266}
267
268
269//________________________________________________________________________
270void AliAnalysisTaskSEITSsaSpectra::UserCreateOutputObjects(){
4d668cf7 271 // Create a TList with histograms and a TNtuple
272 // Called once
029e1796 273
36be14b3 274 fOutput = new TList();
275 fOutput->SetOwner();
276 fOutput->SetName("Spiderman");
277
c841b89f 278 fHistNEvents = new TH1F("fHistNEvents", "Number of processed events",7,-0.5,6.5);
36be14b3 279 fHistNEvents->Sumw2();
280 fHistNEvents->SetMinimum(0);
281 fOutput->Add(fHistNEvents);
282
1f75ac7d 283 fHistMult = new TH1F("fHistMult", "Event Multiplicity",3000,-0.5,2999.5);
c841b89f 284 fHistMult->Sumw2();
285 fHistMult->SetMinimum(0);
286 fOutput->Add(fHistMult);
287
1f75ac7d 288 fHistCen = new TH1F("fHistCen", "Event Centrality",101,-0.5,100.5);
289 fHistCen->Sumw2();
290 fHistCen->SetMinimum(0);
291 fOutput->Add(fHistCen);
292
36be14b3 293 fHistNTracks = new TH1F("fHistNTracks", "Number of ITSsa tracks",20,0.5,20.5);
294 fHistNTracks->Sumw2();
295 fHistNTracks->SetMinimum(0);
296 fOutput->Add(fHistNTracks);
297
c841b89f 298 fHistNTracksPos = new TH1F("fHistNTracksPos", "Number of positive ITSsa tracks",20,0.5,20.5);
299 fHistNTracksPos->Sumw2();
300 fHistNTracksPos->SetMinimum(0);
301 fOutput->Add(fHistNTracksPos);
302
303 fHistNTracksNeg = new TH1F("fHistNTracksNeg", "Number of negative ITSsa tracks",20,0.5,20.5);
304 fHistNTracksNeg->Sumw2();
305 fHistNTracksNeg->SetMinimum(0);
306 fOutput->Add(fHistNTracksNeg);
307
36be14b3 308 //binning for the histogram
309 const Int_t hnbins=400;
310 Double_t hxmin = 0.01;
311 Double_t hxmax = 10;
312 Double_t hlogxmin = TMath::Log10(hxmin);
313 Double_t hlogxmax = TMath::Log10(hxmax);
314 Double_t hbinwidth = (hlogxmax-hlogxmin)/hnbins;
315 Double_t hxbins[hnbins+1];
316 hxbins[0] = 0.01;
317 for (Int_t i=1;i<=hnbins;i++) {
318 hxbins[i] = hxmin + TMath::Power(10,hlogxmin+i*hbinwidth);
319 }
320
321 fHistDEDX = new TH2F("fHistDEDX","",hnbins,hxbins,900,0,1000);
322 fOutput->Add(fHistDEDX);
323
324 fHistDEDXdouble = new TH2F("fHistDEDXdouble","",500,-5,5,900,0,1000);
325 fOutput->Add(fHistDEDXdouble);
326
327
328 fHistBeforeEvSel = new TH1F("fHistBeforeEvSel","fHistBeforeEvSel",kNbins,fPtBinLimits);
329 fHistAfterEvSel = new TH1F("fHistAfterEvSel","fHistAfterEvSel",kNbins,fPtBinLimits);
330 fOutput->Add(fHistBeforeEvSel);
331 fOutput->Add(fHistAfterEvSel);
332
333
334
335 for(Int_t j=0;j<3;j++){
10ba520e 336 fHistPrimMCpos[j] = new TH1F(Form("fHistPrimMCpos%d",j),Form("fHistPrimMCpos%d",j),kNbins,fPtBinLimits);
337 fHistPrimMCneg[j] = new TH1F(Form("fHistPrimMCneg%d",j),Form("fHistPrimMCneg%d",j),kNbins,fPtBinLimits);
338 fOutput->Add(fHistPrimMCneg[j]);
339 fOutput->Add(fHistPrimMCpos[j]);
340 fHistSecStrMCpos[j] = new TH1F(Form("fHistSecStrMCpos%d",j),Form("fHistSecStrMCpos%d",j),kNbins,fPtBinLimits);
341 fHistSecStrMCneg[j] = new TH1F(Form("fHistSecStrMCneg%d",j),Form("fHistSecStrMCneg%d",j),kNbins,fPtBinLimits);
342 fOutput->Add(fHistSecStrMCneg[j]);
343 fOutput->Add(fHistSecStrMCpos[j]);
344 fHistSecMatMCpos[j] = new TH1F(Form("fHistSecMatMCpos%d",j),Form("fHistSecMatMCpos%d",j),kNbins,fPtBinLimits);
345 fHistSecMatMCneg[j] = new TH1F(Form("fHistSecMatMCneg%d",j),Form("fHistSecMatMCneg%d",j),kNbins,fPtBinLimits);
346 fOutput->Add(fHistSecMatMCneg[j]);
347 fOutput->Add(fHistSecMatMCpos[j]);
348 //
349 fHistPrimMCposBefEvSel[j] = new TH1F(Form("fHistPrimMCposBefEvSel%d",j),Form("fHistPrimMCposBefEvSel%d",j),kNbins,fPtBinLimits);
350 fHistPrimMCnegBefEvSel[j] = new TH1F(Form("fHistPrimMCnegBefEvSel%d",j),Form("fHistPrimMCnegBefEvSel%d",j),kNbins,fPtBinLimits);
351 fOutput->Add(fHistPrimMCnegBefEvSel[j]);
352 fOutput->Add(fHistPrimMCposBefEvSel[j]);
353 fHistSecStrMCposBefEvSel[j] = new TH1F(Form("fHistSecStrMCposBefEvSel%d",j),Form("fHistSecStrMCposBefEvSel%d",j),kNbins,fPtBinLimits);
354 fHistSecStrMCnegBefEvSel[j] = new TH1F(Form("fHistSecStrMCnegBefEvSel%d",j),Form("fHistSecStrMCnegBefEvSel%d",j),kNbins,fPtBinLimits);
355 fOutput->Add(fHistSecStrMCnegBefEvSel[j]);
356 fOutput->Add(fHistSecStrMCposBefEvSel[j]);
357 fHistSecMatMCposBefEvSel[j] = new TH1F(Form("fHistSecMatMCposBefEvSel%d",j),Form("fHistSecMatMCposBefEvSel%d",j),kNbins,fPtBinLimits);
358 fHistSecMatMCnegBefEvSel[j] = new TH1F(Form("fHistSecMatMCnegBefEvSel%d",j),Form("fHistSecMatMCnegBefEvSel%d",j),kNbins,fPtBinLimits);
359 fOutput->Add(fHistSecMatMCnegBefEvSel[j]);
360 fOutput->Add(fHistSecMatMCposBefEvSel[j]);
361 //
362 fHistPrimMCposReco[j] = new TH1F(Form("fHistPrimMCposReco%d",j),Form("fHistPrimMCposReco%d",j),kNbins,fPtBinLimits);
363 fHistPrimMCnegReco[j] = new TH1F(Form("fHistPrimMCnegReco%d",j),Form("fHistPrimMCnegReco%d",j),kNbins,fPtBinLimits);
364 fOutput->Add(fHistPrimMCnegReco[j]);
365 fOutput->Add(fHistPrimMCposReco[j]);
366 fHistSecStrMCposReco[j] = new TH1F(Form("fHistSecStrMCposReco%d",j),Form("fHistSecStrMCposReco%d",j),kNbins,fPtBinLimits);
367 fHistSecStrMCnegReco[j] = new TH1F(Form("fHistSecStrMCnegReco%d",j),Form("fHistSecStrMCnegReco%d",j),kNbins,fPtBinLimits);
368 fOutput->Add(fHistSecStrMCnegReco[j]);
369 fOutput->Add(fHistSecStrMCposReco[j]);
370 fHistSecMatMCposReco[j] = new TH1F(Form("fHistSecMatMCposReco%d",j),Form("fHistSecMatMCposReco%d",j),kNbins,fPtBinLimits);
371 fHistSecMatMCnegReco[j] = new TH1F(Form("fHistSecMatMCnegReco%d",j),Form("fHistSecMatMCnegReco%d",j),kNbins,fPtBinLimits);
372 fOutput->Add(fHistSecMatMCnegReco[j]);
373 fOutput->Add(fHistSecMatMCposReco[j]);
36be14b3 374
36be14b3 375 }
376
377 for(Int_t i=0; i<4; i++){
378 fHistCharge[i] = new TH1F(Form("fHistChargeLay%d",i),Form("fHistChargeLay%d",i),100,0,300);
379 fOutput->Add(fHistCharge[i]);
380 }
381
382 for(Int_t i=0; i<kNbins; i++){
383 fHistPosPi[i] = new TH1F(Form("fHistPosPi%d",i),Form("fHistPosPi%d",i),175,-3.5,3.5);
384 fHistPosK[i] = new TH1F(Form("fHistPosK%d",i),Form("fHistPosK%d",i),175,-3.5,3.5);
385 fHistPosP[i] = new TH1F(Form("fHistPosP%d",i),Form("fHistPosP%d",i),175,-3.5,3.5);
386 fHistNegPi[i] = new TH1F(Form("fHistNegPi%d",i),Form("fHistNegPi%d",i),175,-3.5,3.5);
387 fHistNegK[i] = new TH1F(Form("fHistNegK%d",i),Form("fHistNegK%d",i),175,-3.5,3.5);
388 fHistNegP[i] = new TH1F(Form("fHistNegP%d",i),Form("fHistNegP%d",i),175,-3.5,3.5);
389
10ba520e 390 fHistDCAPosPi[i] = new TH1F(Form("fHistDCAPosPi%d",i),Form("fHistDCAPosPi%d",i),2000,-1,1); //DCA distr. with NSigma PID
36be14b3 391 fHistDCAPosK[i] = new TH1F(Form("fHistDCAPosK%d",i),Form("fHistDCAPosK%d",i),2000,-1,1);
392 fHistDCAPosP[i] = new TH1F(Form("fHistDCAPosP%d",i),Form("fHistDCAPosP%d",i),2000,-1,1);
393 fHistDCANegPi[i] = new TH1F(Form("fHistDCANegPi%d",i),Form("fHistDCANegPi%d",i),2000,-1,1);
394 fHistDCANegK[i] = new TH1F(Form("fHistDCANegK%d",i),Form("fHistDCANegK%d",i),2000,-1,1);
395 fHistDCANegP[i] = new TH1F(Form("fHistDCANegP%d",i),Form("fHistDCANegP%d",i),2000,-1,1);
396
10ba520e 397 fHistMCPrimDCAPosPi[i] = new TH1F(Form("fHistMCPrimDCAPosPi%d",i),Form("fHistMCPrimDCAPosPi%d",i),2000,-1,1); //DCA distr. with MC truth
398 fHistMCPrimDCAPosK[i] = new TH1F(Form("fHistMCPrimDCAPosK%d",i),Form("fHistMCPrimDCAPosK%d",i),2000,-1,1);
399 fHistMCPrimDCAPosP[i] = new TH1F(Form("fHistMCPrimDCAPosP%d",i),Form("fHistMCPrimDCAPosP%d",i),2000,-1,1);
400 fHistMCPrimDCANegPi[i] = new TH1F(Form("fHistMCPrimDCANegPi%d",i),Form("fHistMCPrimDCANegPi%d",i),2000,-1,1);
401 fHistMCPrimDCANegK[i] = new TH1F(Form("fHistMCPrimDCANegK%d",i),Form("fHistMCPrimDCANegK%d",i),2000,-1,1);
402 fHistMCPrimDCANegP[i] = new TH1F(Form("fHistMCPrimDCANegP%d",i),Form("fHistMCPrimDCANegP%d",i),2000,-1,1);
403
404 fHistMCSecStDCAPosPi[i] = new TH1F(Form("fHistMCSecStDCAPosPi%d",i),Form("fHistMCSecStDCAPosPi%d",i),2000,-1,1); //DCA distr. with MC truth
405 fHistMCSecStDCAPosK[i] = new TH1F(Form("fHistMCSecStDCAPosK%d",i),Form("fHistMCSecStDCAPosK%d",i),2000,-1,1);
406 fHistMCSecStDCAPosP[i] = new TH1F(Form("fHistMCSecStDCAPosP%d",i),Form("fHistMCSecStDCAPosP%d",i),2000,-1,1);
407 fHistMCSecStDCANegPi[i] = new TH1F(Form("fHistMCSecStDCANegPi%d",i),Form("fHistMCSecStDCANegPi%d",i),2000,-1,1);
408 fHistMCSecStDCANegK[i] = new TH1F(Form("fHistMCSecStDCANegK%d",i),Form("fHistMCSecStDCANegK%d",i),2000,-1,1);
409 fHistMCSecStDCANegP[i] = new TH1F(Form("fHistMCSecStDCANegP%d",i),Form("fHistMCSecStDCANegP%d",i),2000,-1,1);
410
411 fHistMCSecMatDCAPosPi[i] = new TH1F(Form("fHistMCSecMatDCAPosPi%d",i),Form("fHistMCSecMatDCAPosPi%d",i),2000,-1,1); //DCA distr. with MC truth
412 fHistMCSecMatDCAPosK[i] = new TH1F(Form("fHistMCSecMatDCAPosK%d",i),Form("fHistMCSecMatDCAPosK%d",i),2000,-1,1);
413 fHistMCSecMatDCAPosP[i] = new TH1F(Form("fHistMCSecMatDCAPosP%d",i),Form("fHistMCSecMatDCAPosP%d",i),2000,-1,1);
414 fHistMCSecMatDCANegPi[i] = new TH1F(Form("fHistMCSecMatDCANegPi%d",i),Form("fHistMCSecMatDCANegPi%d",i),2000,-1,1);
415 fHistMCSecMatDCANegK[i] = new TH1F(Form("fHistMCSecMatDCANegK%d",i),Form("fHistMCSecMatDCANegK%d",i),2000,-1,1);
416 fHistMCSecMatDCANegP[i] = new TH1F(Form("fHistMCSecMatDCANegP%d",i),Form("fHistMCSecMatDCANegP%d",i),2000,-1,1);
417
36be14b3 418 fHistMCPosPi[i] = new TH1F(Form("fHistMCPosPi%d",i),Form("fHistMCPosPi%d",i),175,-3.5,3.5); //MC truth
419 fHistMCPosK[i] = new TH1F(Form("fHistMCPosK%d",i),Form("fHistMCPosK%d",i),175,-3.5,3.5);
420 fHistMCPosP[i] = new TH1F(Form("fHistMCPosP%d",i),Form("fHistMCPosP%d",i),175,-3.5,3.5);
421 fHistMCNegPi[i] = new TH1F(Form("fHistMCNegPi%d",i),Form("fHistMCNegPi%d",i),175,-3.5,3.5);
422 fHistMCNegK[i] = new TH1F(Form("fHistMCNegK%d",i),Form("fHistMCNegK%d",i),175,-3.5,3.5);
423 fHistMCNegP[i] = new TH1F(Form("fHistMCNegP%d",i),Form("fHistMCNegP%d",i),175,-3.5,3.5);
424 fOutput->Add(fHistPosPi[i]);
425 fOutput->Add(fHistPosK[i]);
426 fOutput->Add(fHistPosP[i]);
427 fOutput->Add(fHistNegPi[i]);
428 fOutput->Add(fHistNegK[i]);
429 fOutput->Add(fHistNegP[i]);
430
6b77d2c0 431 fOutput->Add(fHistDCAPosPi[i]); //DCA distr
36be14b3 432 fOutput->Add(fHistDCAPosK[i]);
433 fOutput->Add(fHistDCAPosP[i]);
434 fOutput->Add(fHistDCANegPi[i]);
435 fOutput->Add(fHistDCANegK[i]);
436 fOutput->Add(fHistDCANegP[i]);
437
10ba520e 438 fOutput->Add(fHistMCPrimDCAPosPi[i]);//DCA distr.
439 fOutput->Add(fHistMCPrimDCAPosK[i]);
440 fOutput->Add(fHistMCPrimDCAPosP[i]);
441 fOutput->Add(fHistMCPrimDCANegPi[i]);
442 fOutput->Add(fHistMCPrimDCANegK[i]);
443 fOutput->Add(fHistMCPrimDCANegP[i]);
444
445 fOutput->Add(fHistMCSecStDCAPosPi[i]);//DCA distr.
446 fOutput->Add(fHistMCSecStDCAPosK[i]);
447 fOutput->Add(fHistMCSecStDCAPosP[i]);
448 fOutput->Add(fHistMCSecStDCANegPi[i]);
449 fOutput->Add(fHistMCSecStDCANegK[i]);
450 fOutput->Add(fHistMCSecStDCANegP[i]);
451
452 fOutput->Add(fHistMCSecMatDCAPosPi[i]);//DCA distr.
453 fOutput->Add(fHistMCSecMatDCAPosK[i]);
454 fOutput->Add(fHistMCSecMatDCAPosP[i]);
455 fOutput->Add(fHistMCSecMatDCANegPi[i]);
456 fOutput->Add(fHistMCSecMatDCANegK[i]);
457 fOutput->Add(fHistMCSecMatDCANegP[i]);
458
36be14b3 459 fOutput->Add(fHistMCPosPi[i]);//MC truth
460 fOutput->Add(fHistMCPosK[i]);
461 fOutput->Add(fHistMCPosP[i]);
462 fOutput->Add(fHistMCNegPi[i]);
463 fOutput->Add(fHistMCNegK[i]);
464 fOutput->Add(fHistMCNegP[i]);
465 }
466
467 //NSigma Histos
468 for(Int_t j=0;j<3;j++){
10ba520e 469
470 fHistPosNSigmaMean[j] = new TH1F(Form("hHistPosNSigmaMean%d",j),Form("hHistPosNSigmaMean%d",j),kNbins,fPtBinLimits);
471 fHistNegNSigmaMean[j] = new TH1F(Form("hHistNegNSigmaMean%d",j),Form("hHistNegNSigmaMean%d",j),kNbins,fPtBinLimits);
472 fHistPosNSigmaMCMean[j] = new TH1F(Form("hHistPosNSigmaMCMean%d",j),Form("hHistPosNSigmaMCMean%d",j),kNbins,fPtBinLimits);
473 fHistNegNSigmaMCMean[j] = new TH1F(Form("hHistNegNSigmaMCMean%d",j),Form("hHistNegNSigmaMCMean%d",j),kNbins,fPtBinLimits);
474 fHistPosNSigmaPrimMean[j] = new TH1F(Form("hHistPosNSigmaPrimMean%d",j),Form("hHistPosNSigmaPrimMean%d",j),kNbins,fPtBinLimits);
475 fHistNegNSigmaPrimMean[j] = new TH1F(Form("hHistNegNSigmaPrimMean%d",j),Form("hHistNegNSigmaPrimMean%d",j),kNbins,fPtBinLimits);
476 fHistPosNSigmaPrimMCMean[j] = new TH1F(Form("hHistPosNSigmaPrimMCMean%d",j),Form("hHistPosNSigmaPrimMCMean%d",j),kNbins,fPtBinLimits);
477 fHistNegNSigmaPrimMCMean[j] = new TH1F(Form("hHistNegNSigmaPrimMCMean%d",j),Form("hHistNegNSigmaPrimMCMean%d",j),kNbins,fPtBinLimits);
478 fOutput->Add(fHistPosNSigmaMean[j]);
479 fOutput->Add(fHistNegNSigmaMean[j]);
480 fOutput->Add(fHistPosNSigmaMCMean[j]);
481 fOutput->Add(fHistNegNSigmaMCMean[j]);
482 fOutput->Add(fHistPosNSigmaPrimMean[j]);
483 fOutput->Add(fHistNegNSigmaPrimMean[j]);
484 fOutput->Add(fHistPosNSigmaPrimMCMean[j]);
485 fOutput->Add(fHistNegNSigmaPrimMCMean[j]);
486
36be14b3 487 fHistPosNSigma[j] = new TH1F(Form("hHistPosNSigma%d",j),Form("hHistPosNSigma%d",j),kNbins,fPtBinLimits);
488 fHistNegNSigma[j] = new TH1F(Form("hHistNegNSigma%d",j),Form("hHistNegNSigma%d",j),kNbins,fPtBinLimits);
10ba520e 489 fHistPosNSigmaMC[j] = new TH1F(Form("hHistPosNSigmaMC%d",j),Form("hHistPosNSigmaMC%d",j),kNbins,fPtBinLimits);
490 fHistNegNSigmaMC[j] = new TH1F(Form("hHistNegNSigmaMC%d",j),Form("hHistNegNSigmaMC%d",j),kNbins,fPtBinLimits);
36be14b3 491 fHistPosNSigmaPrim[j] = new TH1F(Form("hHistPosNSigmaPrim%d",j),Form("hHistPosNSigmaPrim%d",j),kNbins,fPtBinLimits);
492 fHistNegNSigmaPrim[j] = new TH1F(Form("hHistNegNSigmaPrim%d",j),Form("hHistNegNSigmaPrim%d",j),kNbins,fPtBinLimits);
493 fHistPosNSigmaPrimMC[j] = new TH1F(Form("hHistPosNSigmaPrimMC%d",j),Form("hHistPosNSigmaPrimMC%d",j),kNbins,fPtBinLimits);
494 fHistNegNSigmaPrimMC[j] = new TH1F(Form("hHistNegNSigmaPrimMC%d",j),Form("hHistNegNSigmaPrimMC%d",j),kNbins,fPtBinLimits);
495 fOutput->Add(fHistPosNSigma[j]);
496 fOutput->Add(fHistNegNSigma[j]);
10ba520e 497 fOutput->Add(fHistPosNSigmaMC[j]);
498 fOutput->Add(fHistNegNSigmaMC[j]);
36be14b3 499 fOutput->Add(fHistPosNSigmaPrim[j]);
500 fOutput->Add(fHistNegNSigmaPrim[j]);
501 fOutput->Add(fHistPosNSigmaPrimMC[j]);
502 fOutput->Add(fHistNegNSigmaPrimMC[j]);
10ba520e 503
36be14b3 504 }
505
c841b89f 506 fNtupleNSigma = new TNtuple("fNtupleNSigma","fNtupleNSigma","p:pt:dedx:ncls:sign:run:eta:impactXY:impactZ:isph:pdgcode:mfl:chi2ncls");
36be14b3 507 fOutput->Add(fNtupleNSigma);
508 fNtupleMC = new TNtuple("fNtupleMC","fNtupleMC","ptMC:pdgcode:signMC:etaMC:yMC:isph:evSel:run");
509 fOutput->Add(fNtupleMC);
7f8f788e 510
511 // Post output data.
512 PostData(1,fOutput);
513
36be14b3 514 Printf("end of CreateOutputObjects");
515}
516
517//________________________________________________________________________
518void AliAnalysisTaskSEITSsaSpectra::UserExec(Option_t *){
4d668cf7 519 // Main loop
520 // Called for each event
36be14b3 521
522 fESD=(AliESDEvent*)InputEvent();
523 if(!fESD) {
524 printf("AliAnalysisTaskSDDRP::Exec(): bad ESD\n");
525 return;
526 }
c841b89f 527 fHistNEvents->Fill(0);
7f2f1e04 528
1f75ac7d 529 ///////////////////////////////////////
36be14b3 530 //variables
10ba520e 531 Float_t pdgmass[4]={0.13957,0.493677,0.938272,1.8756}; //mass for pi, K, P (Gev/c^2)
36be14b3 532 Int_t listcode[3]={211,321,2212};//code for pi, K, P (Gev/c^2)
533 Double_t s[4];
10ba520e 534 Float_t ptMC=-999,yMC=-999;
535 Int_t code=-999, signMC=-999,isph=-999,mfl=-999;
36be14b3 536 Float_t impactXY=-999, impactZ=-999;
537 Int_t evSel=1;
538 AliESDtrack* track;
539 UInt_t status;
540 AliStack* stack=0;
541 TParticle *part=0;
542 TParticlePDG *pdgPart=0;
543
544 //Nsigma Method
545 Float_t resodedx[4];
546 if(fMC){
7d5daef5 547 resodedx[0]=0.0001; //default value
548 resodedx[1]=0.0001; //default value
549 resodedx[2]=0.108;
550 resodedx[3]=0.0998;
36be14b3 551 }else{
7d5daef5 552 resodedx[0]=0.0001; //default value
553 resodedx[1]=0.0001; //default value
554 resodedx[2]=0.116;
555 resodedx[3]=0.103;
36be14b3 556 }
557 /////////////////////
558 if(fMC){
559 AliMCEventHandler* eventHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
560 if (!eventHandler) {
561 Printf("ERROR: Could not retrieve MC event handler");
562 return;
563 }
564 AliMCEvent* mcEvent = eventHandler->MCEvent();
565 if (!mcEvent) {
566 Printf("ERROR: Could not retrieve MC event");
567 return;
568 }
569 stack = mcEvent->Stack();
570 if (!stack) {
571 printf("ERROR: stack not available\n");
572 return;
573 }
574 }
575 //flags for MC
576 Int_t nTrackMC=0;
577 if(stack) nTrackMC = stack->GetNtrack();
578 const AliESDVertex *vtx = fESD->GetPrimaryVertexSPD();
1f75ac7d 579
580///////////selection of the centrality or multiplicity bin
581
582 //selection on the event centrality
583 if(fHImode){
584 if(!(fLowCentrality<0.0)&&fUpCentrality>0.0)
585 {
586 AliCentrality *centrality = fESD->GetCentrality();
587 if(!centrality->IsEventInCentralityClass(fLowCentrality,fUpCentrality,"V0M"))
588 return;
589 Printf("Centrality of the event: %.1f",centrality->GetCentralityPercentile("V0M"));
590 Printf("Centrality cut: %.1f to %.1f",fLowCentrality,fUpCentrality);
591 fHistCen->Fill(centrality->GetCentralityPercentile("V0M"));
592 }
593 }
594
595 //selection on the event multiplicity based on global tracks
596 Int_t multiplicity = fesdTrackCutsMult->CountAcceptedTracks(fESD);
597 if(!fSPD){
598 if(fLowMult>-1)
599 {
600 if(multiplicity<fLowMult)
601 return;
602 }
603 if(fUpMult>-1)
604 {
605 if(multiplicity>fUpMult)
606 return;
607 }
608
609 Printf("Multiplicity of the event (global tracks) : %i",multiplicity);
610 Printf("Multiplicity cut (global tracks) : %i to %i",fLowMult,fUpMult);
611 fHistMult->Fill(multiplicity);
612 }
613
614 //multipicity selection based on SPD
615 if(fSPD){
616 Float_t spdCorr=-1.0;
617 const AliMultiplicity *mult = fESD->GetMultiplicity();
618 Float_t nClusters[6]={0.0,0.0,0.0,0.0,0.0,0.0};
619 for(Int_t ilay=0; ilay<6; ilay++)
620 {
621 nClusters[ilay] = (Float_t)mult->GetNumberOfITSClusters(ilay);
622 }
623 spdCorr = AliESDUtils::GetCorrSPD2(nClusters[1],vtx->GetZ());
624 {
625 if(fLowMult>-1)
626 {
627 if(((Int_t)spdCorr)<fLowMult)
628 {
629 return;
630 }
631 }
632 if(fUpMult>-1)
633 {
634 if(((Int_t)spdCorr)>fUpMult)
635 {
636 return;
637 }
638 }
639 }
640
641 Printf("Multiplicity of the event (SPD) : %i",(Int_t)spdCorr);
642 Printf("Multiplicity cut (SPD) : %i to %i",fLowMult,fUpMult);
643 fHistMult->Fill(spdCorr);
644 }
645
646
36be14b3 647
648 //event selection
649 fHistNEvents->Fill(1);
650 if(!vtx)evSel=0;
651 else{
652 fHistNEvents->Fill(2);
653 if(vtx->GetNContributors()<0) evSel=0;
654 else{
655 fHistNEvents->Fill(3);
656 if(TMath::Abs(vtx->GetZv())>10) evSel=0;
657 else{
658 fHistNEvents->Fill(4);
659 if(vtx->GetZRes()>0.5) evSel=0;
660 else{
661 fHistNEvents->Fill(5);
662 if(vtx->IsFromVertexerZ() && vtx->GetDispersion()>0.03) evSel=0;
663 else fHistNEvents->Fill(6);
664 }
665 }
666 }
667 }
1f75ac7d 668/////first loop on stack, before event selection, filling MC ntuple
36be14b3 669
670 for(Int_t imc=0; imc<nTrackMC; imc++){
671 part = stack->Particle(imc);
10ba520e 672 isph=1;
673 if(!stack->IsPhysicalPrimary(imc)) isph=0;
36be14b3 674 pdgPart = part->GetPDG();
10ba520e 675 if(!pdgPart)continue;
36be14b3 676 if(pdgPart->Charge()==0) continue; //no neutral particles
36be14b3 677 if(part->Energy() != TMath::Abs(part->Pz())) yMC = 0.5*TMath::Log((part->Energy()+part->Pz())/(part->Energy()-part->Pz()));
678 if(TMath::Abs(yMC) > fMaxY) continue; //rapidity cut
36be14b3 679 if(pdgPart->Charge()>0) signMC=1;
680 else signMC=-1;
681 ptMC=part->Pt();
682 code=pdgPart->PdgCode();
10ba520e 683 Int_t jpart=-1;
684 for(Int_t j=0; j<3; j++){
685 if(TMath::Abs(code)==listcode[j]){
686 jpart=j;
687 break;
688 }
689 }
690 Int_t indexMoth=part->GetFirstMother();
691 if(indexMoth>=0){
692 TParticle* moth = stack->Particle(indexMoth);
693 Float_t codemoth = TMath::Abs(moth->GetPdgCode());
694 mfl = Int_t (codemoth/ TMath::Power(10, Int_t(TMath::Log10(codemoth))));
695 }
36be14b3 696
36be14b3 697 //filling MC ntuple
698 if(TMath::Abs(code)==211 || TMath::Abs(code)==321 || TMath::Abs(code)==2212){
699 Float_t xntMC[8];
700 Int_t indexMC=0;
701 xntMC[indexMC++]=(Float_t)ptMC;
702 xntMC[indexMC++]=(Float_t)code;
703 xntMC[indexMC++]=(Float_t)signMC;
704 xntMC[indexMC++]=(Float_t)part->Eta();
705 xntMC[indexMC++]=(Float_t)yMC;
706 xntMC[indexMC++]=(Float_t)isph;
707 xntMC[indexMC++]=(Float_t)evSel;
708 xntMC[indexMC++]=(Float_t)fESD->GetRunNumber();
10ba520e 709
36be14b3 710 if(fFillNtuple) fNtupleMC->Fill(xntMC);
711 }
712
10ba520e 713 if(jpart>=0){
714 if(stack->IsPhysicalPrimary(imc)){
715 if(signMC>0) fHistPrimMCposBefEvSel[jpart]->Fill(TMath::Abs(ptMC));
716 else fHistPrimMCnegBefEvSel[jpart]->Fill(TMath::Abs(ptMC));
717 if(evSel==1){
718 if(signMC>0) fHistPrimMCpos[jpart]->Fill(TMath::Abs(ptMC));
719 else fHistPrimMCneg[jpart]->Fill(TMath::Abs(ptMC));
720 }
721 }else{
722 if(mfl==3){ // strangeness
723 if(signMC>0) fHistSecStrMCposBefEvSel[jpart]->Fill(TMath::Abs(ptMC));
724 else fHistSecStrMCnegBefEvSel[jpart]->Fill(TMath::Abs(ptMC));
725 if(evSel==1){
726 if(signMC>0) fHistSecStrMCpos[jpart]->Fill(TMath::Abs(ptMC));
727 else fHistSecStrMCneg[jpart]->Fill(TMath::Abs(ptMC));
728 }
729 }else{
730 if(signMC>0) fHistSecMatMCposBefEvSel[jpart]->Fill(TMath::Abs(ptMC));
731 else fHistSecMatMCnegBefEvSel[jpart]->Fill(TMath::Abs(ptMC));
732 if(evSel==1){
733 if(signMC>0) fHistSecMatMCpos[jpart]->Fill(TMath::Abs(ptMC));
734 else fHistSecMatMCneg[jpart]->Fill(TMath::Abs(ptMC));
735 }
36be14b3 736 }
737 }
10ba520e 738 }
739 }
740
36be14b3 741
10ba520e 742 if(evSel==0)return; //event selection
1f75ac7d 743
744
36be14b3 745 //loop on tracks
746 for (Int_t iTrack=0; iTrack<fESD->GetNumberOfTracks(); iTrack++) {
10ba520e 747 isph=-999;
36be14b3 748 code=-999;
749 mfl=-999;
10ba520e 750
36be14b3 751 track = (AliESDtrack*)fESD->GetTrack(iTrack);
752 if (!track) continue;
753
754 //track selection
c841b89f 755 Int_t countBinTrk=1;
756 TString label;
757
758 label="no selection";
759 fHistNTracks->Fill(countBinTrk);
760 if(track->GetSign()>0)fHistNTracksPos->Fill(countBinTrk);
761 if(track->GetSign()<0)fHistNTracksNeg->Fill(countBinTrk);
762 fHistNTracks->GetXaxis()->SetBinLabel(fHistNTracks->FindBin(countBinTrk),label.Data());
763 fHistNTracksPos->GetXaxis()->SetBinLabel(fHistNTracksPos->FindBin(countBinTrk),label.Data());
764 fHistNTracksNeg->GetXaxis()->SetBinLabel(fHistNTracksNeg->FindBin(countBinTrk),label.Data());
765 countBinTrk++;
766
36be14b3 767 status=track->GetStatus();
768 if((status&AliESDtrack::kITSpureSA)==0) continue; //its standalone
c841b89f 769
770 label="ITSsa";
771 fHistNTracks->Fill(countBinTrk);
772 if(track->GetSign()>0)fHistNTracksPos->Fill(countBinTrk);
773 if(track->GetSign()<0)fHistNTracksNeg->Fill(countBinTrk);
774 fHistNTracks->GetXaxis()->SetBinLabel(fHistNTracks->FindBin(countBinTrk),label.Data());
775 fHistNTracksPos->GetXaxis()->SetBinLabel(fHistNTracksPos->FindBin(countBinTrk),label.Data());
776 fHistNTracksNeg->GetXaxis()->SetBinLabel(fHistNTracksNeg->FindBin(countBinTrk),label.Data());
777 countBinTrk++;
778
36be14b3 779 if((status&AliESDtrack::kITSrefit)==0) continue; //its refit
c841b89f 780
781 label="ITSrefit";
782 fHistNTracks->Fill(countBinTrk);
783 if(track->GetSign()>0)fHistNTracksPos->Fill(countBinTrk);
784 if(track->GetSign()<0)fHistNTracksNeg->Fill(countBinTrk);
785 fHistNTracks->GetXaxis()->SetBinLabel(fHistNTracks->FindBin(countBinTrk),label.Data());
786 fHistNTracksPos->GetXaxis()->SetBinLabel(fHistNTracksPos->FindBin(countBinTrk),label.Data());
787 fHistNTracksNeg->GetXaxis()->SetBinLabel(fHistNTracksNeg->FindBin(countBinTrk),label.Data());
788 countBinTrk++;
789
029e1796 790 if(TMath::Abs(track->GetSign())<0.0001) continue; //no neutral particles
c841b89f 791
792 label="neutral particle";
793 fHistNTracks->Fill(countBinTrk);
794 if(track->GetSign()>0)fHistNTracksPos->Fill(countBinTrk);
795 if(track->GetSign()<0)fHistNTracksNeg->Fill(countBinTrk);
796 fHistNTracks->GetXaxis()->SetBinLabel(fHistNTracks->FindBin(countBinTrk),label.Data());
797 fHistNTracksPos->GetXaxis()->SetBinLabel(fHistNTracksPos->FindBin(countBinTrk),label.Data());
798 fHistNTracksNeg->GetXaxis()->SetBinLabel(fHistNTracksNeg->FindBin(countBinTrk),label.Data());
799 countBinTrk++;
10ba520e 800
36be14b3 801 //cluster in ITS
802 UInt_t clumap = track->GetITSClusterMap();
803 Int_t nSPD=0;
804 for(Int_t il=0; il<2; il++) if(TESTBIT(clumap,il)) nSPD++;
805 if(nSPD<fMinSPDPts) continue;
c841b89f 806
807 label="SPDcls";
808 fHistNTracks->Fill(countBinTrk);
809 if(track->GetSign()>0)fHistNTracksPos->Fill(countBinTrk);
810 if(track->GetSign()<0)fHistNTracksNeg->Fill(countBinTrk);
811 fHistNTracks->GetXaxis()->SetBinLabel(fHistNTracks->FindBin(countBinTrk),label.Data());
812 fHistNTracksPos->GetXaxis()->SetBinLabel(fHistNTracksPos->FindBin(countBinTrk),label.Data());
813 fHistNTracksNeg->GetXaxis()->SetBinLabel(fHistNTracksNeg->FindBin(countBinTrk),label.Data());
814 countBinTrk++;
815
36be14b3 816 Int_t count=0;
817 for(Int_t j=2;j<6;j++) if(TESTBIT(clumap,j)) count++;
818 if(count<fMinNdEdxSamples) continue; //at least 3 points on SSD/SDD
c841b89f 819
820 label="SDD+SSD cls";
821 fHistNTracks->Fill(countBinTrk);
822 if(track->GetSign()>0)fHistNTracksPos->Fill(countBinTrk);
823 if(track->GetSign()<0)fHistNTracksNeg->Fill(countBinTrk);
824 fHistNTracks->GetXaxis()->SetBinLabel(fHistNTracks->FindBin(countBinTrk),label.Data());
825 fHistNTracksPos->GetXaxis()->SetBinLabel(fHistNTracksPos->FindBin(countBinTrk),label.Data());
826 fHistNTracksNeg->GetXaxis()->SetBinLabel(fHistNTracksNeg->FindBin(countBinTrk),label.Data());
827 countBinTrk++;
828
36be14b3 829 //chisquare/nclusters
830 Int_t nclu=nSPD+count;
831 if(track->GetITSchi2()/nclu > fMaxChi2Clu) continue;
c841b89f 832
833 label="chi2/ncls";
834 fHistNTracks->Fill(countBinTrk);
835 if(track->GetSign()>0)fHistNTracksPos->Fill(countBinTrk);
836 if(track->GetSign()<0)fHistNTracksNeg->Fill(countBinTrk);
837 fHistNTracks->GetXaxis()->SetBinLabel(fHistNTracks->FindBin(countBinTrk),label.Data());
838 fHistNTracksPos->GetXaxis()->SetBinLabel(fHistNTracksPos->FindBin(countBinTrk),label.Data());
839 fHistNTracksNeg->GetXaxis()->SetBinLabel(fHistNTracksNeg->FindBin(countBinTrk),label.Data());
840 countBinTrk++;
841
36be14b3 842 //pseudorapidity and rapidity
dac69cb0 843 if(TMath::Abs(track->Eta()) > fEtaRange) continue;
c841b89f 844
845 label="eta";
846 fHistNTracks->Fill(countBinTrk);
847 if(track->GetSign()>0)fHistNTracksPos->Fill(countBinTrk);
848 if(track->GetSign()<0)fHistNTracksNeg->Fill(countBinTrk);
849 fHistNTracks->GetXaxis()->SetBinLabel(fHistNTracks->FindBin(countBinTrk),label.Data());
850 fHistNTracksPos->GetXaxis()->SetBinLabel(fHistNTracksPos->FindBin(countBinTrk),label.Data());
851 fHistNTracksNeg->GetXaxis()->SetBinLabel(fHistNTracksNeg->FindBin(countBinTrk),label.Data());
852 countBinTrk++;
853
36be14b3 854 //truncated mean
855 //if(fMC) for(Int_t j=0;j<2;j++) s[j]*=3.34/5.43;//correction for SDD miscalibration of the MCpass4
856 track->GetITSdEdxSamples(s);
857 Double_t dedx = CookdEdx(s);
858 if(dedx<0) continue;
36be14b3 859
c841b89f 860 label="de/dx<0";
861 fHistNTracks->Fill(countBinTrk);
862 if(track->GetSign()>0)fHistNTracksPos->Fill(countBinTrk);
863 if(track->GetSign()<0)fHistNTracksNeg->Fill(countBinTrk);
864 fHistNTracks->GetXaxis()->SetBinLabel(fHistNTracks->FindBin(countBinTrk),label.Data());
865 fHistNTracksPos->GetXaxis()->SetBinLabel(fHistNTracksPos->FindBin(countBinTrk),label.Data());
866 fHistNTracksNeg->GetXaxis()->SetBinLabel(fHistNTracksNeg->FindBin(countBinTrk),label.Data());
867 countBinTrk++;
868
36be14b3 869 Float_t pt = track->Pt();
870 Int_t theBin=-1;
871 for(Int_t m=0; m<kNbins; m++){
872 if(TMath::Abs(pt) > fPtBinLimits[m] && TMath::Abs(pt) < fPtBinLimits[m+1]){
873 theBin=m;
874 break;
875 }
876 }
877 track->GetImpactParameters(impactXY, impactZ);
10ba520e 878
36be14b3 879 //Filling Ntuple
880 //information from the MC kinematics
881 if(fMC){
10ba520e 882 if(track->GetLabel()<0)isph=-1;
36be14b3 883 if(track->GetLabel()>=0){
884 part = (TParticle*)stack->Particle(track->GetLabel());
885 pdgPart = part->GetPDG();
886 code = pdgPart->PdgCode();
10ba520e 887 if(stack->IsPhysicalPrimary(track->GetLabel())) isph=1;
36be14b3 888 else{
10ba520e 889 isph=0;
890 Int_t indexMoth=part->GetFirstMother();
891 if(indexMoth>=0){
892 TParticle* moth = stack->Particle(indexMoth);
893 Float_t codemoth = TMath::Abs(moth->GetPdgCode());
894 mfl = Int_t (codemoth/ TMath::Power(10, Int_t(TMath::Log10(codemoth))));
895 }
896 }
897
898 //Filling DCA distribution with MC truth
899
900 if(theBin>=0 && theBin<kNbins){
901 if(isph==1){//primaries in MC
902 if(track->GetSign()>0){
903 if(TMath::Abs(code)==listcode[0]) fHistMCPrimDCAPosPi[theBin]->Fill(impactXY);
904 if(TMath::Abs(code)==listcode[1]) fHistMCPrimDCAPosK[theBin]->Fill(impactXY);
905 if(TMath::Abs(code)==listcode[2]) fHistMCPrimDCAPosP[theBin]->Fill(impactXY);
906 }else{
907 if(TMath::Abs(code)==listcode[0]) fHistMCPrimDCANegPi[theBin]->Fill(impactXY);
908 if(TMath::Abs(code)==listcode[1]) fHistMCPrimDCANegK[theBin]->Fill(impactXY);
909 if(TMath::Abs(code)==listcode[2]) fHistMCPrimDCANegP[theBin]->Fill(impactXY);
910 }
911 }
912
913 if(isph==0){//primaries in MC
914 if(mfl==3){
915 if(track->GetSign()>0){
916 if(TMath::Abs(code)==listcode[0]) fHistMCSecStDCAPosPi[theBin]->Fill(impactXY);
917 if(TMath::Abs(code)==listcode[1]) fHistMCSecStDCAPosK[theBin]->Fill(impactXY);
918 if(TMath::Abs(code)==listcode[2]) fHistMCSecStDCAPosP[theBin]->Fill(impactXY);
919 }else{
920 if(TMath::Abs(code)==listcode[0]) fHistMCSecStDCANegPi[theBin]->Fill(impactXY);
921 if(TMath::Abs(code)==listcode[1]) fHistMCSecStDCANegK[theBin]->Fill(impactXY);
922 if(TMath::Abs(code)==listcode[2]) fHistMCSecStDCANegP[theBin]->Fill(impactXY);
923 }
924 }else{
925 if(track->GetSign()>0){
926 if(TMath::Abs(code)==listcode[0]) fHistMCSecMatDCAPosPi[theBin]->Fill(impactXY);
927 if(TMath::Abs(code)==listcode[1]) fHistMCSecMatDCAPosK[theBin]->Fill(impactXY);
928 if(TMath::Abs(code)==listcode[2]) fHistMCSecMatDCAPosP[theBin]->Fill(impactXY);
929 }else{
930 if(TMath::Abs(code)==listcode[0]) fHistMCSecMatDCANegPi[theBin]->Fill(impactXY);
931 if(TMath::Abs(code)==listcode[1]) fHistMCSecMatDCANegK[theBin]->Fill(impactXY);
932 if(TMath::Abs(code)==listcode[2]) fHistMCSecMatDCANegP[theBin]->Fill(impactXY);
933 }
934 }
935 }
36be14b3 936 }
937 }
938 }
c841b89f 939 Float_t xnt[13];
36be14b3 940 Int_t index=0;
941 xnt[index++]=(Float_t)track->GetP();
942 xnt[index++]=(Float_t)track->Pt();
943 xnt[index++]=(Float_t)dedx;
944 xnt[index++]=(Float_t)count;
945 xnt[index++]=(Float_t)track->GetSign();
946 xnt[index++]=(Float_t)fESD->GetRunNumber();
947 xnt[index++]=(Float_t)track->Eta();
948 xnt[index++]=(Float_t)impactXY;
949 xnt[index++]=(Float_t)impactZ;
950 xnt[index++]=(Float_t)isph;
951 xnt[index++]=(Float_t)code;
c841b89f 952 xnt[index++]=(Float_t)mfl;
953 xnt[index]=(Float_t)track->GetITSchi2()/nclu;
36be14b3 954
955 if(fFillNtuple) fNtupleNSigma->Fill(xnt);
956
957
958
959 //Compute y and bb
10ba520e 960 Double_t y[4],bbtheo[4],logdiff[4];
36be14b3 961 Float_t p=track->GetP();
962 if(fMC && fSmearMC){
963 dedx=fRandGener->Gaus(dedx,fSmeardEdx*dedx);
964 p=fRandGener->Gaus(p,fSmearP*p);
965 }
966
10ba520e 967 for(Int_t i=0;i<4;i++){
36be14b3 968 y[i] = Eta2y(pt,pdgmass[i],track->Eta());
969 bbtheo[i]=BetheBloch(p/pdgmass[i],fMC);
970 logdiff[i]=TMath::Log(dedx) - TMath::Log(bbtheo[i]);
971 }
36be14b3 972 Int_t resocls=(Int_t)count-1;
10ba520e 973
974 //NSigma Method, with asymmetric bands
975
976 Int_t minPosMean=-1;
977 for(Int_t isp=0; isp<3; isp++){
978 if(dedx<bbtheo[0])continue;
979 Double_t bb=bbtheo[isp];
980 if(dedx<bb){
981 Double_t bbdistance=TMath::Abs((bbtheo[isp]-bbtheo[isp-1])/2);
982 Double_t nsigma=TMath::Abs((dedx-bb)/bbdistance);
983 if(nsigma<1.)minPosMean=isp;
984 }
985 else{
986 Double_t bbdistance=TMath::Abs((bbtheo[isp]-bbtheo[isp+1])/2);
987 Double_t nsigma=TMath::Abs((dedx-bb)/bbdistance);
988 if(nsigma<1.)minPosMean=isp;
989 }
990 }
991 if(dedx<bbtheo[0] && TMath::Abs((dedx-bbtheo[0])/(resodedx[resocls]*bbtheo[0]))<2)minPosMean=0;
992
993 //NSigma method with simmetric bands
36be14b3 994
995 Double_t nsigmas[3];
996 Double_t min=999999.;
997 Int_t minPos=-1;
998 for(Int_t isp=0; isp<3; isp++){
999 Double_t bb=bbtheo[isp];
1000 nsigmas[isp]=TMath::Abs((dedx-bb)/(resodedx[resocls]*bb));
1001 if(nsigmas[isp]<min){
1002 min=nsigmas[isp];
1003 minPos=isp;
1004 }
1005 }
10ba520e 1006
1007 // y calculation
1008 Double_t yPartMean=y[minPosMean];
36be14b3 1009 Double_t yPart=y[minPos];
1010
10ba520e 1011 if(yPartMean<fMaxY){
1012 //DCA distributions, before the DCA cuts, based on asymmetrinc nsigma approach
36be14b3 1013 if(theBin>=0 && theBin<kNbins){
1014 if(track->GetSign()>0){
10ba520e 1015 if(minPosMean==0) fHistDCAPosPi[theBin]->Fill(impactXY);
1016 else if(minPosMean==1) fHistDCAPosK[theBin]->Fill(impactXY);
1017 else if(minPosMean==2) fHistDCAPosP[theBin]->Fill(impactXY);
36be14b3 1018 }else{
10ba520e 1019 if(minPosMean==0) fHistDCANegPi[theBin]->Fill(impactXY);
1020 else if(minPosMean==1) fHistDCANegK[theBin]->Fill(impactXY);
1021 else if(minPosMean==2) fHistDCANegP[theBin]->Fill(impactXY);
36be14b3 1022 }
1023 }
1024 }
1025
1026 //DCA cut on xy and z
1027 if(!DCAcut(impactXY,impactZ,pt,fMC)) continue;
c841b89f 1028
1029 label="DCA";
1030 fHistNTracks->Fill(countBinTrk);
1031 if(track->GetSign()>0)fHistNTracksPos->Fill(countBinTrk);
1032 if(track->GetSign()<0)fHistNTracksNeg->Fill(countBinTrk);
1033 fHistNTracks->GetXaxis()->SetBinLabel(fHistNTracks->FindBin(countBinTrk),label.Data());
1034 fHistNTracksPos->GetXaxis()->SetBinLabel(fHistNTracksPos->FindBin(countBinTrk),label.Data());
1035 fHistNTracksNeg->GetXaxis()->SetBinLabel(fHistNTracksNeg->FindBin(countBinTrk),label.Data());
1036 countBinTrk++;
10ba520e 1037
1038 //Filling Histos for Reco Efficiency
1039 //information from the MC kinematics
1040 if(fMC){
1041 if(track->GetLabel()<0)isph=-1;
1042 if(track->GetLabel()>=0){
1043 part = (TParticle*)stack->Particle(track->GetLabel());
1044 pdgPart = part->GetPDG();
1045 code = pdgPart->PdgCode();
1046 Int_t jpart=-1;
1047 for(Int_t j=0; j<3; j++){
1048 if(TMath::Abs(code)==listcode[j]){
1049 jpart=j;
1050 break;
1051 }
1052 }
1053 if(jpart<0) continue;
1054 if(pdgPart->Charge()>0) signMC=1;
1055 else signMC=-1;
1056 ptMC=part->Pt();
1057 if(stack->IsPhysicalPrimary(track->GetLabel())){
1058 if(signMC>0) fHistPrimMCposReco[jpart]->Fill(TMath::Abs(ptMC));
1059 else fHistPrimMCnegReco[jpart]->Fill(TMath::Abs(ptMC));
1060 }else{
1061 Int_t indexMoth=part->GetFirstMother();
1062 if(indexMoth>=0){
1063 TParticle* moth = stack->Particle(indexMoth);
1064 Float_t codemoth = TMath::Abs(moth->GetPdgCode());
1065 mfl = Int_t (codemoth/ TMath::Power(10, Int_t(TMath::Log10(codemoth))));
1066 }
1067 if(mfl==3){ // strangeness
1068 if(signMC>0) fHistSecStrMCposReco[jpart]->Fill(TMath::Abs(ptMC));
1069 else fHistSecStrMCnegReco[jpart]->Fill(TMath::Abs(ptMC));
1070 }else{
1071 if(signMC>0) fHistSecMatMCposReco[jpart]->Fill(TMath::Abs(ptMC));
1072 else fHistSecMatMCnegReco[jpart]->Fill(TMath::Abs(ptMC));
1073 }
1074 }
1075 }
1076 }
1077
1078 //Nsigma histos with MC truth
1079
1080 //asymmetric approach
1081 if(yPartMean<fMaxY){
1082 //nsigma histos
1083 if(track->GetSign()>0) fHistPosNSigmaMean[minPosMean]->Fill(pt);
1084 else fHistNegNSigmaMean[minPosMean]->Fill(pt);
1085 if(fMC){
1086 //nsigma histos with MC truth on PID
1087 if(TMath::Abs(code)==listcode[minPosMean]){
1088 if(track->GetSign()>0) fHistPosNSigmaMCMean[minPosMean]->Fill(pt);
1089 else fHistNegNSigmaMCMean[minPosMean]->Fill(pt);
1090 }
1091 //nsigma histos with MC truth on IsPhysicalPrimary
1092 if(isph==1){
1093 if(track->GetSign()>0) fHistPosNSigmaPrimMean[minPosMean]->Fill(pt);
1094 else fHistNegNSigmaPrimMean[minPosMean]->Fill(pt);
1095 //nsigma histos with MC truth on IsPhysicalPrimary and PID
1096 if(TMath::Abs(code)==listcode[minPosMean]){
1097 if(track->GetSign()>0) fHistPosNSigmaPrimMCMean[minPosMean]->Fill(pt);
1098 else fHistNegNSigmaPrimMCMean[minPosMean]->Fill(pt);
1099 }
1100 }
1101 }
1102 }
36be14b3 1103
10ba520e 1104 //symmetric bands
36be14b3 1105 if(min<fMinNSigma && yPart<fMaxY){
10ba520e 1106 //nsigma histos
36be14b3 1107 if(track->GetSign()>0) fHistPosNSigma[minPos]->Fill(pt);
1108 else fHistNegNSigma[minPos]->Fill(pt);
1109 if(fMC){
10ba520e 1110 //nsigma histos with MC truth on PID
1111 if(TMath::Abs(code)==listcode[minPos]){
1112 if(track->GetSign()>0) fHistPosNSigmaMC[minPos]->Fill(pt);
1113 else fHistNegNSigmaMC[minPos]->Fill(pt);
1114 }
1115 //nsigma histos with MC truth on IsPhysicalPrimary
36be14b3 1116 if(isph==1){
1117 if(track->GetSign()>0) fHistPosNSigmaPrim[minPos]->Fill(pt);
1118 else fHistNegNSigmaPrim[minPos]->Fill(pt);
10ba520e 1119 //nsigma histos with MC truth on IsPhysicalPrimary and PID
36be14b3 1120 if(TMath::Abs(code)==listcode[minPos]){
1121 if(track->GetSign()>0) fHistPosNSigmaPrimMC[minPos]->Fill(pt);
1122 else fHistNegNSigmaPrimMC[minPos]->Fill(pt);
1123 }
1124 }
1125 }
1126 }
1127
10ba520e 1128
1129 //integral approach histograms
36be14b3 1130 if(theBin>=0 && theBin<kNbins){
1131 if(track->GetSign()>0){
1132 if(TMath::Abs(y[0]) < fMaxY)fHistPosPi[theBin]->Fill(logdiff[0]);
1133 if(TMath::Abs(y[1]) < fMaxY)fHistPosK[theBin]->Fill(logdiff[1]);
1134 if(TMath::Abs(y[2]) < fMaxY)fHistPosP[theBin]->Fill(logdiff[2]);
1135 if(fMC){
1136 if((TMath::Abs(y[0])<fMaxY) && (TMath::Abs(code)==211)) fHistMCPosPi[theBin]->Fill(logdiff[0]);
1137 if((TMath::Abs(y[1])<fMaxY) && (TMath::Abs(code)==321)) fHistMCPosK[theBin]->Fill(logdiff[1]);
1138 if((TMath::Abs(y[2])<fMaxY) && (TMath::Abs(code)==2212)) fHistMCPosP[theBin]->Fill(logdiff[2]);
1139 }
1140 }else{
1141 if(TMath::Abs(y[0]) < fMaxY)fHistNegPi[theBin]->Fill(logdiff[0]);
1142 if(TMath::Abs(y[1]) < fMaxY)fHistNegK[theBin]->Fill(logdiff[1]);
1143 if(TMath::Abs(y[2]) < fMaxY)fHistNegP[theBin]->Fill(logdiff[2]);
1144 if(fMC){
1145 if((TMath::Abs(y[0])<fMaxY) && (TMath::Abs(code)==211)) fHistMCNegPi[theBin]->Fill(logdiff[0]);
1146 if((TMath::Abs(y[1])<fMaxY) && (TMath::Abs(code)==321)) fHistMCNegK[theBin]->Fill(logdiff[1]);
1147 if((TMath::Abs(y[2])<fMaxY) && (TMath::Abs(code)==2212)) fHistMCNegP[theBin]->Fill(logdiff[2]);
1148 }
1149 }
1150 }
1151
1152
1153 //fill propaganda plot with dedx
1154 fHistDEDX->Fill(track->GetP(),dedx);
1155 fHistDEDXdouble->Fill(track->GetP()*track->GetSign(),dedx);
1156
1157 //fill charge distribution histo to check the calibration
1158 for(Int_t j=0;j<4;j++){
1159 if(s[j]<5) continue;
1160 fHistCharge[j]->Fill(s[j]);
1161 }
1162 }
1163
1164 // Post output data.
1165 PostData(1,fOutput);
1166 Printf("............. end of Exec");
1167}
1168
1169//________________________________________________________________________
1170void AliAnalysisTaskSEITSsaSpectra::Terminate(Option_t *) {
4d668cf7 1171 // Merge output
1172 // Called once at the end of the query
36be14b3 1173
1174 fOutput = dynamic_cast<TList*>(GetOutputData(1));
1175 if (!fOutput) {
1176 printf("ERROR: fOutput not available\n");
1177 return;
1178 }
1179 fHistNEvents = dynamic_cast<TH1F*>(fOutput->FindObject("fHistNEvents"));
c841b89f 1180 fHistMult = dynamic_cast<TH1F*>(fOutput->FindObject("fHistMult"));
1f75ac7d 1181 fHistCen = dynamic_cast<TH1F*>(fOutput->FindObject("fHistCen"));
36be14b3 1182 fHistNTracks = dynamic_cast<TH1F*>(fOutput->FindObject("fHistNTracks"));
c841b89f 1183 fHistNTracksPos = dynamic_cast<TH1F*>(fOutput->FindObject("fHistNTracksPos"));
1184 fHistNTracksNeg = dynamic_cast<TH1F*>(fOutput->FindObject("fHistNTracksNeg"));
36be14b3 1185 fHistDEDX = dynamic_cast<TH2F*>(fOutput->FindObject("fHistDEDX"));
1186 fHistDEDXdouble = dynamic_cast<TH2F*>(fOutput->FindObject("fHistDEDXdouble"));
1187
1188 fHistBeforeEvSel = dynamic_cast<TH1F*>(fOutput->FindObject("fHistBeforeEvSel"));
1189 fHistAfterEvSel = dynamic_cast<TH1F*>(fOutput->FindObject("fHistAfterEvSel"));
1190
1191
1192 for(Int_t j=0;j<3;j++){
10ba520e 1193 if(fMC){
1194 fHistPrimMCpos[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistPrimMCpos%d",j)));
1195 fHistPrimMCneg[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistPrimMCneg%d",j)));
1196 fHistSecStrMCpos[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistSecStrMCpos%d",j)));
1197 fHistSecStrMCneg[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistSecStrMCneg%d",j)));
1198 fHistSecMatMCpos[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistSecMatMCpos%d",j)));
1199 fHistSecMatMCneg[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistSecMatMCneg%d",j)));
1200 //
1201 fHistPrimMCposBefEvSel[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistPrimMCposBefEvSel%d",j)));
1202 fHistPrimMCnegBefEvSel[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistPrimMCnegBefEvSel%d",j)));
1203 fHistSecStrMCposBefEvSel[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistSecStrMCposBefEvSel%d",j)));
1204 fHistSecStrMCnegBefEvSel[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistSecStrMCnegBefEvSel%d",j)));
1205 fHistSecMatMCposBefEvSel[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistSecMatMCposBefEvSel%d",j)));
1206 fHistSecMatMCnegBefEvSel[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistSecMatMCnegBefEvSel%d",j)));
1207 //
1208 fHistPrimMCposReco[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistPrimMCposReco%d",j)));
1209 fHistPrimMCnegReco[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistPrimMCnegReco%d",j)));
1210 fHistSecStrMCposReco[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistSecStrMCposReco%d",j)));
1211 fHistSecStrMCnegReco[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistSecStrMCnegReco%d",j)));
1212 fHistSecMatMCposReco[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistSecMatMCposReco%d",j)));
1213 fHistSecMatMCnegReco[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistSecMatMCnegReco%d",j)));
1214 }
36be14b3 1215 }
36be14b3 1216
36be14b3 1217 for(Int_t i=0; i<4; i++){
1218 fHistCharge[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistChargeLay%d",i)));
1219 }
1220
1221 for(Int_t i=0; i<kNbins; i++){
1222 fHistPosPi[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistPosPi%d",i)));
1223 fHistPosK[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistPosK%d",i)));
1224 fHistPosP[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistPosP%d",i)));
1225 fHistNegPi[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistNegPi%d",i)));
1226 fHistNegK[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistNegK%d",i)));
1227 fHistNegP[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistNegP%d",i)));
1228
1229 fHistDCAPosPi[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistDCAPosPi%d",i)));
1230 fHistDCAPosK[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistDCAPosK%d",i)));
1231 fHistDCAPosP[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistDCAPosP%d",i)));
1232 fHistDCANegPi[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistDCANegPi%d",i)));
1233 fHistDCANegK[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistDCANegK%d",i)));
1234 fHistDCANegP[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistDCANegP%d",i)));
1235
1236
1237 if(fMC){
10ba520e 1238
1239 fHistMCPrimDCAPosPi[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCPrimDCAPosPi%d",i)));
1240 fHistMCPrimDCAPosK[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCPrimDCAPosK%d",i)));
1241 fHistMCPrimDCAPosP[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCPrimDCAPosP%d",i)));
1242 fHistMCPrimDCANegPi[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCPrimDCANegPi%d",i)));
1243 fHistMCPrimDCANegK[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCPrimDCANegK%d",i)));
1244 fHistMCPrimDCANegP[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCPrimDCANegP%d",i)));
1245
1246 fHistMCSecStDCAPosPi[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCSecStDCAPosPi%d",i)));
1247 fHistMCSecStDCAPosK[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCSecStDCAPosK%d",i)));
1248 fHistMCSecStDCAPosP[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCSecStDCAPosP%d",i)));
1249 fHistMCSecStDCANegPi[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCSecStDCANegPi%d",i)));
1250 fHistMCSecStDCANegK[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCSecStDCANegK%d",i)));
1251 fHistMCSecStDCANegP[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCSecStDCANegP%d",i)));
1252
1253 fHistMCSecMatDCAPosPi[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCSecMatDCAPosPi%d",i)));
1254 fHistMCSecMatDCAPosK[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCSecMatDCAPosK%d",i)));
1255 fHistMCSecMatDCAPosP[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCSecMatDCAPosP%d",i)));
1256 fHistMCSecMatDCANegPi[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCSecMatDCANegPi%d",i)));
1257 fHistMCSecMatDCANegK[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCSecMatDCANegK%d",i)));
1258 fHistMCSecMatDCANegP[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCSecMatDCANegP%d",i)));
1259
36be14b3 1260 fHistMCPosPi[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCPosPi%d",i)));
1261 fHistMCPosK[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCPosK%d",i)));
1262 fHistMCPosP[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCPosP%d",i)));
1263 fHistMCNegPi[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCNegPi%d",i)));
1264 fHistMCNegK[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCNegK%d",i)));
1265 fHistMCNegP[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCNegP%d",i)));
1266 }
1267 }
10ba520e 1268
36be14b3 1269 for(Int_t j=0;j<3;j++){
10ba520e 1270 fHistPosNSigmaMean[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("hHistPosNSigmaMean%d",j)));
1271 fHistNegNSigmaMean[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("hHistNegNSigmaMean%d",j)));
1272 fHistPosNSigmaMCMean[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("hHistPosNSigmaMCMean%d",j)));
1273 fHistNegNSigmaMCMean[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("hHistNegNSigmaMCMean%d",j)));
1274 fHistPosNSigmaPrimMean[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("hHistPosNSigmaPrimMean%d",j)));
1275 fHistNegNSigmaPrimMean[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("hHistNegNSigmaPrimMean%d",j)));
1276 fHistPosNSigmaPrimMCMean[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("hHistPosNSigmaPrimMCMean%d",j)));
1277 fHistNegNSigmaPrimMCMean[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("hHistNegNSigmaPrimMCMean%d",j)));
1278
36be14b3 1279 fHistPosNSigma[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("hHistPosNSigma%d",j)));
1280 fHistNegNSigma[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("hHistNegNSigma%d",j)));
10ba520e 1281 fHistPosNSigmaMC[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("hHistPosNSigmaMC%d",j)));
1282 fHistNegNSigmaMC[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("hHistNegNSigmaMC%d",j)));
36be14b3 1283 fHistPosNSigmaPrim[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("hHistPosNSigmaPrim%d",j)));
1284 fHistNegNSigmaPrim[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("hHistNegNSigmaPrim%d",j)));
1285 fHistPosNSigmaPrimMC[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("hHistPosNSigmaPrimMC%d",j)));
1286 fHistNegNSigmaPrimMC[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("hHistNegNSigmaPrimMC%d",j)));
10ba520e 1287
36be14b3 1288 }
1289
1290 fNtupleNSigma = dynamic_cast<TNtuple*>(fOutput->FindObject("fNtupleNSigma"));
1291 fNtupleMC = dynamic_cast<TNtuple*>(fOutput->FindObject("fNtupleMC"));
1292
1293 Printf("end of Terminate");
1294 return;
1295}