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 const Int_t gRebin = 1;
5 void drawBalanceFunction2DPsi(const char* filename = "AnalysisResultsPsi.root",
8 const char* gCentralityEstimator = 0x0,
9 Bool_t kShowShuffled = kFALSE,
10 Bool_t kShowMixed = kTRUE,
11 Double_t psiMin = -0.5, Double_t psiMax = 0.5,
12 Double_t ptTriggerMin = -1.,
13 Double_t ptTriggerMax = -1.,
14 Double_t ptAssociatedMin = -1.,
15 Double_t ptAssociatedMax = -1.,
16 Bool_t k2pMethod = kFALSE) {
17 //Macro that draws the BF distributions for each centrality bin
18 //for reaction plane dependent analysis
19 //Author: Panos.Christakoglou@nikhef.nl
20 //Load the PWG2ebye library
21 gSystem->Load("libANALYSIS.so");
22 gSystem->Load("libANALYSISalice.so");
23 gSystem->Load("libEventMixing.so");
24 gSystem->Load("libCORRFW.so");
25 gSystem->Load("libPWGTools.so");
26 gSystem->Load("libPWGCFebye.so");
28 gROOT->LoadMacro("~/SetPlotStyle.C");
30 gStyle->SetPalette(1,0);
32 //Prepare the objects and return them
33 TList *listBF = GetListOfObjects(filename,gCentrality,gBit,gCentralityEstimator,0);
34 TList *listBFShuffled = NULL;
35 if(kShowShuffled) listBFShuffled = GetListOfObjects(filename,gCentrality,gBit,gCentralityEstimator,1);
36 TList *listBFMixed = NULL;
37 if(kShowMixed) listBFMixed = GetListOfObjects(filename,gCentrality,gBit,gCentralityEstimator,2);
39 Printf("The TList object was not created");
43 draw(listBF,listBFShuffled,listBFMixed,gCentrality,psiMin,psiMax,
44 ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax,
48 //______________________________________________________//
49 TList *GetListOfObjects(const char* filename,
52 const char *gCentralityEstimator,
54 //Get the TList objects (QA, bf, bf shuffled)
58 TFile *f = TFile::Open(filename,"UPDATE");
59 if((!f)||(!f->IsOpen())) {
60 Printf("The file %s is not found. Aborting...",filename);
65 TDirectoryFile *dir = dynamic_cast<TDirectoryFile *>(f->Get("PWGCFEbyE.outputBalanceFunctionPsiAnalysis"));
67 Printf("The TDirectoryFile is not found. Aborting...",filename);
74 //cout<<"no shuffling - no mixing"<<endl;
75 listBFName = "listBFPsi_";
78 //cout<<"shuffling - no mixing"<<endl;
79 listBFName = "listBFPsiShuffled_";
82 //cout<<"no shuffling - mixing"<<endl;
83 listBFName = "listBFPsiMixed_";
85 listBFName += centralityArray[gCentrality-1];
87 listBFName += "_Bit"; listBFName += gBit; }
88 if(gCentralityEstimator) {
89 listBFName += "_"; listBFName += gCentralityEstimator;}
91 // histograms were already retrieved (in first iteration)
92 if(dir->Get(Form("%s_histograms",listBFName.Data()))){
93 //cout<<"second iteration"<<endl;
94 listBF = dynamic_cast<TList *>(dir->Get(Form("%s_histograms",listBFName.Data())));
97 // histograms were not yet retrieved (this is the first iteration)
100 listBF = dynamic_cast<TList *>(dir->Get(listBFName.Data()));
101 cout<<"======================================================="<<endl;
102 cout<<"List name (check): "<<listBFName.Data()<<endl;
103 cout<<"List name: "<<listBF->GetName()<<endl;
109 histoName = "fHistPV0M";
111 histoName = "fHistP_shuffleV0M";
113 histoName = "fHistPV0M";
114 AliTHn *fHistP = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
116 Printf("fHistP %s not found!!!",histoName.Data());
119 fHistP->FillParent(); fHistP->DeleteContainers();
122 histoName = "fHistNV0M";
124 histoName = "fHistN_shuffleV0M";
126 histoName = "fHistNV0M";
127 AliTHn *fHistN = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
129 Printf("fHistN %s not found!!!",histoName.Data());
132 fHistN->FillParent(); fHistN->DeleteContainers();
135 histoName = "fHistPNV0M";
137 histoName = "fHistPN_shuffleV0M";
139 histoName = "fHistPNV0M";
140 AliTHn *fHistPN = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
142 Printf("fHistPN %s not found!!!",histoName.Data());
145 fHistPN->FillParent(); fHistPN->DeleteContainers();
148 histoName = "fHistNPV0M";
150 histoName = "fHistNP_shuffleV0M";
152 histoName = "fHistNPV0M";
153 AliTHn *fHistNP = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
155 Printf("fHistNP %s not found!!!",histoName.Data());
158 fHistNP->FillParent(); fHistNP->DeleteContainers();
161 histoName = "fHistPPV0M";
163 histoName = "fHistPP_shuffleV0M";
165 histoName = "fHistPPV0M";
166 AliTHn *fHistPP = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
168 Printf("fHistPP %s not found!!!",histoName.Data());
171 fHistPP->FillParent(); fHistPP->DeleteContainers();
174 histoName = "fHistNNV0M";
176 histoName = "fHistNN_shuffleV0M";
178 histoName = "fHistNNV0M";
179 AliTHn *fHistNN = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
181 Printf("fHistNN %s not found!!!",histoName.Data());
184 fHistNN->FillParent(); fHistNN->DeleteContainers();
187 listBF->Write(Form("%s_histograms",listBFName.Data()), TObject::kSingleKey);
196 //______________________________________________________//
197 void draw(TList *listBF, TList *listBFShuffled, TList *listBFMixed,
198 Int_t gCentrality, Double_t psiMin, Double_t psiMax,
199 Double_t ptTriggerMin, Double_t ptTriggerMax,
200 Double_t ptAssociatedMin, Double_t ptAssociatedMax,
201 Bool_t k2pMethod = kFALSE) {
210 //Printf("=================");
211 hP = (AliTHn*) listBF->FindObject("fHistPV0M");
212 hP->SetName("gHistP");
213 hN = (AliTHn*) listBF->FindObject("fHistNV0M");
214 hN->SetName("gHistN");
215 hPN = (AliTHn*) listBF->FindObject("fHistPNV0M");
216 hPN->SetName("gHistPN");
217 hNP = (AliTHn*) listBF->FindObject("fHistNPV0M");
218 hNP->SetName("gHistNP");
219 hPP = (AliTHn*) listBF->FindObject("fHistPPV0M");
220 hPP->SetName("gHistPP");
221 hNN = (AliTHn*) listBF->FindObject("fHistNNV0M");
222 hNN->SetName("gHistNN");
224 AliBalancePsi *b = new AliBalancePsi();
232 //balance function shuffling
233 AliTHn *hPShuffled = NULL;
234 AliTHn *hNShuffled = NULL;
235 AliTHn *hPNShuffled = NULL;
236 AliTHn *hNPShuffled = NULL;
237 AliTHn *hPPShuffled = NULL;
238 AliTHn *hNNShuffled = NULL;
240 //listBFShuffled->ls();
242 hPShuffled = (AliTHn*) listBFShuffled->FindObject("fHistP_shuffleV0M");
243 hPShuffled->SetName("gHistPShuffled");
244 hNShuffled = (AliTHn*) listBFShuffled->FindObject("fHistN_shuffleV0M");
245 hNShuffled->SetName("gHistNShuffled");
246 hPNShuffled = (AliTHn*) listBFShuffled->FindObject("fHistPN_shuffleV0M");
247 hPNShuffled->SetName("gHistPNShuffled");
248 hNPShuffled = (AliTHn*) listBFShuffled->FindObject("fHistNP_shuffleV0M");
249 hNPShuffled->SetName("gHistNPShuffled");
250 hPPShuffled = (AliTHn*) listBFShuffled->FindObject("fHistPP_shuffleV0M");
251 hPPShuffled->SetName("gHistPPShuffled");
252 hNNShuffled = (AliTHn*) listBFShuffled->FindObject("fHistNN_shuffleV0M");
253 hNNShuffled->SetName("gHistNNShuffled");
255 AliBalancePsi *bShuffled = new AliBalancePsi();
256 bShuffled->SetHistNp(hPShuffled);
257 bShuffled->SetHistNn(hNShuffled);
258 bShuffled->SetHistNpn(hPNShuffled);
259 bShuffled->SetHistNnp(hNPShuffled);
260 bShuffled->SetHistNpp(hPPShuffled);
261 bShuffled->SetHistNnn(hNNShuffled);
264 //balance function mixing
265 AliTHn *hPMixed = NULL;
266 AliTHn *hNMixed = NULL;
267 AliTHn *hPNMixed = NULL;
268 AliTHn *hNPMixed = NULL;
269 AliTHn *hPPMixed = NULL;
270 AliTHn *hNNMixed = NULL;
275 hPMixed = (AliTHn*) listBFMixed->FindObject("fHistPV0M");
276 hPMixed->SetName("gHistPMixed");
277 hNMixed = (AliTHn*) listBFMixed->FindObject("fHistNV0M");
278 hNMixed->SetName("gHistNMixed");
279 hPNMixed = (AliTHn*) listBFMixed->FindObject("fHistPNV0M");
280 hPNMixed->SetName("gHistPNMixed");
281 hNPMixed = (AliTHn*) listBFMixed->FindObject("fHistNPV0M");
282 hNPMixed->SetName("gHistNPMixed");
283 hPPMixed = (AliTHn*) listBFMixed->FindObject("fHistPPV0M");
284 hPPMixed->SetName("gHistPPMixed");
285 hNNMixed = (AliTHn*) listBFMixed->FindObject("fHistNNV0M");
286 hNNMixed->SetName("gHistNNMixed");
288 AliBalancePsi *bMixed = new AliBalancePsi();
289 bMixed->SetHistNp(hPMixed);
290 bMixed->SetHistNn(hNMixed);
291 bMixed->SetHistNpn(hPNMixed);
292 bMixed->SetHistNnp(hNPMixed);
293 bMixed->SetHistNpp(hPPMixed);
294 bMixed->SetHistNnn(hNNMixed);
297 TH2D *gHistBalanceFunction;
298 TH2D *gHistBalanceFunctionSubtracted;
299 TH2D *gHistBalanceFunctionShuffled;
300 TH2D *gHistBalanceFunctionMixed;
301 TString histoTitle, pngName;
303 histoTitle = "Centrality: ";
304 histoTitle += centralityArray[gCentrality-1];
306 if((psiMin == -0.5)&&(psiMax == 0.5))
307 histoTitle += " (-7.5^{o} < #varphi - #Psi_{2} < 7.5^{o})";
308 else if((psiMin == 0.5)&&(psiMax == 1.5))
309 histoTitle += " (37.5^{o} < #varphi - #Psi_{2} < 52.5^{o})";
310 else if((psiMin == 1.5)&&(psiMax == 2.5))
311 histoTitle += " (82.5^{o} < #varphi - #Psi_{2} < 97.5^{o})";
313 histoTitle += " (0^{o} < #varphi - #Psi_{2} < 180^{o})";
317 gHistBalanceFunction = b->GetBalanceFunctionDeltaEtaDeltaPhi2pMethod(psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax,bMixed);
319 cerr<<"NO MIXED BF BUT REQUESTED CORRECTING WITH IT! --> FAIL"<<endl;
323 gHistBalanceFunction = b->GetBalanceFunctionDeltaEtaDeltaPhi(psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
324 gHistBalanceFunction->SetTitle(histoTitle.Data());
325 gHistBalanceFunction->GetYaxis()->SetTitleOffset(1.3);
326 gHistBalanceFunction->SetName("gHistBalanceFunction");
332 gHistBalanceFunctionShuffled = bShuffled->GetBalanceFunctionDeltaEtaDeltaPhi2pMethod(psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax,bMixed);
334 cerr<<"NO MIXED BF BUT REQUESTED CORRECTING WITH IT! --> FAIL"<<endl;
338 gHistBalanceFunctionShuffled = bShuffled->GetBalanceFunctionDeltaEtaDeltaPhi(psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
339 gHistBalanceFunctionShuffled->SetTitle(histoTitle.Data());
340 gHistBalanceFunctionShuffled->GetYaxis()->SetTitleOffset(1.3);
341 gHistBalanceFunctionShuffled->SetName("gHistBalanceFunctionShuffled");
347 gHistBalanceFunctionMixed = bMixed->GetBalanceFunctionDeltaEtaDeltaPhi2pMethod(psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax,bMixed);
349 cerr<<"NO MIXED BF BUT REQUESTED CORRECTING WITH IT! --> FAIL"<<endl;
353 gHistBalanceFunctionMixed = bMixed->GetBalanceFunctionDeltaEtaDeltaPhi(psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
354 gHistBalanceFunctionMixed->SetTitle(histoTitle.Data());
355 gHistBalanceFunctionMixed->GetYaxis()->SetTitleOffset(1.3);
356 gHistBalanceFunctionMixed->SetName("gHistBalanceFunctionMixed");
358 gHistBalanceFunctionSubtracted = dynamic_cast<TH2D *>(gHistBalanceFunction->Clone());
359 gHistBalanceFunctionSubtracted->Add(gHistBalanceFunctionMixed,-1);
360 gHistBalanceFunctionSubtracted->SetTitle(histoTitle.Data());
361 gHistBalanceFunctionSubtracted->GetYaxis()->SetTitleOffset(1.3);
362 gHistBalanceFunctionSubtracted->SetName("gHistBalanceFunctionSubtracted");
366 TCanvas *c1 = new TCanvas("c1","",0,0,600,500);
367 c1->SetFillColor(10);
368 c1->SetHighLightColor(10);
369 c1->SetLeftMargin(0.15);
370 gHistBalanceFunction->DrawCopy("lego2");
371 gPad->SetTheta(30); // default is 30
372 gPad->SetPhi(-60); // default is 30
374 TCanvas *c1a = new TCanvas("c1a","",600,0,600,500);
375 c1a->SetFillColor(10);
376 c1a->SetHighLightColor(10);
377 c1a->SetLeftMargin(0.15);
378 gHistBalanceFunction->DrawCopy("colz");
381 TCanvas *c2 = new TCanvas("c2","",100,100,600,500);
382 c2->SetFillColor(10);
383 c2->SetHighLightColor(10);
384 c2->SetLeftMargin(0.15);
385 gHistBalanceFunctionShuffled->DrawCopy("lego2");
386 gPad->SetTheta(30); // default is 30
387 gPad->SetPhi(-60); // default is 30
389 TCanvas *c2a = new TCanvas("c2a","",700,100,600,500);
390 c2a->SetFillColor(10);
391 c2a->SetHighLightColor(10);
392 c2a->SetLeftMargin(0.15);
393 gHistBalanceFunctionShuffled->DrawCopy("colz");
397 TCanvas *c3 = new TCanvas("c3","",200,200,600,500);
398 c3->SetFillColor(10);
399 c3->SetHighLightColor(10);
400 c3->SetLeftMargin(0.15);
401 gHistBalanceFunctionMixed->DrawCopy("lego2");
402 gPad->SetTheta(30); // default is 30
403 gPad->SetPhi(-60); // default is 30
405 TCanvas *c3a = new TCanvas("c3a","",800,200,600,500);
406 c3a->SetFillColor(10);
407 c3a->SetHighLightColor(10);
408 c3a->SetLeftMargin(0.15);
409 gHistBalanceFunctionMixed->DrawCopy("colz");
411 TCanvas *c4 = new TCanvas("c4","",300,300,600,500);
412 c4->SetFillColor(10);
413 c4->SetHighLightColor(10);
414 c4->SetLeftMargin(0.15);
415 gHistBalanceFunctionSubtracted->DrawCopy("lego2");
416 gPad->SetTheta(30); // default is 30
417 gPad->SetPhi(-60); // default is 30
419 TCanvas *c4a = new TCanvas("c4a","",900,300,600,500);
420 c4a->SetFillColor(10);
421 c4a->SetHighLightColor(10);
422 c4a->SetLeftMargin(0.15);
423 gHistBalanceFunctionSubtracted->DrawCopy("colz");
425 fitbalanceFunction(gCentrality, psiMin = -0.5, psiMax,
426 ptTriggerMin, ptTriggerMax,
427 ptAssociatedMin, ptAssociatedMax,
428 gHistBalanceFunctionSubtracted);
431 TString newFileName = "balanceFunction2D.Centrality";
432 newFileName += gCentrality; newFileName += ".Psi";
433 if((psiMin == -0.5)&&(psiMax == 0.5)) newFileName += "InPlane.Ptt";
434 else if((psiMin == 0.5)&&(psiMax == 1.5)) newFileName += "Intermediate.Ptt";
435 else if((psiMin == 1.5)&&(psiMax == 2.5)) newFileName += "OutOfPlane.Ptt";
436 else if((psiMin == 2.5)&&(psiMax == 3.5)) newFileName += "Rest.PttFrom";
437 else newFileName += "All.PttFrom";
438 newFileName += Form("%.1f",ptTriggerMin); newFileName += "To";
439 newFileName += Form("%.1f",ptTriggerMax); newFileName += "PtaFrom";
440 newFileName += Form("%.1f",ptAssociatedMin); newFileName += "To";
441 newFileName += Form("%.1f",ptAssociatedMax);
442 if(k2pMethod) newFileName += "_2pMethod";
443 newFileName += ".root";
445 TFile *fOutput = new TFile(newFileName.Data(),"recreate");
447 /*hP->Write(); hN->Write();
448 hPN->Write(); hNP->Write();
449 hPP->Write(); hNN->Write();
450 hPShuffled->Write(); hNShuffled->Write();
451 hPNShuffled->Write(); hNPShuffled->Write();
452 hPPShuffled->Write(); hNNShuffled->Write();
453 hPMixed->Write(); hNMixed->Write();
454 hPNMixed->Write(); hNPMixed->Write();
455 hPPMixed->Write(); hNNMixed->Write();*/
456 gHistBalanceFunction->Write();
457 if(listBFShuffled) gHistBalanceFunctionShuffled->Write();
459 gHistBalanceFunctionMixed->Write();
460 gHistBalanceFunctionSubtracted->Write();
465 //____________________________________________________________//
466 void fitbalanceFunction(Int_t gCentrality = 1,
467 Double_t psiMin = -0.5, Double_t psiMax = 3.5,
468 Double_t ptTriggerMin = -1.,
469 Double_t ptTriggerMax = -1.,
470 Double_t ptAssociatedMin = -1.,
471 Double_t ptAssociatedMax = -1.,
474 cout<<"FITTING FUNCTION"<<endl;
476 //near side peak: [1]*TMath::Exp(-TMath::Power((0.5*TMath::Power((x/[2]),2)+0.5*TMath::Power((y/[3]),2)),[4]))
477 //away side ridge: [5]*TMath::Exp(-TMath::Power((0.5*TMath::Power(((y-TMath::Pi())/[6]),2)),[7]))
478 //longitudinal ridge: [8]*TMath::Exp(-TMath::Power((0.5*TMath::Power((x/[9]),2)),[10]))
479 //wing structures: [11]*TMath::Power(x,2)
480 //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))
481 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+[17])/[9]),2)),[10]))+[8]*TMath::Exp(-TMath::Power((0.5*TMath::Power(((x-[17])/[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.);
482 gFitFunction->SetName("gFitFunction");
486 gFitFunction->SetParName(0,"N1");
488 gFitFunction->SetParName(1,"N_{near side}");gFitFunction->FixParameter(1,0.0);
489 gFitFunction->SetParName(2,"Sigma_{near side}(delta eta)"); gFitFunction->FixParameter(2,0.0);
490 gFitFunction->SetParName(3,"Sigma_{near side}(delta phi)"); gFitFunction->FixParameter(3,0.0);
491 gFitFunction->SetParName(4,"Exponent_{near side}"); gFitFunction->FixParameter(4,1.0);
493 gFitFunction->SetParName(5,"N_{away side}"); gFitFunction->FixParameter(5,0.0);
494 gFitFunction->SetParName(6,"Sigma_{away side}(delta phi)"); gFitFunction->FixParameter(6,0.0);
495 gFitFunction->SetParName(7,"Exponent_{away side}"); gFitFunction->FixParameter(7,1.0);
497 gFitFunction->SetParName(8,"N_{long. ridge}"); gFitFunction->SetParameter(8,0.05);//
498 gFitFunction->FixParameter(8,0.0);
499 gFitFunction->SetParName(9,"Sigma_{long. ridge}(delta eta)"); gFitFunction->FixParameter(9,0.0);
500 gFitFunction->SetParName(10,"Exponent_{long. ridge}"); gFitFunction->FixParameter(10,1.0);
502 gFitFunction->SetParName(11,"N_{wing}"); gFitFunction->FixParameter(11,0.0);
504 gFitFunction->SetParName(12,"N_{flow}"); gFitFunction->FixParameter(12,0.0);
505 gFitFunction->SetParName(13,"V1"); gFitFunction->FixParameter(13,0.0);
506 gFitFunction->SetParName(14,"V2"); gFitFunction->FixParameter(14,0.0);
507 gFitFunction->SetParName(15,"V3"); gFitFunction->FixParameter(15,0.0);
508 gFitFunction->SetParName(16,"V4"); gFitFunction->FixParameter(16,0.0);
509 gFitFunction->SetParName(17,"Offset"); gFitFunction->FixParameter(17,0.0);
511 // just release the near side peak
512 gFitFunction->ReleaseParameter(0);gFitFunction->SetParameter(0,1.0);
513 gFitFunction->ReleaseParameter(1);gFitFunction->SetParameter(1,0.3);gFitFunction->SetParLimits(1,0.0,100);
514 gFitFunction->ReleaseParameter(2);gFitFunction->SetParameter(2,0.3);gFitFunction->SetParLimits(2,0.05,0.7);
515 gFitFunction->ReleaseParameter(3);gFitFunction->SetParameter(3,0.3);gFitFunction->SetParLimits(3,0.05,1.7);
517 gHist->Fit("gFitFunction","nm");
520 //Cloning the histogram
521 TString histoName = gHist->GetName(); histoName += "Fit";
522 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());
523 TH2D *gHistResidual = dynamic_cast<TH2D *>(gHist->Clone());
524 gHistResidual->SetName("gHistResidual");
525 gHistResidual->Sumw2();
527 for(Int_t iBinDeltaEta = 1; iBinDeltaEta <= gHist->GetNbinsX(); iBinDeltaEta++) {
528 for(Int_t iBinDeltaPhi = 1; iBinDeltaPhi <= gHist->GetNbinsY(); iBinDeltaPhi++) {
529 gHistFit->SetBinContent(iBinDeltaEta,iBinDeltaPhi,gFitFunction->Eval(gHist->GetXaxis()->GetBinCenter(iBinDeltaEta),gHist->GetYaxis()->GetBinCenter(iBinDeltaPhi)));
532 gHistResidual->Add(gHistFit,-1);
534 //Write to output file
535 TString newFileName = "balanceFunctionFit";
536 newFileName += ".Centrality";
537 newFileName += gCentrality; newFileName += ".Psi";
538 if((psiMin == -0.5)&&(psiMax == 0.5)) newFileName += "InPlane.Ptt";
539 else if((psiMin == 0.5)&&(psiMax == 1.5)) newFileName += "Intermediate.Ptt";
540 else if((psiMin == 1.5)&&(psiMax == 2.5)) newFileName += "OutOfPlane.Ptt";
541 else if((psiMin == 2.5)&&(psiMax == 3.5)) newFileName += "Rest.PttFrom";
542 else newFileName += "All.PttFrom";
543 newFileName += Form("%.1f",ptTriggerMin); newFileName += "To";
544 newFileName += Form("%.1f",ptTriggerMax); newFileName += "PtaFrom";
545 newFileName += Form("%.1f",ptAssociatedMin); newFileName += "To";
546 newFileName += Form("%.1f",ptAssociatedMax);
547 newFileName += ".root";
548 TFile *newFile = TFile::Open(newFileName.Data(),"recreate");
551 gHistResidual->Write();
552 gFitFunction->Write();
560 //____________________________________________________________//
561 void drawBFPsi2D(const char* lhcPeriod = "LHC11h",
562 Int_t gCentrality = 1,
563 Bool_t kShowShuffled = kFALSE,
564 Bool_t kShowMixed = kFALSE,
565 Double_t psiMin = -0.5, Double_t psiMax = 0.5,
566 Double_t ptTriggerMin = -1.,
567 Double_t ptTriggerMax = -1.,
568 Double_t ptAssociatedMin = -1.,
569 Double_t ptAssociatedMax = -1.) {
570 //Macro that draws the BF distributions for each centrality bin
571 //for reaction plane dependent analysis
572 //Author: Panos.Christakoglou@nikhef.nl
573 TGaxis::SetMaxDigits(3);
576 TString filename = lhcPeriod; filename +="/PttFrom";
577 filename += Form("%.1f",ptTriggerMin); filename += "To";
578 filename += Form("%.1f",ptTriggerMax); filename += "PtaFrom";
579 filename += Form("%.1f",ptAssociatedMin); filename += "To";
580 filename += Form("%.1f",ptAssociatedMax);
581 filename += "/balanceFunction2D.Centrality";
582 filename += gCentrality; filename += ".Psi";
583 if((psiMin == -0.5)&&(psiMax == 0.5)) filename += "InPlane.Ptt";
584 else if((psiMin == 0.5)&&(psiMax == 1.5)) filename += "Intermediate.Ptt";
585 else if((psiMin == 1.5)&&(psiMax == 2.5)) filename += "OutOfPlane.Ptt";
586 else if((psiMin == 2.5)&&(psiMax == 3.5)) filename += "Rest.Ptt";
587 else filename += "All.PttFrom";
588 filename += Form("%.1f",ptTriggerMin); filename += "To";
589 filename += Form("%.1f",ptTriggerMax); filename += "PtaFrom";
590 filename += Form("%.1f",ptAssociatedMin); filename += "To";
591 filename += Form("%.1f",ptAssociatedMax); filename += ".root";
594 TFile *f = TFile::Open(filename.Data());
595 if((!f)||(!f->IsOpen())) {
596 Printf("The file %s is not found. Aborting...",filename);
601 //Raw balance function
602 TH1D *gHistBalanceFunction = dynamic_cast<TH1D *>(f->Get("gHistBalanceFunction"));
603 gHistBalanceFunction->SetStats(kFALSE);
604 gHistBalanceFunction->GetXaxis()->SetNdivisions(10);
605 gHistBalanceFunction->GetYaxis()->SetNdivisions(10);
606 gHistBalanceFunction->GetZaxis()->SetNdivisions(10);
607 gHistBalanceFunction->GetXaxis()->SetTitleOffset(1.3);
608 gHistBalanceFunction->GetYaxis()->SetTitleOffset(1.3);
609 gHistBalanceFunction->GetZaxis()->SetTitleOffset(1.3);
610 gHistBalanceFunction->GetXaxis()->SetTitle("#Delta #eta");
611 gHistBalanceFunction->GetYaxis()->SetTitle("#Delta #varphi (deg.)");
612 gHistBalanceFunction->GetZaxis()->SetTitle("B(#Delta #eta, #Delta #varphi)");
614 //Shuffled balance function
616 TH1D *gHistBalanceFunctionShuffled = dynamic_cast<TH1D *>(f->Get("gHistBalanceFunctionShuffled"));
617 gHistBalanceFunctionShuffled->SetStats(kFALSE);
618 gHistBalanceFunctionShuffled->GetXaxis()->SetNdivisions(10);
619 gHistBalanceFunctionShuffled->GetYaxis()->SetNdivisions(10);
620 gHistBalanceFunctionShuffled->GetZaxis()->SetNdivisions(10);
621 gHistBalanceFunctionShuffled->GetXaxis()->SetTitleOffset(1.3);
622 gHistBalanceFunctionShuffled->GetYaxis()->SetTitleOffset(1.3);
623 gHistBalanceFunctionShuffled->GetZaxis()->SetTitleOffset(1.3);
624 gHistBalanceFunctionShuffled->GetXaxis()->SetTitle("#Delta #eta");
625 gHistBalanceFunctionShuffled->GetYaxis()->SetTitle("#Delta #varphi (deg.)");
626 gHistBalanceFunctionShuffled->GetZaxis()->SetTitle("B(#Delta #eta, #Delta #varphi)");
629 //Mixed balance function
631 TH1D *gHistBalanceFunctionMixed = dynamic_cast<TH1D *>(f->Get("gHistBalanceFunctionMixed"));
632 gHistBalanceFunctionMixed->SetStats(kFALSE);
633 gHistBalanceFunctionMixed->GetXaxis()->SetNdivisions(10);
634 gHistBalanceFunctionMixed->GetYaxis()->SetNdivisions(10);
635 gHistBalanceFunctionMixed->GetZaxis()->SetNdivisions(10);
636 gHistBalanceFunctionMixed->GetXaxis()->SetTitleOffset(1.3);
637 gHistBalanceFunctionMixed->GetYaxis()->SetTitleOffset(1.3);
638 gHistBalanceFunctionMixed->GetZaxis()->SetTitleOffset(1.3);
639 gHistBalanceFunctionMixed->GetXaxis()->SetTitle("#Delta #eta");
640 gHistBalanceFunctionMixed->GetYaxis()->SetTitle("#Delta #varphi (deg.)");
641 gHistBalanceFunctionMixed->GetZaxis()->SetTitle("B(#Delta #eta, #Delta #varphi)");
644 //Subtracted balance function
646 TH1D *gHistBalanceFunctionSubtracted = dynamic_cast<TH1D *>(f->Get("gHistBalanceFunctionSubtracted"));
647 gHistBalanceFunctionSubtracted->SetStats(kFALSE);
648 gHistBalanceFunctionSubtracted->GetXaxis()->SetNdivisions(10);
649 gHistBalanceFunctionSubtracted->GetYaxis()->SetNdivisions(10);
650 gHistBalanceFunctionSubtracted->GetZaxis()->SetNdivisions(10);
651 gHistBalanceFunctionSubtracted->GetXaxis()->SetTitleOffset(1.3);
652 gHistBalanceFunctionSubtracted->GetYaxis()->SetTitleOffset(1.3);
653 gHistBalanceFunctionSubtracted->GetZaxis()->SetTitleOffset(1.3);
654 gHistBalanceFunctionSubtracted->GetXaxis()->SetTitle("#Delta #eta");
655 gHistBalanceFunctionSubtracted->GetYaxis()->SetTitle("#Delta #varphi (deg.)");
656 gHistBalanceFunctionSubtracted->GetZaxis()->SetTitle("B(#Delta #eta, #Delta #varphi)");
661 TString centralityLatex = "Centrality: ";
662 centralityLatex += centralityArray[gCentrality-1];
663 centralityLatex += "%";
666 if((psiMin == -0.5)&&(psiMax == 0.5))
667 psiLatex = " -7.5^{o} < #varphi - #Psi_{2} < 7.5^{o}";
668 else if((psiMin == 0.5)&&(psiMax == 1.5))
669 psiLatex = " 37.5^{o} < #varphi - #Psi_{2} < 52.5^{o}";
670 else if((psiMin == 1.5)&&(psiMax == 2.5))
671 psiLatex = " 82.5^{o} < #varphi - #Psi_{2} < 97.5^{o}";
673 psiLatex = " 0^{o} < #varphi - #Psi_{2} < 180^{o}";
675 TString pttLatex = Form("%.1f",ptTriggerMin);
676 pttLatex += " < p_{T}^{t} < "; pttLatex += Form("%.1f",ptTriggerMax);
677 pttLatex += " GeV/c";
679 TString ptaLatex = Form("%.1f",ptAssociatedMin);
680 ptaLatex += " < p_{T}^{a} < "; ptaLatex += Form("%.1f",ptAssociatedMax);
681 ptaLatex += " GeV/c";
683 TLatex *latexInfo1 = new TLatex();
684 latexInfo1->SetNDC();
685 latexInfo1->SetTextSize(0.045);
686 latexInfo1->SetTextColor(1);
689 TCanvas *c1 = new TCanvas("c1","Raw balance function 2D",0,0,600,500);
690 c1->SetFillColor(10); c1->SetHighLightColor(10);
691 c1->SetLeftMargin(0.17); c1->SetTopMargin(0.05);
692 gHistBalanceFunction->SetTitle("Raw balance function");
693 gHistBalanceFunction->GetYaxis()->SetTitleOffset(1.4);
694 gHistBalanceFunction->GetYaxis()->SetNdivisions(10);
695 gHistBalanceFunction->GetXaxis()->SetNdivisions(10);
696 gHistBalanceFunction->DrawCopy("lego2");
697 gPad->SetTheta(30); // default is 30
698 gPad->SetPhi(-60); // default is 30
701 latexInfo1->DrawLatex(0.64,0.88,centralityLatex.Data());
702 latexInfo1->DrawLatex(0.64,0.82,psiLatex.Data());
703 latexInfo1->DrawLatex(0.64,0.76,pttLatex.Data());
704 latexInfo1->DrawLatex(0.64,0.70,ptaLatex.Data());
707 TCanvas *c2 = new TCanvas("c2","Shuffled balance function 2D",100,100,600,500);
708 c2->SetFillColor(10); c2->SetHighLightColor(10);
709 c2->SetLeftMargin(0.17); c2->SetTopMargin(0.05);
710 gHistBalanceFunctionShuffled->SetTitle("Shuffled events");
711 gHistBalanceFunctionShuffled->GetYaxis()->SetTitleOffset(1.4);
712 gHistBalanceFunctionShuffled->GetYaxis()->SetNdivisions(10);
713 gHistBalanceFunctionShuffled->GetXaxis()->SetNdivisions(10);
714 gHistBalanceFunctionShuffled->DrawCopy("lego2");
715 gPad->SetTheta(30); // default is 30
716 gPad->SetPhi(-60); // default is 30
719 latexInfo1->DrawLatex(0.64,0.88,centralityLatex.Data());
720 latexInfo1->DrawLatex(0.64,0.82,psiLatex.Data());
721 latexInfo1->DrawLatex(0.64,0.76,pttLatex.Data());
722 latexInfo1->DrawLatex(0.64,0.70,ptaLatex.Data());
726 TCanvas *c3 = new TCanvas("c3","Mixed balance function 2D",200,200,600,500);
727 c3->SetFillColor(10); c3->SetHighLightColor(10);
728 c3->SetLeftMargin(0.17); c3->SetTopMargin(0.05);
729 gHistBalanceFunctionMixed->SetTitle("Mixed events");
730 gHistBalanceFunctionMixed->GetYaxis()->SetTitleOffset(1.4);
731 gHistBalanceFunctionMixed->GetYaxis()->SetNdivisions(10);
732 gHistBalanceFunctionMixed->GetXaxis()->SetNdivisions(10);
733 gHistBalanceFunctionMixed->DrawCopy("lego2");
734 gPad->SetTheta(30); // default is 30
735 gPad->SetPhi(-60); // default is 30
738 latexInfo1->DrawLatex(0.64,0.88,centralityLatex.Data());
739 latexInfo1->DrawLatex(0.64,0.82,psiLatex.Data());
740 latexInfo1->DrawLatex(0.64,0.76,pttLatex.Data());
741 latexInfo1->DrawLatex(0.64,0.70,ptaLatex.Data());
745 TCanvas *c4 = new TCanvas("c4","Subtracted balance function 2D",300,300,600,500);
746 c4->SetFillColor(10); c4->SetHighLightColor(10);
747 c4->SetLeftMargin(0.17); c4->SetTopMargin(0.05);
748 gHistBalanceFunctionSubtracted->SetTitle("Subtracted balance function");
749 gHistBalanceFunctionSubtracted->GetYaxis()->SetTitleOffset(1.4);
750 gHistBalanceFunctionSubtracted->GetYaxis()->SetNdivisions(10);
751 gHistBalanceFunctionSubtracted->GetXaxis()->SetNdivisions(10);
752 gHistBalanceFunctionSubtracted->DrawCopy("lego2");
753 gPad->SetTheta(30); // default is 30
754 gPad->SetPhi(-60); // default is 30
757 latexInfo1->DrawLatex(0.64,0.88,centralityLatex.Data());
758 latexInfo1->DrawLatex(0.64,0.82,psiLatex.Data());
759 latexInfo1->DrawLatex(0.64,0.76,pttLatex.Data());
760 latexInfo1->DrawLatex(0.64,0.70,ptaLatex.Data());