]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG2/SPECTRA/AliAnalysisTaskSEITSsaSpectra.cxx
Added option for initialization of peak position (Francesco) + new parameterization...
[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];
5c112001 534 Float_t ptMC=-999;
10ba520e 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
5c112001 677 Float_t yMC=-999.;
36be14b3 678 if(part->Energy() != TMath::Abs(part->Pz())) yMC = 0.5*TMath::Log((part->Energy()+part->Pz())/(part->Energy()-part->Pz()));
679 if(TMath::Abs(yMC) > fMaxY) continue; //rapidity cut
36be14b3 680 if(pdgPart->Charge()>0) signMC=1;
681 else signMC=-1;
682 ptMC=part->Pt();
683 code=pdgPart->PdgCode();
10ba520e 684 Int_t jpart=-1;
685 for(Int_t j=0; j<3; j++){
686 if(TMath::Abs(code)==listcode[j]){
687 jpart=j;
688 break;
689 }
690 }
691 Int_t indexMoth=part->GetFirstMother();
692 if(indexMoth>=0){
693 TParticle* moth = stack->Particle(indexMoth);
694 Float_t codemoth = TMath::Abs(moth->GetPdgCode());
695 mfl = Int_t (codemoth/ TMath::Power(10, Int_t(TMath::Log10(codemoth))));
696 }
36be14b3 697
36be14b3 698 //filling MC ntuple
699 if(TMath::Abs(code)==211 || TMath::Abs(code)==321 || TMath::Abs(code)==2212){
700 Float_t xntMC[8];
701 Int_t indexMC=0;
702 xntMC[indexMC++]=(Float_t)ptMC;
703 xntMC[indexMC++]=(Float_t)code;
704 xntMC[indexMC++]=(Float_t)signMC;
705 xntMC[indexMC++]=(Float_t)part->Eta();
706 xntMC[indexMC++]=(Float_t)yMC;
707 xntMC[indexMC++]=(Float_t)isph;
708 xntMC[indexMC++]=(Float_t)evSel;
709 xntMC[indexMC++]=(Float_t)fESD->GetRunNumber();
10ba520e 710
36be14b3 711 if(fFillNtuple) fNtupleMC->Fill(xntMC);
712 }
713
10ba520e 714 if(jpart>=0){
715 if(stack->IsPhysicalPrimary(imc)){
716 if(signMC>0) fHistPrimMCposBefEvSel[jpart]->Fill(TMath::Abs(ptMC));
717 else fHistPrimMCnegBefEvSel[jpart]->Fill(TMath::Abs(ptMC));
718 if(evSel==1){
719 if(signMC>0) fHistPrimMCpos[jpart]->Fill(TMath::Abs(ptMC));
720 else fHistPrimMCneg[jpart]->Fill(TMath::Abs(ptMC));
721 }
722 }else{
723 if(mfl==3){ // strangeness
724 if(signMC>0) fHistSecStrMCposBefEvSel[jpart]->Fill(TMath::Abs(ptMC));
725 else fHistSecStrMCnegBefEvSel[jpart]->Fill(TMath::Abs(ptMC));
726 if(evSel==1){
727 if(signMC>0) fHistSecStrMCpos[jpart]->Fill(TMath::Abs(ptMC));
728 else fHistSecStrMCneg[jpart]->Fill(TMath::Abs(ptMC));
729 }
730 }else{
731 if(signMC>0) fHistSecMatMCposBefEvSel[jpart]->Fill(TMath::Abs(ptMC));
732 else fHistSecMatMCnegBefEvSel[jpart]->Fill(TMath::Abs(ptMC));
733 if(evSel==1){
734 if(signMC>0) fHistSecMatMCpos[jpart]->Fill(TMath::Abs(ptMC));
735 else fHistSecMatMCneg[jpart]->Fill(TMath::Abs(ptMC));
736 }
36be14b3 737 }
738 }
10ba520e 739 }
740 }
741
36be14b3 742
10ba520e 743 if(evSel==0)return; //event selection
1f75ac7d 744
745
36be14b3 746 //loop on tracks
747 for (Int_t iTrack=0; iTrack<fESD->GetNumberOfTracks(); iTrack++) {
10ba520e 748 isph=-999;
36be14b3 749 code=-999;
750 mfl=-999;
10ba520e 751
36be14b3 752 track = (AliESDtrack*)fESD->GetTrack(iTrack);
753 if (!track) continue;
754
755 //track selection
c841b89f 756 Int_t countBinTrk=1;
757 TString label;
758
759 label="no selection";
760 fHistNTracks->Fill(countBinTrk);
761 if(track->GetSign()>0)fHistNTracksPos->Fill(countBinTrk);
762 if(track->GetSign()<0)fHistNTracksNeg->Fill(countBinTrk);
763 fHistNTracks->GetXaxis()->SetBinLabel(fHistNTracks->FindBin(countBinTrk),label.Data());
764 fHistNTracksPos->GetXaxis()->SetBinLabel(fHistNTracksPos->FindBin(countBinTrk),label.Data());
765 fHistNTracksNeg->GetXaxis()->SetBinLabel(fHistNTracksNeg->FindBin(countBinTrk),label.Data());
766 countBinTrk++;
767
36be14b3 768 status=track->GetStatus();
769 if((status&AliESDtrack::kITSpureSA)==0) continue; //its standalone
c841b89f 770
771 label="ITSsa";
772 fHistNTracks->Fill(countBinTrk);
773 if(track->GetSign()>0)fHistNTracksPos->Fill(countBinTrk);
774 if(track->GetSign()<0)fHistNTracksNeg->Fill(countBinTrk);
775 fHistNTracks->GetXaxis()->SetBinLabel(fHistNTracks->FindBin(countBinTrk),label.Data());
776 fHistNTracksPos->GetXaxis()->SetBinLabel(fHistNTracksPos->FindBin(countBinTrk),label.Data());
777 fHistNTracksNeg->GetXaxis()->SetBinLabel(fHistNTracksNeg->FindBin(countBinTrk),label.Data());
778 countBinTrk++;
779
36be14b3 780 if((status&AliESDtrack::kITSrefit)==0) continue; //its refit
c841b89f 781
782 label="ITSrefit";
783 fHistNTracks->Fill(countBinTrk);
784 if(track->GetSign()>0)fHistNTracksPos->Fill(countBinTrk);
785 if(track->GetSign()<0)fHistNTracksNeg->Fill(countBinTrk);
786 fHistNTracks->GetXaxis()->SetBinLabel(fHistNTracks->FindBin(countBinTrk),label.Data());
787 fHistNTracksPos->GetXaxis()->SetBinLabel(fHistNTracksPos->FindBin(countBinTrk),label.Data());
788 fHistNTracksNeg->GetXaxis()->SetBinLabel(fHistNTracksNeg->FindBin(countBinTrk),label.Data());
789 countBinTrk++;
790
029e1796 791 if(TMath::Abs(track->GetSign())<0.0001) continue; //no neutral particles
c841b89f 792
793 label="neutral particle";
794 fHistNTracks->Fill(countBinTrk);
795 if(track->GetSign()>0)fHistNTracksPos->Fill(countBinTrk);
796 if(track->GetSign()<0)fHistNTracksNeg->Fill(countBinTrk);
797 fHistNTracks->GetXaxis()->SetBinLabel(fHistNTracks->FindBin(countBinTrk),label.Data());
798 fHistNTracksPos->GetXaxis()->SetBinLabel(fHistNTracksPos->FindBin(countBinTrk),label.Data());
799 fHistNTracksNeg->GetXaxis()->SetBinLabel(fHistNTracksNeg->FindBin(countBinTrk),label.Data());
800 countBinTrk++;
10ba520e 801
36be14b3 802 //cluster in ITS
803 UInt_t clumap = track->GetITSClusterMap();
804 Int_t nSPD=0;
805 for(Int_t il=0; il<2; il++) if(TESTBIT(clumap,il)) nSPD++;
806 if(nSPD<fMinSPDPts) continue;
c841b89f 807
808 label="SPDcls";
809 fHistNTracks->Fill(countBinTrk);
810 if(track->GetSign()>0)fHistNTracksPos->Fill(countBinTrk);
811 if(track->GetSign()<0)fHistNTracksNeg->Fill(countBinTrk);
812 fHistNTracks->GetXaxis()->SetBinLabel(fHistNTracks->FindBin(countBinTrk),label.Data());
813 fHistNTracksPos->GetXaxis()->SetBinLabel(fHistNTracksPos->FindBin(countBinTrk),label.Data());
814 fHistNTracksNeg->GetXaxis()->SetBinLabel(fHistNTracksNeg->FindBin(countBinTrk),label.Data());
815 countBinTrk++;
816
36be14b3 817 Int_t count=0;
818 for(Int_t j=2;j<6;j++) if(TESTBIT(clumap,j)) count++;
819 if(count<fMinNdEdxSamples) continue; //at least 3 points on SSD/SDD
c841b89f 820
821 label="SDD+SSD cls";
822 fHistNTracks->Fill(countBinTrk);
823 if(track->GetSign()>0)fHistNTracksPos->Fill(countBinTrk);
824 if(track->GetSign()<0)fHistNTracksNeg->Fill(countBinTrk);
825 fHistNTracks->GetXaxis()->SetBinLabel(fHistNTracks->FindBin(countBinTrk),label.Data());
826 fHistNTracksPos->GetXaxis()->SetBinLabel(fHistNTracksPos->FindBin(countBinTrk),label.Data());
827 fHistNTracksNeg->GetXaxis()->SetBinLabel(fHistNTracksNeg->FindBin(countBinTrk),label.Data());
828 countBinTrk++;
829
36be14b3 830 //chisquare/nclusters
831 Int_t nclu=nSPD+count;
832 if(track->GetITSchi2()/nclu > fMaxChi2Clu) continue;
c841b89f 833
834 label="chi2/ncls";
835 fHistNTracks->Fill(countBinTrk);
836 if(track->GetSign()>0)fHistNTracksPos->Fill(countBinTrk);
837 if(track->GetSign()<0)fHistNTracksNeg->Fill(countBinTrk);
838 fHistNTracks->GetXaxis()->SetBinLabel(fHistNTracks->FindBin(countBinTrk),label.Data());
839 fHistNTracksPos->GetXaxis()->SetBinLabel(fHistNTracksPos->FindBin(countBinTrk),label.Data());
840 fHistNTracksNeg->GetXaxis()->SetBinLabel(fHistNTracksNeg->FindBin(countBinTrk),label.Data());
841 countBinTrk++;
842
36be14b3 843 //pseudorapidity and rapidity
dac69cb0 844 if(TMath::Abs(track->Eta()) > fEtaRange) continue;
c841b89f 845
846 label="eta";
847 fHistNTracks->Fill(countBinTrk);
848 if(track->GetSign()>0)fHistNTracksPos->Fill(countBinTrk);
849 if(track->GetSign()<0)fHistNTracksNeg->Fill(countBinTrk);
850 fHistNTracks->GetXaxis()->SetBinLabel(fHistNTracks->FindBin(countBinTrk),label.Data());
851 fHistNTracksPos->GetXaxis()->SetBinLabel(fHistNTracksPos->FindBin(countBinTrk),label.Data());
852 fHistNTracksNeg->GetXaxis()->SetBinLabel(fHistNTracksNeg->FindBin(countBinTrk),label.Data());
853 countBinTrk++;
854
36be14b3 855 //truncated mean
856 //if(fMC) for(Int_t j=0;j<2;j++) s[j]*=3.34/5.43;//correction for SDD miscalibration of the MCpass4
857 track->GetITSdEdxSamples(s);
858 Double_t dedx = CookdEdx(s);
859 if(dedx<0) continue;
36be14b3 860
c841b89f 861 label="de/dx<0";
862 fHistNTracks->Fill(countBinTrk);
863 if(track->GetSign()>0)fHistNTracksPos->Fill(countBinTrk);
864 if(track->GetSign()<0)fHistNTracksNeg->Fill(countBinTrk);
865 fHistNTracks->GetXaxis()->SetBinLabel(fHistNTracks->FindBin(countBinTrk),label.Data());
866 fHistNTracksPos->GetXaxis()->SetBinLabel(fHistNTracksPos->FindBin(countBinTrk),label.Data());
867 fHistNTracksNeg->GetXaxis()->SetBinLabel(fHistNTracksNeg->FindBin(countBinTrk),label.Data());
868 countBinTrk++;
869
36be14b3 870 Float_t pt = track->Pt();
871 Int_t theBin=-1;
872 for(Int_t m=0; m<kNbins; m++){
873 if(TMath::Abs(pt) > fPtBinLimits[m] && TMath::Abs(pt) < fPtBinLimits[m+1]){
874 theBin=m;
875 break;
876 }
877 }
878 track->GetImpactParameters(impactXY, impactZ);
10ba520e 879
36be14b3 880 //Filling Ntuple
881 //information from the MC kinematics
882 if(fMC){
10ba520e 883 if(track->GetLabel()<0)isph=-1;
36be14b3 884 if(track->GetLabel()>=0){
885 part = (TParticle*)stack->Particle(track->GetLabel());
886 pdgPart = part->GetPDG();
887 code = pdgPart->PdgCode();
10ba520e 888 if(stack->IsPhysicalPrimary(track->GetLabel())) isph=1;
36be14b3 889 else{
10ba520e 890 isph=0;
891 Int_t indexMoth=part->GetFirstMother();
892 if(indexMoth>=0){
893 TParticle* moth = stack->Particle(indexMoth);
894 Float_t codemoth = TMath::Abs(moth->GetPdgCode());
895 mfl = Int_t (codemoth/ TMath::Power(10, Int_t(TMath::Log10(codemoth))));
896 }
897 }
898
899 //Filling DCA distribution with MC truth
900
901 if(theBin>=0 && theBin<kNbins){
902 if(isph==1){//primaries in MC
903 if(track->GetSign()>0){
904 if(TMath::Abs(code)==listcode[0]) fHistMCPrimDCAPosPi[theBin]->Fill(impactXY);
905 if(TMath::Abs(code)==listcode[1]) fHistMCPrimDCAPosK[theBin]->Fill(impactXY);
906 if(TMath::Abs(code)==listcode[2]) fHistMCPrimDCAPosP[theBin]->Fill(impactXY);
907 }else{
908 if(TMath::Abs(code)==listcode[0]) fHistMCPrimDCANegPi[theBin]->Fill(impactXY);
909 if(TMath::Abs(code)==listcode[1]) fHistMCPrimDCANegK[theBin]->Fill(impactXY);
910 if(TMath::Abs(code)==listcode[2]) fHistMCPrimDCANegP[theBin]->Fill(impactXY);
911 }
912 }
913
914 if(isph==0){//primaries in MC
915 if(mfl==3){
916 if(track->GetSign()>0){
917 if(TMath::Abs(code)==listcode[0]) fHistMCSecStDCAPosPi[theBin]->Fill(impactXY);
918 if(TMath::Abs(code)==listcode[1]) fHistMCSecStDCAPosK[theBin]->Fill(impactXY);
919 if(TMath::Abs(code)==listcode[2]) fHistMCSecStDCAPosP[theBin]->Fill(impactXY);
920 }else{
921 if(TMath::Abs(code)==listcode[0]) fHistMCSecStDCANegPi[theBin]->Fill(impactXY);
922 if(TMath::Abs(code)==listcode[1]) fHistMCSecStDCANegK[theBin]->Fill(impactXY);
923 if(TMath::Abs(code)==listcode[2]) fHistMCSecStDCANegP[theBin]->Fill(impactXY);
924 }
925 }else{
926 if(track->GetSign()>0){
927 if(TMath::Abs(code)==listcode[0]) fHistMCSecMatDCAPosPi[theBin]->Fill(impactXY);
928 if(TMath::Abs(code)==listcode[1]) fHistMCSecMatDCAPosK[theBin]->Fill(impactXY);
929 if(TMath::Abs(code)==listcode[2]) fHistMCSecMatDCAPosP[theBin]->Fill(impactXY);
930 }else{
931 if(TMath::Abs(code)==listcode[0]) fHistMCSecMatDCANegPi[theBin]->Fill(impactXY);
932 if(TMath::Abs(code)==listcode[1]) fHistMCSecMatDCANegK[theBin]->Fill(impactXY);
933 if(TMath::Abs(code)==listcode[2]) fHistMCSecMatDCANegP[theBin]->Fill(impactXY);
934 }
935 }
936 }
36be14b3 937 }
938 }
939 }
c841b89f 940 Float_t xnt[13];
36be14b3 941 Int_t index=0;
942 xnt[index++]=(Float_t)track->GetP();
943 xnt[index++]=(Float_t)track->Pt();
944 xnt[index++]=(Float_t)dedx;
945 xnt[index++]=(Float_t)count;
946 xnt[index++]=(Float_t)track->GetSign();
947 xnt[index++]=(Float_t)fESD->GetRunNumber();
948 xnt[index++]=(Float_t)track->Eta();
949 xnt[index++]=(Float_t)impactXY;
950 xnt[index++]=(Float_t)impactZ;
951 xnt[index++]=(Float_t)isph;
952 xnt[index++]=(Float_t)code;
c841b89f 953 xnt[index++]=(Float_t)mfl;
954 xnt[index]=(Float_t)track->GetITSchi2()/nclu;
36be14b3 955
956 if(fFillNtuple) fNtupleNSigma->Fill(xnt);
957
958
959
960 //Compute y and bb
10ba520e 961 Double_t y[4],bbtheo[4],logdiff[4];
36be14b3 962 Float_t p=track->GetP();
963 if(fMC && fSmearMC){
964 dedx=fRandGener->Gaus(dedx,fSmeardEdx*dedx);
965 p=fRandGener->Gaus(p,fSmearP*p);
966 }
967
10ba520e 968 for(Int_t i=0;i<4;i++){
36be14b3 969 y[i] = Eta2y(pt,pdgmass[i],track->Eta());
970 bbtheo[i]=BetheBloch(p/pdgmass[i],fMC);
971 logdiff[i]=TMath::Log(dedx) - TMath::Log(bbtheo[i]);
972 }
36be14b3 973 Int_t resocls=(Int_t)count-1;
10ba520e 974
975 //NSigma Method, with asymmetric bands
976
977 Int_t minPosMean=-1;
978 for(Int_t isp=0; isp<3; isp++){
979 if(dedx<bbtheo[0])continue;
980 Double_t bb=bbtheo[isp];
981 if(dedx<bb){
982 Double_t bbdistance=TMath::Abs((bbtheo[isp]-bbtheo[isp-1])/2);
983 Double_t nsigma=TMath::Abs((dedx-bb)/bbdistance);
984 if(nsigma<1.)minPosMean=isp;
985 }
986 else{
987 Double_t bbdistance=TMath::Abs((bbtheo[isp]-bbtheo[isp+1])/2);
988 Double_t nsigma=TMath::Abs((dedx-bb)/bbdistance);
989 if(nsigma<1.)minPosMean=isp;
990 }
991 }
992 if(dedx<bbtheo[0] && TMath::Abs((dedx-bbtheo[0])/(resodedx[resocls]*bbtheo[0]))<2)minPosMean=0;
993
994 //NSigma method with simmetric bands
36be14b3 995
996 Double_t nsigmas[3];
997 Double_t min=999999.;
998 Int_t minPos=-1;
999 for(Int_t isp=0; isp<3; isp++){
1000 Double_t bb=bbtheo[isp];
1001 nsigmas[isp]=TMath::Abs((dedx-bb)/(resodedx[resocls]*bb));
1002 if(nsigmas[isp]<min){
1003 min=nsigmas[isp];
1004 minPos=isp;
1005 }
1006 }
10ba520e 1007
1008 // y calculation
1009 Double_t yPartMean=y[minPosMean];
36be14b3 1010 Double_t yPart=y[minPos];
1011
10ba520e 1012 if(yPartMean<fMaxY){
1013 //DCA distributions, before the DCA cuts, based on asymmetrinc nsigma approach
36be14b3 1014 if(theBin>=0 && theBin<kNbins){
1015 if(track->GetSign()>0){
10ba520e 1016 if(minPosMean==0) fHistDCAPosPi[theBin]->Fill(impactXY);
1017 else if(minPosMean==1) fHistDCAPosK[theBin]->Fill(impactXY);
1018 else if(minPosMean==2) fHistDCAPosP[theBin]->Fill(impactXY);
36be14b3 1019 }else{
10ba520e 1020 if(minPosMean==0) fHistDCANegPi[theBin]->Fill(impactXY);
1021 else if(minPosMean==1) fHistDCANegK[theBin]->Fill(impactXY);
1022 else if(minPosMean==2) fHistDCANegP[theBin]->Fill(impactXY);
36be14b3 1023 }
1024 }
1025 }
1026
1027 //DCA cut on xy and z
1028 if(!DCAcut(impactXY,impactZ,pt,fMC)) continue;
c841b89f 1029
1030 label="DCA";
1031 fHistNTracks->Fill(countBinTrk);
1032 if(track->GetSign()>0)fHistNTracksPos->Fill(countBinTrk);
1033 if(track->GetSign()<0)fHistNTracksNeg->Fill(countBinTrk);
1034 fHistNTracks->GetXaxis()->SetBinLabel(fHistNTracks->FindBin(countBinTrk),label.Data());
1035 fHistNTracksPos->GetXaxis()->SetBinLabel(fHistNTracksPos->FindBin(countBinTrk),label.Data());
1036 fHistNTracksNeg->GetXaxis()->SetBinLabel(fHistNTracksNeg->FindBin(countBinTrk),label.Data());
1037 countBinTrk++;
10ba520e 1038
1039 //Filling Histos for Reco Efficiency
1040 //information from the MC kinematics
1041 if(fMC){
1042 if(track->GetLabel()<0)isph=-1;
1043 if(track->GetLabel()>=0){
1044 part = (TParticle*)stack->Particle(track->GetLabel());
1045 pdgPart = part->GetPDG();
1046 code = pdgPart->PdgCode();
1047 Int_t jpart=-1;
1048 for(Int_t j=0; j<3; j++){
1049 if(TMath::Abs(code)==listcode[j]){
1050 jpart=j;
1051 break;
1052 }
1053 }
1054 if(jpart<0) continue;
1055 if(pdgPart->Charge()>0) signMC=1;
1056 else signMC=-1;
1057 ptMC=part->Pt();
1058 if(stack->IsPhysicalPrimary(track->GetLabel())){
1059 if(signMC>0) fHistPrimMCposReco[jpart]->Fill(TMath::Abs(ptMC));
1060 else fHistPrimMCnegReco[jpart]->Fill(TMath::Abs(ptMC));
1061 }else{
1062 Int_t indexMoth=part->GetFirstMother();
1063 if(indexMoth>=0){
1064 TParticle* moth = stack->Particle(indexMoth);
1065 Float_t codemoth = TMath::Abs(moth->GetPdgCode());
1066 mfl = Int_t (codemoth/ TMath::Power(10, Int_t(TMath::Log10(codemoth))));
1067 }
1068 if(mfl==3){ // strangeness
1069 if(signMC>0) fHistSecStrMCposReco[jpart]->Fill(TMath::Abs(ptMC));
1070 else fHistSecStrMCnegReco[jpart]->Fill(TMath::Abs(ptMC));
1071 }else{
1072 if(signMC>0) fHistSecMatMCposReco[jpart]->Fill(TMath::Abs(ptMC));
1073 else fHistSecMatMCnegReco[jpart]->Fill(TMath::Abs(ptMC));
1074 }
1075 }
1076 }
1077 }
1078
1079 //Nsigma histos with MC truth
1080
1081 //asymmetric approach
1082 if(yPartMean<fMaxY){
1083 //nsigma histos
1084 if(track->GetSign()>0) fHistPosNSigmaMean[minPosMean]->Fill(pt);
1085 else fHistNegNSigmaMean[minPosMean]->Fill(pt);
1086 if(fMC){
1087 //nsigma histos with MC truth on PID
1088 if(TMath::Abs(code)==listcode[minPosMean]){
1089 if(track->GetSign()>0) fHistPosNSigmaMCMean[minPosMean]->Fill(pt);
1090 else fHistNegNSigmaMCMean[minPosMean]->Fill(pt);
1091 }
1092 //nsigma histos with MC truth on IsPhysicalPrimary
1093 if(isph==1){
1094 if(track->GetSign()>0) fHistPosNSigmaPrimMean[minPosMean]->Fill(pt);
1095 else fHistNegNSigmaPrimMean[minPosMean]->Fill(pt);
1096 //nsigma histos with MC truth on IsPhysicalPrimary and PID
1097 if(TMath::Abs(code)==listcode[minPosMean]){
1098 if(track->GetSign()>0) fHistPosNSigmaPrimMCMean[minPosMean]->Fill(pt);
1099 else fHistNegNSigmaPrimMCMean[minPosMean]->Fill(pt);
1100 }
1101 }
1102 }
1103 }
36be14b3 1104
10ba520e 1105 //symmetric bands
36be14b3 1106 if(min<fMinNSigma && yPart<fMaxY){
10ba520e 1107 //nsigma histos
36be14b3 1108 if(track->GetSign()>0) fHistPosNSigma[minPos]->Fill(pt);
1109 else fHistNegNSigma[minPos]->Fill(pt);
1110 if(fMC){
10ba520e 1111 //nsigma histos with MC truth on PID
1112 if(TMath::Abs(code)==listcode[minPos]){
1113 if(track->GetSign()>0) fHistPosNSigmaMC[minPos]->Fill(pt);
1114 else fHistNegNSigmaMC[minPos]->Fill(pt);
1115 }
1116 //nsigma histos with MC truth on IsPhysicalPrimary
36be14b3 1117 if(isph==1){
1118 if(track->GetSign()>0) fHistPosNSigmaPrim[minPos]->Fill(pt);
1119 else fHistNegNSigmaPrim[minPos]->Fill(pt);
10ba520e 1120 //nsigma histos with MC truth on IsPhysicalPrimary and PID
36be14b3 1121 if(TMath::Abs(code)==listcode[minPos]){
1122 if(track->GetSign()>0) fHistPosNSigmaPrimMC[minPos]->Fill(pt);
1123 else fHistNegNSigmaPrimMC[minPos]->Fill(pt);
1124 }
1125 }
1126 }
1127 }
1128
10ba520e 1129
1130 //integral approach histograms
36be14b3 1131 if(theBin>=0 && theBin<kNbins){
1132 if(track->GetSign()>0){
1133 if(TMath::Abs(y[0]) < fMaxY)fHistPosPi[theBin]->Fill(logdiff[0]);
1134 if(TMath::Abs(y[1]) < fMaxY)fHistPosK[theBin]->Fill(logdiff[1]);
1135 if(TMath::Abs(y[2]) < fMaxY)fHistPosP[theBin]->Fill(logdiff[2]);
1136 if(fMC){
1137 if((TMath::Abs(y[0])<fMaxY) && (TMath::Abs(code)==211)) fHistMCPosPi[theBin]->Fill(logdiff[0]);
1138 if((TMath::Abs(y[1])<fMaxY) && (TMath::Abs(code)==321)) fHistMCPosK[theBin]->Fill(logdiff[1]);
1139 if((TMath::Abs(y[2])<fMaxY) && (TMath::Abs(code)==2212)) fHistMCPosP[theBin]->Fill(logdiff[2]);
1140 }
1141 }else{
1142 if(TMath::Abs(y[0]) < fMaxY)fHistNegPi[theBin]->Fill(logdiff[0]);
1143 if(TMath::Abs(y[1]) < fMaxY)fHistNegK[theBin]->Fill(logdiff[1]);
1144 if(TMath::Abs(y[2]) < fMaxY)fHistNegP[theBin]->Fill(logdiff[2]);
1145 if(fMC){
1146 if((TMath::Abs(y[0])<fMaxY) && (TMath::Abs(code)==211)) fHistMCNegPi[theBin]->Fill(logdiff[0]);
1147 if((TMath::Abs(y[1])<fMaxY) && (TMath::Abs(code)==321)) fHistMCNegK[theBin]->Fill(logdiff[1]);
1148 if((TMath::Abs(y[2])<fMaxY) && (TMath::Abs(code)==2212)) fHistMCNegP[theBin]->Fill(logdiff[2]);
1149 }
1150 }
1151 }
1152
1153
1154 //fill propaganda plot with dedx
1155 fHistDEDX->Fill(track->GetP(),dedx);
1156 fHistDEDXdouble->Fill(track->GetP()*track->GetSign(),dedx);
1157
1158 //fill charge distribution histo to check the calibration
1159 for(Int_t j=0;j<4;j++){
1160 if(s[j]<5) continue;
1161 fHistCharge[j]->Fill(s[j]);
1162 }
1163 }
1164
1165 // Post output data.
1166 PostData(1,fOutput);
1167 Printf("............. end of Exec");
1168}
1169
1170//________________________________________________________________________
1171void AliAnalysisTaskSEITSsaSpectra::Terminate(Option_t *) {
4d668cf7 1172 // Merge output
1173 // Called once at the end of the query
36be14b3 1174
1175 fOutput = dynamic_cast<TList*>(GetOutputData(1));
1176 if (!fOutput) {
1177 printf("ERROR: fOutput not available\n");
1178 return;
1179 }
1180 fHistNEvents = dynamic_cast<TH1F*>(fOutput->FindObject("fHistNEvents"));
c841b89f 1181 fHistMult = dynamic_cast<TH1F*>(fOutput->FindObject("fHistMult"));
1f75ac7d 1182 fHistCen = dynamic_cast<TH1F*>(fOutput->FindObject("fHistCen"));
36be14b3 1183 fHistNTracks = dynamic_cast<TH1F*>(fOutput->FindObject("fHistNTracks"));
c841b89f 1184 fHistNTracksPos = dynamic_cast<TH1F*>(fOutput->FindObject("fHistNTracksPos"));
1185 fHistNTracksNeg = dynamic_cast<TH1F*>(fOutput->FindObject("fHistNTracksNeg"));
36be14b3 1186 fHistDEDX = dynamic_cast<TH2F*>(fOutput->FindObject("fHistDEDX"));
1187 fHistDEDXdouble = dynamic_cast<TH2F*>(fOutput->FindObject("fHistDEDXdouble"));
1188
1189 fHistBeforeEvSel = dynamic_cast<TH1F*>(fOutput->FindObject("fHistBeforeEvSel"));
1190 fHistAfterEvSel = dynamic_cast<TH1F*>(fOutput->FindObject("fHistAfterEvSel"));
1191
1192
1193 for(Int_t j=0;j<3;j++){
10ba520e 1194 if(fMC){
1195 fHistPrimMCpos[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistPrimMCpos%d",j)));
1196 fHistPrimMCneg[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistPrimMCneg%d",j)));
1197 fHistSecStrMCpos[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistSecStrMCpos%d",j)));
1198 fHistSecStrMCneg[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistSecStrMCneg%d",j)));
1199 fHistSecMatMCpos[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistSecMatMCpos%d",j)));
1200 fHistSecMatMCneg[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistSecMatMCneg%d",j)));
1201 //
1202 fHistPrimMCposBefEvSel[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistPrimMCposBefEvSel%d",j)));
1203 fHistPrimMCnegBefEvSel[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistPrimMCnegBefEvSel%d",j)));
1204 fHistSecStrMCposBefEvSel[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistSecStrMCposBefEvSel%d",j)));
1205 fHistSecStrMCnegBefEvSel[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistSecStrMCnegBefEvSel%d",j)));
1206 fHistSecMatMCposBefEvSel[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistSecMatMCposBefEvSel%d",j)));
1207 fHistSecMatMCnegBefEvSel[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistSecMatMCnegBefEvSel%d",j)));
1208 //
1209 fHistPrimMCposReco[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistPrimMCposReco%d",j)));
1210 fHistPrimMCnegReco[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistPrimMCnegReco%d",j)));
1211 fHistSecStrMCposReco[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistSecStrMCposReco%d",j)));
1212 fHistSecStrMCnegReco[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistSecStrMCnegReco%d",j)));
1213 fHistSecMatMCposReco[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistSecMatMCposReco%d",j)));
1214 fHistSecMatMCnegReco[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistSecMatMCnegReco%d",j)));
1215 }
36be14b3 1216 }
36be14b3 1217
36be14b3 1218 for(Int_t i=0; i<4; i++){
1219 fHistCharge[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistChargeLay%d",i)));
1220 }
1221
1222 for(Int_t i=0; i<kNbins; i++){
1223 fHistPosPi[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistPosPi%d",i)));
1224 fHistPosK[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistPosK%d",i)));
1225 fHistPosP[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistPosP%d",i)));
1226 fHistNegPi[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistNegPi%d",i)));
1227 fHistNegK[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistNegK%d",i)));
1228 fHistNegP[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistNegP%d",i)));
1229
1230 fHistDCAPosPi[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistDCAPosPi%d",i)));
1231 fHistDCAPosK[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistDCAPosK%d",i)));
1232 fHistDCAPosP[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistDCAPosP%d",i)));
1233 fHistDCANegPi[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistDCANegPi%d",i)));
1234 fHistDCANegK[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistDCANegK%d",i)));
1235 fHistDCANegP[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistDCANegP%d",i)));
1236
1237
1238 if(fMC){
10ba520e 1239
1240 fHistMCPrimDCAPosPi[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCPrimDCAPosPi%d",i)));
1241 fHistMCPrimDCAPosK[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCPrimDCAPosK%d",i)));
1242 fHistMCPrimDCAPosP[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCPrimDCAPosP%d",i)));
1243 fHistMCPrimDCANegPi[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCPrimDCANegPi%d",i)));
1244 fHistMCPrimDCANegK[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCPrimDCANegK%d",i)));
1245 fHistMCPrimDCANegP[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCPrimDCANegP%d",i)));
1246
1247 fHistMCSecStDCAPosPi[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCSecStDCAPosPi%d",i)));
1248 fHistMCSecStDCAPosK[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCSecStDCAPosK%d",i)));
1249 fHistMCSecStDCAPosP[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCSecStDCAPosP%d",i)));
1250 fHistMCSecStDCANegPi[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCSecStDCANegPi%d",i)));
1251 fHistMCSecStDCANegK[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCSecStDCANegK%d",i)));
1252 fHistMCSecStDCANegP[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCSecStDCANegP%d",i)));
1253
1254 fHistMCSecMatDCAPosPi[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCSecMatDCAPosPi%d",i)));
1255 fHistMCSecMatDCAPosK[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCSecMatDCAPosK%d",i)));
1256 fHistMCSecMatDCAPosP[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCSecMatDCAPosP%d",i)));
1257 fHistMCSecMatDCANegPi[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCSecMatDCANegPi%d",i)));
1258 fHistMCSecMatDCANegK[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCSecMatDCANegK%d",i)));
1259 fHistMCSecMatDCANegP[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCSecMatDCANegP%d",i)));
1260
36be14b3 1261 fHistMCPosPi[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCPosPi%d",i)));
1262 fHistMCPosK[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCPosK%d",i)));
1263 fHistMCPosP[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCPosP%d",i)));
1264 fHistMCNegPi[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCNegPi%d",i)));
1265 fHistMCNegK[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCNegK%d",i)));
1266 fHistMCNegP[i] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("fHistMCNegP%d",i)));
1267 }
1268 }
10ba520e 1269
36be14b3 1270 for(Int_t j=0;j<3;j++){
10ba520e 1271 fHistPosNSigmaMean[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("hHistPosNSigmaMean%d",j)));
1272 fHistNegNSigmaMean[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("hHistNegNSigmaMean%d",j)));
1273 fHistPosNSigmaMCMean[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("hHistPosNSigmaMCMean%d",j)));
1274 fHistNegNSigmaMCMean[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("hHistNegNSigmaMCMean%d",j)));
1275 fHistPosNSigmaPrimMean[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("hHistPosNSigmaPrimMean%d",j)));
1276 fHistNegNSigmaPrimMean[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("hHistNegNSigmaPrimMean%d",j)));
1277 fHistPosNSigmaPrimMCMean[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("hHistPosNSigmaPrimMCMean%d",j)));
1278 fHistNegNSigmaPrimMCMean[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("hHistNegNSigmaPrimMCMean%d",j)));
1279
36be14b3 1280 fHistPosNSigma[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("hHistPosNSigma%d",j)));
1281 fHistNegNSigma[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("hHistNegNSigma%d",j)));
10ba520e 1282 fHistPosNSigmaMC[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("hHistPosNSigmaMC%d",j)));
1283 fHistNegNSigmaMC[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("hHistNegNSigmaMC%d",j)));
36be14b3 1284 fHistPosNSigmaPrim[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("hHistPosNSigmaPrim%d",j)));
1285 fHistNegNSigmaPrim[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("hHistNegNSigmaPrim%d",j)));
1286 fHistPosNSigmaPrimMC[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("hHistPosNSigmaPrimMC%d",j)));
1287 fHistNegNSigmaPrimMC[j] = dynamic_cast<TH1F*>(fOutput->FindObject(Form("hHistNegNSigmaPrimMC%d",j)));
10ba520e 1288
36be14b3 1289 }
1290
1291 fNtupleNSigma = dynamic_cast<TNtuple*>(fOutput->FindObject("fNtupleNSigma"));
1292 fNtupleMC = dynamic_cast<TNtuple*>(fOutput->FindObject("fNtupleMC"));
1293
1294 Printf("end of Terminate");
1295 return;
1296}