1 /*************************************************************************
2 * Copyright(c) 1998-2008, 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 **************************************************************************/
16 #define AliFlowAnalysisWithMCEventPlane_cxx
18 #include "Riostream.h" //needed as include
19 #include "TFile.h" //needed as include
20 #include "TProfile.h" //needed as include
21 #include "TComplex.h" //needed as include
26 #include "AliFlowCommonConstants.h" //needed as include
27 #include "AliFlowEventSimple.h"
28 #include "AliFlowTrackSimple.h"
29 #include "AliFlowCommonHist.h"
30 #include "AliFlowCommonHistResults.h"
31 #include "AliFlowAnalysisWithMCEventPlane.h"
35 // AliFlowAnalysisWithMCEventPlane:
36 // Description: Maker to analyze Flow from the generated MC reaction plane.
37 // This class is used to get the real value of the flow
38 // to compare the other methods to when analysing simulated events
39 // author: N. van der Kolk (kolk@nikhef.nl)
41 ClassImp(AliFlowAnalysisWithMCEventPlane)
43 //-----------------------------------------------------------------------
45 AliFlowAnalysisWithMCEventPlane::AliFlowAnalysisWithMCEventPlane():
52 fCommonHistsRes(NULL),
55 fHistProIntFlowRP(NULL),
56 fHistProDiffFlowPtRP(NULL),
57 fHistProDiffFlowEtaRP(NULL),
58 fHistProDiffFlowPtPOI(NULL),
59 fHistProDiffFlowEtaPOI(NULL)
63 fHistList = new TList();
65 fQsum = new TVector2; // flow vector sum
69 //-----------------------------------------------------------------------
72 AliFlowAnalysisWithMCEventPlane::~AliFlowAnalysisWithMCEventPlane()
79 //-----------------------------------------------------------------------
81 void AliFlowAnalysisWithMCEventPlane::WriteHistograms(TString* outputFileName)
83 //store the final results in output .root file
84 TFile *output = new TFile(outputFileName->Data(),"RECREATE");
85 output->WriteObject(fHistList, "cobjMCEP","SingleKey");
89 //-----------------------------------------------------------------------
91 void AliFlowAnalysisWithMCEventPlane::WriteHistograms(TString outputFileName)
93 //store the final results in output .root file
94 TFile *output = new TFile(outputFileName.Data(),"RECREATE");
95 output->WriteObject(fHistList, "cobjMCEP","SingleKey");
99 //-----------------------------------------------------------------------
100 void AliFlowAnalysisWithMCEventPlane::Init() {
102 //Define all histograms
103 cout<<"---Analysis with the real MC Event Plane---"<<endl;
105 Int_t iNbinsPt = AliFlowCommonConstants::GetNbinsPt();
106 Double_t dPtMin = AliFlowCommonConstants::GetPtMin();
107 Double_t dPtMax = AliFlowCommonConstants::GetPtMax();
109 Int_t iNbinsEta = AliFlowCommonConstants::GetNbinsEta();
110 Double_t dEtaMin = AliFlowCommonConstants::GetEtaMin();
111 Double_t dEtaMax = AliFlowCommonConstants::GetEtaMax();
113 fCommonHists = new AliFlowCommonHist("AliFlowCommonHistMCEP");
114 fHistList->Add(fCommonHists);
115 fCommonHistsRes = new AliFlowCommonHistResults("AliFlowCommonHistResultsMCEP");
116 fHistList->Add(fCommonHistsRes);
118 fHistProFlow = new TProfile("FlowPro_VPt_MCEP","FlowPro_VPt_MCEP",iNbinsPt,dPtMin,dPtMax);
119 fHistProFlow->SetXTitle("P_{t}");
120 fHistProFlow->SetYTitle("v_{2}");
121 fHistList->Add(fHistProFlow);
123 fHistRP = new TH1F("Flow_RP_MCEP","Flow_RP_MCEP",100,0.,3.14);
124 fHistRP->SetXTitle("Reaction Plane Angle");
125 fHistRP->SetYTitle("Counts");
126 fHistList->Add(fHistRP);
128 fHistProIntFlowRP = new TProfile("fHistProIntFlowRP","fHistProIntFlowRP",1,0.,1.);
129 fHistProIntFlowRP->SetLabelSize(0.06);
130 (fHistProIntFlowRP->GetXaxis())->SetBinLabel(1,"v_{n}{2}");
131 fHistProIntFlowRP->SetYTitle("");
132 fHistList->Add(fHistProIntFlowRP);
134 fHistProDiffFlowPtRP = new TProfile("fHistProDiffFlowPtRP","fHistProDiffFlowPtRP",iNbinsPt,dPtMin,dPtMax);
135 fHistProDiffFlowPtRP->SetXTitle("P_{t}");
136 fHistProDiffFlowPtRP->SetYTitle("");
137 fHistList->Add(fHistProDiffFlowPtRP);
139 fHistProDiffFlowEtaRP = new TProfile("fHistProDiffFlowEtaRP","fHistProDiffFlowEtaRP",iNbinsEta,dEtaMin,dEtaMax);
140 fHistProDiffFlowEtaRP->SetXTitle("#eta");
141 fHistProDiffFlowEtaRP->SetYTitle("");
142 fHistList->Add(fHistProDiffFlowEtaRP);
144 fHistProDiffFlowPtPOI = new TProfile("fHistProDiffFlowPtPOI","fHistProDiffFlowPtPOI",iNbinsPt,dPtMin,dPtMax);
145 fHistProDiffFlowPtPOI->SetXTitle("P_{t}");
146 fHistProDiffFlowPtPOI->SetYTitle("");
147 fHistList->Add(fHistProDiffFlowPtPOI);
149 fHistProDiffFlowEtaPOI = new TProfile("fHistProDiffFlowEtaPOI","fHistProDiffFlowEtaPOI",iNbinsEta,dEtaMin,dEtaMax);
150 fHistProDiffFlowEtaPOI->SetXTitle("#eta");
151 fHistProDiffFlowEtaPOI->SetYTitle("");
152 fHistList->Add(fHistProDiffFlowEtaPOI);
154 fEventNumber = 0; //set number of events to zero
158 //-----------------------------------------------------------------------
160 void AliFlowAnalysisWithMCEventPlane::Make(AliFlowEventSimple* anEvent, Double_t aRP) {
162 //Calculate v2 from the MC reaction plane
165 //fill control histograms
166 fCommonHists->FillControlHistograms(anEvent);
168 //get the Q vector from the FlowEvent
169 AliFlowVector vQ = anEvent->GetQ();
170 //cout<<"vQ.Mod() = " << vQ.Mod() << endl;
171 //for chi calculation:
173 //cout<<"fQsum.Mod() = "<<fQsum.Mod()<<endl;
175 //cout<<"fQ2sum = "<<fQ2sum<<endl;
183 //Double_t dPi = TMath::Pi();
186 //loop over the tracks of the event
187 Int_t iNumberOfTracks = anEvent->NumberOfTracks();
188 for (Int_t i=0;i<iNumberOfTracks;i++)
190 AliFlowTrackSimple* pTrack = anEvent->GetTrack(i) ;
192 if (pTrack->UseForIntegratedFlow()){
193 dPhi = pTrack->Phi();
194 dv2 = TMath::Cos(2*(dPhi-aRP));
196 dEta = pTrack->Eta();
197 //integrated flow (RP):
198 fHistProIntFlowRP->Fill(0.,dv2);
199 //differential flow (Pt, RP):
200 fHistProDiffFlowPtRP->Fill(dPt,dv2,1.);
201 //differential flow (Eta, RP):
202 fHistProDiffFlowEtaRP->Fill(dEta,dv2,1.);
204 if (pTrack->UseForDifferentialFlow()) {
205 dPhi = pTrack->Phi();
206 //if (dPhi<0.) dPhi+=2*TMath::Pi();
208 dv2 = TMath::Cos(2*(dPhi-aRP));
210 dEta = pTrack->Eta();
212 fHistProFlow->Fill(dPt,dv2);//to be removed
213 //differential flow (Pt, POI):
214 fHistProDiffFlowPtPOI->Fill(dPt,dv2,1.);
215 //differential flow (Eta, POI):
216 fHistProDiffFlowEtaPOI->Fill(dEta,dv2,1.);
222 // cout<<"@@@@@ "<<fEventNumber<<" events processed"<<endl;
226 //--------------------------------------------------------------------
227 void AliFlowAnalysisWithMCEventPlane::Finish() {
229 //*************make histograms etc.
230 if (fDebug) cout<<"AliFlowAnalysisWithMCEventPlane::Terminate()"<<endl;
232 Int_t iNbinsPt = AliFlowCommonConstants::GetNbinsPt();
233 Int_t iNbinsEta = AliFlowCommonConstants::GetNbinsEta();
237 Double_t dVRP = fHistProIntFlowRP->GetBinContent(1);
238 Double_t dErrVRP = fHistProIntFlowRP->GetBinError(1);//to be improved (treatment of errors for non-Gaussian distribution needed!)
239 fCommonHistsRes->FillIntegratedFlowRP(dVRP,dErrVRP);
241 //differential flow (Pt):
242 Double_t dvPtRP = 0.;
243 Double_t dErrvPtRP = 0.;
244 for(Int_t b=0;b<iNbinsPt;b++)
246 dvPtRP = fHistProDiffFlowPtRP->GetBinContent(b+1);
247 dErrvPtRP = fHistProDiffFlowPtRP->GetBinError(b+1);//to be improved (treatment of errors for non-Gaussian distribution needed!)
248 fCommonHistsRes->FillDifferentialFlowPtRP(b+1, dvPtRP, dErrvPtRP);
251 //differential flow (Eta):
252 Double_t dvEtaRP = 0.;
253 Double_t dErrvEtaRP = 0.;
254 for(Int_t b=0;b<iNbinsEta;b++)
256 dvEtaRP = fHistProDiffFlowEtaRP->GetBinContent(b+1);
257 dErrvEtaRP = fHistProDiffFlowEtaRP->GetBinError(b+1);//to be improved (treatment of errors for non-Gaussian distribution needed!)
258 fCommonHistsRes->FillDifferentialFlowEtaRP(b+1, dvEtaRP, dErrvEtaRP);
262 TH1F* fHistPtDiff = fCommonHists->GetHistPtDiff();
263 Double_t dYieldPt = 0.;
267 Double_t dv2proPt = 0.;
268 Double_t dErrdifcombPt = 0.;
269 Double_t dv2proEta = 0.;
270 Double_t dErrdifcombEta = 0.;
272 if(fHistProFlow && fHistProDiffFlowPtPOI) {//to be removed (fHistProFlow)
273 for(Int_t b=0;b<iNbinsPt;b++){
274 //dv2pro = fHistProFlow->GetBinContent(b);//to be removed
275 //dErrdifcomb = fHistProFlow->GetBinError(b);//to be removed
276 dv2proPt = fHistProDiffFlowPtPOI->GetBinContent(b+1);
277 dErrdifcombPt = fHistProDiffFlowPtPOI->GetBinError(b+1);//to be improved (treatment of errors for non-Gaussian distribution needed!)
279 fCommonHistsRes->FillDifferentialFlow(b+1, dv2proPt, dErrdifcombPt);//to be removed
280 fCommonHistsRes->FillDifferentialFlowPtPOI(b+1, dv2proPt, dErrdifcombPt);
282 //integrated flow (POI)
283 dYieldPt = fHistPtDiff->GetBinContent(b+1);
284 dV += dv2proPt*dYieldPt;
286 //error on integrated flow
287 dErrV += dYieldPt*dYieldPt*dErrdifcombPt*dErrdifcombPt;
289 }//end of for(Int_t b=0;b<iNbinsPt;b++)
290 } else { cout<<"fHistProFlow is NULL"<<endl; }
292 dV /= dSum; //because pt distribution should be normalised
293 dErrV /= (dSum*dSum);
294 dErrV = TMath::Sqrt(dErrV);
296 cout<<"dV{MC} is "<<dV<<" +- "<<dErrV<<endl;
297 fCommonHistsRes->FillIntegratedFlow(dV,dErrV);//to be removed
298 fCommonHistsRes->FillIntegratedFlowPOI(dV,dErrV);
301 if(fHistProDiffFlowEtaPOI)
303 for(Int_t b=0;b<iNbinsEta;b++)
305 dv2proEta = fHistProDiffFlowEtaPOI->GetBinContent(b+1);
306 dErrdifcombEta = fHistProDiffFlowEtaPOI->GetBinError(b+1);//to be improved (treatment of errors for non-Gaussian distribution needed!)
307 //fill common hist results:
308 fCommonHistsRes->FillDifferentialFlowEtaPOI(b+1, dv2proEta, dErrdifcombEta);
312 cout<<".....finished"<<endl;