]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGLF/SPECTRA/PiKaPr/TestAOD/AliAnalysisTaskV2AllChAOD.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Merge branch 'feature-movesplit'
[u/mrichter/AliRoot.git] / PWGLF / SPECTRA / PiKaPr / TestAOD / AliAnalysisTaskV2AllChAOD.cxx
CommitLineData
10a99a07 1/**************************************************************************
2 * Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16//-----------------------------------------------------------------
17// AliAnalysisTaskV2AllChAOD class
18//-----------------------------------------------------------------
19
20#include "TChain.h"
21#include "TTree.h"
22#include "TLegend.h"
23#include "TH1F.h"
d0761d58 24#include "TH1D.h"
10a99a07 25#include "TH2F.h"
26#include "THnSparse.h"
27#include "TProfile.h"
28#include "TCanvas.h"
29#include "AliAnalysisTask.h"
30#include "AliAODTrack.h"
31#include "AliAODMCParticle.h"
32#include "AliVParticle.h"
33#include "AliAODEvent.h"
34#include "AliAODInputHandler.h"
35#include "AliAnalysisTaskV2AllChAOD.h"
36#include "AliAnalysisTaskESDfilter.h"
37#include "AliAnalysisDataContainer.h"
38#include "AliSpectraAODTrackCuts.h"
39#include "AliSpectraAODEventCuts.h"
40#include "AliPIDCombined.h"
41#include "AliCentrality.h"
42#include "TProof.h"
43#include "AliVEvent.h"
44#include "AliStack.h"
45#include <TMCProcess.h>
bc9bd5b3 46#include <TRandom.h>
10a99a07 47
48#include <iostream>
49
50using namespace std;
51
52ClassImp(AliAnalysisTaskV2AllChAOD)
53
54//________________________________________________________________________
55AliAnalysisTaskV2AllChAOD::AliAnalysisTaskV2AllChAOD(const char *name) : AliAnalysisTaskSE(name),
8aa04ea0 56fAOD(0x0),
57fTrackCuts(0x0),
58fEventCuts(0x0),
59fIsMC(0),
60fCharge(0),
61fVZEROside(0),
62fOutput(0x0),
63fOutput_lq(0x0),
64fOutput_sq(0x0),
65fnCentBins(20),
66fnQvecBins(100),
67fQvecUpperLim(100),
68fCutLargeQperc(90.),
69fCutSmallQperc(10.),
70fEtaGapMin(-0.5),
71fEtaGapMax(0.5),
72fTrkBit(128),
73fEtaCut(0.8),
74fMinPt(0),
75fMaxPt(20.0),
76fMinTPCNcls(70),
77fFillTHn(kFALSE),
78fCentrality(0),
79fQvector(0),
80fQvector_lq(0),
81fQvector_sq(0),
82fResSP(0),
83fResSP_vs_Cent(0),
84f2partCumQA_vs_Cent(0),
85f2partCumQB_vs_Cent(0),
86fEta_vs_Phi_bef(0),
87fEta_vs_PhiA(0),
88fEta_vs_PhiB(0),
89fResSP_lq(0),
90fResSP_vs_Cent_lq(0),
91f2partCumQA_vs_Cent_lq(0),
92f2partCumQB_vs_Cent_lq(0),
93fResSP_sq(0),
94fResSP_vs_Cent_sq(0),
95f2partCumQA_vs_Cent_sq(0),
96f2partCumQB_vs_Cent_sq(0),
97fResSP_inclusive(0),
98fv2SPGap1A_inclusive_mb(0),
99fv2SPGap1B_inclusive_mb(0),
100fv2SPGap1A_inclusive_lq(0),
101fv2SPGap1B_inclusive_lq(0),
102fv2SPGap1A_inclusive_sq(0),
103fv2SPGap1B_inclusive_sq(0),
104fResSPmc_inclusive(0),
105fv2SPGap1Amc_inclusive_mb(0),
106fv2SPGap1Bmc_inclusive_mb(0),
107fv2SPGap1Amc_inclusive_lq(0),
108fv2SPGap1Bmc_inclusive_lq(0),
109fv2SPGap1Amc_inclusive_sq(0),
110fv2SPGap1Bmc_inclusive_sq(0),
111fResGap1w(0),
112fV2IntGap1w(0),
113fIsRecoEff(0),
114fRecoEffList(0),
115fQvecGen(0),
116fQgenType(0),
117fnNchBins(400),
118fDoCentrSystCentrality(0)
10a99a07 119{
8aa04ea0 120
10a99a07 121 for (Int_t i = 0; i< 9; i++){
8aa04ea0 122
10a99a07 123 fv2SPGap1A[i] = 0;
10a99a07 124 fv2SPGap1B[i] = 0;
bc76879c 125
126 fSinGap1Aq[i] = 0;
127 fCosGap1Aq[i] = 0;
128 fSinGap1Bq[i] = 0;
129 fCosGap1Bq[i] = 0;
10a99a07 130
131 fSinGap1A[i] = 0;
132 fCosGap1A[i] = 0;
133 fSinGap1B[i] = 0;
134 fCosGap1B[i] = 0;
8aa04ea0 135
10a99a07 136 //large q
137 fv2SPGap1A_lq[i] = 0;
138 fv2SPGap1B_lq[i] = 0;
8aa04ea0 139
bc76879c 140 fSinGap1Aq_lq[i] = 0;
141 fCosGap1Aq_lq[i] = 0;
142 fSinGap1Bq_lq[i] = 0;
143 fCosGap1Bq_lq[i] = 0;
8aa04ea0 144
10a99a07 145 fSinGap1A_lq[i] = 0;
146 fCosGap1A_lq[i] = 0;
147 fSinGap1B_lq[i] = 0;
148 fCosGap1B_lq[i] = 0;
8aa04ea0 149
10a99a07 150 //small q
151 fv2SPGap1A_sq[i] = 0;
152 fv2SPGap1B_sq[i] = 0;
8aa04ea0 153
bc76879c 154 fSinGap1Aq_sq[i] = 0;
155 fCosGap1Aq_sq[i] = 0;
156 fSinGap1Bq_sq[i] = 0;
157 fCosGap1Bq_sq[i] = 0;
8aa04ea0 158
10a99a07 159 fSinGap1A_sq[i] = 0;
160 fCosGap1A_sq[i] = 0;
161 fSinGap1B_sq[i] = 0;
162 fCosGap1B_sq[i] = 0;
8aa04ea0 163
52bfa7e6 164 fResSP_vs_Qvec[i] = 0;
165 fV2IntGap1wq[i] = 0;
8aa04ea0 166
10a99a07 167 }
8aa04ea0 168
d0761d58 169 fRecoEffList=new TList();
170 fRecoEffList->SetOwner();
171 fRecoEffList->SetName("fRecoEffList");
8aa04ea0 172
10a99a07 173 // Default constructor
174 DefineInput(0, TChain::Class());
175 DefineOutput(1, TList::Class());
176 DefineOutput(2, AliSpectraAODEventCuts::Class());
177 DefineOutput(3, AliSpectraAODTrackCuts::Class());
178 DefineOutput(4, TList::Class());
179 DefineOutput(5, TList::Class());
180}
181
182//________________________________________________________________________
183void AliAnalysisTaskV2AllChAOD::UserCreateOutputObjects()
184{
185 // create output objects
186 fOutput=new TList();
187 fOutput->SetOwner();
188 fOutput->SetName("fOutput");
8aa04ea0 189
10a99a07 190 fOutput_lq=new TList();
191 fOutput_lq->SetOwner();
192 fOutput_lq->SetName("fOutput_lq");
8aa04ea0 193
10a99a07 194 fOutput_sq=new TList();
195 fOutput_sq->SetOwner();
196 fOutput_sq->SetName("fOutput_sq");
8aa04ea0 197
10a99a07 198 if (!fTrackCuts) AliFatal("Track Cuts should be set in the steering macro");
199 if (!fEventCuts) AliFatal("Event Cuts should be set in the steering macro");
8aa04ea0 200
962f3d11 201 if( fFillTHn ){
202 //dimensions of THnSparse for Q vector checks
34f4a5fc 203 const Int_t nvarev=6;
1ba2a775 204 // cent q-rec_perc qvec_v0a q-rec_v0c qvec-gen_tpc Nch
34f4a5fc 205 Int_t binsHistRealEv[nvarev] = { fnCentBins, 100, fnQvecBins, fnQvecBins, fnQvecBins, fnNchBins};
206 Double_t xminHistRealEv[nvarev] = { 0., 0., 0., 0., 0., 0.};
207 Double_t xmaxHistRealEv[nvarev] = { 100., 100., fQvecUpperLim, fQvecUpperLim, fQvecUpperLim, 2000.};
8aa04ea0 208
962f3d11 209 THnSparseF* NSparseHistEv = new THnSparseF("NSparseHistEv","NSparseHistEv",nvarev,binsHistRealEv,xminHistRealEv,xmaxHistRealEv);
210 NSparseHistEv->GetAxis(0)->SetTitle(Form("%s cent",fEventCuts->GetCentralityMethod().Data()));
211 NSparseHistEv->GetAxis(0)->SetName(Form("%s_cent",fEventCuts->GetCentralityMethod().Data()));
8aa04ea0 212
e3bcd147 213 NSparseHistEv->GetAxis(1)->SetTitle("q-vec rec percentile");
214 NSparseHistEv->GetAxis(1)->SetName("Qrec_perc");
8aa04ea0 215
1ba2a775 216 NSparseHistEv->GetAxis(2)->SetTitle("q-vec V0A");
217 NSparseHistEv->GetAxis(2)->SetName("Qrec_V0A");
8aa04ea0 218
1ba2a775 219 NSparseHistEv->GetAxis(3)->SetTitle("q-vec V0C");
220 NSparseHistEv->GetAxis(3)->SetName("Qrec_V0C");
8aa04ea0 221
1ba2a775 222 NSparseHistEv->GetAxis(4)->SetTitle("q-vec TPC");
223 NSparseHistEv->GetAxis(4)->SetName("Qgen_TPC");
8aa04ea0 224
e3bcd147 225 NSparseHistEv->GetAxis(5)->SetTitle("Ncharged");
226 NSparseHistEv->GetAxis(5)->SetName("Nch");
962f3d11 227 fOutput->Add(NSparseHistEv);
228 }
8aa04ea0 229
29fbaf8a 230 fCentrality = new TH1D("fCentrality", "centrality distribution; centrality", 200, 0., 100);
231 fOutput->Add(fCentrality);
8aa04ea0 232
86f0713c 233 fQvector = new TH1D("fQvector", "q-vector distribution; q-vector", fnQvecBins, 0., fQvecUpperLim);
234 fOutput->Add(fQvector);
8aa04ea0 235
86f0713c 236 fQvector_lq = new TH1D("fQvector_lq", "q-vector distribution; q-vector", fnQvecBins, 0., fQvecUpperLim);
237 fOutput_lq->Add(fQvector_lq);
8aa04ea0 238
86f0713c 239 fQvector_sq = new TH1D("fQvector_sq", "q-vector distribution; q-vector", fnQvecBins, 0., fQvecUpperLim);
240 fOutput_sq->Add(fQvector_sq);
8aa04ea0 241
10a99a07 242 // binning common to all the THn
243 //change it according to your needs + move it to global variables -> setter/getter
8aa04ea0 244 // Double_t ptBins[] = {0., 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0, 1.2, 1.4, 1.6, 1.8, 2.0, 2.4, 2.8, 3.2, 3.6, 4.0, 4.5, 5.0, 5.5, 6.0, 7.0, 8.0, 9.0, 10.0, 12.0, 14.0, 16.0, 20.0};
245 // const Int_t nptBins = 31;
10a99a07 246 Double_t ptBins[] = {0., 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0, 1.2, 1.4, 1.6, 1.8, 2.0, 2.2, 2.4, 2.6, 3.0, 3.4, 3.8, 4.2, 4.6, 5.0, 5.5, 6.0, 7.0, 8.0, 9.0, 10.0, 12.0, 14.0, 16.0, 20.0};
247 const Int_t nptBins = 33;
8aa04ea0 248
10a99a07 249 fResSP = new TProfile("fResSP", "Resolution; centrality; Resolution", 9, -0.5, 8.5);
250 fOutput->Add(fResSP);
8aa04ea0 251
30e20fad 252 fResSP_vs_Cent = new TProfile("fResSP_vs_Cent", "Resolution; centrality; Resolution", 20., 0., 100.);
af960b8a 253 fOutput->Add(fResSP_vs_Cent);
8aa04ea0 254
30e20fad 255 f2partCumQA_vs_Cent = new TProfile("f2partCumQA_vs_Cent", "Resolution; centrality; Resolution", 100., 0., 100.);
af960b8a 256 fOutput->Add(f2partCumQA_vs_Cent);
8aa04ea0 257
30e20fad 258 f2partCumQB_vs_Cent = new TProfile("f2partCumQB_vs_Cent", "Resolution; centrality; Resolution", 100., 0., 100.);
af960b8a 259 fOutput->Add(f2partCumQB_vs_Cent);
260
30e20fad 261 fEta_vs_Phi_bef = new TH2D("fEta_vs_Phi_bef","eta vs phi distribution before eta gap;#eta;#phi",200.,-1.,1.,175.,0.,7.);
bc76879c 262 fOutput->Add(fEta_vs_Phi_bef);
8aa04ea0 263
30e20fad 264 fEta_vs_PhiA = new TH2D("fEta_vs_PhiA","eta vs phi distribution;#eta;#phi",200.,-1.,1.,175.,0.,7.);
2ba0e068 265 fOutput->Add(fEta_vs_PhiA);
8aa04ea0 266
30e20fad 267 fEta_vs_PhiB = new TH2D("fEta_vs_PhiB","eta vs phi distribution;#eta;#phi",200.,-1.,1.,175.,0.,7.);
2ba0e068 268 fOutput->Add(fEta_vs_PhiB);
8aa04ea0 269
270 // MC closure test
271 fResSP_inclusive = new TProfile("fResSP_inclusive", "Resolution; ese; Resolution", 3, 0., 3.);
272 fOutput->Add(fResSP_inclusive);
273
274 fv2SPGap1A_inclusive_mb = new TProfile("fv2SPGap1A_inclusive_mb", "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
275 fOutput->Add(fv2SPGap1A_inclusive_mb);
276
277 fv2SPGap1B_inclusive_mb = new TProfile("fv2SPGap1B_inclusive_mb", "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
278 fOutput->Add(fv2SPGap1B_inclusive_mb);
279
280
29fbaf8a 281 //large q
10a99a07 282 fResSP_lq = new TProfile("fResSP_lq", "Resolution; centrality; Resolution", 9, -0.5, 8.5);
283 fOutput_lq->Add(fResSP_lq);
8aa04ea0 284
30e20fad 285 fResSP_vs_Cent_lq = new TProfile("fResSP_vs_Cent_lq", "Resolution; centrality; Resolution", 20., 0., 100.);
af960b8a 286 fOutput_lq->Add(fResSP_vs_Cent_lq);
8aa04ea0 287
30e20fad 288 f2partCumQA_vs_Cent_lq = new TProfile("f2partCumQA_vs_Cent_lq", "Resolution; centrality; Resolution", 100., 0., 100.);
af960b8a 289 fOutput_lq->Add(f2partCumQA_vs_Cent_lq);
8aa04ea0 290
30e20fad 291 f2partCumQB_vs_Cent_lq = new TProfile("f2partCumQB_vs_Cent_lq", "Resolution; centrality; Resolution", 100., 0., 100.);
af960b8a 292 fOutput_lq->Add(f2partCumQB_vs_Cent_lq);
8aa04ea0 293
294 // MC closure test
295 fv2SPGap1A_inclusive_lq = new TProfile("fv2SPGap1A_inclusive_lq", "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
296 fOutput_lq->Add(fv2SPGap1A_inclusive_lq);
297
298 fv2SPGap1B_inclusive_lq = new TProfile("fv2SPGap1B_inclusive_lq", "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
299 fOutput_lq->Add(fv2SPGap1B_inclusive_lq);
300
10a99a07 301 //small q resolution
302 fResSP_sq = new TProfile("fResSP_sq", "Resolution; centrality; Resolution", 9, -0.5, 8.5);
303 fOutput_sq->Add(fResSP_sq);
af960b8a 304
30e20fad 305 fResSP_vs_Cent_sq = new TProfile("fResSP_vs_Cent_sq", "Resolution; centrality; Resolution", 20., 0., 100.);
af960b8a 306 fOutput_sq->Add(fResSP_vs_Cent_sq);
8aa04ea0 307
30e20fad 308 f2partCumQA_vs_Cent_sq = new TProfile("f2partCumQA_vs_Cent_sq", "Resolution; centrality; Resolution", 100., 0., 100.);
af960b8a 309 fOutput_sq->Add(f2partCumQA_vs_Cent_sq);
8aa04ea0 310
30e20fad 311 f2partCumQB_vs_Cent_sq = new TProfile("f2partCumQB_vs_Cent_sq", "Resolution; centrality; Resolution", 100., 0., 100.);
af960b8a 312 fOutput_sq->Add(f2partCumQB_vs_Cent_sq);
8aa04ea0 313
314 // MC closure test
315 fv2SPGap1A_inclusive_sq = new TProfile("fv2SPGap1A_inclusive_sq", "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
316 fOutput_sq->Add(fv2SPGap1A_inclusive_sq);
317
318 fv2SPGap1B_inclusive_sq = new TProfile("fv2SPGap1B_inclusive_sq", "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
319 fOutput_sq->Add(fv2SPGap1B_inclusive_sq);
320
10a99a07 321 for (Int_t iC = 0; iC < 9; iC++){
322
323 fv2SPGap1A[iC] = new TProfile(Form("fv2SPGap1A_%d", iC), "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
324 fOutput->Add(fv2SPGap1A[iC]);
325
10a99a07 326 fv2SPGap1B[iC] = new TProfile(Form("fv2SPGap1B_%d", iC), "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
327 fOutput->Add(fv2SPGap1B[iC]);
328
bc76879c 329 fSinGap1Aq[iC] = new TProfile(Form("fSinGap1Aq_%d", iC), ";p_{T} (GeV/c);#LT sin(2*#phi) #GT", nptBins, ptBins);
330 fOutput->Add(fSinGap1Aq[iC]);
8aa04ea0 331
bc76879c 332 fCosGap1Aq[iC] = new TProfile(Form("fCosGap1Aq_%d", iC), ";p_{T} (GeV/c);#LT cos(2*#phi) #GT", nptBins, ptBins);
333 fOutput->Add(fCosGap1Aq[iC]);
8aa04ea0 334
bc76879c 335 fSinGap1Bq[iC] = new TProfile(Form("fSinGap1Bq_%d", iC), ";p_{T} (GeV/c);#LT sin(2*#phi) #GT", nptBins, ptBins);
336 fOutput->Add(fSinGap1Bq[iC]);
8aa04ea0 337
bc76879c 338 fCosGap1Bq[iC] = new TProfile(Form("fCosGap1Bq_%d", iC), ";p_{T} (GeV/c);#LT cos(2*#phi) #GT", nptBins, ptBins);
339 fOutput->Add(fCosGap1Bq[iC]);
340
10a99a07 341 fSinGap1A[iC] = new TProfile(Form("fSinGap1A_%d", iC), ";p_{T} (GeV/c);#LT sin(2*#phi) #GT", nptBins, ptBins);
342 fOutput->Add(fSinGap1A[iC]);
8aa04ea0 343
10a99a07 344 fCosGap1A[iC] = new TProfile(Form("fCosGap1A_%d", iC), ";p_{T} (GeV/c);#LT cos(2*#phi) #GT", nptBins, ptBins);
345 fOutput->Add(fCosGap1A[iC]);
8aa04ea0 346
10a99a07 347 fSinGap1B[iC] = new TProfile(Form("fSinGap1B_%d", iC), ";p_{T} (GeV/c);#LT sin(2*#phi) #GT", nptBins, ptBins);
348 fOutput->Add(fSinGap1B[iC]);
8aa04ea0 349
10a99a07 350 fCosGap1B[iC] = new TProfile(Form("fCosGap1B_%d", iC), ";p_{T} (GeV/c);#LT cos(2*#phi) #GT", nptBins, ptBins);
351 fOutput->Add(fCosGap1B[iC]);
8aa04ea0 352
10a99a07 353 //large q
354 fv2SPGap1A_lq[iC] = new TProfile(Form("fv2SPGap1A_lq_%d", iC), "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
355 fOutput_lq->Add(fv2SPGap1A_lq[iC]);
356
357 fv2SPGap1B_lq[iC] = new TProfile(Form("fv2SPGap1B_lq_%d", iC), "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
358 fOutput_lq->Add(fv2SPGap1B_lq[iC]);
359
bc76879c 360 fSinGap1Aq_lq[iC] = new TProfile(Form("fSinGap1Aq_lq_%d", iC), ";p_{T} (GeV/c);#LT sin(2*#phi) #GT", nptBins, ptBins);
361 fOutput_lq->Add(fSinGap1Aq_lq[iC]);
8aa04ea0 362
bc76879c 363 fCosGap1Aq_lq[iC] = new TProfile(Form("fCosGap1Aq_lq_%d", iC), ";p_{T} (GeV/c);#LT cos(2*#phi) #GT", nptBins, ptBins);
364 fOutput_lq->Add(fCosGap1Aq_lq[iC]);
8aa04ea0 365
bc76879c 366 fSinGap1Bq_lq[iC] = new TProfile(Form("fSinGap1Bq_lq_%d", iC), ";p_{T} (GeV/c);#LT sin(2*#phi) #GT", nptBins, ptBins);
367 fOutput_lq->Add(fSinGap1Bq_lq[iC]);
8aa04ea0 368
bc76879c 369 fCosGap1Bq_lq[iC] = new TProfile(Form("fCosGap1Bq_lq_%d", iC), ";p_{T} (GeV/c);#LT cos(2*#phi) #GT", nptBins, ptBins);
370 fOutput_lq->Add(fCosGap1Bq_lq[iC]);
371
10a99a07 372 fSinGap1A_lq[iC] = new TProfile(Form("fSinGap1A_lq_%d", iC), ";p_{T} (GeV/c);#LT sin(2*#phi) #GT", nptBins, ptBins);
373 fOutput_lq->Add(fSinGap1A_lq[iC]);
8aa04ea0 374
10a99a07 375 fCosGap1A_lq[iC] = new TProfile(Form("fCosGap1A_lq_%d", iC), ";p_{T} (GeV/c);#LT cos(2*#phi) #GT", nptBins, ptBins);
376 fOutput_lq->Add(fCosGap1A_lq[iC]);
8aa04ea0 377
10a99a07 378 fSinGap1B_lq[iC] = new TProfile(Form("fSinGap1B_lq_%d", iC), ";p_{T} (GeV/c);#LT sin(2*#phi) #GT", nptBins, ptBins);
379 fOutput_lq->Add(fSinGap1B_lq[iC]);
8aa04ea0 380
10a99a07 381 fCosGap1B_lq[iC] = new TProfile(Form("fCosGap1B_lq_%d", iC), ";p_{T} (GeV/c);#LT cos(2*#phi) #GT", nptBins, ptBins);
382 fOutput_lq->Add(fCosGap1B_lq[iC]);
8aa04ea0 383
10a99a07 384 //small q
385 fv2SPGap1A_sq[iC] = new TProfile(Form("fv2SPGap1A_sq_%d", iC), "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
386 fOutput_sq->Add(fv2SPGap1A_sq[iC]);
387
388 fv2SPGap1B_sq[iC] = new TProfile(Form("fv2SPGap1B_sq_%d", iC), "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
389 fOutput_sq->Add(fv2SPGap1B_sq[iC]);
390
bc76879c 391 fSinGap1Aq_sq[iC] = new TProfile(Form("fSinGap1Aq_sq_%d", iC), ";p_{T} (GeV/c);#LT sin(2*#phi) #GT", nptBins, ptBins);
392 fOutput_sq->Add(fSinGap1Aq_sq[iC]);
8aa04ea0 393
bc76879c 394 fCosGap1Aq_sq[iC] = new TProfile(Form("fCosGap1Aq_sq_%d", iC), ";p_{T} (GeV/c);#LT cos(2*#phi) #GT", nptBins, ptBins);
395 fOutput_sq->Add(fCosGap1Aq_sq[iC]);
8aa04ea0 396
bc76879c 397 fSinGap1Bq_sq[iC] = new TProfile(Form("fSinGap1Bq_sq_%d", iC), ";p_{T} (GeV/c);#LT sin(2*#phi) #GT", nptBins, ptBins);
398 fOutput_sq->Add(fSinGap1Bq_sq[iC]);
8aa04ea0 399
bc76879c 400 fCosGap1Bq_sq[iC] = new TProfile(Form("fCosGap1Bq_sq_%d", iC), "p_{T} (GeV/c);#LT cos(2*#phi) #GT", nptBins, ptBins);
401 fOutput_sq->Add(fCosGap1Bq_sq[iC]);
402
10a99a07 403 fSinGap1A_sq[iC] = new TProfile(Form("fSinGap1A_sq_%d", iC), ";p_{T} (GeV/c);#LT sin(2*#phi) #GT", nptBins, ptBins);
404 fOutput_sq->Add(fSinGap1A_sq[iC]);
8aa04ea0 405
10a99a07 406 fCosGap1A_sq[iC] = new TProfile(Form("fCosGap1A_sq_%d", iC), ";p_{T} (GeV/c);#LT cos(2*#phi) #GT", nptBins, ptBins);
407 fOutput_sq->Add(fCosGap1A_sq[iC]);
8aa04ea0 408
10a99a07 409 fSinGap1B_sq[iC] = new TProfile(Form("fSinGap1B_sq_%d", iC), ";p_{T} (GeV/c);#LT sin(2*#phi) #GT", nptBins, ptBins);
410 fOutput_sq->Add(fSinGap1B_sq[iC]);
8aa04ea0 411
10a99a07 412 fCosGap1B_sq[iC] = new TProfile(Form("fCosGap1B_sq_%d", iC), "p_{T} (GeV/c);#LT cos(2*#phi) #GT", nptBins, ptBins);
413 fOutput_sq->Add(fCosGap1B_sq[iC]);
8aa04ea0 414
52bfa7e6 415 //v2 vs qvec...
416 fResSP_vs_Qvec[iC] = new TProfile(Form("fResSP_vs_Qvec_%d", iC), "Resolution; Qvec (V0A); Resolution", 20., 0., 100.);
417 fResSP_vs_Qvec[iC]->Sumw2();
418 fOutput->Add(fResSP_vs_Qvec[iC]);
8aa04ea0 419
52bfa7e6 420 fV2IntGap1wq[iC] = new TProfile(Form("fV2IntGap1wq_%d", iC), "integrated v_{2} vs q-vector; Qvec (V0A); v_{2}", 20., 0., 100.);
421 fV2IntGap1wq[iC]->Sumw2();
422 fOutput->Add(fV2IntGap1wq[iC]);
10a99a07 423 };
8aa04ea0 424
52bfa7e6 425 fResGap1w = new TProfile("fResGap1w", "Resolution; centrality; Resolution", 9, -0.5, 8.5);
426 fResGap1w->Sumw2();
427 fOutput->Add(fResGap1w);
8aa04ea0 428
52bfa7e6 429 fV2IntGap1w = new TProfile("fV2IntGap1w", "; centrality; v_{2}", 9, -0.5, 8.5);
430 fV2IntGap1w->Sumw2();
431 fOutput->Add(fV2IntGap1w);
8aa04ea0 432
29fbaf8a 433 if(fIsMC){
434 fResSPmc_inclusive = new TProfile("fResSPmc_inclusive", "Resolution; ese; Resolution", 3, 0., 3.);
435 fOutput->Add(fResSPmc_inclusive);
436
437 fv2SPGap1Amc_inclusive_mb = new TProfile("fv2SPGap1Amc_inclusive_mb", "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
438 fOutput->Add(fv2SPGap1Amc_inclusive_mb);
d0761d58 439
29fbaf8a 440 fv2SPGap1Bmc_inclusive_mb = new TProfile("fv2SPGap1Bmc_inclusive_mb", "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
441 fOutput->Add(fv2SPGap1Bmc_inclusive_mb);
8aa04ea0 442
29fbaf8a 443 //large-q
444 fv2SPGap1Amc_inclusive_lq = new TProfile("fv2SPGap1Amc_inclusive_lq", "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
445 fOutput_lq->Add(fv2SPGap1Amc_inclusive_lq);
446
447 fv2SPGap1Bmc_inclusive_lq = new TProfile("fv2SPGap1Bmc_inclusive_lq", "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
448 fOutput_lq->Add(fv2SPGap1Bmc_inclusive_lq);
8aa04ea0 449
29fbaf8a 450 //small-q
451 fv2SPGap1Amc_inclusive_sq = new TProfile("fv2SPGap1Amc_inclusive_sq", "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
452 fOutput_sq->Add(fv2SPGap1Amc_inclusive_sq);
d0761d58 453
29fbaf8a 454 fv2SPGap1Bmc_inclusive_sq = new TProfile("fv2SPGap1Bmc_inclusive_sq", "v_{2}{2} vs p_{T}; p_{T} (GeV/c); v_{2}{2}", nptBins, ptBins);
455 fOutput_sq->Add(fv2SPGap1Bmc_inclusive_sq);
d0761d58 456 }
8aa04ea0 457
10a99a07 458 PostData(1, fOutput );
459 PostData(2, fEventCuts);
460 PostData(3, fTrackCuts);
461 PostData(4, fOutput_lq );
462 PostData(5, fOutput_sq );
463}
464
465//________________________________________________________________________
466
467void AliAnalysisTaskV2AllChAOD::UserExec(Option_t *)
468{
469 //Printf("An event");
470 // main event loop
471 fAOD = dynamic_cast<AliAODEvent*>(fInputEvent);
472 if (!fAOD) {
473 AliWarning("ERROR: AliAODEvent not available \n");
474 return;
475 }
8aa04ea0 476
10a99a07 477 if (strcmp(fAOD->ClassName(), "AliAODEvent"))
8aa04ea0 478 {
479 AliFatal("Not processing AODs");
480 }
481
10a99a07 482 if(!fEventCuts->IsSelected(fAOD,fTrackCuts))return;//event selection
8aa04ea0 483
29fbaf8a 484 //Get q-vector percentile.
29fbaf8a 485 Double_t Qvec=0.;
34f4a5fc 486 if(fIsMC && fQvecGen) Qvec = fEventCuts->GetQvecPercentileMC(fVZEROside, fQgenType);
487 else Qvec = fEventCuts->GetQvecPercentile(fVZEROside);
29fbaf8a 488
86f0713c 489 fQvector->Fill(Qvec);
490 if (Qvec > fCutLargeQperc && Qvec < 100.) fQvector_lq->Fill(Qvec);
491 if (Qvec > 0. && Qvec < fCutSmallQperc) fQvector_sq->Fill(Qvec);
492
509a25f9 493 Double_t Cent=(fDoCentrSystCentrality)?1.01*fEventCuts->GetCent():fEventCuts->GetCent();
29fbaf8a 494 fCentrality->Fill(Cent);
8aa04ea0 495
30e20fad 496 Int_t centV0 = -1;
10a99a07 497 if ((Cent > 0) && (Cent <= 5.0))
8aa04ea0 498 centV0 = 0;
499 else if ((Cent > 5.0) && (Cent <= 10.0))
500 centV0 = 1;
501 else if ((Cent > 10.0) && (Cent <= 20.0))
502 centV0 = 2;
503 else if ((Cent > 20.0) && (Cent <= 30.0))
504 centV0 = 3;
505 else if ((Cent > 30.0) && (Cent <= 40.0))
506 centV0 = 4;
507 else if ((Cent > 40.0) && (Cent <= 50.0))
508 centV0 = 5;
509 else if ((Cent > 50.0) && (Cent <= 60.0))
510 centV0 = 6;
511 else if ((Cent > 60.0) && (Cent <= 70.0))
512 centV0 = 7;
513 else if ((Cent > 70.0) && (Cent <= 80.0))
514 centV0 = 8;
515
52bfa7e6 516 if(centV0==-1)return; // FIXME if the centrality is not defined or >80%... return!!!
8aa04ea0 517
29fbaf8a 518 if(fIsMC) MCclosure(Qvec); // fill mc histograms for montecarlo closure
10a99a07 519
520 Double_t QxGap1A = 0., QyGap1A = 0.;
521 Double_t QxGap1B = 0., QyGap1B = 0.;
e4d189de 522 Double_t multGap1A = 0, multGap1B = 0;
8aa04ea0 523
d57dd83f 524 for (Int_t loop = 0; loop < 2; loop++){
10a99a07 525
d57dd83f 526 //main loop on tracks
527 for (Int_t iTracks = 0; iTracks < fAOD->GetNumberOfTracks(); iTracks++) {
f15c1f69 528 AliAODTrack* track = dynamic_cast<AliAODTrack*>(fAOD->GetTrack(iTracks));
529 if(!track) AliFatal("Not a standard AOD");
d57dd83f 530 if(fCharge != 0 && track->Charge() != fCharge) continue;//if fCharge != 0 only select fCharge
531 if (!fTrackCuts->IsSelected(track,kTRUE)) continue; //track selection (rapidity selection NOT in the standard cuts)
8aa04ea0 532
bc76879c 533 fEta_vs_Phi_bef->Fill( track->Eta(), track->Phi() );
8aa04ea0 534
535 // if (fIsRecoEff){
536 //
537 // // 2) reject randomly tracks at high pT until the reconstruction efficiency becomes flat (add the following before the loop == 0 part): (mail by Alexandru)
538 //
539 // Double_t recoEff = GetRecoEff(track->Pt(), centV0);
540 // if (recoEff < 0){
541 // cout<<"No reconstruction efficiency!"<<endl;
542 // continue;
543 // }
544 //
545 // Double_t rndPt = gRandom->Rndm();
546 // // cout<<"rndPt: "<<rndPt<<endl;
547
548 // Double_t minRecPt = GetRecoEff(0.200001, centV0);
549 // // cout<<"minRecPt: "<<minRecPt<<endl;
550
551 // if (rndPt > minRecPt/recoEff){
552 // // cout<<"Track rejected: "<<iTracks<<" from "<<fAOD->GetNumberOfTracks()<<endl;
553 // continue;
554 // }
555 //
556 // } // end fIsRecoEff
557
d57dd83f 558 if (loop == 0) {
8aa04ea0 559
560 if (track->Eta() > fEtaGapMax){
561
562 if(fIsRecoEff){
563 Double_t receff = GetRecoEff(track->Pt(), centV0);
564 QxGap1A += (TMath::Cos(2.*track->Phi()))/receff;
e4d189de 565 QyGap1A += (TMath::Sin(2.*track->Phi()))/receff;
8aa04ea0 566 multGap1A+=1./receff;
567 } else {
568 QxGap1A += TMath::Cos(2.*track->Phi());
e4d189de 569 QyGap1A += TMath::Sin(2.*track->Phi());
8aa04ea0 570 multGap1A++;
571 }
10a99a07 572
bc76879c 573 fSinGap1Aq[centV0]->Fill(track->Pt(), TMath::Sin(2.*track->Phi()));
574 fCosGap1Aq[centV0]->Fill(track->Pt(), TMath::Cos(2.*track->Phi()));
8aa04ea0 575
576 fEta_vs_PhiA->Fill( track->Eta(), track->Phi() );
577
962f3d11 578 if (Qvec > fCutLargeQperc && Qvec < 100.){
8aa04ea0 579 fSinGap1Aq_lq[centV0]->Fill(track->Pt(), TMath::Sin(2.*track->Phi()));
580 fCosGap1Aq_lq[centV0]->Fill(track->Pt(), TMath::Cos(2.*track->Phi()));
581 }
582
962f3d11 583 if (Qvec > 0. && Qvec < fCutSmallQperc){
8aa04ea0 584 fSinGap1Aq_sq[centV0]->Fill(track->Pt(), TMath::Sin(2.*track->Phi()));
585 fCosGap1Aq_sq[centV0]->Fill(track->Pt(), TMath::Cos(2.*track->Phi()));
586 }
587
d57dd83f 588 }
8aa04ea0 589
590 if (track->Eta() < fEtaGapMin){
591
592 if(fIsRecoEff){
593 Double_t receff = GetRecoEff(track->Pt(), centV0);
594 QxGap1B += (TMath::Cos(2.*track->Phi()))/receff;
595 QyGap1B += (TMath::Sin(2.*track->Phi()))/receff;
596 multGap1B+=1./receff;
597 } else {
598 QxGap1B += TMath::Cos(2.*track->Phi());
599 QyGap1B += TMath::Sin(2.*track->Phi());
600 multGap1B++;
601 }
602
603 fCosGap1Bq[centV0]->Fill(track->Pt(), TMath::Cos(2.*track->Phi()));
604 fSinGap1Bq[centV0]->Fill(track->Pt(), TMath::Sin(2.*track->Phi()));
605
606 fEta_vs_PhiB->Fill( track->Eta(), track->Phi() );
607
608 if (Qvec > fCutLargeQperc && Qvec < 100.){
609 fSinGap1Bq_lq[centV0]->Fill(track->Pt(), TMath::Sin(2.*track->Phi()));
610 fCosGap1Bq_lq[centV0]->Fill(track->Pt(), TMath::Cos(2.*track->Phi()));
611 }
612
613 if (Qvec > 0. && Qvec < fCutSmallQperc){
614 fSinGap1Bq_sq[centV0]->Fill(track->Pt(), TMath::Sin(2.*track->Phi()));
615 fCosGap1Bq_sq[centV0]->Fill(track->Pt(), TMath::Cos(2.*track->Phi()));
616 }
617
d57dd83f 618 }
8aa04ea0 619
620 } else {
621
d57dd83f 622 //eval v2 scalar product
623 if (track->Eta() < fEtaGapMin && multGap1A > 0){
8aa04ea0 624
625 Double_t v2SPGap1A = 0.;
626 if(fIsRecoEff){
627 Double_t receff = GetRecoEff(track->Pt(), centV0);
628 Double_t uxGap1A = (TMath::Cos(2.*track->Phi()))/receff;
629 Double_t uyGap1A = (TMath::Sin(2.*track->Phi()))/receff;
630 // multGap1A = multGap1A/receff;
e4d189de 631 v2SPGap1A = (uxGap1A*QxGap1A + uyGap1A*QyGap1A)/(Double_t)multGap1A;
8aa04ea0 632 } else v2SPGap1A = (TMath::Cos(2.*track->Phi())*QxGap1A + TMath::Sin(2.*track->Phi())*QyGap1A)/(Double_t)multGap1A;
633
634 // Double_t v2SPGap1A = (TMath::Cos(2.*track->Phi())*QxGap1A + TMath::Sin(2.*track->Phi())*QyGap1A)/(Double_t)multGap1A;
635
636 fv2SPGap1A[centV0]->Fill(track->Pt(), v2SPGap1A);
637
638 fv2SPGap1A_inclusive_mb->Fill(track->Pt(), v2SPGap1A); //mb v2 for mc closure
639
640 fSinGap1A[centV0]->Fill(track->Pt(), TMath::Sin(2.*track->Phi()));
bc76879c 641 fCosGap1A[centV0]->Fill(track->Pt(), TMath::Cos(2.*track->Phi()));
8aa04ea0 642
962f3d11 643 if (Qvec > fCutLargeQperc && Qvec < 100.){
8aa04ea0 644 fv2SPGap1A_lq[centV0]->Fill(track->Pt(), v2SPGap1A);
645 fSinGap1A_lq[centV0]->Fill(track->Pt(), TMath::Sin(2.*track->Phi()));
646 fCosGap1A_lq[centV0]->Fill(track->Pt(), TMath::Cos(2.*track->Phi()));
647
648 fv2SPGap1A_inclusive_lq->Fill(track->Pt(), v2SPGap1A); //lq v2 for mc closure
649 }
650
962f3d11 651 if (Qvec > 0. && Qvec < fCutSmallQperc){
8aa04ea0 652 fv2SPGap1A_sq[centV0]->Fill(track->Pt(), v2SPGap1A);
653 fSinGap1A_sq[centV0]->Fill(track->Pt(), TMath::Sin(2.*track->Phi()));
654 fCosGap1A_sq[centV0]->Fill(track->Pt(), TMath::Cos(2.*track->Phi()));
655
656 fv2SPGap1A_inclusive_sq->Fill(track->Pt(), v2SPGap1A); //sq v2 for mc closure
657 }
bc9bd5b3 658
d57dd83f 659 }
8aa04ea0 660
d57dd83f 661 if (track->Eta() > fEtaGapMax && multGap1B > 0){
8aa04ea0 662
663 Double_t v2SPGap1B = 0;
664 if(fIsRecoEff){
665 Double_t receff = GetRecoEff(track->Pt(), centV0);
666 Double_t uxGap1B = (TMath::Cos(2.*track->Phi()))/receff;
667 Double_t uyGap1B = (TMath::Sin(2.*track->Phi()))/receff;
668// multGap1B = multGap1B/receff;
669 v2SPGap1B = (uxGap1B*QxGap1B + uyGap1B*QyGap1B)/(Double_t)multGap1B;
670 }
671 else v2SPGap1B = (TMath::Cos(2.*track->Phi())*QxGap1B + TMath::Sin(2.*track->Phi())*QyGap1B)/(Double_t)multGap1B;
672
673 // Double_t v2SPGap1B = (TMath::Cos(2.*track->Phi())*QxGap1B + TMath::Sin(2.*track->Phi())*QyGap1B)/(Double_t)multGap1B;
674
675 fv2SPGap1B[centV0]->Fill(track->Pt(), v2SPGap1B);
676
677 fv2SPGap1B_inclusive_mb->Fill(track->Pt(), v2SPGap1B); //mb v2 for mc closure
678
bc76879c 679 fCosGap1B[centV0]->Fill(track->Pt(), TMath::Cos(2.*track->Phi()));
680 fSinGap1B[centV0]->Fill(track->Pt(), TMath::Sin(2.*track->Phi()));
8aa04ea0 681
962f3d11 682 if (Qvec > fCutLargeQperc && Qvec < 100.){
8aa04ea0 683 fv2SPGap1B_lq[centV0]->Fill(track->Pt(), v2SPGap1B);
684 fSinGap1B_lq[centV0]->Fill(track->Pt(), TMath::Sin(2.*track->Phi()));
685 fCosGap1B_lq[centV0]->Fill(track->Pt(), TMath::Cos(2.*track->Phi()));
686
687 fv2SPGap1B_inclusive_lq->Fill(track->Pt(), v2SPGap1B); //lq v2 for mc closure
688 }
689
962f3d11 690 if (Qvec > 0. && Qvec < fCutSmallQperc){
8aa04ea0 691 fv2SPGap1B_sq[centV0]->Fill(track->Pt(), v2SPGap1B);
692 fSinGap1B_sq[centV0]->Fill(track->Pt(), TMath::Sin(2.*track->Phi()));
693 fCosGap1B_sq[centV0]->Fill(track->Pt(), TMath::Cos(2.*track->Phi()));
694
695 fv2SPGap1B_inclusive_sq->Fill(track->Pt(), v2SPGap1B); //sq v2 for mc closure
696 }
697
d57dd83f 698 }
699 }// end else
700 } // end loop on tracks
701 } // end loop
8aa04ea0 702
703
10a99a07 704 if (multGap1A > 0 && multGap1B > 0){
705 Double_t res = (QxGap1A*QxGap1B + QyGap1A*QyGap1B)/(Double_t)multGap1A/(Double_t)multGap1B;
706 fResSP->Fill((Double_t)centV0, res);
8aa04ea0 707
29fbaf8a 708 fResSP_inclusive->Fill(0., res); //mb v2 for mc closure
30e20fad 709 fResSP_vs_Cent->Fill(Cent, res);
8aa04ea0 710
af960b8a 711 Double_t f2partCumQA = -999.;
712 if(multGap1A>1)
713 f2partCumQA = ( ( (QxGap1A*QxGap1A + QyGap1A*QyGap1A) - (Double_t)multGap1A ) / ((Double_t)multGap1A*((Double_t)multGap1A-1)) );
30e20fad 714 if(f2partCumQA>0)f2partCumQA_vs_Cent->Fill((Double_t)Cent,f2partCumQA);
8aa04ea0 715
af960b8a 716 Double_t f2partCumQB = -999.;
717 if(multGap1B>1)
718 f2partCumQB = ( ( (QxGap1B*QxGap1B + QyGap1B*QyGap1B) - (Double_t)multGap1B ) / ((Double_t)multGap1B*((Double_t)multGap1B-1)) );
30e20fad 719 if(f2partCumQB>0)f2partCumQB_vs_Cent->Fill((Double_t)Cent,f2partCumQB);
8aa04ea0 720
962f3d11 721 if (Qvec > fCutLargeQperc && Qvec < 100.){
10a99a07 722 fResSP_lq->Fill((Double_t)centV0, res);
30e20fad 723 fResSP_vs_Cent_lq->Fill(Cent, res);
724 if(f2partCumQA>0)f2partCumQA_vs_Cent_lq->Fill((Double_t)Cent,f2partCumQA);
725 if(f2partCumQB>0)f2partCumQB_vs_Cent_lq->Fill((Double_t)Cent,f2partCumQB);
8aa04ea0 726
29fbaf8a 727 fResSP_inclusive->Fill(1., res); //lq v2 for mc closure
962f3d11 728 }
8aa04ea0 729
962f3d11 730 if (Qvec > 0. && Qvec < fCutSmallQperc){
10a99a07 731 fResSP_sq->Fill((Double_t)centV0, res);
30e20fad 732 fResSP_vs_Cent_sq->Fill(Cent, res);
733 if(f2partCumQA>0)f2partCumQA_vs_Cent_sq->Fill((Double_t)Cent,f2partCumQA);
734 if(f2partCumQB>0)f2partCumQB_vs_Cent_sq->Fill((Double_t)Cent,f2partCumQB);
8aa04ea0 735
29fbaf8a 736 fResSP_inclusive->Fill(2., res); //sq v2 for mc closure
962f3d11 737 }
e3bcd147 738 }// end multiplicity if
bc9bd5b3 739
52bfa7e6 740 //v2 vs qvec
741 if ((multGap1A > 0) && (multGap1B > 0)){
8aa04ea0 742
52bfa7e6 743 fResGap1w->Fill(Double_t(centV0), (QxGap1A*QxGap1B + QyGap1A*QyGap1B)/multGap1A/multGap1B, (Double_t)(multGap1A*multGap1B));
8aa04ea0 744
52bfa7e6 745 Double_t nGap1 = multGap1A*multGap1B;
746 Double_t qqGap1 = (QxGap1A*QxGap1B + QyGap1A*QyGap1B)/nGap1;
8aa04ea0 747
52bfa7e6 748 fV2IntGap1w->Fill(Double_t(centV0), qqGap1, nGap1);
8aa04ea0 749
52bfa7e6 750 fResSP_vs_Qvec[centV0]->Fill(Double_t(Qvec), (QxGap1A*QxGap1B + QyGap1A*QyGap1B)/multGap1A/multGap1B, (Double_t)(multGap1A*multGap1B));
751 fV2IntGap1wq[centV0]->Fill(Double_t(Qvec), qqGap1, nGap1);
8aa04ea0 752
10a99a07 753 }
8aa04ea0 754
e4d189de 755 if( fFillTHn ){
756
757 Double_t varEv[6];
758 varEv[0]=Cent;
8aa04ea0 759
e4d189de 760 varEv[1]=(Double_t)Qvec; // qvec_rec_perc
8aa04ea0 761
e4d189de 762 varEv[2]=(Double_t)fEventCuts->GetqV0A();
763
764 varEv[3]=(Double_t)fEventCuts->GetqV0C();
765
766 varEv[4]=(Double_t)fEventCuts->GetqTPC();
8aa04ea0 767
e4d189de 768 varEv[5]=(Double_t)fEventCuts->GetNch(); // Nch
8aa04ea0 769
e4d189de 770 ((THnSparseF*)fOutput->FindObject("NSparseHistEv"))->Fill(varEv);//event loop
771
772 }
e3bcd147 773
8aa04ea0 774
10a99a07 775 PostData(1, fOutput );
776 PostData(2, fEventCuts);
777 PostData(3, fTrackCuts);
778 PostData(4, fOutput_lq );
779 PostData(5, fOutput_sq );
780}
781
782//_________________________________________________________________
783Bool_t AliAnalysisTaskV2AllChAOD::GetDCA(const AliAODTrack* trk, Double_t * p){
8aa04ea0 784
10a99a07 785 //AliAODTrack::DCA(): for newest AOD fTrack->DCA() always gives -999. This should fix.
786 //FIXME should update EventCuts?
787 //FIXME add track->GetXYZ(p) method
8aa04ea0 788
10a99a07 789 double xyz[3],cov[3];
8aa04ea0 790
10a99a07 791 if (!trk->GetXYZ(xyz)) { // dca is not stored
792 AliExternalTrackParam etp;
793 etp.CopyFromVTrack(trk);
794 AliVEvent* ev = (AliVEvent*)trk->GetEvent();
795 if (!ev) {/*printf("Event is not connected to the track\n");*/ return kFALSE;}
796 if (!etp.PropagateToDCA(ev->GetPrimaryVertex(), ev->GetMagneticField(),999,xyz,cov)) return kFALSE; // failed, track is too far from vertex
797 }
798 p[0] = xyz[0];
799 p[1] = xyz[1];
800 return kTRUE;
801
802}
803
a7abb826 804//_________________________________________________________________
29fbaf8a 805void AliAnalysisTaskV2AllChAOD::MCclosure(Double_t qvec){
a7abb826 806 // First do MC to fill up the MC particle array
8aa04ea0 807
a7abb826 808 TClonesArray *arrayMC = 0;
809 if (fIsMC)
8aa04ea0 810 {
811 arrayMC = (TClonesArray*) fAOD->GetList()->FindObject(AliAODMCParticle::StdBranchName());
812 if (!arrayMC) {
813 AliFatal("Error: MC particles branch not found!\n");
814 }
815
816 Double_t QxGap1Amc = 0., QyGap1Amc = 0.;
817 Double_t QxGap1Bmc = 0., QyGap1Bmc = 0.;
818 Int_t multGap1Amc = 0, multGap1Bmc = 0;
819
820 for (Int_t loop = 0; loop < 2; loop++){
821
822 Int_t nMC = arrayMC->GetEntries();
823
824 for (Int_t iMC = 0; iMC < nMC; iMC++)
825 {
826 AliAODMCParticle *partMC = (AliAODMCParticle*) arrayMC->At(iMC);
827 if(!partMC->Charge()) continue;//Skip neutrals
828 if(fCharge != 0 && partMC->Charge()*fCharge < 0.) continue;//if fCharge != 0 only select fCharge
829
830 if (!(partMC->IsPhysicalPrimary()))
831 continue;
832
833 if(partMC->Eta()<fTrackCuts->GetEtaMin() || partMC->Eta()>fTrackCuts->GetEtaMax()) continue;
834
835 //Printf("a particle");
836
837
838 if (loop == 0) {
839
840 if (partMC->Eta() > fEtaGapMax){
841 QxGap1Amc += TMath::Cos(2.*partMC->Phi());
842 QyGap1Amc += TMath::Sin(2.*partMC->Phi());
843 multGap1Amc++;
844 }
845
846 if (partMC->Eta() < fEtaGapMin){
847 QxGap1Bmc += TMath::Cos(2.*partMC->Phi());
848 QyGap1Bmc += TMath::Sin(2.*partMC->Phi());
849 multGap1Bmc++;
850 }
851
852 } else {
853
854 //eval v2 scalar product
855 if (partMC->Eta() < fEtaGapMin && multGap1Amc > 0){
856 Double_t v2SPGap1Amc = (TMath::Cos(2.*partMC->Phi())*QxGap1Amc + TMath::Sin(2.*partMC->Phi())*QyGap1Amc)/(Double_t)multGap1Amc;
857 fv2SPGap1Amc_inclusive_mb->Fill(partMC->Pt(), v2SPGap1Amc);
858
859 if (qvec > fCutLargeQperc && qvec < 100.){
860 fv2SPGap1Amc_inclusive_lq->Fill(partMC->Pt(), v2SPGap1Amc);
861 }
862
863 if (qvec > 0. && qvec < fCutSmallQperc){
864 fv2SPGap1Amc_inclusive_sq->Fill(partMC->Pt(), v2SPGap1Amc);
865 }
866
867 }
868
869 if (partMC->Eta() > fEtaGapMax && multGap1Bmc > 0){
870 Double_t v2SPGap1Bmc = (TMath::Cos(2.*partMC->Phi())*QxGap1Bmc + TMath::Sin(2.*partMC->Phi())*QyGap1Bmc)/(Double_t)multGap1Bmc;
871 fv2SPGap1Bmc_inclusive_mb->Fill(partMC->Pt(), v2SPGap1Bmc);
872
873 if (qvec > fCutLargeQperc && qvec < 100.){
874 fv2SPGap1Bmc_inclusive_lq->Fill(partMC->Pt(), v2SPGap1Bmc);
875 }
876
877 if (qvec > 0. && qvec < fCutSmallQperc){
878 fv2SPGap1Bmc_inclusive_sq->Fill(partMC->Pt(), v2SPGap1Bmc);
879 }
880
881 }
882
883
884 }// end else
885 } // end loop on partMCs
886 } // end loop
887
888 if (multGap1Amc > 0 && multGap1Bmc > 0){
889 Double_t resmc = (QxGap1Amc*QxGap1Bmc + QyGap1Amc*QyGap1Bmc)/(Double_t)multGap1Amc/(Double_t)multGap1Bmc;
890 fResSPmc_inclusive->Fill(0.,resmc);
891
892 if (qvec > fCutLargeQperc && qvec < 100.){
893 fResSPmc_inclusive->Fill(1.,resmc);
a7abb826 894 }
8aa04ea0 895
896 if (qvec > 0. && qvec < fCutSmallQperc){
897 fResSPmc_inclusive->Fill(2.,resmc);
a7abb826 898 }
8aa04ea0 899
900 }
901
902 }// end if MC
a7abb826 903}
8aa04ea0 904
d0761d58 905//_________________________________________________________________
906Double_t AliAnalysisTaskV2AllChAOD::GetRecoEff(Double_t pt, Int_t iC){
907
908 if(iC>8) return 1.;
909
910 if(pt<0.2 || pt>100.) return 1.;
911
8aa04ea0 912 // // // //spectra ese binning
913 // // // //const Double_t ptBins[] = {0.20,0.30,0.4,0.5,0.6,0.7,0.8,0.9,1.0,1.2,1.4,1.6,1.8,2.0,2.4,2.8,3.2,3.6,4.0,5.0,6.0,7.0,8.0,9.0,10.,12.,15.,20.,25.,30.,35.,40.,50.,75.,100.};
914 // // // //const Int_t nptBins=34;
915
d0761d58 916 const Double_t fEpsilon=0.000001;
8aa04ea0 917 // cout<<"list: "<<endl;
918 // fRecoEffList->ls();
e4d189de 919 TH1D *h = (TH1D*)fRecoEffList->At(iC);
8aa04ea0 920
d0761d58 921 Int_t bin = h->FindBin(pt);
8aa04ea0 922
d0761d58 923 Double_t lowlim = h->GetBinLowEdge(bin);
924 Double_t uplim = h->GetBinLowEdge(bin) + h->GetBinWidth(bin);
8aa04ea0 925
e4d189de 926 Double_t eff = 1.;
8aa04ea0 927
e4d189de 928 if( pt>lowlim && pt<uplim ) eff = h->GetBinContent(bin);
929 if( pt == lowlim ) eff = h->GetBinContent( h->FindBin(pt+fEpsilon) );
930 if( pt == uplim ) eff = h->GetBinContent( h->FindBin(pt-fEpsilon) );
8aa04ea0 931
e4d189de 932 return eff;
8aa04ea0 933
d0761d58 934}
10a99a07 935//_________________________________________________________________
936void AliAnalysisTaskV2AllChAOD::Terminate(Option_t *)
937{
938 // Terminate
939}