]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGLF/SPECTRA/PiKaPr/TestAOD/AliSpectraAODEventCuts.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
cxx file
[u/mrichter/AliRoot.git] / PWGLF / SPECTRA / PiKaPr / TestAOD / AliSpectraAODEventCuts.cxx
CommitLineData
c88234ad 1
2/**************************************************************************
3 * Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
4 * *
5 * Author: The ALICE Off-line Project. *
6 * Contributors are mentioned in the code where appropriate. *
7 * *
8 * Permission to use, copy, modify and distribute this software and its *
9 * documentation strictly for non-commercial purposes is hereby granted *
10 * without fee, provided that the above copyright notice appears in all *
11 * copies and that both the copyright notice and this permission notice *
12 * appear in the supporting documentation. The authors make no claims *
13 * about the suitability of this software for any purpose. It is *
14 * provided "as is" without express or implied warranty. *
15 **************************************************************************/
16
17//-----------------------------------------------------------------
18// AliSpectraAODEventCuts class
19//-----------------------------------------------------------------
20
21#include "TChain.h"
22#include "TTree.h"
23#include "TLegend.h"
24#include "TH1F.h"
829b5a81 25#include "TH1I.h"
c88234ad 26#include "TH2F.h"
27#include "TCanvas.h"
93db93de 28#include "TProfile.h"
7b364a00 29#include "TSpline.h"
c88234ad 30#include "AliAnalysisTask.h"
31#include "AliAnalysisManager.h"
32#include "AliAODTrack.h"
33#include "AliAODMCParticle.h"
34#include "AliAODEvent.h"
35#include "AliAODInputHandler.h"
36#include "AliAnalysisTaskESDfilter.h"
37#include "AliAnalysisDataContainer.h"
93db93de 38#include "AliESDVZERO.h"
39#include "AliAODVZERO.h"
c88234ad 40#include "AliSpectraAODEventCuts.h"
f6a38178 41#include "AliSpectraAODTrackCuts.h"
c88234ad 42#include <iostream>
43
44using namespace std;
45
46ClassImp(AliSpectraAODEventCuts)
47
93db93de 48AliSpectraAODEventCuts::AliSpectraAODEventCuts(const char *name) :
49 TNamed(name, "AOD Event Cuts"),
50 fAOD(0),
51 fSelectBit(AliVEvent::kMB),
52 fCentralityMethod("V0M"),
53 fTrackBits(1),
54 fIsMC(0),
48c4d28e 55 fIsLHC10h(1),
93db93de 56 fTrackCuts(0),
57 fIsSelected(0),
58 fCentralityCutMin(0.),
59 fCentralityCutMax(999),
60 fQVectorCutMin(-999.),
61 fQVectorCutMax(999.),
a162b3e2 62 fVertexCutMin(-10.),
63 fVertexCutMax(10.),
93db93de 64 fMultiplicityCutMin(-999.),
65 fMultiplicityCutMax(99999.),
48c4d28e 66 fqV0C(-999.),
67 fqV0A(-999.),
4fd2f0ad
LM		 68 fqV0Cx(-999.),
69 fqV0Ax(-999.),
70 fqV0Cy(-999.),
71 fqV0Ay(-999.),
fdbff5a1 72 fPsiV0C(-999.),
73 fPsiV0A(-999.),
48c4d28e 74 fCent(-999.),
93db93de 75 fOutput(0),
76 fCalib(0),
77 fRun(-1),
78 fMultV0(0),
79 fV0Cpol1(-1),
80 fV0Cpol2(-1),
81 fV0Cpol3(-1),
82 fV0Cpol4(-1),
83 fV0Apol1(-1),
84 fV0Apol2(-1),
85 fV0Apol3(-1),
7b364a00 86 fV0Apol4(-1),
87 fQvecIntList(0),
88 fQvecIntegral(0),
89 fSplineArrayV0A(0),
90 fSplineArrayV0C(0)
c88234ad 91{
f6a38178 92 // Constructor
93db93de 93 fOutput=new TList();
94 fOutput->SetOwner();
95 fOutput->SetName("fOutput");
96
97 fCalib=new TList();
98 fCalib->SetOwner();
99 fCalib->SetName("fCalib");
100
7b364a00 101 fQvecIntList=new TList();
102 fQvecIntList->SetOwner();
103 fQvecIntList->SetName("fQvecIntList");
104
93db93de 105 TH1I *fHistoCuts = new TH1I("fHistoCuts", "Event Cuts", kNVtxCuts, -0.5, kNVtxCuts - 0.5);
23f51e78 106 TH1F *fHistoVtxBefSel = new TH1F("fHistoVtxBefSel", "Vtx distr before event selection;z (cm)",500,-15,15);
107 TH1F *fHistoVtxAftSel = new TH1F("fHistoVtxAftSel", "Vtx distr after event selection;z (cm)",500,-15,15);
108 TH1F *fHistoEtaBefSel = new TH1F("fHistoEtaBefSel", "Eta distr before event selection;eta",500,-2,2);
109 TH1F *fHistoEtaAftSel = new TH1F("fHistoEtaAftSel", "Eta distr after event selection;eta",500,-2,2);
110 TH1F *fHistoNChAftSel = new TH1F("fHistoNChAftSel", "NCh distr after event selection;Nch",2000,-0.5,1999.5);
111 TH2F *fHistoQVector = new TH2F("fHistoQVector", "QVector with VZERO distribution;centrality;Q vector from EP task",20,0,100,100,0,10);
112 TH2F *fHistoEP = new TH2F("fHistoEP", "EP with VZERO distribution;centrality;Psi_{EP} from EP task",20,0,100,100,-2,2);
113 TH2F *fPsiACor = new TH2F("fPsiACor", "EP with VZERO A distribution;centrality;Psi_{EP} VZERO-A",20,0,100,100,-2,2);
114 TH2F *fPsiCCor = new TH2F("fPsiCCor", "EP with VZERO C distribution;centrality;Psi_{EP} VZERO-C",20,0,100,100,-2,2);
115 TH2F *fQVecACor = new TH2F("fQVecACor", "QVec VZERO A;centrality;Qvector VZERO-A",20,0,100,100,0,10);
116 TH2F *fQVecCCor = new TH2F("fQVecCCor", "QVec VZERO C;centrality;Qvector VZERO-C",20,0,100,100,0,10);
117 TH2F *fV0M = new TH2F("fV0M", "V0 Multiplicity, before correction;V0 sector",64,-.5,63.5,500,0,1000);
118 TH2F *fV0MCor = new TH2F("fV0MCor", "V0 Multiplicity, after correction;V0 sector",64,-.5,63.5,500,0,1000);
93db93de 119
7b364a00 120 fSplineArrayV0A = new TObjArray();
121 fSplineArrayV0A->SetOwner();
122 fSplineArrayV0C = new TObjArray();
123 fSplineArrayV0C->SetOwner();
124
93db93de 125 fOutput->Add(fHistoCuts);
126 fOutput->Add(fHistoVtxBefSel);
127 fOutput->Add(fHistoVtxAftSel);
128 fOutput->Add(fHistoEtaBefSel);
129 fOutput->Add(fHistoEtaAftSel);
130 fOutput->Add(fHistoNChAftSel);
131 fOutput->Add(fHistoQVector);
132 fOutput->Add(fHistoEP);
133 fOutput->Add(fPsiACor);
134 fOutput->Add(fPsiCCor);
135 fOutput->Add(fQVecACor);
136 fOutput->Add(fQVecCCor);
23f51e78 137 fOutput->Add(fV0M);
138 fOutput->Add(fV0MCor);
93db93de 139
140 for (Int_t i = 0; i<10; i++){
141 fMeanQxa2[i] = -1;
142 fMeanQya2[i] = -1;
143 fMeanQxc2[i] = -1;
144 fMeanQyc2[i] = -1;
145 }
c88234ad 146}
147
148//______________________________________________________
f6a38178 149Bool_t AliSpectraAODEventCuts::IsSelected(AliAODEvent * aod,AliSpectraAODTrackCuts *trackcuts)
c88234ad 150{
f6a38178 151 // Returns true if Event Cuts are selected and applied
152 fAOD = aod;
31579941 153 fTrackCuts = trackcuts; // FIXME: if track cuts is 0, do not set (use the pre-set member). Do we need to pass this here at all??
154 // FIXME: all those references by name are slow.
93db93de 155 ((TH1I*)fOutput->FindObject("fHistoCuts"))->Fill(kProcessedEvents);
156 Bool_t IsPhysSel = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & fSelectBit);
decf69d9 157 if(!IsPhysSel)return IsPhysSel;
93db93de 158 ((TH1I*)fOutput->FindObject("fHistoCuts"))->Fill(kPhysSelEvents);
f6a38178 159 //loop on tracks, before event selection, filling QA histos
160 AliAODVertex * vertex = fAOD->GetPrimaryVertex();//FIXME vertex is recreated
93db93de 161 if(vertex)((TH1F*)fOutput->FindObject("fHistoVtxBefSel"))->Fill(vertex->GetZ());
decf69d9 162 fIsSelected =kFALSE;
10a8ccbe 163 if(CheckVtxRange() && CheckCentralityCut() && CheckMultiplicityCut()){ //selection on vertex and Centrality
00493191 164 fIsSelected=CheckQVectorCut(); // QVector is calculated only if the centrality and vertex are correct (performance)
decf69d9 165 }
f6a38178 166 if(fIsSelected){
93db93de 167 ((TH1I*)fOutput->FindObject("fHistoCuts"))->Fill(kAcceptedEvents);
168 if(vertex)((TH1F*)fOutput->FindObject("fHistoVtxAftSel"))->Fill(vertex->GetZ());
f6a38178 169 }
ae0fdd7d 170 Int_t Nch=0;
f6a38178 171 for (Int_t iTracks = 0; iTracks < fAOD->GetNumberOfTracks(); iTracks++) {
172 AliAODTrack* track = fAOD->GetTrack(iTracks);
decf69d9 173 if (!fTrackCuts->IsSelected(track,kFALSE)) continue;
93db93de 174 ((TH1F*)fOutput->FindObject("fHistoEtaBefSel"))->Fill(track->Eta());
ae0fdd7d 175 if(fIsSelected){
93db93de 176 ((TH1F*)fOutput->FindObject("fHistoEtaAftSel"))->Fill(track->Eta());
ae0fdd7d 177 Nch++;
178 }
f6a38178 179 }
93db93de 180 if(fIsSelected)((TH1F*)fOutput->FindObject("fHistoNChAftSel"))->Fill(Nch);
f6a38178 181 return fIsSelected;
c88234ad 182}
183
184//______________________________________________________
185Bool_t AliSpectraAODEventCuts::CheckVtxRange()
186{
187 // reject events outside of range
f6a38178 188 AliAODVertex * vertex = fAOD->GetPrimaryVertex();
264de30e 189 //when moving to 2011 wìone has to add a cut using SPD vertex.
190 //The point is that for events with |z|>20 the vertexer tracks is not working (only 2011!). One has to put a safety cut using SPD vertex large e.g. 15cm
f6a38178 191 if (!vertex)
192 {
93db93de 193 ((TH1I*)fOutput->FindObject("fHistoCuts"))->Fill(kVtxNoEvent);
c88234ad 194 return kFALSE;
f6a38178 195 }
10a8ccbe 196 if (vertex->GetZ() > fVertexCutMin && vertex->GetZ() < fVertexCutMax)
f6a38178 197 {
c88234ad 198 return kTRUE;
f6a38178 199 }
93db93de 200 ((TH1I*)fOutput->FindObject("fHistoCuts"))->Fill(kVtxRange);
f6a38178 201 return kFALSE;
c88234ad 202}
203
204//______________________________________________________
205Bool_t AliSpectraAODEventCuts::CheckCentralityCut()
206{
f6a38178 207 // Check centrality cut
48c4d28e 208 fCent=-999.;
209 fCent=fAOD->GetCentrality()->GetCentralityPercentile(fCentralityMethod.Data());
210 if ( (fCent <= fCentralityCutMax) && (fCent >= fCentralityCutMin) ) return kTRUE;
93db93de 211 ((TH1I*)fOutput->FindObject("fHistoCuts"))->Fill(kVtxCentral);
f6a38178 212 return kFALSE;
c88234ad 213}
214
10a8ccbe 215//______________________________________________________
216Bool_t AliSpectraAODEventCuts::CheckMultiplicityCut()
217{
218 // Check multiplicity cut
31579941 219 // FIXME: why this is not tracket in the track stats histos?
10a8ccbe 220 Int_t Ncharged=0;
221 for (Int_t iTracks = 0; iTracks < fAOD->GetNumberOfTracks(); iTracks++){
222 AliAODTrack* track = fAOD->GetTrack(iTracks);
223 if (!fTrackCuts->IsSelected(track,kFALSE)) continue;
224 Ncharged++;
225 }
10a8ccbe 226 if(Ncharged>fMultiplicityCutMin && Ncharged<fMultiplicityCutMax)return kTRUE;
227
228 return kFALSE;
229}
230
decf69d9 231//______________________________________________________
232Bool_t AliSpectraAODEventCuts::CheckQVectorCut()
23f51e78 233{
48c4d28e 234 Double_t qxEPVZERO = -999., qyEPVZERO = -999.;
fdbff5a1 235 Double_t qVZERO = -999.;
236 Double_t psi=-999.;
237
238 if(fIsLHC10h){
239 qVZERO=CalculateQVectorLHC10h();
240 psi=fPsiV0A;
241 }else{
242 psi=fAOD->GetEventplane()->CalculateVZEROEventPlane(fAOD,10,2,qxEPVZERO,qyEPVZERO);//FIXME we can a flag for 2010 and 2011
243 qVZERO= TMath::Sqrt(qxEPVZERO*qxEPVZERO + qyEPVZERO*qyEPVZERO);
244 }
245
246 //cut on q vector
16ee3540 247 if(qVZERO<fQVectorCutMin || qVZERO>fQVectorCutMax)return kFALSE;
23f51e78 248 Double_t cent=fAOD->GetCentrality()->GetCentralityPercentile(fCentralityMethod.Data());
249 ((TH2F*)fOutput->FindObject("fHistoQVector"))->Fill(cent,qVZERO);
250 ((TH2F*)fOutput->FindObject("fHistoEP"))->Fill(cent,psi);
93db93de 251 ((TH1I*)fOutput->FindObject("fHistoCuts"))->Fill(kQVector);
23f51e78 252
fdbff5a1 253
decf69d9 254 return kTRUE;
decf69d9 255}
256
93db93de 257//______________________________________________________
93db93de 258Double_t AliSpectraAODEventCuts::CalculateQVectorLHC10h(){
259
260 Int_t run = fAOD->GetRunNumber();
261 if(run != fRun){
262 // Load the calibrations run dependent
263 if(OpenInfoCalbration(run))fRun=run;
48c4d28e 264 else{
265 fqV0C=-999.;
266 fqV0A=-999.;
267 return -999.;
268 }
93db93de 269 }
270
271 //V0 info
272 Double_t Qxa2 = 0, Qya2 = 0;
273 Double_t Qxc2 = 0, Qyc2 = 0;
274 Double_t sumMc = 0, sumMa = 0;
275
276 AliAODVZERO* aodV0 = fAOD->GetVZEROData();
277
278 for (Int_t iv0 = 0; iv0 < 64; iv0++) {
279
280 Double_t phiV0 = TMath::PiOver4()*(0.5 + iv0 % 8);
281
282 Float_t multv0 = aodV0->GetMultiplicity(iv0);
23f51e78 283 ((TH2F*)fOutput->FindObject("fV0M"))->Fill(iv0,multv0);
284
93db93de 285 if (iv0 < 32){
286
287 Double_t multCorC = -10;
288
289 if (iv0 < 8)
290 multCorC = multv0*fV0Cpol1/fMultV0->GetBinContent(iv0+1);
291 else if (iv0 >= 8 && iv0 < 16)
292 multCorC = multv0*fV0Cpol2/fMultV0->GetBinContent(iv0+1);
293 else if (iv0 >= 16 && iv0 < 24)
294 multCorC = multv0*fV0Cpol3/fMultV0->GetBinContent(iv0+1);
295 else if (iv0 >= 24 && iv0 < 32)
296 multCorC = multv0*fV0Cpol4/fMultV0->GetBinContent(iv0+1);
297
298 if (multCorC < 0){
299 cout<<"Problem with multiplicity in V0C"<<endl;
48c4d28e 300 fqV0C=-999.;
301 fqV0A=-999.;
302 return -999.;
93db93de 303 }
304
305 Qxc2 += TMath::Cos(2*phiV0) * multCorC;
306 Qyc2 += TMath::Sin(2*phiV0) * multCorC;
307
308 sumMc += multCorC;
23f51e78 309 ((TH2F*)fOutput->FindObject("fV0MCor"))->Fill(iv0,multCorC);
93db93de 310
311 } else {
312
313 Double_t multCorA = -10;
314
315 if (iv0 >= 32 && iv0 < 40)
316 multCorA = multv0*fV0Apol1/fMultV0->GetBinContent(iv0+1);
317 else if (iv0 >= 40 && iv0 < 48)
318 multCorA = multv0*fV0Apol2/fMultV0->GetBinContent(iv0+1);
319 else if (iv0 >= 48 && iv0 < 56)
320 multCorA = multv0*fV0Apol3/fMultV0->GetBinContent(iv0+1);
321 else if (iv0 >= 56 && iv0 < 64)
322 multCorA = multv0*fV0Apol4/fMultV0->GetBinContent(iv0+1);
323
324 if (multCorA < 0){
325 cout<<"Problem with multiplicity in V0A"<<endl;
48c4d28e 326 fqV0C=-999.;
327 fqV0A=-999.;
328 return -999.;
93db93de 329 }
330
331 Qxa2 += TMath::Cos(2*phiV0) * multCorA;
332 Qya2 += TMath::Sin(2*phiV0) * multCorA;
333
334 sumMa += multCorA;
23f51e78 335 ((TH2F*)fOutput->FindObject("fV0MCor"))->Fill(iv0,multCorA);
93db93de 336 }
337 }
338
339 Short_t centrV0 = GetCentrCode(fAOD);
340
6a6d2363 341 Double_t Qxamean2 = 0.;
342 Double_t Qyamean2 = 0.;
343 Double_t Qxcmean2 = 0.;
344 Double_t Qycmean2 = 0.;
345
346 if(centrV0!=-1){
347 Qxamean2 = fMeanQxa2[centrV0];
348 Qyamean2 = fMeanQya2[centrV0];
349 Qxcmean2 = fMeanQxc2[centrV0];
350 Qycmean2 = fMeanQyc2[centrV0];
351 }
93db93de 352
353 Double_t QxaCor2 = Qxa2 - Qxamean2*sumMa;
354 Double_t QyaCor2 = Qya2 - Qyamean2*sumMa;
355 Double_t QxcCor2 = Qxc2 - Qxcmean2*sumMc;
356 Double_t QycCor2 = Qyc2 - Qycmean2*sumMc;
357
fdbff5a1 358 fPsiV0A = TMath::ATan2(QyaCor2, QxaCor2)/2.;
359 fPsiV0C = TMath::ATan2(QycCor2, QxcCor2)/2.;
93db93de 360
fdbff5a1 361 ((TH2F*)fOutput->FindObject("fPsiACor"))->Fill((Float_t)fAOD->GetCentrality()->GetCentralityPercentile("V0M"), fPsiV0A);
362 ((TH2F*)fOutput->FindObject("fPsiCCor"))->Fill((Float_t)fAOD->GetCentrality()->GetCentralityPercentile("V0M"), fPsiV0C);
93db93de 363
48c4d28e 364 fqV0A = TMath::Sqrt((QxaCor2*QxaCor2 + QyaCor2*QyaCor2)/sumMa);
fdbff5a1 365 fqV0C = TMath::Sqrt((QxcCor2*QxcCor2 + QycCor2*QycCor2)/sumMc);
4fd2f0ad
LM		 366 fqV0Ax = QxaCor2*TMath::Sqrt(1./sumMa);
367 fqV0Cx = QxcCor2*TMath::Sqrt(1./sumMc);
368 fqV0Ay = QyaCor2*TMath::Sqrt(1./sumMa);
369 fqV0Cy = QycCor2*TMath::Sqrt(1./sumMc);
48c4d28e 370
371 ((TH2F*)fOutput->FindObject("fQVecACor"))->Fill((Float_t)fAOD->GetCentrality()->GetCentralityPercentile("V0M"), fqV0A);
372 ((TH2F*)fOutput->FindObject("fQVecCCor"))->Fill((Float_t)fAOD->GetCentrality()->GetCentralityPercentile("V0M"), fqV0C);
93db93de 373
fdbff5a1 374 return fqV0A; //FIXME we have to combine VZERO-A and C
93db93de 375}
376
fdbff5a1 377//______________________________________________________
93db93de 378Short_t AliSpectraAODEventCuts::GetCentrCode(AliVEvent* ev)
379{
380
381 Short_t centrCode = -1;
382
383 AliCentrality* centrality = 0;
384 AliAODEvent* aod = (AliAODEvent*)ev;
385 centrality = aod->GetHeader()->GetCentralityP();
386
387 Float_t centV0 = centrality->GetCentralityPercentile("V0M");
388 Float_t centTrk = centrality->GetCentralityPercentile("TRK");
389
390
391 if (TMath::Abs(centV0 - centTrk) < 5.0 && centV0 <= 80 && centV0 > 0){
392
393 if ((centV0 > 0) && (centV0 <= 5.0))
394 centrCode = 0;
395 else if ((centV0 > 5.0) && (centV0 <= 10.0))
396 centrCode = 1;
397 else if ((centV0 > 10.0) && (centV0 <= 20.0))
398 centrCode = 2;
399 else if ((centV0 > 20.0) && (centV0 <= 30.0))
400 centrCode = 3;
401 else if ((centV0 > 30.0) && (centV0 <= 40.0))
402 centrCode = 4;
403 else if ((centV0 > 40.0) && (centV0 <= 50.0))
404 centrCode = 5;
405 else if ((centV0 > 50.0) && (centV0 <= 60.0))
406 centrCode = 6;
407 else if ((centV0 > 60.0) && (centV0 <= 70.0))
408 centrCode = 7;
409 else if ((centV0 > 70.0) && (centV0 <= 80.0))
410 centrCode = 8;
411 }
412
413 return centrCode;
414
415}
416
c88234ad 417//______________________________________________________
418void AliSpectraAODEventCuts::PrintCuts()
419{
f6a38178 420 // print info about event cuts
421 cout << "Event Stats" << endl;
93db93de 422 cout << " > Trigger Selection: " << fSelectBit << endl;
423 cout << " > Centrality estimator: " << fCentralityMethod << endl;
424 cout << " > Number of accepted events: " << NumberOfEvents() << endl;
425 cout << " > Number of processed events: " << NumberOfProcessedEvents() << endl;
426 cout << " > Number of PhysSel events: " << NumberOfPhysSelEvents() << endl;
427 cout << " > Vertex out of range: " << ((TH1I*)fOutput->FindObject("fHistoCuts"))->GetBinContent(kVtxRange + 1) << endl;
428 cout << " > Events cut by centrality: " << ((TH1I*)fOutput->FindObject("fHistoCuts"))->GetBinContent(kVtxCentral + 1) << endl;
429 cout << " > Events without vertex: " << ((TH1I*)fOutput->FindObject("fHistoCuts"))->GetBinContent(kVtxNoEvent + 1) << endl;
430 cout << " > QVector cut: " << ((TH1I*)fOutput->FindObject("fHistoCuts"))->GetBinContent(kQVector + 1) << endl;
829b5a81 431 cout << " > Track type used for the QVector calculation: " << fTrackBits << endl;
16ee3540 432 cout << " > QRange: [" << fQVectorCutMin <<"," <<fQVectorCutMax<<"]"<< endl;
10a8ccbe 433 cout << " > Vertex: [" << fVertexCutMin <<"," <<fVertexCutMax<<"]"<< endl;
434 cout << " > Multiplicity: [" << fMultiplicityCutMin <<"," <<fMultiplicityCutMax<<"]"<< endl;
829b5a81 435 cout << " > Centrality: [" << fCentralityCutMin <<"," <<fCentralityCutMax<<"]"<< endl;
f6a38178 436}
c88234ad 437
93db93de 438//_____________________________________________________________________________
439Bool_t AliSpectraAODEventCuts::OpenInfoCalbration(Int_t run)
972a21ad 440{
93db93de 441
442 AliOADBContainer* cont = (AliOADBContainer*) fCalib->FindObject("hMultV0BefCorr");
443 if(!cont){
444 printf("OADB object hMultV0BefCorr is not available in the file\n");
445 return 0;
446 }
447
448 if(!(cont->GetObject(run))){
449 printf("OADB object hMultV0BefCorr is not available for run %i\n",run);
450 return 0;
451 }
452 fMultV0 = ((TH2F*) cont->GetObject(run))->ProfileX();
453
454 TF1* fpolc1 = new TF1("fpolc1","pol0", 0, 7);
455 fMultV0->Fit(fpolc1, "RN");
456 fV0Cpol1 = fpolc1->GetParameter(0);
457
458 TF1* fpolc2 = new TF1("fpolc2","pol0", 8, 15);
459 fMultV0->Fit(fpolc2, "RN");
460 fV0Cpol2 = fpolc2->GetParameter(0);
461
462 TF1* fpolc3 = new TF1("fpolc3","pol0", 16, 23);
463 fMultV0->Fit(fpolc3, "RN");
464 fV0Cpol3 = fpolc3->GetParameter(0);
465
466 TF1* fpolc4 = new TF1("fpolc4","pol0", 24, 31);
467 fMultV0->Fit(fpolc4, "RN");
468 fV0Cpol4 = fpolc4->GetParameter(0);
469
470 TF1* fpola1 = new TF1("fpola1","pol0", 32, 39);
6a6d2363 471 fMultV0->Fit(fpola1, "RN");
93db93de 472 fV0Apol1 = fpola1->GetParameter(0);
473
474 TF1* fpola2 = new TF1("fpola2","pol0", 40, 47);
475 fMultV0->Fit(fpola2, "RN");
476 fV0Apol2 = fpola2->GetParameter(0);
a162b3e2 477
93db93de 478 TF1* fpola3 = new TF1("fpola3","pol0", 48, 55);
479 fMultV0->Fit(fpola3, "RN");
480 fV0Apol3 = fpola3->GetParameter(0);
481
482 TF1* fpola4 = new TF1("fpola4","pol0", 56, 63);
483 fMultV0->Fit(fpola4, "RN");
484 fV0Apol4 = fpola4->GetParameter(0);
485
93db93de 486 for(Int_t i=0; i < 10; i++){
487
488 char nameQxa2[100];
6a6d2363 489 snprintf(nameQxa2,100, "hQxa2m_%i", i);
93db93de 490
491 char nameQya2[100];
6a6d2363 492 snprintf(nameQya2,100, "hQya2m_%i", i);
93db93de 493
494 char nameQxc2[100];
6a6d2363 495 snprintf(nameQxc2,100, "hQxc2m_%i", i);
93db93de 496
497 char nameQyc2[100];
6a6d2363 498 snprintf(nameQyc2,100, "hQyc2m_%i", i);
93db93de 499
500 AliOADBContainer* contQxa2 = (AliOADBContainer*) fCalib->FindObject(nameQxa2);
501 if(!contQxa2){
502 printf("OADB object %s is not available in the file\n", nameQxa2);
503 return 0;
504 }
505
506 if(!(contQxa2->GetObject(run))){
507 printf("OADB object %s is not available for run %i\n", nameQxa2, run);
508 return 0;
509 }
510
511 fMeanQxa2[i] = ((TH1F*) contQxa2->GetObject(run))->GetMean();
512
513
514 AliOADBContainer* contQya2 = (AliOADBContainer*) fCalib->FindObject(nameQya2);
515 if(!contQya2){
516 printf("OADB object %s is not available in the file\n", nameQya2);
517 return 0;
518 }
519
520 if(!(contQya2->GetObject(run))){
521 printf("OADB object %s is not available for run %i\n", nameQya2, run);
522 return 0;
523 }
524
525 fMeanQya2[i] = ((TH1F*) contQya2->GetObject(run))->GetMean();
526
527
528 AliOADBContainer* contQxc2 = (AliOADBContainer*) fCalib->FindObject(nameQxc2);
529 if(!contQxc2){
530 printf("OADB object %s is not available in the file\n", nameQxc2);
531 return 0;
532 }
533
534 if(!(contQxc2->GetObject(run))){
535 printf("OADB object %s is not available for run %i\n", nameQxc2, run);
536 return 0;
537 }
538
539 fMeanQxc2[i] = ((TH1F*) contQxc2->GetObject(run))->GetMean();
540
541
542 AliOADBContainer* contQyc2 = (AliOADBContainer*) fCalib->FindObject(nameQyc2);
543 if(!contQyc2){
544 printf("OADB object %s is not available in the file\n", nameQyc2);
545 return 0;
546 }
547
548 if(!(contQyc2->GetObject(run))){
549 printf("OADB object %s is not available for run %i\n", nameQyc2, run);
550 return 0;
551 }
552
553 fMeanQyc2[i] = ((TH1F*) contQyc2->GetObject(run))->GetMean();
554
555 }
556 return 1;
972a21ad 557}
558
559//______________________________________________________
560
561
c88234ad 562Long64_t AliSpectraAODEventCuts::Merge(TCollection* list)
563{
564 // Merge a list of AliSpectraAODEventCuts objects with this.
565 // Returns the number of merged objects (including this).
93db93de 566
567 AliInfo("Merging");
568
c88234ad 569 if (!list)
570 return 0;
571
572 if (list->IsEmpty())
573 return 1;
93db93de 574
c88234ad 575 TIterator* iter = list->MakeIterator();
576 TObject* obj;
93db93de 577
c88234ad 578 // collections of all histograms
93db93de 579 TList collections;
580
c88234ad 581 Int_t count = 0;
582
583 while ((obj = iter->Next())) {
584 AliSpectraAODEventCuts* entry = dynamic_cast<AliSpectraAODEventCuts*> (obj);
585 if (entry == 0)
586 continue;
587
93db93de 588 TList * l = entry->GetOutputList();
589 collections.Add(l);
c88234ad 590 count++;
591 }
592
93db93de 593 fOutput->Merge(&collections);
c88234ad 594
595 delete iter;
596
597 return count+1;
598}
599
7b364a00 600//______________________________________________________
601Double_t AliSpectraAODEventCuts::GetQvecPercentile(Int_t v0side){
602
603 //check Qvec and Centrality consistency
604 if(fCent>90.) return -999.;
605 if(v0side==0/*V0A*/){ if(fqV0A== -999.) return -999.; }
606 if(v0side==1/*V0C*/){ if(fqV0C== -999.) return -999.; }
607
608 fQvecIntegral = 0x0;
609
610 if(v0side==0/*V0A*/){ fQvecIntegral = (TH2D*)fQvecIntList->FindObject("VZEROA"); }
611 if(v0side==1/*V0C*/){ fQvecIntegral = (TH2D*)fQvecIntList->FindObject("VZEROC"); }
612
613 Double_t ic = (Int_t)fCent; //fQvecIntegral: 1% centrality bin
614
615 TH1D *h1D = (TH1D*)fQvecIntegral->ProjectionY("h1D",ic+1,ic+1);
616
617 TSpline *spline = 0x0;
618
619 if(v0side==0/*V0A*/){
620 if( CheckSplineArray(fSplineArrayV0A) ) {
621 spline = (TSpline*)fSplineArrayV0A->At(ic);
622 cout<<"Qvec V0A - ic: "<<ic<<" - Found TSpline..."<<endl;
623 } else {
624 spline = new TSpline3(h1D,"sp3");
625 fSplineArrayV0A->AddAtAndExpand(spline,ic);
626 cout<<"Qvec V0A - ic: "<<ic<<" - TSpline not found. Creating..."<<endl;
627 }
628 }
629 else if(v0side==1/*V0C*/){
630 if( CheckSplineArray(fSplineArrayV0C) ) {
631 spline = (TSpline*)fSplineArrayV0C->At(ic);
632 } else {
633 spline = new TSpline3(h1D,"sp3");
634 fSplineArrayV0C->AddAtAndExpand(spline,ic);
635 }
636 }
637
638 Double_t percentile = -999.;
639 if(v0side==0/*V0A*/){ percentile = 100*spline->Eval(fqV0A); }
640 if(v0side==1/*V0C*/){ percentile = 100*spline->Eval(fqV0C); }
641
642 if(percentile>100. || percentile<0.) return -999.;
643
644 return percentile;
645}
646
647//______________________________________________________
648Bool_t AliSpectraAODEventCuts::CheckSplineArray(TObjArray * splarr){
649 //check TSpline array consistency
650
651 Int_t n = (Int_t)fCent;
652
653 if(splarr->GetSize()<n) return kFALSE;
654
655 if(!splarr->At(n)) return kFALSE;
656
657 return kTRUE;
658
659}