]>
Commit | Line | Data |
---|---|---|
84bb3fcb | 1 | /************************************************************************** |
2 | * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * | |
3 | * * | |
4 | * Author: Satyajit Jena. * | |
5 | * Contributors are mentioned in the code where appropriate. * | |
6 | * * | |
7 | * Permission to use, copy, modify and distribute this software and its * | |
8 | * documentation strictly for non-commercial purposes is hereby granted * | |
9 | * without fee, provided that the above copyright notice appears in all * | |
10 | * copies and that both the copyright notice and this permission notice * | |
11 | * appear in the supporting documentation. The authors make no claims * | |
12 | * about the suitability of this software for any purpose. It is * | |
13 | * provided "as is" without express or implied warranty. * | |
14 | **************************************************************************/ | |
15 | ||
16 | ||
17 | //=========================================================================// | |
18 | // // | |
19 | // Analysis Task for Net-Charge Higher Moment Analysis // | |
20 | // Author: Satyajit Jena || Nirbhay K. Behera // | |
21 | // sjena@cern.ch || nbehera@cern.ch // | |
22 | // // | |
23 | //=========================================================================// | |
24 | #include "TChain.h" | |
25 | #include "TList.h" | |
26 | #include "TFile.h" | |
27 | #include "TTree.h" | |
28 | #include "TH1D.h" | |
29 | #include "TH2D.h" | |
30 | #include "TH3D.h" | |
31 | #include "THnSparse.h" | |
32 | #include "TCanvas.h" | |
33 | ||
34 | #include "AliAnalysisTaskSE.h" | |
35 | #include "AliAnalysisManager.h" | |
36 | #include "AliAODEvent.h" | |
37 | #include "AliVTrack.h" | |
38 | #include "AliAODTrack.h" | |
39 | #include "AliAODInputHandler.h" | |
40 | #include "AliAODEvent.h" | |
41 | #include "AliStack.h" | |
42 | #include "AliAODMCHeader.h" | |
43 | #include "AliAODMCParticle.h" | |
44 | #include "AliMCEventHandler.h" | |
45 | #include "AliMCEvent.h" | |
46 | #include "AliStack.h" | |
47 | #include "AliGenHijingEventHeader.h" | |
48 | #include "AliGenEventHeader.h" | |
49 | #include "AliPID.h" | |
50 | #include "AliAODPid.h" | |
51 | #include "AliPIDResponse.h" | |
52 | #include "AliAODpidUtil.h" | |
53 | #include "AliPIDCombined.h" | |
54 | #include "AliHelperPID.h" | |
55 | ||
56 | #include "AliEbyEHigherMomentsTaskPID.h" | |
57 | ||
58 | using std::cout; | |
59 | using std::endl; | |
60 | ||
61 | ClassImp(AliEbyEHigherMomentsTaskPID) | |
62 | ||
63 | ||
64 | //----------------------------------------------------------------------- | |
65 | AliEbyEHigherMomentsTaskPID::AliEbyEHigherMomentsTaskPID( const char *name ) | |
66 | : AliAnalysisTaskSE( name ), | |
67 | fListOfHistos(0), | |
68 | fArrayMC(0), | |
69 | fPIDResponse(0), | |
70 | fParticleSpecies(AliPID::kProton), | |
71 | fAnalysisType(0), | |
72 | fCentralityEstimator("V0M"), | |
73 | fCentrality(0), | |
74 | fVxMax(3.), | |
75 | fVyMax(3.), | |
76 | fVzMax(10.), | |
77 | fPtLowerLimit(0.3), | |
78 | fPtHigherLimit(1.5), | |
79 | fNptBins(0), | |
80 | fBin(0), | |
81 | fEtaLowerLimit(-0.8), | |
82 | fEtaHigherLimit(0.8), | |
83 | fRapidityCut(0.5), | |
84 | fNSigmaCut(3.), | |
85 | fAODtrackCutBit(128), | |
86 | fHelperPID(0), | |
87 | fEventCounter(0), | |
88 | fHistVxVy(0), | |
89 | fTHnCentNplusNminusPid(0), | |
90 | fTHnCentNplusNminusPidTruth(0), | |
91 | fPtBinNplusNminusPid(0), | |
92 | fPtBinNplusNminusPidTruth(0) | |
93 | { | |
94 | ||
95 | for ( Int_t i = 0; i < 4; i++) { | |
96 | fHistQA[i] = NULL; | |
97 | } | |
98 | ||
99 | DefineOutput(1, TList::Class()); // Outpput... | |
100 | } | |
101 | ||
102 | AliEbyEHigherMomentsTaskPID::~AliEbyEHigherMomentsTaskPID(){ | |
103 | ||
104 | if(fListOfHistos) delete fListOfHistos; | |
105 | if(fHelperPID) delete fHelperPID; | |
106 | } | |
107 | ||
108 | //--------------------------------------------------------------------------------- | |
109 | void AliEbyEHigherMomentsTaskPID::UserCreateOutputObjects() { | |
110 | ||
111 | fListOfHistos = new TList(); | |
112 | fListOfHistos->SetOwner(kTRUE); | |
113 | ||
114 | fEventCounter = new TH1D("fEventCounter","EventCounter", 10, 0.5,10.5); | |
115 | fEventCounter->GetXaxis()->SetBinLabel(1,"Event Accesed"); | |
116 | fEventCounter->GetXaxis()->SetBinLabel(2,"Within 0-80% centrality"); | |
117 | fEventCounter->GetXaxis()->SetBinLabel(5,"Have a vertex"); | |
118 | fEventCounter->GetXaxis()->SetBinLabel(6,"After vertex Cut"); | |
119 | fEventCounter->GetXaxis()->SetBinLabel(7,"Event Analyzed"); | |
120 | fListOfHistos->Add(fEventCounter); | |
121 | ||
122 | //For QA-Histograms | |
123 | fHistQA[0] = new TH1D("fHistQAvz", "Histo Vz After Cut", 500, -25., 25.); | |
124 | fHistQA[1] = new TH1D("fHistQAPt","p_{T} distribution",600,0,6); | |
125 | fHistQA[2] = new TH1D("fHistQAEta","#eta distribution",240,-1.2,1.2); | |
126 | fHistQA[3] = new TH1D("fHistQAPhi","#phi distribution",340,0,6.8); | |
127 | ||
128 | for(Int_t i = 0; i < 4; i++) | |
129 | { | |
130 | fListOfHistos->Add(fHistQA[i]); | |
131 | } | |
132 | ||
133 | fHistVxVy = new TH2D("fHistVxVy","Vertex-x Vs Vertex-y", 200, -1., 1., 200, -1., 1.); | |
134 | fListOfHistos->Add(fHistVxVy); | |
135 | ||
136 | ||
137 | const Int_t nDim = 3; | |
138 | const Int_t nPid = 5; | |
139 | TString Species[nPid] = {"Electron","Muon","Pion","Kaon","Proton"}; | |
140 | ||
141 | Int_t fBins[nPid][nDim]; | |
142 | Double_t fMin[nPid][nDim]; | |
143 | Double_t fMax[nPid][nDim]; | |
144 | ||
145 | for( Int_t i = 0; i < nPid; i++ ){ | |
146 | for( Int_t j = 0; j < nDim; j++ ){ | |
147 | fBins[i][j] = 0; | |
148 | fMin[i][j] = 0.; | |
149 | fMax[i][j] = 0.; | |
150 | } | |
151 | } | |
152 | ||
153 | ||
154 | const Int_t dim = fNptBins*2 + 1; | |
155 | Int_t bin[nPid][dim]; | |
156 | Double_t min[nPid][dim]; | |
157 | Double_t max[nPid][dim]; | |
158 | ||
159 | bin[2][0] = 100; min[2][0] = -0.5; max[2][0] = 99.5; | |
160 | bin[3][0] = 100; min[3][0] = -0.5; max[3][0] = 99.5; | |
161 | bin[4][0] = 100; min[4][0] = -0.5; max[4][0] = 99.5; | |
162 | ||
163 | for(Int_t i = 1; i < dim; i++){ | |
164 | ||
165 | bin[2][i] = 900; | |
166 | min[2][i] = -0.5; | |
167 | max[2][i] = 899.5; | |
168 | ||
169 | bin[3][i] = 700; | |
170 | min[3][i] = -0.5; | |
171 | max[3][i] = 699.5; | |
172 | ||
173 | bin[4][i] = 400; | |
174 | min[4][i] = -0.5; | |
175 | max[4][i] = 399.5; | |
176 | ||
177 | ||
178 | } | |
179 | ||
180 | ||
181 | ||
182 | ||
183 | ||
184 | TString hname1, hname11; | |
185 | TString htitle1, axisTitle1, axisTitle2; | |
186 | ||
187 | Int_t gPid = (Int_t)fParticleSpecies; | |
188 | ||
189 | //Pion---- | |
190 | fBins[2][0] = 100; fBins[2][1] = 1000; fBins[2][2] = 1000; | |
191 | fMin[2][0] = -0.5; fMin[2][1] = -0.5; fMin[2][2] = -0.5; | |
192 | fMax[2][0] = 99.5; fMax[2][1] = 999.5; fMax[2][2] = 999.5; | |
193 | //Kaon------ | |
194 | fBins[3][0] = 100; fBins[3][1] = 700; fBins[3][2] = 700; | |
195 | fMin[3][0] = -0.5; fMin[3][1] = -0.5; fMin[3][2] = -0.5; | |
196 | fMax[3][0] = 99.5; fMax[3][1] = 699.5; fMax[3][2] = 699.5; | |
197 | //Proton----- | |
198 | fBins[4][0] = 100; fBins[4][1] = 400; fBins[4][2] = 400; | |
199 | fMin[4][0] = -0.5; fMin[4][1] = -0.5; fMin[4][2] = -0.5; | |
200 | fMax[4][0] = 99.5; fMax[4][1] = 399.5; fMax[4][2] = 399.5; | |
201 | ||
202 | hname1 = "fCentNplusNminusPid"; hname1 +=gPid; | |
203 | htitle1 = Species[gPid]; htitle1 +=" And Neg-"; htitle1 +=Species[gPid]; | |
204 | axisTitle1 ="Pos-"; axisTitle1 += Species[gPid]; | |
205 | axisTitle2 ="Neg-"; axisTitle2 += Species[gPid]; | |
206 | ||
207 | fTHnCentNplusNminusPid = new THnSparseD(hname1.Data(),htitle1.Data(), nDim, fBins[gPid], fMin[gPid],fMax[gPid]); | |
208 | ||
209 | fTHnCentNplusNminusPid->GetAxis(0)->SetTitle("Centrality"); | |
210 | fTHnCentNplusNminusPid->GetAxis(1)->SetTitle(axisTitle1.Data()); | |
211 | fTHnCentNplusNminusPid->GetAxis(2)->SetTitle(axisTitle2.Data()); | |
212 | fListOfHistos->Add(fTHnCentNplusNminusPid); | |
213 | ||
214 | TString hname2 = "fPtBinNplusNminusPid"; | |
215 | TString htitle2 = "cent-nplus-nminus-ptbinwise"; | |
216 | hname2 += gPid; | |
217 | htitle2 += gPid; | |
218 | fPtBinNplusNminusPid = new THnSparseI(hname2.Data(),htitle2.Data(), dim, bin[gPid], min[gPid], max[gPid]); | |
219 | fListOfHistos->Add(fPtBinNplusNminusPid); | |
220 | ||
221 | ||
222 | ||
223 | ||
224 | if( fAnalysisType == "MCAOD" ){ | |
225 | ||
226 | hname11 = "fCentNplusNminusPidTruth"; hname11 +=gPid; | |
227 | fTHnCentNplusNminusPidTruth = new THnSparseD(hname11.Data(),htitle1.Data(), nDim, fBins[gPid], fMin[gPid],fMax[gPid]); | |
228 | ||
229 | fTHnCentNplusNminusPidTruth->GetAxis(0)->SetTitle("Centrality"); | |
230 | fTHnCentNplusNminusPidTruth->GetAxis(1)->SetTitle(axisTitle1.Data()); | |
231 | fTHnCentNplusNminusPidTruth->GetAxis(2)->SetTitle(axisTitle2.Data()); | |
232 | fListOfHistos->Add(fTHnCentNplusNminusPidTruth); | |
233 | ||
234 | TString hname22 = "fPtBinNplusNminusPidTruth"; | |
235 | TString htitle22 = "cent-nplus-nminus-ptbinwise"; | |
236 | hname22 += gPid; | |
237 | htitle22 += gPid; | |
238 | ||
239 | fPtBinNplusNminusPidTruth = new THnSparseI(hname22.Data(),htitle22.Data(), dim, bin[gPid], min[gPid], max[gPid]); | |
240 | fListOfHistos->Add(fPtBinNplusNminusPidTruth); | |
241 | ||
242 | }//fUsePid------- | |
243 | ||
244 | PostData(1, fListOfHistos); | |
245 | ||
246 | ||
247 | } | |
248 | ||
249 | ||
250 | //---------------------------------------------------------------------------------- | |
251 | void AliEbyEHigherMomentsTaskPID::UserExec( Option_t * ){ | |
252 | ||
253 | fEventCounter->Fill(1); | |
254 | ||
255 | if(fAnalysisType == "AOD") { | |
256 | ||
257 | doAODEvent(); | |
258 | ||
259 | }//AOD--analysis----- | |
260 | ||
261 | else if(fAnalysisType == "MCAOD") { | |
262 | ||
263 | doMCAODEvent(); | |
264 | ||
265 | } | |
266 | ||
267 | else return; | |
268 | ||
269 | ||
270 | ||
271 | } | |
272 | ||
273 | //-------------------------------------------------------------------------------------- | |
274 | void AliEbyEHigherMomentsTaskPID::doAODEvent(){ | |
275 | ||
276 | Double_t posPidSum = 0.; | |
277 | Double_t negPidSum = 0.; | |
278 | ||
279 | const Int_t dim = fNptBins*2; | |
280 | Double_t PtCh[dim]; | |
281 | ||
282 | for(Int_t idx = 0; idx < dim; idx++){ | |
283 | PtCh[idx] = 0.; | |
284 | } | |
285 | ||
286 | ||
287 | Int_t gPid = (Int_t)fParticleSpecies; | |
288 | ||
289 | AliAnalysisManager* manager = AliAnalysisManager::GetAnalysisManager(); | |
290 | if (!manager) { | |
291 | cout<<"ERROR: Analysis manager not found."<<endl; | |
292 | return; | |
293 | } | |
294 | //coneect to the inputHandler------------ | |
295 | AliAODInputHandler* inputHandler = dynamic_cast<AliAODInputHandler*> (manager->GetInputEventHandler()); | |
296 | if (!inputHandler) { | |
297 | cout<<"ERROR: Input handler not found."<<endl; | |
298 | return; | |
299 | } | |
300 | ||
301 | AliAODEvent *fAOD = dynamic_cast<AliAODEvent*>(InputEvent()); | |
302 | ||
303 | if (!fAOD) | |
304 | { | |
305 | cout<< "ERROR: AOD not found " <<endl; | |
306 | return; | |
307 | } | |
308 | ||
309 | fPIDResponse = inputHandler->GetPIDResponse(); | |
310 | ||
311 | ||
312 | if (!fPIDResponse){ | |
313 | AliFatal("This Task needs the PID response attached to the inputHandler"); | |
314 | return; | |
315 | } | |
316 | ||
317 | if(!ProperVertex(fAOD)) return; | |
318 | ||
319 | AliAODHeader *aodHeader = fAOD->GetHeader(); | |
320 | ||
321 | fCentrality = (Int_t)aodHeader->GetCentralityP()->GetCentralityPercentile(fCentralityEstimator.Data()); | |
322 | ||
323 | if(fCentrality < 0 || fCentrality >= 81) return; | |
324 | ||
325 | fEventCounter->Fill(2); | |
326 | ||
327 | Int_t nTracks = fAOD->GetNumberOfTracks(); | |
328 | ||
329 | for(Int_t i = 0; i < nTracks; i++) { | |
330 | ||
331 | AliAODTrack* aodTrack = dynamic_cast<AliAODTrack *>(fAOD->GetTrack(i)); | |
332 | ||
333 | if(!aodTrack) { | |
334 | AliError(Form("ERROR: Could not retrieve AODtrack %d",i)); | |
335 | continue; | |
336 | } | |
337 | ||
338 | if(!AcceptTrack(aodTrack) ) continue; | |
339 | ||
340 | ||
341 | fBin = GetPtBin(aodTrack->Pt()); | |
342 | ||
343 | Short_t gCharge = aodTrack->Charge(); | |
344 | ||
345 | ||
346 | Double_t rap = aodTrack->Y(AliPID::ParticleMass(fParticleSpecies)); | |
347 | if( fabs(rap ) > fRapidityCut ) continue;//Rapidity cut | |
348 | ||
349 | Int_t PID = fHelperPID->GetParticleSpecies((AliVTrack*)aodTrack, kFALSE); | |
350 | ||
351 | if( (PID+2) == gPid ){ | |
352 | ||
353 | if(gCharge > 0){ | |
354 | posPidSum += 1.; | |
355 | PtCh[fBin] += 1; | |
356 | } | |
357 | if(gCharge < 0){ | |
358 | negPidSum += 1.; | |
359 | PtCh[fNptBins+fBin] += 1; | |
360 | } | |
361 | ||
362 | }//chek the PID | |
363 | ||
364 | ||
365 | }//--------- Track Loop to select with filterbit | |
366 | ||
367 | ||
2942f542 | 368 | Double_t fContainerPid[3] = { static_cast<Double_t>(fCentrality), posPidSum, negPidSum}; |
84bb3fcb | 369 | fTHnCentNplusNminusPid->Fill(fContainerPid); |
370 | ||
371 | //cout << fCentrality <<" "<< posPidSum <<" " << negPidSum << endl; | |
372 | ||
373 | Double_t ptContainer[dim+1]; | |
374 | ||
375 | ptContainer[0] = fCentrality; | |
376 | ||
377 | for(Int_t i = 1; i <= dim; i++){ | |
378 | ptContainer[i] = PtCh[i-1]; | |
379 | //cout << PtCh[i-1] <<" ,"; | |
380 | } | |
381 | //cout << endl; | |
382 | ||
383 | fPtBinNplusNminusPid->Fill(ptContainer); | |
384 | ||
385 | ||
386 | fEventCounter->Fill(7); | |
387 | PostData(1, fListOfHistos); | |
388 | return; | |
389 | ||
390 | } | |
391 | //-------------------------------------------------------------------------------------- | |
392 | //-------------------------------------------------------------------------------------- | |
393 | void AliEbyEHigherMomentsTaskPID::doMCAODEvent(){ | |
394 | ||
395 | ||
396 | //--------- | |
397 | Double_t posPidSumMCRec = 0.; | |
398 | Double_t negPidSumMCRec = 0.; | |
399 | ||
400 | Double_t posPidSumMCTruth = 0.; | |
401 | Double_t negPidSumMCTruth = 0.; | |
402 | ||
403 | const Int_t dim = fNptBins*2; | |
404 | Double_t PtCh[dim]; | |
405 | Double_t ptChMC[dim]; | |
406 | ||
407 | for(Int_t idx = 0; idx < dim; idx++){ | |
408 | PtCh[idx] = 0; | |
409 | ptChMC[idx] = 0; | |
410 | } | |
411 | ||
412 | ||
413 | Int_t gPid = (Int_t)fParticleSpecies; | |
414 | ||
415 | Int_t gPdgCode = AliPID::ParticleCode(fParticleSpecies); | |
416 | ||
417 | //Connect to Anlaysis manager------ | |
418 | AliAnalysisManager* manager = AliAnalysisManager::GetAnalysisManager(); | |
419 | if (!manager) { | |
420 | cout<<"ERROR: Analysis manager not found."<<endl; | |
421 | return; | |
422 | } | |
423 | ||
424 | AliAODInputHandler* inputHandler = dynamic_cast<AliAODInputHandler*> (manager->GetInputEventHandler()); | |
425 | if (!inputHandler) { | |
426 | cout<<"ERROR: Input handler not found."<<endl; | |
427 | return; | |
428 | } | |
429 | ||
430 | //AOD event------ | |
431 | AliAODEvent *fAOD = dynamic_cast<AliAODEvent*>(InputEvent()); | |
432 | ||
433 | if (!fAOD) | |
434 | { | |
435 | cout<< "ERROR: AOD not found " <<endl; | |
436 | return; | |
437 | } | |
438 | ||
439 | fPIDResponse =inputHandler->GetPIDResponse(); | |
440 | ||
441 | ||
442 | if(!fPIDResponse){ | |
443 | AliFatal("This Task needs the PID response attached to the inputHandler"); | |
444 | return; | |
445 | } | |
446 | ||
447 | // -- Get MC header | |
448 | // ------------------------------------------------------------------ | |
449 | ||
450 | fArrayMC = dynamic_cast<TClonesArray*>(fAOD->FindListObject(AliAODMCParticle::StdBranchName())); | |
451 | if (!fArrayMC) { | |
452 | AliFatal("Error: MC particles branch not found!\n"); | |
453 | return; | |
454 | } | |
455 | ||
456 | AliAODMCHeader *mcHdr=NULL; | |
457 | mcHdr=(AliAODMCHeader*)fAOD->GetList()->FindObject(AliAODMCHeader::StdBranchName()); | |
458 | if(!mcHdr) { | |
459 | printf("MC header branch not found!\n"); | |
460 | return; | |
461 | } | |
462 | ||
463 | if(!ProperVertex(fAOD)) return; | |
464 | ||
465 | AliAODHeader *aodHeader = fAOD->GetHeader(); | |
466 | ||
467 | fCentrality = (Int_t)aodHeader->GetCentralityP()->GetCentralityPercentile(fCentralityEstimator.Data()); | |
468 | ||
469 | if( fCentrality < 0 || fCentrality >= 81) return; | |
470 | ||
471 | fEventCounter->Fill(2); | |
472 | ||
473 | Int_t nTracks = fAOD->GetNumberOfTracks(); | |
474 | ||
475 | for(Int_t i = 0; i < nTracks; i++) { | |
476 | ||
477 | AliAODTrack* aodTrack = dynamic_cast<AliAODTrack *>(fAOD->GetTrack(i)); | |
478 | ||
479 | if(!aodTrack) { | |
480 | AliError(Form("ERROR: Could not retrieve AODtrack %d",i)); | |
481 | continue; | |
482 | } | |
483 | ||
484 | if(!AcceptTrack(aodTrack) ) continue; | |
485 | ||
486 | fBin = GetPtBin(aodTrack->Pt()); | |
487 | ||
488 | //cout << "Pt Bin " << fBin << endl; | |
489 | ||
490 | Short_t gCharge = aodTrack->Charge(); | |
491 | ||
492 | Double_t rap = aodTrack->Y(AliPID::ParticleMass(fParticleSpecies)); | |
493 | if( fabs(rap ) > fRapidityCut ) continue;//Rapidity cut | |
494 | ||
495 | Int_t PID = fHelperPID->GetParticleSpecies((AliVTrack*)aodTrack, kFALSE); | |
496 | ||
497 | if( (PID+2) == gPid ){ | |
498 | if(gCharge > 0){ | |
499 | posPidSumMCRec += 1; | |
500 | PtCh[fBin] += 1; | |
501 | } | |
502 | if( gCharge < 0 ){ | |
503 | negPidSumMCRec += 1.; | |
504 | PtCh[fNptBins+fBin] += 1; | |
505 | } | |
506 | }//nSigmaCut----- | |
507 | ||
508 | }//--------- Track Loop to select with filterbit | |
509 | ||
510 | //=============================================================================== | |
511 | //---------------------MC Truth-------------------------------------------------- | |
512 | //=============================================================================== | |
513 | ||
514 | Int_t nMCTrack = fArrayMC->GetEntriesFast(); | |
515 | ||
516 | for (Int_t iMC = 0; iMC < nMCTrack; iMC++) { | |
517 | ||
518 | AliAODMCParticle *partMC = (AliAODMCParticle*) fArrayMC->At(iMC); | |
519 | ||
520 | if(!partMC){ | |
521 | AliError(Form("ERROR: Could not retrieve AODMCtrack %d",iMC)); | |
522 | continue; | |
523 | } | |
524 | ||
525 | if(partMC->Charge() == 0) continue; | |
526 | if(!partMC->IsPhysicalPrimary()) continue; | |
527 | ||
528 | if (partMC->Eta() < fEtaLowerLimit || partMC->Eta() > fEtaHigherLimit) continue; | |
529 | if (partMC->Pt() < fPtLowerLimit || partMC->Pt() > fPtHigherLimit) continue; | |
530 | ||
531 | Short_t gCharge = partMC->Charge(); | |
532 | Int_t mcbin = GetPtBin(partMC->Pt()); | |
533 | ||
534 | Double_t rap = partMC->Y(); | |
535 | if( fabs(rap ) > fRapidityCut ) continue;//Rapidity cut | |
536 | if(TMath::Abs(partMC->GetPdgCode()) != gPdgCode) continue; | |
537 | ||
538 | if(gCharge > 0){ | |
539 | posPidSumMCTruth += 1.; | |
540 | ptChMC[mcbin] += 1; | |
541 | } | |
542 | if(gCharge < 0){ | |
543 | negPidSumMCTruth += 1.; | |
544 | ptChMC[fNptBins+mcbin] += 1; | |
545 | } | |
546 | ||
547 | ||
548 | }//MC-Truth Track loop-- | |
549 | ||
550 | ||
551 | //cout << fCentrality << " " << positiveSumMCRec << " " << negativeSumMCRec <<endl; | |
552 | //cout <<" " << positiveSumMCTruth << " " << negativeSumMCTruth << endl; | |
553 | //cout <<fCentrality<<" " << posPidSumMCRec << " " << negPidSumMCRec << endl; | |
554 | //cout <<fCentrality<<" " << posPidSumMCTruth << " " << negPidSumMCTruth << endl; | |
555 | //cout <<"---------------------------------" << endl; | |
556 | ||
2942f542 | 557 | Double_t fContainerPid[3] = { static_cast<Double_t>(fCentrality), posPidSumMCRec, negPidSumMCRec};//Reco. |
558 | Double_t fContainerPidTruth[3] = { static_cast<Double_t>(fCentrality), posPidSumMCTruth, negPidSumMCTruth}; | |
84bb3fcb | 559 | |
560 | fTHnCentNplusNminusPid->Fill(fContainerPid);//Fill the rec. pid tracks | |
561 | fTHnCentNplusNminusPidTruth->Fill(fContainerPidTruth);//MC-Truth pid | |
562 | ||
563 | Double_t ptContainer[dim+1]; | |
564 | Double_t ptContainerMC[dim+1]; | |
565 | ptContainer[0] = fCentrality; | |
566 | ptContainerMC[0] = fCentrality; | |
567 | ||
568 | for(Int_t i = 1; i <= dim; i++){ | |
569 | ptContainer[i] = PtCh[i-1]; | |
570 | ptContainerMC[i] = ptChMC[i-1]; | |
571 | //cout <<" "<< PtCh[i-1]<<endl; | |
572 | //cout<< " Rec=" << PtCh[i-1]; | |
573 | } | |
574 | ||
575 | //cout << endl; | |
576 | ||
577 | fPtBinNplusNminusPid->Fill(ptContainer); | |
578 | fPtBinNplusNminusPidTruth->Fill(ptContainerMC); | |
579 | ||
580 | fEventCounter->Fill(7); | |
581 | PostData(1, fListOfHistos); | |
582 | return; | |
583 | ||
584 | } | |
585 | ||
586 | //--------------------------------------------------------------------------------------- | |
587 | Bool_t AliEbyEHigherMomentsTaskPID::ProperVertex(AliAODEvent *fAOD) const{ | |
588 | ||
589 | Bool_t IsVtx = kFALSE; | |
590 | ||
591 | const AliAODVertex *vertex = fAOD->GetPrimaryVertex(); | |
592 | ||
593 | if(vertex) { | |
594 | Double32_t fCov[6]; | |
595 | vertex->GetCovarianceMatrix(fCov); | |
596 | if(vertex->GetNContributors() > 0) { | |
597 | if(fCov[5] != 0) { | |
598 | ||
599 | Double_t lvx = vertex->GetX(); | |
600 | Double_t lvy = vertex->GetY(); | |
601 | Double_t lvz = vertex->GetZ(); | |
602 | ||
603 | fEventCounter->Fill(5); | |
604 | ||
605 | if(TMath::Abs(lvx) < fVxMax) { | |
606 | if(TMath::Abs(lvy) < fVyMax) { | |
607 | if(TMath::Abs(lvz) < fVzMax) { | |
608 | ||
609 | fEventCounter->Fill(6); | |
610 | fHistQA[0]->Fill(lvz); | |
611 | fHistVxVy->Fill(lvx,lvy); | |
612 | IsVtx = kTRUE; | |
613 | ||
614 | }//Z-Vertex cut--- | |
615 | }//Y-vertex cut-- | |
616 | }//X-vertex cut--- | |
617 | }//Covariance------ | |
618 | }//Contributors check--- | |
619 | }//If vertex----------- | |
620 | ||
621 | AliCentrality *centrality = fAOD->GetCentrality(); | |
622 | if (centrality->GetQuality() != 0) IsVtx = kFALSE; | |
623 | ||
624 | return IsVtx; | |
625 | } | |
626 | //--------------------------------------------------------------------------------------- | |
627 | ||
628 | //--------------------------------------------------------------------------------------- | |
629 | Bool_t AliEbyEHigherMomentsTaskPID::AcceptTrack(AliAODTrack* track) const{ | |
630 | ||
631 | if(!track) return kFALSE; | |
632 | if( track->Charge() == 0 ) return kFALSE; | |
633 | ||
634 | Double_t pt = track->Pt(); | |
635 | Double_t eta = track->Eta(); | |
636 | if(!track->TestFilterBit(fAODtrackCutBit)) return kFALSE; | |
637 | if( pt < fPtLowerLimit || pt > fPtHigherLimit ) return kFALSE; | |
638 | if( eta < fEtaLowerLimit || eta > fEtaHigherLimit ) return kFALSE; | |
639 | ||
640 | fHistQA[1]->Fill(pt); | |
641 | fHistQA[2]->Fill(eta); | |
642 | fHistQA[3]->Fill(track->Phi()); | |
643 | ||
644 | return kTRUE; | |
645 | } | |
646 | //------------------------------------------------------------------------ | |
647 | ||
648 | //------------------------------------------------------------------------ | |
649 | Int_t AliEbyEHigherMomentsTaskPID::GetPtBin(Double_t pt){ | |
650 | ||
651 | Int_t bin = 0; | |
652 | ||
653 | Double_t BinSize = (fPtHigherLimit - fPtLowerLimit)/fNptBins; | |
654 | ||
655 | for(Int_t iBin = 0; iBin < fNptBins; iBin++){ | |
656 | ||
657 | Double_t xLow = fPtLowerLimit + iBin*BinSize; | |
658 | Double_t xHigh = fPtLowerLimit + (iBin+1)*BinSize; | |
659 | ||
660 | if( pt >= xLow && pt < xHigh){ | |
661 | bin = iBin; | |
662 | //cout << "Pt Bin #" << bin <<"pt=" << pt << endl; | |
663 | break; | |
664 | } | |
665 | ||
666 | }//for | |
667 | ||
668 | ||
669 | return bin; | |
670 | } | |
671 | //------------------------------------------------------------------------ | |
672 | //------------------------------------------------------------------------ | |
673 | void AliEbyEHigherMomentsTaskPID::Terminate( Option_t * ){ | |
674 | ||
675 | Info("AliEbyEHigherMomentTask"," Task Successfully finished"); | |
676 | ||
677 | } | |
678 | ||
679 | //------------------------------------------------------------------------ |