Compilation with Root6
[u/mrichter/AliRoot.git] / PWGCF / Correlations / DPhi / AliAnalysisTaskDptDptQA.cxx
CommitLineData
c5704d5c 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:$ */
349d46ce 17
18#include <TROOT.h>
19#include <TChain.h>
20#include <TFile.h>
21#include <TList.h>
22#include <TMath.h>
23#include <TTree.h>
24#include <TH1F.h>
25#include <TH2F.h>
26#include <TH3F.h>
27#include <TProfile.h>
28#include <TH1D.h>
29#include <TH2D.h>
30#include <TH3D.h>
31#include <TRandom.h>
32#include "AliAnalysisManager.h"
33
34#include "AliAODHandler.h"
35#include "AliAODInputHandler.h"
c5704d5c 36//#include "AliAODMCParticle.h"
349d46ce 37#include "AliInputEventHandler.h"
38#include "AliLog.h"
39#include "AliESDEvent.h"
40#include "AliESDInputHandler.h"
41#include "AliMultiplicity.h"
42#include "AliCentrality.h"
43#include "AliAnalysisTaskDptDptQA.h"
44
e742b9b0 45using std::cout;
46using std::endl;
47
349d46ce 48ClassImp(AliAnalysisTaskDptDptQA)
49
50AliAnalysisTaskDptDptQA::AliAnalysisTaskDptDptQA()
51: AliAnalysisTaskSE(),
52fAODEvent(0),
53fESDEvent(0), //! ESD Event
54fInputHandler(0),
349d46ce 55_outputHistoList(0),
56_twoPi ( 2.0 * 3.1415927),
57_eventCount ( 0),
58_debugLevel ( 0),
59_singlesOnly ( 0),
60_useWeights ( 0),
61_sameFilter ( false),
62_rejectPileup ( 1),
63_rejectPairConversion ( 0),
64_vertexZMin ( -10),
65_vertexZMax ( 10),
66_vertexXYMin ( -10),
67_vertexXYMax ( 10),
68_centralityMethod ( 4),
69_centralityMin ( 0.),
70_centralityMax ( 0.),
71_requestedCharge_1 ( 1),
72_requestedCharge_2 ( -1),
73_dcaZMin ( -3),
74_dcaZMax ( 3.),
c5704d5c 75_dcaXYMin ( -3.),
76_dcaXYMax ( 3.),
349d46ce 77_dedxMin ( 0),
78_dedxMax ( 100000),
c5704d5c 79_nClusterMin ( 70),
80_trackFilterBit ( 128),
349d46ce 81_field ( 1.),
82_nTracks ( 0 ),
83_mult0 ( 0 ),
84_mult1 ( 0 ),
85_mult2 ( 0 ),
86_mult3 ( 0 ),
87_mult4 ( 0 ),
349d46ce 88_mult5 ( 0 ),
89_mult6 ( 0 ),
90arraySize ( 2000),
91_id_1(0),
92_charge_1(0),
93_iEtaPhi_1(0),
94_iPt_1(0),
95_pt_1(0),
96_px_1(0),
97_py_1(0),
98_pz_1(0),
99_correction_1(0),
100_dedx_1(0),
101_id_2(0),
102_charge_2(0),
103_iEtaPhi_2(0),
104_iPt_2(0),
105_pt_2(0),
106_px_2(0),
107_py_2(0),
108_pz_2(0),
109_correction_2(0),
110_dedx_2(0),
111_correctionWeight_1(0),
112_correctionWeight_2(0),
113_nBins_M0(500), _min_M0(0), _max_M0(10000), _width_M0(20),
114_nBins_M1(500), _min_M1(0), _max_M1(10000), _width_M1(20),
115_nBins_M2(500), _min_M2(0), _max_M2(10000), _width_M2(20),
116_nBins_M3(500), _min_M3(0), _max_M3(10000), _width_M3(20),
117_nBins_M4(100), _min_M4(0), _max_M4(1), _width_M4(0.01),
118_nBins_M5(100), _min_M5(0), _max_M5(1), _width_M5(0.01),
119_nBins_M6(100), _min_M6(0), _max_M6(1), _width_M6(0.01),
120_nBins_vertexZ(40), _min_vertexZ(-10), _max_vertexZ(10), _width_vertexZ(0.5),
121
122_nBins_pt_1(18), _min_pt_1(0.2), _max_pt_1(2.0), _width_pt_1(0.1),
123_nBins_phi_1(72), _min_phi_1(0), _max_phi_1(2.*3.1415927),_width_phi_1(2.*3.1415927/72.),
c5704d5c 124_nBins_eta_1(0), _min_eta_1(0), _max_eta_1(0), _width_eta_1(0.1),
349d46ce 125
126_nBins_etaPhi_1(0),
127_nBins_etaPhiPt_1(0),
128_nBins_zEtaPhiPt_1(0),
129
130_nBins_pt_2(18), _min_pt_2(0.2), _max_pt_2(2.0), _width_pt_2(0.1),
131_nBins_phi_2(72), _min_phi_2(0), _max_phi_2(2.*3.1415927),_width_phi_2(2.*3.1415927/72),
c5704d5c 132_nBins_eta_2(0), _min_eta_2(0), _max_eta_2(0), _width_eta_2(0.1),
349d46ce 133
134_nBins_etaPhi_2(0),
135_nBins_etaPhiPt_2(0),
136_nBins_zEtaPhiPt_2(0),
137_nBins_etaPhi_12(0),
138__n1_1(0),
139__n1_2(0),
140__n2_12(0),
141__s1pt_1(0),
142__s1pt_2(0),
143__s2ptpt_12(0),
144__s2NPt_12(0),
145__s2PtN_12(0),
146__n1Nw_1(0),
147__n1Nw_2(0),
148__n2Nw_12(0),
149__s1ptNw_1(0),
150__s1ptNw_2(0),
151__s2ptptNw_12(0),
152__s2NPtNw_12(0),
153__s2PtNNw_12(0),
154__n1_1_vsPt(0),
155__n1_1_vsEtaPhi(0),
156__s1pt_1_vsEtaPhi(0),
157__n1_1_vsZEtaPhiPt(0),
158__n1_2_vsPt(0),
159__n1_2_vsEtaPhi(0),
160__s1pt_2_vsEtaPhi(0),
161__n1_2_vsZEtaPhiPt(0),
162__n2_12_vsPtPt(0),
163__n2_12_vsEtaPhi(0),
164__s2ptpt_12_vsEtaPhi(0),
165__s2PtN_12_vsEtaPhi(0),
166__s2NPt_12_vsEtaPhi(0),
167_weight_1 ( 0 ),
168_weight_2 ( 0 ),
169_eventAccounting ( 0),
170_m0 ( 0),
171_m1 ( 0),
172_m2 ( 0),
173_m3 ( 0),
174_m4 ( 0),
175_m5 ( 0),
176_m6 ( 0),
c5704d5c 177_vertexZ ( 0),
349d46ce 178_etadis ( 0),
179_phidis ( 0),
180_dcaz ( 0),
181_dcaxy ( 0),
349d46ce 182_n1_1_vsPt ( 0),
183_n1_1_vsEtaVsPhi ( 0),
184_s1pt_1_vsEtaVsPhi ( 0),
185_n1_1_vsZVsEtaVsPhiVsPt ( 0),
c5704d5c 186_n1_1_vsM ( 0), // w/ weight
349d46ce 187_s1pt_1_vsM ( 0),
c5704d5c 188_n1Nw_1_vsM ( 0), // w/o weight
349d46ce 189_s1ptNw_1_vsM ( 0),
190_dedxVsP_1 ( 0),
191_corrDedxVsP_1 ( 0),
192_betaVsP_1 ( 0),
193_n1_2_vsPt ( 0),
194_n1_2_vsEtaVsPhi ( 0),
195_s1pt_2_vsEtaVsPhi ( 0),
196_n1_2_vsZVsEtaVsPhiVsPt ( 0),
197_n1_2_vsM ( 0),
198_s1pt_2_vsM ( 0),
199_n1Nw_2_vsM ( 0),
200_s1ptNw_2_vsM ( 0),
201_dedxVsP_2 ( 0),
202_corrDedxVsP_2 ( 0),
203_betaVsP_2 ( 0),
204_n2_12_vsEtaPhi ( 0),
205_n2_12_vsPtVsPt ( 0),
206_s2PtPt_12_vsEtaPhi( 0),
207_s2PtN_12_vsEtaPhi ( 0),
208_s2NPt_12_vsEtaPhi ( 0),
209_n2_12_vsM ( 0),
210_s2PtPt_12_vsM ( 0),
211_s2PtN_12_vsM ( 0),
212_s2NPt_12_vsM ( 0),
213_n2Nw_12_vsM ( 0),
214_s2PtPtNw_12_vsM ( 0),
215_s2PtNNw_12_vsM ( 0),
216_s2NPtNw_12_vsM ( 0),
217_invMass ( 0),
c5704d5c 218_invMassElec ( 0),
349d46ce 219n1Name("NA"),
220n1NwName("NA"),
221n2Name("NA"),
222n2NwName("NA"),
223n3Name("NA"),
224n1n1Name("NA"),
225n1n1n1Name("NA"),
226n2n1Name("NA"),
227r1Name("NA"),
228r2Name("NA"),
229r3Name("NA"),
230r2r1Name("NA"),
231c2Name("NA"),
232c3Name("NA"),
233d3Name("NA"),
234p3Name("NA"),
235cName("NA"),
236
237intR2Name("NA"),
238binCorrName("NA"),
239intBinCorrName("NA"),
240
241countsName("NA"),
242part_1_Name("NA"),
243part_2_Name("NA"),
244part_3_Name("NA"),
245pair_12_Name("NA"),
246pair_13_Name("NA"),
247pair_23_Name("NA"),
248tripletName("NA"),
249
250avg("NA"),
251avgName("NA"),
252sumName("NA"),
253s1ptName("NA"),
254s1ptNwName("NA"),
255s1DptName("NA"),
256
257s2PtPtName("NA"),
258s2NPtName("NA"),
259s2PtNName("NA"),
260s2DptDptName("NA"),
261
262s2PtPtNwName("NA"),
263s2NPtNwName("NA"),
264s2PtNNwName("NA"),
265
266ptName("NA"),
267ptptName("NA"),
268pt1pt1Name("NA"),
269DptName("NA"),
270DptDptName("NA"),
271RDptDptName("NA"),
272nPtName("NA"),
273ptNName("NA"),
274seanName("NA"),
275
276_title_counts("NA"),
277
278_title_m0("NA"),
279_title_m1("NA"),
280_title_m2("NA"),
281_title_m3("NA"),
282_title_m4("NA"),
283_title_m5("NA"),
284_title_m6("NA"),
285
286_title_eta_1("NA"),
287_title_phi_1("NA"),
288_title_pt_1("NA"),
289_title_etaPhi_1("NA"),
290_title_n_1("NA"),
291_title_SumPt_1("NA"),
292_title_AvgPt_1("NA"),
293_title_AvgN_1("NA"),
294_title_AvgSumPt_1("NA"),
295
296_title_eta_2("NA"),
297_title_phi_2("NA"),
298_title_pt_2("NA"),
299_title_etaPhi_2("NA"),
300_title_n_2("NA"),
301_title_SumPt_2("NA"),
302_title_AvgPt_2("NA"),
303_title_AvgN_2("NA"),
304_title_AvgSumPt_2("NA"),
305
306_title_etaPhi_12("NA"),
307
308_title_AvgN2_12("NA"),
309_title_AvgSumPtPt_12("NA"),
310_title_AvgSumPtN_12("NA"),
311_title_AvgNSumPt_12("NA"),
312
313vsZ("NA"),
314vsM("NA"),
315vsPt("NA"),
316vsPhi("NA"),
317vsEta("NA"),
318vsEtaPhi("NA"),
319vsPtVsPt("NA")
320{
321 printf("Default constructor called \n");
322
323 printf("passed \n ");
324
325}
326
327AliAnalysisTaskDptDptQA::AliAnalysisTaskDptDptQA(const TString & name)
328: AliAnalysisTaskSE(name),
329fAODEvent(0),
330fESDEvent(0),
331fInputHandler(0),
349d46ce 332_outputHistoList(0),
333_twoPi ( 2.0 * 3.1415927),
334_eventCount ( 0),
335_debugLevel ( 0),
336_singlesOnly ( 0),
337_useWeights ( 0),
338_sameFilter ( false),
339_rejectPileup ( 1),
340_rejectPairConversion ( 0),
341_vertexZMin ( -10.),
342_vertexZMax ( 10.),
343_vertexXYMin ( -10.),
344_vertexXYMax ( 10.),
345_centralityMethod ( 4),
346_centralityMin ( 0.),
347_centralityMax ( 1.),
348_requestedCharge_1 ( 1),
349_requestedCharge_2 ( -1),
350_dcaZMin ( -3),
351_dcaZMax ( 3.),
c5704d5c 352_dcaXYMin ( -3.),
353_dcaXYMax ( 3.),
349d46ce 354_dedxMin ( 0),
355_dedxMax ( 100000),
c5704d5c 356_nClusterMin ( 70),
349d46ce 357_trackFilterBit ( 0),
349d46ce 358_field ( 1.),
359_nTracks ( 0 ),
360_mult0 ( 0 ),
361_mult1 ( 0 ),
362_mult2 ( 0 ),
363_mult3 ( 0 ),
364_mult4 ( 0 ),
349d46ce 365_mult5 ( 0 ),
366_mult6 ( 0 ),
367arraySize ( 2000),
368_id_1(0),
369_charge_1(0),
370_iEtaPhi_1(0),
371_iPt_1(0),
372_pt_1(0),
373_px_1(0),
374_py_1(0),
375_pz_1(0),
376_correction_1(0),
377_dedx_1(0),
378_id_2(0),
379_charge_2(0),
380_iEtaPhi_2(0),
381_iPt_2(0),
382_pt_2(0),
383_px_2(0),
384_py_2(0),
385_pz_2(0),
386_correction_2(0),
387_dedx_2(0),
388_correctionWeight_1(0),
389_correctionWeight_2(0),
390_nBins_M0(500), _min_M0(0), _max_M0(10000), _width_M0(20),
391_nBins_M1(500), _min_M1(0), _max_M1(10000), _width_M1(20),
392_nBins_M2(500), _min_M2(0), _max_M2(10000), _width_M2(20),
393_nBins_M3(500), _min_M3(0), _max_M3(10000), _width_M3(20),
394_nBins_M4(100), _min_M4(0), _max_M4(1), _width_M4(0.01),
395_nBins_M5(100), _min_M5(0), _max_M5(1), _width_M5(0.01),
396_nBins_M6(100), _min_M6(0), _max_M6(1), _width_M6(0.01),
397_nBins_vertexZ(40), _min_vertexZ(-10), _max_vertexZ(10), _width_vertexZ(0.5),
398
399_nBins_pt_1(18), _min_pt_1(0.2), _max_pt_1(2.0), _width_pt_1(0.1),
400_nBins_phi_1(72), _min_phi_1(0), _max_phi_1(2.*3.1415927),_width_phi_1(2.*3.1415927/72.),
c5704d5c 401_nBins_eta_1(0), _min_eta_1(0), _max_eta_1(0), _width_eta_1(0.1),
349d46ce 402
403_nBins_etaPhi_1(0),
404_nBins_etaPhiPt_1(0),
405_nBins_zEtaPhiPt_1(0),
406
407_nBins_pt_2(18), _min_pt_2(0.2), _max_pt_2(2.0), _width_pt_2(0.1),
408_nBins_phi_2(72), _min_phi_2(0), _max_phi_2(2.*3.1415927),_width_phi_2(2.*3.1415927/72),
c5704d5c 409_nBins_eta_2(0), _min_eta_2(0), _max_eta_2(0), _width_eta_2(0.1),
349d46ce 410
411_nBins_etaPhi_2(0),
412_nBins_etaPhiPt_2(0),
413_nBins_zEtaPhiPt_2(0),
414_nBins_etaPhi_12(0),
415__n1_1(0),
416__n1_2(0),
417__n2_12(0),
418__s1pt_1(0),
419__s1pt_2(0),
420__s2ptpt_12(0),
421__s2NPt_12(0),
422__s2PtN_12(0),
423__n1Nw_1(0),
424__n1Nw_2(0),
425__n2Nw_12(0),
426__s1ptNw_1(0),
427__s1ptNw_2(0),
428__s2ptptNw_12(0),
429__s2NPtNw_12(0),
430__s2PtNNw_12(0),
431__n1_1_vsPt(0),
432__n1_1_vsEtaPhi(0),
433__s1pt_1_vsEtaPhi(0),
434__n1_1_vsZEtaPhiPt(0),
435__n1_2_vsPt(0),
436__n1_2_vsEtaPhi(0),
437__s1pt_2_vsEtaPhi(0),
438__n1_2_vsZEtaPhiPt(0),
439__n2_12_vsPtPt(0),
440__n2_12_vsEtaPhi(0),
441__s2ptpt_12_vsEtaPhi(0),
442__s2PtN_12_vsEtaPhi(0),
443__s2NPt_12_vsEtaPhi(0),
444_weight_1 ( 0 ),
445_weight_2 ( 0 ),
446_eventAccounting ( 0),
447_m0 ( 0),
448_m1 ( 0),
449_m2 ( 0),
450_m3 ( 0),
451_m4 ( 0),
452_m5 ( 0),
453_m6 ( 0),
c5704d5c 454_vertexZ ( 0),
349d46ce 455_etadis ( 0),
456_phidis ( 0),
c5704d5c 457_dcaz ( 0),
458_dcaxy ( 0),
349d46ce 459_n1_1_vsPt ( 0),
460_n1_1_vsEtaVsPhi ( 0),
461_s1pt_1_vsEtaVsPhi ( 0),
462_n1_1_vsZVsEtaVsPhiVsPt ( 0),
c5704d5c 463_n1_1_vsM ( 0), // w/ weight
349d46ce 464_s1pt_1_vsM ( 0),
c5704d5c 465_n1Nw_1_vsM ( 0), // w/o weight
349d46ce 466_s1ptNw_1_vsM ( 0),
467_dedxVsP_1 ( 0),
468_corrDedxVsP_1 ( 0),
469_betaVsP_1 ( 0),
470_n1_2_vsPt ( 0),
471_n1_2_vsEtaVsPhi ( 0),
472_s1pt_2_vsEtaVsPhi ( 0),
473_n1_2_vsZVsEtaVsPhiVsPt ( 0),
474_n1_2_vsM ( 0),
475_s1pt_2_vsM ( 0),
476_n1Nw_2_vsM ( 0),
477_s1ptNw_2_vsM ( 0),
478_dedxVsP_2 ( 0),
479_corrDedxVsP_2 ( 0),
480_betaVsP_2 ( 0),
481_n2_12_vsEtaPhi ( 0),
482_n2_12_vsPtVsPt ( 0),
483_s2PtPt_12_vsEtaPhi( 0),
484_s2PtN_12_vsEtaPhi ( 0),
485_s2NPt_12_vsEtaPhi ( 0),
486_n2_12_vsM ( 0),
487_s2PtPt_12_vsM ( 0),
488_s2PtN_12_vsM ( 0),
489_s2NPt_12_vsM ( 0),
490_n2Nw_12_vsM ( 0),
491_s2PtPtNw_12_vsM ( 0),
492_s2PtNNw_12_vsM ( 0),
493_s2NPtNw_12_vsM ( 0),
494_invMass ( 0),
c5704d5c 495_invMassElec ( 0),
349d46ce 496n1Name("NA"),
497n1NwName("NA"),
498n2Name("NA"),
499n2NwName("NA"),
500n3Name("NA"),
501n1n1Name("NA"),
502n1n1n1Name("NA"),
503n2n1Name("NA"),
504r1Name("NA"),
505r2Name("NA"),
506r3Name("NA"),
507r2r1Name("NA"),
508c2Name("NA"),
509c3Name("NA"),
510d3Name("NA"),
511p3Name("NA"),
512cName("NA"),
513
514intR2Name("NA"),
515binCorrName("NA"),
516intBinCorrName("NA"),
517
518countsName("NA"),
519part_1_Name("NA"),
520part_2_Name("NA"),
521part_3_Name("NA"),
522pair_12_Name("NA"),
523pair_13_Name("NA"),
524pair_23_Name("NA"),
525tripletName("NA"),
526
527avg("NA"),
528avgName("NA"),
529sumName("NA"),
530s1ptName("NA"),
531s1ptNwName("NA"),
532s1DptName("NA"),
533
534s2PtPtName("NA"),
535s2NPtName("NA"),
536s2PtNName("NA"),
537s2DptDptName("NA"),
538
539s2PtPtNwName("NA"),
540s2NPtNwName("NA"),
541s2PtNNwName("NA"),
542
543ptName("NA"),
544ptptName("NA"),
545pt1pt1Name("NA"),
546DptName("NA"),
547DptDptName("NA"),
548RDptDptName("NA"),
549nPtName("NA"),
550ptNName("NA"),
551seanName("NA"),
552
553_title_counts("NA"),
554
555_title_m0("NA"),
556_title_m1("NA"),
557_title_m2("NA"),
558_title_m3("NA"),
559_title_m4("NA"),
560_title_m5("NA"),
561_title_m6("NA"),
562
563_title_eta_1("NA"),
564_title_phi_1("NA"),
565_title_pt_1("NA"),
566_title_etaPhi_1("NA"),
567_title_n_1("NA"),
568_title_SumPt_1("NA"),
569_title_AvgPt_1("NA"),
570_title_AvgN_1("NA"),
571_title_AvgSumPt_1("NA"),
572
573_title_eta_2("NA"),
574_title_phi_2("NA"),
575_title_pt_2("NA"),
576_title_etaPhi_2("NA"),
577_title_n_2("NA"),
578_title_SumPt_2("NA"),
579_title_AvgPt_2("NA"),
580_title_AvgN_2("NA"),
581_title_AvgSumPt_2("NA"),
582
583_title_etaPhi_12("NA"),
584
585_title_AvgN2_12("NA"),
586_title_AvgSumPtPt_12("NA"),
587_title_AvgSumPtN_12("NA"),
588_title_AvgNSumPt_12("NA"),
589
590vsZ("NA"),
591vsM("NA"),
592vsPt("NA"),
593vsPhi("NA"),
594vsEta("NA"),
595vsEtaPhi("NA"),
596vsPtVsPt("NA")
597{
598 printf("2nd constructor called ");
599
600 DefineOutput(0, TList::Class());
601
602 printf("passed ");
603
604}
605
606AliAnalysisTaskDptDptQA::~AliAnalysisTaskDptDptQA()
607{
c5704d5c 608
349d46ce 609}
610
611void AliAnalysisTaskDptDptQA::UserCreateOutputObjects()
612{
c5704d5c 613 cout<< "AliAnalysisTaskDptDptQA::CreateOutputObjects() Starting " << endl;
349d46ce 614 OpenFile(0);
615 _outputHistoList = new TList();
616 _outputHistoList->SetOwner();
617
c5704d5c 618 //Setup the parameters of the histograms
349d46ce 619 _nBins_M0 = 500; _min_M0 = 0.; _max_M0 = 5000.; _width_M0 = (_max_M0-_min_M0)/_nBins_M0;
620 _nBins_M1 = 500; _min_M1 = 0.; _max_M1 = 5000.; _width_M1 = (_max_M1-_min_M1)/_nBins_M1;
621 _nBins_M2 = 500; _min_M2 = 0.; _max_M2 = 5000.; _width_M2 = (_max_M2-_min_M2)/_nBins_M2;
622 _nBins_M3 = 500; _min_M3 = 0.; _max_M3 = 5000.; _width_M3 = (_max_M3-_min_M3)/_nBins_M3;
623 _nBins_M4 = 100; _min_M4 = 0.; _max_M4 = 100.; _width_M4 = (_max_M4-_min_M4)/_nBins_M4;
624 _nBins_M5 = 100; _min_M5 = 0.; _max_M5 = 100.; _width_M5 = (_max_M5-_min_M5)/_nBins_M5;
625 _nBins_M6 = 100; _min_M6 = 0.; _max_M6 = 100.; _width_M6 = (_max_M6-_min_M6)/_nBins_M6;
626
627 _min_vertexZ = _vertexZMin;
628 _max_vertexZ = _vertexZMax;
629 _width_vertexZ = 0.5;
630 _nBins_vertexZ = int(0.5+ (_max_vertexZ - _min_vertexZ)/_width_vertexZ);
631 _nBins_pt_1 = int(0.5+ (_max_pt_1 -_min_pt_1 )/_width_pt_1);
632 _nBins_eta_1 = int(0.5+ (_max_eta_1-_min_eta_1)/_width_eta_1);
633 _width_phi_1 = (_max_phi_1 - _min_phi_1) /_nBins_phi_1;
634 _nBins_etaPhi_1 = _nBins_phi_1 * _nBins_eta_1;
635 _nBins_etaPhiPt_1 = _nBins_etaPhi_1 * _nBins_pt_1;
636 _nBins_zEtaPhiPt_1 = _nBins_vertexZ * _nBins_etaPhiPt_1;
637
638 _nBins_pt_2 = int(0.5+ (_max_pt_2 -_min_pt_2 )/_width_pt_2);
639 _nBins_eta_2 = int(0.5+ (_max_eta_2-_min_eta_2)/_width_eta_2);
640 _width_phi_2 = (_max_phi_2 - _min_phi_2) /_nBins_phi_2;
641 _nBins_etaPhi_2 = _nBins_phi_2 * _nBins_eta_2;
642 _nBins_etaPhiPt_2 = _nBins_etaPhi_2 * _nBins_pt_2;
643 _nBins_zEtaPhiPt_2 = _nBins_vertexZ * _nBins_etaPhiPt_2;
644 _nBins_etaPhi_12 = _nBins_etaPhi_1 * _nBins_etaPhi_2;
c5704d5c 645
646 //setup the work arrays
647
349d46ce 648 _id_1 = new int[arraySize];
649 _charge_1 = new int[arraySize];
c5704d5c 650 //_iPhi_1 = new int[arraySize];
651 //_iEta_1 = new int[arraySize];
349d46ce 652 _iEtaPhi_1 = new int[arraySize];
653 _iPt_1 = new int[arraySize];
654 _pt_1 = new float[arraySize];
655 _px_1 = new float[arraySize];
656 _py_1 = new float[arraySize];
657 _pz_1 = new float[arraySize];
c5704d5c 658 //_phi_1 = new float[arraySize];
659 //_eta_1 = new float[arraySize];
349d46ce 660 _correction_1 = new float[arraySize];
c5704d5c 661 _dedx_1 = new float[arraySize];
662
349d46ce 663 __n1_1_vsPt = getDoubleArray(_nBins_pt_1, 0.);
664 __n1_1_vsEtaPhi = getDoubleArray(_nBins_etaPhi_1, 0.);
665 __s1pt_1_vsEtaPhi = getDoubleArray(_nBins_etaPhi_1, 0.);
666 __n1_1_vsZEtaPhiPt = getFloatArray(_nBins_zEtaPhiPt_1, 0.);
667
c5704d5c 668 cout << "==========================================================================================" << endl;
669 cout << "=============== Booking for particle 1 done." << endl;
670 cout << "_requestedCharge_1: " << _requestedCharge_1 << endl;
671 cout << "_requestedCharge_2: " << _requestedCharge_2 << endl;
672
349d46ce 673 if (_requestedCharge_2!=_requestedCharge_1)
674 {
c5704d5c 675 cout << " creating arrays for particle 2 with size: " << arraySize << endl;
676 _sameFilter = 0;
349d46ce 677 //particle 2
678 _id_2 = new int[arraySize];
679 _charge_2 = new int[arraySize];
c5704d5c 680 //_iPhi_2 = new int[arraySize];
681 //_iEta_2 = new int[arraySize];
349d46ce 682 _iEtaPhi_2 = new int[arraySize];
683 _iPt_2 = new int[arraySize];
684 _pt_2 = new float[arraySize];
685 _px_2 = new float[arraySize];
686 _py_2 = new float[arraySize];
687 _pz_2 = new float[arraySize];
c5704d5c 688 //_phi_2 = new float[arraySize];
689 //_eta_2 = new float[arraySize];
349d46ce 690 _correction_2 = new float[arraySize];
c5704d5c 691 _dedx_2 = new float[arraySize];
692
349d46ce 693 __n1_2_vsPt = getDoubleArray(_nBins_pt_2, 0.);
694 __n1_2_vsEtaPhi = getDoubleArray(_nBins_etaPhi_2, 0.);
695 __s1pt_2_vsEtaPhi = getDoubleArray(_nBins_etaPhi_2, 0.);
696 __n1_2_vsZEtaPhiPt = getFloatArray(_nBins_zEtaPhiPt_2, 0.);
697
698 }
699
700 __n2_12_vsPtPt = getDoubleArray(_nBins_pt_1*_nBins_pt_2,0.);
701 __n2_12_vsEtaPhi = getFloatArray(_nBins_etaPhi_12, 0.);
702 __s2ptpt_12_vsEtaPhi = getFloatArray(_nBins_etaPhi_12, 0.);
703 __s2PtN_12_vsEtaPhi = getFloatArray(_nBins_etaPhi_12, 0.);
704 __s2NPt_12_vsEtaPhi = getFloatArray(_nBins_etaPhi_12, 0.);
705
706 // Setup all the labels needed.
707
708 part_1_Name = "_1";
709 part_2_Name = "_2";
710 pair_12_Name = "_12";
711
712 n1Name = "n1";
713 n2Name = "n2";
714 n1NwName = "n1Nw";
715 n2NwName = "n2Nw";
716 r1Name = "r1";
717 r2Name = "r2";
718 r3Name = "r3";
719 r2r1Name = "r2r1";
720 c2Name = "c2";
721 c3Name = "c3";
722 d3Name = "d3";
723 p3Name = "p3";
724 cName = "sean";
725
726 intR2Name = "intR2";
727 binCorrName = "binCorr";
728 intBinCorrName = "intBinCorr";
729
730 avgName = "avg";
731 sumName = "sum";
732 s1ptName = "sumPt";
733 s1ptNwName = "sumPtNw";
734 s1DptName = "sumDpt";
735 s2PtPtName = "sumPtPt";
736 s2PtPtNwName = "sumPtPtNw";
737 s2DptDptName = "sumDptDpt";
738 s2NPtName = "sumNPt";
739 s2NPtNwName = "sumNPtNw";
740 s2PtNName = "sumPtN";
741 s2NPtNwName = "sumNPtNw";
742 s2PtNNwName = "sumPtNNw";
743 ptName = "avgPt";
744 ptptName = "avgPtPt";
745 pt1pt1Name = "avgPtavgPt";
746 DptName = "avgDpt";
747 DptDptName = "avgDptDpt";
748 RDptDptName = "relDptDpt"; // ratio of avgDptDpt by avgPt*avgPt
749 nPtName = "avgNpt";
750 ptNName = "avgPtN";
751 seanName = "seanC";
752
753 _title_counts = "yield";
754
755 _title_m0 = "M_{0}";
756 _title_m1 = "M_{1}";
757 _title_m2 = "M_{2}";
758 _title_m3 = "M_{3}";
759 _title_m4 = "V0Centrality";
760 _title_m5 = "TrkCentrality";
761 _title_m6 = "SpdCentrality";
762
763 _title_eta_1 = "#eta_{1}";
764 _title_phi_1 = "#varphi_{1} (radian)";
765 _title_etaPhi_1 = "#eta_{1}#times#varphi_{1}";
766 _title_pt_1 = "p_{t,1} (GeV/c)";
767 _title_n_1 = "n_{1}";
768 _title_SumPt_1 = "#Sigma p_{t,1} (GeV/c)";
769 _title_AvgPt_1 = "#LT p_{t,1} #GT (GeV/c)";
770 _title_AvgN_1 = "#LT n_{1} #GT";
771 _title_AvgSumPt_1 = "#LT #Sigma p_{t,1} #GT (GeV/c)";
772
773 _title_eta_2 = "#eta_{2}";
774 _title_phi_2 = "#varphi_{2} (radian)";
775 _title_etaPhi_2 = "#eta_{2}#times#varphi_{2}";
776 _title_pt_2 = "p_{t,2} (GeV/c)";
777 _title_n_2 = "n_{2}";
778 _title_SumPt_2 = "#Sigma p_{t,1} (GeV/c)";
779 _title_AvgPt_2 = "#LT p_{t,2} #GT (GeV/c)";
780 _title_AvgN_2 = "#LT n_{2} #GT";
781 _title_AvgSumPt_2 = "#LT #Sigma p_{t,2} #GT (GeV/c)";
782
783 _title_etaPhi_12 = "#eta_{1}#times#varphi_{1}#times#eta_{2}#times#varphi_{2}";
784
785 _title_AvgN2_12 = "#LT n_{2} #GT";;
786 _title_AvgSumPtPt_12 = "#LT #Sigma p_{t,1}p_{t,2} #GT";;
787 _title_AvgSumPtN_12 = "#LT #Sigma p_{t,1}N #GT";;
788 _title_AvgNSumPt_12 = "#LT N#Sigma p_{t,2} #GT";;
789
790
791 vsZ = "_vsZ";
792 vsM = "_vsM";
793 vsPt = "_vsPt";
794 vsPhi = "_vsPhi";
795 vsEta = "_vsEta";
796 vsEtaPhi = "_vsEtaPhi";
797 vsPtVsPt = "_vsPtVsPt";
798
799
800 if (_useWeights)
801 {
80d059c6 802 int iZ, iEtaPhi, iPt;
803 int iZ1,iEtaPhi1,iPt1;
804 int a, b;
805 if (_weight_1)
806 {
807
808 cout<<"Prabhat Sanu weeding is here *************"<<" "<<endl;
809
810 _correctionWeight_1 = new float[_nBins_vertexZ*_nBins_etaPhi_1*_nBins_pt_1];
811 a = _nBins_pt_1;
812 b = _nBins_etaPhi_1*_nBins_pt_1;
813 for (iZ=0,iZ1=1; iZ<_nBins_vertexZ; iZ++, iZ1++)
814 {
815 for (iEtaPhi=0,iEtaPhi1=1; iEtaPhi<_nBins_etaPhi_1; iEtaPhi++, iEtaPhi1++)
816 {
817 for (iPt=0,iPt1=1; iPt<_nBins_pt_1; iPt++, iPt1++)
818 {
819 _correctionWeight_1[iZ*b+iEtaPhi*a+iPt] = _weight_1->GetBinContent(iZ1,iEtaPhi1,iPt1);
820 }
821 }
822 }
823 } // _weight_1
349d46ce 824 else
80d059c6 825 {
826 AliError("AliAnalysisTaskDptDptQA:: _weight_1 is a null pointer.");
827 return;
828 }
829 if (!_sameFilter)
830 {
831 if (_weight_2)
832 {
833 _correctionWeight_2 = new float[_nBins_vertexZ*_nBins_etaPhi_2*_nBins_pt_2];
834 a = _nBins_pt_2;
835 b = _nBins_etaPhi_2*_nBins_pt_2;
836 for (iZ=0,iZ1=1; iZ<_nBins_vertexZ; iZ++, iZ1++)
837 {
838 for (iEtaPhi=0,iEtaPhi1=1; iEtaPhi<_nBins_etaPhi_2; iEtaPhi++, iEtaPhi1++)
839 {
840 for (iPt=0,iPt1=1; iPt<_nBins_pt_2; iPt++, iPt1++)
841 {
842 _correctionWeight_2[iZ*b+iEtaPhi*a+iPt] = _weight_2->GetBinContent(iZ1,iEtaPhi1,iPt1);
843 }
844 }
845 }
846 } // _weight_2
847 else
848 {
849 AliError("AliAnalysisTaskDptDptQA:: _weight_1 is a null pointer.");
850 return;
851 }
852 }
349d46ce 853 }
854
855 createHistograms();
856 PostData(0,_outputHistoList);
857
c5704d5c 858 cout<< "AliAnalysisTaskDptDptQA::CreateOutputObjects() DONE " << endl;
349d46ce 859
860}
861
862void AliAnalysisTaskDptDptQA::createHistograms()
863{
864 AliInfo(" AliAnalysisTaskDptDptQA::createHistoHistograms() Creating Event Histos");
865 TString name;
866
867 name = "eventAccounting";
868
c5704d5c 869 // bin index : what it is...
870 // 0 : number of event submitted
871 // 1 : number accepted by centrality cut
872 // 2 : number accepted by centrality cut and z cut
873 // 3 : total number of particles that satisfied filter 1
874 // 4 : total number of particles that satisfied filter 2
875 _eventAccounting = createHisto1D(name,name,10, -0.5, 9.5, "event Code", _title_counts);
349d46ce 876
877 name = "m0"; _m0 = createHisto1D(name,name,_nBins_M1, _min_M1, _max_M1, _title_m0, _title_counts);
878 name = "m1"; _m1 = createHisto1D(name,name,_nBins_M1, _min_M1, _max_M1, _title_m1, _title_counts);
879 name = "m2"; _m2 = createHisto1D(name,name,_nBins_M2, _min_M2, _max_M2, _title_m2, _title_counts);
880 name = "m3"; _m3 = createHisto1D(name,name,_nBins_M3, _min_M3, _max_M3, _title_m3, _title_counts);
881 name = "m4"; _m4 = createHisto1D(name,name,_nBins_M4, _min_M4, _max_M4, _title_m4, _title_counts);
882 name = "m5"; _m5 = createHisto1D(name,name,_nBins_M5, _min_M5, _max_M5, _title_m5, _title_counts);
883 name = "m6"; _m6 = createHisto1D(name,name,_nBins_M6, _min_M6, _max_M6, _title_m6, _title_counts);
c5704d5c 884 name = "zV"; _vertexZ = createHisto1D(name,name,_nBins_vertexZ, _min_vertexZ, _max_vertexZ, "z-Vertex (cm)", _title_counts);
349d46ce 885
349d46ce 886 name = "Eta"; _etadis = createHisto1F(name,name, 200, -1.0, 1.0, "#eta","counts");
887 name = "Phi"; _phidis = createHisto1F(name,name, 360, 0.0, 6.4, "#phi","counts");
80d059c6 888
349d46ce 889 if (_singlesOnly)
890 {
80d059c6 891 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);
892 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);
893 //name = "dedxVsP_1"; _dedxVsP_1 = createHisto2F(name,name,400,-2.,2.,120,0.,120.,"p (GeV/c)", "dedx", "counts");
894 //name = "corrDedxVsP_1"; _corrDedxVsP_1 = createHisto2F(name,name,400,-2.,2.,120,0.,120.,"p (GeV/c)", "dedx", "counts");
895 //name = "betaVsP_1"; _betaVsP_1 = createHisto2F(name,name,400,-2.,2.,120,0.5,1.1,"p (GeV/c)", "beta", "counts");
896
897 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);
898 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);
899 //name = "dedxVsP_2"; _dedxVsP_2 = createHisto2F(name,name,400,-2.,2.,120,0.,120.,"p (GeV/c)", "dedx", "counts");
900 //name = "corrDedxVsP_2"; _corrDedxVsP_2 = createHisto2F(name,name,400,-2.,2.,120,0.,120.,"p (GeV/c)", "dedx", "counts");
901 //name = "betaVsP_2"; _betaVsP_2 = createHisto2F(name,name,400,-2.,2.,120,0.5,1.1,"p (GeV/c)", "beta", "counts");
349d46ce 902 }
903 else
904 {
80d059c6 905 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);
906 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);
907 name = n1Name+part_1_Name+vsM; _n1_1_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgN_1);
908 name = s1ptName+part_1_Name+vsM; _s1pt_1_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgSumPt_1);
909 name = n1NwName+part_1_Name+vsM; _n1Nw_1_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgN_1);
910 name = s1ptNwName+part_1_Name+vsM; _s1ptNw_1_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgSumPt_1);
911
912 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);
913 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);
914 name = n1Name+part_2_Name + vsM; _n1_2_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgN_2);
915 name = s1ptName+part_2_Name + vsM; _s1pt_2_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgSumPt_2);
916 name = n1NwName+part_2_Name+vsM; _n1Nw_2_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgN_1);
917 name = s1ptNwName+part_2_Name+vsM; _s1ptNw_2_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgSumPt_1);
918
919 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);
349d46ce 920 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);
921 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);
922 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);
923 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);
924
925 name = n2Name+pair_12_Name + vsM; _n2_12_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgN2_12);
926 name = s2PtPtName+pair_12_Name + vsM; _s2PtPt_12_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgSumPtPt_12);
927 name = s2PtNName+pair_12_Name + vsM; _s2PtN_12_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgSumPtN_12);
928 name = s2NPtName+pair_12_Name + vsM; _s2NPt_12_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgNSumPt_12);
929
930 name = n2NwName+pair_12_Name + vsM; _n2Nw_12_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgN2_12);
931 name = s2PtPtNwName+pair_12_Name + vsM; _s2PtPtNw_12_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgSumPtPt_12);
932 name = s2PtNNwName+pair_12_Name + vsM; _s2PtNNw_12_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgSumPtN_12);
933 name = s2NPtNwName+pair_12_Name + vsM; _s2NPtNw_12_vsM = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgNSumPt_12);
934
c5704d5c 935 name = "mInv"; _invMass = createHisto1F(name,name, 500, 0., 1.000, "M_{inv}","counts");
936 name = "mInvElec"; _invMassElec = createHisto1F(name,name, 500, 0., 1.000, "M_{inv}","counts");
349d46ce 937 }
938
939 AliInfo(" AliAnalysisTaskDptDptQA::createHistoHistograms() All Done");
940}
941//-----------------------//
942
943void AliAnalysisTaskDptDptQA::finalizeHistograms()
944{
945
946 AliInfo("AliAnalysisTaskDptDptQA::finalizeHistograms() starting");
947 AliInfo(Form("CorrelationAnalyzers::finalizeHistograms() _eventCount : %d",int(_eventCount)));
c5704d5c 948 if (_singlesOnly)
949 {
950 if (_sameFilter)
951 {
952 fillHistoWithArray(_n1_1_vsPt, __n1_1_vsPt, _nBins_pt_1);
953 fillHistoWithArray(_n1_1_vsZVsEtaVsPhiVsPt, __n1_1_vsZEtaPhiPt, _nBins_vertexZ, _nBins_etaPhi_1, _nBins_pt_1);
954 fillHistoWithArray(_n1_2_vsPt, __n1_1_vsPt, _nBins_pt_1);
955 fillHistoWithArray(_n1_2_vsZVsEtaVsPhiVsPt, __n1_1_vsZEtaPhiPt, _nBins_vertexZ, _nBins_etaPhi_1, _nBins_pt_1);
956 }
957 else
958 {
959 fillHistoWithArray(_n1_1_vsPt, __n1_1_vsPt, _nBins_pt_1);
960 fillHistoWithArray(_n1_1_vsZVsEtaVsPhiVsPt, __n1_1_vsZEtaPhiPt, _nBins_vertexZ, _nBins_etaPhi_1, _nBins_pt_1);
961 fillHistoWithArray(_n1_2_vsPt, __n1_2_vsPt, _nBins_pt_2);
962 fillHistoWithArray(_n1_2_vsZVsEtaVsPhiVsPt, __n1_2_vsZEtaPhiPt, _nBins_vertexZ, _nBins_etaPhi_2, _nBins_pt_2);
963 }
964 }
965 else
966 {
967 if (_sameFilter)
968 {
969 fillHistoWithArray(_n1_1_vsEtaVsPhi, __n1_1_vsEtaPhi, _nBins_eta_1, _nBins_phi_1);
970 fillHistoWithArray(_s1pt_1_vsEtaVsPhi, __s1pt_1_vsEtaPhi, _nBins_eta_1, _nBins_phi_1);
971 fillHistoWithArray(_n1_2_vsEtaVsPhi, __n1_1_vsEtaPhi, _nBins_eta_1, _nBins_phi_1);
972 fillHistoWithArray(_s1pt_2_vsEtaVsPhi, __s1pt_1_vsEtaPhi, _nBins_eta_1, _nBins_phi_1);
973 }
974 else
975 {
976 fillHistoWithArray(_n1_1_vsEtaVsPhi, __n1_1_vsEtaPhi, _nBins_eta_1, _nBins_phi_1);
977 fillHistoWithArray(_s1pt_1_vsEtaVsPhi, __s1pt_1_vsEtaPhi, _nBins_eta_1, _nBins_phi_1);
978 fillHistoWithArray(_n1_2_vsEtaVsPhi, __n1_2_vsEtaPhi, _nBins_eta_2, _nBins_phi_2);
979 fillHistoWithArray(_s1pt_2_vsEtaVsPhi, __s1pt_2_vsEtaPhi, _nBins_eta_2, _nBins_phi_2);
980 }
981 fillHistoWithArray(_n2_12_vsEtaPhi, __n2_12_vsEtaPhi, _nBins_etaPhi_12);
982 fillHistoWithArray(_s2PtPt_12_vsEtaPhi, __s2ptpt_12_vsEtaPhi, _nBins_etaPhi_12);
983 fillHistoWithArray(_s2PtN_12_vsEtaPhi, __s2PtN_12_vsEtaPhi, _nBins_etaPhi_12);
984 fillHistoWithArray(_s2NPt_12_vsEtaPhi, __s2NPt_12_vsEtaPhi, _nBins_etaPhi_12);
985 fillHistoWithArray(_n2_12_vsPtVsPt, __n2_12_vsPtPt, _nBins_pt_1, _nBins_pt_2);
986
987 }
349d46ce 988 AliInfo("AliAnalysisTaskDptDptQA::finalizeHistograms() Done ");
989}
990//--------------//
991
992
993void AliAnalysisTaskDptDptQA::UserExec(Option_t */*option*/)
994{
995
996 int k1,k2;
c5704d5c 997 int iPhi, iEta, iEtaPhi, iPt, charge;
998 float q, p, phi, pt, eta, corr, corrPt, dedx, px, py, pz;
999 int ij;
1000 int id_1, q_1, iEtaPhi_1, iPt_1;
1001 float pt_1, px_1, py_1, pz_1, corr_1, dedx_1;
1002 int id_2, q_2, iEtaPhi_2, iPt_2;
1003 float pt_2, px_2, py_2, pz_2, corr_2, dedx_2;
1004 float ptpt;
349d46ce 1005 int iVertex, iVertexP1, iVertexP2;
c5704d5c 1006 int iZEtaPhiPt;
1007 float massElecSq = 2.5e-7;
349d46ce 1008 const AliAODVertex* vertex;
c5704d5c 1009 int nClus;
349d46ce 1010 bool bitOK;
c5704d5c 1011
349d46ce 1012 // count all events looked at here
1013 _eventCount++;
c5704d5c 1014
1015 // A. Adare - Fix from above to avoid coverity complaint about null dereference
349d46ce 1016 if (_eventAccounting)
1017 {
1018 _eventAccounting->Fill(0);// count all calls to this function
1019 }
1020 else
1021 {
c5704d5c 1022 cout << "AliAnalysisTaskDptDptQA::UserExec(Option_t *option) - !_eventAccounting" << endl;
1023 cout << "AliAnalysisTaskDptDptQA::UserExec(Option_t *option) - Skip this task" << endl;
1024 return;
1025 }
1026
1027 //cout << "AliAnalysisTaskDptDptQA::UserExec(Option_t *option) - 1" << endl;
1028
1029 //Get pointer to current event
1030 //fESDEvent = dynamic_cast<AliESDEvent*> (InputEvent());
1031 fAODEvent = dynamic_cast<AliAODEvent*> (InputEvent());
1032
1033 // if(!fAODEvent && !fESDEvent)
1034 if(!fAODEvent)
1035 {
1036 AliError("AliAnalysisTaskDptDptQA::Exec(Option_t *option) !fAODEvent");
1037 return;
349d46ce 1038 }
1039
1040 _eventAccounting->Fill(1);// count all calls to this function with a valid pointer
1041 //reset single particle counters
1042 k1 = k2 = 0;
1043 __n1_1 = __n1_2 = __s1pt_1 = __s1pt_2 = __n1Nw_1 = __n1Nw_2 = __s1ptNw_1 = __s1ptNw_2 = 0;
1044
1045 float v0Centr = -999.;
349d46ce 1046 float trkCentr = -999.;
1047 float spdCentr = -999.;
349d46ce 1048 float vertexX = -999;
1049 float vertexY = -999;
1050 float vertexZ = -999;
c5704d5c 1051 float vertexXY = -999;
349d46ce 1052 float centrality = -999;
1053
349d46ce 1054 if(fAODEvent)
1055 {
c5704d5c 1056 //cout << "AliAnalysisTaskDptDptQA::UserExec(Option_t *option) - 5" << endl;
1057
1058 //Centrality
1059 AliCentrality* centralityObject = fAODEvent->GetHeader()->GetCentralityP();
1060 if (centralityObject)
1061 {
1062 //cout << "AliAnalysisTaskDptDptQA::UserExec(Option_t *option) - 6" << endl;
349d46ce 1063
c5704d5c 1064 v0Centr = centralityObject->GetCentralityPercentile("V0M");
1065 trkCentr = centralityObject->GetCentralityPercentile("TRK");
1066 spdCentr = centralityObject->GetCentralityPercentile("CL1");
1067 //cout << "AliAnalysisTaskDptDptQA::UserExec(Option_t *option) - 7" << endl;
349d46ce 1068
c5704d5c 1069 }
1070 //cout << "AliAnalysisTaskDptDptQA::UserExec(Option_t *option) - 8" << endl;
1071
1072 _nTracks = fAODEvent->GetNTracks();
1073 //_mult0 = 0;
1074 //_mult1 = 0;
1075 //_mult2 = 0;
1076 _mult3 = _nTracks;
1077 _mult4 = v0Centr;
1078 _mult5 = trkCentr;
1079 _mult6 = spdCentr;
1080 _field = fAODEvent->GetMagneticField();
1081
1082 //cout << "AliAnalysisTaskDptDptQA::UserExec(Option_t *option) - _field:" << _field << endl;
1083
1084 //_centralityMethod
1085 switch (_centralityMethod)
1086 {
1087 case 0: centrality = _mult0; break;
1088 case 1: centrality = _mult1; break;
1089 case 2: centrality = _mult2; break;
1090 case 3: centrality = _mult3; break;
1091 case 4: centrality = _mult4; break;
1092 case 5: centrality = _mult5; break;
1093 case 6: centrality = _mult6; break;
1094 }
1095
1096 //cout << "AliAnalysisTaskDptDptQA::UserExec(Option_t *option) - 10" << endl;
1097
1098 //filter on centrality
1099 // require the v0 and trk centralities to agree within 5%.
1100
1101 //if ((Float_t(GetTPCMult(fAOD)) > (-40.3+1.22*GetGlobalMult(fAOD))) && (Float_t(GetTPCMult(fAOD)) < (32.1+1.59*GetGlobalMult(fAOD)))) {//3 sigma
1102
1103
1104 if ( centrality < _centralityMin ||
1105 centrality > _centralityMax ||
1106 fabs(v0Centr-trkCentr)>5.0)
1107 {
1108 //cout << "AliAnalysisTaskDptDptQA::UserExec(Option_t *option) - 11" << endl;
349d46ce 1109
c5704d5c 1110 // we dont analyze this event here...
1111 return;
1112 }
1113 _eventAccounting->Fill(2);// count all events with right centrality
1114 //cout << "AliAnalysisTaskDptDptQA::UserExec(Option_t *option) - 12" << endl;
1115
1116 // filter on z and xy vertex
1117 vertex = (AliAODVertex*) fAODEvent->GetPrimaryVertexSPD();
1118 if (!vertex || vertex->GetNContributors()<1)
1119 {
1120 vertexZ = -999;
1121 vertexXY = -999;
1122 cout << "AliAnalysisTaskDptDptQA::UserExec(Option_t *option) - No valid vertex object or poor vertex" << endl;
1123 }
1124 else
1125 {
1126 vertexX = vertex->GetX();
1127 vertexY = vertex->GetY();
1128 vertexZ = vertex->GetZ();
1129 vertexXY = sqrt(vertexX*vertexX+vertexY*vertexY);
1130 }
1131 if (!vertex ||
1132 vertexZ < _vertexZMin ||
1133 vertexZ > _vertexZMax ||
1134 vertexXY < _vertexXYMin ||
1135 vertexXY > _vertexXYMax)
1136 return;
1137
1138 iVertex = int((vertexZ-_min_vertexZ)/_width_vertexZ);
1139 iVertexP1 = iVertex*_nBins_etaPhiPt_1;
1140 iVertexP2 = iVertex*_nBins_etaPhiPt_2;
1141 if (iVertex<0 || iVertex>=_nBins_vertexZ)
1142 {
1143 AliError("AliAnalysisTaskDptDptQA::Exec(Option_t *option) iVertex<0 || iVertex>=_nBins_vertexZ ");
1144 return;
1145 }
1146 _eventAccounting->Fill(3);// count all calls to this function with a valid pointer
1147
1148 // loop over all particles for singles
1149 //debug() << "_nTracks:"<< _nTracks << endl;
1150 for (int iTrack=0; iTrack< _nTracks; iTrack++)
1151 {
1152 //cout << "AliAnalysisTaskDptDptQA::UserExec(Option_t *option) - 16 iTrack:" << iTrack << endl;
349d46ce 1153
c5704d5c 1154 AliAODTrack * t = (AliAODTrack *) fAODEvent->GetTrack(iTrack);
1155 if (!t)
1156 {
1157 AliError(Form("AliAnalysisTaskDptDptQA::Exec(Option_t *option) No track ofr iTrack=%d", iTrack));
1158 continue;
1159 }
349d46ce 1160
c5704d5c 1161 q = t->Charge();
1162 charge = int(q);
1163 phi = t->Phi();
1164 p = t->P();
1165 pt = t->Pt();
1166 px = t->Px();
1167 py = t->Py();
1168 pz = t->Pz();
1169 eta = t->Eta();
1170 dedx = t->GetTPCsignal();
1171 nClus = t->GetTPCNcls();
1172 bitOK = t->TestFilterBit(_trackFilterBit);
349d46ce 1173
c5704d5c 1174 //cout << "_trackFilterBit:" << _trackFilterBit << " Track returns:" << bitOK << endl;
1175 //cout << " q:" << q << " _requestedCharge_1:" << _requestedCharge_1 << endl;
1176 //cout << " pt:" << pt << " _min_pt_1:" << _min_pt_1 << " _max_pt_1:" << _max_pt_1<< endl;
1177 //cout << " phi:" << phi << endl;
1178 //cout << " eta:" << eta << " _min_eta_1:" << _min_eta_1 << " _max_eta_1:" << _max_eta_1<< endl;
1179 //cout << " dedx:" << dedx << " _dedxMin:" << _dedxMin << " _dedxMax:" << _dedxMax << endl;
1180 //cout << "nclus:" << nClus<< " _nClusterMin:" << _nClusterMin << endl;
349d46ce 1181
c5704d5c 1182 //common cuts first
349d46ce 1183
c5704d5c 1184 if (!bitOK || dedx<_dedxMin || dedx>_dedxMax || nClus<_nClusterMin) continue;
1185 // Kinematics cuts used
1186
1187 if( pt < _min_pt_1 || pt > _max_pt_1) continue;
1188 if( eta < _min_eta_1 || eta > _max_eta_1) continue;
1189
1190 _etadis->Fill(eta);
1191 _phidis->Fill(phi);
349d46ce 1192
80d059c6 1193 //Particle 1
1194 if (_requestedCharge_1 == charge &&
1195 eta >= 0.2 &&
1196 eta < 0.3)
c5704d5c 1197 {
1198 iPhi = int( phi/_width_phi_1);
349d46ce 1199
c5704d5c 1200 if (iPhi<0 || iPhi>=_nBins_phi_1 )
1201 {
1202 AliWarning("AliAnalysisTaskDptDptQA::analyze() iPhi<0 || iPhi>=_nBins_phi_1");
1203 return;
1204 }
349d46ce 1205
c5704d5c 1206 iEta = int((eta-_min_eta_1)/_width_eta_1);
1207 if (iEta<0 || iEta>=_nBins_eta_1)
1208 {
1209 AliWarning(Form("AliAnalysisTaskDptDptQA::analyze(AliceEvent * event) Mismatched iEta: %d", iEta));
1210 continue;
1211 }
1212 iPt = int((pt -_min_pt_1 )/_width_pt_1 );
1213 if (iPt<0 || iPt >=_nBins_pt_1)
1214 {
1215 AliWarning(Form("AliAnalysisTaskDptDptQA::analyze(AliceEvent * event) Mismatched iPt: %d",iPt));
1216 continue;
1217 }
1218 iEtaPhi = iEta*_nBins_phi_1+iPhi;
1219 iZEtaPhiPt = iVertexP1 + iEtaPhi*_nBins_pt_1 + iPt;
349d46ce 1220
c5704d5c 1221 if (_correctionWeight_1)
1222 corr = _correctionWeight_1[iZEtaPhiPt];
1223 else
1224 corr = 1;
1225 if (iZEtaPhiPt<0 || iZEtaPhiPt>=_nBins_zEtaPhiPt_1)
1226 {
1227 AliWarning("AliAnalysisTaskDptDptQA::analyze(AliceEvent * event) iZEtaPhiPt<0 || iZEtaPhiPt>=_nBins_zEtaPhiPt_1");
1228 continue;
1229 }
1230 //cout << "all good; process track:" << endl;
349d46ce 1231
c5704d5c 1232 if (_singlesOnly)
1233 {
1234 //_betaVsP_1->Fill(p*q,_trackFilter_1->getBeta() );
1235 //_dedxVsP_1->Fill(p*q,dedx);
1236 //_corrDedxVsP_1->Fill(p*q,_trackFilter_1->getCorrectedDedx() );
1237 __n1_1_vsPt[iPt] += corr; //cout << "step 15" << endl;
1238 __n1_1_vsZEtaPhiPt[iZEtaPhiPt] += corr; //cout << "step 12" << endl;
1239
1240 }
1241 else
1242 {
1243 corrPt = corr*pt;
1244 _id_1[k1] = iTrack;
1245 _charge_1[k1] = charge;
1246 _iEtaPhi_1[k1] = iEtaPhi;
1247 _iPt_1[k1] = iPt;
1248 _pt_1[k1] = pt;
1249 _px_1[k1] = px;
1250 _py_1[k1] = py;
1251 _pz_1[k1] = pz;
1252 _correction_1[k1] = corr;
1253 __n1_1 += corr;
1254 __n1_1_vsEtaPhi[iEtaPhi] += corr;
1255 __s1pt_1 += corrPt;
1256 __s1pt_1_vsEtaPhi[iEtaPhi] += corrPt;
1257 __n1Nw_1 += 1;
1258 __s1ptNw_1 += pt;
1259 ++k1;
1260 if (k1>=arraySize)
1261 {
1262 AliError(Form("AliAnalysisTaskDptDptQA::analyze(AliceEvent * event) k1 >=arraySize; arraySize: %d",arraySize));
1263 return;
1264 }
1265 }
1266 }
1267
80d059c6 1268 if (!_sameFilter && _requestedCharge_2 == charge &&
1269 eta >= -0.8 &&
1270 eta < 0.8)
1271 {
c5704d5c 1272 iPhi = int( phi/_width_phi_2);
1273 if (iPhi<0 || iPhi>=_nBins_phi_2 )
1274 {
1275 AliWarning("AliAnalysisTaskDptDptQA::analyze() iPhi<0 || iPhi>=_nBins_phi_1");
1276 return;
1277 }
349d46ce 1278
c5704d5c 1279 iEta = int((eta-_min_eta_2)/_width_eta_2);
1280 if (iEta<0 || iEta>=_nBins_eta_2)
1281 {
1282 AliWarning(Form("AliAnalysisTaskDptDptQA::analyze(AliceEvent * event) Mismatched iEta: %d", iEta));
1283 continue;
1284 }
1285 iPt = int((pt -_min_pt_2 )/_width_pt_2 );
1286 if (iPt<0 || iPt >=_nBins_pt_2)
1287 {
1288 AliWarning(Form("AliAnalysisTaskDptDptQA::analyze(AliceEvent * event) Mismatched iPt: %d",iPt));
1289 continue;
1290 }
349d46ce 1291
c5704d5c 1292 iEtaPhi = iEta*_nBins_phi_2+iPhi;
1293 iZEtaPhiPt = iVertexP2 + iEtaPhi*_nBins_pt_2 + iPt;
1294 if (iZEtaPhiPt<0 || iZEtaPhiPt>=_nBins_zEtaPhiPt_2)
349d46ce 1295 {
c5704d5c 1296 AliWarning("AliAnalysisTaskDptDptQA::analyze(AliceEvent * event) iZEtaPhiPt<0 || iZEtaPhiPt>=_nBins_zEtaPhiPt_2");
1297 continue;
1298 }
1299
1300 //cout << " iEtaPhi:" << iEtaPhi << " _nBins_etaPhi_1: "<< _nBins_etaPhi_1<< " _nBins_etaPhi_2: "<< _nBins_etaPhi_2<< endl;
1301 //if (_useEbyECorrections) corr = (charge>0) ? _correction_2p[iEtaPhi] : _correction_2m[iEtaPhi]; else
1302 //cout << "_correctionWeight_2:" << _correctionWeight_2 << endl;
1303 if (_correctionWeight_2)
1304 corr = _correctionWeight_2[iZEtaPhiPt];
1305 else
1306 corr = 1;
1307 //dpt = pt - (charge>0) ? _avgPt_vsEtaPhi_2p[iEtaPhi] : _avgPt_vsEtaPhi_2m[iEtaPhi];
1308
1309 if (_singlesOnly)
1310 {
1311 //_dedxVsP_2->Fill(p*q,dedx);
1312 __n1_2_vsPt[iPt] += corr; //cout << "step 15" << endl;
1313 __n1_2_vsZEtaPhiPt[iZEtaPhiPt] += corr; //cout << "step 12" << endl;
1314 }
1315 else
1316 {
1317 corrPt = corr*pt;
1318 _id_2[k2] = iTrack; //cout << "step 1" << endl;
1319 _charge_2[k2] = charge; //cout << "step 2" << endl;
1320 _iEtaPhi_2[k2] = iEtaPhi; //cout << "step 3" << endl;
1321 _iPt_2[k2] = iPt; //cout << "step 4" << endl;
1322 _pt_2[k2] = pt; //cout << "step 5" << endl;
1323 _px_2[k2] = px; //cout << "step 6" << endl;
1324 _py_2[k2] = py; //cout << "step 7" << endl;
1325 _pz_2[k2] = pz; //cout << "step 8" << endl;
1326 _correction_2[k2] = corr; //cout << "step 9" << endl;
1327 __n1_2 += corr; //cout << "step 10" << endl;
1328 __s1pt_2 += corrPt; //cout << "step 13" << endl;
1329 __n1Nw_2 += 1;
1330 __n1_2_vsEtaPhi[iEtaPhi] += corr; //cout << "step 11" << endl;
1331 __s1pt_2_vsEtaPhi[iEtaPhi] += corrPt; //cout << "step 14" << endl;
1332 __s1ptNw_2 += pt;
1333 ++k2;
1334 if (k2>=arraySize)
349d46ce 1335 {
c5704d5c 1336 AliWarning(Form("-W- k2 >=arraySize; arraySize: %d",arraySize));
1337 return;
349d46ce 1338 }
1339 }
c5704d5c 1340
1341 //cout << "done with track" << endl;
1342 } //iTrack
1343 } //aod
1344 }
1345
1346
1347 //cout << "Filling histograms now" << endl;
1348 _m0->Fill(_mult0);
1349 _m1->Fill(_mult1);
1350 _m2->Fill(_mult2);
1351 _m3->Fill(_mult3);
1352 _m4->Fill(_mult4);
1353 _m5->Fill(_mult5);
1354 _m6->Fill(_mult6);
1355 _vertexZ->Fill(vertexZ);
1356
1357 if (_singlesOnly)
1358 {
1359 // nothing to do here.
1360 }
1361 else
1362 {
1363 if (_sameFilter)
1364 {
1365 _n1_1_vsM->Fill(centrality, __n1_1);
1366 _s1pt_1_vsM->Fill(centrality, __s1pt_1);
1367 _n1Nw_1_vsM->Fill(centrality, __n1Nw_1);
1368 _s1ptNw_1_vsM->Fill(centrality, __s1ptNw_1);
1369 _n1_2_vsM->Fill(centrality, __n1_1);
1370 _s1pt_2_vsM->Fill(centrality, __s1pt_1);
1371 _n1Nw_2_vsM->Fill(centrality, __n1Nw_1);
1372 _s1ptNw_2_vsM->Fill(centrality, __s1ptNw_1);
1373 // reset pair counters
1374 __n2_12 = __s2ptpt_12 = __s2NPt_12 = __s2PtN_12 = 0;
1375 __n2Nw_12 = __s2ptptNw_12 = __s2NPtNw_12 = __s2PtNNw_12 = 0;
1376 if (_field>0)
1377 {
1378 for (int i1=0; i1<k1; i1++)
1379 {
1380 ////cout << " i1:" << i1 << endl;
1381 id_1 = _id_1[i1]; ////cout << " id_1:" << id_1 << endl;
1382 q_1 = _charge_1[i1]; ////cout << " q_1:" << q_1 << endl;
1383 iEtaPhi_1 = _iEtaPhi_1[i1]; ////cout << " iEtaPhi_1:" << iEtaPhi_1 << endl;
1384 iPt_1 = _iPt_1[i1]; ////cout << " iPt_1:" << iPt_1 << endl;
1385 corr_1 = _correction_1[i1]; ////cout << " corr_1:" << corr_1 << endl;
1386 pt_1 = _pt_1[i1]; ////cout << " pt_1:" << pt_1 << endl;
1387 dedx_1 = _dedx_1[i1]; ////cout << " dedx_1:" << dedx_1 << endl;
1388 //1 and 2
1389 for (int i2=i1+1; i2<k1; i2++)
1390 {
1391 ////cout << " i2:" << i2 << endl;
1392 id_2 = _id_1[i2]; ////cout << " id_2:" << id_2 << endl;
1393 if (id_1!=id_2)
1394 {
1395 q_2 = _charge_1[i2]; ////cout << " q_1:" << q_1 << endl;
1396 iEtaPhi_2 = _iEtaPhi_1[i2]; ////cout << " iEtaPhi_1:" << iEtaPhi_1 << endl;
1397 iPt_2 = _iPt_1[i2]; ////cout << " iPt_1:" << iPt_1 << endl;
1398 corr_2 = _correction_1[i2]; ////cout << " corr_1:" << corr_1 << endl;
1399 pt_2 = _pt_1[i2]; ////cout << " pt_1:" << pt_1 << endl;
1400 dedx_2 = _dedx_1[i2]; ////cout << " dedx_2:" << dedx_2 << endl;
1401 corr = corr_1*corr_2;
1402 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))
1403 {
1404 ij = iEtaPhi_1*_nBins_etaPhi_1 + iEtaPhi_2; ////cout << " ij:" << ij<< endl;
1405 }
1406 else // swap particles
1407 {
1408 ij = iEtaPhi_2*_nBins_etaPhi_1 + iEtaPhi_1; ////cout << " ij:" << ij<< endl;
1409 }
1410
1411 __n2_12 += corr;
1412 __n2_12_vsEtaPhi[ij] += corr;
1413 ptpt = pt_1*pt_2;
1414 __s2ptpt_12 += corr*ptpt;
1415 __s2PtN_12 += corr*pt_1;
1416 __s2NPt_12 += corr*pt_2;
1417 __s2ptpt_12_vsEtaPhi[ij] += corr*ptpt;
1418 __s2PtN_12_vsEtaPhi[ij] += corr*pt_1;
1419 __s2NPt_12_vsEtaPhi[ij] += corr*pt_2;
1420 __n2_12_vsPtPt[iPt_1*_nBins_pt_2 + iPt_2] += corr;
1421
1422 __n2Nw_12 += 1;
1423 __s2ptptNw_12 += ptpt;
1424 __s2PtNNw_12 += pt_1;
1425 __s2NPtNw_12 += pt_2;
1426
1427 }
1428 } //i2
1429 } //i1
1430 }
1431 else // field<0
1432 {
1433 for (int i1=0; i1<k1; i1++)
1434 {
1435 ////cout << " i1:" << i1 << endl;
1436 id_1 = _id_1[i1]; ////cout << " id_1:" << id_1 << endl;
1437 q_1 = _charge_1[i1]; ////cout << " q_1:" << q_1 << endl;
1438 iEtaPhi_1 = _iEtaPhi_1[i1]; ////cout << " iEtaPhi_1:" << iEtaPhi_1 << endl;
1439 iPt_1 = _iPt_1[i1]; ////cout << " iPt_1:" << iPt_1 << endl;
1440 corr_1 = _correction_1[i1]; ////cout << " corr_1:" << corr_1 << endl;
1441 pt_1 = _pt_1[i1]; ////cout << " pt_1:" << pt_1 << endl;
1442 dedx_1 = _dedx_1[i1]; ////cout << " dedx_1:" << dedx_1 << endl;
1443 //1 and 2
1444 for (int i2=i1+1; i2<k1; i2++)
1445 {
1446 ////cout << " i2:" << i2 << endl;
1447 id_2 = _id_1[i2]; ////cout << " id_2:" << id_2 << endl;
1448 if (id_1!=id_2)
1449 {
1450 q_2 = _charge_1[i2]; ////cout << " q_2:" << q_2 << endl;
1451 iEtaPhi_2 = _iEtaPhi_1[i2]; ////cout << " iEtaPhi_2:" << iEtaPhi_2 << endl;
1452 iPt_2 = _iPt_1[i2]; ////cout << " iPt_2:" << iPt_2 << endl;
1453 corr_2 = _correction_1[i2]; ////cout << " corr_2:" << corr_2 << endl;
1454 pt_2 = _pt_1[i2]; ////cout << " pt_2:" << pt_2 << endl;
1455 dedx_2 = _dedx_1[i2]; ////cout << " dedx_2:" << dedx_2 << endl;
1456 corr = corr_1*corr_2;
1457 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))
1458 {
1459 ij = iEtaPhi_1*_nBins_etaPhi_1 + iEtaPhi_2; ////cout << " ij:" << ij<< endl;
1460 }
1461 else // swap particles
1462 {
1463 ij = iEtaPhi_2*_nBins_etaPhi_1 + iEtaPhi_1; ////cout << " ij:" << ij<< endl;
1464 }
1465
1466 __n2_12 += corr;
1467 __n2_12_vsEtaPhi[ij] += corr;
1468 ptpt = pt_1*pt_2;
1469 __s2ptpt_12 += corr*ptpt;
1470 __s2PtN_12 += corr*pt_1;
1471 __s2NPt_12 += corr*pt_2;
1472 __s2ptpt_12_vsEtaPhi[ij] += corr*ptpt;
1473 __s2PtN_12_vsEtaPhi[ij] += corr*pt_1;
1474 __s2NPt_12_vsEtaPhi[ij] += corr*pt_2;
1475 __n2_12_vsPtPt[iPt_1*_nBins_pt_2 + iPt_2] += corr;
1476
1477 __n2Nw_12 += 1;
1478 __s2ptptNw_12 += ptpt;
1479 __s2PtNNw_12 += pt_1;
1480 __s2NPtNw_12 += pt_2;
1481
1482 }
1483 } //i2
1484 } //i1
1485 }
1486 }
1487 else // filter 1 and 2 are different -- must do all particle pairs...
1488 {
1489 _n1_1_vsM->Fill(centrality, __n1_1);
1490 _s1pt_1_vsM->Fill(centrality, __s1pt_1);
1491 _n1Nw_1_vsM->Fill(centrality, __n1Nw_1);
1492 _s1ptNw_1_vsM->Fill(centrality, __s1ptNw_1);
1493 _n1_2_vsM->Fill(centrality, __n1_2);
1494 _s1pt_2_vsM->Fill(centrality, __s1pt_2);
1495 _n1Nw_2_vsM->Fill(centrality, __n1Nw_2);
1496 _s1ptNw_2_vsM->Fill(centrality, __s1ptNw_2);
1497 // reset pair counters
1498 __n2_12 = __s2ptpt_12 = __s2NPt_12 = __s2PtN_12 = 0;
1499 __n2Nw_12 = __s2ptptNw_12 = __s2NPtNw_12 = __s2PtNNw_12 = 0;
1500 for (int i1=0; i1<k1; i1++)
1501 {
1502 ////cout << " i1:" << i1 << endl;
1503 id_1 = _id_1[i1]; ////cout << " id_1:" << id_1 << endl;
1504 q_1 = _charge_1[i1]; ////cout << " q_1:" << q_1 << endl;
1505 iEtaPhi_1 = _iEtaPhi_1[i1]; ////cout << " iEtaPhi_1:" << iEtaPhi_1 << endl;
1506 iPt_1 = _iPt_1[i1]; ////cout << " iPt_1:" << iPt_1 << endl;
1507 corr_1 = _correction_1[i1]; ////cout << " corr_1:" << corr_1 << endl;
1508 pt_1 = _pt_1[i1]; ////cout << " pt_1:" << pt_1 << endl;
1509 px_1 = _px_1[i1]; ////cout << " px_1:" << px_1 << endl;
1510 py_1 = _py_1[i1]; ////cout << " py_1:" << py_1 << endl;
1511 pz_1 = _pz_1[i1]; ////cout << " pz_1:" << pz_1 << endl;
1512 dedx_1 = _dedx_1[i1]; ////cout << " dedx_1:" << dedx_1 << endl;
1513
1514 //1 and 2
1515 for (int i2=0; i2<k2; i2++)
1516 {
1517 ////cout << " i2:" << i2 << endl;
1518 id_2 = _id_2[i2]; ////cout << " id_2:" << id_2 << endl;
1519 if (id_1!=id_2) // exclude auto correlation
1520 {
1521 q_2 = _charge_2[i2]; ////cout << " q_2:" << q_2 << endl;
1522 iEtaPhi_2 = _iEtaPhi_2[i2]; ////cout << " iEtaPhi_2:" << iEtaPhi_2 << endl;
1523 iPt_2 = _iPt_2[i2]; ////cout << " iPt_2:" << iPt_2 << endl;
1524 corr_2 = _correction_2[i2]; ////cout << " corr_2:" << corr_2 << endl;
1525 pt_2 = _pt_2[i2]; ////cout << " pt_2:" << pt_2 << endl;
1526 px_2 = _px_2[i2]; ////cout << " px_2:" << px_2 << endl;
1527 py_2 = _py_2[i2]; ////cout << " py_2:" << py_2 << endl;
1528 pz_2 = _pz_2[i2]; ////cout << " pz_2:" << pz_2 << endl;
1529 dedx_2 = _dedx_2[i2]; ////cout << " dedx_2:" << dedx_2 << endl;
1530 if (_rejectPairConversion)
1531 {
1532 float e1Sq = massElecSq + pt_1*pt_1 + pz_1*pz_1;
1533 float e2Sq = massElecSq + pt_2*pt_2 + pz_2*pz_2;
1534 float mInvSq = 2*(massElecSq + sqrt(e1Sq*e2Sq) - px_1*px_2 - py_1*py_2 - pz_1*pz_2 );
1535 float mInv = sqrt(mInvSq);
1536 _invMass->Fill(mInv);
1537 if (mInv<0.5)
1538 {
1539 if (dedx_1>75. && dedx_2>75.)
1540 {
1541 _invMassElec->Fill(mInv);
1542 if (mInv<0.05) continue;
1543 }
1544 }
1545 }
1546 corr = corr_1*corr_2;
1547 ij = iEtaPhi_1*_nBins_etaPhi_1 + iEtaPhi_2; ////cout << " ij:" << ij<< endl;
1548 __n2_12 += corr;
1549 __n2_12_vsEtaPhi[ij] += corr;
1550 ptpt = pt_1*pt_2;
1551 __s2ptpt_12 += corr*ptpt;
1552 __s2PtN_12 += corr*pt_1;
1553 __s2NPt_12 += corr*pt_2;
1554 __s2ptpt_12_vsEtaPhi[ij] += corr*ptpt;
1555 __s2PtN_12_vsEtaPhi[ij] += corr*pt_1;
1556 __s2NPt_12_vsEtaPhi[ij] += corr*pt_2;
1557 __n2_12_vsPtPt[iPt_1*_nBins_pt_2 + iPt_2] += corr;
1558 __n2Nw_12 += 1;
1559 __s2ptptNw_12 += ptpt;
1560 __s2PtNNw_12 += pt_1;
1561 __s2NPtNw_12 += pt_2;
1562
1563 }
1564 } //i2
1565 } //i1
1566 }
1567
1568 _n2_12_vsM->Fill(centrality, __n2_12);
1569 _s2PtPt_12_vsM->Fill(centrality, __s2ptpt_12);
1570 _s2PtN_12_vsM->Fill(centrality, __s2NPt_12);
1571 _s2NPt_12_vsM->Fill(centrality, __s2PtN_12);
1572
1573 _n2Nw_12_vsM->Fill(centrality, __n2Nw_12);
1574 _s2PtPtNw_12_vsM->Fill(centrality, __s2ptptNw_12);
1575 _s2PtNNw_12_vsM->Fill(centrality, __s2NPtNw_12);
1576 _s2NPtNw_12_vsM->Fill(centrality, __s2PtNNw_12);
1577
1578 }
1579
1580
1581 cout << "AliAnalysisTaskDptDptQA::UserExec() -----------------Event Done " << endl;
1582 PostData(0,_outputHistoList);
349d46ce 1583
349d46ce 1584}
1585
1586void AliAnalysisTaskDptDptQA::FinishTaskOutput()
1587{
c5704d5c 1588 cout << "AliAnalysisTaskDptDptQA::FinishTaskOutput() Starting." << endl;
1589 cout << "= 0 ====================================================================" << endl;
349d46ce 1590 finalizeHistograms();
c5704d5c 1591 cout << "= 1 ====================================================================" << endl;
349d46ce 1592 PostData(0,_outputHistoList);
c5704d5c 1593 cout << "= 2 ====================================================================" << endl;
1594 cout << "AliAnalysisTaskDptDptQA::FinishTaskOutput() Done." << endl;
349d46ce 1595}
1596
1597void AliAnalysisTaskDptDptQA::Terminate(Option_t* /*option*/)
1598{
c5704d5c 1599 cout << "AliAnalysisTaskDptDptQA::Terminate() Starting/Done." << endl;
349d46ce 1600}
1601
1602
1603//Tools
1604//===================================================================================================
1605void AliAnalysisTaskDptDptQA::fillHistoWithArray(TH1 * h, float * array, int size)
1606{
1607 int i, i1;
1608 float v1, ev1, v2, ev2, sum, esum;
1609 for (i=0, i1=1; i<size; ++i,++i1)
1610 {
1611 v1 = array[i]; ev1 = sqrt(v1);
1612 v2 = h->GetBinContent(i1);
1613 ev2 = h->GetBinError(i1);
1614 sum = v1 + v2;
1615 esum = sqrt(ev1*ev1+ev2*ev2);
1616 h->SetBinContent(i1,sum);
1617 h->SetBinError(i1,esum);
1618 }
1619}
1620
1621void AliAnalysisTaskDptDptQA::fillHistoWithArray(TH2 * h, float * array, int size1, int size2)
1622{
1623 int i, i1;
1624 int j, j1;
1625 float v1, ev1, v2, ev2, sum, esum;
1626 for (i=0, i1=1; i<size1; ++i,++i1)
1627 {
1628 for (j=0, j1=1; j<size2; ++j,++j1)
1629 {
1630 v1 = array[i*size2+j]; ev1 = sqrt(v1);
1631 v2 = h->GetBinContent(i1,j1);
1632 ev2 = h->GetBinError(i1,j1);
1633 sum = v1 + v2;
1634 esum = sqrt(ev1*ev1+ev2*ev2);
1635 h->SetBinContent(i1,j1,sum);
1636 h->SetBinError(i1,j1,esum);
1637 }
1638 }
1639}
1640
1641void AliAnalysisTaskDptDptQA::fillHistoWithArray(TH3 * h, float * array, int size1, int size2, int size3)
1642{
1643 int i, i1;
1644 int j, j1;
1645 int k, k1;
1646 float v1, ev1, v2, ev2, sum, esum;
1647 int size23 = size2*size3;
1648 for (i=0, i1=1; i<size1; ++i,++i1)
1649 {
1650 for (j=0, j1=1; j<size2; ++j,++j1)
1651 {
1652 for (k=0, k1=1; k<size3; ++k,++k1)
1653 {
1654 v1 = array[i*size23+j*size3+k]; ev1 = sqrt(v1);
1655 v2 = h->GetBinContent(i1,j1,k1);
1656 ev2 = h->GetBinError(i1,j1,k1);
1657 sum = v1 + v2;
1658 esum = sqrt(ev1*ev1+ev2*ev2);
1659 h->SetBinContent(i1,j1,k1,sum);
1660 h->SetBinError(i1,j1,k1,esum);
1661 }
1662 }
1663 }
1664}
1665
1666void AliAnalysisTaskDptDptQA::fillHistoWithArray(TH1 * h, double * array, int size)
1667{
1668 int i, i1;
1669 double v1, ev1, v2, ev2, sum, esum;
1670 for (i=0, i1=1; i<size; ++i,++i1)
1671 {
1672 v1 = array[i]; ev1 = sqrt(v1);
1673 v2 = h->GetBinContent(i1);
1674 ev2 = h->GetBinError(i1);
1675 sum = v1 + v2;
1676 esum = sqrt(ev1*ev1+ev2*ev2);
1677 h->SetBinContent(i1,sum);
1678 h->SetBinError(i1,esum);
1679 }
1680}
1681
1682void AliAnalysisTaskDptDptQA::fillHistoWithArray(TH2 * h, double * array, int size1, int size2)
1683{
1684 int i, i1;
1685 int j, j1;
1686 double v1, ev1, v2, ev2, sum, esum;
1687 for (i=0, i1=1; i<size1; ++i,++i1)
1688 {
1689 for (j=0, j1=1; j<size2; ++j,++j1)
1690 {
1691 v1 = array[i*size2+j]; ev1 = sqrt(v1);
1692 v2 = h->GetBinContent(i1,j1);
1693 ev2 = h->GetBinError(i1,j1);
1694 sum = v1 + v2;
1695 esum = sqrt(ev1*ev1+ev2*ev2);
1696 h->SetBinContent(i1,j1,sum);
1697 h->SetBinError(i1,j1,esum);
1698 }
1699 }
1700}
1701
1702void AliAnalysisTaskDptDptQA::fillHistoWithArray(TH3 * h, double * array, int size1, int size2, int size3)
1703{
1704 int i, i1;
1705 int j, j1;
1706 int k, k1;
1707 double v1, ev1, v2, ev2, sum, esum;
1708 int size23 = size2*size3;
1709 for (i=0, i1=1; i<size1; ++i,++i1)
1710 {
1711 for (j=0, j1=1; j<size2; ++j,++j1)
1712 {
1713 for (k=0, k1=1; k<size3; ++k,++k1)
1714 {
1715 v1 = array[i*size23+j*size3+k]; ev1 = sqrt(v1);
1716 v2 = h->GetBinContent(i1,j1,k1);
1717 ev2 = h->GetBinError(i1,j1,k1);
1718 sum = v1 + v2;
1719 esum = sqrt(ev1*ev1+ev2*ev2);
1720 h->SetBinContent(i1,j1,k1,sum);
1721 h->SetBinError(i1,j1,k1,esum);
1722 }
1723 }
1724 }
1725}
1726
1727//________________________________________________________________________
1728double * AliAnalysisTaskDptDptQA::getDoubleArray(int size, double v)
1729{
1730 /// Allocate an array of type double with n values
1731 /// Initialize the array to the given value
1732 double * array = new double [size];
1733 for (int i=0;i<size;++i) array[i]=v;
1734 return array;
1735}
1736
1737//________________________________________________________________________
1738float * AliAnalysisTaskDptDptQA::getFloatArray(int size, float v)
1739{
1740 /// Allocate an array of type float with n values
1741 /// Initialize the array to the given value
1742 float * array = new float [size];
1743 for (int i=0;i<size;++i) array[i]=v;
1744 return array;
1745}
1746
1747
1748//________________________________________________________________________
1749TH1D * AliAnalysisTaskDptDptQA::createHisto1D(const TString & name, const TString & title,
1750 int n, double xMin, double xMax,
1751 const TString & xTitle, const TString & yTitle)
1752{
1753 //CreateHisto new 1D historgram
1754 AliInfo(Form("createHisto 1D histo %s nBins: %d xMin: %f xMax: %f",name.Data(),n,xMin,xMax));
1755 TH1D * h = new TH1D(name,title,n,xMin,xMax);
1756 h->GetXaxis()->SetTitle(xTitle);
1757 h->GetYaxis()->SetTitle(yTitle);
1758 addToList(h);
1759 return h;
1760}
1761
1762
1763//________________________________________________________________________
1764TH1D * AliAnalysisTaskDptDptQA::createHisto1D(const TString & name, const TString & title,
1765 int n, double * bins,
1766 const TString & xTitle, const TString & yTitle)
1767{
1768 AliInfo(Form("createHisto 1D histo %s with %d non uniform nBins",name.Data(),n));
1769 TH1D * h = new TH1D(name,title,n,bins);
1770 h->GetXaxis()->SetTitle(xTitle);
1771 h->GetYaxis()->SetTitle(yTitle);
1772 addToList(h);
1773 return h;
1774}
1775
1776
1777//________________________________________________________________________
1778TH2D * AliAnalysisTaskDptDptQA::createHisto2D(const TString & name, const TString & title,
1779 int nx, double xMin, double xMax, int ny, double yMin, double yMax,
1780 const TString & xTitle, const TString & yTitle, const TString & zTitle)
1781{
1782 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));
1783 TH2D * h = new TH2D(name,title,nx,xMin,xMax,ny,yMin,yMax);
1784 h->GetXaxis()->SetTitle(xTitle);
1785 h->GetYaxis()->SetTitle(yTitle);
1786 h->GetZaxis()->SetTitle(zTitle);
1787 addToList(h);
1788 return h;
1789}
1790
1791//________________________________________________________________________
1792TH2D * AliAnalysisTaskDptDptQA::createHisto2D(const TString & name, const TString & title,
1793 int nx, double* xbins, int ny, double yMin, double yMax,
1794 const TString & xTitle, const TString & yTitle, const TString & zTitle)
1795{
1796 AliInfo(Form("createHisto 2D histo %s with %d non uniform nBins",name.Data(),nx));
1797 TH2D * h;
1798 h = new TH2D(name,title,nx,xbins,ny,yMin,yMax);
1799 h->GetXaxis()->SetTitle(xTitle);
1800 h->GetYaxis()->SetTitle(yTitle);
1801 h->GetZaxis()->SetTitle(zTitle);
1802 addToList(h);
1803 return h;
1804}
1805
1806//// F /////
1807//________________________________________________________________________
1808TH1F * AliAnalysisTaskDptDptQA::createHisto1F(const TString & name, const TString & title,
1809 int n, double xMin, double xMax,
1810 const TString & xTitle, const TString & yTitle)
1811{
1812 //CreateHisto new 1D historgram
1813 AliInfo(Form("createHisto 1D histo %s nBins: %d xMin: %f xMax: %f",name.Data(),n,xMin,xMax));
1814 TH1F * h = new TH1F(name,title,n,xMin,xMax);
1815 h->GetXaxis()->SetTitle(xTitle);
1816 h->GetYaxis()->SetTitle(yTitle);
1817 addToList(h);
1818 return h;
1819}
1820
1821
1822//________________________________________________________________________
1823TH1F * AliAnalysisTaskDptDptQA::createHisto1F(const TString & name, const TString & title,
1824 int n, double * bins,
1825 const TString & xTitle, const TString & yTitle)
1826{
1827 AliInfo(Form("createHisto 1D histo %s with %d non uniform nBins",name.Data(),n));
1828 TH1F * h = new TH1F(name,title,n,bins);
1829 h->GetXaxis()->SetTitle(xTitle);
1830 h->GetYaxis()->SetTitle(yTitle);
1831 addToList(h);
1832 return h;
1833}
1834
1835
1836//________________________________________________________________________
1837TH2F * AliAnalysisTaskDptDptQA::createHisto2F(const TString & name, const TString & title,
1838 int nx, double xMin, double xMax, int ny, double yMin, double yMax,
1839 const TString & xTitle, const TString & yTitle, const TString & zTitle)
1840{
1841 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));
1842 TH2F * h = new TH2F(name,title,nx,xMin,xMax,ny,yMin,yMax);
1843 h->GetXaxis()->SetTitle(xTitle);
1844 h->GetYaxis()->SetTitle(yTitle);
1845 h->GetZaxis()->SetTitle(zTitle);
1846 addToList(h);
1847 return h;
1848}
1849
1850//________________________________________________________________________
1851TH2F * AliAnalysisTaskDptDptQA::createHisto2F(const TString & name, const TString & title,
1852 int nx, double* xbins, int ny, double yMin, double yMax,
1853 const TString & xTitle, const TString & yTitle, const TString & zTitle)
1854{
1855 AliInfo(Form("createHisto 2D histo %s with %d non uniform nBins",name.Data(),nx));
1856 TH2F * h;
1857 h = new TH2F(name,title,nx,xbins,ny,yMin,yMax);
1858 h->GetXaxis()->SetTitle(xTitle);
1859 h->GetYaxis()->SetTitle(yTitle);
1860 h->GetZaxis()->SetTitle(zTitle);
1861 addToList(h);
1862 return h;
1863}
1864
1865
1866//________________________________________________________________________
1867TH3F * AliAnalysisTaskDptDptQA::createHisto3F(const TString & name, const TString & title,
1868 int nx, double xMin, double xMax,
1869 int ny, double yMin, double yMax,
1870 int nz, double zMin, double zMax,
1871 const TString & xTitle, const TString & yTitle, const TString & zTitle)
1872{
1873 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));
1874 TH3F * h = new TH3F(name,title,nx,xMin,xMax,ny,yMin,yMax,nz,zMin,zMax);
1875 h->GetXaxis()->SetTitle(xTitle);
1876 h->GetYaxis()->SetTitle(yTitle);
1877 h->GetZaxis()->SetTitle(zTitle);
1878 addToList(h);
1879 return h;
1880}
1881
1882
1883//________________________________________________________________________
1884TProfile * AliAnalysisTaskDptDptQA::createProfile(const TString & name, const TString & description,
1885 int nx,double xMin,double xMax,
1886 const TString & xTitle, const TString & yTitle)
1887{
1888 AliInfo(Form("createHisto 1D profile %s nBins: %d xMin: %f10.4 xMax: %f10.4",name.Data(),nx,xMin,xMax));
1889 TProfile * h = new TProfile(name,description,nx,xMin,xMax);
1890 h->GetXaxis()->SetTitle(xTitle);
1891 h->GetYaxis()->SetTitle(yTitle);
1892 addToList(h);
1893 return h;
1894}
1895
1896//________________________________________________________________________
1897TProfile * AliAnalysisTaskDptDptQA::createProfile(const TString & name,const TString & description,
1898 int nx, double* bins,
1899 const TString & xTitle, const TString & yTitle)
1900{
1901 AliInfo(Form("createHisto 1D profile %s with %d non uniform bins",name.Data(),nx));
1902 TProfile * h = new TProfile(name,description,nx,bins);
1903 h->GetXaxis()->SetTitle(xTitle);
1904 h->GetYaxis()->SetTitle(yTitle);
1905 addToList(h);
1906 return h;
1907}
1908
1909
1910void AliAnalysisTaskDptDptQA::addToList(TH1 *h)
1911{
1912 if (_outputHistoList)
1913 {
1914 _outputHistoList->Add(h);
1915 }
1916 else
c5704d5c 1917 cout << "addToList(TH1 *h) _outputHistoList is null!!!!! Shoudl abort ship" << endl;
349d46ce 1918
1919}
1920
1921
1922