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 **************************************************************************/
16 #include "Riostream.h" //needed as include
17 #include "AliFlowCommonConstants.h" //needed as include
18 #include "AliFlowCommonHist.h"
19 #include "AliFlowEventSimple.h"
20 #include "AliFlowTrackSimple.h"
24 #include "TMath.h" //needed as include
26 #include "AliFlowVector.h"
34 // Description: Class to organise common histograms for Flow Analysis
36 // authors: N. van der Kolk (kolk@nikhef.nl), A. Bilandzic (anteb@nikhef.nl), RS
39 ClassImp(AliFlowCommonHist)
41 //-----------------------------------------------------------------------
43 AliFlowCommonHist::AliFlowCommonHist():TNamed(),
53 fHistProMeanPtperBin(NULL),
62 AliFlowCommonHist::AliFlowCommonHist(const AliFlowCommonHist& a):
64 fHistMultOrig(new TH1F(*a.fHistMultOrig)),
65 fHistMultInt(new TH1F(*a.fHistMultInt)),
66 fHistMultDiff(new TH1F(*a.fHistMultDiff)),
67 fHistPtInt(new TH1F(*a.fHistPtInt)),
68 fHistPtDiff(new TH1F(*a.fHistPtDiff)),
69 fHistPhiInt(new TH1F(*a.fHistPhiInt)),
70 fHistPhiDiff(new TH1F(*a.fHistPhiDiff)),
71 fHistEtaInt(new TH1F(*a.fHistEtaInt)),
72 fHistEtaDiff(new TH1F(*a.fHistEtaDiff)),
73 fHistProMeanPtperBin(new TProfile(*a.fHistProMeanPtperBin)),
74 fHistQ(new TH1F(*a.fHistQ)),
79 fHistList = new TList();
80 fHistList-> Add(fHistMultOrig);
81 fHistList-> Add(fHistMultInt);
82 fHistList-> Add(fHistMultDiff);
83 fHistList-> Add(fHistPtInt);
84 fHistList-> Add(fHistPtDiff);
85 fHistList-> Add(fHistPhiInt);
86 fHistList-> Add(fHistPhiDiff);
87 fHistList-> Add(fHistEtaInt);
88 fHistList-> Add(fHistEtaDiff);
89 fHistList-> Add(fHistProMeanPtperBin);
90 fHistList-> Add(fHistQ);
91 // TListIter next = TListIter(a.fHistList);
95 // AliFlowCommonHist& AliFlowCommonHist::operator=(const AliFlowCommonHist& a)
97 // *fHistMultOrig = *a.fHistMultOrig;
98 // *fHistMultInt = *a.fHistMultInt;
99 // *fHistMultDiff = *a.fHistMultDiff;
100 // *fHistPtInt = *a.fHistPtInt;
101 // *fHistPtDiff = *a.fHistPtDiff;
102 // *fHistPhiInt = *a.fHistPhiInt;
103 // *fHistPhiDiff = *a.fHistPhiDiff;
104 // *fHistEtaInt = *a.fHistEtaInt;
105 // *fHistEtaDiff = *a.fHistEtaDiff;
106 // *fHistProMeanPtperBin = *a.fHistProMeanPtperBin;
107 // *fHistQ = *a.fHistQ;
108 // // *fHistList = *a.fHistList;
114 //-----------------------------------------------------------------------
116 AliFlowCommonHist::AliFlowCommonHist(const char *anInput,const char *title):TNamed(anInput,title),
126 fHistProMeanPtperBin(NULL),
131 //constructor creating histograms
132 Int_t iNbinsMult = AliFlowCommonConstants::GetNbinsMult();
133 Int_t iNbinsPt = AliFlowCommonConstants::GetNbinsPt();
134 Int_t iNbinsPhi = AliFlowCommonConstants::GetNbinsPhi();
135 Int_t iNbinsEta = AliFlowCommonConstants::GetNbinsEta();
136 Int_t iNbinsQ = AliFlowCommonConstants::GetNbinsQ();
139 Double_t dMultMin = AliFlowCommonConstants::GetMultMin();
140 Double_t dMultMax = AliFlowCommonConstants::GetMultMax();
141 Double_t dPtMin = AliFlowCommonConstants::GetPtMin();
142 Double_t dPtMax = AliFlowCommonConstants::GetPtMax();
143 Double_t dPhiMin = AliFlowCommonConstants::GetPhiMin();
144 Double_t dPhiMax = AliFlowCommonConstants::GetPhiMax();
145 Double_t dEtaMin = AliFlowCommonConstants::GetEtaMin();
146 Double_t dEtaMax = AliFlowCommonConstants::GetEtaMax();
147 Double_t dQMin = AliFlowCommonConstants::GetQMin();
148 Double_t dQMax = AliFlowCommonConstants::GetQMax();
150 cout<<"The settings for the common histograms are as follows:"<<endl;
151 cout<<"Multiplicity: "<<iNbinsMult<<" bins between "<<dMultMin<<" and "<<dMultMax<<endl;
152 cout<<"Pt: "<<iNbinsPt<<" bins between "<<dPtMin<<" and "<<dPtMax<<endl;
153 cout<<"Phi: "<<iNbinsPhi<<" bins between "<<dPhiMin<<" and "<<dPhiMax<<endl;
154 cout<<"Eta: "<<iNbinsEta<<" bins between "<<dEtaMin<<" and "<<dEtaMax<<endl;
155 cout<<"Q: "<<iNbinsQ<<" bins between "<<dQMin<<" and "<<dQMax<<endl;
158 sName = "Control_Flow_OrigMult_";
160 fHistMultOrig = new TH1F(sName.Data(), sName.Data(),iNbinsMult, dMultMin, dMultMax);
161 fHistMultOrig ->SetXTitle("Original Multiplicity");
162 fHistMultOrig ->SetYTitle("Counts");
164 sName = "Control_Flow_MultInt_";
166 fHistMultInt = new TH1F(sName.Data(), sName.Data(),iNbinsMult, dMultMin, dMultMax);
167 fHistMultInt ->SetXTitle("Multiplicity for integrated flow");
168 fHistMultInt ->SetYTitle("Counts");
170 sName = "Control_Flow_MultDiff_";
172 fHistMultDiff = new TH1F(sName.Data(), sName.Data(),iNbinsMult, dMultMin, dMultMax);
173 fHistMultDiff ->SetXTitle("Multiplicity for differential flow");
174 fHistMultDiff ->SetYTitle("Counts");
177 sName = "Control_Flow_PtInt_";
179 fHistPtInt = new TH1F(sName.Data(), sName.Data(),iNbinsPt, dPtMin, dPtMax);
180 fHistPtInt ->SetXTitle("P_{t} (GeV/c) for integrated flow");
181 fHistPtInt ->SetYTitle("Counts");
183 sName = "Control_Flow_PtDiff_";
185 fHistPtDiff = new TH1F(sName.Data(), sName.Data(),iNbinsPt, dPtMin, dPtMax);
186 //binning has to be the same as for fHistProVPt! use to get Nprime!
187 fHistPtDiff ->SetXTitle("P_{t} (GeV/c) for differential flow");
188 fHistPtDiff ->SetYTitle("Counts");
191 sName = "Control_Flow_PhiInt_";
193 fHistPhiInt = new TH1F(sName.Data(), sName.Data(),iNbinsPhi, dPhiMin, dPhiMax);
194 fHistPhiInt ->SetXTitle("#phi for integrated flow");
195 fHistPhiInt ->SetYTitle("Counts");
197 sName = "Control_Flow_PhiDiff_";
199 fHistPhiDiff = new TH1F(sName.Data(), sName.Data(),iNbinsPhi, dPhiMin, dPhiMax);
200 fHistPhiDiff ->SetXTitle("#phi for differential flow");
201 fHistPhiDiff ->SetYTitle("Counts");
204 sName = "Control_Flow_EtaInt_";
206 fHistEtaInt = new TH1F(sName.Data(), sName.Data(),iNbinsEta, dEtaMin, dEtaMax);
207 fHistEtaInt ->SetXTitle("#eta for integrated flow");
208 fHistEtaInt ->SetYTitle("Counts");
210 sName = "Control_Flow_EtaDiff_";
212 fHistEtaDiff = new TH1F(sName.Data(), sName.Data(),iNbinsEta, dEtaMin, dEtaMax);
213 fHistEtaDiff ->SetXTitle("#eta for differential flow");
214 fHistEtaDiff ->SetYTitle("Counts");
217 sName = "Control_FlowPro_MeanPtperBin_";
219 fHistProMeanPtperBin = new TProfile(sName.Data(), sName.Data(),iNbinsPt,dPtMin,dPtMax);
220 fHistProMeanPtperBin ->SetXTitle("P_{t} (GeV/c) ");
221 fHistProMeanPtperBin ->SetYTitle("<P_{t}> (GeV/c) ");
224 sName = "Control_Flow_Q_";
226 fHistQ = new TH1F(sName.Data(), sName.Data(),iNbinsQ, dQMin, dQMax);
227 fHistQ ->SetXTitle("Q_{vector}/Mult");
228 fHistQ ->SetYTitle("Counts");
230 //list of histograms if added here also add in copy constructor
231 fHistList = new TList();
232 fHistList-> Add(fHistMultOrig);
233 fHistList-> Add(fHistMultInt);
234 fHistList-> Add(fHistMultDiff);
235 fHistList-> Add(fHistPtInt);
236 fHistList-> Add(fHistPtDiff);
237 fHistList-> Add(fHistPhiInt);
238 fHistList-> Add(fHistPhiDiff);
239 fHistList-> Add(fHistEtaInt);
240 fHistList-> Add(fHistEtaDiff);
241 fHistList-> Add(fHistProMeanPtperBin);
242 fHistList-> Add(fHistQ);
249 //-----------------------------------------------------------------------
251 AliFlowCommonHist::~AliFlowCommonHist()
254 delete fHistMultOrig;
256 delete fHistMultDiff;
263 delete fHistProMeanPtperBin;
269 //-----------------------------------------------------------------------
271 Bool_t AliFlowCommonHist::FillControlHistograms(AliFlowEventSimple* anEvent)
273 //Fills the control histograms
275 cout<<"##### FillControlHistograms: FlowEvent pointer null"<<endl;
279 Double_t dPt, dPhi, dEta;
282 //fill the histograms
283 Int_t iNumberOfTracks = anEvent->NumberOfTracks();
284 fHistMultOrig->Fill(iNumberOfTracks);
286 AliFlowVector vQ = anEvent->GetQ();
287 //weight by the Multiplicity
290 if (vQ.GetMult()!=0) {
291 dQX = vQ.X()/vQ.GetMult();
292 dQY = vQ.Y()/vQ.GetMult();
295 fHistQ->Fill(vQ.Mod());
300 AliFlowTrackSimple* pTrack = NULL;
302 for (Int_t i=0;i<iNumberOfTracks;i++) {
303 pTrack = anEvent->GetTrack(i);
305 if (pTrack->UseForIntegratedFlow()){
307 fHistPtInt->Fill(dPt);
308 dPhi = pTrack->Phi();
309 if (dPhi<0.) dPhi+=2*TMath::Pi();
310 fHistPhiInt->Fill(dPhi);
311 dEta = pTrack->Eta();
312 fHistEtaInt->Fill(dEta);
315 if (pTrack->UseForDifferentialFlow()){
317 fHistPtDiff->Fill(dPt);
318 dPhi = pTrack->Phi();
319 if (dPhi<0.) dPhi+=2*TMath::Pi();
320 fHistPhiDiff->Fill(dPhi);
321 dEta = pTrack->Eta();
322 fHistEtaDiff->Fill(dEta);
323 fHistProMeanPtperBin->Fill(dPt,dPt);
329 fHistMultInt->Fill(iMultInt);
330 fHistMultDiff->Fill(iMultDiff);
335 //-----------------------------------------------------------------------
337 Double_t AliFlowCommonHist::GetEntriesInPtBinRP(Int_t aBin)
339 //get entries in bin aBin from fHistPtDiff
340 Double_t dEntries = fHistPtInt->GetBinContent(aBin);
346 //-----------------------------------------------------------------------
348 Double_t AliFlowCommonHist::GetEntriesInPtBinPOI(Int_t aBin)
350 //get entries in bin aBin from fHistPtDiff
351 Double_t dEntries = fHistPtDiff->GetBinContent(aBin);
357 //-----------------------------------------------------------------------
359 Double_t AliFlowCommonHist::GetEntriesInEtaBinRP(Int_t aBin)
361 //get entries in bin aBin from fHistPtDiff
362 Double_t dEntries = fHistEtaInt->GetBinContent(aBin);
368 //-----------------------------------------------------------------------
370 Double_t AliFlowCommonHist::GetEntriesInEtaBinPOI(Int_t aBin)
372 //get entries in bin aBin from fHistPtDiff
373 Double_t dEntries = fHistEtaDiff->GetBinContent(aBin);
379 //-----------------------------------------------------------------------
381 Double_t AliFlowCommonHist::GetMeanPt(Int_t aBin)
383 //Get entry from bin aBin from fHistProMeanPtperBin
384 Double_t dMeanPt = fHistProMeanPtperBin->GetBinContent(aBin);
391 //-----------------------------------------------------------------------
392 Double_t AliFlowCommonHist::Merge(TCollection *aList)
395 cout<<"entering merge function"<<endl;
396 if (!aList) return 0;
397 if (aList->IsEmpty()) return 0; //no merging is needed
400 TIter next(aList); // list is supposed to contain only objects of the same type as this
401 AliFlowCommonHist *toMerge;
402 // make a temporary list
403 TList *pTemp = new TList();
404 while ((toMerge=(AliFlowCommonHist*)next())) {
405 pTemp->Add(toMerge->GetHistList());
408 // Now call merge for fHistList providing temp list
409 fHistList->Merge(pTemp);
413 cout<<"Merged"<<endl;
414 return (double)iCount;
418 void AliFlowCommonHist::Print(Option_t *option) const
420 // -*-*-*-*-*Print some global quantities for this histogram collection class *-*-*-*-*-*-*-*
421 // ===============================================
422 // printf( "TH1.Print Name = %s, Entries= %d, Total sum= %g\n",GetName(),Int_t(fEntries),GetSumOfWeights());
423 printf( "Class.Print Name = %s, Histogram list:\n",GetName());
426 fHistList->Print(option);
430 printf( "Empty histogram list \n");
434 //-----------------------------------------------------------------------
435 void AliFlowCommonHist::Browse(TBrowser *b)
439 if (fHistList) b->Add(fHistList,"AliFlowCommonHistList");