]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGCF/EBYE/macros/drawCorrelationFunctionPsi.C
new task for higher moments efficiency and contamination (Nirbhay Kumar Behera <nirbh...
[u/mrichter/AliRoot.git] / PWGCF / EBYE / macros / drawCorrelationFunctionPsi.C
CommitLineData
52daf7b2 1const Int_t numberOfCentralityBins = 8;
2TString centralityArray[numberOfCentralityBins] = {"0-10","10-20","20-30","30-40","40-50","50-60","60-70","70-80"};
a38fd7f3 3
4const Int_t gRebin = 1;
07d0a35c 5
6//____________________________________________________________//
7void drawCorrelationFunctionPsiAllPtCombinations(const char* filename = "AnalysisResults.root",
8 Int_t gCentrality = 1,
9 Int_t gBit = -1,
10 const char* gCentralityEstimator = 0x0,
11 Bool_t kShowShuffled = kFALSE,
12 Bool_t kShowMixed = kTRUE,
13 Double_t psiMin = -0.5,
14 Double_t psiMax = 3.5){
15
16 // this could also be retrieved directly from AliBalancePsi
17 const Int_t kNPtBins = 16;
18 Double_t ptBins[kNPtBins+1] = {0.2,0.6,1.0,1.5,2.0,2.5,3.0,3.5,4.0,5.0,6.0,7.0,8.0,10.,12.,15.,20.};
19
20 cout<<"You have chosen to do all pT combinations --> this could take some time."<<endl;
21
22 for(Int_t iTrig = 0; iTrig < 2/*kNPtBins*/; iTrig++){
23 for(Int_t iAssoc = 0; iAssoc < 2/*kNPtBins*/; iAssoc++){
24 cout<<"================================================================="<<endl;
25 cout<<"PROCESS NOW: "<<endl;
26 cout<<" -> "<< ptBins[iTrig]<<" < pTtrig < "<<ptBins[iTrig+1]<<" "<<ptBins[iAssoc]<<" < pTassoc < "<<ptBins[iAssoc+1]<<endl;
27 drawCorrelationFunctionPsi(filename,gCentrality,gBit,gCentralityEstimator,kShowShuffled,kShowMixed,psiMin,psiMax,ptBins[iTrig],ptBins[iTrig+1],ptBins[iAssoc],ptBins[iAssoc+1]);
28 cout<<"================================================================="<<endl;
29 cout<<endl;
30 }
31 }
32
33}
34
35//____________________________________________________________//
36void drawCorrelationFunctionsAllPtCombinations(const char* lhcPeriod = "LHC11h",
5ea3b761 37 Int_t gTrainID = 208,
07d0a35c 38 Int_t gCentrality = 1,
39 Double_t psiMin = -0.5, Double_t psiMax = 3.5)
40{
41
42 // this could also be retrieved directly from AliBalancePsi
43 const Int_t kNPtBins = 16;
44 Double_t ptBins[kNPtBins+1] = {0.2,0.6,1.0,1.5,2.0,2.5,3.0,3.5,4.0,5.0,6.0,7.0,8.0,10.,12.,15.,20.};
45
46 cout<<"You have chosen to do all pT combinations --> this could take some time."<<endl;
47
60f1c352 48 for(Int_t iTrig = 0; iTrig < kNPtBins; iTrig++){
49 for(Int_t iAssoc = 0; iAssoc < kNPtBins; iAssoc++){
07d0a35c 50 cout<<"================================================================="<<endl;
51 cout<<"FIT NOW: "<<endl;
52 cout<<" -> "<< ptBins[iTrig]<<" < pTtrig < "<<ptBins[iTrig+1]<<" "<<ptBins[iAssoc]<<" < pTassoc < "<<ptBins[iAssoc+1]<<endl;
53 drawCorrelationFunctions(lhcPeriod,gTrainID,gCentrality,psiMin,psiMax,ptBins[iTrig],ptBins[iTrig+1],ptBins[iAssoc],ptBins[iAssoc+1]);
54 cout<<"================================================================="<<endl;
55 cout<<endl;
56 }
57 }
58
59}
60
61
a38fd7f3 62void drawCorrelationFunctionPsi(const char* filename = "AnalysisResults.root",
52daf7b2 63 Int_t gCentrality = 1,
5de9ad1a 64 Int_t gBit = -1,
65 const char* gCentralityEstimator = 0x0,
52daf7b2 66 Bool_t kShowShuffled = kFALSE,
67 Bool_t kShowMixed = kTRUE,
6acdbcb2 68 Double_t psiMin = -0.5,
5de9ad1a 69 Double_t psiMax = 3.5,
6acdbcb2 70 Double_t ptTriggerMin = -1.,
71 Double_t ptTriggerMax = -1.,
72 Double_t ptAssociatedMin = -1.,
73 Double_t ptAssociatedMax = -1.) {
a38fd7f3 74 //Macro that draws the correlation functions from the balance function
75 //analysis vs the reaction plane
76 //Author: Panos.Christakoglou@nikhef.nl
52daf7b2 77 gROOT->LoadMacro("~/SetPlotStyle.C");
78 SetPlotStyle();
79 gStyle->SetPalette(1,0);
80
a38fd7f3 81 //Load the PWG2ebye library
82 gSystem->Load("libANALYSIS.so");
83 gSystem->Load("libANALYSISalice.so");
84 gSystem->Load("libEventMixing.so");
85 gSystem->Load("libCORRFW.so");
86 gSystem->Load("libPWGTools.so");
87 gSystem->Load("libPWGCFebye.so");
88
89 //Prepare the objects and return them
5de9ad1a 90 TList *list = GetListOfObjects(filename,gCentrality,gBit,gCentralityEstimator,0);
52daf7b2 91 TList *listShuffled = NULL;
5de9ad1a 92 if(kShowShuffled) listShuffled = GetListOfObjects(filename,gCentrality,gBit,gCentralityEstimator,1);
52daf7b2 93 TList *listMixed = NULL;
5de9ad1a 94 if(kShowMixed) listMixed = GetListOfObjects(filename,gCentrality,gBit,gCentralityEstimator,2);
52daf7b2 95
a38fd7f3 96 if(!list) {
97 Printf("The TList object was not created");
98 return;
99 }
100 else
52daf7b2 101 draw(list,listShuffled,listMixed,gCentrality,psiMin,psiMax,
6acdbcb2 102 ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
a38fd7f3 103}
104
105//______________________________________________________//
52daf7b2 106TList *GetListOfObjects(const char* filename,
107 Int_t gCentrality,
5de9ad1a 108 Int_t gBit,
109 const char *gCentralityEstimator,
52daf7b2 110 Int_t kData = 1) {
a38fd7f3 111 //Get the TList objects (QA, bf, bf shuffled)
a38fd7f3 112 TList *listBF = 0x0;
a38fd7f3 113
114 //Open the file
ff0805a7 115 TFile *f = TFile::Open(filename,"UPDATE");
a38fd7f3 116 if((!f)||(!f->IsOpen())) {
117 Printf("The file %s is not found. Aborting...",filename);
118 return listBF;
119 }
120 //f->ls();
121
6acdbcb2 122 TDirectoryFile *dir = dynamic_cast<TDirectoryFile *>(f->Get("PWGCFEbyE.outputBalanceFunctionPsiAnalysis"));
a38fd7f3 123 if(!dir) {
124 Printf("The TDirectoryFile is not found. Aborting...",filename);
125 return listBF;
126 }
127 //dir->ls();
128
52daf7b2 129 TString listBFName;
130 if(kData == 0) {
131 //cout<<"no shuffling - no mixing"<<endl;
132 listBFName = "listBFPsi_";
133 }
134 else if(kData == 1) {
135 //cout<<"shuffling - no mixing"<<endl;
136 listBFName = "listBFPsiShuffled_";
137 }
138 else if(kData == 2) {
139 //cout<<"no shuffling - mixing"<<endl;
140 listBFName = "listBFPsiMixed_";
141 }
142 listBFName += centralityArray[gCentrality-1];
5de9ad1a 143 if(gBit > -1) {
144 listBFName += "_Bit"; listBFName += gBit; }
145 if(gCentralityEstimator) {
146 listBFName += "_"; listBFName += gCentralityEstimator;}
a38fd7f3 147
5365d1d7 148 // histograms were already retrieved (in first iteration)
149 if(dir->Get(Form("%s_histograms",listBFName.Data()))){
150 listBF = dynamic_cast<TList *>(dir->Get(Form("%s_histograms",listBFName.Data())));
a38fd7f3 151 }
a38fd7f3 152
5365d1d7 153 // histograms were not yet retrieved (this is the first iteration)
154 else{
155
156 listBF = dynamic_cast<TList *>(dir->Get(listBFName.Data()));
157 cout<<"======================================================="<<endl;
158 cout<<"List name: "<<listBF->GetName()<<endl;
159 //listBF->ls();
a38fd7f3 160
5365d1d7 161 //Get the histograms
162 TString histoName;
163 if(kData == 0)
164 histoName = "fHistPV0M";
165 else if(kData == 1)
166 histoName = "fHistP_shuffleV0M";
167 else if(kData == 2)
168 histoName = "fHistPV0M";
169 AliTHn *fHistP = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
170 if(!fHistP) {
171 Printf("fHistP %s not found!!!",histoName.Data());
172 break;
173 }
174 fHistP->FillParent(); fHistP->DeleteContainers();
175
176 if(kData == 0)
177 histoName = "fHistNV0M";
178 if(kData == 1)
179 histoName = "fHistN_shuffleV0M";
180 if(kData == 2)
181 histoName = "fHistNV0M";
182 AliTHn *fHistN = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
183 if(!fHistN) {
184 Printf("fHistN %s not found!!!",histoName.Data());
185 break;
186 }
187 fHistN->FillParent(); fHistN->DeleteContainers();
188
189 if(kData == 0)
190 histoName = "fHistPNV0M";
191 if(kData == 1)
192 histoName = "fHistPN_shuffleV0M";
193 if(kData == 2)
194 histoName = "fHistPNV0M";
195 AliTHn *fHistPN = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
196 if(!fHistPN) {
197 Printf("fHistPN %s not found!!!",histoName.Data());
198 break;
199 }
200 fHistPN->FillParent(); fHistPN->DeleteContainers();
201
202 if(kData == 0)
203 histoName = "fHistNPV0M";
204 if(kData == 1)
205 histoName = "fHistNP_shuffleV0M";
206 if(kData == 2)
207 histoName = "fHistNPV0M";
208 AliTHn *fHistNP = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
209 if(!fHistNP) {
210 Printf("fHistNP %s not found!!!",histoName.Data());
211 break;
212 }
213 fHistNP->FillParent(); fHistNP->DeleteContainers();
214
215 if(kData == 0)
216 histoName = "fHistPPV0M";
217 if(kData == 1)
218 histoName = "fHistPP_shuffleV0M";
219 if(kData == 2)
220 histoName = "fHistPPV0M";
221 AliTHn *fHistPP = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
222 if(!fHistPP) {
223 Printf("fHistPP %s not found!!!",histoName.Data());
224 break;
225 }
226 fHistPP->FillParent(); fHistPP->DeleteContainers();
227
228 if(kData == 0)
229 histoName = "fHistNNV0M";
230 if(kData == 1)
231 histoName = "fHistNN_shuffleV0M";
232 if(kData == 2)
233 histoName = "fHistNNV0M";
234 AliTHn *fHistNN = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
235 if(!fHistNN) {
236 Printf("fHistNN %s not found!!!",histoName.Data());
237 break;
238 }
239 fHistNN->FillParent(); fHistNN->DeleteContainers();
240
241 dir->cd();
242 listBF->Write(Form("%s_histograms",listBFName.Data()), TObject::kSingleKey);
243
244 }// first iteration
a38fd7f3 245
ff0805a7 246 f->Close();
5365d1d7 247
a38fd7f3 248 return listBF;
249}
250
251//______________________________________________________//
52daf7b2 252void draw(TList *list, TList *listBFShuffled, TList *listBFMixed,
253 Int_t gCentrality, Double_t psiMin, Double_t psiMax,
6acdbcb2 254 Double_t ptTriggerMin, Double_t ptTriggerMax,
255 Double_t ptAssociatedMin, Double_t ptAssociatedMax) {
a38fd7f3 256 //Draws the correlation functions for every centrality bin
52daf7b2 257 //(+-), (-+), (++), (--)
a38fd7f3 258 AliTHn *hP = NULL;
259 AliTHn *hN = NULL;
260 AliTHn *hPN = NULL;
261 AliTHn *hNP = NULL;
262 AliTHn *hPP = NULL;
263 AliTHn *hNN = NULL;
264
265 hP = (AliTHn*) list->FindObject("fHistPV0M");
266 hN = (AliTHn*) list->FindObject("fHistNV0M");
267 hPN = (AliTHn*) list->FindObject("fHistPNV0M");
268 hNP = (AliTHn*) list->FindObject("fHistNPV0M");
269 hPP = (AliTHn*) list->FindObject("fHistPPV0M");
270 hNN = (AliTHn*) list->FindObject("fHistNNV0M");
271
272 //Create the AliBalancePsi object and fill it with the AliTHn objects
273 AliBalancePsi *b = new AliBalancePsi();
274 b->SetHistNp(hP);
275 b->SetHistNn(hN);
276 b->SetHistNpn(hPN);
277 b->SetHistNnp(hNP);
278 b->SetHistNpp(hPP);
279 b->SetHistNnn(hNN);
52daf7b2 280
281 //balance function shuffling
282 AliTHn *hPShuffled = NULL;
283 AliTHn *hNShuffled = NULL;
284 AliTHn *hPNShuffled = NULL;
285 AliTHn *hNPShuffled = NULL;
286 AliTHn *hPPShuffled = NULL;
287 AliTHn *hNNShuffled = NULL;
288 if(listBFShuffled) {
289 //listBFShuffled->ls();
290
291 hPShuffled = (AliTHn*) listBFShuffled->FindObject("fHistP_shuffleV0M");
292 hPShuffled->SetName("gHistPShuffled");
293 hNShuffled = (AliTHn*) listBFShuffled->FindObject("fHistN_shuffleV0M");
294 hNShuffled->SetName("gHistNShuffled");
295 hPNShuffled = (AliTHn*) listBFShuffled->FindObject("fHistPN_shuffleV0M");
296 hPNShuffled->SetName("gHistPNShuffled");
297 hNPShuffled = (AliTHn*) listBFShuffled->FindObject("fHistNP_shuffleV0M");
298 hNPShuffled->SetName("gHistNPShuffled");
299 hPPShuffled = (AliTHn*) listBFShuffled->FindObject("fHistPP_shuffleV0M");
300 hPPShuffled->SetName("gHistPPShuffled");
301 hNNShuffled = (AliTHn*) listBFShuffled->FindObject("fHistNN_shuffleV0M");
302 hNNShuffled->SetName("gHistNNShuffled");
303
304 AliBalancePsi *bShuffled = new AliBalancePsi();
305 bShuffled->SetHistNp(hPShuffled);
306 bShuffled->SetHistNn(hNShuffled);
307 bShuffled->SetHistNpn(hPNShuffled);
308 bShuffled->SetHistNnp(hNPShuffled);
309 bShuffled->SetHistNpp(hPPShuffled);
310 bShuffled->SetHistNnn(hNNShuffled);
311 }
312
313 //balance function mixing
314 AliTHn *hPMixed = NULL;
315 AliTHn *hNMixed = NULL;
316 AliTHn *hPNMixed = NULL;
317 AliTHn *hNPMixed = NULL;
318 AliTHn *hPPMixed = NULL;
319 AliTHn *hNNMixed = NULL;
320
321 if(listBFMixed) {
322 //listBFMixed->ls();
323
324 hPMixed = (AliTHn*) listBFMixed->FindObject("fHistPV0M");
325 hPMixed->SetName("gHistPMixed");
326 hNMixed = (AliTHn*) listBFMixed->FindObject("fHistNV0M");
327 hNMixed->SetName("gHistNMixed");
328 hPNMixed = (AliTHn*) listBFMixed->FindObject("fHistPNV0M");
329 hPNMixed->SetName("gHistPNMixed");
330 hNPMixed = (AliTHn*) listBFMixed->FindObject("fHistNPV0M");
331 hNPMixed->SetName("gHistNPMixed");
332 hPPMixed = (AliTHn*) listBFMixed->FindObject("fHistPPV0M");
333 hPPMixed->SetName("gHistPPMixed");
334 hNNMixed = (AliTHn*) listBFMixed->FindObject("fHistNNV0M");
335 hNNMixed->SetName("gHistNNMixed");
336
337 AliBalancePsi *bMixed = new AliBalancePsi();
338 bMixed->SetHistNp(hPMixed);
339 bMixed->SetHistNn(hNMixed);
340 bMixed->SetHistNpn(hPNMixed);
341 bMixed->SetHistNnp(hNPMixed);
342 bMixed->SetHistNpp(hPPMixed);
343 bMixed->SetHistNnn(hNNMixed);
344 }
345
77b9ff18 346 TH2D *gHistPN[4];
347 TH2D *gHistNP[4];
348 TH2D *gHistPP[4];
349 TH2D *gHistNN[4];
a38fd7f3 350
77b9ff18 351 TCanvas *cPN[4];
352 TCanvas *cNP[4];
353 TCanvas *cPP[4];
354 TCanvas *cNN[4];
a38fd7f3 355 TString histoTitle, pngName;
6acdbcb2 356
52daf7b2 357 //(+-)
358 histoTitle = "(+-) | Centrality: ";
359 histoTitle += centralityArray[gCentrality-1];
6acdbcb2 360 histoTitle += "%";
52daf7b2 361 if((psiMin == -0.5)&&(psiMax == 0.5))
362 histoTitle += " (-7.5^{o} < #varphi - #Psi_{2} < 7.5^{o})";
363 else if((psiMin == 0.5)&&(psiMax == 1.5))
364 histoTitle += " (37.5^{o} < #varphi - #Psi_{2} < 52.5^{o})";
365 else if((psiMin == 1.5)&&(psiMax == 2.5))
366 histoTitle += " (82.5^{o} < #varphi - #Psi_{2} < 97.5^{o})";
367 else
368 histoTitle += " (0^{o} < #varphi - #Psi_{2} < 180^{o})";
369
370 gHistPN[0] = b->GetCorrelationFunctionPN(psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
371 gHistPN[0]->GetYaxis()->SetTitleOffset(1.5);
77b9ff18 372 gHistPN[0]->GetYaxis()->SetTitle("#Delta #varphi (rad)");
52daf7b2 373 gHistPN[0]->SetTitle(histoTitle.Data());
374 cPN[0] = new TCanvas("cPN0","",0,0,600,500);
375 cPN[0]->SetFillColor(10); cPN[0]->SetHighLightColor(10);
5de9ad1a 376 gHistPN[0]->DrawCopy("surf1fb");
377 gPad->SetTheta(30); // default is 30
378 //gPad->SetPhi(130); // default is 30
379 gPad->SetPhi(-60); // default is 30
380 gPad->Update();
6acdbcb2 381 pngName = "DeltaPhiDeltaEta.Centrality";
52daf7b2 382 pngName += centralityArray[gCentrality-1];
6acdbcb2 383 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
384 pngName += ".PositiveNegative.png";
5c9b4733 385 //cPN[0]->SaveAs(pngName.Data());
6acdbcb2 386
52daf7b2 387 if(listBFShuffled) {
388 histoTitle = "(+-) shuffled | Centrality: ";
389 histoTitle += centralityArray[gCentrality-1];
390 histoTitle += "%";
391 if((psiMin == -0.5)&&(psiMax == 0.5))
392 histoTitle += " (-7.5^{o} < #varphi - #Psi_{2} < 7.5^{o})";
393 else if((psiMin == 0.5)&&(psiMax == 1.5))
394 histoTitle += " (37.5^{o} < #varphi - #Psi_{2} < 52.5^{o})";
395 else if((psiMin == 1.5)&&(psiMax == 2.5))
396 histoTitle += " (82.5^{o} < #varphi - #Psi_{2} < 97.5^{o})";
397 else
398 histoTitle += " (0^{o} < #varphi - #Psi_{2} < 180^{o})";
399
400 gHistPN[1] = bShuffled->GetCorrelationFunctionPN(psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
401 gHistPN[1]->GetYaxis()->SetTitleOffset(1.5);
77b9ff18 402 gHistPN[1]->GetYaxis()->SetTitle("#Delta #varphi (rad)");
52daf7b2 403 gHistPN[1]->SetTitle(histoTitle.Data());
404 cPN[1] = new TCanvas("cPN1","",0,100,600,500);
405 cPN[1]->SetFillColor(10);
406 cPN[1]->SetHighLightColor(10);
5de9ad1a 407 gHistPN[1]->DrawCopy("surf1fb");
408 gPad->SetTheta(30); // default is 30
409 //gPad->SetPhi(130); // default is 30
410 gPad->SetPhi(-60); // default is 30
411 gPad->Update();
52daf7b2 412 pngName = "DeltaPhiDeltaEtaShuffled.Centrality";
413 pngName += centralityArray[gCentrality-1];
414 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
415 pngName += ".PositiveNegative.png";
5c9b4733 416 //cPN[1]->SaveAs(pngName.Data());
52daf7b2 417 }
418
419 if(listBFMixed) {
420 histoTitle = "(+-) mixed | Centrality: ";
421 histoTitle += centralityArray[gCentrality-1];
422 histoTitle += "%";
423 if((psiMin == -0.5)&&(psiMax == 0.5))
424 histoTitle += " (-7.5^{o} < #varphi - #Psi_{2} < 7.5^{o})";
425 else if((psiMin == 0.5)&&(psiMax == 1.5))
426 histoTitle += " (37.5^{o} < #varphi - #Psi_{2} < 52.5^{o})";
427 else if((psiMin == 1.5)&&(psiMax == 2.5))
428 histoTitle += " (82.5^{o} < #varphi - #Psi_{2} < 97.5^{o})";
429 else
430 histoTitle += " (0^{o} < #varphi - #Psi_{2} < 180^{o})";
431
432 gHistPN[2] = bMixed->GetCorrelationFunctionPN(psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
433 gHistPN[2]->GetYaxis()->SetTitleOffset(1.5);
77b9ff18 434 gHistPN[2]->GetYaxis()->SetTitle("#Delta #varphi (rad)");
52daf7b2 435 gHistPN[2]->SetTitle(histoTitle.Data());
436 cPN[2] = new TCanvas("cPN2","",0,200,600,500);
437 cPN[2]->SetFillColor(10);
438 cPN[2]->SetHighLightColor(10);
5de9ad1a 439 gHistPN[2]->DrawCopy("surf1fb");
440 gPad->SetTheta(30); // default is 30
441 //gPad->SetPhi(130); // default is 30
442 gPad->SetPhi(-60); // default is 30
443 gPad->Update();
52daf7b2 444 pngName = "DeltaPhiDeltaEtaMixed.Centrality";
445 pngName += centralityArray[gCentrality-1];
446 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
447 pngName += ".PositiveNegative.png";
5c9b4733 448 //cPN[2]->SaveAs(pngName.Data());
77b9ff18 449
450 //Correlation function (+-)
451 gHistPN[3] = dynamic_cast<TH2D *>(gHistPN[0]->Clone());
452 gHistPN[3]->Divide(gHistPN[2]);
453 gHistPN[3]->GetXaxis()->SetRangeUser(-1.5,1.5);
454 gHistPN[3]->GetZaxis()->SetTitle("C_{+-}(#Delta#eta,#Delta#varphi)");
455 cPN[3] = new TCanvas("cPN3","",0,300,600,500);
456 cPN[3]->SetFillColor(10);
457 cPN[3]->SetHighLightColor(10);
458 gHistPN[3]->DrawCopy("surf1fb");
459 gPad->SetTheta(30); // default is 30
460 //gPad->SetPhi(130); // default is 30
461 gPad->SetPhi(-60); // default is 30
462 gPad->Update();
463 pngName = "CorrelationFunction.Centrality";
464 pngName += centralityArray[gCentrality-1];
465 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
466 pngName += ".PositiveNegative.png";
5c9b4733 467 //cPN[3]->SaveAs(pngName.Data());
52daf7b2 468 }
469
470 //(-+)
6acdbcb2 471 histoTitle = "(-+) | Centrality: ";
52daf7b2 472 histoTitle += centralityArray[gCentrality-1];
6acdbcb2 473 histoTitle += "%";
52daf7b2 474 if((psiMin == -0.5)&&(psiMax == 0.5))
475 histoTitle += " (-7.5^{o} < #varphi - #Psi_{2} < 7.5^{o})";
476 else if((psiMin == 0.5)&&(psiMax == 1.5))
477 histoTitle += " (37.5^{o} < #varphi - #Psi_{2} < 52.5^{o})";
478 else if((psiMin == 1.5)&&(psiMax == 2.5))
479 histoTitle += " (82.5^{o} < #varphi - #Psi_{2} < 97.5^{o})";
480 else
481 histoTitle += " (0^{o} < #varphi - #Psi_{2} < 180^{o})";
482
483 gHistNP[0] = b->GetCorrelationFunctionNP(psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
484 gHistNP[0]->GetYaxis()->SetTitleOffset(1.5);
77b9ff18 485 gHistNP[0]->GetYaxis()->SetTitle("#Delta #varphi (rad)");
52daf7b2 486 gHistNP[0]->SetTitle(histoTitle.Data());
487 cNP[0] = new TCanvas("cNP0","",100,0,600,500);
488 cNP[0]->SetFillColor(10);
489 cNP[0]->SetHighLightColor(10);
5de9ad1a 490 gHistNP[0]->DrawCopy("surf1fb");
491 gPad->SetTheta(30); // default is 30
492 //gPad->SetPhi(130); // default is 30
493 gPad->SetPhi(-60); // default is 30
494 gPad->Update();
6acdbcb2 495 pngName = "DeltaPhiDeltaEta.Centrality";
52daf7b2 496 pngName += centralityArray[gCentrality-1];
6acdbcb2 497 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
498 pngName += ".NegativePositive.png";
5c9b4733 499 //cNP[0]->SaveAs(pngName.Data());
52daf7b2 500
501 if(listBFShuffled) {
502 histoTitle = "(-+) shuffled | Centrality: ";
503 histoTitle += centralityArray[gCentrality-1];
504 histoTitle += "%";
505 if((psiMin == -0.5)&&(psiMax == 0.5))
506 histoTitle += " (-7.5^{o} < #varphi - #Psi_{2} < 7.5^{o})";
507 else if((psiMin == 0.5)&&(psiMax == 1.5))
508 histoTitle += " (37.5^{o} < #varphi - #Psi_{2} < 52.5^{o})";
509 else if((psiMin == 1.5)&&(psiMax == 2.5))
510 histoTitle += " (82.5^{o} < #varphi - #Psi_{2} < 97.5^{o})";
511 else
512 histoTitle += " (0^{o} < #varphi - #Psi_{2} < 180^{o})";
513
514 gHistNP[1] = bShuffled->GetCorrelationFunctionNP(psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
515 gHistNP[1]->GetYaxis()->SetTitleOffset(1.5);
77b9ff18 516 gHistNP[1]->GetYaxis()->SetTitle("#Delta #varphi (rad)");
52daf7b2 517 gHistNP[1]->SetTitle(histoTitle.Data());
518 cNP[1] = new TCanvas("cNP1","",100,100,600,500);
519 cNP[1]->SetFillColor(10);
520 cNP[1]->SetHighLightColor(10);
5de9ad1a 521 gHistNP[1]->DrawCopy("surf1fb");
522 gPad->SetTheta(30); // default is 30
523 //gPad->SetPhi(130); // default is 30
524 gPad->SetPhi(-60); // default is 30
525 gPad->Update();
52daf7b2 526 pngName = "DeltaPhiDeltaEtaShuffled.Centrality";
527 pngName += centralityArray[gCentrality-1];
528 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
529 pngName += ".NegativePositive.png";
5c9b4733 530 //cNP[1]->SaveAs(pngName.Data());
52daf7b2 531 }
532
533 if(listBFMixed) {
534 histoTitle = "(-+) mixed | Centrality: ";
535 histoTitle += centralityArray[gCentrality-1];
536 histoTitle += "%";
537 if((psiMin == -0.5)&&(psiMax == 0.5))
538 histoTitle += " (-7.5^{o} < #varphi - #Psi_{2} < 7.5^{o})";
539 else if((psiMin == 0.5)&&(psiMax == 1.5))
540 histoTitle += " (37.5^{o} < #varphi - #Psi_{2} < 52.5^{o})";
541 else if((psiMin == 1.5)&&(psiMax == 2.5))
542 histoTitle += " (82.5^{o} < #varphi - #Psi_{2} < 97.5^{o})";
543 else
544 histoTitle += " (0^{o} < #varphi - #Psi_{2} < 180^{o})";
545
546 gHistNP[2] = bMixed->GetCorrelationFunctionNP(psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
547 gHistNP[2]->GetYaxis()->SetTitleOffset(1.5);
77b9ff18 548 gHistNP[2]->GetYaxis()->SetTitle("#Delta #varphi (rad)");
52daf7b2 549 gHistNP[2]->SetTitle(histoTitle.Data());
550 cNP[2] = new TCanvas("cNP2","",100,200,600,500);
551 cNP[2]->SetFillColor(10);
552 cNP[2]->SetHighLightColor(10);
5de9ad1a 553 gHistNP[2]->DrawCopy("surf1fb");
554 gPad->SetTheta(30); // default is 30
555 //gPad->SetPhi(130); // default is 30
556 gPad->SetPhi(-60); // default is 30
557 gPad->Update();
52daf7b2 558 pngName = "DeltaPhiDeltaEtaMixed.Centrality";
559 pngName += centralityArray[gCentrality-1];
560 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
561 pngName += ".NegativePositive.png";
5c9b4733 562 //cNP[2]->SaveAs(pngName.Data());
77b9ff18 563
564 //Correlation function (-+)
565 gHistNP[3] = dynamic_cast<TH2D *>(gHistNP[0]->Clone());
566 gHistNP[3]->Divide(gHistNP[2]);
567 gHistNP[3]->GetXaxis()->SetRangeUser(-1.5,1.5);
568 gHistNP[3]->GetZaxis()->SetTitle("C_{-+}(#Delta#eta,#Delta#varphi)");
569 cNP[3] = new TCanvas("cNP3","",100,300,600,500);
570 cNP[3]->SetFillColor(10);
571 cNP[3]->SetHighLightColor(10);
572 gHistNP[3]->DrawCopy("surf1fb");
573 gPad->SetTheta(30); // default is 30
574 //gPad->SetPhi(130); // default is 30
575 gPad->SetPhi(-60); // default is 30
576 gPad->Update();
577 pngName = "CorrelationFunction.Centrality";
578 pngName += centralityArray[gCentrality-1];
579 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
580 pngName += ".NegativePositive.png";
5c9b4733 581 //cNP[3]->SaveAs(pngName.Data());
52daf7b2 582 }
6acdbcb2 583
52daf7b2 584 //(++)
6acdbcb2 585 histoTitle = "(++) | Centrality: ";
52daf7b2 586 histoTitle += centralityArray[gCentrality-1];
6acdbcb2 587 histoTitle += "%";
52daf7b2 588 if((psiMin == -0.5)&&(psiMax == 0.5))
589 histoTitle += " (-7.5^{o} < #varphi - #Psi_{2} < 7.5^{o})";
590 else if((psiMin == 0.5)&&(psiMax == 1.5))
591 histoTitle += " (37.5^{o} < #varphi - #Psi_{2} < 52.5^{o})";
592 else if((psiMin == 1.5)&&(psiMax == 2.5))
593 histoTitle += " (82.5^{o} < #varphi - #Psi_{2} < 97.5^{o})";
594 else
595 histoTitle += " (0^{o} < #varphi - #Psi_{2} < 180^{o})";
596
597 gHistPP[0] = b->GetCorrelationFunctionPP(psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
598 gHistPP[0]->GetYaxis()->SetTitleOffset(1.5);
77b9ff18 599 gHistPP[0]->GetYaxis()->SetTitle("#Delta #varphi (rad)");
52daf7b2 600 gHistPP[0]->SetTitle(histoTitle.Data());
601 cPP[0] = new TCanvas("cPP0","",200,0,600,500);
602 cPP[0]->SetFillColor(10);
603 cPP[0]->SetHighLightColor(10);
5de9ad1a 604 gHistPP[0]->DrawCopy("surf1fb");
605 gPad->SetTheta(30); // default is 30
606 //gPad->SetPhi(130); // default is 30
607 gPad->SetPhi(-60); // default is 30
608 gPad->Update();
6acdbcb2 609 pngName = "DeltaPhiDeltaEta.Centrality";
52daf7b2 610 pngName += centralityArray[gCentrality-1];
6acdbcb2 611 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
612 pngName += ".PositivePositive.png";
5c9b4733 613 //cPP[0]->SaveAs(pngName.Data());
6acdbcb2 614
52daf7b2 615 if(listBFShuffled) {
616 histoTitle = "(++) shuffled | Centrality: ";
617 histoTitle += centralityArray[gCentrality-1];
618 histoTitle += "%";
619 if((psiMin == -0.5)&&(psiMax == 0.5))
620 histoTitle += " (-7.5^{o} < #varphi - #Psi_{2} < 7.5^{o})";
621 else if((psiMin == 0.5)&&(psiMax == 1.5))
622 histoTitle += " (37.5^{o} < #varphi - #Psi_{2} < 52.5^{o})";
623 else if((psiMin == 1.5)&&(psiMax == 2.5))
624 histoTitle += " (82.5^{o} < #varphi - #Psi_{2} < 97.5^{o})";
625 else
626 histoTitle += " (0^{o} < #varphi - #Psi_{2} < 180^{o})";
627
628 gHistPP[1] = bShuffled->GetCorrelationFunctionPP(psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
629 gHistPP[1]->GetYaxis()->SetTitleOffset(1.5);
77b9ff18 630 gHistPP[1]->GetYaxis()->SetTitle("#Delta #varphi (rad)");
52daf7b2 631 gHistPP[1]->SetTitle(histoTitle.Data());
632 cPP[1] = new TCanvas("cPP1","",200,100,600,500);
633 cPP[1]->SetFillColor(10);
634 cPP[1]->SetHighLightColor(10);
5de9ad1a 635 gHistPP[1]->DrawCopy("surf1fb");
636 gPad->SetTheta(30); // default is 30
637 //gPad->SetPhi(130); // default is 30
638 gPad->SetPhi(-60); // default is 30
639 gPad->Update();
52daf7b2 640 pngName = "DeltaPhiDeltaEtaShuffled.Centrality";
641 pngName += centralityArray[gCentrality-1];
642 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
643 pngName += ".PositivePositive.png";
5c9b4733 644 //cPP[1]->SaveAs(pngName.Data());
52daf7b2 645 }
646
647 if(listBFMixed) {
648 histoTitle = "(++) mixed | Centrality: ";
649 histoTitle += centralityArray[gCentrality-1];
650 histoTitle += "%";
651 if((psiMin == -0.5)&&(psiMax == 0.5))
652 histoTitle += " (-7.5^{o} < #varphi - #Psi_{2} < 7.5^{o})";
653 else if((psiMin == 0.5)&&(psiMax == 1.5))
654 histoTitle += " (37.5^{o} < #varphi - #Psi_{2} < 52.5^{o})";
655 else if((psiMin == 1.5)&&(psiMax == 2.5))
656 histoTitle += " (82.5^{o} < #varphi - #Psi_{2} < 97.5^{o})";
657 else
658 histoTitle += " (0^{o} < #varphi - #Psi_{2} < 180^{o})";
659
660 gHistPP[2] = bMixed->GetCorrelationFunctionPP(psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
661 gHistPP[2]->GetYaxis()->SetTitleOffset(1.5);
77b9ff18 662 gHistPP[2]->GetYaxis()->SetTitle("#Delta #varphi (rad)");
52daf7b2 663 gHistPP[2]->SetTitle(histoTitle.Data());
664 cPP[2] = new TCanvas("cPP2","",200,200,600,500);
665 cPP[2]->SetFillColor(10);
666 cPP[2]->SetHighLightColor(10);
5de9ad1a 667 gHistPP[2]->DrawCopy("surf1fb");
668 gPad->SetTheta(30); // default is 30
669 //gPad->SetPhi(130); // default is 30
670 gPad->SetPhi(-60); // default is 30
671 gPad->Update();
52daf7b2 672 pngName = "DeltaPhiDeltaEtaMixed.Centrality";
673 pngName += centralityArray[gCentrality-1];
674 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
675 pngName += ".PositivePositive.png";
5c9b4733 676 //cPP[2]->SaveAs(pngName.Data());
77b9ff18 677
678 //Correlation function (++)
679 gHistPP[3] = dynamic_cast<TH2D *>(gHistPP[0]->Clone());
680 gHistPP[3]->Divide(gHistPP[2]);
681 gHistPP[3]->GetXaxis()->SetRangeUser(-1.5,1.5);
682 gHistPP[3]->GetZaxis()->SetTitle("C_{++}(#Delta#eta,#Delta#varphi)");
683 cPP[3] = new TCanvas("cPP3","",200,300,600,500);
684 cPP[3]->SetFillColor(10);
685 cPP[3]->SetHighLightColor(10);
686 gHistPP[3]->DrawCopy("surf1fb");
687 gPad->SetTheta(30); // default is 30
688 //gPad->SetPhi(130); // default is 30
689 gPad->SetPhi(-60); // default is 30
690 gPad->Update();
691 pngName = "CorrelationFunction.Centrality";
692 pngName += centralityArray[gCentrality-1];
693 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
694 pngName += ".PositivePositive.png";
5c9b4733 695 //cPP[3]->SaveAs(pngName.Data());
52daf7b2 696 }
697
698 //(--)
6acdbcb2 699 histoTitle = "(--) | Centrality: ";
52daf7b2 700 histoTitle += centralityArray[gCentrality-1];
6acdbcb2 701 histoTitle += "%";
52daf7b2 702 if((psiMin == -0.5)&&(psiMax == 0.5))
703 histoTitle += " (-7.5^{o} < #varphi - #Psi_{2} < 7.5^{o})";
704 else if((psiMin == 0.5)&&(psiMax == 1.5))
705 histoTitle += " (37.5^{o} < #varphi - #Psi_{2} < 52.5^{o})";
706 else if((psiMin == 1.5)&&(psiMax == 2.5))
707 histoTitle += " (82.5^{o} < #varphi - #Psi_{2} < 97.5^{o})";
708 else
709 histoTitle += " (0^{o} < #varphi - #Psi_{2} < 180^{o})";
710
711 gHistNN[0] = b->GetCorrelationFunctionNN(psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
712 gHistNN[0]->GetYaxis()->SetTitleOffset(1.5);
77b9ff18 713 gHistNN[0]->GetYaxis()->SetTitle("#Delta #varphi (rad)");
52daf7b2 714 gHistNN[0]->SetTitle(histoTitle.Data());
715 cNN[0] = new TCanvas("cNN0","",300,0,600,500);
716 cNN[0]->SetFillColor(10);
717 cNN[0]->SetHighLightColor(10);
5de9ad1a 718 gHistNN[0]->DrawCopy("surf1fb");
719 gPad->SetTheta(30); // default is 30
720 gPad->SetPhi(-60); // default is 30
721 //gPad->SetPhi(-60); // default is 30
722 gPad->Update();
6acdbcb2 723 pngName = "DeltaPhiDeltaEta.Centrality";
52daf7b2 724 pngName += centralityArray[gCentrality-1];
6acdbcb2 725 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
726 pngName += ".NegativeNegative.png";
5c9b4733 727 //cNN[0]->SaveAs(pngName.Data());
52daf7b2 728
729 if(listBFShuffled) {
730 histoTitle = "(--) shuffled | Centrality: ";
731 histoTitle += centralityArray[gCentrality-1];
732 histoTitle += "%";
733 if((psiMin == -0.5)&&(psiMax == 0.5))
734 histoTitle += " (-7.5^{o} < #varphi - #Psi_{2} < 7.5^{o})";
735 else if((psiMin == 0.5)&&(psiMax == 1.5))
736 histoTitle += " (37.5^{o} < #varphi - #Psi_{2} < 52.5^{o})";
737 else if((psiMin == 1.5)&&(psiMax == 2.5))
738 histoTitle += " (82.5^{o} < #varphi - #Psi_{2} < 97.5^{o})";
739 else
740 histoTitle += " (0^{o} < #varphi - #Psi_{2} < 180^{o})";
741
742 gHistNN[1] = bShuffled->GetCorrelationFunctionNN(psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
743 gHistNN[1]->GetYaxis()->SetTitleOffset(1.5);
77b9ff18 744 gHistNN[1]->GetYaxis()->SetTitle("#Delta #varphi (rad)");
52daf7b2 745 gHistNN[1]->SetTitle(histoTitle.Data());
746 cNN[1] = new TCanvas("cNN1","",300,100,600,500);
747 cNN[1]->SetFillColor(10);
748 cNN[1]->SetHighLightColor(10);
5de9ad1a 749 gHistNN[1]->DrawCopy("surf1fb");
750 gPad->SetTheta(30); // default is 30
751 //gPad->SetPhi(130); // default is 30
752 gPad->SetPhi(-60); // default is 30
753 gPad->Update();
52daf7b2 754 pngName = "DeltaPhiDeltaEtaShuffled.Centrality";
755 pngName += centralityArray[gCentrality-1];
756 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
757 pngName += ".NegativeNegative.png";
5c9b4733 758 //cNN[1]->SaveAs(pngName.Data());
52daf7b2 759 }
760
761 if(listBFMixed) {
762 histoTitle = "(--) mixed | Centrality: ";
763 histoTitle += centralityArray[gCentrality-1];
764 histoTitle += "%";
765 if((psiMin == -0.5)&&(psiMax == 0.5))
766 histoTitle += " (-7.5^{o} < #varphi - #Psi_{2} < 7.5^{o})";
767 else if((psiMin == 0.5)&&(psiMax == 1.5))
768 histoTitle += " (37.5^{o} < #varphi - #Psi_{2} < 52.5^{o})";
769 else if((psiMin == 1.5)&&(psiMax == 2.5))
770 histoTitle += " (82.5^{o} < #varphi - #Psi_{2} < 97.5^{o})";
771 else
772 histoTitle += " (0^{o} < #varphi - #Psi_{2} < 180^{o})";
773
774 gHistNN[2] = bMixed->GetCorrelationFunctionNN(psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
775 gHistNN[2]->GetYaxis()->SetTitleOffset(1.5);
77b9ff18 776 gHistNN[2]->GetYaxis()->SetTitle("#Delta #varphi (rad)");
52daf7b2 777 gHistNN[2]->SetTitle(histoTitle.Data());
778 cNN[2] = new TCanvas("cNN2","",300,200,600,500);
779 cNN[2]->SetFillColor(10);
780 cNN[2]->SetHighLightColor(10);
5de9ad1a 781 gHistNN[2]->DrawCopy("surf1fb");
782 gPad->SetTheta(30); // default is 30
783 //gPad->SetPhi(130); // default is 30
784 gPad->SetPhi(-60); // default is 30
785 gPad->Update();
52daf7b2 786 pngName = "DeltaPhiDeltaEtaMixed.Centrality";
787 pngName += centralityArray[gCentrality-1];
788 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
789 pngName += ".NegativeNegative.png";
5c9b4733 790 //cNN[2]->SaveAs(pngName.Data());
77b9ff18 791
792 //Correlation function (--)
793 gHistNN[3] = dynamic_cast<TH2D *>(gHistNN[0]->Clone());
794 gHistNN[3]->Divide(gHistNN[2]);
795 gHistNN[3]->GetXaxis()->SetRangeUser(-1.5,1.5);
796 gHistNN[3]->GetZaxis()->SetTitle("C_{--}(#Delta#eta,#Delta#varphi)");
797 cNN[3] = new TCanvas("cNN3","",300,300,600,500);
798 cNN[3]->SetFillColor(10);
799 cNN[3]->SetHighLightColor(10);
800 gHistNN[3]->DrawCopy("surf1fb");
801 gPad->SetTheta(30); // default is 30
802 //gPad->SetPhi(130); // default is 30
803 gPad->SetPhi(-60); // default is 30
804 gPad->Update();
805 pngName = "CorrelationFunction.Centrality";
806 pngName += centralityArray[gCentrality-1];
807 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
808 pngName += ".NegativeNegative.png";
5c9b4733 809 //cNN[3]->SaveAs(pngName.Data());
52daf7b2 810 }
5c9b4733 811
812 //Write to output file
813 TString newFileName = "correlationFunction.Centrality";
814 newFileName += gCentrality; newFileName += ".Psi";
815 if((psiMin == -0.5)&&(psiMax == 0.5)) newFileName += "InPlane.Ptt";
816 else if((psiMin == 0.5)&&(psiMax == 1.5)) newFileName += "Intermediate.Ptt";
817 else if((psiMin == 1.5)&&(psiMax == 2.5)) newFileName += "OutOfPlane.Ptt";
818 else if((psiMin == 2.5)&&(psiMax == 3.5)) newFileName += "Rest.PttFrom";
819 else newFileName += "All.PttFrom";
820 newFileName += Form("%.1f",ptTriggerMin); newFileName += "To";
821 newFileName += Form("%.1f",ptTriggerMax); newFileName += "PtaFrom";
822 newFileName += Form("%.1f",ptAssociatedMin); newFileName += "To";
823 newFileName += Form("%.1f",ptAssociatedMax);
824 newFileName += ".root";
825 TFile *newFile = TFile::Open(newFileName.Data(),"recreate");
826 gHistPN[0]->SetName("gHistPNRaw"); gHistPN[0]->Write();
827 gHistNP[0]->SetName("gHistNPRaw"); gHistNP[0]->Write();
828 gHistPP[0]->SetName("gHistPPRaw"); gHistPP[0]->Write();
829 gHistNN[0]->SetName("gHistNNRaw"); gHistNN[0]->Write();
830 if(listBFShuffled) {
831 gHistPN[1]->SetName("gHistPNShuffled"); gHistPN[1]->Write();
832 gHistNP[1]->SetName("gHistNPShuffled"); gHistNP[1]->Write();
833 gHistPP[1]->SetName("gHistPPShuffled"); gHistPP[1]->Write();
834 gHistNN[1]->SetName("gHistNNShuffled"); gHistNN[1]->Write();
835 }
836 if(listBFMixed) {
837 gHistPN[2]->SetName("gHistPNMixed"); gHistPN[2]->Write();
838 gHistNP[2]->SetName("gHistNPMixed"); gHistNP[2]->Write();
839 gHistPP[2]->SetName("gHistPPMixed"); gHistPP[2]->Write();
840 gHistNN[2]->SetName("gHistNNMixed"); gHistNN[2]->Write();
841
842 gHistPN[3]->SetName("gHistPNCorrelationFunctions"); gHistPN[3]->Write();
843 gHistNP[3]->SetName("gHistNPCorrelationFunctions"); gHistNP[3]->Write();
844 gHistPP[3]->SetName("gHistPPCorrelationFunctions"); gHistPP[3]->Write();
845 gHistNN[3]->SetName("gHistNNCorrelationFunctions"); gHistNN[3]->Write();
846 }
847 newFile->Close();
07d0a35c 848
849 // some cleaning
850 for(Int_t i = 0; i < 4; i++){
851
852 if(!listBFShuffled && i == 1) continue;
853 if(!listBFMixed && (i == 2 || i == 3)) continue;
854
855 if(gHistPP[i]) delete gHistPP[i];
856 if(gHistPN[i]) delete gHistPN[i];
857 if(gHistNP[i]) delete gHistNP[i];
858 if(gHistNN[i]) delete gHistNN[i];
859
860 if(cPN[i]) delete cPN[i];
861 if(cNP[i]) delete cNP[i];
862 if(cPP[i]) delete cPP[i];
863 if(cNN[i]) delete cNN[i];
864 }
865
866 delete hP;
867 delete hN;
868 delete hPP;
869 delete hPN;
870 delete hNP;
871 delete hNN;
872
873 delete hPMixed;
874 delete hNMixed;
875 delete hPPMixed;
876 delete hPNMixed;
877 delete hNPMixed;
878 delete hNNMixed;
879
880 delete hPShuffled;
881 delete hNShuffled;
882 delete hPPShuffled;
883 delete hPNShuffled;
884 delete hNPShuffled;
885 delete hNNShuffled;
886
5c9b4733 887}
888
889//____________________________________________________________//
890void drawCorrelationFunctions(const char* lhcPeriod = "LHC11h",
5ea3b761 891 Int_t gTrainID = 208,
5c9b4733 892 Int_t gCentrality = 1,
893 Double_t psiMin = -0.5, Double_t psiMax = 3.5,
894 Double_t ptTriggerMin = -1.,
895 Double_t ptTriggerMax = -1.,
896 Double_t ptAssociatedMin = -1.,
897 Double_t ptAssociatedMax = -1.) {
898 //Macro that draws the charge dependent correlation functions
899 //for each centrality bin for the different pT of trigger and
900 //associated particles
901 //Author: Panos.Christakoglou@nikhef.nl
902 TGaxis::SetMaxDigits(3);
903
904 //Get the input file
905 TString filename = "PbPb/"; filename += lhcPeriod;
07d0a35c 906 filename +="/Train"; filename += gTrainID;
907 filename +="/Centrality"; filename += gCentrality;
5c9b4733 908 filename += "/correlationFunction.Centrality";
909 filename += gCentrality; filename += ".Psi";
910 if((psiMin == -0.5)&&(psiMax == 0.5)) filename += "InPlane.Ptt";
911 else if((psiMin == 0.5)&&(psiMax == 1.5)) filename += "Intermediate.Ptt";
912 else if((psiMin == 1.5)&&(psiMax == 2.5)) filename += "OutOfPlane.Ptt";
913 else if((psiMin == 2.5)&&(psiMax == 3.5)) filename += "Rest.Ptt";
914 else filename += "All.PttFrom";
915 filename += Form("%.1f",ptTriggerMin); filename += "To";
916 filename += Form("%.1f",ptTriggerMax); filename += "PtaFrom";
917 filename += Form("%.1f",ptAssociatedMin); filename += "To";
918 filename += Form("%.1f",ptAssociatedMax);
919 filename += ".root";
920
921 //Open the file
922 TFile *f = TFile::Open(filename.Data());
923 if((!f)||(!f->IsOpen())) {
924 Printf("The file %s is not found. Aborting...",filename);
925 return listBF;
926 }
927 //f->ls();
928
929 //Latex
930 TString centralityLatex = "Centrality: ";
931 centralityLatex += centralityArray[gCentrality-1];
932 centralityLatex += "%";
933
934 TString psiLatex;
935 if((psiMin == -0.5)&&(psiMax == 0.5))
936 psiLatex = " -7.5^{o} < #varphi - #Psi_{2} < 7.5^{o}";
937 else if((psiMin == 0.5)&&(psiMax == 1.5))
938 psiLatex = " 37.5^{o} < #varphi - #Psi_{2} < 52.5^{o}";
939 else if((psiMin == 1.5)&&(psiMax == 2.5))
940 psiLatex = " 82.5^{o} < #varphi - #Psi_{2} < 97.5^{o}";
941 else
942 psiLatex = " 0^{o} < #varphi - #Psi_{2} < 180^{o}";
943
944 TString pttLatex = Form("%.1f",ptTriggerMin);
945 pttLatex += " < p_{T,trig} < "; pttLatex += Form("%.1f",ptTriggerMax);
946 pttLatex += " GeV/c";
947
948 TString ptaLatex = Form("%.1f",ptAssociatedMin);
949 ptaLatex += " < p_{T,assoc} < "; ptaLatex += Form("%.1f",ptAssociatedMax);
950 ptaLatex += " GeV/c";
951
952 TLatex *latexInfo1 = new TLatex();
953 latexInfo1->SetNDC();
954 latexInfo1->SetTextSize(0.045);
955 latexInfo1->SetTextColor(1);
956
957 TString pngName;
958
959 //============================================================//
960 //Get the +- correlation function
961 TH2D *gHistPN = dynamic_cast<TH2D *>(f->Get("gHistPNCorrelationFunctions"));
962 gHistPN->SetStats(kFALSE);
963 gHistPN->SetTitle("");
964 gHistPN->GetXaxis()->SetRangeUser(-1.45,1.45);
965 gHistPN->GetXaxis()->CenterTitle();
966 gHistPN->GetXaxis()->SetTitleOffset(1.2);
967 gHistPN->GetYaxis()->CenterTitle();
968 gHistPN->GetYaxis()->SetTitleOffset(1.2);
969 gHistPN->GetZaxis()->SetTitleOffset(1.2);
970 TCanvas *cPN = new TCanvas("cPN","",0,0,600,500);
971 cPN->SetFillColor(10); cPN->SetHighLightColor(10);
972 cPN->SetLeftMargin(0.15);
973 gHistPN->DrawCopy("surf1fb");
974 gPad->SetTheta(30); // default is 30
975 gPad->SetPhi(-60); // default is 30
976 gPad->Update();
977
978 latexInfo1->DrawLatex(0.44,0.88,centralityLatex.Data());
979 //latexInfo1->DrawLatex(0.44,0.82,psiLatex.Data());
980 latexInfo1->DrawLatex(0.44,0.82,pttLatex.Data());
981 latexInfo1->DrawLatex(0.44,0.76,ptaLatex.Data());
982
983 pngName = "CorrelationFunction.Centrality";
984 pngName += centralityArray[gCentrality-1];
985 pngName += ".Psi";
986 if((psiMin == -0.5)&&(psiMax == 0.5)) pngName += "InPlane.Ptt";
987 else if((psiMin == 0.5)&&(psiMax == 1.5)) pngName += "Intermediate.Ptt";
988 else if((psiMin == 1.5)&&(psiMax == 2.5)) pngName += "OutOfPlane.Ptt";
989 else if((psiMin == 2.5)&&(psiMax == 3.5)) pngName += "Rest.Ptt";
990 else pngName += "All.PttFrom";
991 pngName += Form("%.1f",ptTriggerMin); pngName += "To";
992 pngName += Form("%.1f",ptTriggerMax); pngName += "PtaFrom";
993 pngName += Form("%.1f",ptAssociatedMin); pngName += "To";
994 pngName += Form("%.1f",ptAssociatedMax);
995 pngName += ".PositiveNegative.png";
996 cPN->SaveAs(pngName.Data());
997 fitCorrelationFunctions(gCentrality, psiMin, psiMax,
998 ptTriggerMin,ptTriggerMax,
999 ptAssociatedMin, ptAssociatedMax,gHistPN);
1000 //============================================================//
1001 //Get the -+ correlation function
1002 TH2D *gHistNP = dynamic_cast<TH2D *>(f->Get("gHistNPCorrelationFunctions"));
1003 gHistNP->SetStats(kFALSE);
1004 gHistNP->SetTitle("");
1005 gHistNP->GetXaxis()->SetRangeUser(-1.45,1.45);
1006 gHistNP->GetXaxis()->CenterTitle();
1007 gHistNP->GetXaxis()->SetTitleOffset(1.2);
1008 gHistNP->GetYaxis()->CenterTitle();
1009 gHistNP->GetYaxis()->SetTitleOffset(1.2);
1010 gHistNP->GetZaxis()->SetTitleOffset(1.2);
1011 TCanvas *cNP = new TCanvas("cNP","",50,50,600,500);
1012 cNP->SetFillColor(10); cNP->SetHighLightColor(10);
1013 cNP->SetLeftMargin(0.15);
1014 gHistNP->DrawCopy("surf1fb");
1015 gPad->SetTheta(30); // default is 30
1016 gPad->SetPhi(-60); // default is 30
1017 gPad->Update();
1018
1019 latexInfo1->DrawLatex(0.44,0.88,centralityLatex.Data());
1020 //latexInfo1->DrawLatex(0.44,0.82,psiLatex.Data());
1021 latexInfo1->DrawLatex(0.44,0.82,pttLatex.Data());
1022 latexInfo1->DrawLatex(0.44,0.76,ptaLatex.Data());
1023
1024 pngName = "CorrelationFunction.Centrality";
1025 pngName += centralityArray[gCentrality-1];
1026 pngName += ".Psi";
1027 if((psiMin == -0.5)&&(psiMax == 0.5)) pngName += "InPlane.Ptt";
1028 else if((psiMin == 0.5)&&(psiMax == 1.5)) pngName += "Intermediate.Ptt";
1029 else if((psiMin == 1.5)&&(psiMax == 2.5)) pngName += "OutOfPlane.Ptt";
1030 else if((psiMin == 2.5)&&(psiMax == 3.5)) pngName += "Rest.Ptt";
1031 else pngName += "All.PttFrom";
1032 pngName += Form("%.1f",ptTriggerMin); pngName += "To";
1033 pngName += Form("%.1f",ptTriggerMax); pngName += "PtaFrom";
1034 pngName += Form("%.1f",ptAssociatedMin); pngName += "To";
1035 pngName += Form("%.1f",ptAssociatedMax);
1036 pngName += ".NegativePositive.png";
1037 cNP->SaveAs(pngName.Data());
1038
1039 //============================================================//
1040 //Get the ++ correlation function
1041 TH2D *gHistPP = dynamic_cast<TH2D *>(f->Get("gHistPPCorrelationFunctions"));
1042 gHistPP->SetStats(kFALSE);
1043 gHistPP->SetTitle("");
1044 gHistPP->GetXaxis()->SetRangeUser(-1.45,1.45);
1045 gHistPP->GetXaxis()->CenterTitle();
1046 gHistPP->GetXaxis()->SetTitleOffset(1.2);
1047 gHistPP->GetYaxis()->CenterTitle();
1048 gHistPP->GetYaxis()->SetTitleOffset(1.2);
1049 gHistPP->GetZaxis()->SetTitleOffset(1.2);
1050 TCanvas *cPP = new TCanvas("cPP","",100,100,600,500);
1051 cPP->SetFillColor(10); cPP->SetHighLightColor(10);
1052 cPP->SetLeftMargin(0.15);
1053 gHistPP->DrawCopy("surf1fb");
1054 gPad->SetTheta(30); // default is 30
1055 gPad->SetPhi(-60); // default is 30
1056 gPad->Update();
1057
1058 latexInfo1->DrawLatex(0.44,0.88,centralityLatex.Data());
1059 //latexInfo1->DrawLatex(0.44,0.82,psiLatex.Data());
1060 latexInfo1->DrawLatex(0.44,0.82,pttLatex.Data());
1061 latexInfo1->DrawLatex(0.44,0.76,ptaLatex.Data());
1062
1063 pngName = "CorrelationFunction.Centrality";
1064 pngName += centralityArray[gCentrality-1];
1065 pngName += ".Psi";
1066 if((psiMin == -0.5)&&(psiMax == 0.5)) pngName += "InPlane.Ptt";
1067 else if((psiMin == 0.5)&&(psiMax == 1.5)) pngName += "Intermediate.Ptt";
1068 else if((psiMin == 1.5)&&(psiMax == 2.5)) pngName += "OutOfPlane.Ptt";
1069 else if((psiMin == 2.5)&&(psiMax == 3.5)) pngName += "Rest.Ptt";
1070 else pngName += "All.PttFrom";
1071 pngName += Form("%.1f",ptTriggerMin); pngName += "To";
1072 pngName += Form("%.1f",ptTriggerMax); pngName += "PtaFrom";
1073 pngName += Form("%.1f",ptAssociatedMin); pngName += "To";
1074 pngName += Form("%.1f",ptAssociatedMax);
1075 pngName += ".PositivePositive.png";
1076 cPP->SaveAs(pngName.Data());
1077
1078 //============================================================//
1079 //Get the -- correlation function
1080 TH2D *gHistNN = dynamic_cast<TH2D *>(f->Get("gHistNNCorrelationFunctions"));
1081 gHistNN->SetStats(kFALSE);
1082 gHistNN->SetTitle("");
1083 gHistNN->GetXaxis()->SetRangeUser(-1.45,1.45);
1084 gHistNN->GetXaxis()->CenterTitle();
1085 gHistNN->GetXaxis()->SetTitleOffset(1.2);
1086 gHistNN->GetYaxis()->CenterTitle();
1087 gHistNN->GetYaxis()->SetTitleOffset(1.2);
1088 gHistNN->GetZaxis()->SetTitleOffset(1.2);
1089 TCanvas *cNN = new TCanvas("cNN","",150,150,600,500);
1090 cNN->SetFillColor(10); cNN->SetHighLightColor(10);
1091 cNN->SetLeftMargin(0.15);
1092 gHistNN->DrawCopy("surf1fb");
1093 gPad->SetTheta(30); // default is 30
1094 gPad->SetPhi(-60); // default is 30
1095 gPad->Update();
1096
1097 latexInfo1->DrawLatex(0.44,0.88,centralityLatex.Data());
1098 //latexInfo1->DrawLatex(0.44,0.82,psiLatex.Data());
1099 latexInfo1->DrawLatex(0.44,0.82,pttLatex.Data());
1100 latexInfo1->DrawLatex(0.44,0.76,ptaLatex.Data());
1101
1102 pngName = "CorrelationFunction.Centrality";
1103 pngName += centralityArray[gCentrality-1];
1104 pngName += ".Psi";
1105 if((psiMin == -0.5)&&(psiMax == 0.5)) pngName += "InPlane.Ptt";
1106 else if((psiMin == 0.5)&&(psiMax == 1.5)) pngName += "Intermediate.Ptt";
1107 else if((psiMin == 1.5)&&(psiMax == 2.5)) pngName += "OutOfPlane.Ptt";
1108 else if((psiMin == 2.5)&&(psiMax == 3.5)) pngName += "Rest.Ptt";
1109 else pngName += "All.PttFrom";
1110 pngName += Form("%.1f",ptTriggerMin); pngName += "To";
1111 pngName += Form("%.1f",ptTriggerMax); pngName += "PtaFrom";
1112 pngName += Form("%.1f",ptAssociatedMin); pngName += "To";
1113 pngName += Form("%.1f",ptAssociatedMax);
1114 pngName += ".NegativeNegative.png";
1115 cNN->SaveAs(pngName.Data());
a38fd7f3 1116}
1117
5c9b4733 1118//____________________________________________________________//
1119void fitCorrelationFunctions(Int_t gCentrality = 1,
1120 Double_t psiMin = -0.5, Double_t psiMax = 3.5,
1121 Double_t ptTriggerMin = -1.,
1122 Double_t ptTriggerMax = -1.,
1123 Double_t ptAssociatedMin = -1.,
1124 Double_t ptAssociatedMax = -1.,
1125 TH2D *gHist) {
07d0a35c 1126
1127 cout<<"FITTING FUNCTION"<<endl;
1128
5c9b4733 1129 //near side peak: [1]*TMath::Exp(-TMath::Power((0.5*TMath::Power((x/[2]),2)+0.5*TMath::Power((y/[3]),2)),[4]))
1130 //away side ridge: [5]*TMath::Exp(-TMath::Power((0.5*TMath::Power(((y-TMath::Pi())/[6]),2)),[7]))
1131 //longitudinal ridge: [8]*TMath::Exp(-TMath::Power((0.5*TMath::Power((x/[9]),2)),[10]))
1132 //wing structures: [11]*TMath::Power(x,2)
1133 //flow contribution (v1 up to v4): 2.*([12]*TMath::Cos(y) + [13]*TMath::Cos(2.*y) + [14]*TMath::Cos(3.*y) + [15]*TMath::Cos(4.*y))
1134 TF2 *gFitFunction = new TF2("gFitFunction","[0]+[1]*TMath::Exp(-TMath::Power((0.5*TMath::Power((x/[2]),2)+0.5*TMath::Power((y/[3]),2)),[4]))+[5]*TMath::Exp(-TMath::Power((0.5*TMath::Power(((y-TMath::Pi())/[6]),2)),[7]))+[8]*TMath::Exp(-TMath::Power((0.5*TMath::Power((x/[9]),2)),[10]))+[11]*TMath::Power(x,2)+2.*[12]*([13]*TMath::Cos(y) + [14]*TMath::Cos(2.*y) + [15]*TMath::Cos(3.*y) + [16]*TMath::Cos(4.*y))",-2.0,2.0,-TMath::Pi()/2.,3.*TMath::Pi()/2.);
1135 gFitFunction->SetName("gFitFunction");
1136 //Normalization
1137 gFitFunction->SetParName(0,"N1"); gFitFunction->SetParameter(0,1.0);
1138 //near side peak
1139 gFitFunction->SetParName(1,"N_{near side}"); gFitFunction->SetParameter(1,0.3);
1140 gFitFunction->SetParName(2,"Sigma_{near side}(delta eta)"); gFitFunction->SetParameter(2,0.3);
1141 gFitFunction->SetParName(3,"Sigma_{near side}(delta phi)"); gFitFunction->SetParameter(3,0.1);
1142 gFitFunction->SetParName(4,"Exponent_{near side}"); gFitFunction->SetParameter(4,1.1);
1143 //away side ridge
1144 gFitFunction->SetParName(5,"N_{away side}"); gFitFunction->SetParameter(5,0.1);
1145 gFitFunction->SetParName(6,"Sigma_{away side}(delta phi)"); gFitFunction->SetParameter(6,1.1);
1146 gFitFunction->SetParName(7,"Exponent_{away side}"); gFitFunction->SetParameter(7,1.0);
1147 //longitudianl ridge
1148 gFitFunction->SetParName(8,"N_{long. ridge}"); gFitFunction->SetParameter(8,0.05);
1149 gFitFunction->SetParName(9,"Sigma_{long. ridge}(delta eta)"); gFitFunction->SetParameter(9,0.6);
1150 gFitFunction->SetParName(10,"Exponent_{long. ridge}"); gFitFunction->SetParameter(10,1.0);
1151 //wing structures
1152 gFitFunction->SetParName(11,"N_{wing}"); gFitFunction->SetParameter(11,0.01);
1153 //flow contribution
1154 gFitFunction->SetParName(12,"N_{flow}"); gFitFunction->SetParameter(12,0.2);
1155 gFitFunction->SetParName(13,"V1"); gFitFunction->SetParameter(13,0.005);
1156 gFitFunction->SetParName(14,"V2"); gFitFunction->SetParameter(14,0.1);
1157 gFitFunction->SetParName(15,"V3"); gFitFunction->SetParameter(15,0.05);
1158 gFitFunction->SetParName(16,"V4"); gFitFunction->SetParameter(16,0.005);
1159
1160 //Fitting the correlation function
1161 gHist->Fit("gFitFunction","nm");
1162
1163 //Cloning the histogram
1164 TString histoName = gHist->GetName(); histoName += "Fit";
1165 TH2D *gHistFit = new TH2D(histoName.Data(),";#Delta#eta;#Delta#varphi (rad);C(#Delta#eta,#Delta#varphi)",gHist->GetNbinsX(),gHist->GetXaxis()->GetXmin(),gHist->GetXaxis()->GetXmax(),gHist->GetNbinsY(),gHist->GetYaxis()->GetXmin(),gHist->GetYaxis()->GetXmax());
1166 TH2D *gHistResidual = dynamic_cast<TH2D *>(gHist->Clone());
1167 gHistResidual->SetName("gHistResidual");
1168 gHistResidual->Sumw2();
1169
1170 for(Int_t iBinDeltaEta = 1; iBinDeltaEta <= gHist->GetNbinsX(); iBinDeltaEta++) {
1171 for(Int_t iBinDeltaPhi = 1; iBinDeltaPhi <= gHist->GetNbinsY(); iBinDeltaPhi++) {
1172 gHistFit->SetBinContent(iBinDeltaEta,iBinDeltaPhi,gFitFunction->Eval(gHist->GetXaxis()->GetBinCenter(iBinDeltaEta),gHist->GetYaxis()->GetBinCenter(iBinDeltaPhi)));
1173 }
1174 }
1175 gHistResidual->Add(gHistFit,-1);
1176
1177 //Write to output file
1178 TString newFileName = "correlationFunctionFit";
1179 if(histoName.Contains("PN")) newFileName += "PN";
1180 else if(histoName.Contains("NP")) newFileName += "NP";
1181 else if(histoName.Contains("PP")) newFileName += "PP";
1182 else if(histoName.Contains("NN")) newFileName += "NN";
1183 newFileName += ".Centrality";
1184 newFileName += gCentrality; newFileName += ".Psi";
1185 if((psiMin == -0.5)&&(psiMax == 0.5)) newFileName += "InPlane.Ptt";
1186 else if((psiMin == 0.5)&&(psiMax == 1.5)) newFileName += "Intermediate.Ptt";
1187 else if((psiMin == 1.5)&&(psiMax == 2.5)) newFileName += "OutOfPlane.Ptt";
1188 else if((psiMin == 2.5)&&(psiMax == 3.5)) newFileName += "Rest.PttFrom";
1189 else newFileName += "All.PttFrom";
1190 newFileName += Form("%.1f",ptTriggerMin); newFileName += "To";
1191 newFileName += Form("%.1f",ptTriggerMax); newFileName += "PtaFrom";
1192 newFileName += Form("%.1f",ptAssociatedMin); newFileName += "To";
1193 newFileName += Form("%.1f",ptAssociatedMax);
1194 newFileName += ".root";
1195 TFile *newFile = TFile::Open(newFileName.Data(),"recreate");
5ea3b761 1196 gHist->Write();
5c9b4733 1197 gHistFit->Write();
1198 gHistResidual->Write();
1199 gFitFunction->Write();
1200 newFile->Close();
1201
1202
1203}