]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGCF/Correlations/DPhi/AliDptDptInMC.cxx
Fixed dependencies
[u/mrichter/AliRoot.git] / PWGCF / Correlations / DPhi / AliDptDptInMC.cxx
CommitLineData
1486c432
MW
1/**************************************************************************
2 * Copyright(c) 1998-1999, 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/* $Id:$ */
17
18#include "TChain.h"
19#include "TList.h"
20#include "TFile.h"
21#include "TTree.h"
22#include "TH1D.h"
23#include "TH2D.h"
24#include "TH3D.h"
25#include "THnSparse.h"
26#include "TCanvas.h"
27
28#include <TROOT.h>
29#include <TChain.h>
30#include <TFile.h>
31#include <TList.h>
32#include <TMath.h>
33#include <TTree.h>
34#include <TH1F.h>
35#include <TH2F.h>
36#include <TH3F.h>
37#include <TProfile.h>
38#include <TH1D.h>
39#include <TH2D.h>
40#include <TH3D.h>
41#include <TRandom.h>
42#include "AliAnalysisManager.h"
43
44#include "AliAODHandler.h"
45#include "AliAODInputHandler.h"
46#include "AliInputEventHandler.h"
47#include "AliLog.h"
48#include "AliESDEvent.h"
49#include "AliESDInputHandler.h"
50#include "AliMultiplicity.h"
51#include "AliCentrality.h"
52#include "AliDptDptInMC.h"
53
54#include "AliPID.h"
947cb887 55#include "AliPIDResponse.h"
1486c432
MW
56
57#include "AliESDVertex.h"
58#include "AliESDEvent.h"
59#include "AliESDInputHandler.h"
60#include "AliAODEvent.h"
61#include "AliAODTrack.h"
62#include "AliAODInputHandler.h"
63#include "AliESD.h"
64#include "AliESDEvent.h"
65#include "AliAODEvent.h"
66#include "AliStack.h"
67#include "AliESDtrackCuts.h"
68#include "AliAODMCHeader.h"
69#include "AliAODMCParticle.h"
70#include "AliAODMCParticle.h"
71#include "AliAODMCHeader.h"
72#include "AliAODHeader.h"
73
74#include "AliGenHijingEventHeader.h"
75#include "AliGenEventHeader.h"
76#include "AliPID.h"
77#include "AliAODPid.h"
947cb887 78#include "AliPIDResponse.h"
1486c432
MW
79#include "AliAODpidUtil.h"
80#include "AliPIDCombined.h"
81
82using std::cout;
83using std::endl;
84
85ClassImp(AliDptDptInMC)
86
87AliDptDptInMC::AliDptDptInMC()
88: AliAnalysisTaskSE(),
89fAODEvent(0),
90fESDEvent(0), //! ESD Event
91fInputHandler(0),
947cb887 92fPIDResponse(0x0),
1486c432
MW
93_outputHistoList(0),
94fArrayMC (0),
95_twoPi ( 2.0 * 3.1415927),
96_eventCount ( 0),
97_debugLevel ( 0),
98_singlesOnly ( 0),
99_useWeights ( 0),
100_sameFilter ( false),
101_rejectPileup ( 1),
102_rejectPairConversion ( 0),
103_vertexZMin ( -10),
104_vertexZMax ( 10),
105_vertexXYMin ( -10),
106_vertexXYMax ( 10),
107_centralityMethod ( 4),
108_centralityMin ( 0.),
109_centralityMax ( 0.),
110_requestedCharge_1 ( 1),
111_requestedCharge_2 ( -1),
112_dcaZMin ( -3),
113_dcaZMax ( 3.),
114_dcaXYMin ( -2.4),
115_dcaXYMax ( 2.4),
116_dedxMin ( 0),
117_dedxMax ( 100000),
118_nClusterMin ( 80),
119_trackFilterBit (0),
120fNSigmaCut (3.),
729278ef 121fAnalysisType ("MCAOD"),
1486c432
MW
122fExcludeResonancesInMC (kFALSE),
123fExcludeElectronsInMC (kFALSE),
124_field ( 1.),
125_nTracks ( 0 ),
126_mult0 ( 0 ),
127_mult1 ( 0 ),
128_mult2 ( 0 ),
129_mult3 ( 0 ),
130_mult4 ( 0 ),
131_mult4a ( 0 ),
132_mult5 ( 0 ),
133_mult6 ( 0 ),
134arraySize ( 2000),
135_id_1(0),
136_charge_1(0),
137_iEtaPhi_1(0),
138_iPt_1(0),
139_pt_1(0),
140_px_1(0),
141_py_1(0),
142_pz_1(0),
143_correction_1(0),
144_dedx_1(0),
145_id_2(0),
146_charge_2(0),
147_iEtaPhi_2(0),
148_iPt_2(0),
149_pt_2(0),
150_px_2(0),
151_py_2(0),
152_pz_2(0),
153_correction_2(0),
154_dedx_2(0),
155_correctionWeight_1(0),
156_correctionWeight_2(0),
157_nBins_M0(500), _min_M0(0), _max_M0(10000), _width_M0(20),
158_nBins_M1(500), _min_M1(0), _max_M1(10000), _width_M1(20),
159_nBins_M2(500), _min_M2(0), _max_M2(10000), _width_M2(20),
160_nBins_M3(500), _min_M3(0), _max_M3(10000), _width_M3(20),
161_nBins_M4(100), _min_M4(0), _max_M4(1), _width_M4(0.01),
162_nBins_M5(100), _min_M5(0), _max_M5(1), _width_M5(0.01),
163_nBins_M6(100), _min_M6(0), _max_M6(1), _width_M6(0.01),
164_nBins_vertexZ(40), _min_vertexZ(-10), _max_vertexZ(10), _width_vertexZ(0.5),
165
166_nBins_pt_1(18), _min_pt_1(0.2), _max_pt_1(2.0), _width_pt_1(0.1),
167_nBins_phi_1(72), _min_phi_1(0), _max_phi_1(2.*3.1415927),_width_phi_1(2.*3.1415927/72.),
168//_nBins_eta_1(18),
169//_min_eta_1(-0.9), _max_eta_1(0.9), _width_eta_1(0.1),
170 _nBins_eta_1(0),
171 _min_eta_1(0), _max_eta_1(0), _width_eta_1(0.1),
172
173_nBins_etaPhi_1(0),
174_nBins_etaPhiPt_1(0),
175_nBins_zEtaPhiPt_1(0),
176
177_nBins_pt_2(18), _min_pt_2(0.2), _max_pt_2(2.0), _width_pt_2(0.1),
178_nBins_phi_2(72), _min_phi_2(0), _max_phi_2(2.*3.1415927),_width_phi_2(2.*3.1415927/72),
179//_nBins_eta_2(18),
180//_min_eta_2(-0.9), _max_eta_2(0.9), _width_eta_2(0.1),
181 _nBins_eta_2(0),
182 _min_eta_2(0), _max_eta_2(0), _width_eta_2(0.1),
183
184
185
186_nBins_etaPhi_2(0),
187_nBins_etaPhiPt_2(0),
188_nBins_zEtaPhiPt_2(0),
189_nBins_etaPhi_12(0),
190__n1_1(0),
191__n1_2(0),
192__n2_12(0),
193__s1pt_1(0),
194__s1pt_2(0),
195__s2ptpt_12(0),
196__s2NPt_12(0),
197__s2PtN_12(0),
198__n1Nw_1(0),
199__n1Nw_2(0),
200__n2Nw_12(0),
201__s1ptNw_1(0),
202__s1ptNw_2(0),
203__s2ptptNw_12(0),
204__s2NPtNw_12(0),
205__s2PtNNw_12(0),
206__n1_1_vsPt(0),
207__n1_1_vsEtaPhi(0),
208__s1pt_1_vsEtaPhi(0),
209__n1_1_vsZEtaPhiPt(0),
210__n1_2_vsPt(0),
211__n1_2_vsEtaPhi(0),
212__s1pt_2_vsEtaPhi(0),
213__n1_2_vsZEtaPhiPt(0),
214__n2_12_vsPtPt(0),
215__n2_12_vsEtaPhi(0),
216__s2ptpt_12_vsEtaPhi(0),
217__s2PtN_12_vsEtaPhi(0),
218__s2NPt_12_vsEtaPhi(0),
219_weight_1 ( 0 ),
220_weight_2 ( 0 ),
221_eventAccounting ( 0),
222_m0 ( 0),
223_m1 ( 0),
224_m2 ( 0),
225_m3 ( 0),
226_m4 ( 0),
227_m5 ( 0),
228_m6 ( 0),
229_vertexZ ( 0),
230
231_Ncluster1 ( 0),
232_Ncluster2 ( 0),
233_etadis ( 0),
234_phidis ( 0),
235_dcaz ( 0),
236_dcaxy ( 0),
0f50230b 237_etadis2 ( 0),
238_phidis2 ( 0),
239_dcaz2 ( 0),
240_dcaxy2 ( 0),
241_etadis3 ( 0),
242_phidis3 ( 0),
243_dcaz3 ( 0),
244_dcaxy3 ( 0),
1486c432
MW
245_n1_1_vsPt ( 0),
246_n1_1_vsEtaVsPhi ( 0),
247_s1pt_1_vsEtaVsPhi ( 0),
248_n1_1_vsZVsEtaVsPhiVsPt ( 0),
249_n1_1_vsM ( 0), // w/ weight
250_s1pt_1_vsM ( 0),
251_n1Nw_1_vsM ( 0), // w/o weight
252_s1ptNw_1_vsM ( 0),
253_dedxVsP_1 ( 0),
254_corrDedxVsP_1 ( 0),
255_betaVsP_1 ( 0),
256_n1_2_vsPt ( 0),
257_n1_2_vsEtaVsPhi ( 0),
258_s1pt_2_vsEtaVsPhi ( 0),
259_n1_2_vsZVsEtaVsPhiVsPt ( 0),
260_n1_2_vsM ( 0),
261_s1pt_2_vsM ( 0),
262_n1Nw_2_vsM ( 0),
263_s1ptNw_2_vsM ( 0),
264_dedxVsP_2 ( 0),
265_corrDedxVsP_2 ( 0),
266_betaVsP_2 ( 0),
267_n2_12_vsEtaPhi ( 0),
268_n2_12_vsPtVsPt ( 0),
269_s2PtPt_12_vsEtaPhi( 0),
270_s2PtN_12_vsEtaPhi ( 0),
271_s2NPt_12_vsEtaPhi ( 0),
272_n2_12_vsM ( 0),
273_s2PtPt_12_vsM ( 0),
274_s2PtN_12_vsM ( 0),
275_s2NPt_12_vsM ( 0),
276_n2Nw_12_vsM ( 0),
277_s2PtPtNw_12_vsM ( 0),
278_s2PtNNw_12_vsM ( 0),
279_s2NPtNw_12_vsM ( 0),
280_invMass ( 0),
281_invMassElec ( 0),
282n1Name("NA"),
283n1NwName("NA"),
284n2Name("NA"),
285n2NwName("NA"),
286n3Name("NA"),
287n1n1Name("NA"),
288n1n1n1Name("NA"),
289n2n1Name("NA"),
290r1Name("NA"),
291r2Name("NA"),
292r3Name("NA"),
293r2r1Name("NA"),
294c2Name("NA"),
295c3Name("NA"),
296d3Name("NA"),
297p3Name("NA"),
298cName("NA"),
299
300intR2Name("NA"),
301binCorrName("NA"),
302intBinCorrName("NA"),
303
304countsName("NA"),
305part_1_Name("NA"),
306part_2_Name("NA"),
307part_3_Name("NA"),
308pair_12_Name("NA"),
309pair_13_Name("NA"),
310pair_23_Name("NA"),
311tripletName("NA"),
312
313avg("NA"),
314avgName("NA"),
315sumName("NA"),
316s1ptName("NA"),
317s1ptNwName("NA"),
318s1DptName("NA"),
319
320s2PtPtName("NA"),
321s2NPtName("NA"),
322s2PtNName("NA"),
323s2DptDptName("NA"),
324
325s2PtPtNwName("NA"),
326s2NPtNwName("NA"),
327s2PtNNwName("NA"),
328
329ptName("NA"),
330ptptName("NA"),
331pt1pt1Name("NA"),
332DptName("NA"),
333DptDptName("NA"),
334RDptDptName("NA"),
335nPtName("NA"),
336ptNName("NA"),
337seanName("NA"),
338
339_title_counts("NA"),
340
341_title_m0("NA"),
342_title_m1("NA"),
343_title_m2("NA"),
344_title_m3("NA"),
345_title_m4("NA"),
346_title_m5("NA"),
347_title_m6("NA"),
348
349_title_eta_1("NA"),
350_title_phi_1("NA"),
351_title_pt_1("NA"),
352_title_etaPhi_1("NA"),
353_title_n_1("NA"),
354_title_SumPt_1("NA"),
355_title_AvgPt_1("NA"),
356_title_AvgN_1("NA"),
357_title_AvgSumPt_1("NA"),
358
359_title_eta_2("NA"),
360_title_phi_2("NA"),
361_title_pt_2("NA"),
362_title_etaPhi_2("NA"),
363_title_n_2("NA"),
364_title_SumPt_2("NA"),
365_title_AvgPt_2("NA"),
366_title_AvgN_2("NA"),
367_title_AvgSumPt_2("NA"),
368
369_title_etaPhi_12("NA"),
370
371_title_AvgN2_12("NA"),
372_title_AvgSumPtPt_12("NA"),
373_title_AvgSumPtN_12("NA"),
374_title_AvgNSumPt_12("NA"),
375
376vsZ("NA"),
377vsM("NA"),
378vsPt("NA"),
379vsPhi("NA"),
380vsEta("NA"),
381vsEtaPhi("NA"),
382vsPtVsPt("NA")
383{
384 printf("Default constructor called \n");
385
386 printf("passed \n ");
387
388}
389
390AliDptDptInMC::AliDptDptInMC(const TString & name)
391: AliAnalysisTaskSE(name),
392fAODEvent(0),
393fESDEvent(0),
394fInputHandler(0),
947cb887 395fPIDResponse(0),
1486c432 396_outputHistoList(0),
729278ef 397fArrayMC (0),
1486c432
MW
398_twoPi ( 2.0 * 3.1415927),
399_eventCount ( 0),
400_debugLevel ( 0),
401_singlesOnly ( 0),
402_useWeights ( 0),
403_sameFilter ( false),
404_rejectPileup ( 1),
405_rejectPairConversion ( 0),
406_vertexZMin ( -10.),
407_vertexZMax ( 10.),
408_vertexXYMin ( -10.),
409_vertexXYMax ( 10.),
410_centralityMethod ( 4),
411_centralityMin ( 0.),
412_centralityMax ( 1.),
413_requestedCharge_1 ( 1),
414_requestedCharge_2 ( -1),
415_dcaZMin ( -3),
416_dcaZMax ( 3.),
417_dcaXYMin ( -2.4),
418_dcaXYMax ( 2.4),
419_dedxMin ( 0),
420_dedxMax ( 100000),
421_nClusterMin ( 80),
422_trackFilterBit ( 0),
423fNSigmaCut ( 3.),
729278ef 424
425fAnalysisType ("MCAOD"),
426fExcludeResonancesInMC (kFALSE),
427fExcludeElectronsInMC (kFALSE),
428
1486c432
MW
429_field ( 1.),
430_nTracks ( 0 ),
431_mult0 ( 0 ),
432_mult1 ( 0 ),
433_mult2 ( 0 ),
434_mult3 ( 0 ),
435_mult4 ( 0 ),
436_mult4a ( 0 ),
437_mult5 ( 0 ),
438_mult6 ( 0 ),
439arraySize ( 2000),
440_id_1(0),
441_charge_1(0),
442_iEtaPhi_1(0),
443_iPt_1(0),
444_pt_1(0),
445_px_1(0),
446_py_1(0),
447_pz_1(0),
448_correction_1(0),
449_dedx_1(0),
450_id_2(0),
451_charge_2(0),
452_iEtaPhi_2(0),
453_iPt_2(0),
454_pt_2(0),
455_px_2(0),
456_py_2(0),
457_pz_2(0),
458_correction_2(0),
459_dedx_2(0),
460_correctionWeight_1(0),
461_correctionWeight_2(0),
462_nBins_M0(500), _min_M0(0), _max_M0(10000), _width_M0(20),
463_nBins_M1(500), _min_M1(0), _max_M1(10000), _width_M1(20),
464_nBins_M2(500), _min_M2(0), _max_M2(10000), _width_M2(20),
465_nBins_M3(500), _min_M3(0), _max_M3(10000), _width_M3(20),
466_nBins_M4(100), _min_M4(0), _max_M4(1), _width_M4(0.01),
467_nBins_M5(100), _min_M5(0), _max_M5(1), _width_M5(0.01),
468_nBins_M6(100), _min_M6(0), _max_M6(1), _width_M6(0.01),
469_nBins_vertexZ(40), _min_vertexZ(-10), _max_vertexZ(10), _width_vertexZ(0.5),
470
471_nBins_pt_1(18), _min_pt_1(0.2), _max_pt_1(2.0), _width_pt_1(0.1),
472_nBins_phi_1(72), _min_phi_1(0), _max_phi_1(2.*3.1415927),_width_phi_1(2.*3.1415927/72.),
473 //_nBins_eta_1(18),
474 //_min_eta_1(-0.9), _max_eta_1(0.9), _width_eta_1(0.1),
475 _nBins_eta_1(0),
476 _min_eta_1(0), _max_eta_1(0), _width_eta_1(0.1),
477
478_nBins_etaPhi_1(0),
479_nBins_etaPhiPt_1(0),
480_nBins_zEtaPhiPt_1(0),
481
482_nBins_pt_2(18), _min_pt_2(0.2), _max_pt_2(2.0), _width_pt_2(0.1),
483_nBins_phi_2(72), _min_phi_2(0), _max_phi_2(2.*3.1415927),_width_phi_2(2.*3.1415927/72),
484 //_nBins_eta_2(18),
485 //_min_eta_2(-0.9), _max_eta_2(0.9), _width_eta_2(0.1),
486 _nBins_eta_2(0),
487 _min_eta_2(0), _max_eta_2(0), _width_eta_2(0.1),
488
489
490_nBins_etaPhi_2(0),
491_nBins_etaPhiPt_2(0),
492_nBins_zEtaPhiPt_2(0),
493_nBins_etaPhi_12(0),
494__n1_1(0),
495__n1_2(0),
496__n2_12(0),
497__s1pt_1(0),
498__s1pt_2(0),
499__s2ptpt_12(0),
500__s2NPt_12(0),
501__s2PtN_12(0),
502__n1Nw_1(0),
503__n1Nw_2(0),
504__n2Nw_12(0),
505__s1ptNw_1(0),
506__s1ptNw_2(0),
507__s2ptptNw_12(0),
508__s2NPtNw_12(0),
509__s2PtNNw_12(0),
510__n1_1_vsPt(0),
511__n1_1_vsEtaPhi(0),
512__s1pt_1_vsEtaPhi(0),
513__n1_1_vsZEtaPhiPt(0),
514__n1_2_vsPt(0),
515__n1_2_vsEtaPhi(0),
516__s1pt_2_vsEtaPhi(0),
517__n1_2_vsZEtaPhiPt(0),
518__n2_12_vsPtPt(0),
519__n2_12_vsEtaPhi(0),
520__s2ptpt_12_vsEtaPhi(0),
521__s2PtN_12_vsEtaPhi(0),
522__s2NPt_12_vsEtaPhi(0),
523_weight_1 ( 0 ),
524_weight_2 ( 0 ),
525_eventAccounting ( 0),
526_m0 ( 0),
527_m1 ( 0),
528_m2 ( 0),
529_m3 ( 0),
530_m4 ( 0),
531_m5 ( 0),
532_m6 ( 0),
533_vertexZ ( 0),
534_Ncluster1 ( 0),
535_Ncluster2 ( 0),
536_etadis ( 0),
537_phidis ( 0),
538
539_dcaz ( 0),
540_dcaxy ( 0),
0f50230b 541_etadis2 ( 0),
542_phidis2 ( 0),
543_dcaz2 ( 0),
544_dcaxy2 ( 0),
545_etadis3 ( 0),
546_phidis3 ( 0),
547_dcaz3 ( 0),
548_dcaxy3 ( 0),
1486c432
MW
549_n1_1_vsPt ( 0),
550_n1_1_vsEtaVsPhi ( 0),
551_s1pt_1_vsEtaVsPhi ( 0),
552_n1_1_vsZVsEtaVsPhiVsPt ( 0),
553_n1_1_vsM ( 0), // w/ weight
554_s1pt_1_vsM ( 0),
555_n1Nw_1_vsM ( 0), // w/o weight
556_s1ptNw_1_vsM ( 0),
557_dedxVsP_1 ( 0),
558_corrDedxVsP_1 ( 0),
559_betaVsP_1 ( 0),
560_n1_2_vsPt ( 0),
561_n1_2_vsEtaVsPhi ( 0),
562_s1pt_2_vsEtaVsPhi ( 0),
563_n1_2_vsZVsEtaVsPhiVsPt ( 0),
564_n1_2_vsM ( 0),
565_s1pt_2_vsM ( 0),
566_n1Nw_2_vsM ( 0),
567_s1ptNw_2_vsM ( 0),
568_dedxVsP_2 ( 0),
569_corrDedxVsP_2 ( 0),
570_betaVsP_2 ( 0),
571_n2_12_vsEtaPhi ( 0),
572_n2_12_vsPtVsPt ( 0),
573_s2PtPt_12_vsEtaPhi( 0),
574_s2PtN_12_vsEtaPhi ( 0),
575_s2NPt_12_vsEtaPhi ( 0),
576_n2_12_vsM ( 0),
577_s2PtPt_12_vsM ( 0),
578_s2PtN_12_vsM ( 0),
579_s2NPt_12_vsM ( 0),
580_n2Nw_12_vsM ( 0),
581_s2PtPtNw_12_vsM ( 0),
582_s2PtNNw_12_vsM ( 0),
583_s2NPtNw_12_vsM ( 0),
584_invMass ( 0),
585_invMassElec ( 0),
586n1Name("NA"),
587n1NwName("NA"),
588n2Name("NA"),
589n2NwName("NA"),
590n3Name("NA"),
591n1n1Name("NA"),
592n1n1n1Name("NA"),
593n2n1Name("NA"),
594r1Name("NA"),
595r2Name("NA"),
596r3Name("NA"),
597r2r1Name("NA"),
598c2Name("NA"),
599c3Name("NA"),
600d3Name("NA"),
601p3Name("NA"),
602cName("NA"),
603
604intR2Name("NA"),
605binCorrName("NA"),
606intBinCorrName("NA"),
607
608countsName("NA"),
609part_1_Name("NA"),
610part_2_Name("NA"),
611part_3_Name("NA"),
612pair_12_Name("NA"),
613pair_13_Name("NA"),
614pair_23_Name("NA"),
615tripletName("NA"),
616
617avg("NA"),
618avgName("NA"),
619sumName("NA"),
620s1ptName("NA"),
621s1ptNwName("NA"),
622s1DptName("NA"),
623
624s2PtPtName("NA"),
625s2NPtName("NA"),
626s2PtNName("NA"),
627s2DptDptName("NA"),
628
629s2PtPtNwName("NA"),
630s2NPtNwName("NA"),
631s2PtNNwName("NA"),
632
633ptName("NA"),
634ptptName("NA"),
635pt1pt1Name("NA"),
636DptName("NA"),
637DptDptName("NA"),
638RDptDptName("NA"),
639nPtName("NA"),
640ptNName("NA"),
641seanName("NA"),
642
643_title_counts("NA"),
644
645_title_m0("NA"),
646_title_m1("NA"),
647_title_m2("NA"),
648_title_m3("NA"),
649_title_m4("NA"),
650_title_m5("NA"),
651_title_m6("NA"),
652
653_title_eta_1("NA"),
654_title_phi_1("NA"),
655_title_pt_1("NA"),
656_title_etaPhi_1("NA"),
657_title_n_1("NA"),
658_title_SumPt_1("NA"),
659_title_AvgPt_1("NA"),
660_title_AvgN_1("NA"),
661_title_AvgSumPt_1("NA"),
662
663_title_eta_2("NA"),
664_title_phi_2("NA"),
665_title_pt_2("NA"),
666_title_etaPhi_2("NA"),
667_title_n_2("NA"),
668_title_SumPt_2("NA"),
669_title_AvgPt_2("NA"),
670_title_AvgN_2("NA"),
671_title_AvgSumPt_2("NA"),
672
673_title_etaPhi_12("NA"),
674
675_title_AvgN2_12("NA"),
676_title_AvgSumPtPt_12("NA"),
677_title_AvgSumPtN_12("NA"),
678_title_AvgNSumPt_12("NA"),
679
680vsZ("NA"),
681vsM("NA"),
682vsPt("NA"),
683vsPhi("NA"),
684vsEta("NA"),
685vsEtaPhi("NA"),
686vsPtVsPt("NA")
687{
688 printf("2nd constructor called ");
689
690 DefineOutput(0, TList::Class());
691
692 printf("passed ");
693
694}
695
696AliDptDptInMC::~AliDptDptInMC()
697{
698 /*
699 delete _id_1;
700 delete _charge_1;
701 delete _iEtaPhi_1;
702 delete _iPt_1;
703 delete _pt_1;
704 delete _px_1;
705 delete _py_1;
706 delete _pz_1;
707 delete _correction_1;
708 delete _dedx_1;
709 delete __n1_1_vsPt;
710 delete __n1_1_vsEtaPhi;
711 delete __s1pt_1_vsEtaPhi;
712 delete __n1_1_vsZEtaPhiPt;
713 if (_correctionWeight_1) delete _correctionWeight_1;
714
715 if (!_sameFilter)
716 {
717 delete _id_2;
718 delete _charge_2;
719 delete _iEtaPhi_2;
720 delete _iPt_2;
721 delete _pt_2;
722 delete _px_2;
723 delete _py_2;
724 delete _pz_2;
725 delete _correction_2;
726 delete _dedx_2;
727 delete __n1_2_vsPt;
728 delete __n1_2_vsEtaPhi;
729 delete __s1pt_2_vsEtaPhi;
730 delete __n1_2_vsZEtaPhiPt;
731 if (_correctionWeight_2) delete _correctionWeight_2;
732 }
733
734 if (!_singlesOnly)
735 {
736 delete __n2_12_vsPtPt;
737 delete __n2_12_vsEtaPhi;
738 delete __s2ptpt_12_vsEtaPhi;
739 delete __s2PtN_12_vsEtaPhi;
740 delete __s2NPt_12_vsEtaPhi;
741 }
742 */
743}
744
745void AliDptDptInMC::UserCreateOutputObjects()
746{
747 //cout<< "AliDptDptInMC::CreateOutputObjects() Starting " << endl;
748 OpenFile(0);
749 _outputHistoList = new TList();
750 _outputHistoList->SetOwner();
751
752 //if (_useWeights) DefineInput(2, TList::Class());
753
754 //Setup the parameters of the histograms
755 _nBins_M0 = 500; _min_M0 = 0.; _max_M0 = 5000.; _width_M0 = (_max_M0-_min_M0)/_nBins_M0;
756 _nBins_M1 = 500; _min_M1 = 0.; _max_M1 = 5000.; _width_M1 = (_max_M1-_min_M1)/_nBins_M1;
757 _nBins_M2 = 500; _min_M2 = 0.; _max_M2 = 5000.; _width_M2 = (_max_M2-_min_M2)/_nBins_M2;
758 _nBins_M3 = 500; _min_M3 = 0.; _max_M3 = 5000.; _width_M3 = (_max_M3-_min_M3)/_nBins_M3;
759 _nBins_M4 = 100; _min_M4 = 0.; _max_M4 = 100.; _width_M4 = (_max_M4-_min_M4)/_nBins_M4;
760 _nBins_M5 = 100; _min_M5 = 0.; _max_M5 = 100.; _width_M5 = (_max_M5-_min_M5)/_nBins_M5;
761 _nBins_M6 = 100; _min_M6 = 0.; _max_M6 = 100.; _width_M6 = (_max_M6-_min_M6)/_nBins_M6;
762
763 _min_vertexZ = _vertexZMin;
764 _max_vertexZ = _vertexZMax;
765 _width_vertexZ = 0.5;
766 _nBins_vertexZ = int(0.5+ (_max_vertexZ - _min_vertexZ)/_width_vertexZ);
767 _nBins_pt_1 = int(0.5+ (_max_pt_1 -_min_pt_1 )/_width_pt_1);
768 _nBins_eta_1 = int(0.5+ (_max_eta_1-_min_eta_1)/_width_eta_1);
769 _width_phi_1 = (_max_phi_1 - _min_phi_1) /_nBins_phi_1;
770 _nBins_etaPhi_1 = _nBins_phi_1 * _nBins_eta_1;
771 _nBins_etaPhiPt_1 = _nBins_etaPhi_1 * _nBins_pt_1;
772 _nBins_zEtaPhiPt_1 = _nBins_vertexZ * _nBins_etaPhiPt_1;
773
774 _nBins_pt_2 = int(0.5+ (_max_pt_2 -_min_pt_2 )/_width_pt_2);
775 _nBins_eta_2 = int(0.5+ (_max_eta_2-_min_eta_2)/_width_eta_2);
776 _width_phi_2 = (_max_phi_2 - _min_phi_2) /_nBins_phi_2;
777 _nBins_etaPhi_2 = _nBins_phi_2 * _nBins_eta_2;
778 _nBins_etaPhiPt_2 = _nBins_etaPhi_2 * _nBins_pt_2;
779 _nBins_zEtaPhiPt_2 = _nBins_vertexZ * _nBins_etaPhiPt_2;
780 _nBins_etaPhi_12 = _nBins_etaPhi_1 * _nBins_etaPhi_2;
781
782 //setup the work arrays
783
784 _id_1 = new int[arraySize];
785 _charge_1 = new int[arraySize];
786 //_iPhi_1 = new int[arraySize];
787 //_iEta_1 = new int[arraySize];
788 _iEtaPhi_1 = new int[arraySize];
789 _iPt_1 = new int[arraySize];
790 _pt_1 = new float[arraySize];
791 _px_1 = new float[arraySize];
792 _py_1 = new float[arraySize];
793 _pz_1 = new float[arraySize];
794 //_phi_1 = new float[arraySize];
795 //_eta_1 = new float[arraySize];
796 _correction_1 = new float[arraySize];
797 //_dedx_1 = new float[arraySize];
798
799 __n1_1_vsPt = getDoubleArray(_nBins_pt_1, 0.);
947cb887 800
1486c432
MW
801 __n1_1_vsEtaPhi = getDoubleArray(_nBins_etaPhi_1, 0.);
802 __s1pt_1_vsEtaPhi = getDoubleArray(_nBins_etaPhi_1, 0.);
803 __n1_1_vsZEtaPhiPt = getFloatArray(_nBins_zEtaPhiPt_1, 0.);
804
805 //cout << "==========================================================================================" << endl;
806 //cout << "=============== Booking for particle 1 done." << endl;
807 //cout << "_requestedCharge_1: " << _requestedCharge_1 << endl;
808 //cout << "_requestedCharge_2: " << _requestedCharge_2 << endl;
809
810 if (_requestedCharge_2!=_requestedCharge_1)
811 {
812 //cout << " creating arrays for particle 2 with size: " << arraySize << endl;
813 _sameFilter = 0;
814 //particle 2
815 _id_2 = new int[arraySize];
816 _charge_2 = new int[arraySize];
817 //_iPhi_2 = new int[arraySize];
818 //_iEta_2 = new int[arraySize];
819 _iEtaPhi_2 = new int[arraySize];
820 _iPt_2 = new int[arraySize];
821 _pt_2 = new float[arraySize];
822 _px_2 = new float[arraySize];
823 _py_2 = new float[arraySize];
824 _pz_2 = new float[arraySize];
825 //_phi_2 = new float[arraySize];
826 //_eta_2 = new float[arraySize];
827 _correction_2 = new float[arraySize];
828 //_dedx_2 = new float[arraySize];
829
830 __n1_2_vsPt = getDoubleArray(_nBins_pt_2, 0.);
831 __n1_2_vsEtaPhi = getDoubleArray(_nBins_etaPhi_2, 0.);
832 __s1pt_2_vsEtaPhi = getDoubleArray(_nBins_etaPhi_2, 0.);
833 __n1_2_vsZEtaPhiPt = getFloatArray(_nBins_zEtaPhiPt_2, 0.);
834
835 }
836
837 __n2_12_vsPtPt = getDoubleArray(_nBins_pt_1*_nBins_pt_2,0.);
838 __n2_12_vsEtaPhi = getFloatArray(_nBins_etaPhi_12, 0.);
839 __s2ptpt_12_vsEtaPhi = getFloatArray(_nBins_etaPhi_12, 0.);
840 __s2PtN_12_vsEtaPhi = getFloatArray(_nBins_etaPhi_12, 0.);
841 __s2NPt_12_vsEtaPhi = getFloatArray(_nBins_etaPhi_12, 0.);
842
843 // Setup all the labels needed.
844
845 part_1_Name = "_1";
846 part_2_Name = "_2";
847 pair_12_Name = "_12";
848
849 n1Name = "n1";
850 n2Name = "n2";
851 n1NwName = "n1Nw";
852 n2NwName = "n2Nw";
853 r1Name = "r1";
854 r2Name = "r2";
855 r3Name = "r3";
856 r2r1Name = "r2r1";
857 c2Name = "c2";
858 c3Name = "c3";
859 d3Name = "d3";
860 p3Name = "p3";
861 cName = "sean";
862
863 intR2Name = "intR2";
864 binCorrName = "binCorr";
865 intBinCorrName = "intBinCorr";
866
867 avgName = "avg";
868 sumName = "sum";
869 s1ptName = "sumPt";
870 s1ptNwName = "sumPtNw";
871 s1DptName = "sumDpt";
872 s2PtPtName = "sumPtPt";
873 s2PtPtNwName = "sumPtPtNw";
874 s2DptDptName = "sumDptDpt";
875 s2NPtName = "sumNPt";
876 s2NPtNwName = "sumNPtNw";
877 s2PtNName = "sumPtN";
878 s2NPtNwName = "sumNPtNw";
879 s2PtNNwName = "sumPtNNw";
880 ptName = "avgPt";
881 ptptName = "avgPtPt";
882 pt1pt1Name = "avgPtavgPt";
883 DptName = "avgDpt";
884 DptDptName = "avgDptDpt";
885 RDptDptName = "relDptDpt"; // ratio of avgDptDpt by avgPt*avgPt
886 nPtName = "avgNpt";
887 ptNName = "avgPtN";
888 seanName = "seanC";
889
890 _title_counts = "yield";
891
892 _title_m0 = "M_{0}";
893 _title_m1 = "M_{1}";
894 _title_m2 = "M_{2}";
895 _title_m3 = "M_{3}";
896 _title_m4 = "V0Centrality";
897 _title_m5 = "TrkCentrality";
898 _title_m6 = "SpdCentrality";
899
900 _title_eta_1 = "#eta_{1}";
901 _title_phi_1 = "#varphi_{1} (radian)";
902 _title_etaPhi_1 = "#eta_{1}#times#varphi_{1}";
903 _title_pt_1 = "p_{t,1} (GeV/c)";
904 _title_n_1 = "n_{1}";
905 _title_SumPt_1 = "#Sigma p_{t,1} (GeV/c)";
906 _title_AvgPt_1 = "#LT p_{t,1} #GT (GeV/c)";
907 _title_AvgN_1 = "#LT n_{1} #GT";
908 _title_AvgSumPt_1 = "#LT #Sigma p_{t,1} #GT (GeV/c)";
909
910 _title_eta_2 = "#eta_{2}";
911 _title_phi_2 = "#varphi_{2} (radian)";
912 _title_etaPhi_2 = "#eta_{2}#times#varphi_{2}";
913 _title_pt_2 = "p_{t,2} (GeV/c)";
914 _title_n_2 = "n_{2}";
915 _title_SumPt_2 = "#Sigma p_{t,1} (GeV/c)";
916 _title_AvgPt_2 = "#LT p_{t,2} #GT (GeV/c)";
917 _title_AvgN_2 = "#LT n_{2} #GT";
918 _title_AvgSumPt_2 = "#LT #Sigma p_{t,2} #GT (GeV/c)";
919
920 _title_etaPhi_12 = "#eta_{1}#times#varphi_{1}#times#eta_{2}#times#varphi_{2}";
921
922 _title_AvgN2_12 = "#LT n_{2} #GT";;
923 _title_AvgSumPtPt_12 = "#LT #Sigma p_{t,1}p_{t,2} #GT";;
924 _title_AvgSumPtN_12 = "#LT #Sigma p_{t,1}N #GT";;
925 _title_AvgNSumPt_12 = "#LT N#Sigma p_{t,2} #GT";;
926
927
928 vsZ = "_vsZ";
929 vsM = "_vsM";
930 vsPt = "_vsPt";
931 vsPhi = "_vsPhi";
932 vsEta = "_vsEta";
933 vsEtaPhi = "_vsEtaPhi";
934 vsPtVsPt = "_vsPtVsPt";
935
936
937 if (_useWeights)
938 {
939 int iZ, iEtaPhi, iPt;
940 int iZ1,iEtaPhi1,iPt1;
941 int a, b;
942 if (_weight_1)
943 {
944 _correctionWeight_1 = new float[_nBins_vertexZ*_nBins_etaPhi_1*_nBins_pt_1];
945 a = _nBins_pt_1;
946 b = _nBins_etaPhi_1*_nBins_pt_1;
947 for (iZ=0,iZ1=1; iZ<_nBins_vertexZ; iZ++, iZ1++)
948 {
949 for (iEtaPhi=0,iEtaPhi1=1; iEtaPhi<_nBins_etaPhi_1; iEtaPhi++, iEtaPhi1++)
950 {
951 for (iPt=0,iPt1=1; iPt<_nBins_pt_1; iPt++, iPt1++)
952 {
953 _correctionWeight_1[iZ*b+iEtaPhi*a+iPt] = _weight_1->GetBinContent(iZ1,iEtaPhi1,iPt1);
954 }
955 }
956 }
957 } // _weight_1
958 else
959 {
960 AliError("AliDptDptInMC:: _weight_1 is a null pointer.");
961 return;
962 }
963 if (!_sameFilter)
964 {
965 if (_weight_2)
966 {
967 _correctionWeight_2 = new float[_nBins_vertexZ*_nBins_etaPhi_2*_nBins_pt_2];
968 a = _nBins_pt_2;
969 b = _nBins_etaPhi_2*_nBins_pt_2;
970 for (iZ=0,iZ1=1; iZ<_nBins_vertexZ; iZ++, iZ1++)
971 {
972 for (iEtaPhi=0,iEtaPhi1=1; iEtaPhi<_nBins_etaPhi_2; iEtaPhi++, iEtaPhi1++)
973 {
974 for (iPt=0,iPt1=1; iPt<_nBins_pt_2; iPt++, iPt1++)
975 {
976 _correctionWeight_2[iZ*b+iEtaPhi*a+iPt] = _weight_2->GetBinContent(iZ1,iEtaPhi1,iPt1);
977 }
978 }
979 }
980 } // _weight_2
981 else
982 {
983 AliError("AliDptDptInMC:: _weight_1 is a null pointer.");
984 return;
985 }
986 }
987 }
988
989 createHistograms();
990 PostData(0,_outputHistoList);
991
992 //cout<< "AliDptDptInMC::CreateOutputObjects() DONE " << endl;
993
994}
995
996void AliDptDptInMC::createHistograms()
997{
998 AliInfo(" AliDptDptInMC::createHistoHistograms() Creating Event Histos");
999 TString name;
1000
1001 name = "eventAccounting";
1002
1003 // bin index : what it is...
1004 // 0 : number of event submitted
1005 // 1 : number accepted by centrality cut
1006 // 2 : number accepted by centrality cut and z cut
1007 // 3 : total number of particles that satisfied filter 1
1008 // 4 : total number of particles that satisfied filter 2
1009 _eventAccounting = createHisto1D(name,name,10, -0.5, 9.5, "event Code", _title_counts);
1010
1011 name = "m0"; _m0 = createHisto1D(name,name,_nBins_M1, _min_M1, _max_M1, _title_m0, _title_counts);
1012 name = "m1"; _m1 = createHisto1D(name,name,_nBins_M1, _min_M1, _max_M1, _title_m1, _title_counts);
1013 name = "m2"; _m2 = createHisto1D(name,name,_nBins_M2, _min_M2, _max_M2, _title_m2, _title_counts);
1014 name = "m3"; _m3 = createHisto1D(name,name,_nBins_M3, _min_M3, _max_M3, _title_m3, _title_counts);
1015 name = "m4"; _m4 = createHisto1D(name,name,_nBins_M4, _min_M4, _max_M4, _title_m4, _title_counts);
1016 name = "m5"; _m5 = createHisto1D(name,name,_nBins_M5, _min_M5, _max_M5, _title_m5, _title_counts);
1017 name = "m6"; _m6 = createHisto1D(name,name,_nBins_M6, _min_M6, _max_M6, _title_m6, _title_counts);
1018 name = "zV"; _vertexZ = createHisto1D(name,name,_nBins_vertexZ, _min_vertexZ, _max_vertexZ, "z-Vertex (cm)", _title_counts);
1019
1020
0f50230b 1021 name = "Eta"; _etadis = createHisto1F(name,name, 200, -1.0, 1.0, "#eta","counts");
1022 name = "Phi"; _phidis = createHisto1F(name,name, 360, 0.0, 6.4, "#phi","counts");
1023 name = "DCAz"; _dcaz = createHisto1F(name,name, 340, -3.3, 3.3, "dcaZ","counts");
1024 name = "DCAxy"; _dcaxy = createHisto1F(name,name, 100, -0.1, 2.5, "dcaXY","counts");
947cb887 1025
0f50230b 1026 name = "Eta2"; _etadis2 = createHisto1F(name,name, 200, -1.0, 1.0, "#eta","counts");
1027 name = "Phi2"; _phidis2 = createHisto1F(name,name, 360, 0.0, 6.4, "#phi","counts");
1028 name = "DCAz2"; _dcaz2 = createHisto1F(name,name, 340, -3.3, 3.3, "dcaZ","counts");
1029 name = "DCAxy2"; _dcaxy2 = createHisto1F(name,name, 100, -0.1, 2.5, "dcaXY","counts");
1030
1031 name = "Eta3"; _etadis3 = createHisto1F(name,name, 200, -1.0, 1.0, "#eta","counts");
1032 name = "Phi3"; _phidis3 = createHisto1F(name,name, 360, 0.0, 6.4, "#phi","counts");
1033 name = "DCAz3"; _dcaz3 = createHisto1F(name,name, 340, -3.3, 3.3, "dcaZ","counts");
1034 name = "DCAxy3"; _dcaxy3 = createHisto1F(name,name, 100, -0.1, 2.5, "dcaXY","counts");
1035
1036 //name = "Eta"; _etadis = createHisto1F(name,name, 250, 0.0, 2.5, "#eta","counts"); //temporaryly
1037 //name = "Phi"; _phidis = createHisto1F(name,name, 250, 0.0, 2.5, "#phi","counts");
1038 //name = "DCAz"; _dcaz = createHisto1F(name,name, 250, 0.0, 2.5, "dcaZ","counts");
564a0b8f 1039
1040 //name = "Nclus1"; _Ncluster1 = createHisto1F(name,name, 200, 0, 200, "Ncluster1","counts");
1041 //name = "Nclus2"; _Ncluster2 = createHisto1F(name,name, 200, 0, 200, "Ncluster2","counts");
1486c432
MW
1042
1043 if (_singlesOnly)
1044 {
1045 name = n1Name+part_1_Name+vsPt; _n1_1_vsPt = createHisto1F(name,name, _nBins_pt_1, _min_pt_1, _max_pt_1, _title_pt_1, _title_AvgN_1);
947cb887 1046
1486c432
MW
1047 name = n1Name+part_1_Name+vsZ+vsEtaPhi+vsPt; _n1_1_vsZVsEtaVsPhiVsPt = createHisto3F(name,name, _nBins_vertexZ,_min_vertexZ,_max_vertexZ, _nBins_etaPhi_1, 0., double(_nBins_etaPhi_1), _nBins_pt_1, _min_pt_1, _max_pt_1, "zVertex", _title_etaPhi_1, _title_pt_1);
1048 //name = "dedxVsP_1"; _dedxVsP_1 = createHisto2F(name,name,400,-2.,2.,120,0.,120.,"p (GeV/c)", "dedx", "counts");
1049 //name = "corrDedxVsP_1"; _corrDedxVsP_1 = createHisto2F(name,name,400,-2.,2.,120,0.,120.,"p (GeV/c)", "dedx", "counts");
1050 //name = "betaVsP_1"; _betaVsP_1 = createHisto2F(name,name,400,-2.,2.,120,0.5,1.1,"p (GeV/c)", "beta", "counts");
1051
1052 name = n1Name+part_2_Name+vsPt; _n1_2_vsPt = createHisto1F(name,name, _nBins_pt_2, _min_pt_2, _max_pt_2, _title_pt_2, _title_AvgN_2);
1053 name = n1Name+part_2_Name+vsZ+vsEtaPhi+vsPt; _n1_2_vsZVsEtaVsPhiVsPt = createHisto3F(name,name, _nBins_vertexZ,_min_vertexZ,_max_vertexZ, _nBins_etaPhi_2, 0., double(_nBins_etaPhi_2), _nBins_pt_2, _min_pt_2, _max_pt_2, "zVertex", _title_etaPhi_2, _title_pt_2);
1054 //name = "dedxVsP_2"; _dedxVsP_2 = createHisto2F(name,name,400,-2.,2.,120,0.,120.,"p (GeV/c)", "dedx", "counts");
1055 //name = "corrDedxVsP_2"; _corrDedxVsP_2 = createHisto2F(name,name,400,-2.,2.,120,0.,120.,"p (GeV/c)", "dedx", "counts");
1056 //name = "betaVsP_2"; _betaVsP_2 = createHisto2F(name,name,400,-2.,2.,120,0.5,1.1,"p (GeV/c)", "beta", "counts");
1057
1058 }
1059 else
1060 {
1061 name = n1Name+part_1_Name+vsEtaPhi; _n1_1_vsEtaVsPhi = createHisto2F(name,name, _nBins_eta_1, _min_eta_1, _max_eta_1, _nBins_phi_1, _min_phi_1, _max_phi_1, _title_eta_1, _title_phi_1, _title_AvgN_1);
1062 name = s1ptName+part_1_Name+vsEtaPhi; _s1pt_1_vsEtaVsPhi = createHisto2F(name,name, _nBins_eta_1, _min_eta_1, _max_eta_1, _nBins_phi_1, _min_phi_1, _max_phi_1, _title_eta_1, _title_phi_1, _title_AvgSumPt_1);
1063 name = n1Name+part_1_Name+vsM; _n1_1_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgN_1);
1064 name = s1ptName+part_1_Name+vsM; _s1pt_1_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgSumPt_1);
1065 name = n1NwName+part_1_Name+vsM; _n1Nw_1_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgN_1);
1066 name = s1ptNwName+part_1_Name+vsM; _s1ptNw_1_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgSumPt_1);
1067
1068 name = n1Name+part_2_Name+vsEtaPhi; _n1_2_vsEtaVsPhi = createHisto2F(name,name, _nBins_eta_2, _min_eta_2, _max_eta_2, _nBins_phi_2, _min_phi_2, _max_phi_2, _title_eta_2, _title_phi_2, _title_AvgN_2);
1069 name = s1ptName+part_2_Name+vsEtaPhi; _s1pt_2_vsEtaVsPhi = createHisto2F(name,name, _nBins_eta_2, _min_eta_2, _max_eta_2, _nBins_phi_2, _min_phi_2, _max_phi_2, _title_eta_2, _title_phi_2, _title_AvgSumPt_2);
1070 name = n1Name+part_2_Name + vsM; _n1_2_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgN_2);
1071 name = s1ptName+part_2_Name + vsM; _s1pt_2_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgSumPt_2);
1072 name = n1NwName+part_2_Name+vsM; _n1Nw_2_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgN_1);
1073 name = s1ptNwName+part_2_Name+vsM; _s1ptNw_2_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgSumPt_1);
1074
1075 name = n2Name+pair_12_Name+vsEtaPhi; _n2_12_vsEtaPhi = createHisto1F(name,name, _nBins_etaPhi_12, 0., double(_nBins_etaPhi_12), _title_etaPhi_12, _title_AvgN2_12);
1076 name = s2PtPtName+pair_12_Name + vsEtaPhi;_s2PtPt_12_vsEtaPhi = createHisto1F(name,name, _nBins_etaPhi_12, 0., double(_nBins_etaPhi_12), _title_etaPhi_12, _title_AvgSumPtPt_12);
1077 name = s2PtNName+pair_12_Name + vsEtaPhi; _s2PtN_12_vsEtaPhi = createHisto1F(name,name, _nBins_etaPhi_12, 0., double(_nBins_etaPhi_12), _title_etaPhi_12, _title_AvgSumPtN_12);
1078 name = s2NPtName+pair_12_Name + vsEtaPhi; _s2NPt_12_vsEtaPhi = createHisto1F(name,name, _nBins_etaPhi_12, 0., double(_nBins_etaPhi_12), _title_etaPhi_12, _title_AvgNSumPt_12);
1079 name = n2Name+pair_12_Name+vsPtVsPt; _n2_12_vsPtVsPt = createHisto2F(name,name, _nBins_pt_1, _min_pt_1, _max_pt_1, _nBins_pt_2, _min_pt_2, _max_pt_2, _title_pt_1, _title_pt_2, _title_AvgN2_12);
1080
1081 name = n2Name+pair_12_Name + vsM; _n2_12_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgN2_12);
1082 name = s2PtPtName+pair_12_Name + vsM; _s2PtPt_12_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgSumPtPt_12);
1083 name = s2PtNName+pair_12_Name + vsM; _s2PtN_12_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgSumPtN_12);
1084 name = s2NPtName+pair_12_Name + vsM; _s2NPt_12_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgNSumPt_12);
1085
1086 name = n2NwName+pair_12_Name + vsM; _n2Nw_12_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgN2_12);
1087 name = s2PtPtNwName+pair_12_Name + vsM; _s2PtPtNw_12_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgSumPtPt_12);
1088 name = s2PtNNwName+pair_12_Name + vsM; _s2PtNNw_12_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgSumPtN_12);
1089 name = s2NPtNwName+pair_12_Name + vsM; _s2NPtNw_12_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgNSumPt_12);
1090
1091 name = "mInv"; _invMass = createHisto1F(name,name, 500, 0., 1.000, "M_{inv}","counts");
1092 name = "mInvElec"; _invMassElec = createHisto1F(name,name, 500, 0., 1.000, "M_{inv}","counts");
1093 }
1094
1095 AliInfo(" AliDptDptInMC::createHistoHistograms() All Done");
1096}
1097//-----------------------//
1098
1099void AliDptDptInMC::finalizeHistograms()
1100{
1101
1102 AliInfo("AliDptDptInMC::finalizeHistograms() starting");
1103 AliInfo(Form("CorrelationAnalyzers::finalizeHistograms() _eventCount : %d",int(_eventCount)));
1104 if (_singlesOnly)
1105 {
1106 if (_sameFilter)
1107 {
1108 fillHistoWithArray(_n1_1_vsPt, __n1_1_vsPt, _nBins_pt_1);
1109 fillHistoWithArray(_n1_1_vsZVsEtaVsPhiVsPt, __n1_1_vsZEtaPhiPt, _nBins_vertexZ, _nBins_etaPhi_1, _nBins_pt_1);
1110 fillHistoWithArray(_n1_2_vsPt, __n1_1_vsPt, _nBins_pt_1);
1111 fillHistoWithArray(_n1_2_vsZVsEtaVsPhiVsPt, __n1_1_vsZEtaPhiPt, _nBins_vertexZ, _nBins_etaPhi_1, _nBins_pt_1);
1112 }
1113 else
1114 {
1115 fillHistoWithArray(_n1_1_vsPt, __n1_1_vsPt, _nBins_pt_1);
1116 fillHistoWithArray(_n1_1_vsZVsEtaVsPhiVsPt, __n1_1_vsZEtaPhiPt, _nBins_vertexZ, _nBins_etaPhi_1, _nBins_pt_1);
1117 fillHistoWithArray(_n1_2_vsPt, __n1_2_vsPt, _nBins_pt_2);
1118 fillHistoWithArray(_n1_2_vsZVsEtaVsPhiVsPt, __n1_2_vsZEtaPhiPt, _nBins_vertexZ, _nBins_etaPhi_2, _nBins_pt_2);
1119 }
1120 }
1121 else
1122 {
1123 if (_sameFilter)
1124 {
1125 fillHistoWithArray(_n1_1_vsEtaVsPhi, __n1_1_vsEtaPhi, _nBins_eta_1, _nBins_phi_1);
1126 fillHistoWithArray(_s1pt_1_vsEtaVsPhi, __s1pt_1_vsEtaPhi, _nBins_eta_1, _nBins_phi_1);
1127 fillHistoWithArray(_n1_2_vsEtaVsPhi, __n1_1_vsEtaPhi, _nBins_eta_1, _nBins_phi_1);
1128 fillHistoWithArray(_s1pt_2_vsEtaVsPhi, __s1pt_1_vsEtaPhi, _nBins_eta_1, _nBins_phi_1);
1129 }
1130 else
1131 {
1132 fillHistoWithArray(_n1_1_vsEtaVsPhi, __n1_1_vsEtaPhi, _nBins_eta_1, _nBins_phi_1);
1133 fillHistoWithArray(_s1pt_1_vsEtaVsPhi, __s1pt_1_vsEtaPhi, _nBins_eta_1, _nBins_phi_1);
1134 fillHistoWithArray(_n1_2_vsEtaVsPhi, __n1_2_vsEtaPhi, _nBins_eta_2, _nBins_phi_2);
1135 fillHistoWithArray(_s1pt_2_vsEtaVsPhi, __s1pt_2_vsEtaPhi, _nBins_eta_2, _nBins_phi_2);
1136 }
1137 fillHistoWithArray(_n2_12_vsEtaPhi, __n2_12_vsEtaPhi, _nBins_etaPhi_12);
1138 fillHistoWithArray(_s2PtPt_12_vsEtaPhi, __s2ptpt_12_vsEtaPhi, _nBins_etaPhi_12);
1139 fillHistoWithArray(_s2PtN_12_vsEtaPhi, __s2PtN_12_vsEtaPhi, _nBins_etaPhi_12);
1140 fillHistoWithArray(_s2NPt_12_vsEtaPhi, __s2NPt_12_vsEtaPhi, _nBins_etaPhi_12);
1141 fillHistoWithArray(_n2_12_vsPtVsPt, __n2_12_vsPtPt, _nBins_pt_1, _nBins_pt_2);
1142
1143 }
1144 AliInfo("AliDptDptInMC::finalizeHistograms() Done ");
1145}
1146//--------------//
1147
1148
1149void AliDptDptInMC::UserExec(Option_t */*option*/)
1150{
1151
1152
1153 int k1,k2;
1154 int iPhi, iEta, iEtaPhi, iPt, charge;
1155 float q, phi, pt, eta, corr, corrPt, px, py, pz;
1156 int ij;
1157 int id_1, q_1, iEtaPhi_1, iPt_1;
1158 float pt_1, px_1, py_1, pz_1, corr_1;
1159 int id_2, q_2, iEtaPhi_2, iPt_2;
1160 float pt_2, px_2, py_2, pz_2, corr_2;
1161 float ptpt;
1162 int iVertex, iVertexP1, iVertexP2;
1163 int iZEtaPhiPt;
1164 float massElecSq = 2.5e-7;
1165 const AliAODVertex* vertex;
1166 //int nClus;
1167 bool bitOK;
1168
1169 AliAnalysisManager* manager = AliAnalysisManager::GetAnalysisManager();
1170 if (!manager) {
1171 //cout<<"ERROR: Analysis manager not found."<<endl;
1172 return;
1173 }
1174 //coneect to the inputHandler------------
1175 AliAODInputHandler* inputHandler = dynamic_cast<AliAODInputHandler*> (manager->GetInputEventHandler());
1176 if (!inputHandler) {
1177 return;
1178 }
1179
1180 fAODEvent = dynamic_cast<AliAODEvent*>(InputEvent());
1181 if (!fAODEvent)
1182 {
1183 return;
1184 }
947cb887 1185 fPIDResponse =inputHandler->GetPIDResponse();
1186 if (!fPIDResponse){
1187 AliFatal("This Task needs the PID response attached to the inputHandler");
ddcb10ae 1188 return;
1189 }
1486c432
MW
1190 // count all events looked at here
1191 _eventCount++;
1192
1193 if (_eventAccounting)
1194 {
1195 _eventAccounting->Fill(0);// count all calls to this function
1196 }
1197 else
1198 {
1199 return;
1200 }
1201
1202 _eventAccounting->Fill(1);// count all calls to this function with a valid pointer
1203 //reset single particle counters
1204 k1 = k2 = 0;
1205 __n1_1 = __n1_2 = __s1pt_1 = __s1pt_2 = __n1Nw_1 = __n1Nw_2 = __s1ptNw_1 = __s1ptNw_2 = 0;
1206
1207 float v0Centr = -999.;
1208 float v0ACentr = -999.;
1209 float trkCentr = -999.;
1210 float spdCentr = -999.;
030e4697 1211 //float vertexX = -999;
1212 //float vertexY = -999;
1486c432 1213 float vertexZ = -999;
0da1fab7 1214 //float vertexXY = -999;
1486c432
MW
1215
1216 float centrality = -999;
1217 //Double_t nSigma =-999;
1218 if(fAODEvent)
1219 {
1220
1221 //Centrality
0a918d8d 1222 //AliAODHeader* centralityObject = ((AliVAODHeader*)fAODEvent->GetHeader())->GetCentralityP();
1223 AliCentrality* centralityObject = ((AliVAODHeader*)fAODEvent->GetHeader())->GetCentralityP();
1486c432
MW
1224 if (centralityObject)
1225 {
1226 v0Centr = centralityObject->GetCentralityPercentile("V0M");
1227 v0ACentr = centralityObject->GetCentralityPercentile("V0A");
1228 trkCentr = centralityObject->GetCentralityPercentile("TRK");
1229 spdCentr = centralityObject->GetCentralityPercentile("CL1");
1230
1231 }
1232
1233 _mult4 = v0Centr;
1234 _mult4a = v0ACentr;
1235 _mult5 = trkCentr;
1236 _mult6 = spdCentr;
1237 _field = fAODEvent->GetMagneticField();
1238 //_centralityMethod
1239 switch (_centralityMethod)
1240 {
1241 case 0: centrality = _mult0; break;
1242 case 1: centrality = _mult1; break;
1243 case 2: centrality = _mult2; break;
1244 case 3: centrality = _mult3; break;
1245 case 4: centrality = _mult4; break;
1246 case 5: centrality = _mult5; break;
1247 case 6: centrality = _mult6; break;
1248 case 7: centrality = _mult4a; break;
1249 }
1250
1251 if ( centrality < _centralityMin || centrality > _centralityMax ) return;
1252
1253 _eventAccounting->Fill(2);// count all events with right centrality
1254
1255 vertex = (AliAODVertex*) fAODEvent->GetPrimaryVertex();
1256
1257 if(vertex->GetNContributors() > 0)
1258 {
030e4697 1259 //vertexX = vertex->GetX();
1260 //vertexY = vertex->GetY();
1486c432 1261 vertexZ = vertex->GetZ();
0da1fab7 1262 //vertexXY = sqrt(vertexX*vertexX+vertexY*vertexY);
1486c432
MW
1263 }
1264 if (vertexZ < _vertexZMin || vertexZ > _vertexZMax ) return;
1265
1266 iVertex = int((vertexZ-_min_vertexZ)/_width_vertexZ);
1267 iVertexP1 = iVertex*_nBins_etaPhiPt_1;
1268 iVertexP2 = iVertex*_nBins_etaPhiPt_2;
1269 if (iVertex<0 || iVertex>=_nBins_vertexZ)
1270 {
1271 AliError("AliTaskDptCorrMC::Exec(Option_t *option) iVertex<0 || iVertex>=_nBins_vertexZ ");
1272 return;
1273 }
1274 _eventAccounting->Fill(3);// count all calls to this function with a valid pointer
1275 //======================
9aba9023 1276 //***********************************
1277 //MC AOD Truth
1486c432 1278 if(fAnalysisType == "MCAOD")
9aba9023 1279 {
947cb887 1280 AliMCEvent* mcEvent = MCEvent();
1281 _nTracks = mcEvent->GetNumberOfTracks();
1486c432 1282 _mult3 = _nTracks;
9aba9023 1283 //loop over tracks starts here
1486c432
MW
1284 for (int iTrack=0; iTrack< _nTracks; iTrack++)
1285 {
947cb887 1286 AliAODMCParticle *aodTrack = (AliAODMCParticle*) mcEvent->GetTrack(iTrack);
1486c432 1287
9aba9023 1288 if (!aodTrack)
1486c432
MW
1289 {
1290 AliError(Form("AliTaskDptCorrMC::Exec(Option_t *option) No track ofr iTrack=%d", iTrack));
1291 continue;
9aba9023 1292 }
94694b54 1293
947cb887 1294 if(!aodTrack->IsPhysicalPrimary()) continue;
1486c432 1295
9aba9023 1296 q = aodTrack->Charge();
1486c432 1297 charge = int(q);
9aba9023 1298 phi = aodTrack->Phi();
1299 pt = aodTrack->Pt();
1300 px = aodTrack->Px();
1301 py = aodTrack->Py();
1302 pz = aodTrack->Pz();
1303 eta = aodTrack->Eta();
1304 // Kinematics cuts from ESD track cuts
1305 if( pt < 0.2 || pt > 2.0) continue;
1306 if( eta < -0.8 || eta > 0.8) continue;
947cb887 1307
1308 if(q == 0) continue;
1309 _dcaz->Fill(pt); //AllCh
1310
1311 if(TMath::Abs(aodTrack->GetPdgCode()) == 211)
1312 {
1313 _etadis->Fill(pt); //pion
1314 }
1315 if(TMath::Abs(aodTrack->GetPdgCode()) == 321)
1316 {
1317 _phidis->Fill(pt); //kaon
1318 }
1319
1320 //_etadis->Fill(eta);
1321 //_phidis->Fill(phi);
1486c432 1322
9aba9023 1323 if(fExcludeResonancesInMC)
1324 {
947cb887 1325 //cout<<"***************Prabhat on Weak Decay Particles ************"<<endl;
1326 Int_t gMotherIndex = aodTrack->GetMother();
1327 if(gMotherIndex != -1) {
1328 AliAODMCParticle* motherTrack = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(gMotherIndex));
1329 if(motherTrack) {
1330 Int_t pdgCodeOfMother = motherTrack->GetPdgCode();
1331
1332 if(pdgCodeOfMother == 311 ||
1333 pdgCodeOfMother == -311 ||
1334 pdgCodeOfMother == 310 ||
1335 pdgCodeOfMother == 3122 ||
1336 pdgCodeOfMother == -3122 ||
1337 pdgCodeOfMother == 111 ||
1338 pdgCodeOfMother == 22 ) continue;
1339 }
1340 }
1341 }
9aba9023 1342 //Exclude electrons with PDG
1343
1344 if(fExcludeElectronsInMC) {
1345 if(TMath::Abs(aodTrack->GetPdgCode()) == 11) continue;
1346 }
1347
1348 //Particle loop 1st particle
1349 if (q > 0)
1350 {
1351
1486c432 1352 iPhi = int( phi/_width_phi_1);
9aba9023 1353
1486c432
MW
1354 if (iPhi<0 || iPhi>=_nBins_phi_1 )
1355 {
1356 AliWarning("AliTaskDptCorrMC::analyze() iPhi<0 || iPhi>=_nBins_phi_1");
1357 return;
1358 }
1359
1360 iEta = int((eta-_min_eta_1)/_width_eta_1);
1361 if (iEta<0 || iEta>=_nBins_eta_1)
1362 {
1363 AliWarning(Form("AliTaskDptCorrMC::analyze(AliceEvent * event) Mismatched iEta: %d", iEta));
1364 continue;
1365 }
1366
1367 iPt = int((pt -_min_pt_1 )/_width_pt_1 );
1368 if (iPt<0 || iPt >=_nBins_pt_1)
1369 {
1370 AliWarning(Form("AliTaskDptCorrMC::analyze(AliceEvent * event) Mismatched iPt: %d",iPt));
9aba9023 1371 continue;
1486c432
MW
1372 }
1373 iEtaPhi = iEta*_nBins_phi_1+iPhi;
1374 iZEtaPhiPt = iVertexP1 + iEtaPhi*_nBins_pt_1 + iPt;
1375
1376 if (_correctionWeight_1)
1377 corr = _correctionWeight_1[iZEtaPhiPt];
1378 else
1379 corr = 1;
1380 if (iZEtaPhiPt<0 || iZEtaPhiPt>=_nBins_zEtaPhiPt_1)
1381 {
1382 AliWarning("AliTaskDptCorrMC::analyze(AliceEvent * event) iZEtaPhiPt<0 || iZEtaPhiPt>=_nBins_zEtaPhiPt_1");
1383 continue;
1384 }
1385
1386 if (_singlesOnly)
1387 {
1388
9aba9023 1389 __n1_1_vsPt[iPt] += corr; //cout << "step 15" << endl;
947cb887 1390
9aba9023 1391 __n1_1_vsZEtaPhiPt[iZEtaPhiPt] += corr; //cout << "step 12" << endl;
1486c432
MW
1392
1393 }
1394
9aba9023 1395 else
1486c432
MW
1396 {
1397 corrPt = corr*pt;
1398 _id_1[k1] = iTrack;
1399 _charge_1[k1] = charge;
1400 _iEtaPhi_1[k1] = iEtaPhi;
1401 _iPt_1[k1] = iPt;
1402 _pt_1[k1] = pt;
1403 _px_1[k1] = px;
1404 _py_1[k1] = py;
1405 _pz_1[k1] = pz;
1406 _correction_1[k1] = corr;
1407 __n1_1 += corr;
1408 __n1_1_vsEtaPhi[iEtaPhi] += corr;
1409 __s1pt_1 += corrPt;
1410 __s1pt_1_vsEtaPhi[iEtaPhi] += corrPt;
1411 __n1Nw_1 += 1;
1412 __s1ptNw_1 += pt;
1413 ++k1;
1414 if (k1>=arraySize)
1415 {
1416 AliError(Form("AliTaskDptCorrMC::analyze(AliceEvent * event) k1 >=arraySize; arraySize: %d",arraySize));
1417 return;
1418 }
1419 }
1420 }
1421
9aba9023 1422 if (q < 0)
1486c432
MW
1423 {
1424 iPhi = int( phi/_width_phi_2);
1425
1426 if (iPhi<0 || iPhi>=_nBins_phi_2 )
1427 {
1428 AliWarning("AliTaskDptCorrMC::analyze() iPhi<0 || iPhi>=_nBins_phi_1");
1429 return;
1430 }
1431
1432 iEta = int((eta-_min_eta_2)/_width_eta_2);
1433 if (iEta<0 || iEta>=_nBins_eta_2)
1434 {
1435 AliWarning(Form("AliTaskDptCorrMC::analyze(AliceEvent * event) Mismatched iEta: %d", iEta));
1436 continue;
1437 }
9aba9023 1438
1486c432 1439 iPt = int((pt -_min_pt_2 )/_width_pt_2 );
9aba9023 1440 if (iPt<0 || iPt >=_nBins_pt_2)
1486c432
MW
1441 {
1442 AliWarning(Form("AliTaskDptCorrMC::analyze(AliceEvent * event) Mismatched iPt: %d",iPt));
1443 continue;
1444 }
1445
1446 iEtaPhi = iEta*_nBins_phi_2+iPhi;
1447 iZEtaPhiPt = iVertexP2 + iEtaPhi*_nBins_pt_2 + iPt;
1448 if (iZEtaPhiPt<0 || iZEtaPhiPt>=_nBins_zEtaPhiPt_2)
1449 {
1450 AliWarning("AliTaskDptCorrMC::analyze(AliceEvent * event) iZEtaPhiPt<0 || iZEtaPhiPt>=_nBins_zEtaPhiPt_2");
1451 continue;
1452 }
1453
1486c432
MW
1454 if (_correctionWeight_2)
1455 corr = _correctionWeight_2[iZEtaPhiPt];
1456 else
1457 corr = 1;
9aba9023 1458 //dpt = pt - (charge>0) ? _avgPt_vsEtaPhi_2p[iEtaPhi] : _avgPt_vsEtaPhi_2m[iEtaPhi];
1486c432
MW
1459
1460 if (_singlesOnly)
1461 {
9aba9023 1462 __n1_2_vsPt[iPt] += corr; //cout << "step 15" << endl;
1463 __n1_2_vsZEtaPhiPt[iZEtaPhiPt] += corr; //cout << "step 12" << endl;
1486c432
MW
1464 }
1465 else
1466 {
1467 corrPt = corr*pt;
9aba9023 1468 _id_2[k2] = iTrack; //cout << "step 1" << endl;
1469 _charge_2[k2] = charge; //cout << "step 2" << endl;
1470 _iEtaPhi_2[k2] = iEtaPhi; //cout << "step 3" << endl;
1471 _iPt_2[k2] = iPt; //cout << "step 4" << endl;
1472 _pt_2[k2] = pt; //cout << "step 5" << endl;
1473 _px_2[k2] = px; //cout << "step 6" << endl;
1474 _py_2[k2] = py; //cout << "step 7" << endl;
1475 _pz_2[k2] = pz; //cout << "step 8" << endl;
1476 _correction_2[k2] = corr; //cout << "step 9" << endl;
1477 __n1_2 += corr; //cout << "step 10" << endl;
1478 __s1pt_2 += corrPt; //cout << "step 13" << endl;
1486c432 1479 __n1Nw_2 += 1;
9aba9023 1480 __n1_2_vsEtaPhi[iEtaPhi] += corr; //cout << "step 11" << endl;
1481 __s1pt_2_vsEtaPhi[iEtaPhi] += corrPt; //cout << "step 14" << endl;
1486c432
MW
1482 __s1ptNw_2 += pt;
1483 ++k2;
1484 if (k2>=arraySize)
1485 {
1486 AliWarning(Form("-W- k2 >=arraySize; arraySize: %d",arraySize));
1487 return;
1488 }
1489 }
1490
9aba9023 1491 //cout << "done with track" << endl;
1492 } //iTrack
1493 } //data aod loop
1494 } //MC AOD loop ends
1486c432 1495
9aba9023 1496 //***********************************************
1497 //MC AOD Reconstructed tracks
1498 if(fAnalysisType == "MCAODreco")
1486c432 1499 {
9aba9023 1500
1501 fArrayMC = dynamic_cast<TClonesArray*>(fAODEvent->FindListObject(AliAODMCParticle::StdBranchName()));
1502 if (!fArrayMC) {
1503 AliError("No array of MC particles found !!!");
1504 }
1505
1506 AliMCEvent* mcEvent = MCEvent();
1507 if (!mcEvent) {
1508 AliError("ERROR: Could not retrieve MC event");
1509 }
1510
947cb887 1511 TExMap *trackMap = new TExMap();//Mapping matrix----
1512 //1st loop track for Global tracks
1513 for(Int_t i = 0; i < fAODEvent->GetNumberOfTracks(); i++)
1514 {
1515 AliAODTrack* aodTrack = dynamic_cast<AliAODTrack *>(fAODEvent->GetTrack(i));
1516 if(!aodTrack) {
1517 AliError(Form("ERROR: Could not retrieve AODtrack %d",i));
1518 continue;
1519 }
1520 Int_t gID = aodTrack->GetID();
1521 if (aodTrack->TestFilterBit(1)) trackMap->Add(gID, i);//Global tracks
1522 }
1523
1524 AliAODTrack* newAodTrack;
1525
1526 for (int iTrack=0; iTrack< fAODEvent->GetNumberOfTracks(); iTrack++)
1486c432 1527 {
564a0b8f 1528
1486c432 1529 AliAODTrack *t = dynamic_cast<AliAODTrack *>(fAODEvent->GetTrack(iTrack));
564a0b8f 1530
1486c432
MW
1531 if(!t) {
1532 AliError(Form("ERROR: Could not retrieve AODtrack %d",iTrack));
1533 continue;
1534 }
564a0b8f 1535
1486c432
MW
1536 bitOK = t->TestFilterBit(_trackFilterBit);
1537 if (!bitOK) continue;
947cb887 1538 Int_t gID = t->GetID();
f15c1f69 1539 newAodTrack = gID >= 0 ?t : dynamic_cast<AliAODTrack*>(fAODEvent->GetTrack(trackMap->GetValue(-1-gID)));
1540 if(!newAodTrack) {
1541 AliFatal("Not a standard AOD?");
1542 }
564a0b8f 1543
1486c432
MW
1544 q = t->Charge();
1545 charge = int(q);
1546 phi = t->Phi();
1547 pt = t->Pt();
1548 px = t->Px();
1549 py = t->Py();
1550 pz = t->Pz();
1551 eta = t->Eta();
564a0b8f 1552 //Float_t dcaXY = t->DCA();
fe9d7d9d 1553 //Float_t dcaZ = t->ZAtDCA();
564a0b8f 1554
947cb887 1555 // get the electron nsigma
0f50230b 1556 Double_t nSigma = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(newAodTrack,(AliPID::EParticleType)AliPID::kElectron));
947cb887 1557 Double_t nSigmaPions = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(newAodTrack,(AliPID::EParticleType)AliPID::kPion));
1558 Double_t nSigmaKaons = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(newAodTrack,(AliPID::EParticleType)AliPID::kKaon));
0f50230b 1559 Double_t nSigmaProtons = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(newAodTrack,(AliPID::EParticleType)AliPID::kProton));
947cb887 1560
1561 //Make the decision based on the n-sigma of electrons exclusively
0f50230b 1562 if(nSigma < fNSigmaCut
947cb887 1563 && nSigmaPions > fNSigmaCut
1564 && nSigmaKaons > fNSigmaCut
1565 && nSigmaProtons > fNSigmaCut ) continue;
0f50230b 1566
947cb887 1567
1568 if(q == 0) continue;
564a0b8f 1569 // Kinematics cuts
1570 if( pt < 0.2 || pt > 2.0) continue;
1571 if( eta < _min_eta_1 || eta > _max_eta_1) continue;
1572
0f50230b 1573
1574 /*
1575 _dcaz->Fill(pt);
1576 if (nSigmaPions < fNSigmaCut)
947cb887 1577 {
0f50230b 1578 _etadis->Fill(pt);
947cb887 1579 }
0f50230b 1580 if (nSigmaKaons < fNSigmaCut)
947cb887 1581 {
0f50230b 1582 _phidis->Fill(pt);
1583 }
1584 */ //particle ratio calculation
1585
fe9d7d9d 1586 Double_t pos[3];
947cb887 1587 newAodTrack->GetXYZ(pos);
1588
0f50230b 1589 //Double_t DCAX = pos[0] - vertexX;
1590 //Double_t DCAY = pos[1] - vertexY;
947cb887 1591 Double_t DCAZ = pos[2] - vertexZ;
1592
0f50230b 1593 //Double_t DCAXY = TMath::Sqrt((DCAX*DCAX) + (DCAY*DCAY));
947cb887 1594
0f50230b 1595 /*if (DCAZ < _dcaZMin ||
1596 DCAZ > _dcaZMax ||
1597 DCAXY > _dcaXYMax ) continue;
947cb887 1598 */
947cb887 1599
0f50230b 1600 Int_t label = TMath::Abs(t->GetLabel());
1601 AliAODMCParticle *AODmcTrack = (AliAODMCParticle*) fArrayMC->At(label);
1602 //------------Test of tracks QAs------------
1603 if (AODmcTrack->IsPhysicalPrimary())
1604 {
1605 _dcaz->Fill(DCAZ);
1606 //_dcaxy->Fill(DCAXY);
1607 //_etadis->Fill(eta);
1608 _phidis->Fill(phi);
1609 }
1610 if (AODmcTrack->IsSecondaryFromWeakDecay())
1611 {
1612 _dcaz2->Fill(DCAZ);
1613 //_dcaxy2->Fill(DCAXY);
1614 //_etadis2->Fill(eta);
1615 _phidis2->Fill(phi);
1616 }
fe9d7d9d 1617
0f50230b 1618 if (AODmcTrack->IsSecondaryFromMaterial())
564a0b8f 1619 {
0f50230b 1620 _dcaz3->Fill(DCAZ);
1621 //_dcaxy3->Fill(DCAXY);
1622 //_etadis3->Fill(eta);
1623 _phidis3->Fill(phi);
564a0b8f 1624 }
0f50230b 1625 //---------------------------------
947cb887 1626
0f50230b 1627 /* //W/Wo Secondaries
1628 //if (!AODmcTrack->IsPhysicalPrimary()) continue;
1629 //cout<<"***************Prabhat on Weak Decay Particles ************"<<endl;
1630 if(fExcludeResonancesInMC)
1631 {
947cb887 1632 Int_t label = TMath::Abs(t->GetLabel());
1633 AliAODMCParticle *AODmcTrack = (AliAODMCParticle*) fArrayMC->At(label);
0f50230b 1634
1635 Int_t gMotherIndex = AODmcTrack->GetMother();
1636 if(gMotherIndex != -1) {
1637 AliAODMCParticle* motherTrack = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(gMotherIndex));
1638 if(motherTrack) {
1639 Int_t pdgCodeOfMother = motherTrack->GetPdgCode();
1640
1641 if(pdgCodeOfMother == 311 ||
1642 pdgCodeOfMother == -311 ||
1643 pdgCodeOfMother == 310 ||
1644 pdgCodeOfMother == 3122 ||
1645 pdgCodeOfMother == -3122 ||
1646 pdgCodeOfMother == 111 ||
1647 pdgCodeOfMother == 22 ) continue;
1648 }
1649 }
1650 } */
1651
1652
1653 //Int_t label = TMath::Abs(t->GetLabel());
1654 //AliAODMCParticle *AODmcTrack = (AliAODMCParticle*) fArrayMC->At(label);
947cb887 1655 if (AODmcTrack)
1656 {
1657 if(TMath::Abs(AODmcTrack->GetPdgCode()) == 11) continue;
1658 }
1659
9aba9023 1660 //Particle 1
1661 if (t->Charge() > 0)
1662 {
564a0b8f 1663
1486c432
MW
1664 iPhi = int( phi/_width_phi_1);
1665
1666 if (iPhi<0 || iPhi>=_nBins_phi_1 )
1667 {
1668 AliWarning("AliTaskDptCorrMC::analyze() iPhi<0 || iPhi>=_nBins_phi_1");
1669 return;
1670 }
1671
1672 iEta = int((eta-_min_eta_1)/_width_eta_1);
1673 if (iEta<0 || iEta>=_nBins_eta_1)
1674 {
1675 AliWarning(Form("AliTaskDptCorrMC::analyze(AliceEvent * event) Mismatched iEta: %d", iEta));
1676 continue;
1677 }
1678 iPt = int((pt -_min_pt_1 )/_width_pt_1 );
1679 if (iPt<0 || iPt >=_nBins_pt_1)
1680 {
1681 AliWarning(Form("AliTaskDptCorrMC::analyze(AliceEvent * event) Mismatched iPt: %d",iPt));
1682 continue;
1683 }
1684 iEtaPhi = iEta*_nBins_phi_1+iPhi;
1685 iZEtaPhiPt = iVertexP1 + iEtaPhi*_nBins_pt_1 + iPt;
1686
1687 if (_correctionWeight_1)
1688 corr = _correctionWeight_1[iZEtaPhiPt];
1689 else
1690 corr = 1;
1691 if (iZEtaPhiPt<0 || iZEtaPhiPt>=_nBins_zEtaPhiPt_1)
1692 {
1693 AliWarning("AliTaskDptCorrMC::analyze(AliceEvent * event) iZEtaPhiPt<0 || iZEtaPhiPt>=_nBins_zEtaPhiPt_1");
1694 continue;
1695 }
1696
1697
1698 if (_singlesOnly)
1699 {
1700
1701 __n1_1_vsPt[iPt] += corr; //cout << "step 15" << endl;
947cb887 1702
1486c432
MW
1703 __n1_1_vsZEtaPhiPt[iZEtaPhiPt] += corr; //cout << "step 12" << endl;
1704
1705 }
1706
1707 else
1708 {
1709 corrPt = corr*pt;
1710 _id_1[k1] = iTrack;
1711 _charge_1[k1] = charge;
1712 _iEtaPhi_1[k1] = iEtaPhi;
1713 _iPt_1[k1] = iPt;
1714 _pt_1[k1] = pt;
1715 _px_1[k1] = px;
1716 _py_1[k1] = py;
1717 _pz_1[k1] = pz;
1718 _correction_1[k1] = corr;
1719 __n1_1 += corr;
1720 __n1_1_vsEtaPhi[iEtaPhi] += corr;
1721 __s1pt_1 += corrPt;
1722 __s1pt_1_vsEtaPhi[iEtaPhi] += corrPt;
1723 __n1Nw_1 += 1;
1724 __s1ptNw_1 += pt;
1725 ++k1;
1726 if (k1>=arraySize)
1727 {
1728 AliError(Form("AliTaskDptCorrMC::analyze(AliceEvent * event) k1 >=arraySize; arraySize: %d",arraySize));
1729 return;
1730 }
1731 }
1732 }
1733
9aba9023 1734 if (t->Charge() < 0 )
1486c432
MW
1735 {
1736 //===================================
1737
1738 iPhi = int( phi/_width_phi_2);
1739
1740 if (iPhi<0 || iPhi>=_nBins_phi_2 )
1741 {
1742 AliWarning("AliTaskDptCorrMC::analyze() iPhi<0 || iPhi>=_nBins_phi_1");
1743 return;
1744 }
1745
1746 iEta = int((eta-_min_eta_2)/_width_eta_2);
1747 if (iEta<0 || iEta>=_nBins_eta_2)
1748 {
1749 AliWarning(Form("AliTaskDptCorrMC::analyze(AliceEvent * event) Mismatched iEta: %d", iEta));
1750 continue;
1751 }
1752 iPt = int((pt -_min_pt_2 )/_width_pt_2 );
1753 if (iPt<0 || iPt >=_nBins_pt_2)
1754 {
1755 AliWarning(Form("AliTaskDptCorrMC::analyze(AliceEvent * event) Mismatched iPt: %d",iPt));
1756 continue;
1757 }
1758
1759 iEtaPhi = iEta*_nBins_phi_2+iPhi;
1760 iZEtaPhiPt = iVertexP2 + iEtaPhi*_nBins_pt_2 + iPt;
1761 if (iZEtaPhiPt<0 || iZEtaPhiPt>=_nBins_zEtaPhiPt_2)
1762 {
1763 AliWarning("AliTaskDptCorrMC::analyze(AliceEvent * event) iZEtaPhiPt<0 || iZEtaPhiPt>=_nBins_zEtaPhiPt_2");
1764 continue;
1765 }
1766
1767
1768 if (_correctionWeight_2)
1769 corr = _correctionWeight_2[iZEtaPhiPt];
1770 else
1771 corr = 1;
1772 //dpt = pt - (charge>0) ? _avgPt_vsEtaPhi_2p[iEtaPhi] : _avgPt_vsEtaPhi_2m[iEtaPhi];
1773
1774 if (_singlesOnly)
1775 {
1776 __n1_2_vsPt[iPt] += corr; //cout << "step 15" << endl;
1777 __n1_2_vsZEtaPhiPt[iZEtaPhiPt] += corr; //cout << "step 12" << endl;
1778 }
1779
1780 else
1781 {
1782 corrPt = corr*pt;
1783 _id_2[k2] = iTrack; //cout << "step 1" << endl;
1784 _charge_2[k2] = charge; //cout << "step 2" << endl;
1785 _iEtaPhi_2[k2] = iEtaPhi; //cout << "step 3" << endl;
1786 _iPt_2[k2] = iPt; //cout << "step 4" << endl;
1787 _pt_2[k2] = pt; //cout << "step 5" << endl;
1788 _px_2[k2] = px; //cout << "step 6" << endl;
1789 _py_2[k2] = py; //cout << "step 7" << endl;
1790 _pz_2[k2] = pz; //cout << "step 8" << endl;
1791 _correction_2[k2] = corr; //cout << "step 9" << endl;
1792 __n1_2 += corr; //cout << "step 10" << endl;
1793 __s1pt_2 += corrPt; //cout << "step 13" << endl;
1794 __n1Nw_2 += 1;
1795 __n1_2_vsEtaPhi[iEtaPhi] += corr; //cout << "step 11" << endl;
1796 __s1pt_2_vsEtaPhi[iEtaPhi] += corrPt; //cout << "step 14" << endl;
1797 __s1ptNw_2 += pt;
1798 ++k2;
1799
1800 if (k2>=arraySize)
1801 {
1802 AliWarning(Form("-W- k2 >=arraySize; arraySize: %d",arraySize));
1803 return;
1804 }
1805 }
1806
1807 //cout << "done with track" << endl;
1808 } //iTrack
1809 } //data aod loop
1810 } //MC AOD loop ends
1811
9aba9023 1812 //************************************************
1486c432 1813 } //AOD events
9aba9023 1814
1815
947cb887 1816 // _m0->Fill(_mult0);
1817 //_m1->Fill(_mult1);
1818 //_m2->Fill(_mult2);
1819 //_m3->Fill(_mult3);
1820 //_m4->Fill(_mult4);
1821 //_m5->Fill(_mult5);
1822 //_m6->Fill(_mult6);
1486c432
MW
1823 _vertexZ->Fill(vertexZ);
1824
1825 if (_singlesOnly)
1826 {
1827 // nothing to do here.
1828 }
1829 else
1830 {
1831 if (_sameFilter)
1832 {
1833 _n1_1_vsM->Fill(centrality, __n1_1);
1834 _s1pt_1_vsM->Fill(centrality, __s1pt_1);
1835 _n1Nw_1_vsM->Fill(centrality, __n1Nw_1);
1836 _s1ptNw_1_vsM->Fill(centrality, __s1ptNw_1);
1837 _n1_2_vsM->Fill(centrality, __n1_1);
1838 _s1pt_2_vsM->Fill(centrality, __s1pt_1);
1839 _n1Nw_2_vsM->Fill(centrality, __n1Nw_1);
1840 _s1ptNw_2_vsM->Fill(centrality, __s1ptNw_1);
1841 // reset pair counters
1842 __n2_12 = __s2ptpt_12 = __s2NPt_12 = __s2PtN_12 = 0;
1843 __n2Nw_12 = __s2ptptNw_12 = __s2NPtNw_12 = __s2PtNNw_12 = 0;
1844 if (_field>0)
1845 {
1846 for (int i1=0; i1<k1; i1++)
1847 {
1848 ////cout << " i1:" << i1 << endl;
1849 id_1 = _id_1[i1]; ////cout << " id_1:" << id_1 << endl;
1850 q_1 = _charge_1[i1]; ////cout << " q_1:" << q_1 << endl;
1851 iEtaPhi_1 = _iEtaPhi_1[i1]; ////cout << " iEtaPhi_1:" << iEtaPhi_1 << endl;
1852 iPt_1 = _iPt_1[i1]; ////cout << " iPt_1:" << iPt_1 << endl;
1853 corr_1 = _correction_1[i1]; ////cout << " corr_1:" << corr_1 << endl;
1854 pt_1 = _pt_1[i1]; ////cout << " pt_1:" << pt_1 << endl;
1855
1856
1857 //1 and 2
1858 for (int i2=i1+1; i2<k1; i2++)
1859 {
1860 ////cout << " i2:" << i2 << endl;
1861 id_2 = _id_1[i2]; ////cout << " id_2:" << id_2 << endl;
1862 if (id_1!=id_2)
1863 {
1864 q_2 = _charge_1[i2]; ////cout << " q_1:" << q_1 << endl;
1865 iEtaPhi_2 = _iEtaPhi_1[i2]; ////cout << " iEtaPhi_1:" << iEtaPhi_1 << endl;
1866 iPt_2 = _iPt_1[i2]; ////cout << " iPt_1:" << iPt_1 << endl;
1867 corr_2 = _correction_1[i2]; ////cout << " corr_1:" << corr_1 << endl;
1868 pt_2 = _pt_1[i2]; ////cout << " pt_1:" << pt_1 << endl;
1869
1870 corr = corr_1*corr_2;
1871 if (q_2>q_1 || (q_1>0 && q_2>0 && pt_2<=pt_1) || (q_1<0 && q_2<0 && pt_2>=pt_1))
1872 {
1873 ij = iEtaPhi_1*_nBins_etaPhi_1 + iEtaPhi_2; ////cout << " ij:" << ij<< endl;
1874 }
1875 else // swap particles
1876 {
1877 ij = iEtaPhi_2*_nBins_etaPhi_1 + iEtaPhi_1; ////cout << " ij:" << ij<< endl;
1878 }
1879
1880 __n2_12 += corr;
1881 __n2_12_vsEtaPhi[ij] += corr;
1882 ptpt = pt_1*pt_2;
1883 __s2ptpt_12 += corr*ptpt;
1884 __s2PtN_12 += corr*pt_1;
1885 __s2NPt_12 += corr*pt_2;
1886 __s2ptpt_12_vsEtaPhi[ij] += corr*ptpt;
1887 __s2PtN_12_vsEtaPhi[ij] += corr*pt_1;
1888 __s2NPt_12_vsEtaPhi[ij] += corr*pt_2;
1889 __n2_12_vsPtPt[iPt_1*_nBins_pt_2 + iPt_2] += corr;
1890
1891 __n2Nw_12 += 1;
1892 __s2ptptNw_12 += ptpt;
1893 __s2PtNNw_12 += pt_1;
1894 __s2NPtNw_12 += pt_2;
1895
1896 }
1897 } //i2
1898 } //i1
1899 }
1900
1901 else // field<0
1902 {
1903 for (int i1=0; i1<k1; i1++)
1904 {
1905 ////cout << " i1:" << i1 << endl;
1906 id_1 = _id_1[i1]; ////cout << " id_1:" << id_1 << endl;
1907 q_1 = _charge_1[i1]; ////cout << " q_1:" << q_1 << endl;
1908 iEtaPhi_1 = _iEtaPhi_1[i1]; ////cout << " iEtaPhi_1:" << iEtaPhi_1 << endl;
1909 iPt_1 = _iPt_1[i1]; ////cout << " iPt_1:" << iPt_1 << endl;
1910 corr_1 = _correction_1[i1]; ////cout << " corr_1:" << corr_1 << endl;
1911 pt_1 = _pt_1[i1]; ////cout << " pt_1:" << pt_1 << endl;
1912
1913 //1 and 2
1914 for (int i2=i1+1; i2<k1; i2++)
1915 {
1916 ////cout << " i2:" << i2 << endl;
1917 id_2 = _id_1[i2]; ////cout << " id_2:" << id_2 << endl;
1918 if (id_1!=id_2)
1919 {
1920 q_2 = _charge_1[i2]; ////cout << " q_2:" << q_2 << endl;
1921 iEtaPhi_2 = _iEtaPhi_1[i2]; ////cout << " iEtaPhi_2:" << iEtaPhi_2 << endl;
1922 iPt_2 = _iPt_1[i2]; ////cout << " iPt_2:" << iPt_2 << endl;
1923 corr_2 = _correction_1[i2]; ////cout << " corr_2:" << corr_2 << endl;
1924 pt_2 = _pt_1[i2]; ////cout << " pt_2:" << pt_2 << endl;
1925 corr = corr_1*corr_2;
1926 if ( q_2<q_1 || (q_1>0 && q_2>0 && pt_2>=pt_1) || (q_1<0 && q_2<0 && pt_2<=pt_1))
1927 {
1928 ij = iEtaPhi_1*_nBins_etaPhi_1 + iEtaPhi_2; ////cout << " ij:" << ij<< endl;
1929 }
1930 else // swap particles
1931 {
1932 ij = iEtaPhi_2*_nBins_etaPhi_1 + iEtaPhi_1; ////cout << " ij:" << ij<< endl;
1933 }
1934
1935 __n2_12 += corr;
1936 __n2_12_vsEtaPhi[ij] += corr;
1937 ptpt = pt_1*pt_2;
1938 __s2ptpt_12 += corr*ptpt;
1939 __s2PtN_12 += corr*pt_1;
1940 __s2NPt_12 += corr*pt_2;
1941 __s2ptpt_12_vsEtaPhi[ij] += corr*ptpt;
1942 __s2PtN_12_vsEtaPhi[ij] += corr*pt_1;
1943 __s2NPt_12_vsEtaPhi[ij] += corr*pt_2;
1944 __n2_12_vsPtPt[iPt_1*_nBins_pt_2 + iPt_2] += corr;
1945
1946 __n2Nw_12 += 1;
1947 __s2ptptNw_12 += ptpt;
1948 __s2PtNNw_12 += pt_1;
1949 __s2NPtNw_12 += pt_2;
1950
1951 }
1952 } //i2
1953 } //i1
1954 }
1955 }
1956 else // filter 1 and 2 are different -- must do all particle pairs...
1957 {
1958 _n1_1_vsM->Fill(centrality, __n1_1);
1959 _s1pt_1_vsM->Fill(centrality, __s1pt_1);
1960 _n1Nw_1_vsM->Fill(centrality, __n1Nw_1);
1961 _s1ptNw_1_vsM->Fill(centrality, __s1ptNw_1);
1962 _n1_2_vsM->Fill(centrality, __n1_2);
1963 _s1pt_2_vsM->Fill(centrality, __s1pt_2);
1964 _n1Nw_2_vsM->Fill(centrality, __n1Nw_2);
1965 _s1ptNw_2_vsM->Fill(centrality, __s1ptNw_2);
1966 // reset pair counters
1967 __n2_12 = __s2ptpt_12 = __s2NPt_12 = __s2PtN_12 = 0;
1968 __n2Nw_12 = __s2ptptNw_12 = __s2NPtNw_12 = __s2PtNNw_12 = 0;
1969 for (int i1=0; i1<k1; i1++)
1970 {
1971 ////cout << " i1:" << i1 << endl;
1972 id_1 = _id_1[i1]; ////cout << " id_1:" << id_1 << endl;
1973 q_1 = _charge_1[i1]; ////cout << " q_1:" << q_1 << endl;
1974 iEtaPhi_1 = _iEtaPhi_1[i1]; ////cout << " iEtaPhi_1:" << iEtaPhi_1 << endl;
1975 iPt_1 = _iPt_1[i1]; ////cout << " iPt_1:" << iPt_1 << endl;
1976 corr_1 = _correction_1[i1]; ////cout << " corr_1:" << corr_1 << endl;
1977 pt_1 = _pt_1[i1]; ////cout << " pt_1:" << pt_1 << endl;
1978 px_1 = _px_1[i1]; ////cout << " px_1:" << px_1 << endl;
1979 py_1 = _py_1[i1]; ////cout << " py_1:" << py_1 << endl;
1980 pz_1 = _pz_1[i1]; ////cout << " pz_1:" << pz_1 << endl;
1981 //1 and 2
1982 for (int i2=0; i2<k2; i2++)
1983 {
1984 ////cout << " i2:" << i2 << endl;
1985 id_2 = _id_2[i2]; ////cout << " id_2:" << id_2 << endl;
1986 if (id_1!=id_2) // exclude auto correlation
1987 {
1988 q_2 = _charge_2[i2]; ////cout << " q_2:" << q_2 << endl;
1989 iEtaPhi_2 = _iEtaPhi_2[i2]; ////cout << " iEtaPhi_2:" << iEtaPhi_2 << endl;
1990 iPt_2 = _iPt_2[i2]; ////cout << " iPt_2:" << iPt_2 << endl;
1991 corr_2 = _correction_2[i2]; ////cout << " corr_2:" << corr_2 << endl;
1992 pt_2 = _pt_2[i2]; ////cout << " pt_2:" << pt_2 << endl;
1993 px_2 = _px_2[i2]; ////cout << " px_2:" << px_2 << endl;
1994 py_2 = _py_2[i2]; ////cout << " py_2:" << py_2 << endl;
1995 pz_2 = _pz_2[i2]; ////cout << " pz_2:" << pz_2 << endl;
1996
1997
1998 if (_rejectPairConversion)
1999 {
2000 float e1Sq = massElecSq + pt_1*pt_1 + pz_1*pz_1;
2001 float e2Sq = massElecSq + pt_2*pt_2 + pz_2*pz_2;
2002 float mInvSq = 2*(massElecSq + sqrt(e1Sq*e2Sq) - px_1*px_2 - py_1*py_2 - pz_1*pz_2 );
2003 float mInv = sqrt(mInvSq);
2004 _invMass->Fill(mInv);
2005 }
2006
2007 corr = corr_1*corr_2;
2008 ij = iEtaPhi_1*_nBins_etaPhi_1 + iEtaPhi_2; ////cout << " ij:" << ij<< endl;
2009 __n2_12 += corr;
2010 __n2_12_vsEtaPhi[ij] += corr;
2011 ptpt = pt_1*pt_2;
2012 __s2ptpt_12 += corr*ptpt;
2013 __s2PtN_12 += corr*pt_1;
2014 __s2NPt_12 += corr*pt_2;
2015 __s2ptpt_12_vsEtaPhi[ij] += corr*ptpt;
2016 __s2PtN_12_vsEtaPhi[ij] += corr*pt_1;
2017 __s2NPt_12_vsEtaPhi[ij] += corr*pt_2;
2018 __n2_12_vsPtPt[iPt_1*_nBins_pt_2 + iPt_2] += corr;
2019 __n2Nw_12 += 1;
2020 __s2ptptNw_12 += ptpt;
2021 __s2PtNNw_12 += pt_1;
2022 __s2NPtNw_12 += pt_2;
2023
2024 }
2025 } //i2
2026 } //i1
2027 }
2028
2029 _n2_12_vsM->Fill(centrality, __n2_12);
2030 _s2PtPt_12_vsM->Fill(centrality, __s2ptpt_12);
2031 _s2PtN_12_vsM->Fill(centrality, __s2NPt_12);
2032 _s2NPt_12_vsM->Fill(centrality, __s2PtN_12);
2033
2034 _n2Nw_12_vsM->Fill(centrality, __n2Nw_12);
2035 _s2PtPtNw_12_vsM->Fill(centrality, __s2ptptNw_12);
2036 _s2PtNNw_12_vsM->Fill(centrality, __s2NPtNw_12);
2037 _s2NPtNw_12_vsM->Fill(centrality, __s2PtNNw_12);
2038
2039 }
2040
2041
2042 AliInfo("AliTaskDptCorrMC::UserExec() -----------------Event Done ");
2043 PostData(0,_outputHistoList);
2044
2045} //End of UserExec
2046
2047
2048
2049
2050void AliDptDptInMC::FinishTaskOutput()
2051{
2052 AliInfo("AliDptDptInMC::FinishTaskOutput() Starting.");
2053 Printf("= 0 ====================================================================");
2054 finalizeHistograms();
2055 AliInfo("= 1 ====================================================================");
2056 PostData(0,_outputHistoList);
2057 AliInfo("= 2 ====================================================================");
2058 AliInfo("AliDptDptInMC::FinishTaskOutput() Done.");
2059}
2060
2061void AliDptDptInMC::Terminate(Option_t* /*option*/)
2062{
2063 AliInfo("AliDptDptInMC::Terminate() Starting/Done.");
2064}
2065
2066
2067//Tools
2068//===================================================================================================
2069void AliDptDptInMC::fillHistoWithArray(TH1 * h, float * array, int size)
2070{
2071 int i, i1;
2072 float v1, ev1, v2, ev2, sum, esum;
2073 for (i=0, i1=1; i<size; ++i,++i1)
2074 {
2075 v1 = array[i]; ev1 = sqrt(v1);
2076 v2 = h->GetBinContent(i1);
2077 ev2 = h->GetBinError(i1);
2078 sum = v1 + v2;
2079 esum = sqrt(ev1*ev1+ev2*ev2);
2080 h->SetBinContent(i1,sum);
2081 h->SetBinError(i1,esum);
2082 }
2083}
2084
2085void AliDptDptInMC::fillHistoWithArray(TH2 * h, float * array, int size1, int size2)
2086{
2087 int i, i1;
2088 int j, j1;
2089 float v1, ev1, v2, ev2, sum, esum;
2090 for (i=0, i1=1; i<size1; ++i,++i1)
2091 {
2092 for (j=0, j1=1; j<size2; ++j,++j1)
2093 {
2094 v1 = array[i*size2+j]; ev1 = sqrt(v1);
2095 v2 = h->GetBinContent(i1,j1);
2096 ev2 = h->GetBinError(i1,j1);
2097 sum = v1 + v2;
2098 esum = sqrt(ev1*ev1+ev2*ev2);
2099 h->SetBinContent(i1,j1,sum);
2100 h->SetBinError(i1,j1,esum);
2101 }
2102 }
2103}
2104
2105void AliDptDptInMC::fillHistoWithArray(TH3 * h, float * array, int size1, int size2, int size3)
2106{
2107 int i, i1;
2108 int j, j1;
2109 int k, k1;
2110 float v1, ev1, v2, ev2, sum, esum;
2111 int size23 = size2*size3;
2112 for (i=0, i1=1; i<size1; ++i,++i1)
2113 {
2114 for (j=0, j1=1; j<size2; ++j,++j1)
2115 {
2116 for (k=0, k1=1; k<size3; ++k,++k1)
2117 {
2118 v1 = array[i*size23+j*size3+k]; ev1 = sqrt(v1);
2119 v2 = h->GetBinContent(i1,j1,k1);
2120 ev2 = h->GetBinError(i1,j1,k1);
2121 sum = v1 + v2;
2122 esum = sqrt(ev1*ev1+ev2*ev2);
2123 h->SetBinContent(i1,j1,k1,sum);
2124 h->SetBinError(i1,j1,k1,esum);
2125 }
2126 }
2127 }
2128}
2129
2130void AliDptDptInMC::fillHistoWithArray(TH1 * h, double * array, int size)
2131{
2132 int i, i1;
2133 double v1, ev1, v2, ev2, sum, esum;
2134 for (i=0, i1=1; i<size; ++i,++i1)
2135 {
2136 v1 = array[i]; ev1 = sqrt(v1);
2137 v2 = h->GetBinContent(i1);
2138 ev2 = h->GetBinError(i1);
2139 sum = v1 + v2;
2140 esum = sqrt(ev1*ev1+ev2*ev2);
2141 h->SetBinContent(i1,sum);
2142 h->SetBinError(i1,esum);
2143 }
2144}
2145
2146void AliDptDptInMC::fillHistoWithArray(TH2 * h, double * array, int size1, int size2)
2147{
2148 int i, i1;
2149 int j, j1;
2150 double v1, ev1, v2, ev2, sum, esum;
2151 for (i=0, i1=1; i<size1; ++i,++i1)
2152 {
2153 for (j=0, j1=1; j<size2; ++j,++j1)
2154 {
2155 v1 = array[i*size2+j]; ev1 = sqrt(v1);
2156 v2 = h->GetBinContent(i1,j1);
2157 ev2 = h->GetBinError(i1,j1);
2158 sum = v1 + v2;
2159 esum = sqrt(ev1*ev1+ev2*ev2);
2160 h->SetBinContent(i1,j1,sum);
2161 h->SetBinError(i1,j1,esum);
2162 }
2163 }
2164}
2165
2166void AliDptDptInMC::fillHistoWithArray(TH3 * h, double * array, int size1, int size2, int size3)
2167{
2168 int i, i1;
2169 int j, j1;
2170 int k, k1;
2171 double v1, ev1, v2, ev2, sum, esum;
2172 int size23 = size2*size3;
2173 for (i=0, i1=1; i<size1; ++i,++i1)
2174 {
2175 for (j=0, j1=1; j<size2; ++j,++j1)
2176 {
2177 for (k=0, k1=1; k<size3; ++k,++k1)
2178 {
2179 v1 = array[i*size23+j*size3+k]; ev1 = sqrt(v1);
2180 v2 = h->GetBinContent(i1,j1,k1);
2181 ev2 = h->GetBinError(i1,j1,k1);
2182 sum = v1 + v2;
2183 esum = sqrt(ev1*ev1+ev2*ev2);
2184 h->SetBinContent(i1,j1,k1,sum);
2185 h->SetBinError(i1,j1,k1,esum);
2186 }
2187 }
2188 }
2189}
2190
2191//________________________________________________________________________
2192double * AliDptDptInMC::getDoubleArray(int size, double v)
2193{
2194 /// Allocate an array of type double with n values
2195 /// Initialize the array to the given value
2196 double * array = new double [size];
2197 for (int i=0;i<size;++i) array[i]=v;
2198 return array;
2199}
2200
2201//________________________________________________________________________
2202float * AliDptDptInMC::getFloatArray(int size, float v)
2203{
2204 /// Allocate an array of type float with n values
2205 /// Initialize the array to the given value
2206 float * array = new float [size];
2207 for (int i=0;i<size;++i) array[i]=v;
2208 return array;
2209}
2210
2211
2212//________________________________________________________________________
2213TH1D * AliDptDptInMC::createHisto1D(const TString & name, const TString & title,
2214 int n, double xMin, double xMax,
2215 const TString & xTitle, const TString & yTitle)
2216{
2217 //CreateHisto new 1D historgram
2218 AliInfo(Form("createHisto 1D histo %s nBins: %d xMin: %f xMax: %f",name.Data(),n,xMin,xMax));
2219 TH1D * h = new TH1D(name,title,n,xMin,xMax);
2220 h->GetXaxis()->SetTitle(xTitle);
2221 h->GetYaxis()->SetTitle(yTitle);
2222 addToList(h);
2223 return h;
2224}
2225
2226
2227//________________________________________________________________________
2228TH1D * AliDptDptInMC::createHisto1D(const TString & name, const TString & title,
2229 int n, double * bins,
2230 const TString & xTitle, const TString & yTitle)
2231{
2232 AliInfo(Form("createHisto 1D histo %s with %d non uniform nBins",name.Data(),n));
2233 TH1D * h = new TH1D(name,title,n,bins);
2234 h->GetXaxis()->SetTitle(xTitle);
2235 h->GetYaxis()->SetTitle(yTitle);
2236 addToList(h);
2237 return h;
2238}
2239
2240
2241//________________________________________________________________________
2242TH2D * AliDptDptInMC::createHisto2D(const TString & name, const TString & title,
2243 int nx, double xMin, double xMax, int ny, double yMin, double yMax,
2244 const TString & xTitle, const TString & yTitle, const TString & zTitle)
2245{
2246 AliInfo(Form("createHisto 2D histo %s nx: %d xMin: %f10.4 xMax: %f10.4 ny: %d yMin: %f10.4 yMax: %f10.4",name.Data(),nx,xMin,xMax,ny,yMin,yMax));
2247 TH2D * h = new TH2D(name,title,nx,xMin,xMax,ny,yMin,yMax);
2248 h->GetXaxis()->SetTitle(xTitle);
2249 h->GetYaxis()->SetTitle(yTitle);
2250 h->GetZaxis()->SetTitle(zTitle);
2251 addToList(h);
2252 return h;
2253}
2254
2255//________________________________________________________________________
2256TH2D * AliDptDptInMC::createHisto2D(const TString & name, const TString & title,
2257 int nx, double* xbins, int ny, double yMin, double yMax,
2258 const TString & xTitle, const TString & yTitle, const TString & zTitle)
2259{
2260 AliInfo(Form("createHisto 2D histo %s with %d non uniform nBins",name.Data(),nx));
2261 TH2D * h;
2262 h = new TH2D(name,title,nx,xbins,ny,yMin,yMax);
2263 h->GetXaxis()->SetTitle(xTitle);
2264 h->GetYaxis()->SetTitle(yTitle);
2265 h->GetZaxis()->SetTitle(zTitle);
2266 addToList(h);
2267 return h;
2268}
2269
2270//// F /////
2271//________________________________________________________________________
2272TH1F * AliDptDptInMC::createHisto1F(const TString & name, const TString & title,
2273 int n, double xMin, double xMax,
2274 const TString & xTitle, const TString & yTitle)
2275{
2276 //CreateHisto new 1D historgram
2277 AliInfo(Form("createHisto 1D histo %s nBins: %d xMin: %f xMax: %f",name.Data(),n,xMin,xMax));
2278 TH1F * h = new TH1F(name,title,n,xMin,xMax);
2279 h->GetXaxis()->SetTitle(xTitle);
2280 h->GetYaxis()->SetTitle(yTitle);
2281 addToList(h);
2282 return h;
2283}
2284
2285
2286//________________________________________________________________________
2287TH1F * AliDptDptInMC::createHisto1F(const TString & name, const TString & title,
2288 int n, double * bins,
2289 const TString & xTitle, const TString & yTitle)
2290{
2291 AliInfo(Form("createHisto 1D histo %s with %d non uniform nBins",name.Data(),n));
2292 TH1F * h = new TH1F(name,title,n,bins);
2293 h->GetXaxis()->SetTitle(xTitle);
2294 h->GetYaxis()->SetTitle(yTitle);
2295 addToList(h);
2296 return h;
2297}
2298
2299
2300//________________________________________________________________________
2301TH2F * AliDptDptInMC::createHisto2F(const TString & name, const TString & title,
2302 int nx, double xMin, double xMax, int ny, double yMin, double yMax,
2303 const TString & xTitle, const TString & yTitle, const TString & zTitle)
2304{
2305 AliInfo(Form("createHisto 2D histo %s nx: %d xMin: %f10.4 xMax: %f10.4 ny: %d yMin: %f10.4 yMax: %f10.4",name.Data(),nx,xMin,xMax,ny,yMin,yMax));
2306 TH2F * h = new TH2F(name,title,nx,xMin,xMax,ny,yMin,yMax);
2307 h->GetXaxis()->SetTitle(xTitle);
2308 h->GetYaxis()->SetTitle(yTitle);
2309 h->GetZaxis()->SetTitle(zTitle);
2310 addToList(h);
2311 return h;
2312}
2313
2314//________________________________________________________________________
2315TH2F * AliDptDptInMC::createHisto2F(const TString & name, const TString & title,
2316 int nx, double* xbins, int ny, double yMin, double yMax,
2317 const TString & xTitle, const TString & yTitle, const TString & zTitle)
2318{
2319 AliInfo(Form("createHisto 2D histo %s with %d non uniform nBins",name.Data(),nx));
2320 TH2F * h;
2321 h = new TH2F(name,title,nx,xbins,ny,yMin,yMax);
2322 h->GetXaxis()->SetTitle(xTitle);
2323 h->GetYaxis()->SetTitle(yTitle);
2324 h->GetZaxis()->SetTitle(zTitle);
2325 addToList(h);
2326 return h;
2327}
2328
2329
2330//________________________________________________________________________
2331TH3F * AliDptDptInMC::createHisto3F(const TString & name, const TString & title,
2332 int nx, double xMin, double xMax,
2333 int ny, double yMin, double yMax,
2334 int nz, double zMin, double zMax,
2335 const TString & xTitle, const TString & yTitle, const TString & zTitle)
2336{
2337 AliInfo(Form("createHisto 2D histo %s nx: %d xMin: %f10.4 xMax: %f10.4 ny: %d yMin: %f10.4 yMax: %f10.4 nz: %d zMin: %f10.4 zMax: %f10.4",name.Data(),nx,xMin,xMax,ny,yMin,yMax,nz,zMin,zMax));
2338 TH3F * h = new TH3F(name,title,nx,xMin,xMax,ny,yMin,yMax,nz,zMin,zMax);
2339 h->GetXaxis()->SetTitle(xTitle);
2340 h->GetYaxis()->SetTitle(yTitle);
2341 h->GetZaxis()->SetTitle(zTitle);
2342 addToList(h);
2343 return h;
2344}
2345
2346
2347//________________________________________________________________________
2348TProfile * AliDptDptInMC::createProfile(const TString & name, const TString & description,
2349 int nx,double xMin,double xMax,
2350 const TString & xTitle, const TString & yTitle)
2351{
2352 AliInfo(Form("createHisto 1D profile %s nBins: %d xMin: %f10.4 xMax: %f10.4",name.Data(),nx,xMin,xMax));
2353 TProfile * h = new TProfile(name,description,nx,xMin,xMax);
2354 h->GetXaxis()->SetTitle(xTitle);
2355 h->GetYaxis()->SetTitle(yTitle);
2356 addToList(h);
2357 return h;
2358}
2359
2360//________________________________________________________________________
2361TProfile * AliDptDptInMC::createProfile(const TString & name,const TString & description,
2362 int nx, double* bins,
2363 const TString & xTitle, const TString & yTitle)
2364{
2365 AliInfo(Form("createHisto 1D profile %s with %d non uniform bins",name.Data(),nx));
2366 TProfile * h = new TProfile(name,description,nx,bins);
2367 h->GetXaxis()->SetTitle(xTitle);
2368 h->GetYaxis()->SetTitle(yTitle);
2369 addToList(h);
2370 return h;
2371}
2372
2373
2374void AliDptDptInMC::addToList(TH1 *h)
2375{
2376 if (_outputHistoList)
2377 {
2378 _outputHistoList->Add(h);
2379 }
2380 else
2381 AliInfo("addToList(TH1 *h) _outputHistoList is null!!!!! Should abort ship");
2382
2383}
2384
2385
2386