1 const Int_t numberOfCentralityBins = 10;
2 TString centralityArray[numberOfCentralityBins] = {"0-10","10-20","20-30","30-40","40-50","50-60","60-70","70-80","0-1","1-2"};
4 void drawBalanceFunctionPsi(const char* filename = "AnalysisResultsPsi.root",
6 Int_t gDeltaEtaDeltaPhi = 2,
7 Double_t psiMin = -0.5, Double_t psiMax = 0.5,
8 Double_t ptTriggerMin = -1.,
9 Double_t ptTriggerMax = -1.,
10 Double_t ptAssociatedMin = -1.,
11 Double_t ptAssociatedMax = -1.) {
12 //Macro that draws the BF distributions for each centrality bin
13 //for reaction plane dependent analysis
14 //Author: Panos.Christakoglou@nikhef.nl
15 //Load the PWG2ebye library
16 gSystem->Load("libANALYSIS.so");
17 gSystem->Load("libANALYSISalice.so");
18 gSystem->Load("libEventMixing.so");
19 gSystem->Load("libCORRFW.so");
20 gSystem->Load("libPWGTools.so");
21 gSystem->Load("libPWGCFebye.so");
23 //Prepare the objects and return them
24 TList *listBF = GetListOfObjects(filename,gCentrality,0);
25 TList *listBFShuffled = GetListOfObjects(filename,gCentrality,1);
26 TList *listBFMixed = GetListOfObjects(filename,gCentrality,2);
28 Printf("The TList object was not created");
32 draw(listBF,listBFShuffled,listBFMixed,
33 gCentrality,gDeltaEtaDeltaPhi,
35 ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
38 //______________________________________________________//
39 TList *GetListOfObjects(const char* filename,
40 Int_t gCentrality, Int_t kData = 1) {
41 //Get the TList objects (QA, bf, bf shuffled)
43 //kData == 1: shuffling
49 TFile *f = TFile::Open(filename);
50 if((!f)||(!f->IsOpen())) {
51 Printf("The file %s is not found. Aborting...",filename);
56 TDirectoryFile *dir = dynamic_cast<TDirectoryFile *>(f->Get("PWGCFEbyE.outputBalanceFunctionPsiAnalysis"));
58 Printf("The TDirectoryFile is not found. Aborting...",filename);
65 //cout<<"no shuffling - no mixing"<<endl;
66 listBFName = "listBFPsi_";
69 //cout<<"shuffling - no mixing"<<endl;
70 listBFName = "listBFPsiShuffled_";
73 //cout<<"no shuffling - mixing"<<endl;
74 listBFName = "listBFPsiMixed_";
76 listBFName += centralityArray[gCentrality-1];
77 listBF = dynamic_cast<TList *>(dir->Get(listBFName.Data()));
78 cout<<"======================================================="<<endl;
79 cout<<"List name: "<<listBF->GetName()<<endl;
85 histoName = "fHistPV0M";
87 histoName = "fHistP_shuffleV0M";
89 histoName = "fHistPV0M";
90 AliTHn *fHistP = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
92 Printf("fHistP %s not found!!!",histoName.Data());
95 fHistP->FillParent(); fHistP->DeleteContainers();
98 histoName = "fHistNV0M";
100 histoName = "fHistN_shuffleV0M";
102 histoName = "fHistNV0M";
103 AliTHn *fHistN = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
105 Printf("fHistN %s not found!!!",histoName.Data());
108 fHistN->FillParent(); fHistN->DeleteContainers();
111 histoName = "fHistPNV0M";
113 histoName = "fHistPN_shuffleV0M";
115 histoName = "fHistPNV0M";
116 AliTHn *fHistPN = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
118 Printf("fHistPN %s not found!!!",histoName.Data());
121 fHistPN->FillParent(); fHistPN->DeleteContainers();
124 histoName = "fHistNPV0M";
126 histoName = "fHistNP_shuffleV0M";
128 histoName = "fHistNPV0M";
129 AliTHn *fHistNP = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
131 Printf("fHistNP %s not found!!!",histoName.Data());
134 fHistNP->FillParent(); fHistNP->DeleteContainers();
137 histoName = "fHistPPV0M";
139 histoName = "fHistPP_shuffleV0M";
141 histoName = "fHistPPV0M";
142 AliTHn *fHistPP = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
144 Printf("fHistPP %s not found!!!",histoName.Data());
147 fHistPP->FillParent(); fHistPP->DeleteContainers();
150 histoName = "fHistNNV0M";
152 histoName = "fHistNN_shuffleV0M";
154 histoName = "fHistNNV0M";
155 AliTHn *fHistNN = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
157 Printf("fHistNN %s not found!!!",histoName.Data());
160 fHistNN->FillParent(); fHistNN->DeleteContainers();
165 //______________________________________________________//
166 void draw(TList *listBF, TList *listBFShuffled, TList *listBFMixed,
167 Int_t gCentrality, Int_t gDeltaEtaDeltaPhi,
168 Double_t psiMin, Double_t psiMax,
169 Double_t ptTriggerMin, Double_t ptTriggerMax,
170 Double_t ptAssociatedMin, Double_t ptAssociatedMax) {
171 gROOT->LoadMacro("~/SetPlotStyle.C");
173 gStyle->SetPalette(1,0);
175 const Int_t gRebin = gDeltaEtaDeltaPhi; //rebin by 2 the Delta phi projection
185 //Printf("=================");
186 hP = (AliTHn*) listBF->FindObject("fHistPV0M");
187 hN = (AliTHn*) listBF->FindObject("fHistNV0M");
188 hPN = (AliTHn*) listBF->FindObject("fHistPNV0M");
189 hNP = (AliTHn*) listBF->FindObject("fHistNPV0M");
190 hPP = (AliTHn*) listBF->FindObject("fHistPPV0M");
191 hNN = (AliTHn*) listBF->FindObject("fHistNNV0M");
193 AliBalancePsi *b = new AliBalancePsi();
201 //balance function shuffling
202 AliTHn *hPShuffled = NULL;
203 AliTHn *hNShuffled = NULL;
204 AliTHn *hPNShuffled = NULL;
205 AliTHn *hNPShuffled = NULL;
206 AliTHn *hPPShuffled = NULL;
207 AliTHn *hNNShuffled = NULL;
208 //listBFShuffled->ls();
209 hPShuffled = (AliTHn*) listBFShuffled->FindObject("fHistP_shuffleV0M");
210 hNShuffled = (AliTHn*) listBFShuffled->FindObject("fHistN_shuffleV0M");
211 hPNShuffled = (AliTHn*) listBFShuffled->FindObject("fHistPN_shuffleV0M");
212 hNPShuffled = (AliTHn*) listBFShuffled->FindObject("fHistNP_shuffleV0M");
213 hPPShuffled = (AliTHn*) listBFShuffled->FindObject("fHistPP_shuffleV0M");
214 hNNShuffled = (AliTHn*) listBFShuffled->FindObject("fHistNN_shuffleV0M");
216 AliBalancePsi *bShuffled = new AliBalancePsi();
217 bShuffled->SetHistNp(hPShuffled);
218 bShuffled->SetHistNn(hNShuffled);
219 bShuffled->SetHistNpn(hPNShuffled);
220 bShuffled->SetHistNnp(hNPShuffled);
221 bShuffled->SetHistNpp(hPPShuffled);
222 bShuffled->SetHistNnn(hNNShuffled);
224 //balance function mixing
225 AliTHn *hPMixed = NULL;
226 AliTHn *hNMixed = NULL;
227 AliTHn *hPNMixed = NULL;
228 AliTHn *hNPMixed = NULL;
229 AliTHn *hPPMixed = NULL;
230 AliTHn *hNNMixed = NULL;
232 hPMixed = (AliTHn*) listBFMixed->FindObject("fHistPV0M");
233 hNMixed = (AliTHn*) listBFMixed->FindObject("fHistNV0M");
234 hPNMixed = (AliTHn*) listBFMixed->FindObject("fHistPNV0M");
235 hNPMixed = (AliTHn*) listBFMixed->FindObject("fHistNPV0M");
236 hPPMixed = (AliTHn*) listBFMixed->FindObject("fHistPPV0M");
237 hNNMixed = (AliTHn*) listBFMixed->FindObject("fHistNNV0M");
239 AliBalancePsi *bMixed = new AliBalancePsi();
240 bMixed->SetHistNp(hPMixed);
241 bMixed->SetHistNn(hNMixed);
242 bMixed->SetHistNpn(hPNMixed);
243 bMixed->SetHistNnp(hNPMixed);
244 bMixed->SetHistNpp(hPPMixed);
245 bMixed->SetHistNnn(hNNMixed);
247 TH1D *gHistBalanceFunction;
248 TH1D *gHistBalanceFunctionShuffled;
249 TH1D *gHistBalanceFunctionMixed;
250 TH1D *gHistBalanceFunctionSubtracted;
251 TString histoTitle, pngName;
254 histoTitle = "Centrality: ";
255 histoTitle += centralityArray[gCentrality-1];
257 if((psiMin == -0.5)&&(psiMax == 0.5))
258 histoTitle += " (-7.5^{o} < #phi - #Psi_{2} < 7.5^{o})";
259 else if((psiMin == 0.5)&&(psiMax == 1.5))
260 histoTitle += " (37.5^{o} < #phi - #Psi_{2} < 52.5^{o})";
261 else if((psiMin == 1.5)&&(psiMax == 2.5))
262 histoTitle += " (82.5^{o} < #phi - #Psi_{2} < 97.5^{o})";
264 histoTitle += " (0^{o} < #phi - #Psi_{2} < 180^{o})";
266 //Raw balance function
267 gHistBalanceFunction = b->GetBalanceFunctionHistogram(0,gDeltaEtaDeltaPhi,psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
268 gHistBalanceFunction->SetMarkerStyle(20);
269 gHistBalanceFunction->SetTitle(histoTitle.Data());
270 gHistBalanceFunction->GetYaxis()->SetTitleOffset(1.3);
271 gHistBalanceFunction->SetName("gHistBalanceFunction");
273 //Shuffled balance function
274 gHistBalanceFunctionShuffled = bShuffled->GetBalanceFunctionHistogram(0,gDeltaEtaDeltaPhi,psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
275 gHistBalanceFunctionShuffled->SetMarkerStyle(24);
276 gHistBalanceFunctionShuffled->SetName("gHistBalanceFunctionShuffled");
278 //Mixed balance function
279 gHistBalanceFunctionMixed = bMixed->GetBalanceFunctionHistogram(0,gDeltaEtaDeltaPhi,psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
280 gHistBalanceFunctionMixed->SetMarkerStyle(25);
281 gHistBalanceFunctionMixed->SetName("gHistBalanceFunctionMixed");
283 //Subtracted balance function
284 gHistBalanceFunctionSubtracted = dynamic_cast<TH1D *>(gHistBalanceFunction->Clone());
285 gHistBalanceFunctionSubtracted->Add(gHistBalanceFunctionMixed,-1);
286 gHistBalanceFunctionSubtracted->Rebin(gRebin);
287 gHistBalanceFunctionSubtracted->SetMarkerStyle(20);
288 gHistBalanceFunctionSubtracted->SetTitle(histoTitle.Data());
289 gHistBalanceFunctionSubtracted->GetYaxis()->SetTitleOffset(1.3);
290 gHistBalanceFunctionSubtracted->SetName("gHistBalanceFunctionSubtracted");
292 TCanvas *c1 = new TCanvas("c1","",0,0,600,500);
293 c1->SetFillColor(10);
294 c1->SetHighLightColor(10);
295 c1->SetLeftMargin(0.15);
296 gHistBalanceFunction->DrawCopy("E");
297 gHistBalanceFunctionShuffled->DrawCopy("ESAME");
298 gHistBalanceFunctionMixed->DrawCopy("ESAME");
300 legend = new TLegend(0.18,0.62,0.45,0.82,"","brNDC");
301 legend->SetTextSize(0.045);
302 legend->SetTextFont(42);
303 legend->SetBorderSize(0);
304 legend->SetFillStyle(0);
305 legend->SetFillColor(10);
306 legend->SetMargin(0.25);
307 legend->SetShadowColor(10);
308 legend->AddEntry(gHistBalanceFunction,"Data","lp");
309 legend->AddEntry(gHistBalanceFunctionShuffled,"Shuffled data","lp");
310 legend->AddEntry(gHistBalanceFunctionMixed,"Mixed data","lp");
313 pngName = "BalanceFunctionDeltaEta.Centrality";
314 pngName += centralityArray[gCentrality-1];
315 pngName += ".Psi"; //pngName += psiMin; pngName += "To"; pngName += psiMax;
317 c1->SaveAs(pngName.Data());
319 GetWeightedMean(gHistBalanceFunction);
320 GetWeightedMean(gHistBalanceFunctionShuffled);
322 TString meanLatex, rmsLatex, skewnessLatex, kurtosisLatex;
323 meanLatex = "#mu = ";
324 meanLatex += Form("%.3f",gHistBalanceFunctionSubtracted->GetMean());
325 meanLatex += " #pm ";
326 meanLatex += Form("%.3f",gHistBalanceFunctionSubtracted->GetMeanError());
328 rmsLatex = "#sigma = ";
329 rmsLatex += Form("%.3f",gHistBalanceFunctionSubtracted->GetRMS());
331 rmsLatex += Form("%.3f",gHistBalanceFunctionSubtracted->GetRMSError());
333 skewnessLatex = "S = ";
334 skewnessLatex += Form("%.3f",gHistBalanceFunctionSubtracted->GetSkewness(1));
335 skewnessLatex += " #pm ";
336 skewnessLatex += Form("%.3f",gHistBalanceFunctionSubtracted->GetSkewness(11));
338 kurtosisLatex = "K = ";
339 kurtosisLatex += Form("%.3f",gHistBalanceFunctionSubtracted->GetKurtosis(1));
340 kurtosisLatex += " #pm ";
341 kurtosisLatex += Form("%.3f",gHistBalanceFunctionSubtracted->GetKurtosis(11));
342 Printf("Mean: %lf - Error: %lf",gHistBalanceFunctionSubtracted->GetMean(),gHistBalanceFunctionSubtracted->GetMeanError());
343 Printf("RMS: %lf - Error: %lf",gHistBalanceFunctionSubtracted->GetRMS(),gHistBalanceFunctionSubtracted->GetRMSError());
344 Printf("Skeweness: %lf - Error: %lf",gHistBalanceFunctionSubtracted->GetSkewness(1),gHistBalanceFunctionSubtracted->GetSkewness(11));
345 Printf("Kurtosis: %lf - Error: %lf",gHistBalanceFunctionSubtracted->GetKurtosis(1),gHistBalanceFunctionSubtracted->GetKurtosis(11));
347 TCanvas *c2 = new TCanvas("c2","",600,0,600,500);
348 c2->SetFillColor(10);
349 c2->SetHighLightColor(10);
350 c2->SetLeftMargin(0.15);
351 gHistBalanceFunctionSubtracted->DrawCopy("E");
353 TLatex *latex = new TLatex();
355 latex->SetTextSize(0.035);
356 latex->SetTextColor(1);
357 latex->DrawLatex(0.64,0.85,meanLatex.Data());
358 latex->DrawLatex(0.64,0.81,rmsLatex.Data());
359 latex->DrawLatex(0.64,0.77,skewnessLatex.Data());
360 latex->DrawLatex(0.64,0.73,kurtosisLatex.Data());
362 TString newFileName = "balanceFunction.Centrality";
363 newFileName += gCentrality; newFileName += ".In";
364 if(gDeltaEtaDeltaPhi == 1) newFileName += "DeltaEta.Psi";
365 else if(gDeltaEtaDeltaPhi == 2) newFileName += "DeltaPhi.Psi";
366 if((psiMin == -0.5)&&(psiMax == 0.5)) newFileName += "InPlane.Ptt";
367 else if((psiMin == 0.5)&&(psiMax == 1.5)) newFileName += "Intermediate.Ptt";
368 else if((psiMin == 1.5)&&(psiMax == 2.5)) newFileName += "OutOfPlane.Ptt";
369 else if((psiMin == 2.5)&&(psiMax == 3.5)) newFileName += "Rest.PttFrom";
370 else newFileName += "All.PttFrom";
371 newFileName += Form("%.1f",ptTriggerMin); newFileName += "To";
372 newFileName += Form("%.1f",ptTriggerMax); newFileName += "PtaFrom";
373 newFileName += Form("%.1f",ptAssociatedMin); newFileName += "To";
374 newFileName += Form("%.1f",ptAssociatedMax);
375 newFileName += ".root";
377 TFile *fOutput = new TFile(newFileName.Data(),"recreate");
379 gHistBalanceFunction->Write();
380 gHistBalanceFunctionShuffled->Write();
381 gHistBalanceFunctionMixed->Write();
382 gHistBalanceFunctionSubtracted->Write();
386 //____________________________________________________________________//
387 void GetWeightedMean(TH1D *gHistBalance, Int_t fStartBin = 1) {
388 //Prints the calculated width of the BF and its error
389 Double_t gSumXi = 0.0, gSumBi = 0.0, gSumBiXi = 0.0;
390 Double_t gSumBiXi2 = 0.0, gSumBi2Xi2 = 0.0;
391 Double_t gSumDeltaBi2 = 0.0, gSumXi2DeltaBi2 = 0.0;
392 Double_t deltaBalP2 = 0.0, integral = 0.0;
393 Double_t deltaErrorNew = 0.0;
395 //Retrieve this variables from Histogram
396 Int_t fNumberOfBins = gHistBalance->GetNbinsX();
397 Double_t fP2Step = gHistBalance->GetBinWidth(1); // assume equal binning!
399 cout<<"=================================================="<<endl;
400 cout<<"RECALCULATION OF BF WIDTH (StartBin = "<<fStartBin<<")"<<endl;
401 cout<<"HISTOGRAM has "<<fNumberOfBins<<" bins with bin size of "<<fP2Step<<endl;
402 cout<<"=================================================="<<endl;
403 for(Int_t i = 1; i <= fNumberOfBins; i++) {
404 // this is to simulate |Delta eta| or |Delta phi|
405 if(fNumberOfBins/2 - fStartBin + 1 < i && i < fNumberOfBins/2 + fStartBin ) continue;
407 cout<<"B: "<<gHistBalance->GetBinContent(i)<<"\t Error: "<<gHistBalance->GetBinError(i)<<"\t bin: "<<TMath::Abs(gHistBalance->GetBinCenter(i))<<endl;
409 gSumXi += TMath::Abs(gHistBalance->GetBinCenter(i)); // this is to simulate |Delta eta| or |Delta phi|
410 gSumBi += gHistBalance->GetBinContent(i);
411 gSumBiXi += gHistBalance->GetBinContent(i)*TMath::Abs(gHistBalance->GetBinCenter(i));
412 gSumBiXi2 += gHistBalance->GetBinContent(i)*TMath::Power(TMath::Abs(gHistBalance->GetBinCenter(i)),2);
413 gSumBi2Xi2 += TMath::Power(gHistBalance->GetBinContent(i),2)*TMath::Power(TMath::Abs(gHistBalance->GetBinCenter(i)),2);
414 gSumDeltaBi2 += TMath::Power(gHistBalance->GetBinError(i),2);
415 gSumXi2DeltaBi2 += TMath::Power(TMath::Abs(gHistBalance->GetBinCenter(i)),2) * TMath::Power(gHistBalance->GetBinError(i),2);
417 deltaBalP2 += fP2Step*TMath::Power(gHistBalance->GetBinError(i),2);
418 integral += fP2Step*gHistBalance->GetBinContent(i);
420 for(Int_t i = fStartBin; i < fNumberOfBins; i++)
421 deltaErrorNew += gHistBalance->GetBinError(i)*(TMath::Abs(gHistBalance->GetBinCenter(i))*gSumBi - gSumBiXi)/TMath::Power(gSumBi,2);
423 Double_t integralError = TMath::Sqrt(deltaBalP2);
425 Double_t delta = gSumBiXi / gSumBi;
426 Double_t deltaError = (gSumBiXi / gSumBi) * TMath::Sqrt(TMath::Power((TMath::Sqrt(gSumXi2DeltaBi2)/gSumBiXi),2) + TMath::Power((gSumDeltaBi2/gSumBi),2) );
427 cout<<"=================================================="<<endl;
428 cout<<"Width: "<<delta<<"\t Error: "<<deltaError<<endl;
429 cout<<"New error: "<<deltaErrorNew<<endl;
430 cout<<"Integral: "<<integral<<"\t Error: "<<integralError<<endl;
431 cout<<"=================================================="<<endl;
435 //______________________________________________________//
436 void drawBFPsi(const char* lhcPeriod = "LHC10h",
437 Int_t gCentrality = 1,
438 Int_t gDeltaEtaDeltaPhi = 2,
439 Double_t psiMin = -0.5, Double_t psiMax = 0.5,
440 Double_t ptTriggerMin = -1.,
441 Double_t ptTriggerMax = -1.,
442 Double_t ptAssociatedMin = -1.,
443 Double_t ptAssociatedMax = -1.) {
444 //Macro that draws the BF distributions for each centrality bin
445 //for reaction plane dependent analysis
446 //Author: Panos.Christakoglou@nikhef.nl
447 gROOT->LoadMacro("~/SetPlotStyle.C");
451 TString filename = lhcPeriod; filename +="/PttFrom";
452 filename += Form("%.1f",ptTriggerMin); filename += "To";
453 filename += Form("%.1f",ptTriggerMax); filename += "PtaFrom";
454 filename += Form("%.1f",ptAssociatedMin); filename += "To";
455 filename += Form("%.1f",ptAssociatedMax);
456 filename += "/balanceFunction.Centrality";
457 filename += gCentrality; filename += ".In";
458 if(gDeltaEtaDeltaPhi == 1) filename += "DeltaEta.Psi";
459 else if(gDeltaEtaDeltaPhi == 2) filename += "DeltaPhi.Psi";
460 if((psiMin == -0.5)&&(psiMax == 0.5)) filename += "InPlane.Ptt";
461 else if((psiMin == 0.5)&&(psiMax == 1.5)) filename += "Intermediate.Ptt";
462 else if((psiMin == 1.5)&&(psiMax == 2.5)) filename += "OutOfPlane.Ptt";
463 else if((psiMin == 2.5)&&(psiMax == 3.5)) filename += "Rest.Ptt";
464 else filename += "All.PttFrom";
465 filename += Form("%.1f",ptTriggerMin); filename += "To";
466 filename += Form("%.1f",ptTriggerMax); filename += "PtaFrom";
467 filename += Form("%.1f",ptAssociatedMin); filename += "To";
468 filename += Form("%.1f",ptAssociatedMax); filename += ".root";
471 TFile *f = TFile::Open(filename.Data());
472 if((!f)||(!f->IsOpen())) {
473 Printf("The file %s is not found. Aborting...",filename);
478 //Raw balance function
479 TH1D *gHistBalanceFunction = dynamic_cast<TH1D *>(f->Get("gHistBalanceFunction"));
480 gHistBalanceFunction->SetStats(kFALSE);
481 gHistBalanceFunction->SetMarkerStyle(20);
482 gHistBalanceFunction->GetYaxis()->SetTitleOffset(1.3);
483 if(gDeltaEtaDeltaPhi == 2) {
484 gHistBalanceFunction->GetYaxis()->SetTitle("B(#Delta #varphi)");
485 gHistBalanceFunction->GetXaxis()->SetTitle("#Delta #varphi (deg.)");
488 //Shuffled balance function
489 TH1D *gHistBalanceFunctionShuffled = dynamic_cast<TH1D *>(f->Get("gHistBalanceFunctionShuffled"));
490 gHistBalanceFunction->SetStats(kFALSE);
491 gHistBalanceFunctionShuffled->SetMarkerStyle(24);
493 //Mixed balance function
494 TH1D *gHistBalanceFunctionMixed = dynamic_cast<TH1D *>(f->Get("gHistBalanceFunctionMixed"));
495 gHistBalanceFunction->SetStats(kFALSE);
496 gHistBalanceFunctionMixed->SetMarkerStyle(25);
498 //Subtracted balance function
499 TH1D *gHistBalanceFunctionSubtracted = dynamic_cast<TH1D *>(f->Get("gHistBalanceFunctionSubtracted"));
500 gHistBalanceFunctionSubtracted->SetStats(kFALSE);
501 gHistBalanceFunctionSubtracted->SetMarkerStyle(20);
502 gHistBalanceFunctionSubtracted->GetYaxis()->SetTitleOffset(1.3);
503 if(gDeltaEtaDeltaPhi == 2) {
504 gHistBalanceFunctionSubtracted->GetYaxis()->SetTitle("B(#Delta #varphi)");
505 gHistBalanceFunctionSubtracted->GetXaxis()->SetTitle("#Delta #varphi (deg.)");
511 TString centralityLatex = "Centrality: ";
512 centralityLatex += centralityArray[gCentrality-1];
513 centralityLatex += "%";
516 if((psiMin == -0.5)&&(psiMax == 0.5))
517 psiLatex = " -7.5^{o} < #varphi - #Psi_{2} < 7.5^{o}";
518 else if((psiMin == 0.5)&&(psiMax == 1.5))
519 psiLatex = " 37.5^{o} < #varphi - #Psi_{2} < 52.5^{o}";
520 else if((psiMin == 1.5)&&(psiMax == 2.5))
521 psiLatex = " 82.5^{o} < #varphi - #Psi_{2} < 97.5^{o}";
523 psiLatex = " 0^{o} < #varphi - #Psi_{2} < 180^{o}";
525 TString pttLatex = Form("%.1f",ptTriggerMin);
526 pttLatex += " < p_{T}^{t} < "; pttLatex += Form("%.1f",ptTriggerMax);
527 pttLatex += " GeV/c";
529 TString ptaLatex = Form("%.1f",ptAssociatedMin);
530 ptaLatex += " < p_{T}^{a} < "; ptaLatex += Form("%.1f",ptAssociatedMax);
531 ptaLatex += " GeV/c";
534 TCanvas *c1 = new TCanvas("c1","",0,0,600,500);
535 c1->SetFillColor(10); c1->SetHighLightColor(10);
536 c1->SetLeftMargin(0.17); c1->SetTopMargin(0.05);
537 gHistBalanceFunction->SetTitle("");
538 gHistBalanceFunction->GetYaxis()->SetTitleOffset(1.4);
539 gHistBalanceFunction->GetYaxis()->SetNdivisions(10);
540 gHistBalanceFunction->GetXaxis()->SetNdivisions(10);
541 gHistBalanceFunction->DrawCopy("E");
542 gHistBalanceFunctionShuffled->DrawCopy("ESAME");
543 gHistBalanceFunctionMixed->DrawCopy("ESAME");
545 legend = new TLegend(0.2,0.72,0.45,0.92,"","brNDC");
546 legend->SetTextSize(0.05);
547 legend->SetTextFont(42);
548 legend->SetBorderSize(0);
549 legend->SetFillStyle(0);
550 legend->SetFillColor(10);
551 legend->SetMargin(0.25);
552 legend->SetShadowColor(10);
553 legend->AddEntry(gHistBalanceFunction,"Data","lp");
554 legend->AddEntry(gHistBalanceFunctionShuffled,"Shuffled data","lp");
555 legend->AddEntry(gHistBalanceFunctionMixed,"Mixed data","lp");
558 TLatex *latexInfo1 = new TLatex();
559 latexInfo1->SetNDC();
560 latexInfo1->SetTextSize(0.04);
561 latexInfo1->SetTextColor(1);
562 latexInfo1->DrawLatex(0.58,0.88,centralityLatex.Data());
563 latexInfo1->DrawLatex(0.58,0.82,psiLatex.Data());
564 latexInfo1->DrawLatex(0.58,0.76,pttLatex.Data());
565 latexInfo1->DrawLatex(0.58,0.70,ptaLatex.Data());
567 //pngName = "BalanceFunctionDeltaEta.Centrality";
568 //pngName += centralityArray[gCentrality-1];
569 //pngName += ".Psi"; //pngName += psiMin; pngName += "To"; pngName += psiMax;
571 //c1->SaveAs(pngName.Data());
573 TCanvas *c2 = new TCanvas("c2","",600,0,600,500);
574 c2->SetFillColor(10); c2->SetHighLightColor(10);
575 c2->SetLeftMargin(0.17); c2->SetTopMargin(0.05);
576 gHistBalanceFunctionSubtracted->SetTitle("");
577 gHistBalanceFunctionSubtracted->GetYaxis()->SetTitleOffset(1.4);
578 gHistBalanceFunctionSubtracted->GetYaxis()->SetNdivisions(10);
579 gHistBalanceFunctionSubtracted->GetXaxis()->SetNdivisions(10);
580 gHistBalanceFunctionSubtracted->DrawCopy("E");
582 //Opening the output ascii files
583 TString filenameMean = "meanIn";
584 TString filenameSigma = "sigmaIn";
585 TString filenameSkewness = "skewnessIn";
586 TString filenameKurtosis = "kurtosisIn";
587 if(gDeltaEtaDeltaPhi == 1) {
588 filenameMean += "DeltaEta.Psi"; filenameSigma += "DeltaEta.Psi";
589 filenameSkewness += "DeltaEta.Psi"; filenameKurtosis += "DeltaEta.Psi";
591 else if(gDeltaEtaDeltaPhi == 2) {
592 filenameMean += "DeltaPhi.Psi"; filenameSigma += "DeltaPhi.Psi";
593 filenameSkewness += "DeltaPhi.Psi"; filenameKurtosis += "DeltaPhi.Psi";
595 if((psiMin == -0.5)&&(psiMax == 0.5)) {
596 filenameMean += "InPlane.Ptt"; filenameSigma += "InPlane.Ptt";
597 filenameSkewness += "InPlane.Ptt"; filenameKurtosis += "InPlane.Ptt";
599 else if((psiMin == 0.5)&&(psiMax == 1.5)) {
600 filenameMean += "Intermediate.Ptt"; filenameSigma += "Intermediate.Ptt";
601 filenameSkewness += "Intermediate.Ptt";
602 filenameKurtosis += "Intermediate.Ptt";
604 else if((psiMin == 1.5)&&(psiMax == 2.5)) {
605 filenameMean += "OutOfPlane.Ptt"; filenameSigma += "OutOfPlane.Ptt";
606 filenameSkewness += "OutOfPlane.Ptt";
607 filenameKurtosis += "OutOfPlane.Ptt";
609 else if((psiMin == 2.5)&&(psiMax == 3.5)) {
610 filenameMean += "Rest.Ptt"; filenameSigma += "Rest.Ptt";
611 filenameSkewness += "Rest.Ptt"; filenameKurtosis += "Rest.Ptt";
614 filenameMean += "All.Ptt"; filenameSigma += "All.Ptt";
615 filenameSkewness += "All.Ptt"; filenameKurtosis += "All.Ptt";
617 filenameMean += Form("%.1f",ptTriggerMin); filenameMean += "To";
618 filenameMean += Form("%.1f",ptTriggerMax); filenameMean += "PtaFrom";
619 filenameMean += Form("%.1f",ptAssociatedMin); filenameMean += "To";
620 filenameMean += Form("%.1f",ptAssociatedMax); filenameMean += ".txt";
621 filenameSigma += Form("%.1f",ptTriggerMin); filenameSigma += "To";
622 filenameSigma += Form("%.1f",ptTriggerMax); filenameSigma += "PtaFrom";
623 filenameSigma += Form("%.1f",ptAssociatedMin); filenameSigma += "To";
624 filenameSigma += Form("%.1f",ptAssociatedMax); filenameSigma += ".txt";
625 filenameSkewness += Form("%.1f",ptTriggerMin); filenameSkewness += "To";
626 filenameSkewness += Form("%.1f",ptTriggerMax); filenameSkewness += "PtaFrom";
627 filenameSkewness += Form("%.1f",ptAssociatedMin); filenameSkewness += "To";
628 filenameSkewness += Form("%.1f",ptAssociatedMax);
629 filenameSkewness += ".txt";
630 filenameKurtosis += Form("%.1f",ptTriggerMin); filenameKurtosis += "To";
631 filenameKurtosis += Form("%.1f",ptTriggerMax); filenameKurtosis += "PtaFrom";
632 filenameKurtosis += Form("%.1f",ptAssociatedMin); filenameKurtosis += "To";
633 filenameKurtosis += Form("%.1f",ptAssociatedMax);
634 filenameKurtosis += ".txt";
636 TString meanLatex, rmsLatex, skewnessLatex, kurtosisLatex;
637 meanLatex = "#mu = ";
638 meanLatex += Form("%.3f",gHistBalanceFunctionSubtracted->GetMean());
639 meanLatex += " #pm ";
640 meanLatex += Form("%.3f",gHistBalanceFunctionSubtracted->GetMeanError());
641 ofstream fileMean(filenameMean.Data(),ios::app);
642 fileMean << gCentrality << " " << gHistBalanceFunctionSubtracted->GetMean() << " " <<gHistBalanceFunctionSubtracted->GetMeanError()<<endl;
645 rmsLatex = "#sigma = ";
646 rmsLatex += Form("%.3f",gHistBalanceFunctionSubtracted->GetRMS());
648 rmsLatex += Form("%.3f",gHistBalanceFunctionSubtracted->GetRMSError());
649 ofstream fileSigma(filenameSigma.Data(),ios::app);
650 fileSigma << gCentrality << " " << gHistBalanceFunctionSubtracted->GetRMS() << " " <<gHistBalanceFunctionSubtracted->GetRMSError()<<endl;
653 skewnessLatex = "S = ";
654 skewnessLatex += Form("%.3f",gHistBalanceFunctionSubtracted->GetSkewness(1));
655 skewnessLatex += " #pm ";
656 skewnessLatex += Form("%.3f",gHistBalanceFunctionSubtracted->GetSkewness(11));
657 ofstream fileSkewness(filenameSkewness.Data(),ios::app);
658 fileSkewness << gCentrality << " " << gHistBalanceFunctionSubtracted->GetSkewness(1) << " " <<gHistBalanceFunctionSubtracted->GetSkewness(11)<<endl;
659 fileSkewness.close();
661 kurtosisLatex = "K = ";
662 kurtosisLatex += Form("%.3f",gHistBalanceFunctionSubtracted->GetKurtosis(1));
663 kurtosisLatex += " #pm ";
664 kurtosisLatex += Form("%.3f",gHistBalanceFunctionSubtracted->GetKurtosis(11));
665 ofstream fileKurtosis(filenameKurtosis.Data(),ios::app);
666 fileKurtosis << gCentrality << " " << gHistBalanceFunctionSubtracted->GetKurtosis(1) << " " <<gHistBalanceFunctionSubtracted->GetKurtosis(11)<<endl;
667 fileKurtosis.close();
669 Printf("Mean: %lf - Error: %lf",gHistBalanceFunctionSubtracted->GetMean(),gHistBalanceFunctionSubtracted->GetMeanError());
670 Printf("RMS: %lf - Error: %lf",gHistBalanceFunctionSubtracted->GetRMS(),gHistBalanceFunctionSubtracted->GetRMSError());
671 Printf("Skeweness: %lf - Error: %lf",gHistBalanceFunctionSubtracted->GetSkewness(1),gHistBalanceFunctionSubtracted->GetSkewness(11));
672 Printf("Kurtosis: %lf - Error: %lf",gHistBalanceFunctionSubtracted->GetKurtosis(1),gHistBalanceFunctionSubtracted->GetKurtosis(11));
674 latexInfo1->DrawLatex(0.18,0.88,centralityLatex.Data());
675 latexInfo1->DrawLatex(0.18,0.82,psiLatex.Data());
676 latexInfo1->DrawLatex(0.18,0.76,pttLatex.Data());
677 latexInfo1->DrawLatex(0.18,0.70,ptaLatex.Data());
679 TLatex *latexResults = new TLatex();
680 latexResults->SetNDC();
681 latexResults->SetTextSize(0.04);
682 latexResults->SetTextColor(1);
683 latexResults->DrawLatex(0.6,0.88,meanLatex.Data());
684 latexResults->DrawLatex(0.6,0.82,rmsLatex.Data());
685 latexResults->DrawLatex(0.6,0.76,skewnessLatex.Data());
686 latexResults->DrawLatex(0.6,0.70,kurtosisLatex.Data());