1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
17 // VZERO event plane task for 2010
18 // Gain equalization + Recentering
22 // Raphaelle Bailhache <R.Bailhache@gsi.de>
26 #include "AliHFEVZEROEventPlane.h"
29 #include "AliAnalysisManager.h"
30 #include "AliInputEventHandler.h"
31 #include "AliVEvent.h"
32 //#include "AliESDEvent.h"
33 #include "AliCentrality.h"
34 //#include "AliESDVZERO.h"
35 #include "AliVVZERO.h"
37 #include "AliOADBContainer.h"
42 #include "THnSparse.h"
45 #include "AliEventplane.h"
53 //_____________________________________________________________________________
54 AliHFEVZEROEventPlane::AliHFEVZEROEventPlane():
56 fEventPlaneV0A(-100.0),
57 fEventPlaneV0C(-100.0),
58 fEventPlaneV0(-100.0),
64 fnamefile(TString("")),
70 // Default constructor (should not be used)
73 for(Int_t k = 0; k < fgknCentrBin; k++) {
74 for(Int_t iside = 0; iside < 2; iside++) {
75 for(Int_t icoord = 0; icoord < 2; icoord++) {
76 fQBefore[k][iside][icoord] = 0x0;
77 fQAfter[k][iside][icoord] = 0x0;
84 //______________________________________________________________________________
85 AliHFEVZEROEventPlane::AliHFEVZEROEventPlane(const char *name, const Char_t *title):
87 fEventPlaneV0A(-100.0),
88 fEventPlaneV0C(-100.0),
89 fEventPlaneV0(-100.0),
95 fnamefile(TString("")),
104 TString namelist("QA_hist_");
106 fOutputList = new TList();
107 fOutputList->SetName((const char*)namelist);
108 fOutputList->SetOwner(kTRUE);
111 fMultV0Before = new TProfile("MultV0Before","",64,0,64);
112 fMultV0Before->Sumw2();
113 fMultV0After = new TProfile("MultV0After","",64,0,64);
114 fMultV0After->Sumw2();
115 fOutputList->Add(fMultV0Before);
116 fOutputList->Add(fMultV0After);
119 for(Int_t k = 0; k < fgknCentrBin; k++) {
120 for(Int_t iside = 0; iside < 2; iside++) {
121 for(Int_t icoord = 0; icoord < 2; icoord++) {
123 TString namecontbefore("");
124 TString namecontafter("");
126 if(iside==0 && icoord==0) {
127 namecontbefore += "hQxc2_";
130 else if(iside==1 && icoord==0){
131 namecontbefore += "hQxa2_";
134 else if(iside==0 && icoord==1){
135 namecontbefore += "hQyc2_";
138 else if(iside==1 && icoord==1) {
139 namecontbefore += "hQya2_";
143 if(iside==0 && icoord==0){
144 namecontafter += "hQxc2_";
146 namecontafter += "_after";
148 else if(iside==1 && icoord==0){
149 namecontafter += "hQxa2_";
151 namecontafter += "_after";
153 else if(iside==0 && icoord==1) {
154 namecontafter += "hQyc2_";
156 namecontafter += "_after";
158 else if(iside==1 && icoord==1) {
159 namecontafter += "hQya2_";
161 namecontafter += "_after";
164 fQBefore[k][iside][icoord] = new TH1F(((const char*)namecontbefore),"",1600,-800.0,800.0);
165 fQBefore[k][iside][icoord]->Sumw2();
166 fOutputList->Add(fQBefore[k][iside][icoord]);
167 fQAfter[k][iside][icoord] = new TH1F(((const char*)namecontafter),"",1600,-800.0,800.0);
168 fQAfter[k][iside][icoord]->Sumw2();
169 fOutputList->Add(fQAfter[k][iside][icoord]);
176 //____________________________________________________________
177 AliHFEVZEROEventPlane::AliHFEVZEROEventPlane(const AliHFEVZEROEventPlane &ref):
179 fEventPlaneV0A(ref.fEventPlaneV0A),
180 fEventPlaneV0C(ref.fEventPlaneV0C),
181 fEventPlaneV0(ref.fEventPlaneV0),
187 fnamefile(ref.fnamefile),
198 //____________________________________________________________
199 AliHFEVZEROEventPlane &AliHFEVZEROEventPlane::operator=(const AliHFEVZEROEventPlane &ref){
201 // Assignment operator
208 //____________________________________________________________
209 void AliHFEVZEROEventPlane::Copy(TObject &o) const {
211 // Copy into object o
213 AliHFEVZEROEventPlane &target = dynamic_cast<AliHFEVZEROEventPlane &>(o);
215 target.fEventPlaneV0A = fEventPlaneV0A;
216 target.fEventPlaneV0C = fEventPlaneV0C;
217 target.fEventPlaneV0 = fEventPlaneV0;
218 target.fnamefile = fnamefile;
222 //__________________________________________________________________
223 AliHFEVZEROEventPlane::~AliHFEVZEROEventPlane(){
228 fOutputList->Delete();
234 //______________________________________________________________________________
235 void AliHFEVZEROEventPlane::ProcessEvent(AliVEvent *event)
242 fEventPlaneV0A = -100.0;
243 fEventPlaneV0C = -100.0;
244 fEventPlaneV0 = -100.0;
247 Int_t run = event->GetRunNumber();
250 if(fMultV0) delete fMultV0;
252 // Load the calibrations run dependent
253 if(!OpenInfoCalbration(run)) ok = kFALSE;
257 if(ok) Analyze(event);
261 //________________________________________________________________________
262 void AliHFEVZEROEventPlane::Analyze(AliVEvent* esdEvent)
265 // Do VZERO calibration + centering
269 //printf("Did not find calibration VZERO\n");
272 //printf("HERE!!!\n");
275 Float_t v0Centr = -10.;
276 AliCentrality* centrality = esdEvent->GetCentrality();
278 v0Centr = centrality->GetCentralityPercentile("V0M");
281 // Analyse only for 0-80% PbPb collisions
283 if (v0Centr >0 && v0Centr < 80){
284 if(v0Centr < 5) iC = 0;
285 else if(v0Centr < 10) iC = 1;
286 else if(v0Centr < 20) iC = 2;
287 else if(v0Centr < 30) iC = 3;
288 else if(v0Centr < 40) iC = 4;
289 else if(v0Centr < 50) iC = 5;
290 else if(v0Centr < 60) iC = 6;
291 else if(v0Centr < 70) iC = 7;
295 Double_t qxa2 = 0, qya2 = 0;
296 Double_t qxc2 = 0, qyc2 = 0;
299 AliVVZERO* esdV0 = esdEvent->GetVZEROData();
301 for (Int_t iv0 = 0; iv0 < 64; iv0++) {
302 Double_t phiV0 = TMath::PiOver4()*(0.5 + iv0 % 8);
303 Float_t multv0 = esdV0->GetMultiplicity(iv0);
305 //printf("Value %f\n",fMultV0->GetBinContent(iv0+1));
306 if(fMultV0->GetBinContent(iv0+1)>0.0) {
307 qxc2 += TMath::Cos(2*phiV0) * multv0*fV0Cpol/fMultV0->GetBinContent(iv0+1);
308 qyc2 += TMath::Sin(2*phiV0) * multv0*fV0Cpol/fMultV0->GetBinContent(iv0+1);
309 fMultV0After->Fill(iv0,multv0*fV0Cpol/fMultV0->GetBinContent(iv0+1));
312 if(fMultV0->GetBinContent(iv0+1)>0.0) {
313 qxa2 += TMath::Cos(2*phiV0) * multv0*fV0Apol/fMultV0->GetBinContent(iv0+1);
314 qya2 += TMath::Sin(2*phiV0) * multv0*fV0Apol/fMultV0->GetBinContent(iv0+1);
315 fMultV0After->Fill(iv0,multv0*fV0Apol/fMultV0->GetBinContent(iv0+1));
318 fMultV0Before->Fill(iv0,multv0);
322 fQBefore[iC][1][0] -> Fill(qxa2);
323 fQBefore[iC][1][1] -> Fill(qya2);
324 fQBefore[iC][0][0] -> Fill(qxc2);
325 fQBefore[iC][0][1] -> Fill(qyc2);
327 //grab for each centrality the proper histo with the Qx and Qy to do the recentering
328 Double_t qxamean2 = fMeanQ[iC][1][0];
329 Double_t qxarms2 = fWidthQ[iC][1][0];
330 Double_t qyamean2 = fMeanQ[iC][1][1];
331 Double_t qyarms2 = fWidthQ[iC][1][1];
333 Double_t qxcmean2 = fMeanQ[iC][0][0];
334 Double_t qxcrms2 = fWidthQ[iC][0][0];
335 Double_t qycmean2 = fMeanQ[iC][0][1];
336 Double_t qycrms2 = fWidthQ[iC][0][1];
338 if((TMath::Abs(qxarms2) < 0.00000001) || (TMath::Abs(qyarms2) < 0.00000001) || (TMath::Abs(qxcrms2) < 0.00000001) || (TMath::Abs(qycrms2) < 0.00000001)) return;
340 Double_t qxaCor2 = (qxa2 - qxamean2)/qxarms2;
341 Double_t qyaCor2 = (qya2 - qyamean2)/qyarms2;
342 Double_t qxcCor2 = (qxc2 - qxcmean2)/qxcrms2;
343 Double_t qycCor2 = (qyc2 - qycmean2)/qycrms2;
345 Double_t qxCor2 = qxaCor2 + qxcCor2;
346 Double_t qyCor2 = qyaCor2 + qycCor2;
349 fQAfter[iC][1][0] -> Fill(qxaCor2);
350 fQAfter[iC][1][1] -> Fill(qyaCor2);
351 fQAfter[iC][0][0] -> Fill(qxcCor2);
352 fQAfter[iC][0][1] -> Fill(qycCor2);
354 Double_t evPlAngV0ACor2 = TVector2::Phi_0_2pi(TMath::ATan2(qyaCor2, qxaCor2))/2.;
355 Double_t evPlAngV0CCor2 = TVector2::Phi_0_2pi(TMath::ATan2(qycCor2, qxcCor2))/2.;
356 Double_t evPlAngV0Cor2 = TVector2::Phi_0_2pi(TMath::ATan2(qyCor2, qxCor2))/2.;
358 fEventPlaneV0A = evPlAngV0ACor2;
359 fEventPlaneV0C = evPlAngV0CCor2;
360 fEventPlaneV0 = evPlAngV0Cor2;
362 //printf("Eventplane V0A %f, V0C %f\n",fEventPlaneV0A,fEventPlaneV0C);
366 //_____________________________________________________________________________
367 Bool_t AliHFEVZEROEventPlane::OpenInfoCalbration(Int_t run)
370 // Take the calibration coefficients
373 //printf("Name of the file %s\n",(const char*)fnamefile);
375 TFile *foadb = TFile::Open(fnamefile.Data());
377 printf("OADB file %s cannot be opened\n",fnamefile.Data());
383 AliOADBContainer *cont = (AliOADBContainer*) foadb->Get("hMultV0BefCorr");
385 printf("OADB object hMultV0BefCorr is not available in the file\n");
389 if(!(cont->GetObject(run))){
390 printf("OADB object hMultV0BefCorr is not available for run %i\n",run);
393 //TProfile *multV0 = ((TH2F *) cont->GetObject(run))->ProfileX();
394 TProfile *multV0 = ((TProfile *) cont->GetObject(run));
395 if(!multV0) return kFALSE;
396 fMultV0 = (TProfile *) multV0->Clone();
397 fMultV0->SetDirectory(0);
399 TF1 *fpol0 = new TF1("fpol0","pol0");
400 fMultV0->Fit(fpol0,"Q0","",0,31);
401 fV0Cpol = fpol0->GetParameter(0);
402 fMultV0->Fit(fpol0,"Q0","",32,64);
403 fV0Apol = fpol0->GetParameter(0);
405 for(Int_t iside=0;iside<2;iside++){
406 for(Int_t icoord=0;icoord<2;icoord++){
407 for(Int_t i=0;i < fgknCentrBin;i++){
409 TString namecont("");
410 if(iside==0 && icoord==0) {
411 namecont += "hQxc2_";
414 else if(iside==1 && icoord==0){
415 namecont += "hQxa2_";
418 else if(iside==0 && icoord==1){
419 namecont += "hQyc2_";
422 else if(iside==1 && icoord==1) {
423 namecont += "hQya2_";
427 cont = (AliOADBContainer*) foadb->Get(namecont);
429 printf("OADB object %s is not available in the file\n",(const char*)namecont);
434 if(!(cont->GetObject(run))){
435 printf("OADB object %s is not available for run %i\n",(const char*)namecont,run);
439 fMeanQ[i][iside][icoord] = ((TH1F *) cont->GetObject(run))->GetMean();
440 fWidthQ[i][iside][icoord] = ((TH1F *) cont->GetObject(run))->GetRMS();