]>
Commit | Line | Data |
---|---|---|
f456b167 | 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 | ||
f456b167 | 16 | |
448e8856 | 17 | //#include "Riostream.h" //in case one wants to make print statements |
f456b167 | 18 | #include "TProfile.h" |
19 | #include "TString.h" | |
8de6876d | 20 | #include "TComplex.h" |
21 | #include "TList.h" | |
448e8856 | 22 | #include "AliFlowLYZConstants.h" |
23 | #include "AliFlowLYZHist1.h" | |
24 | ||
f456b167 | 25 | class TH1D; |
26 | ||
f456b167 | 27 | // Class to organize the histograms in the first run |
28 | // in the Lee Yang Zeros Flow analysis. | |
29 | // Also contains methods to find the first minimum R0 | |
30 | // of the generating function. | |
31 | // author: N. van der Kolk (kolk@nikhef.nl) | |
32 | ||
f456b167 | 33 | ClassImp(AliFlowLYZHist1) |
34 | ||
8de6876d | 35 | //----------------------------------------------------------------------- |
36 | ||
37 | AliFlowLYZHist1::AliFlowLYZHist1(): | |
38 | TObject(), | |
39 | fHistGtheta(0), | |
40 | fHistProReGtheta(0), | |
41 | fHistProImGtheta(0), | |
42 | fHistList(NULL) | |
43 | { | |
44 | //default constructor | |
45 | } | |
f456b167 | 46 | |
47 | //----------------------------------------------------------------------- | |
48 | ||
448e8856 | 49 | AliFlowLYZHist1::AliFlowLYZHist1(Int_t theta): |
8de6876d | 50 | TObject(), |
448e8856 | 51 | fHistGtheta(0), |
52 | fHistProReGtheta(0), | |
8de6876d | 53 | fHistProImGtheta(0), |
54 | fHistList(NULL) | |
f456b167 | 55 | { |
448e8856 | 56 | |
f456b167 | 57 | //constructor creating histograms |
882ffd6a | 58 | Int_t iNbins = AliFlowLYZConstants::kNbins; |
59 | Double_t dMin = AliFlowLYZConstants::fgMin; | |
60 | Double_t dMax = AliFlowLYZConstants::fgMax; | |
f456b167 | 61 | TString title, name; |
62 | ||
63 | ||
64 | //fHistGtheta | |
448e8856 | 65 | name = "First_Flow_Gtheta"; |
f456b167 | 66 | name +=theta; |
448e8856 | 67 | name +="_LYZ"; |
68 | title = "First_Flow_Gtheta"; | |
f456b167 | 69 | title +=theta; |
448e8856 | 70 | title +="_LYZ"; |
882ffd6a | 71 | fHistGtheta = new TH1D(name.Data(),title.Data(),iNbins,dMin,dMax); |
f456b167 | 72 | fHistGtheta->SetXTitle("r"); |
73 | fHistGtheta->SetYTitle("|G^{#theta}(ir)|^{2}"); | |
74 | ||
75 | //fHistProReGtheta | |
448e8856 | 76 | name = "First_FlowPro_ReGtheta"; |
f456b167 | 77 | name +=theta; |
448e8856 | 78 | name +="_LYZ"; |
79 | title = "First_FlowPro_ReGtheta"; | |
f456b167 | 80 | title +=theta; |
448e8856 | 81 | title +="_LYZ"; |
882ffd6a | 82 | fHistProReGtheta = new TProfile(name.Data(),title.Data(),iNbins,dMin,dMax); |
f456b167 | 83 | fHistProReGtheta->SetXTitle("r"); |
84 | fHistProReGtheta->SetYTitle("Re G^{#theta}(ir)"); | |
85 | ||
86 | //fHistProImGtheta | |
448e8856 | 87 | name = "First_FlowPro_ImGtheta"; |
f456b167 | 88 | name +=theta; |
448e8856 | 89 | name +="_LYZ"; |
90 | title = "First_FlowPro_ImGtheta"; | |
f456b167 | 91 | title +=theta; |
448e8856 | 92 | title +="_LYZ"; |
882ffd6a | 93 | fHistProImGtheta = new TProfile(name.Data(),title.Data(),iNbins,dMin,dMax); |
f456b167 | 94 | fHistProImGtheta->SetXTitle("r"); |
95 | fHistProImGtheta->SetYTitle("Im G^{#theta}(ir)"); | |
8de6876d | 96 | |
97 | //list of histograms | |
98 | fHistList = new TList(); | |
99 | fHistList-> Add(fHistGtheta); | |
100 | fHistList-> Add(fHistProReGtheta); | |
101 | fHistList-> Add(fHistProImGtheta); | |
448e8856 | 102 | |
f456b167 | 103 | } |
104 | ||
f456b167 | 105 | //----------------------------------------------------------------------- |
106 | ||
107 | AliFlowLYZHist1::~AliFlowLYZHist1() | |
108 | { | |
109 | //deletes histograms | |
110 | delete fHistGtheta; | |
111 | delete fHistProReGtheta; | |
112 | delete fHistProImGtheta; | |
8de6876d | 113 | delete fHistList; |
f456b167 | 114 | } |
115 | ||
f456b167 | 116 | //----------------------------------------------------------------------- |
117 | ||
882ffd6a | 118 | void AliFlowLYZHist1::Fill(Double_t f, TComplex c) |
f456b167 | 119 | { |
120 | //fill the histograms | |
121 | ||
882ffd6a | 122 | fHistProReGtheta->Fill(f, c.Re()); |
123 | fHistProImGtheta->Fill(f, c.Im()); | |
f456b167 | 124 | } |
125 | ||
f456b167 | 126 | //----------------------------------------------------------------------- |
127 | ||
128 | TH1D* AliFlowLYZHist1::FillGtheta() | |
129 | { | |
130 | //fill the fHistGtheta histograms | |
882ffd6a | 131 | Int_t iNbins = fHistGtheta->GetNbinsX(); |
132 | for (Int_t bin=1;bin<=iNbins;bin++) | |
f456b167 | 133 | { |
134 | //get bincentre of bins in histogram | |
882ffd6a | 135 | Double_t dBin = fHistGtheta->GetBinCenter(bin); |
136 | Double_t dRe = fHistProReGtheta->GetBinContent(bin); | |
137 | Double_t dIm = fHistProImGtheta->GetBinContent(bin); | |
138 | TComplex cGtheta(dRe,dIm); | |
139 | //fill fHistGtheta with the modulus squared of cGtheta | |
140 | fHistGtheta->Fill(dBin,cGtheta.Rho2()); | |
f456b167 | 141 | } |
142 | ||
143 | return fHistGtheta; | |
144 | } | |
145 | ||
146 | //----------------------------------------------------------------------- | |
147 | ||
448e8856 | 148 | Double_t AliFlowLYZHist1::GetR0() |
f456b167 | 149 | { |
150 | //find the first minimum of the square of the modulus of Gtheta | |
151 | ||
882ffd6a | 152 | Int_t iNbins = fHistGtheta->GetNbinsX(); |
153 | Double_t dR0 = 0.; | |
f456b167 | 154 | |
882ffd6a | 155 | for (Int_t b=2;b<iNbins;b++) |
f456b167 | 156 | { |
882ffd6a | 157 | Double_t dG0 = fHistGtheta->GetBinContent(b); |
158 | Double_t dGnext = fHistGtheta->GetBinContent(b+1); | |
159 | Double_t dGnextnext = fHistGtheta->GetBinContent(b+2); | |
448e8856 | 160 | |
882ffd6a | 161 | if (dGnext > dG0 && dGnextnext > dG0) |
f456b167 | 162 | { |
882ffd6a | 163 | Double_t dGlast = fHistGtheta->GetBinContent(b-1); |
164 | Double_t dXlast = fHistGtheta->GetBinCenter(b-1); | |
165 | Double_t dX0 = fHistGtheta->GetBinCenter(b); | |
166 | Double_t dXnext = fHistGtheta->GetBinCenter(b+1); | |
f456b167 | 167 | |
882ffd6a | 168 | dR0 = dX0 - ((dX0-dXlast)*(dX0-dXlast)*(dG0-dGnext) - (dX0-dXnext)*(dX0-dXnext)*(dG0-dGlast))/ |
169 | (2.*((dX0-dXlast)*(dG0-dGnext) - (dX0-dXnext)*(dG0-dGlast))); //interpolated minimum | |
f456b167 | 170 | |
171 | break; //stop loop if minimum is found | |
172 | } //if | |
173 | ||
174 | }//b | |
175 | ||
176 | ||
882ffd6a | 177 | return dR0; |
f456b167 | 178 | } |
179 | ||
180 | //----------------------------------------------------------------------- | |
448e8856 | 181 | Double_t AliFlowLYZHist1::GetBinCenter(Int_t i) |
f456b167 | 182 | { |
183 | //gets bincenter of histogram | |
882ffd6a | 184 | Double_t dR = fHistGtheta->GetBinCenter(i); |
185 | return dR; | |
f456b167 | 186 | } |
187 | ||
188 | //----------------------------------------------------------------------- | |
189 | ||
190 | Int_t AliFlowLYZHist1::GetNBins() | |
191 | { | |
882ffd6a | 192 | //gets iNbins |
193 | Int_t iNbins = fHistGtheta->GetNbinsX(); | |
194 | return iNbins; | |
f456b167 | 195 | } |
196 | ||
8de6876d | 197 | //----------------------------------------------------------------------- |
198 | Double_t AliFlowLYZHist1::Merge(TCollection *aList) | |
199 | { | |
200 | //merge fuction | |
201 | if (!aList) return 0; | |
202 | if (aList->IsEmpty()) return 0; //no merging is needed | |
203 | ||
204 | Int_t iCount = 0; | |
205 | TIter next(aList); // list is supposed to contain only objects of the same type as this | |
206 | AliFlowLYZHist1 *toMerge; | |
207 | // make a temporary list | |
208 | TList *pTemp = new TList(); | |
209 | while ((toMerge=(AliFlowLYZHist1*)next())) { | |
210 | pTemp->Add(toMerge->GetHistList()); | |
211 | iCount++; | |
212 | } | |
213 | // Now call merge for fHistList providing temp list | |
214 | fHistList->Merge(pTemp); | |
215 | // Cleanup | |
216 | delete pTemp; | |
217 | ||
218 | return (double)iCount; | |
219 | ||
220 | } | |
221 |