AliPhysics  9fe175b (9fe175b)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliFlowLYZHist2.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-1999, 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 $Log$
18 */
19 
20 #include "Riostream.h"
21 #include "AliFlowLYZHist2.h"
22 #include "AliFlowCommonConstants.h"
23 #include "TProfile.h"
24 #include "TProfile2D.h"
25 #include "TString.h"
26 #include "TComplex.h"
27 #include "TList.h"
28 #include "TBrowser.h"
29 
30 class TH1D;
31 
32 // Class to organize the histograms in the second run
33 // in the Lee Yang Zeros Flow analysis.
34 // Also contains methods to get values from the histograms
35 // which are called in AliFlowLeeYandZerosMaker::Finish().
36 // author: N. van der Kolk (kolk@nikhef.nl)
37 
39 
40 //-----------------------------------------------------------------------
41 
43  TNamed(),
44  fHistProReNumer(0),
45  fHistProImNumer(0),
46  fHistProReNumerPt(0),
47  fHistProImNumerPt(0),
48  fHistProReNumer2D(0),
49  fHistProImNumer2D(0),
50  fHistList(NULL)
51 {
52  //default constructor
53 }
54 
55 //-----------------------------------------------------------------------
56 
57  AliFlowLYZHist2::AliFlowLYZHist2(Int_t theta, const char* aSelection, const char *anInput, Bool_t useSum):
58  TNamed(anInput,anInput),
59  fHistProReNumer(0),
60  fHistProImNumer(0),
61  fHistProReNumerPt(0),
62  fHistProImNumerPt(0),
63  fHistProReNumer2D(0),
64  fHistProImNumer2D(0),
65  fHistList(NULL)
66 {
67 
68  //constructor creating histograms
69  TString title, name, addlast;
70 
71  if (useSum) { addlast = "LYZSUM"; }
72  else { addlast = "LYZPROD"; }
73 
74  Int_t iNbinsPt = AliFlowCommonConstants::GetMaster()->GetNbinsPt();
75  Int_t iNbinsEta = AliFlowCommonConstants::GetMaster()->GetNbinsEta();
76 
77  Double_t dPtMin = AliFlowCommonConstants::GetMaster()->GetPtMin();
78  Double_t dPtMax = AliFlowCommonConstants::GetMaster()->GetPtMax();
79  Double_t dEtaMin = AliFlowCommonConstants::GetMaster()->GetEtaMin();
80  Double_t dEtaMax = AliFlowCommonConstants::GetMaster()->GetEtaMax();
81 
82  //fHistProReNumer
83  name = "Second_FlowPro_ReNumer";
84  name +=theta;
85  name +=aSelection;
86  name +=addlast;
87  title = name;
88  fHistProReNumer = new TProfile(name.Data(),title.Data(),iNbinsEta,dEtaMin,dEtaMax);
89  fHistProReNumer->SetXTitle("eta");
90  fHistProReNumer->SetYTitle("v (%)");
91 
92  //fHistProImNumer
93  name = "Second_FlowPro_ImNumer";
94  name +=theta;
95  name +=aSelection;
96  name +=addlast;
97  title = name;
98  fHistProImNumer = new TProfile(name.Data(),title.Data(),iNbinsEta,dEtaMin,dEtaMax);
99  fHistProImNumer->SetXTitle("eta");
100  fHistProImNumer->SetYTitle("v (%)");
101 
102  //fHistProReNumerPt
103  name = "Second_FlowPro_ReNumerPt";
104  name +=theta;
105  name +=aSelection;
106  name +=addlast;
107  title = name;
108  fHistProReNumerPt = new TProfile(name.Data(),title.Data(),iNbinsPt,dPtMin,dPtMax);
109  fHistProReNumerPt->SetXTitle("Pt");
110  fHistProReNumerPt->SetYTitle("v (%)");
111 
112  //fHistProImNumerPt
113  name = "Second_FlowPro_ImNumerPt";
114  name +=theta;
115  name +=aSelection;
116  name +=addlast;
117  title = name;
118  fHistProImNumerPt = new TProfile(name.Data(),title.Data(),iNbinsPt,dPtMin,dPtMax);
119  fHistProImNumerPt->SetXTitle("Pt");
120  fHistProImNumerPt->SetYTitle("v (%)");
121 
122  //fHistProReNumer2D
123  name = "Second_FlowPro_ReNumer2D";
124  name +=theta;
125  name +=aSelection;
126  name +=addlast;
127  title = name;
128  fHistProReNumer2D = new TProfile2D(name.Data(),title.Data(),iNbinsEta,dEtaMin,dEtaMax,iNbinsPt,dPtMin,dPtMax);
129  fHistProReNumer2D->SetXTitle("eta");
130  fHistProReNumer2D->SetYTitle("Pt (GeV/c)");
131 
132  //fHistProImNumer2D
133  name = "Second_FlowPro_ImNumer2D";
134  name +=theta;
135  name +=aSelection;
136  name +=addlast;
137  title = name;
138  fHistProImNumer2D = new TProfile2D(name.Data(),title.Data(),iNbinsEta,dEtaMin,dEtaMax,iNbinsPt,dPtMin,dPtMax);
139  fHistProImNumer2D->SetXTitle("eta");
140  fHistProImNumer2D->SetYTitle("Pt (GeV/c)");
141 
142  //list of histograms
143  fHistList = new TList();
144  fHistList-> Add(fHistProReNumer);
145  fHistList-> Add(fHistProImNumer);
150 
151 }
152 
153 //-----------------------------------------------------------------------
154 
156 {
157  //deletes histograms
158  delete fHistProReNumer;
159  delete fHistProImNumer;
160  delete fHistProReNumerPt;
161  delete fHistProImNumerPt;
162  delete fHistProReNumer2D;
163  delete fHistProImNumer2D;
164  delete fHistList;
165 }
166 
167 //-----------------------------------------------------------------------
168 void AliFlowLYZHist2::Fill(Double_t d1, Double_t d2, TComplex c)
169 {
170  //fill the real and imaginary part of fNumer
171 
172  fHistProReNumer->Fill(d1, c.Re());
173  fHistProImNumer->Fill(d1, c.Im());
174 
175  fHistProReNumerPt->Fill(d2, c.Re());
176  fHistProImNumerPt->Fill(d2, c.Im());
177 
178  fHistProReNumer2D->Fill(d1, d2, c.Re());
179  fHistProImNumer2D->Fill(d1, d2, c.Im());
180 }
181 
182 //-----------------------------------------------------------------------
184 {
185  //get the real and imaginary part of fNumer
186  Double_t dReNumer = fHistProReNumer->GetBinContent(i);
187  Double_t dImNumer = fHistProImNumer->GetBinContent(i);
188  TComplex cNumer(dReNumer,dImNumer);
189  //if (dNumer.Rho()==0) {cerr<<"modulus of dNumer is zero in AliFlowLYZHist2::GetNumer(Int_t i)"<<endl;}
190  return cNumer;
191 }
192 
193 //-----------------------------------------------------------------------
195 {
196  //get the real and imaginary part of fNumer
197  Double_t dReNumer = fHistProReNumerPt->GetBinContent(i);
198  Double_t dImNumer = fHistProImNumerPt->GetBinContent(i);
199  TComplex cNumer(dReNumer,dImNumer);
200  return cNumer;
201 }
202 
203 //-----------------------------------------------------------------------
204  Double_t AliFlowLYZHist2::Merge(TCollection *aList)
205 {
206  //merge fuction
207  if (!aList) return 0;
208  if (aList->IsEmpty()) return 0; //no merging is needed
209 
210  Int_t iCount = 0;
211  TIter next(aList); // list is supposed to contain only objects of the same type as this
212  AliFlowLYZHist2 *toMerge;
213  // make a temporary list
214  TList *pTemp = new TList();
215  while ((toMerge=(AliFlowLYZHist2*)next())) {
216  pTemp->Add(toMerge->GetHistList());
217  iCount++;
218  }
219  // Now call merge for fHistList providing temp list
220  fHistList->Merge(pTemp);
221  // Cleanup
222  delete pTemp;
223 
224  return (double)iCount;
225 
226 }
227 
228 //-----------------------------------------------------------------------
229 void AliFlowLYZHist2::Print(Option_t *option) const
230 {
231  // -*-*-*-*-*Print some global quantities for this histogram collection class *-*-*-*-*-*-*-*
232  // ===============================================
233  // printf( "TH1.Print Name = %s, Entries= %d, Total sum= %g\n",GetName(),Int_t(fEntries),GetSumOfWeights());
234  printf( "Class.Print Name = %s, Histogram list:\n",GetName());
235 
236  if (fHistList) {
237  fHistList->Print(option);
238  }
239  else
240  {
241  printf( "Empty histogram list \n");
242  }
243 }
244 
245 //-----------------------------------------------------------------------
246  void AliFlowLYZHist2::Browse(TBrowser *b)
247 {
248 
249  if (!b) return;
250  if (fHistList) b->Add(fHistList,"AliFlowLYZHist2List");
251 }
TProfile * fHistProReNumerPt
const char * title
Definition: MakeQAPdf.C:26
TProfile * fHistProReNumer
void Browse(TBrowser *b)
TList * GetHistList()
virtual Double_t Merge(TCollection *aList)
TProfile * fHistProImNumerPt
TComplex GetNumerPt(Int_t i)
ClassImp(AliFlowLYZHist2) AliFlowLYZHist2
void Fill(Double_t d1, Double_t d2, TComplex c)
static AliFlowCommonConstants * GetMaster()
TProfile2D * fHistProImNumer2D
TComplex GetNumerEta(Int_t i)
TProfile2D * fHistProReNumer2D
TProfile * fHistProImNumer
void Print(Option_t *option="") const
virtual ~AliFlowLYZHist2()