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