]>
Commit | Line | Data |
---|---|---|
e4391c02 | 1 | const Int_t numberOfCentralityBins = 12; |
a9f20288 | 2 | TString centralityArray[numberOfCentralityBins] = {"0-100","10-20","20-30","30-40","40-50","50-60","60-70","70-80","0-100","0-1","1-2","2-3"}; |
e4391c02 | 3 | |
6acdbcb2 | 4 | |
5 | const Int_t gRebin = 1; | |
6 | void drawBalanceFunction2DPsi(const char* filename = "AnalysisResultsPsi.root", | |
eb63b883 | 7 | Int_t gCentrality = 1, |
5de9ad1a | 8 | Int_t gBit = -1, |
9 | const char* gCentralityEstimator = 0x0, | |
52daf7b2 | 10 | Bool_t kShowShuffled = kFALSE, |
11 | Bool_t kShowMixed = kTRUE, | |
eb63b883 | 12 | Double_t psiMin = -0.5, Double_t psiMax = 0.5, |
34616eda | 13 | Double_t vertexZMin = -10., |
14 | Double_t vertexZMax = 10., | |
6acdbcb2 | 15 | Double_t ptTriggerMin = -1., |
16 | Double_t ptTriggerMax = -1., | |
17 | Double_t ptAssociatedMin = -1., | |
f0c5040c | 18 | Double_t ptAssociatedMax = -1., |
832b21d8 | 19 | Bool_t kUseVzBinning = kTRUE, |
20 | Bool_t k2pMethod = kTRUE, | |
d67eae53 | 21 | TString eventClass = "EventPlane") //Can be "EventPlane", "Centrality", "Multiplicity" |
22 | { | |
6acdbcb2 | 23 | //Macro that draws the BF distributions for each centrality bin |
24 | //for reaction plane dependent analysis | |
25 | //Author: Panos.Christakoglou@nikhef.nl | |
26 | //Load the PWG2ebye library | |
27 | gSystem->Load("libANALYSIS.so"); | |
28 | gSystem->Load("libANALYSISalice.so"); | |
29 | gSystem->Load("libEventMixing.so"); | |
30 | gSystem->Load("libCORRFW.so"); | |
31 | gSystem->Load("libPWGTools.so"); | |
32 | gSystem->Load("libPWGCFebye.so"); | |
33 | ||
9a0495b0 | 34 | //gROOT->LoadMacro("~/SetPlotStyle.C"); |
35 | //SetPlotStyle(); | |
8795e993 | 36 | gStyle->SetPalette(1,0); |
15dd45a0 | 37 | |
6acdbcb2 | 38 | //Prepare the objects and return them |
5de9ad1a | 39 | TList *listBF = GetListOfObjects(filename,gCentrality,gBit,gCentralityEstimator,0); |
52daf7b2 | 40 | TList *listBFShuffled = NULL; |
5de9ad1a | 41 | if(kShowShuffled) listBFShuffled = GetListOfObjects(filename,gCentrality,gBit,gCentralityEstimator,1); |
52daf7b2 | 42 | TList *listBFMixed = NULL; |
5de9ad1a | 43 | if(kShowMixed) listBFMixed = GetListOfObjects(filename,gCentrality,gBit,gCentralityEstimator,2); |
6acdbcb2 | 44 | if(!listBF) { |
45 | Printf("The TList object was not created"); | |
46 | return; | |
47 | } | |
48 | else | |
27634c13 | 49 | draw(listBF,listBFShuffled,listBFMixed,gCentrality,gCentralityEstimator, |
34616eda | 50 | psiMin,psiMax,vertexZMin,vertexZMax, |
f0c5040c | 51 | ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax, |
832b21d8 | 52 | kUseVzBinning,k2pMethod,eventClass); |
6acdbcb2 | 53 | } |
54 | ||
55 | //______________________________________________________// | |
5365d1d7 | 56 | TList *GetListOfObjects(const char* filename, |
57 | Int_t gCentrality, | |
58 | Int_t gBit, | |
59 | const char *gCentralityEstimator, | |
5de9ad1a | 60 | Int_t kData = 1) { |
6acdbcb2 | 61 | //Get the TList objects (QA, bf, bf shuffled) |
6acdbcb2 | 62 | TList *listBF = 0x0; |
6acdbcb2 | 63 | |
64 | //Open the file | |
5365d1d7 | 65 | TFile *f = TFile::Open(filename,"UPDATE"); |
6acdbcb2 | 66 | if((!f)||(!f->IsOpen())) { |
67 | Printf("The file %s is not found. Aborting...",filename); | |
68 | return listBF; | |
69 | } | |
70 | //f->ls(); | |
71 | ||
72 | TDirectoryFile *dir = dynamic_cast<TDirectoryFile *>(f->Get("PWGCFEbyE.outputBalanceFunctionPsiAnalysis")); | |
73 | if(!dir) { | |
74 | Printf("The TDirectoryFile is not found. Aborting...",filename); | |
75 | return listBF; | |
76 | } | |
77 | //dir->ls(); | |
78 | ||
79 | TString listBFName; | |
eb63b883 | 80 | if(kData == 0) { |
81 | //cout<<"no shuffling - no mixing"<<endl; | |
82 | listBFName = "listBFPsi_"; | |
83 | } | |
84 | else if(kData == 1) { | |
85 | //cout<<"shuffling - no mixing"<<endl; | |
86 | listBFName = "listBFPsiShuffled_"; | |
87 | } | |
88 | else if(kData == 2) { | |
89 | //cout<<"no shuffling - mixing"<<endl; | |
90 | listBFName = "listBFPsiMixed_"; | |
91 | } | |
92 | listBFName += centralityArray[gCentrality-1]; | |
5de9ad1a | 93 | if(gBit > -1) { |
94 | listBFName += "_Bit"; listBFName += gBit; } | |
95 | if(gCentralityEstimator) { | |
96 | listBFName += "_"; listBFName += gCentralityEstimator;} | |
6acdbcb2 | 97 | |
5365d1d7 | 98 | // histograms were already retrieved (in first iteration) |
99 | if(dir->Get(Form("%s_histograms",listBFName.Data()))){ | |
100 | //cout<<"second iteration"<<endl; | |
101 | listBF = dynamic_cast<TList *>(dir->Get(Form("%s_histograms",listBFName.Data()))); | |
6acdbcb2 | 102 | } |
6acdbcb2 | 103 | |
5365d1d7 | 104 | // histograms were not yet retrieved (this is the first iteration) |
105 | else{ | |
106 | ||
107 | listBF = dynamic_cast<TList *>(dir->Get(listBFName.Data())); | |
108 | cout<<"======================================================="<<endl; | |
5365d1d7 | 109 | cout<<"List name: "<<listBF->GetName()<<endl; |
110 | //listBF->ls(); | |
6acdbcb2 | 111 | |
5365d1d7 | 112 | //Get the histograms |
113 | TString histoName; | |
114 | if(kData == 0) | |
27634c13 | 115 | histoName = "fHistP"; |
5365d1d7 | 116 | else if(kData == 1) |
27634c13 | 117 | histoName = "fHistP_shuffle"; |
5365d1d7 | 118 | else if(kData == 2) |
27634c13 | 119 | histoName = "fHistP"; |
120 | if(gCentralityEstimator) histoName += gCentralityEstimator; | |
5365d1d7 | 121 | AliTHn *fHistP = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data())); |
122 | if(!fHistP) { | |
123 | Printf("fHistP %s not found!!!",histoName.Data()); | |
124 | break; | |
125 | } | |
126 | fHistP->FillParent(); fHistP->DeleteContainers(); | |
127 | ||
128 | if(kData == 0) | |
27634c13 | 129 | histoName = "fHistN"; |
5365d1d7 | 130 | if(kData == 1) |
27634c13 | 131 | histoName = "fHistN_shuffle"; |
5365d1d7 | 132 | if(kData == 2) |
27634c13 | 133 | histoName = "fHistN"; |
134 | if(gCentralityEstimator) histoName += gCentralityEstimator; | |
5365d1d7 | 135 | AliTHn *fHistN = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data())); |
136 | if(!fHistN) { | |
137 | Printf("fHistN %s not found!!!",histoName.Data()); | |
138 | break; | |
139 | } | |
140 | fHistN->FillParent(); fHistN->DeleteContainers(); | |
141 | ||
142 | if(kData == 0) | |
27634c13 | 143 | histoName = "fHistPN"; |
5365d1d7 | 144 | if(kData == 1) |
27634c13 | 145 | histoName = "fHistPN_shuffle"; |
5365d1d7 | 146 | if(kData == 2) |
27634c13 | 147 | histoName = "fHistPN"; |
148 | if(gCentralityEstimator) histoName += gCentralityEstimator; | |
5365d1d7 | 149 | AliTHn *fHistPN = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data())); |
150 | if(!fHistPN) { | |
151 | Printf("fHistPN %s not found!!!",histoName.Data()); | |
152 | break; | |
153 | } | |
154 | fHistPN->FillParent(); fHistPN->DeleteContainers(); | |
155 | ||
156 | if(kData == 0) | |
27634c13 | 157 | histoName = "fHistNP"; |
5365d1d7 | 158 | if(kData == 1) |
27634c13 | 159 | histoName = "fHistNP_shuffle"; |
5365d1d7 | 160 | if(kData == 2) |
27634c13 | 161 | histoName = "fHistNP"; |
162 | if(gCentralityEstimator) histoName += gCentralityEstimator; | |
5365d1d7 | 163 | AliTHn *fHistNP = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data())); |
164 | if(!fHistNP) { | |
165 | Printf("fHistNP %s not found!!!",histoName.Data()); | |
166 | break; | |
167 | } | |
168 | fHistNP->FillParent(); fHistNP->DeleteContainers(); | |
169 | ||
170 | if(kData == 0) | |
27634c13 | 171 | histoName = "fHistPP"; |
5365d1d7 | 172 | if(kData == 1) |
27634c13 | 173 | histoName = "fHistPP_shuffle"; |
5365d1d7 | 174 | if(kData == 2) |
27634c13 | 175 | histoName = "fHistPP"; |
176 | if(gCentralityEstimator) histoName += gCentralityEstimator; | |
5365d1d7 | 177 | AliTHn *fHistPP = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data())); |
178 | if(!fHistPP) { | |
179 | Printf("fHistPP %s not found!!!",histoName.Data()); | |
180 | break; | |
181 | } | |
182 | fHistPP->FillParent(); fHistPP->DeleteContainers(); | |
183 | ||
184 | if(kData == 0) | |
27634c13 | 185 | histoName = "fHistNN"; |
5365d1d7 | 186 | if(kData == 1) |
27634c13 | 187 | histoName = "fHistNN_shuffle"; |
5365d1d7 | 188 | if(kData == 2) |
27634c13 | 189 | histoName = "fHistNN"; |
190 | if(gCentralityEstimator) histoName += gCentralityEstimator; | |
5365d1d7 | 191 | AliTHn *fHistNN = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data())); |
192 | if(!fHistNN) { | |
193 | Printf("fHistNN %s not found!!!",histoName.Data()); | |
194 | break; | |
195 | } | |
196 | fHistNN->FillParent(); fHistNN->DeleteContainers(); | |
197 | ||
198 | dir->cd(); | |
199 | listBF->Write(Form("%s_histograms",listBFName.Data()), TObject::kSingleKey); | |
200 | ||
201 | }// first iteration | |
6acdbcb2 | 202 | |
5365d1d7 | 203 | f->Close(); |
6acdbcb2 | 204 | |
205 | return listBF; | |
206 | } | |
207 | ||
208 | //______________________________________________________// | |
eb63b883 | 209 | void draw(TList *listBF, TList *listBFShuffled, TList *listBFMixed, |
27634c13 | 210 | Int_t gCentrality, const char* gCentralityEstimator, |
211 | Double_t psiMin, Double_t psiMax, | |
34616eda | 212 | Double_t vertexZMin, |
213 | Double_t vertexZMax, | |
6acdbcb2 | 214 | Double_t ptTriggerMin, Double_t ptTriggerMax, |
f0c5040c | 215 | Double_t ptAssociatedMin, Double_t ptAssociatedMax, |
832b21d8 | 216 | Bool_t kUseVzBinning=kFALSE, |
d67eae53 | 217 | Bool_t k2pMethod = kFALSE, TString eventClass) { |
6acdbcb2 | 218 | //balance function |
219 | AliTHn *hP = NULL; | |
220 | AliTHn *hN = NULL; | |
221 | AliTHn *hPN = NULL; | |
222 | AliTHn *hNP = NULL; | |
223 | AliTHn *hPP = NULL; | |
224 | AliTHn *hNN = NULL; | |
225 | //listBF->ls(); | |
226 | //Printf("================="); | |
27634c13 | 227 | TString histoName = "fHistP"; |
228 | if(gCentralityEstimator) histoName += gCentralityEstimator; | |
229 | hP = (AliTHn*) listBF->FindObject(histoName.Data()); | |
52daf7b2 | 230 | hP->SetName("gHistP"); |
27634c13 | 231 | histoName = "fHistN"; |
232 | if(gCentralityEstimator) histoName += gCentralityEstimator; | |
233 | hN = (AliTHn*) listBF->FindObject(histoName.Data()); | |
52daf7b2 | 234 | hN->SetName("gHistN"); |
27634c13 | 235 | histoName = "fHistPN"; |
236 | if(gCentralityEstimator) histoName += gCentralityEstimator; | |
237 | hPN = (AliTHn*) listBF->FindObject(histoName.Data()); | |
52daf7b2 | 238 | hPN->SetName("gHistPN"); |
27634c13 | 239 | histoName = "fHistNP"; |
240 | if(gCentralityEstimator) histoName += gCentralityEstimator; | |
241 | hNP = (AliTHn*) listBF->FindObject(histoName.Data()); | |
52daf7b2 | 242 | hNP->SetName("gHistNP"); |
27634c13 | 243 | histoName = "fHistPP"; |
244 | if(gCentralityEstimator) histoName += gCentralityEstimator; | |
245 | hPP = (AliTHn*) listBF->FindObject(histoName.Data()); | |
52daf7b2 | 246 | hPP->SetName("gHistPP"); |
27634c13 | 247 | histoName = "fHistNN"; |
248 | if(gCentralityEstimator) histoName += gCentralityEstimator; | |
249 | hNN = (AliTHn*) listBF->FindObject(histoName.Data()); | |
52daf7b2 | 250 | hNN->SetName("gHistNN"); |
6acdbcb2 | 251 | |
252 | AliBalancePsi *b = new AliBalancePsi(); | |
d67eae53 | 253 | b->SetEventClass(eventClass); |
6acdbcb2 | 254 | b->SetHistNp(hP); |
255 | b->SetHistNn(hN); | |
256 | b->SetHistNpn(hPN); | |
257 | b->SetHistNnp(hNP); | |
258 | b->SetHistNpp(hPP); | |
259 | b->SetHistNnn(hNN); | |
832b21d8 | 260 | if(kUseVzBinning) b->SetVertexZBinning(kTRUE); |
261 | ||
6acdbcb2 | 262 | |
263 | //balance function shuffling | |
264 | AliTHn *hPShuffled = NULL; | |
265 | AliTHn *hNShuffled = NULL; | |
266 | AliTHn *hPNShuffled = NULL; | |
267 | AliTHn *hNPShuffled = NULL; | |
268 | AliTHn *hPPShuffled = NULL; | |
269 | AliTHn *hNNShuffled = NULL; | |
52daf7b2 | 270 | if(listBFShuffled) { |
271 | //listBFShuffled->ls(); | |
abf0c200 | 272 | histoName = "fHistP_shuffle"; |
27634c13 | 273 | if(gCentralityEstimator) histoName += gCentralityEstimator; |
274 | hPShuffled = (AliTHn*) listBFShuffled->FindObject(histoName.Data()); | |
52daf7b2 | 275 | hPShuffled->SetName("gHistPShuffled"); |
abf0c200 | 276 | histoName = "fHistN_shuffle"; |
27634c13 | 277 | if(gCentralityEstimator) histoName += gCentralityEstimator; |
278 | hNShuffled = (AliTHn*) listBFShuffled->FindObject(histoName.Data()); | |
52daf7b2 | 279 | hNShuffled->SetName("gHistNShuffled"); |
abf0c200 | 280 | histoName = "fHistPN_shuffle"; |
27634c13 | 281 | if(gCentralityEstimator) histoName += gCentralityEstimator; |
282 | hPNShuffled = (AliTHn*) listBFShuffled->FindObject(histoName.Data()); | |
52daf7b2 | 283 | hPNShuffled->SetName("gHistPNShuffled"); |
abf0c200 | 284 | histoName = "fHistNP_shuffle"; |
27634c13 | 285 | if(gCentralityEstimator) histoName += gCentralityEstimator; |
286 | hNPShuffled = (AliTHn*) listBFShuffled->FindObject(histoName.Data()); | |
52daf7b2 | 287 | hNPShuffled->SetName("gHistNPShuffled"); |
abf0c200 | 288 | histoName = "fHistPP_shuffle"; |
27634c13 | 289 | if(gCentralityEstimator) histoName += gCentralityEstimator; |
290 | hPPShuffled = (AliTHn*) listBFShuffled->FindObject(histoName.Data()); | |
52daf7b2 | 291 | hPPShuffled->SetName("gHistPPShuffled"); |
abf0c200 | 292 | histoName = "fHistNN_shuffle"; |
27634c13 | 293 | if(gCentralityEstimator) histoName += gCentralityEstimator; |
294 | hNNShuffled = (AliTHn*) listBFShuffled->FindObject(histoName.Data()); | |
52daf7b2 | 295 | hNNShuffled->SetName("gHistNNShuffled"); |
296 | ||
297 | AliBalancePsi *bShuffled = new AliBalancePsi(); | |
d67eae53 | 298 | bShuffled->SetEventClass(eventClass); |
52daf7b2 | 299 | bShuffled->SetHistNp(hPShuffled); |
300 | bShuffled->SetHistNn(hNShuffled); | |
301 | bShuffled->SetHistNpn(hPNShuffled); | |
302 | bShuffled->SetHistNnp(hNPShuffled); | |
303 | bShuffled->SetHistNpp(hPPShuffled); | |
304 | bShuffled->SetHistNnn(hNNShuffled); | |
832b21d8 | 305 | if(kUseVzBinning) bShuffled->SetVertexZBinning(kTRUE); |
306 | ||
52daf7b2 | 307 | } |
6acdbcb2 | 308 | |
eb63b883 | 309 | //balance function mixing |
310 | AliTHn *hPMixed = NULL; | |
311 | AliTHn *hNMixed = NULL; | |
312 | AliTHn *hPNMixed = NULL; | |
313 | AliTHn *hNPMixed = NULL; | |
314 | AliTHn *hPPMixed = NULL; | |
315 | AliTHn *hNNMixed = NULL; | |
eb63b883 | 316 | |
52daf7b2 | 317 | if(listBFMixed) { |
318 | //listBFMixed->ls(); | |
27634c13 | 319 | histoName = "fHistP"; |
320 | if(gCentralityEstimator) histoName += gCentralityEstimator; | |
321 | hPMixed = (AliTHn*) listBFMixed->FindObject(histoName.Data()); | |
52daf7b2 | 322 | hPMixed->SetName("gHistPMixed"); |
27634c13 | 323 | histoName = "fHistN"; |
324 | if(gCentralityEstimator) histoName += gCentralityEstimator; | |
325 | hNMixed = (AliTHn*) listBFMixed->FindObject(histoName.Data()); | |
52daf7b2 | 326 | hNMixed->SetName("gHistNMixed"); |
27634c13 | 327 | histoName = "fHistPN"; |
328 | if(gCentralityEstimator) histoName += gCentralityEstimator; | |
329 | hPNMixed = (AliTHn*) listBFMixed->FindObject(histoName.Data()); | |
52daf7b2 | 330 | hPNMixed->SetName("gHistPNMixed"); |
27634c13 | 331 | histoName = "fHistNP"; |
332 | if(gCentralityEstimator) histoName += gCentralityEstimator; | |
333 | hNPMixed = (AliTHn*) listBFMixed->FindObject(histoName.Data()); | |
334 | histoName = "fHistNP"; | |
335 | if(gCentralityEstimator) histoName += gCentralityEstimator; | |
52daf7b2 | 336 | hNPMixed->SetName("gHistNPMixed"); |
27634c13 | 337 | histoName = "fHistPP"; |
338 | if(gCentralityEstimator) histoName += gCentralityEstimator; | |
339 | hPPMixed = (AliTHn*) listBFMixed->FindObject(histoName.Data()); | |
52daf7b2 | 340 | hPPMixed->SetName("gHistPPMixed"); |
27634c13 | 341 | histoName = "fHistNN"; |
342 | if(gCentralityEstimator) histoName += gCentralityEstimator; | |
343 | hNNMixed = (AliTHn*) listBFMixed->FindObject(histoName.Data()); | |
52daf7b2 | 344 | hNNMixed->SetName("gHistNNMixed"); |
345 | ||
346 | AliBalancePsi *bMixed = new AliBalancePsi(); | |
d67eae53 | 347 | bMixed->SetEventClass(eventClass); |
52daf7b2 | 348 | bMixed->SetHistNp(hPMixed); |
349 | bMixed->SetHistNn(hNMixed); | |
350 | bMixed->SetHistNpn(hPNMixed); | |
351 | bMixed->SetHistNnp(hNPMixed); | |
352 | bMixed->SetHistNpp(hPPMixed); | |
353 | bMixed->SetHistNnn(hNNMixed); | |
832b21d8 | 354 | if(kUseVzBinning) bMixed->SetVertexZBinning(kTRUE); |
355 | ||
52daf7b2 | 356 | } |
eb63b883 | 357 | |
6acdbcb2 | 358 | TH2D *gHistBalanceFunction; |
eb63b883 | 359 | TH2D *gHistBalanceFunctionSubtracted; |
6acdbcb2 | 360 | TH2D *gHistBalanceFunctionShuffled; |
eb63b883 | 361 | TH2D *gHistBalanceFunctionMixed; |
6acdbcb2 | 362 | TString histoTitle, pngName; |
6acdbcb2 | 363 | |
d67eae53 | 364 | if(eventClass == "Centrality"){ |
365 | histoTitle = "Centrality: "; | |
366 | histoTitle += psiMin; | |
367 | histoTitle += " - "; | |
368 | histoTitle += psiMax; | |
369 | histoTitle += " % "; | |
52daf7b2 | 370 | histoTitle += " (0^{o} < #varphi - #Psi_{2} < 180^{o})"; |
d67eae53 | 371 | } |
372 | else if(eventClass == "Multiplicity"){ | |
373 | histoTitle = "Multiplicity: "; | |
374 | histoTitle += psiMin; | |
375 | histoTitle += " - "; | |
376 | histoTitle += psiMax; | |
377 | histoTitle += " tracks"; | |
378 | histoTitle += " (0^{o} < #varphi - #Psi_{2} < 180^{o})"; | |
379 | } | |
380 | else{ // "EventPlane" (default) | |
381 | histoTitle = "Centrality: "; | |
382 | histoTitle += centralityArray[gCentrality-1]; | |
383 | histoTitle += "%"; | |
384 | if((psiMin == -0.5)&&(psiMax == 0.5)) | |
385 | histoTitle += " (-7.5^{o} < #varphi - #Psi_{2} < 7.5^{o})"; | |
386 | else if((psiMin == 0.5)&&(psiMax == 1.5)) | |
387 | histoTitle += " (37.5^{o} < #varphi - #Psi_{2} < 52.5^{o})"; | |
388 | else if((psiMin == 1.5)&&(psiMax == 2.5)) | |
389 | histoTitle += " (82.5^{o} < #varphi - #Psi_{2} < 97.5^{o})"; | |
390 | else | |
391 | histoTitle += " (0^{o} < #varphi - #Psi_{2} < 180^{o})"; | |
392 | } | |
393 | ||
f0c5040c | 394 | if(k2pMethod) |
395 | if(bMixed) | |
34616eda | 396 | gHistBalanceFunction = b->GetBalanceFunctionDeltaEtaDeltaPhi2pMethod(psiMin,psiMax,vertexZMin,vertexZMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax,bMixed); |
f0c5040c | 397 | else{ |
398 | cerr<<"NO MIXED BF BUT REQUESTED CORRECTING WITH IT! --> FAIL"<<endl; | |
399 | return; | |
400 | } | |
401 | else | |
34616eda | 402 | gHistBalanceFunction = b->GetBalanceFunctionDeltaEtaDeltaPhi(psiMin,psiMax,vertexZMin,vertexZMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax); |
6acdbcb2 | 403 | gHistBalanceFunction->SetTitle(histoTitle.Data()); |
404 | gHistBalanceFunction->GetYaxis()->SetTitleOffset(1.3); | |
eb63b883 | 405 | gHistBalanceFunction->SetName("gHistBalanceFunction"); |
406 | ||
52daf7b2 | 407 | if(listBFShuffled) { |
f0c5040c | 408 | |
409 | if(k2pMethod) | |
410 | if(bMixed) | |
34616eda | 411 | gHistBalanceFunctionShuffled = bShuffled->GetBalanceFunctionDeltaEtaDeltaPhi2pMethod(psiMin,psiMax,vertexZMin,vertexZMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax,bMixed); |
f0c5040c | 412 | else{ |
413 | cerr<<"NO MIXED BF BUT REQUESTED CORRECTING WITH IT! --> FAIL"<<endl; | |
414 | return; | |
415 | } | |
416 | else | |
34616eda | 417 | gHistBalanceFunctionShuffled = bShuffled->GetBalanceFunctionDeltaEtaDeltaPhi(psiMin,psiMax,vertexZMin,vertexZMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax); |
52daf7b2 | 418 | gHistBalanceFunctionShuffled->SetTitle(histoTitle.Data()); |
419 | gHistBalanceFunctionShuffled->GetYaxis()->SetTitleOffset(1.3); | |
420 | gHistBalanceFunctionShuffled->SetName("gHistBalanceFunctionShuffled"); | |
421 | } | |
eb63b883 | 422 | |
52daf7b2 | 423 | if(listBFMixed) { |
f0c5040c | 424 | if(k2pMethod) |
425 | if(bMixed) | |
34616eda | 426 | gHistBalanceFunctionMixed = bMixed->GetBalanceFunctionDeltaEtaDeltaPhi2pMethod(psiMin,psiMax,vertexZMin,vertexZMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax,bMixed); |
f0c5040c | 427 | else{ |
428 | cerr<<"NO MIXED BF BUT REQUESTED CORRECTING WITH IT! --> FAIL"<<endl; | |
429 | return; | |
430 | } | |
431 | else | |
34616eda | 432 | gHistBalanceFunctionMixed = bMixed->GetBalanceFunctionDeltaEtaDeltaPhi(psiMin,psiMax,vertexZMin,vertexZMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax); |
52daf7b2 | 433 | gHistBalanceFunctionMixed->SetTitle(histoTitle.Data()); |
434 | gHistBalanceFunctionMixed->GetYaxis()->SetTitleOffset(1.3); | |
435 | gHistBalanceFunctionMixed->SetName("gHistBalanceFunctionMixed"); | |
436 | ||
437 | gHistBalanceFunctionSubtracted = dynamic_cast<TH2D *>(gHistBalanceFunction->Clone()); | |
438 | gHistBalanceFunctionSubtracted->Add(gHistBalanceFunctionMixed,-1); | |
439 | gHistBalanceFunctionSubtracted->SetTitle(histoTitle.Data()); | |
440 | gHistBalanceFunctionSubtracted->GetYaxis()->SetTitleOffset(1.3); | |
441 | gHistBalanceFunctionSubtracted->SetName("gHistBalanceFunctionSubtracted"); | |
442 | } | |
eb63b883 | 443 | |
444 | //Draw the results | |
445 | TCanvas *c1 = new TCanvas("c1","",0,0,600,500); | |
6acdbcb2 | 446 | c1->SetFillColor(10); |
447 | c1->SetHighLightColor(10); | |
448 | c1->SetLeftMargin(0.15); | |
eb63b883 | 449 | gHistBalanceFunction->DrawCopy("lego2"); |
5de9ad1a | 450 | gPad->SetTheta(30); // default is 30 |
451 | gPad->SetPhi(-60); // default is 30 | |
452 | gPad->Update(); | |
eb63b883 | 453 | TCanvas *c1a = new TCanvas("c1a","",600,0,600,500); |
454 | c1a->SetFillColor(10); | |
455 | c1a->SetHighLightColor(10); | |
456 | c1a->SetLeftMargin(0.15); | |
457 | gHistBalanceFunction->DrawCopy("colz"); | |
458 | ||
52daf7b2 | 459 | if(listBFShuffled) { |
460 | TCanvas *c2 = new TCanvas("c2","",100,100,600,500); | |
461 | c2->SetFillColor(10); | |
462 | c2->SetHighLightColor(10); | |
463 | c2->SetLeftMargin(0.15); | |
464 | gHistBalanceFunctionShuffled->DrawCopy("lego2"); | |
5de9ad1a | 465 | gPad->SetTheta(30); // default is 30 |
466 | gPad->SetPhi(-60); // default is 30 | |
467 | gPad->Update(); | |
52daf7b2 | 468 | TCanvas *c2a = new TCanvas("c2a","",700,100,600,500); |
469 | c2a->SetFillColor(10); | |
470 | c2a->SetHighLightColor(10); | |
471 | c2a->SetLeftMargin(0.15); | |
472 | gHistBalanceFunctionShuffled->DrawCopy("colz"); | |
473 | } | |
eb63b883 | 474 | |
52daf7b2 | 475 | if(listBFMixed) { |
476 | TCanvas *c3 = new TCanvas("c3","",200,200,600,500); | |
477 | c3->SetFillColor(10); | |
478 | c3->SetHighLightColor(10); | |
479 | c3->SetLeftMargin(0.15); | |
480 | gHistBalanceFunctionMixed->DrawCopy("lego2"); | |
5de9ad1a | 481 | gPad->SetTheta(30); // default is 30 |
482 | gPad->SetPhi(-60); // default is 30 | |
483 | gPad->Update(); | |
52daf7b2 | 484 | TCanvas *c3a = new TCanvas("c3a","",800,200,600,500); |
485 | c3a->SetFillColor(10); | |
486 | c3a->SetHighLightColor(10); | |
487 | c3a->SetLeftMargin(0.15); | |
488 | gHistBalanceFunctionMixed->DrawCopy("colz"); | |
eb63b883 | 489 | |
52daf7b2 | 490 | TCanvas *c4 = new TCanvas("c4","",300,300,600,500); |
491 | c4->SetFillColor(10); | |
492 | c4->SetHighLightColor(10); | |
493 | c4->SetLeftMargin(0.15); | |
494 | gHistBalanceFunctionSubtracted->DrawCopy("lego2"); | |
5de9ad1a | 495 | gPad->SetTheta(30); // default is 30 |
496 | gPad->SetPhi(-60); // default is 30 | |
497 | gPad->Update(); | |
52daf7b2 | 498 | TCanvas *c4a = new TCanvas("c4a","",900,300,600,500); |
499 | c4a->SetFillColor(10); | |
500 | c4a->SetHighLightColor(10); | |
501 | c4a->SetLeftMargin(0.15); | |
502 | gHistBalanceFunctionSubtracted->DrawCopy("colz"); | |
742af4bd | 503 | |
bd36d661 | 504 | fitbalanceFunction(gCentrality, psiMin , psiMax, |
742af4bd | 505 | ptTriggerMin, ptTriggerMax, |
506 | ptAssociatedMin, ptAssociatedMax, | |
d67eae53 | 507 | gHistBalanceFunctionSubtracted,k2pMethod, eventClass); |
52daf7b2 | 508 | } |
eb63b883 | 509 | |
d67eae53 | 510 | TString newFileName = "balanceFunction2D."; |
511 | if(eventClass == "Centrality"){ | |
512 | newFileName += Form("Centrality%.1fTo%.1f",psiMin,psiMax); | |
513 | newFileName += ".PsiAll.PttFrom"; | |
514 | } | |
515 | else if(eventClass == "Multiplicity"){ | |
516 | newFileName += Form("Multiplicity%.0fTo%.0f",psiMin,psiMax); | |
517 | newFileName += ".PsiAll.PttFrom"; | |
518 | } | |
519 | else{ // "EventPlane" (default) | |
520 | newFileName += "Centrality"; | |
521 | newFileName += gCentrality; newFileName += ".Psi"; | |
522 | if((psiMin == -0.5)&&(psiMax == 0.5)) newFileName += "InPlane.Ptt"; | |
523 | else if((psiMin == 0.5)&&(psiMax == 1.5)) newFileName += "Intermediate.Ptt"; | |
524 | else if((psiMin == 1.5)&&(psiMax == 2.5)) newFileName += "OutOfPlane.Ptt"; | |
525 | else if((psiMin == 2.5)&&(psiMax == 3.5)) newFileName += "Rest.PttFrom"; | |
526 | else newFileName += "All.PttFrom"; | |
527 | } | |
648f1a5a | 528 | newFileName += Form("%.1f",ptTriggerMin); newFileName += "To"; |
529 | newFileName += Form("%.1f",ptTriggerMax); newFileName += "PtaFrom"; | |
530 | newFileName += Form("%.1f",ptAssociatedMin); newFileName += "To"; | |
531 | newFileName += Form("%.1f",ptAssociatedMax); | |
f0c5040c | 532 | if(k2pMethod) newFileName += "_2pMethod"; |
a9f20288 | 533 | |
a2da950e | 534 | newFileName += "_"; |
535 | newFileName += Form("%.1f",psiMin); | |
536 | newFileName += "-"; | |
537 | newFileName += Form("%.1f",psiMax); | |
538 | newFileName += ".root"; | |
eb63b883 | 539 | |
540 | TFile *fOutput = new TFile(newFileName.Data(),"recreate"); | |
541 | fOutput->cd(); | |
52daf7b2 | 542 | /*hP->Write(); hN->Write(); |
543 | hPN->Write(); hNP->Write(); | |
544 | hPP->Write(); hNN->Write(); | |
545 | hPShuffled->Write(); hNShuffled->Write(); | |
546 | hPNShuffled->Write(); hNPShuffled->Write(); | |
547 | hPPShuffled->Write(); hNNShuffled->Write(); | |
548 | hPMixed->Write(); hNMixed->Write(); | |
549 | hPNMixed->Write(); hNPMixed->Write(); | |
550 | hPPMixed->Write(); hNNMixed->Write();*/ | |
eb63b883 | 551 | gHistBalanceFunction->Write(); |
52daf7b2 | 552 | if(listBFShuffled) gHistBalanceFunctionShuffled->Write(); |
553 | if(listBFMixed) { | |
554 | gHistBalanceFunctionMixed->Write(); | |
555 | gHistBalanceFunctionSubtracted->Write(); | |
556 | } | |
eb63b883 | 557 | fOutput->Close(); |
6acdbcb2 | 558 | } |
559 | ||
742af4bd | 560 | //____________________________________________________________// |
561 | void fitbalanceFunction(Int_t gCentrality = 1, | |
562 | Double_t psiMin = -0.5, Double_t psiMax = 3.5, | |
563 | Double_t ptTriggerMin = -1., | |
564 | Double_t ptTriggerMax = -1., | |
565 | Double_t ptAssociatedMin = -1., | |
566 | Double_t ptAssociatedMax = -1., | |
d67eae53 | 567 | TH2D *gHist, |
568 | Bool_t k2pMethod = kFALSE, | |
569 | TString eventClass="EventPlane") { | |
48f521eb | 570 | //balancing charges: [1]*TMath::Exp(-0.5*TMath::Power(((x - [3])/[2]),2)-0.5*TMath::Power(((y - [5])/[4]),2)) |
571 | //short range correlations: [6]*TMath::Exp(-0.5*TMath::Power(((x - [8])/[7]),2)-0.5*TMath::Power(((y - [10])/[9]),2)) | |
742af4bd | 572 | cout<<"FITTING FUNCTION"<<endl; |
573 | ||
48f521eb | 574 | TF2 *gFitFunction = new TF2("gFitFunction","[0] + [1]*TMath::Exp(-0.5*TMath::Power(((x - [3])/[2]),2)-0.5*TMath::Power(((y - [5])/[4]),2)) + [6]*TMath::Exp(-0.5*TMath::Power(((x - [8])/[7]),2)-0.5*TMath::Power(((y - [10])/[9]),2))",-1.2,1.2,-TMath::Pi()/2.,3.*TMath::Pi()/2.); |
742af4bd | 575 | gFitFunction->SetName("gFitFunction"); |
576 | ||
742af4bd | 577 | //Normalization |
578 | gFitFunction->SetParName(0,"N1"); | |
48f521eb | 579 | gFitFunction->SetParameter(0,1.0); |
580 | ||
581 | //2D balance function | |
582 | gFitFunction->SetParName(1,"N_{BF}"); | |
583 | gFitFunction->SetParameter(1,1.0); | |
584 | gFitFunction->SetParLimits(1, 0., 100.); | |
585 | gFitFunction->SetParName(2,"Sigma_{BF}(delta eta)"); | |
586 | gFitFunction->SetParameter(2,0.6); | |
587 | gFitFunction->SetParLimits(2, 0., 1.); | |
588 | gFitFunction->SetParName(3,"Mean_{BF}(delta eta)"); | |
589 | gFitFunction->SetParameter(3,0.0); | |
590 | gFitFunction->SetParLimits(3, -0.2, 0.2); | |
591 | gFitFunction->SetParName(4,"Sigma_{BF}(delta phi)"); | |
592 | gFitFunction->SetParameter(4,0.6); | |
593 | gFitFunction->SetParLimits(4, 0., 1.); | |
594 | gFitFunction->SetParName(5,"Mean_{BF}(delta phi)"); | |
595 | gFitFunction->SetParameter(5,0.0); | |
596 | gFitFunction->SetParLimits(5, -0.2, 0.2); | |
597 | ||
598 | //short range structure | |
599 | gFitFunction->SetParName(6,"N_{SR}"); | |
600 | gFitFunction->SetParameter(6,5.0); | |
601 | gFitFunction->SetParLimits(6, 0., 100.); | |
602 | gFitFunction->SetParName(7,"Sigma_{SR}(delta eta)"); | |
603 | gFitFunction->SetParameter(7,0.01); | |
604 | gFitFunction->SetParLimits(7, 0.0, 0.1); | |
605 | gFitFunction->SetParName(8,"Mean_{SR}(delta eta)"); | |
606 | gFitFunction->SetParameter(8,0.0); | |
607 | gFitFunction->SetParLimits(8, -0.01, 0.01); | |
608 | gFitFunction->SetParName(9,"Sigma_{SR}(delta phi)"); | |
609 | gFitFunction->SetParameter(9,0.01); | |
610 | gFitFunction->SetParLimits(9, 0.0, 0.1); | |
611 | gFitFunction->SetParName(10,"Mean_{SR}(delta phi)"); | |
612 | gFitFunction->SetParameter(10,0.0); | |
613 | gFitFunction->SetParLimits(10, -0.01, 0.01); | |
742af4bd | 614 | |
615 | ||
616 | //Cloning the histogram | |
742af4bd | 617 | TH2D *gHistResidual = dynamic_cast<TH2D *>(gHist->Clone()); |
618 | gHistResidual->SetName("gHistResidual"); | |
619 | gHistResidual->Sumw2(); | |
620 | ||
48f521eb | 621 | //Fitting the 2D bf |
622 | for(Int_t iAttempt = 0; iAttempt < 10; iAttempt++) { | |
623 | gHist->Fit("gFitFunction","nm"); | |
624 | for(Int_t iParam = 0; iParam < 11; iParam++) | |
625 | gFitFunction->SetParameter(iParam,gFitFunction->GetParameter(iParam)); | |
742af4bd | 626 | } |
48f521eb | 627 | cout<<"======================================================"<<endl; |
628 | cout<<"Fit chi2/ndf: "<<gFitFunction->GetChisquare()/gFitFunction->GetNDF()<<" - chi2: "<<gFitFunction->GetChisquare()<<" - ndf: "<<gFitFunction->GetNDF()<<endl; | |
629 | cout<<"======================================================"<<endl; | |
630 | ||
631 | //Getting the residual | |
632 | gHistResidual->Add(gFitFunction,-1); | |
742af4bd | 633 | |
634 | //Write to output file | |
d67eae53 | 635 | TString newFileName = "balanceFunctionFit2D."; |
636 | if(eventClass == "Centrality"){ | |
637 | newFileName += Form("Centrality%.1fTo%.1f",psiMin,psiMax); | |
638 | newFileName += ".PsiAll.PttFrom"; | |
639 | } | |
640 | else if(eventClass == "Multiplicity"){ | |
641 | newFileName += Form("Multiplicity%.0fTo%.0f",psiMin,psiMax); | |
642 | newFileName += ".PsiAll.PttFrom"; | |
643 | } | |
644 | else{ // "EventPlane" (default) | |
645 | newFileName += "Centrality"; | |
646 | newFileName += gCentrality; newFileName += ".Psi"; | |
647 | if((psiMin == -0.5)&&(psiMax == 0.5)) newFileName += "InPlane.Ptt"; | |
648 | else if((psiMin == 0.5)&&(psiMax == 1.5)) newFileName += "Intermediate.Ptt"; | |
649 | else if((psiMin == 1.5)&&(psiMax == 2.5)) newFileName += "OutOfPlane.Ptt"; | |
650 | else if((psiMin == 2.5)&&(psiMax == 3.5)) newFileName += "Rest.PttFrom"; | |
651 | else newFileName += "All.PttFrom"; | |
652 | } | |
742af4bd | 653 | newFileName += Form("%.1f",ptTriggerMin); newFileName += "To"; |
654 | newFileName += Form("%.1f",ptTriggerMax); newFileName += "PtaFrom"; | |
655 | newFileName += Form("%.1f",ptAssociatedMin); newFileName += "To"; | |
656 | newFileName += Form("%.1f",ptAssociatedMax); | |
d67eae53 | 657 | if(k2pMethod) newFileName += "_2pMethod"; |
a9f20288 | 658 | |
659 | newFileName += "_"; | |
660 | newFileName += Form("%.1f",psiMin); | |
661 | newFileName += "-"; | |
662 | newFileName += Form("%.1f",psiMax); | |
742af4bd | 663 | newFileName += ".root"; |
a9f20288 | 664 | |
742af4bd | 665 | TFile *newFile = TFile::Open(newFileName.Data(),"recreate"); |
666 | gHist->Write(); | |
742af4bd | 667 | gHistResidual->Write(); |
668 | gFitFunction->Write(); | |
669 | newFile->Close(); | |
742af4bd | 670 | } |
671 | ||
eb63b883 | 672 | //____________________________________________________________// |
648f1a5a | 673 | void drawBFPsi2D(const char* lhcPeriod = "LHC11h", |
bd36d661 | 674 | const char* gCentralityEstimator = "V0M", |
675 | Int_t gBit = 128, | |
676 | const char* gEventPlaneEstimator = "VZERO", | |
648f1a5a | 677 | Int_t gCentrality = 1, |
678 | Bool_t kShowShuffled = kFALSE, | |
679 | Bool_t kShowMixed = kFALSE, | |
eb63b883 | 680 | Double_t psiMin = -0.5, Double_t psiMax = 0.5, |
681 | Double_t ptTriggerMin = -1., | |
682 | Double_t ptTriggerMax = -1., | |
683 | Double_t ptAssociatedMin = -1., | |
bd36d661 | 684 | Double_t ptAssociatedMax = -1., |
685 | Bool_t k2pMethod = kTRUE) { | |
eb63b883 | 686 | //Macro that draws the BF distributions for each centrality bin |
687 | //for reaction plane dependent analysis | |
688 | //Author: Panos.Christakoglou@nikhef.nl | |
db7174c0 | 689 | TGaxis::SetMaxDigits(3); |
eb63b883 | 690 | |
691 | //Get the input file | |
bd36d661 | 692 | TString filename = lhcPeriod; |
693 | filename += "/Centrality"; filename += gCentralityEstimator; | |
694 | filename += "_Bit"; filename += gBit; | |
695 | filename += "_"; filename += gEventPlaneEstimator; | |
696 | filename +="/PttFrom"; | |
648f1a5a | 697 | filename += Form("%.1f",ptTriggerMin); filename += "To"; |
698 | filename += Form("%.1f",ptTriggerMax); filename += "PtaFrom"; | |
699 | filename += Form("%.1f",ptAssociatedMin); filename += "To"; | |
700 | filename += Form("%.1f",ptAssociatedMax); | |
701 | filename += "/balanceFunction2D.Centrality"; | |
eb63b883 | 702 | filename += gCentrality; filename += ".Psi"; |
703 | if((psiMin == -0.5)&&(psiMax == 0.5)) filename += "InPlane.Ptt"; | |
704 | else if((psiMin == 0.5)&&(psiMax == 1.5)) filename += "Intermediate.Ptt"; | |
705 | else if((psiMin == 1.5)&&(psiMax == 2.5)) filename += "OutOfPlane.Ptt"; | |
648f1a5a | 706 | else if((psiMin == 2.5)&&(psiMax == 3.5)) filename += "Rest.Ptt"; |
707 | else filename += "All.PttFrom"; | |
708 | filename += Form("%.1f",ptTriggerMin); filename += "To"; | |
db7174c0 | 709 | filename += Form("%.1f",ptTriggerMax); filename += "PtaFrom"; |
648f1a5a | 710 | filename += Form("%.1f",ptAssociatedMin); filename += "To"; |
bd36d661 | 711 | filename += Form("%.1f",ptAssociatedMax); |
712 | if(k2pMethod) filename += "_2pMethod"; | |
a9f20288 | 713 | |
714 | filename += "_"; | |
715 | filename += Form("%.1f",psiMin); | |
716 | filename += "-"; | |
717 | filename += Form("%.1f",psiMax); | |
bd36d661 | 718 | filename += ".root"; |
eb63b883 | 719 | |
720 | //Open the file | |
721 | TFile *f = TFile::Open(filename.Data()); | |
722 | if((!f)||(!f->IsOpen())) { | |
723 | Printf("The file %s is not found. Aborting...",filename); | |
724 | return listBF; | |
6acdbcb2 | 725 | } |
eb63b883 | 726 | //f->ls(); |
727 | ||
728 | //Raw balance function | |
729 | TH1D *gHistBalanceFunction = dynamic_cast<TH1D *>(f->Get("gHistBalanceFunction")); | |
730 | gHistBalanceFunction->SetStats(kFALSE); | |
731 | gHistBalanceFunction->GetXaxis()->SetNdivisions(10); | |
732 | gHistBalanceFunction->GetYaxis()->SetNdivisions(10); | |
733 | gHistBalanceFunction->GetZaxis()->SetNdivisions(10); | |
734 | gHistBalanceFunction->GetXaxis()->SetTitleOffset(1.3); | |
735 | gHistBalanceFunction->GetYaxis()->SetTitleOffset(1.3); | |
736 | gHistBalanceFunction->GetZaxis()->SetTitleOffset(1.3); | |
737 | gHistBalanceFunction->GetXaxis()->SetTitle("#Delta #eta"); | |
27634c13 | 738 | gHistBalanceFunction->GetYaxis()->SetTitle("#Delta #varphi (rad)"); |
15dd45a0 | 739 | gHistBalanceFunction->GetZaxis()->SetTitle("B(#Delta #eta, #Delta #varphi)"); |
eb63b883 | 740 | |
741 | //Shuffled balance function | |
648f1a5a | 742 | if(kShowShuffled) { |
743 | TH1D *gHistBalanceFunctionShuffled = dynamic_cast<TH1D *>(f->Get("gHistBalanceFunctionShuffled")); | |
744 | gHistBalanceFunctionShuffled->SetStats(kFALSE); | |
745 | gHistBalanceFunctionShuffled->GetXaxis()->SetNdivisions(10); | |
746 | gHistBalanceFunctionShuffled->GetYaxis()->SetNdivisions(10); | |
747 | gHistBalanceFunctionShuffled->GetZaxis()->SetNdivisions(10); | |
748 | gHistBalanceFunctionShuffled->GetXaxis()->SetTitleOffset(1.3); | |
749 | gHistBalanceFunctionShuffled->GetYaxis()->SetTitleOffset(1.3); | |
750 | gHistBalanceFunctionShuffled->GetZaxis()->SetTitleOffset(1.3); | |
751 | gHistBalanceFunctionShuffled->GetXaxis()->SetTitle("#Delta #eta"); | |
27634c13 | 752 | gHistBalanceFunctionShuffled->GetYaxis()->SetTitle("#Delta #varphi (rad)"); |
648f1a5a | 753 | gHistBalanceFunctionShuffled->GetZaxis()->SetTitle("B(#Delta #eta, #Delta #varphi)"); |
754 | } | |
eb63b883 | 755 | |
756 | //Mixed balance function | |
648f1a5a | 757 | if(kShowMixed) { |
758 | TH1D *gHistBalanceFunctionMixed = dynamic_cast<TH1D *>(f->Get("gHistBalanceFunctionMixed")); | |
759 | gHistBalanceFunctionMixed->SetStats(kFALSE); | |
760 | gHistBalanceFunctionMixed->GetXaxis()->SetNdivisions(10); | |
761 | gHistBalanceFunctionMixed->GetYaxis()->SetNdivisions(10); | |
762 | gHistBalanceFunctionMixed->GetZaxis()->SetNdivisions(10); | |
763 | gHistBalanceFunctionMixed->GetXaxis()->SetTitleOffset(1.3); | |
764 | gHistBalanceFunctionMixed->GetYaxis()->SetTitleOffset(1.3); | |
765 | gHistBalanceFunctionMixed->GetZaxis()->SetTitleOffset(1.3); | |
766 | gHistBalanceFunctionMixed->GetXaxis()->SetTitle("#Delta #eta"); | |
27634c13 | 767 | gHistBalanceFunctionMixed->GetYaxis()->SetTitle("#Delta #varphi (rad)"); |
648f1a5a | 768 | gHistBalanceFunctionMixed->GetZaxis()->SetTitle("B(#Delta #eta, #Delta #varphi)"); |
769 | } | |
eb63b883 | 770 | |
771 | //Subtracted balance function | |
648f1a5a | 772 | if(kShowMixed) { |
773 | TH1D *gHistBalanceFunctionSubtracted = dynamic_cast<TH1D *>(f->Get("gHistBalanceFunctionSubtracted")); | |
774 | gHistBalanceFunctionSubtracted->SetStats(kFALSE); | |
775 | gHistBalanceFunctionSubtracted->GetXaxis()->SetNdivisions(10); | |
776 | gHistBalanceFunctionSubtracted->GetYaxis()->SetNdivisions(10); | |
777 | gHistBalanceFunctionSubtracted->GetZaxis()->SetNdivisions(10); | |
778 | gHistBalanceFunctionSubtracted->GetXaxis()->SetTitleOffset(1.3); | |
779 | gHistBalanceFunctionSubtracted->GetYaxis()->SetTitleOffset(1.3); | |
780 | gHistBalanceFunctionSubtracted->GetZaxis()->SetTitleOffset(1.3); | |
781 | gHistBalanceFunctionSubtracted->GetXaxis()->SetTitle("#Delta #eta"); | |
27634c13 | 782 | gHistBalanceFunctionSubtracted->GetYaxis()->SetTitle("#Delta #varphi (rad)"); |
648f1a5a | 783 | gHistBalanceFunctionSubtracted->GetZaxis()->SetTitle("B(#Delta #eta, #Delta #varphi)"); |
784 | } | |
785 | ||
eb63b883 | 786 | TString pngName; |
787 | ||
788 | TString centralityLatex = "Centrality: "; | |
789 | centralityLatex += centralityArray[gCentrality-1]; | |
790 | centralityLatex += "%"; | |
791 | ||
792 | TString psiLatex; | |
793 | if((psiMin == -0.5)&&(psiMax == 0.5)) | |
794 | psiLatex = " -7.5^{o} < #varphi - #Psi_{2} < 7.5^{o}"; | |
795 | else if((psiMin == 0.5)&&(psiMax == 1.5)) | |
796 | psiLatex = " 37.5^{o} < #varphi - #Psi_{2} < 52.5^{o}"; | |
797 | else if((psiMin == 1.5)&&(psiMax == 2.5)) | |
798 | psiLatex = " 82.5^{o} < #varphi - #Psi_{2} < 97.5^{o}"; | |
799 | else | |
800 | psiLatex = " 0^{o} < #varphi - #Psi_{2} < 180^{o}"; | |
801 | ||
802 | TString pttLatex = Form("%.1f",ptTriggerMin); | |
803 | pttLatex += " < p_{T}^{t} < "; pttLatex += Form("%.1f",ptTriggerMax); | |
804 | pttLatex += " GeV/c"; | |
805 | ||
806 | TString ptaLatex = Form("%.1f",ptAssociatedMin); | |
807 | ptaLatex += " < p_{T}^{a} < "; ptaLatex += Form("%.1f",ptAssociatedMax); | |
808 | ptaLatex += " GeV/c"; | |
809 | ||
810 | TLatex *latexInfo1 = new TLatex(); | |
811 | latexInfo1->SetNDC(); | |
15dd45a0 | 812 | latexInfo1->SetTextSize(0.045); |
eb63b883 | 813 | latexInfo1->SetTextColor(1); |
814 | ||
815 | //Draw the results | |
816 | TCanvas *c1 = new TCanvas("c1","Raw balance function 2D",0,0,600,500); | |
817 | c1->SetFillColor(10); c1->SetHighLightColor(10); | |
818 | c1->SetLeftMargin(0.17); c1->SetTopMargin(0.05); | |
bd36d661 | 819 | gHistBalanceFunction->SetTitle(""); |
eb63b883 | 820 | gHistBalanceFunction->GetYaxis()->SetTitleOffset(1.4); |
821 | gHistBalanceFunction->GetYaxis()->SetNdivisions(10); | |
27634c13 | 822 | gHistBalanceFunction->GetXaxis()->SetRangeUser(-1.4,1.4); |
eb63b883 | 823 | gHistBalanceFunction->GetXaxis()->SetNdivisions(10); |
27634c13 | 824 | gHistBalanceFunction->GetYaxis()->SetTitle("#Delta #varphi (rad)"); |
eb63b883 | 825 | gHistBalanceFunction->DrawCopy("lego2"); |
5de9ad1a | 826 | gPad->SetTheta(30); // default is 30 |
827 | gPad->SetPhi(-60); // default is 30 | |
828 | gPad->Update(); | |
eb63b883 | 829 | |
15dd45a0 | 830 | latexInfo1->DrawLatex(0.64,0.88,centralityLatex.Data()); |
831 | latexInfo1->DrawLatex(0.64,0.82,psiLatex.Data()); | |
832 | latexInfo1->DrawLatex(0.64,0.76,pttLatex.Data()); | |
833 | latexInfo1->DrawLatex(0.64,0.70,ptaLatex.Data()); | |
eb63b883 | 834 | |
bd36d661 | 835 | TString pngName = "BalanceFunction2D."; |
836 | pngName += "Centrality"; | |
837 | pngName += gCentrality; | |
838 | if((psiMin == -0.5)&&(psiMax == 0.5)) pngName += "InPlane.Ptt"; | |
839 | else if((psiMin == 0.5)&&(psiMax == 1.5)) pngName += "Intermediate.Ptt"; | |
840 | else if((psiMin == 1.5)&&(psiMax == 2.5)) pngName += "OutOfPlane.Ptt"; | |
841 | else if((psiMin == 2.5)&&(psiMax == 3.5)) pngName += "Rest.PttFrom"; | |
842 | else pngName += "All.PttFrom"; | |
843 | pngName += Form("%.1f",ptTriggerMin); pngName += "To"; | |
844 | pngName += Form("%.1f",ptTriggerMax); pngName += "PtaFrom"; | |
845 | pngName += Form("%.1f",ptAssociatedMin); pngName += "To"; | |
846 | pngName += Form("%.1f",ptAssociatedMax); | |
847 | if(k2pMethod) pngName += "_2pMethod"; | |
848 | pngName += ".png"; | |
849 | ||
850 | c1->SaveAs(pngName.Data()); | |
851 | ||
648f1a5a | 852 | if(kShowShuffled) { |
853 | TCanvas *c2 = new TCanvas("c2","Shuffled balance function 2D",100,100,600,500); | |
854 | c2->SetFillColor(10); c2->SetHighLightColor(10); | |
855 | c2->SetLeftMargin(0.17); c2->SetTopMargin(0.05); | |
856 | gHistBalanceFunctionShuffled->SetTitle("Shuffled events"); | |
857 | gHistBalanceFunctionShuffled->GetYaxis()->SetTitleOffset(1.4); | |
858 | gHistBalanceFunctionShuffled->GetYaxis()->SetNdivisions(10); | |
859 | gHistBalanceFunctionShuffled->GetXaxis()->SetNdivisions(10); | |
860 | gHistBalanceFunctionShuffled->DrawCopy("lego2"); | |
5de9ad1a | 861 | gPad->SetTheta(30); // default is 30 |
862 | gPad->SetPhi(-60); // default is 30 | |
863 | gPad->Update(); | |
864 | ||
648f1a5a | 865 | latexInfo1->DrawLatex(0.64,0.88,centralityLatex.Data()); |
866 | latexInfo1->DrawLatex(0.64,0.82,psiLatex.Data()); | |
867 | latexInfo1->DrawLatex(0.64,0.76,pttLatex.Data()); | |
868 | latexInfo1->DrawLatex(0.64,0.70,ptaLatex.Data()); | |
869 | } | |
eb63b883 | 870 | |
648f1a5a | 871 | if(kShowMixed) { |
872 | TCanvas *c3 = new TCanvas("c3","Mixed balance function 2D",200,200,600,500); | |
873 | c3->SetFillColor(10); c3->SetHighLightColor(10); | |
874 | c3->SetLeftMargin(0.17); c3->SetTopMargin(0.05); | |
875 | gHistBalanceFunctionMixed->SetTitle("Mixed events"); | |
876 | gHistBalanceFunctionMixed->GetYaxis()->SetTitleOffset(1.4); | |
877 | gHistBalanceFunctionMixed->GetYaxis()->SetNdivisions(10); | |
878 | gHistBalanceFunctionMixed->GetXaxis()->SetNdivisions(10); | |
879 | gHistBalanceFunctionMixed->DrawCopy("lego2"); | |
5de9ad1a | 880 | gPad->SetTheta(30); // default is 30 |
881 | gPad->SetPhi(-60); // default is 30 | |
882 | gPad->Update(); | |
883 | ||
648f1a5a | 884 | latexInfo1->DrawLatex(0.64,0.88,centralityLatex.Data()); |
885 | latexInfo1->DrawLatex(0.64,0.82,psiLatex.Data()); | |
886 | latexInfo1->DrawLatex(0.64,0.76,pttLatex.Data()); | |
887 | latexInfo1->DrawLatex(0.64,0.70,ptaLatex.Data()); | |
888 | } | |
eb63b883 | 889 | |
648f1a5a | 890 | if(kShowMixed) { |
891 | TCanvas *c4 = new TCanvas("c4","Subtracted balance function 2D",300,300,600,500); | |
892 | c4->SetFillColor(10); c4->SetHighLightColor(10); | |
893 | c4->SetLeftMargin(0.17); c4->SetTopMargin(0.05); | |
894 | gHistBalanceFunctionSubtracted->SetTitle("Subtracted balance function"); | |
895 | gHistBalanceFunctionSubtracted->GetYaxis()->SetTitleOffset(1.4); | |
896 | gHistBalanceFunctionSubtracted->GetYaxis()->SetNdivisions(10); | |
897 | gHistBalanceFunctionSubtracted->GetXaxis()->SetNdivisions(10); | |
898 | gHistBalanceFunctionSubtracted->DrawCopy("lego2"); | |
5de9ad1a | 899 | gPad->SetTheta(30); // default is 30 |
900 | gPad->SetPhi(-60); // default is 30 | |
901 | gPad->Update(); | |
902 | ||
648f1a5a | 903 | latexInfo1->DrawLatex(0.64,0.88,centralityLatex.Data()); |
904 | latexInfo1->DrawLatex(0.64,0.82,psiLatex.Data()); | |
905 | latexInfo1->DrawLatex(0.64,0.76,pttLatex.Data()); | |
906 | latexInfo1->DrawLatex(0.64,0.70,ptaLatex.Data()); | |
907 | } | |
6acdbcb2 | 908 | } |
17c69b4e | 909 | |
910 | //____________________________________________________________// | |
79dc837f | 911 | void drawProjections(TH2D *gHistBalanceFunction2D = 0x0, |
912 | Bool_t kProjectInEta = kFALSE, | |
17c69b4e | 913 | Int_t binMin = 1, |
914 | Int_t binMax = 80, | |
915 | Int_t gCentrality = 1, | |
916 | Double_t psiMin = -0.5, | |
917 | Double_t psiMax = 3.5, | |
17c69b4e | 918 | Double_t ptTriggerMin = -1., |
919 | Double_t ptTriggerMax = -1., | |
920 | Double_t ptAssociatedMin = -1., | |
921 | Double_t ptAssociatedMax = -1., | |
17c69b4e | 922 | Bool_t k2pMethod = kTRUE, |
923 | TString eventClass = "Centrality", | |
924 | Bool_t bRootMoments = kFALSE) { | |
438dab2c | 925 | //Macro that draws the balance functions PROJECTIONS |
17c69b4e | 926 | //for each centrality bin for the different pT of trigger and |
927 | //associated particles | |
928 | TGaxis::SetMaxDigits(3); | |
929 | ||
930 | //first we need some libraries | |
a9f20288 | 931 | gSystem->Load("libTree"); |
932 | gSystem->Load("libGeom"); | |
933 | gSystem->Load("libVMC"); | |
934 | gSystem->Load("libXMLIO"); | |
935 | gSystem->Load("libPhysics"); | |
936 | ||
937 | gSystem->Load("libSTEERBase"); | |
938 | gSystem->Load("libESD"); | |
939 | gSystem->Load("libAOD"); | |
940 | ||
17c69b4e | 941 | gSystem->Load("libANALYSIS.so"); |
942 | gSystem->Load("libANALYSISalice.so"); | |
943 | gSystem->Load("libEventMixing.so"); | |
944 | gSystem->Load("libCORRFW.so"); | |
945 | gSystem->Load("libPWGTools.so"); | |
946 | gSystem->Load("libPWGCFebye.so"); | |
947 | ||
438dab2c | 948 | gStyle->SetOptStat(0); |
949 | ||
17c69b4e | 950 | //Get the input file |
951 | TString filename = "balanceFunction2D."; | |
952 | if(eventClass == "Centrality"){ | |
953 | filename += Form("Centrality%.1fTo%.1f",psiMin,psiMax); | |
954 | filename += ".PsiAll.PttFrom"; | |
955 | } | |
956 | else if(eventClass == "Multiplicity"){ | |
957 | filename += Form("Multiplicity%.0fTo%.0f",psiMin,psiMax); | |
958 | filename += ".PsiAll.PttFrom"; | |
959 | } | |
960 | else{ // "EventPlane" (default) | |
961 | filename += "Centrality"; | |
962 | filename += gCentrality; filename += ".Psi"; | |
963 | if((psiMin == -0.5)&&(psiMax == 0.5)) filename += "InPlane.Ptt"; | |
964 | else if((psiMin == 0.5)&&(psiMax == 1.5)) filename += "Intermediate.Ptt"; | |
965 | else if((psiMin == 1.5)&&(psiMax == 2.5)) filename += "OutOfPlane.Ptt"; | |
966 | else if((psiMin == 2.5)&&(psiMax == 3.5)) filename += "Rest.PttFrom"; | |
967 | else filename += "All.PttFrom"; | |
968 | } | |
969 | filename += Form("%.1f",ptTriggerMin); filename += "To"; | |
970 | filename += Form("%.1f",ptTriggerMax); filename += "PtaFrom"; | |
971 | filename += Form("%.1f",ptAssociatedMin); filename += "To"; | |
972 | filename += Form("%.1f",ptAssociatedMax); | |
973 | if(k2pMethod) filename += "_2pMethod"; | |
a9f20288 | 974 | |
a2da950e | 975 | filename += "_"; |
976 | filename += Form("%.1f",psiMin); | |
977 | filename += "-"; | |
978 | filename += Form("%.1f",psiMax); | |
17c69b4e | 979 | filename += ".root"; |
980 | ||
981 | //Open the file | |
79dc837f | 982 | TFile *f = 0x0; |
983 | if(!gHistBalanceFunction2D) { | |
984 | TFile::Open(filename.Data()); | |
985 | if((!f)||(!f->IsOpen())) { | |
986 | Printf("The file %s is not found. Aborting...",filename); | |
987 | return listBF; | |
988 | } | |
989 | //f->ls(); | |
17c69b4e | 990 | } |
79dc837f | 991 | |
17c69b4e | 992 | //Latex |
993 | TString centralityLatex = "Centrality: "; | |
994 | centralityLatex += centralityArray[gCentrality-1]; | |
995 | centralityLatex += "%"; | |
996 | ||
997 | TString psiLatex; | |
998 | if((psiMin == -0.5)&&(psiMax == 0.5)) | |
999 | psiLatex = " -7.5^{o} < #varphi^{t} - #Psi_{2} < 7.5^{o}"; | |
1000 | else if((psiMin == 0.5)&&(psiMax == 1.5)) | |
1001 | psiLatex = " 37.5^{o} < #varphi^{t} - #Psi_{2} < 52.5^{o}"; | |
1002 | else if((psiMin == 1.5)&&(psiMax == 2.5)) | |
1003 | psiLatex = " 82.5^{o} < #varphi^{t} - #Psi_{2} < 97.5^{o}"; | |
1004 | else | |
1005 | psiLatex = " 0^{o} < #varphi^{t} - #Psi_{2} < 180^{o}"; | |
1006 | ||
1007 | TString pttLatex = Form("%.1f",ptTriggerMin); | |
1008 | pttLatex += " < p_{T,trig} < "; pttLatex += Form("%.1f",ptTriggerMax); | |
1009 | pttLatex += " GeV/c"; | |
1010 | ||
1011 | TString ptaLatex = Form("%.1f",ptAssociatedMin); | |
1012 | ptaLatex += " < p_{T,assoc} < "; ptaLatex += Form("%.1f",ptAssociatedMax); | |
1013 | ptaLatex += " GeV/c"; | |
1014 | ||
1015 | TLatex *latexInfo1 = new TLatex(); | |
1016 | latexInfo1->SetNDC(); | |
1017 | latexInfo1->SetTextSize(0.045); | |
1018 | latexInfo1->SetTextColor(1); | |
1019 | ||
17c69b4e | 1020 | |
1021 | //============================================================// | |
1022 | //Get subtracted and mixed balance function | |
79dc837f | 1023 | TH2D *gHistBalanceFunctionSubtracted2D = 0x0; |
1024 | TH2D *gHistBalanceFunctionMixed2D = 0x0; | |
1025 | if(!gHistBalanceFunction2D) { | |
1026 | gHistBalanceFunctionSubtracted2D = (TH2D*)f->Get("gHistBalanceFunctionSubtracted"); | |
1027 | gHistBalanceFunctionMixed2D = (TH2D*)f->Get("gHistBalanceFunctionMixed"); | |
1028 | } | |
1029 | else { | |
1030 | gHistBalanceFunctionSubtracted2D = dynamic_cast<TH2D*>(gHistBalanceFunction2D->Clone()); | |
1031 | gHistBalanceFunctionMixed2D = dynamic_cast<TH2D*>(gHistBalanceFunction2D->Clone()); | |
1032 | } | |
17c69b4e | 1033 | |
1034 | TH1D *gHistBalanceFunctionSubtracted = NULL; | |
1035 | TH1D *gHistBalanceFunctionMixed = NULL; | |
1036 | ||
a9f20288 | 1037 | TH1D *gHistBalanceFunctionSubtracted_scale = NULL; |
1038 | TH1D *gHistBalanceFunctionMixed_scale = NULL; | |
1039 | ||
17c69b4e | 1040 | if(kProjectInEta){ |
86044267 | 1041 | gHistBalanceFunctionSubtracted = dynamic_cast<TH1D *>(gHistBalanceFunctionSubtracted2D->ProjectionX("gHistBalanceFunctionSubtractedEta",binMin,binMax)); |
a9f20288 | 1042 | gHistBalanceFunctionSubtracted->Scale(gHistBalanceFunctionSubtracted2D->GetYaxis()->GetBinWidth(1)); // to remove normalization to phi bin width |
86044267 | 1043 | gHistBalanceFunctionMixed = dynamic_cast<TH1D *>(gHistBalanceFunctionMixed2D->ProjectionX("gHistBalanceFunctionMixedEta",binMin,binMax)); |
a9f20288 | 1044 | gHistBalanceFunctionMixed->Scale(gHistBalanceFunctionMixed2D->GetYaxis()->GetBinWidth(1)); // to remove normalization to phi bin width |
17c69b4e | 1045 | gHistBalanceFunctionSubtracted->SetTitle("B(#Delta#eta)"); |
1046 | gHistBalanceFunctionMixed->SetTitle("B_{mix}(#Delta#eta)"); | |
1047 | } | |
1048 | else{ | |
86044267 | 1049 | gHistBalanceFunctionSubtracted = dynamic_cast<TH1D *>(gHistBalanceFunctionSubtracted2D->ProjectionY("gHistBalanceFunctionSubtractedPhi",binMin,binMax)); |
a9f20288 | 1050 | gHistBalanceFunctionSubtracted->Scale(gHistBalanceFunctionSubtracted2D->GetXaxis()->GetBinWidth(1)); // to remove normalization to eta bin width |
86044267 | 1051 | gHistBalanceFunctionMixed = dynamic_cast<TH1D *>(gHistBalanceFunctionMixed2D->ProjectionY("gHistBalanceFunctionMixedPhi",binMin,binMax)); |
a9f20288 | 1052 | gHistBalanceFunctionMixed->Scale(gHistBalanceFunctionMixed2D->GetXaxis()->GetBinWidth(1)); // to remove normalization to eta bin width |
17c69b4e | 1053 | gHistBalanceFunctionSubtracted->SetTitle("B(#Delta#varphi)"); |
1054 | gHistBalanceFunctionMixed->SetTitle("B_{mix}(#Delta#varphi)"); | |
1055 | } | |
1056 | ||
1057 | gHistBalanceFunctionSubtracted->SetMarkerStyle(20); | |
1058 | gHistBalanceFunctionSubtracted->GetYaxis()->SetTitleOffset(1.3); | |
17c69b4e | 1059 | |
1060 | gHistBalanceFunctionMixed->SetMarkerStyle(25); | |
17c69b4e | 1061 | |
1062 | TCanvas *c1 = new TCanvas("c1","",0,0,600,500); | |
1063 | c1->SetFillColor(10); | |
1064 | c1->SetHighLightColor(10); | |
1065 | c1->SetLeftMargin(0.15); | |
1066 | gHistBalanceFunctionSubtracted->DrawCopy("E"); | |
a9f20288 | 1067 | gHistBalanceFunctionMixed->DrawCopy("ESAME"); |
17c69b4e | 1068 | |
1069 | legend = new TLegend(0.18,0.62,0.45,0.82,"","brNDC"); | |
1070 | legend->SetTextSize(0.045); | |
1071 | legend->SetTextFont(42); | |
1072 | legend->SetBorderSize(0); | |
1073 | legend->SetFillStyle(0); | |
1074 | legend->SetFillColor(10); | |
1075 | legend->SetMargin(0.25); | |
1076 | legend->SetShadowColor(10); | |
1077 | legend->AddEntry(gHistBalanceFunctionSubtracted,"Data","lp"); | |
1078 | legend->AddEntry(gHistBalanceFunctionMixed,"Mixed data","lp"); | |
1079 | legend->Draw(); | |
1080 | ||
17c69b4e | 1081 | |
17c69b4e | 1082 | TString meanLatex, rmsLatex, skewnessLatex, kurtosisLatex; |
1083 | ||
1084 | if(bRootMoments){ | |
1085 | meanLatex = "#mu = "; | |
1086 | meanLatex += Form("%.3f",gHistBalanceFunctionSubtracted->GetMean()); | |
1087 | meanLatex += " #pm "; | |
1088 | meanLatex += Form("%.3f",gHistBalanceFunctionSubtracted->GetMeanError()); | |
1089 | ||
1090 | rmsLatex = "#sigma = "; | |
1091 | rmsLatex += Form("%.3f",gHistBalanceFunctionSubtracted->GetRMS()); | |
1092 | rmsLatex += " #pm "; | |
1093 | rmsLatex += Form("%.3f",gHistBalanceFunctionSubtracted->GetRMSError()); | |
1094 | ||
1095 | skewnessLatex = "S = "; | |
1096 | skewnessLatex += Form("%.3f",gHistBalanceFunctionSubtracted->GetSkewness(1)); | |
1097 | skewnessLatex += " #pm "; | |
1098 | skewnessLatex += Form("%.3f",gHistBalanceFunctionSubtracted->GetSkewness(11)); | |
1099 | ||
1100 | kurtosisLatex = "K = "; | |
1101 | kurtosisLatex += Form("%.3f",gHistBalanceFunctionSubtracted->GetKurtosis(1)); | |
1102 | kurtosisLatex += " #pm "; | |
1103 | kurtosisLatex += Form("%.3f",gHistBalanceFunctionSubtracted->GetKurtosis(11)); | |
a2da950e | 1104 | |
17c69b4e | 1105 | Printf("Mean: %lf - Error: %lf",gHistBalanceFunctionSubtracted->GetMean(),gHistBalanceFunctionSubtracted->GetMeanError()); |
1106 | Printf("RMS: %lf - Error: %lf",gHistBalanceFunctionSubtracted->GetRMS(),gHistBalanceFunctionSubtracted->GetRMSError()); | |
1107 | Printf("Skeweness: %lf - Error: %lf",gHistBalanceFunctionSubtracted->GetSkewness(1),gHistBalanceFunctionSubtracted->GetSkewness(11)); | |
1108 | Printf("Kurtosis: %lf - Error: %lf",gHistBalanceFunctionSubtracted->GetKurtosis(1),gHistBalanceFunctionSubtracted->GetKurtosis(11)); | |
438dab2c | 1109 | |
1110 | ||
1111 | // store in txt files | |
438dab2c | 1112 | TString meanFileName = filename; |
79dc837f | 1113 | if(kProjectInEta) |
1114 | meanFileName= "deltaEtaProjection_Mean.txt"; | |
1115 | //meanFileName.ReplaceAll(".root","_DeltaEtaProjection_Mean.txt"); | |
1116 | else | |
1117 | meanFileName = "deltaPhiProjection_Mean.txt"; | |
1118 | //meanFileName.ReplaceAll(".root","_DeltaPhiProjection_Mean.txt"); | |
1119 | ofstream fileMean(meanFileName.Data(),ios::app); | |
438dab2c | 1120 | fileMean << " " << gHistBalanceFunctionSubtracted->GetMean() << " " <<gHistBalanceFunctionSubtracted->GetMeanError()<<endl; |
1121 | fileMean.close(); | |
1122 | ||
1123 | TString rmsFileName = filename; | |
79dc837f | 1124 | if(kProjectInEta) |
1125 | rmsFileName = "deltaEtaProjection_Rms.txt"; | |
1126 | //rmsFileName.ReplaceAll(".root","_DeltaEtaProjection_Rms.txt"); | |
1127 | else | |
1128 | rmsFileName = "deltaPhiProjection_Rms.txt"); | |
1129 | //rmsFileName.ReplaceAll(".root","_DeltaPhiProjection_Rms.txt"); | |
1130 | ofstream fileRms(rmsFileName.Data(),ios::app); | |
438dab2c | 1131 | fileRms << " " << gHistBalanceFunctionSubtracted->GetRMS() << " " <<gHistBalanceFunctionSubtracted->GetRMSError()<<endl; |
1132 | fileRms.close(); | |
1133 | ||
1134 | TString skewnessFileName = filename; | |
79dc837f | 1135 | if(kProjectInEta) |
1136 | skewnessFileName = "deltaEtaProjection_Skewness.txt"; | |
1137 | //skewnessFileName.ReplaceAll(".root","_DeltaEtaProjection_Skewness.txt"); | |
1138 | else | |
1139 | skewnessFileName = "deltaPhiProjection_Skewness.txt"; | |
1140 | //skewnessFileName.ReplaceAll(".root","_DeltaPhiProjection_Skewness.txt"); | |
1141 | ofstream fileSkewness(skewnessFileName.Data(),ios::app); | |
438dab2c | 1142 | fileSkewness << " " << gHistBalanceFunctionSubtracted->GetSkewness(1) << " " <<gHistBalanceFunctionSubtracted->GetSkewness(11)<<endl; |
1143 | fileSkewness.close(); | |
1144 | ||
1145 | TString kurtosisFileName = filename; | |
79dc837f | 1146 | if(kProjectInEta) |
1147 | kurtosisFileName = "deltaEtaProjection_Kurtosis.txt"; | |
1148 | //kurtosisFileName.ReplaceAll(".root","_DeltaEtaProjection_Kurtosis.txt"); | |
1149 | else | |
1150 | kurtosisFileName = "deltaPhiProjection_Kurtosis.txt"; | |
1151 | //kurtosisFileName.ReplaceAll(".root","_DeltaPhiProjection_Kurtosis.txt"); | |
1152 | ofstream fileKurtosis(kurtosisFileName.Data(),ios::app); | |
438dab2c | 1153 | fileKurtosis << " " << gHistBalanceFunctionSubtracted->GetKurtosis(1) << " " <<gHistBalanceFunctionSubtracted->GetKurtosis(11)<<endl; |
1154 | fileKurtosis.close(); | |
17c69b4e | 1155 | } |
1156 | // calculate the moments by hand | |
1157 | else{ | |
1158 | ||
1159 | Double_t meanAnalytical, meanAnalyticalError; | |
1160 | Double_t sigmaAnalytical, sigmaAnalyticalError; | |
1161 | Double_t skewnessAnalytical, skewnessAnalyticalError; | |
1162 | Double_t kurtosisAnalytical, kurtosisAnalyticalError; | |
1163 | ||
1164 | Int_t gDeltaEtaPhi = 2; | |
1165 | if(kProjectInEta) gDeltaEtaPhi = 1; | |
1166 | ||
1167 | AliBalancePsi *bHelper = new AliBalancePsi; | |
1168 | bHelper->GetMomentsAnalytical(gDeltaEtaPhi,gHistBalanceFunctionSubtracted,meanAnalytical,meanAnalyticalError,sigmaAnalytical,sigmaAnalyticalError,skewnessAnalytical,skewnessAnalyticalError,kurtosisAnalytical,kurtosisAnalyticalError); | |
1169 | ||
1170 | meanLatex = "#mu = "; | |
1171 | meanLatex += Form("%.3f",meanAnalytical); | |
1172 | meanLatex += " #pm "; | |
1173 | meanLatex += Form("%.3f",meanAnalyticalError); | |
1174 | ||
1175 | rmsLatex = "#sigma = "; | |
1176 | rmsLatex += Form("%.3f",sigmaAnalytical); | |
1177 | rmsLatex += " #pm "; | |
1178 | rmsLatex += Form("%.3f",sigmaAnalyticalError); | |
1179 | ||
1180 | skewnessLatex = "S = "; | |
1181 | skewnessLatex += Form("%.3f",skewnessAnalytical); | |
1182 | skewnessLatex += " #pm "; | |
1183 | skewnessLatex += Form("%.3f",skewnessAnalyticalError); | |
1184 | ||
1185 | kurtosisLatex = "K = "; | |
1186 | kurtosisLatex += Form("%.3f",kurtosisAnalytical); | |
1187 | kurtosisLatex += " #pm "; | |
1188 | kurtosisLatex += Form("%.3f",kurtosisAnalyticalError); | |
1189 | Printf("Mean: %lf - Error: %lf",meanAnalytical, meanAnalyticalError); | |
1190 | Printf("Sigma: %lf - Error: %lf",sigmaAnalytical, sigmaAnalyticalError); | |
1191 | Printf("Skeweness: %lf - Error: %lf",skewnessAnalytical, skewnessAnalyticalError); | |
1192 | Printf("Kurtosis: %lf - Error: %lf",kurtosisAnalytical, kurtosisAnalyticalError); | |
438dab2c | 1193 | |
1194 | // store in txt files | |
1195 | TString meanFileName = filename; | |
79dc837f | 1196 | if(kProjectInEta) |
1197 | meanFileName = "deltaEtaProjection_Mean.txt"; | |
1198 | //meanFileName.ReplaceAll(".root","_DeltaEtaProjection_Mean.txt"); | |
1199 | else | |
1200 | meanFileName = "deltaPhiProjection_Mean.txt"; | |
1201 | //meanFileName.ReplaceAll(".root","_DeltaPhiProjection_Mean.txt"); | |
1202 | ofstream fileMean(meanFileName.Data(),ios::app); | |
438dab2c | 1203 | fileMean << " " << meanAnalytical << " " <<meanAnalyticalError<<endl; |
1204 | fileMean.close(); | |
1205 | ||
1206 | TString rmsFileName = filename; | |
79dc837f | 1207 | if(kProjectInEta) |
1208 | rmsFileName = "deltaEtaProjection_Rms.txt"; | |
1209 | //rmsFileName.ReplaceAll(".root","_DeltaEtaProjection_Rms.txt"); | |
1210 | else | |
1211 | rmsFileName = "deltaPhiProjection_Rms.txt"; | |
1212 | //rmsFileName.ReplaceAll(".root","_DeltaPhiProjection_Rms.txt"); | |
1213 | ofstream fileRms(rmsFileName.Data(),ios::app); | |
438dab2c | 1214 | fileRms << " " << sigmaAnalytical << " " <<sigmaAnalyticalError<<endl; |
1215 | fileRms.close(); | |
1216 | ||
1217 | TString skewnessFileName = filename; | |
79dc837f | 1218 | if(kProjectInEta) |
1219 | skewnessFileName = "deltaEtaProjection_Skewness.txt"; | |
1220 | //skewnessFileName.ReplaceAll(".root","_DeltaEtaProjection_Skewness.txt"); | |
1221 | else | |
1222 | skewnessFileName = "deltaPhiProjection_Skewness.txt"; | |
1223 | //skewnessFileName.ReplaceAll(".root","_DeltaPhiProjection_Skewness.txt"); | |
1224 | ofstream fileSkewness(skewnessFileName.Data(),ios::app); | |
438dab2c | 1225 | fileSkewness << " " << skewnessAnalytical << " " <<skewnessAnalyticalError<<endl; |
1226 | fileSkewness.close(); | |
1227 | ||
1228 | TString kurtosisFileName = filename; | |
79dc837f | 1229 | if(kProjectInEta) |
1230 | kurtosisFileName = "deltaEtaProjection_Kurtosis.txt"; | |
1231 | //kurtosisFileName.ReplaceAll(".root","_DeltaEtaProjection_Kurtosis.txt"); | |
1232 | else | |
1233 | kurtosisFileName = "deltaPhiProjection_Kurtosis.txt"; | |
1234 | //kurtosisFileName.ReplaceAll(".root","_DeltaPhiProjection_Kurtosis.txt"); | |
1235 | ofstream fileKurtosis(kurtosisFileName.Data(),ios::app); | |
438dab2c | 1236 | fileKurtosis << " " << kurtosisAnalytical << " " <<kurtosisAnalyticalError<<endl; |
1237 | fileKurtosis.close(); | |
17c69b4e | 1238 | } |
1239 | ||
438dab2c | 1240 | |
1241 | ||
17c69b4e | 1242 | TCanvas *c2 = new TCanvas("c2","",600,0,600,500); |
1243 | c2->SetFillColor(10); | |
1244 | c2->SetHighLightColor(10); | |
1245 | c2->SetLeftMargin(0.15); | |
1246 | gHistBalanceFunctionSubtracted->DrawCopy("E"); | |
438dab2c | 1247 | |
17c69b4e | 1248 | TLatex *latex = new TLatex(); |
1249 | latex->SetNDC(); | |
1250 | latex->SetTextSize(0.035); | |
1251 | latex->SetTextColor(1); | |
1252 | latex->DrawLatex(0.64,0.85,meanLatex.Data()); | |
1253 | latex->DrawLatex(0.64,0.81,rmsLatex.Data()); | |
1254 | latex->DrawLatex(0.64,0.77,skewnessLatex.Data()); | |
1255 | latex->DrawLatex(0.64,0.73,kurtosisLatex.Data()); | |
1256 | ||
438dab2c | 1257 | TString pngName = filename; |
1258 | if(kProjectInEta) pngName.ReplaceAll(".root","_DeltaEtaProjection.png"); | |
1259 | else pngName.ReplaceAll(".root","_DeltaPhiProjection.png"); | |
1260 | ||
1261 | c2->SaveAs(pngName.Data()); | |
1262 | ||
17c69b4e | 1263 | TString outFileName = filename; |
1264 | if(kProjectInEta) outFileName.ReplaceAll(".root","_DeltaEtaProjection.root"); | |
1265 | else outFileName.ReplaceAll(".root","_DeltaPhiProjection.root"); | |
1266 | TFile *fProjection = TFile::Open(outFileName.Data(),"recreate"); | |
1267 | gHistBalanceFunctionSubtracted->Write(); | |
1268 | gHistBalanceFunctionMixed->Write(); | |
1269 | fProjection->Close(); | |
1270 | } | |
79dc837f | 1271 | |
1272 | //____________________________________________________________// | |
1273 | void drawBFPsi2DFromCorrelationFunctions(const char* lhcPeriod = "LHC10h", | |
1274 | Int_t gTrainNumber = 64, | |
1275 | const char* gCentralityEstimator = "V0M", | |
1276 | Int_t gBit = 128, | |
1277 | const char* gEventPlaneEstimator = "VZERO", | |
1278 | Int_t gCentrality = 1, | |
1279 | Double_t psiMin = -0.5, Double_t psiMax = 3.5, | |
1280 | Double_t vertexZMin = -10., | |
1281 | Double_t vertexZMax = 10., | |
1282 | Double_t ptTriggerMin = -1., | |
1283 | Double_t ptTriggerMax = -1., | |
1284 | Double_t ptAssociatedMin = -1., | |
86044267 | 1285 | Double_t ptAssociatedMax = -1., |
1286 | TString eventClass = "Multiplicity" | |
1287 | ) { | |
79dc837f | 1288 | //Macro that draws the BF distributions for each centrality bin |
1289 | //for reaction plane dependent analysis | |
1290 | //Author: Panos.Christakoglou@nikhef.nl | |
1291 | TGaxis::SetMaxDigits(3); | |
1292 | gStyle->SetPalette(55,0); | |
1293 | ||
1294 | //Get the input file | |
1295 | TString filename = lhcPeriod; | |
86044267 | 1296 | if(lhcPeriod != ""){ |
1297 | //filename += "/Train"; filename += gTrainNumber; | |
1298 | filename +="/PttFrom"; | |
1299 | filename += Form("%.1f",ptTriggerMin); filename += "To"; | |
1300 | filename += Form("%.1f",ptTriggerMax); filename += "PtaFrom"; | |
1301 | filename += Form("%.1f",ptAssociatedMin); filename += "To"; | |
1302 | filename += Form("%.1f",ptAssociatedMax); | |
1303 | filename += "/correlationFunction."; | |
1304 | } | |
1305 | else{ | |
1306 | filename += "correlationFunction."; | |
1307 | } | |
1308 | if(eventClass == "Centrality"){ | |
1309 | filename += Form("Centrality%.1fTo%.1f",psiMin,psiMax); | |
1310 | filename += ".PsiAll.PttFrom"; | |
1311 | } | |
1312 | else if(eventClass == "Multiplicity"){ | |
1313 | filename += Form("Multiplicity%.0fTo%.0f",psiMin,psiMax); | |
1314 | filename += ".PsiAll.PttFrom"; | |
1315 | } | |
1316 | else{ // "EventPlane" (default) | |
1317 | filename += "Centrality"; | |
1318 | filename += gCentrality; filename += ".Psi"; | |
1319 | if((psiMin == -0.5)&&(psiMax == 0.5)) filename += "InPlane.Ptt"; | |
1320 | else if((psiMin == 0.5)&&(psiMax == 1.5)) filename += "Intermediate.Ptt"; | |
1321 | else if((psiMin == 1.5)&&(psiMax == 2.5)) filename += "OutOfPlane.Ptt"; | |
1322 | else if((psiMin == 2.5)&&(psiMax == 3.5)) filename += "Rest.PttFrom"; | |
1323 | else filename += "All.PttFrom"; | |
1324 | } | |
79dc837f | 1325 | filename += Form("%.1f",ptTriggerMin); filename += "To"; |
1326 | filename += Form("%.1f",ptTriggerMax); filename += "PtaFrom"; | |
1327 | filename += Form("%.1f",ptAssociatedMin); filename += "To"; | |
1328 | filename += Form("%.1f",ptAssociatedMax); | |
1329 | filename += "_"; | |
1330 | filename += Form("%.1f",psiMin); | |
1331 | filename += "-"; | |
1332 | filename += Form("%.1f",psiMax); | |
1333 | filename += ".root"; | |
1334 | ||
1335 | //Open the file | |
1336 | TFile *f = TFile::Open(filename.Data()); | |
1337 | if((!f)||(!f->IsOpen())) { | |
1338 | Printf("The file %s is not found. Aborting...",filename); | |
1339 | return listBF; | |
1340 | } | |
1341 | //f->ls(); | |
1342 | ||
1343 | TH2D *gHistPN = dynamic_cast<TH2D *>(f->Get("gHistPNCorrelationFunctions")); | |
1344 | if(!gHistPN) return; | |
1345 | TH2D *gHistNP = dynamic_cast<TH2D *>(f->Get("gHistNPCorrelationFunctions")); | |
1346 | if(!gHistNP) return; | |
1347 | TH2D *gHistPP = dynamic_cast<TH2D *>(f->Get("gHistPPCorrelationFunctions")); | |
1348 | if(!gHistPP) return; | |
1349 | TH2D *gHistNN = dynamic_cast<TH2D *>(f->Get("gHistNNCorrelationFunctions")); | |
1350 | if(!gHistNN) return; | |
1351 | ||
1352 | gHistPN->Sumw2(); | |
1353 | gHistPP->Sumw2(); | |
1354 | gHistPN->Add(gHistPP,-1); | |
1355 | gHistNP->Sumw2(); | |
1356 | gHistNN->Sumw2(); | |
1357 | gHistNP->Add(gHistNN,-1); | |
1358 | gHistPN->Add(gHistNP); | |
1359 | gHistPN->Scale(0.5); | |
1360 | TH2D *gHistBalanceFunction2D = dynamic_cast<TH2D *>(gHistPN->Clone()); | |
1361 | gHistBalanceFunction2D->SetStats(kFALSE); | |
1362 | gHistBalanceFunction2D->GetXaxis()->SetTitle("#Delta#eta"); | |
1363 | gHistBalanceFunction2D->GetYaxis()->SetTitle("#Delta#varphi (rad)"); | |
1364 | gHistBalanceFunction2D->GetZaxis()->SetTitle("B(#Delta#eta,#Delta#varphi)"); | |
1365 | ||
1366 | //Draw the results | |
1367 | TCanvas *c0 = new TCanvas("c0","Balance function 2D",0,0,600,500); | |
1368 | c0->SetFillColor(10); c0->SetHighLightColor(10); | |
1369 | c0->SetLeftMargin(0.17); c0->SetTopMargin(0.05); | |
1370 | gHistBalanceFunction2D->SetTitle(""); | |
1371 | gHistBalanceFunction2D->GetZaxis()->SetTitleOffset(1.4); | |
1372 | gHistBalanceFunction2D->GetZaxis()->SetNdivisions(10); | |
1373 | gHistBalanceFunction2D->GetYaxis()->SetTitleOffset(1.4); | |
1374 | gHistBalanceFunction2D->GetYaxis()->SetNdivisions(10); | |
1375 | gHistBalanceFunction2D->GetXaxis()->SetRangeUser(-1.4,1.4); | |
1376 | gHistBalanceFunction2D->GetXaxis()->SetNdivisions(10); | |
1377 | gHistBalanceFunction2D->DrawCopy("lego2"); | |
1378 | gPad->SetTheta(30); // default is 30 | |
1379 | gPad->SetPhi(-60); // default is 30 | |
1380 | gPad->Update(); | |
1381 | ||
1382 | TString multLatex = Form("Multiplicity: %.1f - %.1f",psiMin,psiMax); | |
1383 | ||
1384 | TString pttLatex = Form("%.1f",ptTriggerMin); | |
1385 | pttLatex += " < p_{T}^{t} < "; pttLatex += Form("%.1f",ptTriggerMax); | |
1386 | pttLatex += " GeV/c"; | |
1387 | ||
1388 | TString ptaLatex = Form("%.1f",ptAssociatedMin); | |
1389 | ptaLatex += " < p_{T}^{a} < "; ptaLatex += Form("%.1f",ptAssociatedMax); | |
1390 | ptaLatex += " GeV/c"; | |
1391 | ||
1392 | TLatex *latexInfo1 = new TLatex(); | |
1393 | latexInfo1->SetNDC(); | |
1394 | latexInfo1->SetTextSize(0.045); | |
1395 | latexInfo1->SetTextColor(1); | |
1396 | latexInfo1->DrawLatex(0.54,0.88,multLatex.Data()); | |
1397 | latexInfo1->DrawLatex(0.54,0.82,pttLatex.Data()); | |
1398 | latexInfo1->DrawLatex(0.54,0.76,ptaLatex.Data()); | |
1399 | ||
1400 | TString pngName = "BalanceFunction2D."; | |
86044267 | 1401 | pngName += Form("%s: %.1f - %.1f",eventClass.Data(),psiMin,psiMax); |
79dc837f | 1402 | pngName += ".PttFrom"; |
1403 | pngName += Form("%.1f",ptTriggerMin); pngName += "To"; | |
1404 | pngName += Form("%.1f",ptTriggerMax); pngName += "PtaFrom"; | |
1405 | pngName += Form("%.1f",ptAssociatedMin); pngName += "To"; | |
1406 | pngName += Form("%.1f",ptAssociatedMax); | |
1407 | pngName += ".png"; | |
1408 | ||
1409 | c0->SaveAs(pngName.Data()); | |
1410 | ||
1411 | drawProjections(gHistBalanceFunction2D, | |
1412 | kTRUE, | |
1413 | 1,36, | |
1414 | gCentrality, | |
1415 | psiMin,psiMax, | |
1416 | ptTriggerMin,ptTriggerMax, | |
1417 | ptAssociatedMin,ptAssociatedMax, | |
1418 | kTRUE, | |
86044267 | 1419 | eventClass, |
79dc837f | 1420 | kFALSE); |
1421 | ||
1422 | drawProjections(gHistBalanceFunction2D, | |
86044267 | 1423 | kFALSE, |
1424 | 1,80, | |
1425 | gCentrality, | |
1426 | psiMin,psiMax, | |
1427 | ptTriggerMin,ptTriggerMax, | |
1428 | ptAssociatedMin,ptAssociatedMax, | |
1429 | kTRUE, | |
1430 | eventClass.Data(), | |
1431 | kFALSE); | |
79dc837f | 1432 | } |