]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGCF/EBYE/macros/drawCorrelationFunctionPsi.C
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
1) Bugfix for correlation function calculation (when averaging over several centralit...
[u/mrichter/AliRoot.git] / PWGCF / EBYE / macros / drawCorrelationFunctionPsi.C
CommitLineData
e5c5c33f 1const Int_t numberOfCentralityBins = 12;
41f2bc59 2//TString centralityArray[numberOfCentralityBins] = {"0-4","4-5","6-14","30-40","40-50","50-60","60-70","70-80","0-100","0-1","1-2","2-3"};
09b14834 3TString centralityArray[numberOfCentralityBins] = {"0-100","10-20","20-30","30-40","40-50","50-60","60-70","70-80","0-100","0-1","1-2","2-3"};
a38fd7f3 4
5const Int_t gRebin = 1;
07d0a35c 6
2b4abbc5 7void drawCorrelationFunctionPsi(const char* filename = "AnalysisResultsPsi_train_06feb.root",
52daf7b2 8 Int_t gCentrality = 1,
5de9ad1a 9 Int_t gBit = -1,
10 const char* gCentralityEstimator = 0x0,
52daf7b2 11 Bool_t kShowShuffled = kFALSE,
12 Bool_t kShowMixed = kTRUE,
6acdbcb2 13 Double_t psiMin = -0.5,
5de9ad1a 14 Double_t psiMax = 3.5,
20006629 15 Double_t vertexZMin = -10.,
16 Double_t vertexZMax = 10.,
6acdbcb2 17 Double_t ptTriggerMin = -1.,
18 Double_t ptTriggerMax = -1.,
19 Double_t ptAssociatedMin = -1.,
742af4bd 20 Double_t ptAssociatedMax = -1.,
2b4abbc5 21 Bool_t normToTrig = kTRUE,
20006629 22 Bool_t kUseVzBinning = kTRUE,
742af4bd 23 Int_t rebinEta = 1,
32e94079 24 Int_t rebinPhi = 1,
dd46a0c3 25 TString eventClass = "EventPlane", //Can be "EventPlane", "Centrality", "Multiplicity"
633996e7 26 Bool_t bToy = kFALSE,
27 TString listNameAdd = "")
32e94079 28{
a38fd7f3 29 //Macro that draws the correlation functions from the balance function
30 //analysis vs the reaction plane
31 //Author: Panos.Christakoglou@nikhef.nl
2b4abbc5 32 //gROOT->LoadMacro("~/SetPlotStyle.C");
33 //SetPlotStyle();
52daf7b2 34 gStyle->SetPalette(1,0);
35
a38fd7f3 36 //Load the PWG2ebye library
2b4abbc5 37 gSystem->Load("libCore.so");
38 gSystem->Load("libGeom.so");
39 gSystem->Load("libVMC.so");
40 gSystem->Load("libPhysics.so");
41 gSystem->Load("libTree.so");
42 gSystem->Load("libSTEERBase.so");
43 gSystem->Load("libESD.so");
44 gSystem->Load("libAOD.so");
45
a38fd7f3 46 gSystem->Load("libANALYSIS.so");
47 gSystem->Load("libANALYSISalice.so");
48 gSystem->Load("libEventMixing.so");
49 gSystem->Load("libCORRFW.so");
50 gSystem->Load("libPWGTools.so");
51 gSystem->Load("libPWGCFebye.so");
52
53 //Prepare the objects and return them
3ce45c12 54 TList *listQA = NULL;
633996e7 55 TList *list = GetListOfObjects(filename,gCentrality,gBit,gCentralityEstimator,0,bToy,listNameAdd);
52daf7b2 56 TList *listShuffled = NULL;
633996e7 57 if(kShowShuffled) listShuffled = GetListOfObjects(filename,gCentrality,gBit,gCentralityEstimator,1,bToy,listNameAdd);
52daf7b2 58 TList *listMixed = NULL;
633996e7 59 if(kShowMixed) listMixed = GetListOfObjects(filename,gCentrality,gBit,gCentralityEstimator,2,bToy,listNameAdd);
52daf7b2 60
3ce45c12 61 // else get the QA histograms (for convolution)
62 else{
63 //Open the file again
64 TFile *f = TFile::Open(filename,"UPDATE");
65 if((!f)||(!f->IsOpen())) {
66 Printf("The file %s is not found.",filename);
67 }
68
2b4abbc5 69
3ce45c12 70 TDirectoryFile *dir = dynamic_cast<TDirectoryFile *>(f->Get("PWGCFEbyE.outputBalanceFunctionPsiAnalysis"));
71 if(!dir) {
72 Printf("The TDirectoryFile is not found.",filename);
73 }
74
75 TString listQAName = "listQAPsi_";
76
77 listQAName += centralityArray[gCentrality-1];
78 if(gBit > -1) {
79 listQAName += "_Bit"; listQAName += gBit; }
80 if(gCentralityEstimator) {
81 listQAName += "_"; listQAName += gCentralityEstimator;}
a9737921 82 listQAName += listNameAdd;
83
3ce45c12 84
85 listQA = (TList*)dir->Get(Form("%s",listQAName.Data()));
86 if(!listQA) {
a9737921 87 Printf("TList %s not found!!!",listQAName.Data());
3ce45c12 88 }
89 }
90
a38fd7f3 91 if(!list) {
92 Printf("The TList object was not created");
93 return;
94 }
95 else
3ce45c12 96 draw(list,listShuffled,listMixed,listQA,
20006629 97 gCentralityEstimator,gCentrality,psiMin,psiMax,vertexZMin,vertexZMax,
98 ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax,normToTrig,kUseVzBinning,rebinEta,rebinPhi,eventClass);
a38fd7f3 99}
100
101//______________________________________________________//
52daf7b2 102TList *GetListOfObjects(const char* filename,
103 Int_t gCentrality,
5de9ad1a 104 Int_t gBit,
105 const char *gCentralityEstimator,
dd46a0c3 106 Int_t kData = 1,
633996e7 107 Bool_t bToy = kFALSE,
108 TString listNameAdd = "") {
a38fd7f3 109 //Get the TList objects (QA, bf, bf shuffled)
a38fd7f3 110 TList *listBF = 0x0;
a38fd7f3 111
112 //Open the file
ff0805a7 113 TFile *f = TFile::Open(filename,"UPDATE");
a38fd7f3 114 if((!f)||(!f->IsOpen())) {
115 Printf("The file %s is not found. Aborting...",filename);
116 return listBF;
117 }
118 //f->ls();
119
6acdbcb2 120 TDirectoryFile *dir = dynamic_cast<TDirectoryFile *>(f->Get("PWGCFEbyE.outputBalanceFunctionPsiAnalysis"));
a38fd7f3 121 if(!dir) {
122 Printf("The TDirectoryFile is not found. Aborting...",filename);
123 return listBF;
124 }
125 //dir->ls();
2b4abbc5 126
52daf7b2 127 TString listBFName;
128 if(kData == 0) {
129 //cout<<"no shuffling - no mixing"<<endl;
dd46a0c3 130 listBFName = "listBFPsi";
52daf7b2 131 }
132 else if(kData == 1) {
133 //cout<<"shuffling - no mixing"<<endl;
dd46a0c3 134 listBFName = "listBFPsiShuffled";
52daf7b2 135 }
136 else if(kData == 2) {
137 //cout<<"no shuffling - mixing"<<endl;
dd46a0c3 138 listBFName = "listBFPsiMixed";
139 }
140
141 // different list names in case of toy model
142 if(!bToy){
143 listBFName += "_";
144 listBFName += centralityArray[gCentrality-1];
145 if(gBit > -1) {
146 listBFName += "_Bit"; listBFName += gBit; }
147 if(gCentralityEstimator) {
148 listBFName += "_"; listBFName += gCentralityEstimator;}
149 }
150 else{
151 listBFName.ReplaceAll("Psi","");
52daf7b2 152 }
633996e7 153 listBFName += listNameAdd;
a38fd7f3 154
5365d1d7 155 // histograms were already retrieved (in first iteration)
156 if(dir->Get(Form("%s_histograms",listBFName.Data()))){
157 listBF = dynamic_cast<TList *>(dir->Get(Form("%s_histograms",listBFName.Data())));
a38fd7f3 158 }
a38fd7f3 159
5365d1d7 160 // histograms were not yet retrieved (this is the first iteration)
161 else{
162
163 listBF = dynamic_cast<TList *>(dir->Get(listBFName.Data()));
633996e7 164 if(!listBF) dir->ls();
5365d1d7 165 cout<<"======================================================="<<endl;
32e94079 166 cout<<"List name (control): "<<listBFName.Data()<<endl;
17d10bdd 167 cout<<"List name: "<<listBF->GetName()<<endl;
5365d1d7 168 //listBF->ls();
a38fd7f3 169
5365d1d7 170 //Get the histograms
171 TString histoName;
17d10bdd 172 if(kData == 0)
173 histoName = "fHistP";
5365d1d7 174 else if(kData == 1)
17d10bdd 175 histoName = "fHistP_shuffle";
5365d1d7 176 else if(kData == 2)
17d10bdd 177 histoName = "fHistP";
178 if(gCentralityEstimator)
179 histoName += gCentralityEstimator;
5365d1d7 180 AliTHn *fHistP = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
181 if(!fHistP) {
182 Printf("fHistP %s not found!!!",histoName.Data());
183 break;
184 }
185 fHistP->FillParent(); fHistP->DeleteContainers();
186
187 if(kData == 0)
17d10bdd 188 histoName = "fHistN";
5365d1d7 189 if(kData == 1)
17d10bdd 190 histoName = "fHistN_shuffle";
5365d1d7 191 if(kData == 2)
17d10bdd 192 histoName = "fHistN";
193 if(gCentralityEstimator)
194 histoName += gCentralityEstimator;
5365d1d7 195 AliTHn *fHistN = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
196 if(!fHistN) {
197 Printf("fHistN %s not found!!!",histoName.Data());
198 break;
199 }
200 fHistN->FillParent(); fHistN->DeleteContainers();
201
202 if(kData == 0)
17d10bdd 203 histoName = "fHistPN";
5365d1d7 204 if(kData == 1)
17d10bdd 205 histoName = "fHistPN_shuffle";
5365d1d7 206 if(kData == 2)
17d10bdd 207 histoName = "fHistPN";
208 if(gCentralityEstimator)
209 histoName += gCentralityEstimator;
5365d1d7 210 AliTHn *fHistPN = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
211 if(!fHistPN) {
212 Printf("fHistPN %s not found!!!",histoName.Data());
213 break;
214 }
215 fHistPN->FillParent(); fHistPN->DeleteContainers();
216
217 if(kData == 0)
17d10bdd 218 histoName = "fHistNP";
5365d1d7 219 if(kData == 1)
17d10bdd 220 histoName = "fHistNP_shuffle";
5365d1d7 221 if(kData == 2)
17d10bdd 222 histoName = "fHistNP";
223 if(gCentralityEstimator)
224 histoName += gCentralityEstimator;
5365d1d7 225 AliTHn *fHistNP = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
226 if(!fHistNP) {
227 Printf("fHistNP %s not found!!!",histoName.Data());
228 break;
229 }
230 fHistNP->FillParent(); fHistNP->DeleteContainers();
231
232 if(kData == 0)
17d10bdd 233 histoName = "fHistPP";
5365d1d7 234 if(kData == 1)
17d10bdd 235 histoName = "fHistPP_shuffle";
5365d1d7 236 if(kData == 2)
17d10bdd 237 histoName = "fHistPP";
238 if(gCentralityEstimator)
239 histoName += gCentralityEstimator;
5365d1d7 240 AliTHn *fHistPP = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
241 if(!fHistPP) {
242 Printf("fHistPP %s not found!!!",histoName.Data());
243 break;
244 }
245 fHistPP->FillParent(); fHistPP->DeleteContainers();
246
247 if(kData == 0)
17d10bdd 248 histoName = "fHistNN";
5365d1d7 249 if(kData == 1)
17d10bdd 250 histoName = "fHistNN_shuffle";
5365d1d7 251 if(kData == 2)
17d10bdd 252 histoName = "fHistNN";
253 if(gCentralityEstimator)
254 histoName += gCentralityEstimator;
5365d1d7 255 AliTHn *fHistNN = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
256 if(!fHistNN) {
257 Printf("fHistNN %s not found!!!",histoName.Data());
258 break;
259 }
260 fHistNN->FillParent(); fHistNN->DeleteContainers();
17d10bdd 261
5365d1d7 262 dir->cd();
263 listBF->Write(Form("%s_histograms",listBFName.Data()), TObject::kSingleKey);
264
265 }// first iteration
a38fd7f3 266
ff0805a7 267 f->Close();
5365d1d7 268
a38fd7f3 269 return listBF;
270}
271
272//______________________________________________________//
20006629 273void draw(TList *list, TList *listBFShuffled, TList *listBFMixed,
274 TList *listQA,
17d10bdd 275 const char *gCentralityEstimator,
20006629 276 Int_t gCentrality, Double_t psiMin, Double_t psiMax,
277 Double_t vertexZMin,Double_t vertexZMax,
6acdbcb2 278 Double_t ptTriggerMin, Double_t ptTriggerMax,
742af4bd 279 Double_t ptAssociatedMin, Double_t ptAssociatedMax,
20006629 280 Bool_t normToTrig,
281 Bool_t kUseVzBinning,
282 Int_t rebinEta, Int_t rebinPhi,TString eventClass) {
a38fd7f3 283 //Draws the correlation functions for every centrality bin
52daf7b2 284 //(+-), (-+), (++), (--)
a38fd7f3 285 AliTHn *hP = NULL;
286 AliTHn *hN = NULL;
287 AliTHn *hPN = NULL;
288 AliTHn *hNP = NULL;
289 AliTHn *hPP = NULL;
290 AliTHn *hNN = NULL;
291
17d10bdd 292 TString gHistPname = "fHistP";
293 if(gCentralityEstimator) gHistPname += gCentralityEstimator;
294 TString gHistNname = "fHistN";
295 if(gCentralityEstimator) gHistNname += gCentralityEstimator;
296 TString gHistPNname = "fHistPN";
297 if(gCentralityEstimator) gHistPNname += gCentralityEstimator;
298 TString gHistNPname = "fHistNP";
299 if(gCentralityEstimator) gHistNPname += gCentralityEstimator;
300 TString gHistPPname = "fHistPP";
301 if(gCentralityEstimator) gHistPPname += gCentralityEstimator;
302 TString gHistNNname = "fHistNN";
303 if(gCentralityEstimator) gHistNNname += gCentralityEstimator;
304
305 hP = (AliTHn*) list->FindObject(gHistPname.Data());
306 hN = (AliTHn*) list->FindObject(gHistNname.Data());
307 hPN = (AliTHn*) list->FindObject(gHistPNname.Data());
308 hNP = (AliTHn*) list->FindObject(gHistNPname.Data());
309 hPP = (AliTHn*) list->FindObject(gHistPPname.Data());
310 hNN = (AliTHn*) list->FindObject(gHistNNname.Data());
32e94079 311 hNN->Print();
312
a38fd7f3 313
314 //Create the AliBalancePsi object and fill it with the AliTHn objects
315 AliBalancePsi *b = new AliBalancePsi();
32e94079 316 b->SetEventClass(eventClass);
a38fd7f3 317 b->SetHistNp(hP);
318 b->SetHistNn(hN);
319 b->SetHistNpn(hPN);
320 b->SetHistNnp(hNP);
321 b->SetHistNpp(hPP);
322 b->SetHistNnn(hNN);
20006629 323 if(kUseVzBinning) b->SetVertexZBinning(kTRUE);
52daf7b2 324
325 //balance function shuffling
326 AliTHn *hPShuffled = NULL;
327 AliTHn *hNShuffled = NULL;
328 AliTHn *hPNShuffled = NULL;
329 AliTHn *hNPShuffled = NULL;
330 AliTHn *hPPShuffled = NULL;
331 AliTHn *hNNShuffled = NULL;
332 if(listBFShuffled) {
333 //listBFShuffled->ls();
334
17d10bdd 335 gHistPname = "fHistP_shuffle";
336 if(gCentralityEstimator) gHistPname += gCentralityEstimator;
337 gHistNname = "fHistN_shuffle";
338 if(gCentralityEstimator) gHistNname += gCentralityEstimator;
339 gHistPNname = "fHistPN_shuffle";
340 if(gCentralityEstimator) gHistPNname += gCentralityEstimator;
341 gHistNPname = "fHistNP_shuffle";
342 if(gCentralityEstimator) gHistNPname += gCentralityEstimator;
343 gHistPPname = "fHistPP_shuffle";
344 if(gCentralityEstimator) gHistPPname += gCentralityEstimator;
345 gHistNNname = "fHistNN_shuffle";
346 if(gCentralityEstimator) gHistNNname += gCentralityEstimator;
347
348 hPShuffled = (AliTHn*) listBFShuffled->FindObject(gHistPname.Data());
52daf7b2 349 hPShuffled->SetName("gHistPShuffled");
17d10bdd 350 hNShuffled = (AliTHn*) listBFShuffled->FindObject(gHistNname.Data());
52daf7b2 351 hNShuffled->SetName("gHistNShuffled");
17d10bdd 352 hPNShuffled = (AliTHn*) listBFShuffled->FindObject(gHistPNname.Data());
52daf7b2 353 hPNShuffled->SetName("gHistPNShuffled");
17d10bdd 354 hNPShuffled = (AliTHn*) listBFShuffled->FindObject(gHistNPname.Data());
52daf7b2 355 hNPShuffled->SetName("gHistNPShuffled");
17d10bdd 356 hPPShuffled = (AliTHn*) listBFShuffled->FindObject(gHistPPname.Data());
52daf7b2 357 hPPShuffled->SetName("gHistPPShuffled");
17d10bdd 358 hNNShuffled = (AliTHn*) listBFShuffled->FindObject(gHistNNname.Data());
52daf7b2 359 hNNShuffled->SetName("gHistNNShuffled");
360
361 AliBalancePsi *bShuffled = new AliBalancePsi();
32e94079 362 bShuffled->SetEventClass(eventClass);
52daf7b2 363 bShuffled->SetHistNp(hPShuffled);
364 bShuffled->SetHistNn(hNShuffled);
365 bShuffled->SetHistNpn(hPNShuffled);
366 bShuffled->SetHistNnp(hNPShuffled);
367 bShuffled->SetHistNpp(hPPShuffled);
368 bShuffled->SetHistNnn(hNNShuffled);
20006629 369 if(kUseVzBinning) bShuffled->SetVertexZBinning(kTRUE);
52daf7b2 370 }
371
372 //balance function mixing
373 AliTHn *hPMixed = NULL;
374 AliTHn *hNMixed = NULL;
375 AliTHn *hPNMixed = NULL;
376 AliTHn *hNPMixed = NULL;
377 AliTHn *hPPMixed = NULL;
378 AliTHn *hNNMixed = NULL;
379
380 if(listBFMixed) {
381 //listBFMixed->ls();
382
17d10bdd 383 gHistPname = "fHistP";
384 if(gCentralityEstimator) gHistPname += gCentralityEstimator;
385 gHistNname = "fHistN";
386 if(gCentralityEstimator) gHistNname += gCentralityEstimator;
387 gHistPNname = "fHistPN";
388 if(gCentralityEstimator) gHistPNname += gCentralityEstimator;
389 gHistNPname = "fHistNP";
390 if(gCentralityEstimator) gHistNPname += gCentralityEstimator;
391 gHistPPname = "fHistPP";
392 if(gCentralityEstimator) gHistPPname += gCentralityEstimator;
393 gHistNNname = "fHistNN";
394 if(gCentralityEstimator) gHistNNname += gCentralityEstimator;
395 hPMixed = (AliTHn*) listBFMixed->FindObject(gHistPname.Data());
52daf7b2 396 hPMixed->SetName("gHistPMixed");
17d10bdd 397 hNMixed = (AliTHn*) listBFMixed->FindObject(gHistNname.Data());
52daf7b2 398 hNMixed->SetName("gHistNMixed");
17d10bdd 399 hPNMixed = (AliTHn*) listBFMixed->FindObject(gHistPNname.Data());
52daf7b2 400 hPNMixed->SetName("gHistPNMixed");
17d10bdd 401 hNPMixed = (AliTHn*) listBFMixed->FindObject(gHistNPname.Data());
52daf7b2 402 hNPMixed->SetName("gHistNPMixed");
17d10bdd 403 hPPMixed = (AliTHn*) listBFMixed->FindObject(gHistPPname.Data());
52daf7b2 404 hPPMixed->SetName("gHistPPMixed");
17d10bdd 405 hNNMixed = (AliTHn*) listBFMixed->FindObject(gHistNNname.Data());
52daf7b2 406 hNNMixed->SetName("gHistNNMixed");
407
408 AliBalancePsi *bMixed = new AliBalancePsi();
32e94079 409 bMixed->SetEventClass(eventClass);
52daf7b2 410 bMixed->SetHistNp(hPMixed);
411 bMixed->SetHistNn(hNMixed);
412 bMixed->SetHistNpn(hPNMixed);
413 bMixed->SetHistNnp(hNPMixed);
414 bMixed->SetHistNpp(hPPMixed);
415 bMixed->SetHistNnn(hNNMixed);
20006629 416 if(kUseVzBinning) bMixed->SetVertexZBinning(kTRUE);
52daf7b2 417 }
418
3ce45c12 419
77b9ff18 420 TH2D *gHistPN[4];
421 TH2D *gHistNP[4];
422 TH2D *gHistPP[4];
423 TH2D *gHistNN[4];
a38fd7f3 424
77b9ff18 425 TCanvas *cPN[4];
426 TCanvas *cNP[4];
427 TCanvas *cPP[4];
428 TCanvas *cNN[4];
a38fd7f3 429 TString histoTitle, pngName;
3ce45c12 430
431 // if no mixing then divide by convoluted histograms
432 if(!listBFMixed && listQA){
3ce45c12 433
e5c5c33f 434 if(!listQA->FindObject("fHistEtaPhiPos") || !listQA->FindObject("fHistEtaPhiNeg")){
3ce45c12 435
e5c5c33f 436 Printf("fHistEtaPhiPos or fHistEtaPhiNeg not found! --> no convolution");
437 listQA = NULL;
3ce45c12 438
e5c5c33f 439 }
440 else{
3ce45c12 441
e5c5c33f 442 TH2D* fHistPos = (TH2D*)((TH3D*)listQA->FindObject("fHistEtaPhiPos"))->Project3D("xy");
443 fHistPos->GetYaxis()->SetRangeUser(-0.79,0.79);
e5c5c33f 444
445 TH2D* fHistNeg = (TH2D*)((TH3D*)listQA->FindObject("fHistEtaPhiNeg"))->Project3D("xy");
446 fHistNeg->GetYaxis()->SetRangeUser(-0.79,0.79);
e5c5c33f 447
448 gHistPN[2] = convolute2D(fHistPos, fHistNeg, "hConvPN");
449 gHistPN[2]->Scale(1./gHistPN[2]->GetBinContent(gHistPN[2]->FindBin(0,0)));
450
451 gHistNP[2] = convolute2D(fHistNeg, fHistPos, "hConvNP");
452 gHistNP[2]->Scale(1./gHistNP[2]->GetBinContent(gHistNP[2]->FindBin(0,0)));
453
454 gHistNN[2] = convolute2D(fHistNeg, fHistNeg, "hConvNN");
455 gHistNN[2]->Scale(1./gHistNN[2]->GetBinContent(gHistNN[2]->FindBin(0,0)));
456
457 gHistPP[2] = convolute2D(fHistPos, fHistPos, "hConvPP");
458 gHistPP[2]->Scale(1./gHistPP[2]->GetBinContent(gHistPP[2]->FindBin(0,0)));
459 }
3ce45c12 460 }
6acdbcb2 461
52daf7b2 462 //(+-)
32e94079 463 if(eventClass == "Centrality"){
464 histoTitle = "(+-) | Centrality: ";
465 histoTitle += psiMin;
466 histoTitle += " - ";
467 histoTitle += psiMax;
468 histoTitle += " % ";
3b8a460b 469 histoTitle += " (0^{o} < #varphi^{t} - #Psi_{2} < 180^{o})";
32e94079 470 }
471 else if(eventClass == "Multiplicity"){
472 histoTitle = "(+-) | Multiplicity: ";
473 histoTitle += psiMin;
474 histoTitle += " - ";
475 histoTitle += psiMax;
476 histoTitle += " tracks";
3b8a460b 477 histoTitle += " (0^{o} < #varphi^{t} - #Psi_{2} < 180^{o})";
32e94079 478 }
479 else{ // "EventPlane" (default)
480 histoTitle = "(+-) | Centrality: ";
481 histoTitle += centralityArray[gCentrality-1];
482 histoTitle += "%";
483 if((psiMin == -0.5)&&(psiMax == 0.5))
3b8a460b 484 histoTitle += " (-7.5^{o} < #varphi^{t} - #Psi_{2} < 7.5^{o})";
32e94079 485 else if((psiMin == 0.5)&&(psiMax == 1.5))
3b8a460b 486 histoTitle += " (37.5^{o} < #varphi^{t} - #Psi_{2} < 52.5^{o})";
32e94079 487 else if((psiMin == 1.5)&&(psiMax == 2.5))
3b8a460b 488 histoTitle += " (82.5^{o} < #varphi^{t} - #Psi_{2} < 97.5^{o})";
32e94079 489 else
3b8a460b 490 histoTitle += " (0^{o} < #varphi^{t} - #Psi_{2} < 180^{o})";
32e94079 491 }
20006629 492 gHistPN[0] = b->GetCorrelationFunctionPN(psiMin,psiMax,vertexZMin,vertexZMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
faa03677 493 if(rebinEta > 1 || rebinPhi > 1){
494 gHistPN[0]->Rebin2D(rebinEta,rebinPhi);
495 gHistPN[0]->Scale(1./(Double_t)(rebinEta*rebinPhi));
496 }
52daf7b2 497 gHistPN[0]->GetYaxis()->SetTitleOffset(1.5);
77b9ff18 498 gHistPN[0]->GetYaxis()->SetTitle("#Delta #varphi (rad)");
52daf7b2 499 gHistPN[0]->SetTitle(histoTitle.Data());
500 cPN[0] = new TCanvas("cPN0","",0,0,600,500);
501 cPN[0]->SetFillColor(10); cPN[0]->SetHighLightColor(10);
5de9ad1a 502 gHistPN[0]->DrawCopy("surf1fb");
503 gPad->SetTheta(30); // default is 30
504 //gPad->SetPhi(130); // default is 30
505 gPad->SetPhi(-60); // default is 30
506 gPad->Update();
6acdbcb2 507 pngName = "DeltaPhiDeltaEta.Centrality";
52daf7b2 508 pngName += centralityArray[gCentrality-1];
6acdbcb2 509 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
510 pngName += ".PositiveNegative.png";
5c9b4733 511 //cPN[0]->SaveAs(pngName.Data());
6acdbcb2 512
52daf7b2 513 if(listBFShuffled) {
34616eda 514
515 histoTitle.ReplaceAll("(+-)","(+-) shuffled");
52daf7b2 516
20006629 517 gHistPN[1] = bShuffled->GetCorrelationFunctionPN(psiMin,psiMax,vertexZMin,vertexZMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
faa03677 518 if(rebinEta > 1 || rebinPhi > 1){
519 gHistPN[1]->Rebin2D(rebinEta,rebinPhi);
520 gHistPN[1]->Scale(1./(Double_t)(rebinEta*rebinPhi));
521 }
52daf7b2 522 gHistPN[1]->GetYaxis()->SetTitleOffset(1.5);
77b9ff18 523 gHistPN[1]->GetYaxis()->SetTitle("#Delta #varphi (rad)");
52daf7b2 524 gHistPN[1]->SetTitle(histoTitle.Data());
525 cPN[1] = new TCanvas("cPN1","",0,100,600,500);
526 cPN[1]->SetFillColor(10);
527 cPN[1]->SetHighLightColor(10);
5de9ad1a 528 gHistPN[1]->DrawCopy("surf1fb");
529 gPad->SetTheta(30); // default is 30
530 //gPad->SetPhi(130); // default is 30
531 gPad->SetPhi(-60); // default is 30
532 gPad->Update();
52daf7b2 533 pngName = "DeltaPhiDeltaEtaShuffled.Centrality";
534 pngName += centralityArray[gCentrality-1];
535 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
536 pngName += ".PositiveNegative.png";
5c9b4733 537 //cPN[1]->SaveAs(pngName.Data());
52daf7b2 538 }
539
540 if(listBFMixed) {
34616eda 541
542 if(histoTitle.Contains("shuffled")) histoTitle.ReplaceAll("(+-) shuffled","(+-) mixed");
543 else histoTitle.ReplaceAll("(+-)","(+-) mixed");
52daf7b2 544
742af4bd 545 // if normalization to trigger then do not divide Event mixing by number of trigger particles
20006629 546 gHistPN[2] = bMixed->GetCorrelationFunctionPN(psiMin,psiMax,vertexZMin,vertexZMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
faa03677 547 if(rebinEta > 1 || rebinPhi > 1){
548 gHistPN[2]->Rebin2D(rebinEta,rebinPhi);
549 gHistPN[2]->Scale(1./(Double_t)(rebinEta*rebinPhi));
550 }
742af4bd 551
552 // normalization to 1 at (0,0) --> Jan Fietes method
553 if(normToTrig){
554 Double_t mixedNorm = gHistPN[2]->Integral(gHistPN[2]->GetXaxis()->FindBin(0-10e-5),gHistPN[2]->GetXaxis()->FindBin(0+10e-5),1,gHistPN[2]->GetNbinsX());
555 mixedNorm /= gHistPN[2]->GetNbinsY()*(gHistPN[2]->GetXaxis()->FindBin(0.01) - gHistPN[2]->GetXaxis()->FindBin(-0.01) + 1);
556 gHistPN[2]->Scale(1./mixedNorm);
557 }
558
52daf7b2 559 gHistPN[2]->GetYaxis()->SetTitleOffset(1.5);
77b9ff18 560 gHistPN[2]->GetYaxis()->SetTitle("#Delta #varphi (rad)");
52daf7b2 561 gHistPN[2]->SetTitle(histoTitle.Data());
562 cPN[2] = new TCanvas("cPN2","",0,200,600,500);
563 cPN[2]->SetFillColor(10);
564 cPN[2]->SetHighLightColor(10);
5de9ad1a 565 gHistPN[2]->DrawCopy("surf1fb");
566 gPad->SetTheta(30); // default is 30
567 //gPad->SetPhi(130); // default is 30
568 gPad->SetPhi(-60); // default is 30
569 gPad->Update();
52daf7b2 570 pngName = "DeltaPhiDeltaEtaMixed.Centrality";
571 pngName += centralityArray[gCentrality-1];
572 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
573 pngName += ".PositiveNegative.png";
5c9b4733 574 //cPN[2]->SaveAs(pngName.Data());
77b9ff18 575
576 //Correlation function (+-)
34616eda 577 gHistPN[3] = b->GetCorrelationFunction("PN",psiMin,psiMax,vertexZMin,vertexZMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax,bMixed);
77b9ff18 578 gHistPN[3]->GetXaxis()->SetRangeUser(-1.5,1.5);
41f2bc59 579 if(normToTrig)
580 gHistPN[3]->GetZaxis()->SetTitle("#frac{1}{N_{trig}}#frac{d^{2}N_{assoc}}{d#Delta#eta#Delta#varphi} (rad^{-1})");
581 else
582 gHistPN[3]->GetZaxis()->SetTitle("C_{+-}(#Delta#eta,#Delta#varphi)");
2b4abbc5 583 //gHistPN[3]->SetTitle("C_{+-}(#Delta#eta,#Delta#varphi)");
77b9ff18 584 cPN[3] = new TCanvas("cPN3","",0,300,600,500);
585 cPN[3]->SetFillColor(10);
586 cPN[3]->SetHighLightColor(10);
587 gHistPN[3]->DrawCopy("surf1fb");
588 gPad->SetTheta(30); // default is 30
589 //gPad->SetPhi(130); // default is 30
590 gPad->SetPhi(-60); // default is 30
591 gPad->Update();
592 pngName = "CorrelationFunction.Centrality";
593 pngName += centralityArray[gCentrality-1];
594 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
595 pngName += ".PositiveNegative.png";
5c9b4733 596 //cPN[3]->SaveAs(pngName.Data());
52daf7b2 597 }
3ce45c12 598 // if no mixing then divide by convoluted histograms
599 else if(listQA){
34616eda 600
601 if(histoTitle.Contains("shuffled")) histoTitle.ReplaceAll("(+-) shuffled","(+-) convoluted");
602 else histoTitle.ReplaceAll("(+-)","(+-) convoluted");
3ce45c12 603
faa03677 604 if(rebinEta > 1 || rebinPhi > 1){
605 gHistPN[2]->Rebin2D(rebinEta,rebinPhi);
606 gHistPN[0]->Scale(1./(Double_t)(rebinEta*rebinPhi));
607 }
608
3ce45c12 609 // normalization to 1 at (0,0) --> Jan Fietes method
610 if(normToTrig){
611 Double_t mixedNorm = gHistPN[2]->Integral(gHistPN[2]->GetXaxis()->FindBin(0-10e-5),gHistPN[2]->GetXaxis()->FindBin(0+10e-5),1,gHistPN[2]->GetNbinsX());
612 mixedNorm /= gHistPN[2]->GetNbinsY()*(gHistPN[2]->GetXaxis()->FindBin(0.01) - gHistPN[2]->GetXaxis()->FindBin(-0.01) + 1);
613 gHistPN[2]->Scale(1./mixedNorm);
614 }
615
616 gHistPN[2]->GetYaxis()->SetTitleOffset(1.5);
617 gHistPN[2]->GetYaxis()->SetTitle("#Delta #varphi (rad)");
618 gHistPN[2]->SetTitle(histoTitle.Data());
619 cPN[2] = new TCanvas("cPN2","",0,200,600,500);
620 cPN[2]->SetFillColor(10);
621 cPN[2]->SetHighLightColor(10);
622 gHistPN[2]->DrawCopy("surf1fb");
623 gPad->SetTheta(30); // default is 30
624 //gPad->SetPhi(130); // default is 30
625 gPad->SetPhi(-60); // default is 30
626 gPad->Update();
627 pngName = "DeltaPhiDeltaEtaMixed.Centrality";
628 pngName += centralityArray[gCentrality-1];
629 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
630 pngName += ".PositiveNegative.png";
631 //cPN[2]->SaveAs(pngName.Data());
632
633 //Correlation function (+-)
634 gHistPN[3] = dynamic_cast<TH2D *>(gHistPN[0]->Clone());
635 gHistPN[3]->Divide(gHistPN[2]);
636 gHistPN[3]->GetXaxis()->SetRangeUser(-1.5,1.5);
41f2bc59 637 if(normToTrig)
638 gHistPN[3]->GetZaxis()->SetTitle("#frac{1}{N_{trig}}#frac{d^{2}N_{assoc}}{d#Delta#eta#Delta#varphi} (rad^{-1})");
639 else
640 gHistPN[3]->GetZaxis()->SetTitle("C_{+-}(#Delta#eta,#Delta#varphi)");
2b4abbc5 641 //gHistPN[3]->SetTitle("C_{+-}(#Delta#eta,#Delta#varphi)");
3ce45c12 642 cPN[3] = new TCanvas("cPN3","",0,300,600,500);
643 cPN[3]->SetFillColor(10);
644 cPN[3]->SetHighLightColor(10);
645 gHistPN[3]->DrawCopy("surf1fb");
646 gPad->SetTheta(30); // default is 30
647 //gPad->SetPhi(130); // default is 30
648 gPad->SetPhi(-60); // default is 30
649 gPad->Update();
650 pngName = "CorrelationFunction.Centrality";
651 pngName += centralityArray[gCentrality-1];
652 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
653 pngName += ".PositiveNegative.png";
654 //cPN[3]->SaveAs(pngName.Data());
655 }
52daf7b2 656
657 //(-+)
32e94079 658 if(eventClass == "Centrality"){
659 histoTitle = "(-+) | Centrality: ";
660 histoTitle += psiMin;
661 histoTitle += " - ";
662 histoTitle += psiMax;
663 histoTitle += " % ";
3b8a460b 664 histoTitle += " (0^{o} < #varphi^{t} - #Psi_{2} < 180^{o})";
32e94079 665 }
666 else if(eventClass == "Multiplicity"){
667 histoTitle = "(-+) | Multiplicity: ";
668 histoTitle += psiMin;
669 histoTitle += " - ";
670 histoTitle += psiMax;
671 histoTitle += " tracks";
3b8a460b 672 histoTitle += " (0^{o} < #varphi^{t} - #Psi_{2} < 180^{o})";
32e94079 673 }
674 else{ // "EventPlane" (default)
675 histoTitle = "(-+) | Centrality: ";
676 histoTitle += centralityArray[gCentrality-1];
677 histoTitle += "%";
678 if((psiMin == -0.5)&&(psiMax == 0.5))
3b8a460b 679 histoTitle += " (-7.5^{o} < #varphi^{t} - #Psi_{2} < 7.5^{o})";
32e94079 680 else if((psiMin == 0.5)&&(psiMax == 1.5))
3b8a460b 681 histoTitle += " (37.5^{o} < #varphi^{t} - #Psi_{2} < 52.5^{o})";
32e94079 682 else if((psiMin == 1.5)&&(psiMax == 2.5))
3b8a460b 683 histoTitle += " (82.5^{o} < #varphi^{t} - #Psi_{2} < 97.5^{o})";
32e94079 684 else
3b8a460b 685 histoTitle += " (0^{o} < #varphi^{t} - #Psi_{2} < 180^{o})";
32e94079 686 }
52daf7b2 687
20006629 688 gHistNP[0] = b->GetCorrelationFunctionNP(psiMin,psiMax,vertexZMin,vertexZMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
faa03677 689 if(rebinEta > 1 || rebinPhi > 1){
690 gHistNP[0]->Rebin2D(rebinEta,rebinPhi);
691 gHistNP[0]->Scale(1./(Double_t)(rebinEta*rebinPhi));
692 }
52daf7b2 693 gHistNP[0]->GetYaxis()->SetTitleOffset(1.5);
77b9ff18 694 gHistNP[0]->GetYaxis()->SetTitle("#Delta #varphi (rad)");
52daf7b2 695 gHistNP[0]->SetTitle(histoTitle.Data());
696 cNP[0] = new TCanvas("cNP0","",100,0,600,500);
697 cNP[0]->SetFillColor(10);
698 cNP[0]->SetHighLightColor(10);
5de9ad1a 699 gHistNP[0]->DrawCopy("surf1fb");
700 gPad->SetTheta(30); // default is 30
701 //gPad->SetPhi(130); // default is 30
702 gPad->SetPhi(-60); // default is 30
703 gPad->Update();
6acdbcb2 704 pngName = "DeltaPhiDeltaEta.Centrality";
52daf7b2 705 pngName += centralityArray[gCentrality-1];
6acdbcb2 706 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
707 pngName += ".NegativePositive.png";
5c9b4733 708 //cNP[0]->SaveAs(pngName.Data());
52daf7b2 709
710 if(listBFShuffled) {
34616eda 711
712 histoTitle.ReplaceAll("(-+)","(-+) shuffled");
52daf7b2 713
20006629 714 gHistNP[1] = bShuffled->GetCorrelationFunctionNP(psiMin,psiMax,vertexZMin,vertexZMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
faa03677 715 if(rebinEta > 1 || rebinPhi > 1){
716 gHistNP[1]->Rebin2D(rebinEta,rebinPhi);
717 gHistNP[1]->Scale(1./(Double_t)(rebinEta*rebinPhi));
718 }
52daf7b2 719 gHistNP[1]->GetYaxis()->SetTitleOffset(1.5);
77b9ff18 720 gHistNP[1]->GetYaxis()->SetTitle("#Delta #varphi (rad)");
52daf7b2 721 gHistNP[1]->SetTitle(histoTitle.Data());
722 cNP[1] = new TCanvas("cNP1","",100,100,600,500);
723 cNP[1]->SetFillColor(10);
724 cNP[1]->SetHighLightColor(10);
5de9ad1a 725 gHistNP[1]->DrawCopy("surf1fb");
726 gPad->SetTheta(30); // default is 30
727 //gPad->SetPhi(130); // default is 30
728 gPad->SetPhi(-60); // default is 30
729 gPad->Update();
52daf7b2 730 pngName = "DeltaPhiDeltaEtaShuffled.Centrality";
731 pngName += centralityArray[gCentrality-1];
732 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
733 pngName += ".NegativePositive.png";
5c9b4733 734 //cNP[1]->SaveAs(pngName.Data());
52daf7b2 735 }
736
737 if(listBFMixed) {
34616eda 738
739 if(histoTitle.Contains("shuffled")) histoTitle.ReplaceAll("(-+) shuffled","(-+) mixed");
740 else histoTitle.ReplaceAll("(-+)","(-+) mixed");
741
742af4bd 742 // if normalization to trigger then do not divide Event mixing by number of trigger particles
20006629 743 gHistNP[2] = bMixed->GetCorrelationFunctionNP(psiMin,psiMax,vertexZMin,vertexZMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
faa03677 744 if(rebinEta > 1 || rebinPhi > 1){
745 gHistNP[2]->Rebin2D(rebinEta,rebinPhi);
746 gHistNP[2]->Scale(1./(Double_t)(rebinEta*rebinPhi));
747 }
742af4bd 748 // normalization to 1 at (0,0) --> Jan Fietes method
749 if(normToTrig){
750 Double_t mixedNorm = gHistNP[2]->Integral(gHistNP[2]->GetXaxis()->FindBin(0-10e-5),gHistNP[2]->GetXaxis()->FindBin(0+10e-5),1,gHistNP[2]->GetNbinsX());
751 mixedNorm /= gHistNP[2]->GetNbinsY()*(gHistNP[2]->GetXaxis()->FindBin(0.01) - gHistNP[2]->GetXaxis()->FindBin(-0.01) + 1);
752 gHistNP[2]->Scale(1./mixedNorm);
753 }
754
52daf7b2 755 gHistNP[2]->GetYaxis()->SetTitleOffset(1.5);
77b9ff18 756 gHistNP[2]->GetYaxis()->SetTitle("#Delta #varphi (rad)");
52daf7b2 757 gHistNP[2]->SetTitle(histoTitle.Data());
758 cNP[2] = new TCanvas("cNP2","",100,200,600,500);
759 cNP[2]->SetFillColor(10);
760 cNP[2]->SetHighLightColor(10);
5de9ad1a 761 gHistNP[2]->DrawCopy("surf1fb");
762 gPad->SetTheta(30); // default is 30
763 //gPad->SetPhi(130); // default is 30
764 gPad->SetPhi(-60); // default is 30
765 gPad->Update();
52daf7b2 766 pngName = "DeltaPhiDeltaEtaMixed.Centrality";
767 pngName += centralityArray[gCentrality-1];
768 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
769 pngName += ".NegativePositive.png";
5c9b4733 770 //cNP[2]->SaveAs(pngName.Data());
77b9ff18 771
772 //Correlation function (-+)
34616eda 773 gHistNP[3] = b->GetCorrelationFunction("NP",psiMin,psiMax,vertexZMin,vertexZMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax,bMixed);
77b9ff18 774 gHistNP[3]->GetXaxis()->SetRangeUser(-1.5,1.5);
41f2bc59 775 if(normToTrig)
776 gHistNP[3]->GetZaxis()->SetTitle("#frac{1}{N_{trig}}#frac{d^{2}N_{assoc}}{d#Delta#eta#Delta#varphi} (rad^{-1})");
777 else
778 gHistNP[3]->GetZaxis()->SetTitle("C_{-+}(#Delta#eta,#Delta#varphi)");
2b4abbc5 779 //gHistNP[3]->SetTitle("C_{+-}(#Delta#eta,#Delta#varphi)");
77b9ff18 780 cNP[3] = new TCanvas("cNP3","",100,300,600,500);
781 cNP[3]->SetFillColor(10);
782 cNP[3]->SetHighLightColor(10);
783 gHistNP[3]->DrawCopy("surf1fb");
784 gPad->SetTheta(30); // default is 30
785 //gPad->SetPhi(130); // default is 30
786 gPad->SetPhi(-60); // default is 30
787 gPad->Update();
788 pngName = "CorrelationFunction.Centrality";
789 pngName += centralityArray[gCentrality-1];
790 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
791 pngName += ".NegativePositive.png";
5c9b4733 792 //cNP[3]->SaveAs(pngName.Data());
52daf7b2 793 }
3ce45c12 794 // if no mixing then divide by convoluted histograms
795 else if(listQA){
34616eda 796
797 if(histoTitle.Contains("shuffled")) histoTitle.ReplaceAll("(-+) shuffled","(-+) convoluted");
798 else histoTitle.ReplaceAll("(-+)","(-+) convoluted");
3ce45c12 799
faa03677 800 if(rebinEta > 1 || rebinPhi > 1){
801 gHistNP[2]->Rebin2D(rebinEta,rebinPhi);
802 gHistNP[2]->Scale(1./(Double_t)(rebinEta*rebinPhi));
803 }
3ce45c12 804
805 // normalization to 1 at (0,0) --> Jan Fietes method
806 if(normToTrig){
807 Double_t mixedNorm = gHistNP[2]->Integral(gHistNP[2]->GetXaxis()->FindBin(0-10e-5),gHistNP[2]->GetXaxis()->FindBin(0+10e-5),1,gHistNP[2]->GetNbinsX());
808 mixedNorm /= gHistNP[2]->GetNbinsY()*(gHistNP[2]->GetXaxis()->FindBin(0.01) - gHistNP[2]->GetXaxis()->FindBin(-0.01) + 1);
809 gHistNP[2]->Scale(1./mixedNorm);
810 }
811
812 gHistNP[2]->GetYaxis()->SetTitleOffset(1.5);
813 gHistNP[2]->GetYaxis()->SetTitle("#Delta #varphi (rad)");
814 gHistNP[2]->SetTitle(histoTitle.Data());
815 cNP[2] = new TCanvas("cNP2","",100,200,600,500);
816 cNP[2]->SetFillColor(10);
817 cNP[2]->SetHighLightColor(10);
818 gHistNP[2]->DrawCopy("surf1fb");
819 gPad->SetTheta(30); // default is 30
820 //gPad->SetPhi(130); // default is 30
821 gPad->SetPhi(-60); // default is 30
822 gPad->Update();
823 pngName = "DeltaPhiDeltaEtaMixed.Centrality";
824 pngName += centralityArray[gCentrality-1];
825 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
826 pngName += ".NegativePositive.png";
827 //cNP[2]->SaveAs(pngName.Data());
828
829 //Correlation function (-+)
830 gHistNP[3] = dynamic_cast<TH2D *>(gHistNP[0]->Clone());
831 gHistNP[3]->Divide(gHistNP[2]);
832 gHistNP[3]->GetXaxis()->SetRangeUser(-1.5,1.5);
41f2bc59 833 if(normToTrig)
834 gHistNP[3]->GetZaxis()->SetTitle("#frac{1}{N_{trig}}#frac{d^{2}N_{assoc}}{d#Delta#eta#Delta#varphi} (rad^{-1})");
835 else
836 gHistNP[3]->GetZaxis()->SetTitle("C_{-+}(#Delta#eta,#Delta#varphi)");
2b4abbc5 837 //gHistNP[3]->SetTitle("C_{-+}(#Delta#eta,#Delta#varphi)");
3ce45c12 838 cNP[3] = new TCanvas("cNP3","",100,300,600,500);
839 cNP[3]->SetFillColor(10);
840 cNP[3]->SetHighLightColor(10);
841 gHistNP[3]->DrawCopy("surf1fb");
842 gPad->SetTheta(30); // default is 30
843 //gPad->SetPhi(130); // default is 30
844 gPad->SetPhi(-60); // default is 30
845 gPad->Update();
846 pngName = "CorrelationFunction.Centrality";
847 pngName += centralityArray[gCentrality-1];
848 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
849 pngName += ".NegativePositive.png";
850 //cNP[3]->SaveAs(pngName.Data());
851 }
852
853
52daf7b2 854 //(++)
32e94079 855 if(eventClass == "Centrality"){
856 histoTitle = "(++) | Centrality: ";
857 histoTitle += psiMin;
858 histoTitle += " - ";
859 histoTitle += psiMax;
860 histoTitle += " % ";
3b8a460b 861 histoTitle += " (0^{o} < #varphi^{t} - #Psi_{2} < 180^{o})";
32e94079 862 }
863 else if(eventClass == "Multiplicity"){
864 histoTitle = "(++) | Multiplicity: ";
865 histoTitle += psiMin;
866 histoTitle += " - ";
867 histoTitle += psiMax;
868 histoTitle += " tracks";
3b8a460b 869 histoTitle += " (0^{o} < #varphi^{t} - #Psi_{2} < 180^{o})";
32e94079 870 }
871 else{ // "EventPlane" (default)
872 histoTitle = "(++) | Centrality: ";
873 histoTitle += centralityArray[gCentrality-1];
874 histoTitle += "%";
875 if((psiMin == -0.5)&&(psiMax == 0.5))
3b8a460b 876 histoTitle += " (-7.5^{o} < #varphi^{t} - #Psi_{2} < 7.5^{o})";
32e94079 877 else if((psiMin == 0.5)&&(psiMax == 1.5))
3b8a460b 878 histoTitle += " (37.5^{o} < #varphi^{t} - #Psi_{2} < 52.5^{o})";
32e94079 879 else if((psiMin == 1.5)&&(psiMax == 2.5))
3b8a460b 880 histoTitle += " (82.5^{o} < #varphi^{t} - #Psi_{2} < 97.5^{o})";
32e94079 881 else
3b8a460b 882 histoTitle += " (0^{o} < #varphi^{t} - #Psi_{2} < 180^{o})";
32e94079 883 }
52daf7b2 884
20006629 885 gHistPP[0] = b->GetCorrelationFunctionPP(psiMin,psiMax,vertexZMin,vertexZMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
faa03677 886 if(rebinEta > 1 || rebinPhi > 1){
887 gHistPP[0]->Rebin2D(rebinEta,rebinPhi);
888 gHistPP[0]->Scale(1./(Double_t)(rebinEta*rebinPhi));
889 }
52daf7b2 890 gHistPP[0]->GetYaxis()->SetTitleOffset(1.5);
77b9ff18 891 gHistPP[0]->GetYaxis()->SetTitle("#Delta #varphi (rad)");
52daf7b2 892 gHistPP[0]->SetTitle(histoTitle.Data());
893 cPP[0] = new TCanvas("cPP0","",200,0,600,500);
894 cPP[0]->SetFillColor(10);
895 cPP[0]->SetHighLightColor(10);
5de9ad1a 896 gHistPP[0]->DrawCopy("surf1fb");
897 gPad->SetTheta(30); // default is 30
898 //gPad->SetPhi(130); // default is 30
899 gPad->SetPhi(-60); // default is 30
900 gPad->Update();
6acdbcb2 901 pngName = "DeltaPhiDeltaEta.Centrality";
52daf7b2 902 pngName += centralityArray[gCentrality-1];
6acdbcb2 903 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
904 pngName += ".PositivePositive.png";
5c9b4733 905 //cPP[0]->SaveAs(pngName.Data());
6acdbcb2 906
52daf7b2 907 if(listBFShuffled) {
34616eda 908
909 histoTitle.ReplaceAll("(++)","(++) shuffled");
52daf7b2 910
20006629 911 gHistPP[1] = bShuffled->GetCorrelationFunctionPP(psiMin,psiMax,vertexZMin,vertexZMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
faa03677 912 if(rebinEta > 1 || rebinPhi > 1){
913 gHistPP[1]->Rebin2D(rebinEta,rebinPhi);
914 gHistPP[1]->Scale(1./(Double_t)(rebinEta*rebinPhi));
915 }
52daf7b2 916 gHistPP[1]->GetYaxis()->SetTitleOffset(1.5);
77b9ff18 917 gHistPP[1]->GetYaxis()->SetTitle("#Delta #varphi (rad)");
52daf7b2 918 gHistPP[1]->SetTitle(histoTitle.Data());
919 cPP[1] = new TCanvas("cPP1","",200,100,600,500);
920 cPP[1]->SetFillColor(10);
921 cPP[1]->SetHighLightColor(10);
5de9ad1a 922 gHistPP[1]->DrawCopy("surf1fb");
923 gPad->SetTheta(30); // default is 30
924 //gPad->SetPhi(130); // default is 30
925 gPad->SetPhi(-60); // default is 30
926 gPad->Update();
52daf7b2 927 pngName = "DeltaPhiDeltaEtaShuffled.Centrality";
928 pngName += centralityArray[gCentrality-1];
929 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
930 pngName += ".PositivePositive.png";
5c9b4733 931 //cPP[1]->SaveAs(pngName.Data());
52daf7b2 932 }
933
934 if(listBFMixed) {
34616eda 935
936 if(histoTitle.Contains("shuffled")) histoTitle.ReplaceAll("(++) shuffled","(++) mixed");
937 else histoTitle.ReplaceAll("(++)","(++) mixed");
52daf7b2 938
742af4bd 939 // if normalization to trigger then do not divide Event mixing by number of trigger particles
20006629 940 gHistPP[2] = bMixed->GetCorrelationFunctionPP(psiMin,psiMax,vertexZMin,vertexZMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
faa03677 941 if(rebinEta > 1 || rebinPhi > 1){
942 gHistPP[2]->Rebin2D(rebinEta,rebinPhi);
943 gHistPP[2]->Scale(1./(Double_t)(rebinEta*rebinPhi));
944 }
742af4bd 945 // normalization to 1 at (0,0) --> Jan Fietes method
946 if(normToTrig){
947 Double_t mixedNorm = gHistPP[2]->Integral(gHistPP[2]->GetXaxis()->FindBin(0-10e-5),gHistPP[2]->GetXaxis()->FindBin(0+10e-5),1,gHistPP[2]->GetNbinsX());
948 mixedNorm /= gHistPP[2]->GetNbinsY()*(gHistPP[2]->GetXaxis()->FindBin(0.01) - gHistPP[2]->GetXaxis()->FindBin(-0.01) + 1);
949 gHistPP[2]->Scale(1./mixedNorm);
950 }
951
52daf7b2 952 gHistPP[2]->GetYaxis()->SetTitleOffset(1.5);
77b9ff18 953 gHistPP[2]->GetYaxis()->SetTitle("#Delta #varphi (rad)");
52daf7b2 954 gHistPP[2]->SetTitle(histoTitle.Data());
955 cPP[2] = new TCanvas("cPP2","",200,200,600,500);
956 cPP[2]->SetFillColor(10);
957 cPP[2]->SetHighLightColor(10);
5de9ad1a 958 gHistPP[2]->DrawCopy("surf1fb");
959 gPad->SetTheta(30); // default is 30
960 //gPad->SetPhi(130); // default is 30
961 gPad->SetPhi(-60); // default is 30
962 gPad->Update();
52daf7b2 963 pngName = "DeltaPhiDeltaEtaMixed.Centrality";
964 pngName += centralityArray[gCentrality-1];
965 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
966 pngName += ".PositivePositive.png";
5c9b4733 967 //cPP[2]->SaveAs(pngName.Data());
77b9ff18 968
969 //Correlation function (++)
34616eda 970 gHistPP[3] = b->GetCorrelationFunction("PP",psiMin,psiMax,vertexZMin,vertexZMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax,bMixed);
77b9ff18 971 gHistPP[3]->GetXaxis()->SetRangeUser(-1.5,1.5);
41f2bc59 972 if(normToTrig)
973 gHistPP[3]->GetZaxis()->SetTitle("#frac{1}{N_{trig}}#frac{d^{2}N_{assoc}}{d#Delta#eta#Delta#varphi} (rad^{-1})");
974 else
975 gHistPP[3]->GetZaxis()->SetTitle("C_{++}(#Delta#eta,#Delta#varphi)");
2b4abbc5 976 // gHistPP[3]->SetTitle("C_{++}(#Delta#eta,#Delta#varphi)");
77b9ff18 977 cPP[3] = new TCanvas("cPP3","",200,300,600,500);
978 cPP[3]->SetFillColor(10);
979 cPP[3]->SetHighLightColor(10);
980 gHistPP[3]->DrawCopy("surf1fb");
981 gPad->SetTheta(30); // default is 30
982 //gPad->SetPhi(130); // default is 30
983 gPad->SetPhi(-60); // default is 30
984 gPad->Update();
985 pngName = "CorrelationFunction.Centrality";
986 pngName += centralityArray[gCentrality-1];
987 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
988 pngName += ".PositivePositive.png";
5c9b4733 989 //cPP[3]->SaveAs(pngName.Data());
52daf7b2 990 }
3ce45c12 991 // if no mixing then divide by convoluted histograms
992 else if(listQA){
34616eda 993
994 if(histoTitle.Contains("shuffled")) histoTitle.ReplaceAll("(++) shuffled","(++) convoluted");
995 else histoTitle.ReplaceAll("(++)","(++) convoluted");
3ce45c12 996
faa03677 997 if(rebinEta > 1 || rebinPhi > 1){
998 gHistPP[2]->Rebin2D(rebinEta,rebinPhi);
999 gHistPP[2]->Scale(1./(Double_t)(rebinEta*rebinPhi));
1000 }
3ce45c12 1001 // normalization to 1 at (0,0) --> Jan Fietes method
1002 if(normToTrig){
1003 Double_t mixedNorm = gHistPP[2]->Integral(gHistPP[2]->GetXaxis()->FindBin(0-10e-5),gHistPP[2]->GetXaxis()->FindBin(0+10e-5),1,gHistPP[2]->GetNbinsX());
1004 mixedNorm /= gHistPP[2]->GetNbinsY()*(gHistPP[2]->GetXaxis()->FindBin(0.01) - gHistPP[2]->GetXaxis()->FindBin(-0.01) + 1);
1005 gHistPP[2]->Scale(1./mixedNorm);
1006 }
1007
1008 gHistPP[2]->GetYaxis()->SetTitleOffset(1.5);
1009 gHistPP[2]->GetYaxis()->SetTitle("#Delta #varphi (rad)");
1010 gHistPP[2]->SetTitle(histoTitle.Data());
1011 cPP[2] = new TCanvas("cPP2","",200,200,600,500);
1012 cPP[2]->SetFillColor(10);
1013 cPP[2]->SetHighLightColor(10);
1014 gHistPP[2]->DrawCopy("surf1fb");
1015 gPad->SetTheta(30); // default is 30
1016 //gPad->SetPhi(130); // default is 30
1017 gPad->SetPhi(-60); // default is 30
1018 gPad->Update();
1019 pngName = "DeltaPhiDeltaEtaMixed.Centrality";
1020 pngName += centralityArray[gCentrality-1];
1021 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
1022 pngName += ".PositivePositive.png";
1023 //cPP[2]->SaveAs(pngName.Data());
1024
1025 //Correlation function (++)
1026 gHistPP[3] = dynamic_cast<TH2D *>(gHistPP[0]->Clone());
1027 gHistPP[3]->Divide(gHistPP[2]);
1028 gHistPP[3]->GetXaxis()->SetRangeUser(-1.5,1.5);
41f2bc59 1029 if(normToTrig)
1030 gHistPP[3]->GetZaxis()->SetTitle("#frac{1}{N_{trig}}#frac{d^{2}N_{assoc}}{d#Delta#eta#Delta#varphi} (rad^{-1})");
1031 else
1032 gHistPP[3]->GetZaxis()->SetTitle("C_{++}(#Delta#eta,#Delta#varphi)");
2b4abbc5 1033 //gHistPP[3]->SetTitle("C_{++}(#Delta#eta,#Delta#varphi)");
3ce45c12 1034 cPP[3] = new TCanvas("cPP3","",200,300,600,500);
1035 cPP[3]->SetFillColor(10);
1036 cPP[3]->SetHighLightColor(10);
1037 gHistPP[3]->DrawCopy("surf1fb");
1038 gPad->SetTheta(30); // default is 30
1039 //gPad->SetPhi(130); // default is 30
1040 gPad->SetPhi(-60); // default is 30
1041 gPad->Update();
1042 pngName = "CorrelationFunction.Centrality";
1043 pngName += centralityArray[gCentrality-1];
1044 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
1045 pngName += ".PositivePositive.png";
1046 //cPP[3]->SaveAs(pngName.Data());
1047 }
52daf7b2 1048
1049 //(--)
32e94079 1050 if(eventClass == "Centrality"){
1051 histoTitle = "(--) | Centrality: ";
1052 histoTitle += psiMin;
1053 histoTitle += " - ";
1054 histoTitle += psiMax;
1055 histoTitle += " % ";
3b8a460b 1056 histoTitle += " (0^{o} < #varphi^{t} - #Psi_{2} < 180^{o})";
32e94079 1057 }
1058 else if(eventClass == "Multiplicity"){
1059 histoTitle = "(--) | Multiplicity: ";
1060 histoTitle += psiMin;
1061 histoTitle += " - ";
1062 histoTitle += psiMax;
1063 histoTitle += " tracks";
3b8a460b 1064 histoTitle += " (0^{o} < #varphi^{t} - #Psi_{2} < 180^{o})";
32e94079 1065 }
1066 else{ // "EventPlane" (default)
1067 histoTitle = "(--) | Centrality: ";
1068 histoTitle += centralityArray[gCentrality-1];
1069 histoTitle += "%";
1070 if((psiMin == -0.5)&&(psiMax == 0.5))
3b8a460b 1071 histoTitle += " (-7.5^{o} < #varphi^{t} - #Psi_{2} < 7.5^{o})";
32e94079 1072 else if((psiMin == 0.5)&&(psiMax == 1.5))
3b8a460b 1073 histoTitle += " (37.5^{o} < #varphi^{t} - #Psi_{2} < 52.5^{o})";
32e94079 1074 else if((psiMin == 1.5)&&(psiMax == 2.5))
3b8a460b 1075 histoTitle += " (82.5^{o} < #varphi^{t} - #Psi_{2} < 97.5^{o})";
32e94079 1076 else
3b8a460b 1077 histoTitle += " (0^{o} < #varphi^{t} - #Psi_{2} < 180^{o})";
32e94079 1078 }
52daf7b2 1079
20006629 1080 gHistNN[0] = b->GetCorrelationFunctionNN(psiMin,psiMax,vertexZMin,vertexZMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
faa03677 1081 if(rebinEta > 1 || rebinPhi > 1){
1082 gHistNN[0]->Rebin2D(rebinEta,rebinPhi);
1083 gHistNN[0]->Scale(1./(Double_t)(rebinEta*rebinPhi));
1084 }
52daf7b2 1085 gHistNN[0]->GetYaxis()->SetTitleOffset(1.5);
77b9ff18 1086 gHistNN[0]->GetYaxis()->SetTitle("#Delta #varphi (rad)");
52daf7b2 1087 gHistNN[0]->SetTitle(histoTitle.Data());
1088 cNN[0] = new TCanvas("cNN0","",300,0,600,500);
1089 cNN[0]->SetFillColor(10);
1090 cNN[0]->SetHighLightColor(10);
5de9ad1a 1091 gHistNN[0]->DrawCopy("surf1fb");
1092 gPad->SetTheta(30); // default is 30
1093 gPad->SetPhi(-60); // default is 30
1094 //gPad->SetPhi(-60); // default is 30
1095 gPad->Update();
6acdbcb2 1096 pngName = "DeltaPhiDeltaEta.Centrality";
52daf7b2 1097 pngName += centralityArray[gCentrality-1];
6acdbcb2 1098 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
1099 pngName += ".NegativeNegative.png";
5c9b4733 1100 //cNN[0]->SaveAs(pngName.Data());
52daf7b2 1101
1102 if(listBFShuffled) {
34616eda 1103
1104 histoTitle.ReplaceAll("(--)","(--) shuffled");
52daf7b2 1105
20006629 1106 gHistNN[1] = bShuffled->GetCorrelationFunctionNN(psiMin,psiMax,vertexZMin,vertexZMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
faa03677 1107 if(rebinEta > 1 || rebinPhi > 1){
1108 gHistNN[1]->Rebin2D(rebinEta,rebinPhi);
1109 gHistNN[1]->Scale(1./(Double_t)(rebinEta*rebinPhi));
1110 }
52daf7b2 1111 gHistNN[1]->GetYaxis()->SetTitleOffset(1.5);
77b9ff18 1112 gHistNN[1]->GetYaxis()->SetTitle("#Delta #varphi (rad)");
52daf7b2 1113 gHistNN[1]->SetTitle(histoTitle.Data());
1114 cNN[1] = new TCanvas("cNN1","",300,100,600,500);
1115 cNN[1]->SetFillColor(10);
1116 cNN[1]->SetHighLightColor(10);
5de9ad1a 1117 gHistNN[1]->DrawCopy("surf1fb");
1118 gPad->SetTheta(30); // default is 30
1119 //gPad->SetPhi(130); // default is 30
1120 gPad->SetPhi(-60); // default is 30
1121 gPad->Update();
52daf7b2 1122 pngName = "DeltaPhiDeltaEtaShuffled.Centrality";
1123 pngName += centralityArray[gCentrality-1];
1124 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
1125 pngName += ".NegativeNegative.png";
5c9b4733 1126 //cNN[1]->SaveAs(pngName.Data());
52daf7b2 1127 }
1128
1129 if(listBFMixed) {
34616eda 1130
1131 if(histoTitle.Contains("shuffled")) histoTitle.ReplaceAll("(--) shuffled","(--) mixed");
1132 else histoTitle.ReplaceAll("(--)","(--) mixed");
52daf7b2 1133
742af4bd 1134 // if normalization to trigger then do not divide Event mixing by number of trigger particles
20006629 1135 gHistNN[2] = bMixed->GetCorrelationFunctionNN(psiMin,psiMax,vertexZMin,vertexZMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
faa03677 1136 if(rebinEta > 1 || rebinPhi > 1){
1137 gHistNN[2]->Rebin2D(rebinEta,rebinPhi);
1138 gHistNN[2]->Scale(1./(Double_t)(rebinEta*rebinPhi));
1139 }
742af4bd 1140 // normalization to 1 at (0,0) --> Jan Fietes method
1141 if(normToTrig){
1142 Double_t mixedNorm = gHistNN[2]->Integral(gHistNN[2]->GetXaxis()->FindBin(0-10e-5),gHistNN[2]->GetXaxis()->FindBin(0+10e-5),1,gHistNN[2]->GetNbinsX());
1143 mixedNorm /= gHistNN[2]->GetNbinsY()*(gHistNN[2]->GetXaxis()->FindBin(0.01) - gHistNN[2]->GetXaxis()->FindBin(-0.01) + 1);
1144 gHistNN[2]->Scale(1./mixedNorm);
1145 }
1146
52daf7b2 1147 gHistNN[2]->GetYaxis()->SetTitleOffset(1.5);
77b9ff18 1148 gHistNN[2]->GetYaxis()->SetTitle("#Delta #varphi (rad)");
52daf7b2 1149 gHistNN[2]->SetTitle(histoTitle.Data());
1150 cNN[2] = new TCanvas("cNN2","",300,200,600,500);
1151 cNN[2]->SetFillColor(10);
1152 cNN[2]->SetHighLightColor(10);
5de9ad1a 1153 gHistNN[2]->DrawCopy("surf1fb");
1154 gPad->SetTheta(30); // default is 30
1155 //gPad->SetPhi(130); // default is 30
1156 gPad->SetPhi(-60); // default is 30
1157 gPad->Update();
52daf7b2 1158 pngName = "DeltaPhiDeltaEtaMixed.Centrality";
1159 pngName += centralityArray[gCentrality-1];
1160 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
1161 pngName += ".NegativeNegative.png";
5c9b4733 1162 //cNN[2]->SaveAs(pngName.Data());
77b9ff18 1163
1164 //Correlation function (--)
34616eda 1165 gHistNN[3] = b->GetCorrelationFunction("NN",psiMin,psiMax,vertexZMin,vertexZMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax,bMixed);
77b9ff18 1166 gHistNN[3]->GetXaxis()->SetRangeUser(-1.5,1.5);
41f2bc59 1167 if(normToTrig)
1168 gHistNN[3]->GetZaxis()->SetTitle("#frac{1}{N_{trig}}#frac{d^{2}N_{assoc}}{d#Delta#eta#Delta#varphi} (rad^{-1})");
1169 else
1170 gHistNN[3]->GetZaxis()->SetTitle("C_{--}(#Delta#eta,#Delta#varphi)");
2b4abbc5 1171 // gHistNN[3]->SetTitle("C_{--}(#Delta#eta,#Delta#varphi)");
77b9ff18 1172 cNN[3] = new TCanvas("cNN3","",300,300,600,500);
1173 cNN[3]->SetFillColor(10);
1174 cNN[3]->SetHighLightColor(10);
1175 gHistNN[3]->DrawCopy("surf1fb");
1176 gPad->SetTheta(30); // default is 30
1177 //gPad->SetPhi(130); // default is 30
1178 gPad->SetPhi(-60); // default is 30
1179 gPad->Update();
1180 pngName = "CorrelationFunction.Centrality";
1181 pngName += centralityArray[gCentrality-1];
1182 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
1183 pngName += ".NegativeNegative.png";
5c9b4733 1184 //cNN[3]->SaveAs(pngName.Data());
52daf7b2 1185 }
3ce45c12 1186 // if no mixing then divide by convoluted histograms
1187 else if(listQA){
34616eda 1188
1189 if(histoTitle.Contains("shuffled")) histoTitle.ReplaceAll("(--) shuffled","(--) convoluted");
1190 else histoTitle.ReplaceAll("(--)","(--) convoluted");
3ce45c12 1191
faa03677 1192 if(rebinEta > 1 || rebinPhi > 1){
1193 gHistNN[2]->Rebin2D(rebinEta,rebinPhi);
1194 gHistNN[2]->Scale(1./(Double_t)(rebinEta*rebinPhi));
1195 }
3ce45c12 1196 // normalization to 1 at (0,0) --> Jan Fietes method
1197 if(normToTrig){
1198 Double_t mixedNorm = gHistNN[2]->Integral(gHistNN[2]->GetXaxis()->FindBin(0-10e-5),gHistNN[2]->GetXaxis()->FindBin(0+10e-5),1,gHistNN[2]->GetNbinsX());
1199 mixedNorm /= gHistNN[2]->GetNbinsY()*(gHistNN[2]->GetXaxis()->FindBin(0.01) - gHistNN[2]->GetXaxis()->FindBin(-0.01) + 1);
1200 gHistNN[2]->Scale(1./mixedNorm);
1201 }
1202
1203 gHistNN[2]->GetYaxis()->SetTitleOffset(1.5);
1204 gHistNN[2]->GetYaxis()->SetTitle("#Delta #varphi (rad)");
1205 gHistNN[2]->SetTitle(histoTitle.Data());
1206 cNN[2] = new TCanvas("cNN2","",300,200,600,500);
1207 cNN[2]->SetFillColor(10);
1208 cNN[2]->SetHighLightColor(10);
1209 gHistNN[2]->DrawCopy("surf1fb");
1210 gPad->SetTheta(30); // default is 30
1211 //gPad->SetPhi(130); // default is 30
1212 gPad->SetPhi(-60); // default is 30
1213 gPad->Update();
1214 pngName = "DeltaPhiDeltaEtaMixed.Centrality";
1215 pngName += centralityArray[gCentrality-1];
1216 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
1217 pngName += ".NegativeNegative.png";
1218 //cNN[2]->SaveAs(pngName.Data());
1219
1220 //Correlation function (--)
1221 gHistNN[3] = dynamic_cast<TH2D *>(gHistNN[0]->Clone());
1222 gHistNN[3]->Divide(gHistNN[2]);
1223 gHistNN[3]->GetXaxis()->SetRangeUser(-1.5,1.5);
41f2bc59 1224 if(normToTrig)
1225 gHistNN[3]->GetZaxis()->SetTitle("#frac{1}{N_{trig}}#frac{d^{2}N_{assoc}}{d#Delta#eta#Delta#varphi} (rad^{-1})");
1226 else
1227 gHistNN[3]->GetZaxis()->SetTitle("C_{--}(#Delta#eta,#Delta#varphi)");
2b4abbc5 1228 //gHistNN[3]->SetTitle("C_{--}(#Delta#eta,#Delta#varphi)");
3ce45c12 1229 cNN[3] = new TCanvas("cNN3","",300,300,600,500);
1230 cNN[3]->SetFillColor(10);
1231 cNN[3]->SetHighLightColor(10);
1232 gHistNN[3]->DrawCopy("surf1fb");
1233 gPad->SetTheta(30); // default is 30
1234 //gPad->SetPhi(130); // default is 30
1235 gPad->SetPhi(-60); // default is 30
1236 gPad->Update();
1237 pngName = "CorrelationFunction.Centrality";
1238 pngName += centralityArray[gCentrality-1];
1239 pngName += ".Psi"; pngName += psiMin; pngName += "To"; pngName += psiMax;
1240 pngName += ".NegativeNegative.png";
1241 //cNN[3]->SaveAs(pngName.Data());
1242 }
5c9b4733 1243
1244 //Write to output file
32e94079 1245 TString newFileName = "correlationFunction.";
1246 if(eventClass == "Centrality"){
1247 newFileName += Form("Centrality%.1fTo%.1f",psiMin,psiMax);
1248 newFileName += ".PsiAll.PttFrom";
1249 }
1250 else if(eventClass == "Multiplicity"){
1251 newFileName += Form("Multiplicity%.0fTo%.0f",psiMin,psiMax);
1252 newFileName += ".PsiAll.PttFrom";
1253 }
1254 else{ // "EventPlane" (default)
1255 newFileName += "Centrality";
1256 newFileName += gCentrality; newFileName += ".Psi";
1257 if((psiMin == -0.5)&&(psiMax == 0.5)) newFileName += "InPlane.Ptt";
1258 else if((psiMin == 0.5)&&(psiMax == 1.5)) newFileName += "Intermediate.Ptt";
1259 else if((psiMin == 1.5)&&(psiMax == 2.5)) newFileName += "OutOfPlane.Ptt";
1260 else if((psiMin == 2.5)&&(psiMax == 3.5)) newFileName += "Rest.PttFrom";
1261 else newFileName += "All.PttFrom";
1262 }
5c9b4733 1263 newFileName += Form("%.1f",ptTriggerMin); newFileName += "To";
1264 newFileName += Form("%.1f",ptTriggerMax); newFileName += "PtaFrom";
1265 newFileName += Form("%.1f",ptAssociatedMin); newFileName += "To";
1266 newFileName += Form("%.1f",ptAssociatedMax);
09b14834 1267
1268 newFileName += "_";
1269 newFileName += Form("%.1f",psiMin);
1270 newFileName += "-";
1271 newFileName += Form("%.1f",psiMax);
5c9b4733 1272 newFileName += ".root";
09b14834 1273
5c9b4733 1274 TFile *newFile = TFile::Open(newFileName.Data(),"recreate");
1275 gHistPN[0]->SetName("gHistPNRaw"); gHistPN[0]->Write();
1276 gHistNP[0]->SetName("gHistNPRaw"); gHistNP[0]->Write();
1277 gHistPP[0]->SetName("gHistPPRaw"); gHistPP[0]->Write();
1278 gHistNN[0]->SetName("gHistNNRaw"); gHistNN[0]->Write();
1279 if(listBFShuffled) {
1280 gHistPN[1]->SetName("gHistPNShuffled"); gHistPN[1]->Write();
1281 gHistNP[1]->SetName("gHistNPShuffled"); gHistNP[1]->Write();
1282 gHistPP[1]->SetName("gHistPPShuffled"); gHistPP[1]->Write();
1283 gHistNN[1]->SetName("gHistNNShuffled"); gHistNN[1]->Write();
1284 }
3ce45c12 1285 if(listBFMixed || (!listBFMixed&&listQA)) {
5c9b4733 1286 gHistPN[2]->SetName("gHistPNMixed"); gHistPN[2]->Write();
1287 gHistNP[2]->SetName("gHistNPMixed"); gHistNP[2]->Write();
1288 gHistPP[2]->SetName("gHistPPMixed"); gHistPP[2]->Write();
1289 gHistNN[2]->SetName("gHistNNMixed"); gHistNN[2]->Write();
1290
1291 gHistPN[3]->SetName("gHistPNCorrelationFunctions"); gHistPN[3]->Write();
1292 gHistNP[3]->SetName("gHistNPCorrelationFunctions"); gHistNP[3]->Write();
1293 gHistPP[3]->SetName("gHistPPCorrelationFunctions"); gHistPP[3]->Write();
1294 gHistNN[3]->SetName("gHistNNCorrelationFunctions"); gHistNN[3]->Write();
1295 }
1296 newFile->Close();
07d0a35c 1297
1298 // some cleaning
1299 for(Int_t i = 0; i < 4; i++){
1300
1301 if(!listBFShuffled && i == 1) continue;
1302 if(!listBFMixed && (i == 2 || i == 3)) continue;
1303
1304 if(gHistPP[i]) delete gHistPP[i];
1305 if(gHistPN[i]) delete gHistPN[i];
1306 if(gHistNP[i]) delete gHistNP[i];
1307 if(gHistNN[i]) delete gHistNN[i];
1308
1309 if(cPN[i]) delete cPN[i];
1310 if(cNP[i]) delete cNP[i];
1311 if(cPP[i]) delete cPP[i];
1312 if(cNN[i]) delete cNN[i];
1313 }
1314
1315 delete hP;
1316 delete hN;
1317 delete hPP;
1318 delete hPN;
1319 delete hNP;
1320 delete hNN;
1321
1322 delete hPMixed;
1323 delete hNMixed;
1324 delete hPPMixed;
1325 delete hPNMixed;
1326 delete hNPMixed;
1327 delete hNNMixed;
1328
1329 delete hPShuffled;
1330 delete hNShuffled;
1331 delete hPPShuffled;
1332 delete hPNShuffled;
1333 delete hNPShuffled;
1334 delete hNNShuffled;
1335
5c9b4733 1336}
1337
1338//____________________________________________________________//
3b8a460b 1339void drawCorrelationFunctions(const char* lhcPeriod = "LHC10h",
1340 const char* gCentralityEstimator = "V0M",
1341 Int_t gBit = 128,
1342 const char* gEventPlaneEstimator = "VZERO",
5c9b4733 1343 Int_t gCentrality = 1,
1344 Double_t psiMin = -0.5, Double_t psiMax = 3.5,
20006629 1345 Double_t vertexZMin = -10.,
1346 Double_t vertexZMax = 10.,
5c9b4733 1347 Double_t ptTriggerMin = -1.,
1348 Double_t ptTriggerMax = -1.,
1349 Double_t ptAssociatedMin = -1.,
3b8a460b 1350 Double_t ptAssociatedMax = -1.,
1351 Bool_t kFit = kFALSE) {
5c9b4733 1352 //Macro that draws the charge dependent correlation functions
1353 //for each centrality bin for the different pT of trigger and
1354 //associated particles
1355 //Author: Panos.Christakoglou@nikhef.nl
1356 TGaxis::SetMaxDigits(3);
1357
1358 //Get the input file
09b14834 1359 /* TString filename = lhcPeriod;
3b8a460b 1360 filename += "/Centrality"; filename += gCentralityEstimator;
1361 filename += "_Bit"; filename += gBit;
1362 filename += "_"; filename += gEventPlaneEstimator;
1363 filename +="/PttFrom";
1364 filename += Form("%.1f",ptTriggerMin); filename += "To";
1365 filename += Form("%.1f",ptTriggerMax); filename += "PtaFrom";
1366 filename += Form("%.1f",ptAssociatedMin); filename += "To";
09b14834 1367 filename += Form("%.1f",ptAssociatedMax); */
1368
1369 TString filename = "correlationFunction.Centrality";
5c9b4733 1370 filename += gCentrality; filename += ".Psi";
1371 if((psiMin == -0.5)&&(psiMax == 0.5)) filename += "InPlane.Ptt";
1372 else if((psiMin == 0.5)&&(psiMax == 1.5)) filename += "Intermediate.Ptt";
1373 else if((psiMin == 1.5)&&(psiMax == 2.5)) filename += "OutOfPlane.Ptt";
1374 else if((psiMin == 2.5)&&(psiMax == 3.5)) filename += "Rest.Ptt";
1375 else filename += "All.PttFrom";
1376 filename += Form("%.1f",ptTriggerMin); filename += "To";
1377 filename += Form("%.1f",ptTriggerMax); filename += "PtaFrom";
1378 filename += Form("%.1f",ptAssociatedMin); filename += "To";
1379 filename += Form("%.1f",ptAssociatedMax);
09b14834 1380
1381 filename += "_";
1382 filename += Form("%.1f",psiMin);
1383 filename += "-";
1384 filename += Form("%.1f",psiMax);
5c9b4733 1385 filename += ".root";
1386
1387 //Open the file
1388 TFile *f = TFile::Open(filename.Data());
1389 if((!f)||(!f->IsOpen())) {
1390 Printf("The file %s is not found. Aborting...",filename);
3b8a460b 1391 return;
5c9b4733 1392 }
1393 //f->ls();
1394
1395 //Latex
1396 TString centralityLatex = "Centrality: ";
1397 centralityLatex += centralityArray[gCentrality-1];
1398 centralityLatex += "%";
1399
1400 TString psiLatex;
1401 if((psiMin == -0.5)&&(psiMax == 0.5))
3b8a460b 1402 psiLatex = " -7.5^{o} < #varphi^{t} - #Psi_{2} < 7.5^{o}";
5c9b4733 1403 else if((psiMin == 0.5)&&(psiMax == 1.5))
3b8a460b 1404 psiLatex = " 37.5^{o} < #varphi^{t} - #Psi_{2} < 52.5^{o}";
5c9b4733 1405 else if((psiMin == 1.5)&&(psiMax == 2.5))
3b8a460b 1406 psiLatex = " 82.5^{o} < #varphi^{t} - #Psi_{2} < 97.5^{o}";
5c9b4733 1407 else
3b8a460b 1408 psiLatex = " 0^{o} < #varphi^{t} - #Psi_{2} < 180^{o}";
5c9b4733 1409
1410 TString pttLatex = Form("%.1f",ptTriggerMin);
1411 pttLatex += " < p_{T,trig} < "; pttLatex += Form("%.1f",ptTriggerMax);
1412 pttLatex += " GeV/c";
1413
1414 TString ptaLatex = Form("%.1f",ptAssociatedMin);
1415 ptaLatex += " < p_{T,assoc} < "; ptaLatex += Form("%.1f",ptAssociatedMax);
1416 ptaLatex += " GeV/c";
1417
1418 TLatex *latexInfo1 = new TLatex();
1419 latexInfo1->SetNDC();
1420 latexInfo1->SetTextSize(0.045);
1421 latexInfo1->SetTextColor(1);
1422
1423 TString pngName;
1424
1425 //============================================================//
1426 //Get the +- correlation function
1427 TH2D *gHistPN = dynamic_cast<TH2D *>(f->Get("gHistPNCorrelationFunctions"));
1428 gHistPN->SetStats(kFALSE);
2b4abbc5 1429 gHistPN->SetTitle("C_{+-}(#Delta#eta,#Delta#varphi)");
3b8a460b 1430 gHistPN->GetXaxis()->SetRangeUser(-1.4,1.4);
5c9b4733 1431 gHistPN->GetXaxis()->CenterTitle();
1432 gHistPN->GetXaxis()->SetTitleOffset(1.2);
1433 gHistPN->GetYaxis()->CenterTitle();
1434 gHistPN->GetYaxis()->SetTitleOffset(1.2);
2b4abbc5 1435 gHistPN->GetZaxis()->SetTitleOffset(1.5);
5c9b4733 1436 TCanvas *cPN = new TCanvas("cPN","",0,0,600,500);
1437 cPN->SetFillColor(10); cPN->SetHighLightColor(10);
1438 cPN->SetLeftMargin(0.15);
1439 gHistPN->DrawCopy("surf1fb");
2b4abbc5 1440
5c9b4733 1441 gPad->SetTheta(30); // default is 30
1442 gPad->SetPhi(-60); // default is 30
1443 gPad->Update();
1444
3b8a460b 1445 latexInfo1->DrawLatex(0.6,0.95,centralityLatex.Data());
1446 latexInfo1->DrawLatex(0.6,0.89,psiLatex.Data());
1447 latexInfo1->DrawLatex(0.6,0.83,pttLatex.Data());
1448 latexInfo1->DrawLatex(0.6,0.77,ptaLatex.Data());
5c9b4733 1449
1450 pngName = "CorrelationFunction.Centrality";
1451 pngName += centralityArray[gCentrality-1];
1452 pngName += ".Psi";
1453 if((psiMin == -0.5)&&(psiMax == 0.5)) pngName += "InPlane.Ptt";
1454 else if((psiMin == 0.5)&&(psiMax == 1.5)) pngName += "Intermediate.Ptt";
1455 else if((psiMin == 1.5)&&(psiMax == 2.5)) pngName += "OutOfPlane.Ptt";
1456 else if((psiMin == 2.5)&&(psiMax == 3.5)) pngName += "Rest.Ptt";
1457 else pngName += "All.PttFrom";
1458 pngName += Form("%.1f",ptTriggerMin); pngName += "To";
1459 pngName += Form("%.1f",ptTriggerMax); pngName += "PtaFrom";
1460 pngName += Form("%.1f",ptAssociatedMin); pngName += "To";
1461 pngName += Form("%.1f",ptAssociatedMax);
1462 pngName += ".PositiveNegative.png";
1463 cPN->SaveAs(pngName.Data());
3b8a460b 1464 if(kFit)
20006629 1465 fitCorrelationFunctions(gCentrality, psiMin, psiMax,vertexZMin, vertexZMax,
3b8a460b 1466 ptTriggerMin,ptTriggerMax,
1467 ptAssociatedMin, ptAssociatedMax,gHistPN);
1468
5c9b4733 1469 //============================================================//
1470 //Get the -+ correlation function
1471 TH2D *gHistNP = dynamic_cast<TH2D *>(f->Get("gHistNPCorrelationFunctions"));
1472 gHistNP->SetStats(kFALSE);
2b4abbc5 1473 gHistNP->SetTitle("C_{-+}(#Delta#eta,#Delta#varphi)");
3b8a460b 1474 gHistNP->GetXaxis()->SetRangeUser(-1.4,1);
5c9b4733 1475 gHistNP->GetXaxis()->CenterTitle();
1476 gHistNP->GetXaxis()->SetTitleOffset(1.2);
1477 gHistNP->GetYaxis()->CenterTitle();
1478 gHistNP->GetYaxis()->SetTitleOffset(1.2);
2b4abbc5 1479 gHistNP->GetZaxis()->SetTitleOffset(1.5);
5c9b4733 1480 TCanvas *cNP = new TCanvas("cNP","",50,50,600,500);
1481 cNP->SetFillColor(10); cNP->SetHighLightColor(10);
1482 cNP->SetLeftMargin(0.15);
1483 gHistNP->DrawCopy("surf1fb");
2b4abbc5 1484
5c9b4733 1485 gPad->SetTheta(30); // default is 30
1486 gPad->SetPhi(-60); // default is 30
1487 gPad->Update();
1488
3b8a460b 1489 latexInfo1->DrawLatex(0.6,0.95,centralityLatex.Data());
1490 latexInfo1->DrawLatex(0.6,0.89,psiLatex.Data());
1491 latexInfo1->DrawLatex(0.6,0.83,pttLatex.Data());
1492 latexInfo1->DrawLatex(0.6,0.77,ptaLatex.Data());
5c9b4733 1493
1494 pngName = "CorrelationFunction.Centrality";
1495 pngName += centralityArray[gCentrality-1];
1496 pngName += ".Psi";
1497 if((psiMin == -0.5)&&(psiMax == 0.5)) pngName += "InPlane.Ptt";
1498 else if((psiMin == 0.5)&&(psiMax == 1.5)) pngName += "Intermediate.Ptt";
1499 else if((psiMin == 1.5)&&(psiMax == 2.5)) pngName += "OutOfPlane.Ptt";
1500 else if((psiMin == 2.5)&&(psiMax == 3.5)) pngName += "Rest.Ptt";
1501 else pngName += "All.PttFrom";
1502 pngName += Form("%.1f",ptTriggerMin); pngName += "To";
1503 pngName += Form("%.1f",ptTriggerMax); pngName += "PtaFrom";
1504 pngName += Form("%.1f",ptAssociatedMin); pngName += "To";
1505 pngName += Form("%.1f",ptAssociatedMax);
1506 pngName += ".NegativePositive.png";
1507 cNP->SaveAs(pngName.Data());
1508
3b8a460b 1509 if(kFit)
20006629 1510 fitCorrelationFunctions(gCentrality, psiMin, psiMax,vertexZMin, vertexZMax,
3b8a460b 1511 ptTriggerMin,ptTriggerMax,
1512 ptAssociatedMin, ptAssociatedMax,gHistNP);
1513
5c9b4733 1514 //============================================================//
1515 //Get the ++ correlation function
1516 TH2D *gHistPP = dynamic_cast<TH2D *>(f->Get("gHistPPCorrelationFunctions"));
1517 gHistPP->SetStats(kFALSE);
2b4abbc5 1518 gHistPP->SetTitle("C_{++}(#Delta#eta,#Delta#varphi)");
3b8a460b 1519 gHistPP->GetXaxis()->SetRangeUser(-1.4,1.4);
5c9b4733 1520 gHistPP->GetXaxis()->CenterTitle();
1521 gHistPP->GetXaxis()->SetTitleOffset(1.2);
1522 gHistPP->GetYaxis()->CenterTitle();
1523 gHistPP->GetYaxis()->SetTitleOffset(1.2);
2b4abbc5 1524 gHistPP->GetZaxis()->SetTitleOffset(1.5);
5c9b4733 1525 TCanvas *cPP = new TCanvas("cPP","",100,100,600,500);
1526 cPP->SetFillColor(10); cPP->SetHighLightColor(10);
1527 cPP->SetLeftMargin(0.15);
1528 gHistPP->DrawCopy("surf1fb");
2b4abbc5 1529
5c9b4733 1530 gPad->SetTheta(30); // default is 30
1531 gPad->SetPhi(-60); // default is 30
1532 gPad->Update();
1533
3b8a460b 1534 latexInfo1->DrawLatex(0.6,0.95,centralityLatex.Data());
1535 latexInfo1->DrawLatex(0.6,0.89,psiLatex.Data());
1536 latexInfo1->DrawLatex(0.6,0.83,pttLatex.Data());
1537 latexInfo1->DrawLatex(0.6,0.77,ptaLatex.Data());
5c9b4733 1538
1539 pngName = "CorrelationFunction.Centrality";
1540 pngName += centralityArray[gCentrality-1];
1541 pngName += ".Psi";
1542 if((psiMin == -0.5)&&(psiMax == 0.5)) pngName += "InPlane.Ptt";
1543 else if((psiMin == 0.5)&&(psiMax == 1.5)) pngName += "Intermediate.Ptt";
1544 else if((psiMin == 1.5)&&(psiMax == 2.5)) pngName += "OutOfPlane.Ptt";
1545 else if((psiMin == 2.5)&&(psiMax == 3.5)) pngName += "Rest.Ptt";
1546 else pngName += "All.PttFrom";
1547 pngName += Form("%.1f",ptTriggerMin); pngName += "To";
1548 pngName += Form("%.1f",ptTriggerMax); pngName += "PtaFrom";
1549 pngName += Form("%.1f",ptAssociatedMin); pngName += "To";
1550 pngName += Form("%.1f",ptAssociatedMax);
1551 pngName += ".PositivePositive.png";
1552 cPP->SaveAs(pngName.Data());
1553
3b8a460b 1554 if(kFit)
20006629 1555 fitCorrelationFunctions(gCentrality, psiMin, psiMax,vertexZMin, vertexZMax,
3b8a460b 1556 ptTriggerMin,ptTriggerMax,
1557 ptAssociatedMin, ptAssociatedMax,gHistPP);
1558
5c9b4733 1559 //============================================================//
1560 //Get the -- correlation function
1561 TH2D *gHistNN = dynamic_cast<TH2D *>(f->Get("gHistNNCorrelationFunctions"));
1562 gHistNN->SetStats(kFALSE);
2b4abbc5 1563 gHistNN->SetTitle("C_{--}(#Delta#eta,#Delta#varphi)");
3b8a460b 1564 gHistNN->GetXaxis()->SetRangeUser(-1.4,1.4);
5c9b4733 1565 gHistNN->GetXaxis()->CenterTitle();
1566 gHistNN->GetXaxis()->SetTitleOffset(1.2);
1567 gHistNN->GetYaxis()->CenterTitle();
1568 gHistNN->GetYaxis()->SetTitleOffset(1.2);
2b4abbc5 1569 gHistNN->GetZaxis()->SetTitleOffset(1.5);
5c9b4733 1570 TCanvas *cNN = new TCanvas("cNN","",150,150,600,500);
1571 cNN->SetFillColor(10); cNN->SetHighLightColor(10);
1572 cNN->SetLeftMargin(0.15);
1573 gHistNN->DrawCopy("surf1fb");
2b4abbc5 1574
5c9b4733 1575 gPad->SetTheta(30); // default is 30
1576 gPad->SetPhi(-60); // default is 30
1577 gPad->Update();
1578
3b8a460b 1579 latexInfo1->DrawLatex(0.6,0.95,centralityLatex.Data());
1580 latexInfo1->DrawLatex(0.6,0.89,psiLatex.Data());
1581 latexInfo1->DrawLatex(0.6,0.83,pttLatex.Data());
1582 latexInfo1->DrawLatex(0.6,0.77,ptaLatex.Data());
5c9b4733 1583
1584 pngName = "CorrelationFunction.Centrality";
1585 pngName += centralityArray[gCentrality-1];
1586 pngName += ".Psi";
1587 if((psiMin == -0.5)&&(psiMax == 0.5)) pngName += "InPlane.Ptt";
1588 else if((psiMin == 0.5)&&(psiMax == 1.5)) pngName += "Intermediate.Ptt";
1589 else if((psiMin == 1.5)&&(psiMax == 2.5)) pngName += "OutOfPlane.Ptt";
1590 else if((psiMin == 2.5)&&(psiMax == 3.5)) pngName += "Rest.Ptt";
1591 else pngName += "All.PttFrom";
1592 pngName += Form("%.1f",ptTriggerMin); pngName += "To";
1593 pngName += Form("%.1f",ptTriggerMax); pngName += "PtaFrom";
1594 pngName += Form("%.1f",ptAssociatedMin); pngName += "To";
1595 pngName += Form("%.1f",ptAssociatedMax);
1596 pngName += ".NegativeNegative.png";
1597 cNN->SaveAs(pngName.Data());
24ca71a6 1598
3b8a460b 1599 if(kFit)
20006629 1600 fitCorrelationFunctions(gCentrality, psiMin, psiMax,vertexZMin, vertexZMax,
3b8a460b 1601 ptTriggerMin,ptTriggerMax,
1602 ptAssociatedMin, ptAssociatedMax,gHistNN);
a38fd7f3 1603}
1604
3b8a460b 1605//____________________________________________________________//
2340240c 1606void drawProjections(Bool_t kProjectInEta = kFALSE,
41f2bc59 1607 Int_t binMin = 1,
1608 Int_t binMax = 80,
3b8a460b 1609 Int_t gCentrality = 1,
1610 Double_t psiMin = -0.5,
1611 Double_t psiMax = 3.5,
1612 Double_t ptTriggerMin = -1.,
1613 Double_t ptTriggerMax = -1.,
1614 Double_t ptAssociatedMin = -1.,
41f2bc59 1615 Double_t ptAssociatedMax = -1.,
2340240c 1616 Bool_t kUseZYAM = kFALSE,
1617 TString eventClass = "Centrality") {
2b4abbc5 1618 //Macro that draws the charge dependent correlation functions PROJECTIONS
1619 //for each centrality bin for the different pT of trigger and
1620 //associated particles
1621 TGaxis::SetMaxDigits(3);
1622
1623 //Get the input file
09b14834 1624 /*TString filename = lhcPeriod;
3b8a460b 1625 filename += "/Centrality"; filename += gCentralityEstimator;
1626 filename += "_Bit"; filename += gBit;
1627 filename += "_"; filename += gEventPlaneEstimator;
1628 filename +="/PttFrom";
1629 filename += Form("%.1f",ptTriggerMin); filename += "To";
1630 filename += Form("%.1f",ptTriggerMax); filename += "PtaFrom";
1631 filename += Form("%.1f",ptAssociatedMin); filename += "To";
09b14834 1632 filename += Form("%.1f",ptAssociatedMax); */
1633
2340240c 1634 TString filename = "correlationFunction.";
1635 if(eventClass == "Centrality"){
1636 filename += Form("Centrality%.1fTo%.1f",psiMin,psiMax);
1637 filename += ".PsiAll.PttFrom";
1638 }
1639 else if(eventClass == "Multiplicity"){
1640 filename += Form("Multiplicity%.0fTo%.0f",psiMin,psiMax);
1641 filename += ".PsiAll.PttFrom";
1642 }
1643 else{ // "EventPlane" (default)
1644 filename += "Centrality";
1645 filename += gCentrality; filename += ".Psi";
1646 if((psiMin == -0.5)&&(psiMax == 0.5)) filename += "InPlane.Ptt";
1647 else if((psiMin == 0.5)&&(psiMax == 1.5)) filename += "Intermediate.Ptt";
1648 else if((psiMin == 1.5)&&(psiMax == 2.5)) filename += "OutOfPlane.Ptt";
1649 else if((psiMin == 2.5)&&(psiMax == 3.5)) filename += "Rest.PttFrom";
1650 else filename += "All.PttFrom";
1651 }
2b4abbc5 1652 filename += Form("%.1f",ptTriggerMin); filename += "To";
1653 filename += Form("%.1f",ptTriggerMax); filename += "PtaFrom";
1654 filename += Form("%.1f",ptAssociatedMin); filename += "To";
1655 filename += Form("%.1f",ptAssociatedMax);
2340240c 1656 //if(k2pMethod) filename += "_2pMethod";
09b14834 1657
1658 filename += "_";
1659 filename += Form("%.1f",psiMin);
1660 filename += "-";
1661 filename += Form("%.1f",psiMax);
2340240c 1662
2b4abbc5 1663 filename += ".root";
1664
1665 //Open the file
1666 TFile *f = TFile::Open(filename.Data());
1667 if((!f)||(!f->IsOpen())) {
1668 Printf("The file %s is not found. Aborting...",filename);
1669 return listBF;
1670 }
1671 //f->ls();
1672
1673 //Latex
1674 TString centralityLatex = "Centrality: ";
2340240c 1675 if(eventClass == "Centrality")
1676 centralityLatex += Form("%.0f - %.0f",psiMin,psiMax);
1677 else
1678 centralityLatex += centralityArray[gCentrality-1];
1679 centralityLatex += " %";
2b4abbc5 1680
1681 TString psiLatex;
1682 if((psiMin == -0.5)&&(psiMax == 0.5))
3b8a460b 1683 psiLatex = " -7.5^{o} < #varphi^{t} - #Psi_{2} < 7.5^{o}";
2b4abbc5 1684 else if((psiMin == 0.5)&&(psiMax == 1.5))
3b8a460b 1685 psiLatex = " 37.5^{o} < #varphi^{t} - #Psi_{2} < 52.5^{o}";
2b4abbc5 1686 else if((psiMin == 1.5)&&(psiMax == 2.5))
3b8a460b 1687 psiLatex = " 82.5^{o} < #varphi^{t} - #Psi_{2} < 97.5^{o}";
2b4abbc5 1688 else
3b8a460b 1689 psiLatex = " 0^{o} < #varphi^{t} - #Psi_{2} < 180^{o}";
2b4abbc5 1690
1691 TString pttLatex = Form("%.1f",ptTriggerMin);
1692 pttLatex += " < p_{T,trig} < "; pttLatex += Form("%.1f",ptTriggerMax);
1693 pttLatex += " GeV/c";
1694
1695 TString ptaLatex = Form("%.1f",ptAssociatedMin);
1696 ptaLatex += " < p_{T,assoc} < "; ptaLatex += Form("%.1f",ptAssociatedMax);
1697 ptaLatex += " GeV/c";
1698
1699 TLatex *latexInfo1 = new TLatex();
1700 latexInfo1->SetNDC();
1701 latexInfo1->SetTextSize(0.045);
1702 latexInfo1->SetTextColor(1);
1703
1704 TString pngName;
1705
1706 //============================================================//
1707 //Get the +- correlation function
1708 TH2D *gHistPN = dynamic_cast<TH2D *>(f->Get("gHistPNCorrelationFunctions"));
1709 gHistPN->SetStats(kFALSE);
1710 gHistPN->SetTitle("C_{+-}(#Delta#eta,#Delta#varphi)");
3b8a460b 1711 gHistPN->GetXaxis()->SetRangeUser(-1.4,1.4);
2b4abbc5 1712 gHistPN->GetXaxis()->CenterTitle();
1713 gHistPN->GetXaxis()->SetTitleOffset(1.2);
1714 gHistPN->GetYaxis()->CenterTitle();
1715 gHistPN->GetYaxis()->SetTitleOffset(1.2);
1716 gHistPN->GetZaxis()->SetTitleOffset(1.5);
1717 TCanvas *cPN = new TCanvas("cPN","",0,0,600,500);
1718 cPN->SetFillColor(10);
1719 cPN->SetHighLightColor(10);
1720 cPN->SetLeftMargin(0.15);
2b4abbc5 1721
1722 //================//
3b8a460b 1723 TH1D* gHistPNprojection = 0x0;
1724 Double_t sum = 0.0;
1725 Double_t gError = 0.0;
1726 Int_t nCounter = 0;
41f2bc59 1727 //projection in delta eta
3b8a460b 1728 if(kProjectInEta) {
1729 gHistPNprojection = new TH1D("gHistPNprojection","",gHistPN->GetNbinsX(),gHistPN->GetXaxis()->GetXmin(),gHistPN->GetXaxis()->GetXmax());
1730 for(Int_t iBinX = 1; iBinX <= gHistPN->GetNbinsX(); iBinX++) {
1731 sum = 0.; gError = 0.0; nCounter = 0;
41f2bc59 1732 for(Int_t iBinY = binMin; iBinY <= binMax; iBinY++) {
3b8a460b 1733 sum += gHistPN->GetBinContent(iBinX,iBinY);
1734 if(gHistPN->GetBinContent(iBinX,iBinY) != 0.) nCounter += 1;
1735 Double_t exy = gHistPN->GetCellError(iBinX,iBinY);
1736 gError += exy*exy;
1737 }
1738 if(nCounter != 0) {
1739 sum /= nCounter;
1740 gError = TMath::Sqrt(gError)/nCounter;
1741 }
1742 gHistPNprojection->SetBinContent(iBinX,sum);
1743 gHistPNprojection->SetBinError(iBinX,gError);
1744 }
1745 gHistPNprojection->GetXaxis()->SetRangeUser(-1.4,1.4);
41f2bc59 1746 if(kUseZYAM)
1747 gHistPNprojection->GetYaxis()->SetTitle("#frac{1}{N_{trig}}#frac{dN_{assoc}}{#Delta#eta} - b_{ZYAM}");
1748 else
1749 gHistPNprojection->GetYaxis()->SetTitle("C_{+-}(#Delta#eta)");
3b8a460b 1750 gHistPNprojection->GetXaxis()->SetTitle("#Delta#eta");
41f2bc59 1751 }//projection in delta eta
1752 //projection in delta phi
3b8a460b 1753 else {
1754 gHistPNprojection = new TH1D("gHistPNprojection","",gHistPN->GetNbinsY(),gHistPN->GetYaxis()->GetXmin(),gHistPN->GetYaxis()->GetXmax());
1755 for(Int_t iBinY = 1; iBinY <= gHistPN->GetNbinsY(); iBinY++) {
1756 sum = 0.; gError = 0.0; nCounter = 0;
09b14834 1757 for(Int_t iBinX = binMin; iBinX <= binMax; iBinX++) {
1758 //for(Int_t iBinX = 1; iBinX <= gHistPN->GetNbinsX(); iBinX++) {
3b8a460b 1759 sum += gHistPN->GetBinContent(iBinX,iBinY);
1760 if(gHistPN->GetBinContent(iBinX,iBinY) != 0.) nCounter += 1;
1761 Double_t exy = gHistPN->GetCellError(iBinX,iBinY);
1762 gError += exy*exy;
1763 }
1764 if(nCounter != 0) {
1765 sum /= nCounter;
1766 gError = TMath::Sqrt(gError)/nCounter;
1767 }
1768 gHistPNprojection->SetBinContent(iBinY,sum);
1769 gHistPNprojection->SetBinError(iBinY,gError);
1770 }
41f2bc59 1771 if(kUseZYAM)
1772 gHistPNprojection->GetYaxis()->SetTitle("#frac{1}{N_{trig}}#frac{dN_{assoc}}{#Delta#varphi} - b_{ZYAM} (rad^{-1})");
1773 else
1774 gHistPNprojection->GetYaxis()->SetTitle("C_{+-}(#Delta#varphi)");
3b8a460b 1775 gHistPNprojection->GetXaxis()->SetTitle("#Delta#varphi (rad)");
2b4abbc5 1776 }
41f2bc59 1777
3b8a460b 1778 //ZYAM
41f2bc59 1779 if(kUseZYAM) {
1780 Double_t reference = gHistPNprojection->GetBinContent(gHistPNprojection->GetMinimumBin());
1781 for(Int_t iBinX = 1; iBinX <= gHistPNprojection->GetNbinsX(); iBinX++)
1782 gHistPNprojection->SetBinContent(iBinX,gHistPNprojection->GetBinContent(iBinX) - reference);
1783 }
3b8a460b 1784
2b4abbc5 1785 gHistPNprojection->GetYaxis()->SetTitleOffset(1.5);
3b8a460b 1786 gHistPNprojection->SetMarkerStyle(20);
1787 gHistPNprojection->SetStats(kFALSE);
1788 gHistPNprojection->DrawCopy("E");
2b4abbc5 1789 //=================//
1790
1791 gPad->SetTheta(30); // default is 30
1792 gPad->SetPhi(-60); // default is 30
1793 gPad->Update();
1794
3b8a460b 1795 latexInfo1->DrawLatex(0.6,0.95,centralityLatex.Data());
1796 latexInfo1->DrawLatex(0.6,0.89,psiLatex.Data());
1797 latexInfo1->DrawLatex(0.6,0.83,pttLatex.Data());
1798 latexInfo1->DrawLatex(0.6,0.77,ptaLatex.Data());
2b4abbc5 1799
2340240c 1800 pngName = filename;
1801 if(kProjectInEta) pngName.ReplaceAll(".root","_DeltaEtaProjection");
1802 else pngName.ReplaceAll(".root","_DeltaPhiProjection");
2b4abbc5 1803 pngName += ".PositiveNegative.png";
1804 cPN->SaveAs(pngName.Data());
1805
1806 //============================================================//
1807 //Get the -+ correlation function
1808 TH2D *gHistNP = dynamic_cast<TH2D *>(f->Get("gHistNPCorrelationFunctions"));
1809 gHistNP->SetStats(kFALSE);
1810 gHistNP->SetTitle("C_{-+}(#Delta#eta,#Delta#varphi)");
3b8a460b 1811 gHistNP->GetXaxis()->SetRangeUser(-1.4,1.4);
2b4abbc5 1812 gHistNP->GetXaxis()->CenterTitle();
1813 gHistNP->GetXaxis()->SetTitleOffset(1.2);
1814 gHistNP->GetYaxis()->CenterTitle();
1815 gHistNP->GetYaxis()->SetTitleOffset(1.2);
1816 gHistNP->GetZaxis()->SetTitleOffset(1.5);
1817 TCanvas *cNP = new TCanvas("cNP","",50,50,600,500);
1818 cNP->SetFillColor(10); cNP->SetHighLightColor(10);
1819 cNP->SetLeftMargin(0.15);
2b4abbc5 1820
1821 //================//
3b8a460b 1822 TH1D* gHistNPprojection = 0x0;
1823 Double_t sum = 0.0;
1824 Double_t gError = 0.0;
1825 Int_t nCounter = 0;
1826 if(kProjectInEta) {
1827 gHistNPprojection = new TH1D("gHistNPprojection","",gHistNP->GetNbinsX(),gHistNP->GetXaxis()->GetXmin(),gHistNP->GetXaxis()->GetXmax());
1828 for(Int_t iBinX = 1; iBinX <= gHistNP->GetNbinsX(); iBinX++) {
1829 sum = 0.; gError = 0.0; nCounter = 0;
09b14834 1830 for(Int_t iBinY = binMin; iBinY <= binMax; iBinY++) {
1831 //for(Int_t iBinY = 1; iBinY <= gHistNP->GetNbinsY(); iBinY++) {
3b8a460b 1832 sum += gHistNP->GetBinContent(iBinX,iBinY);
1833 if(gHistNP->GetBinContent(iBinX,iBinY) != 0.) nCounter += 1;
1834 Double_t exy = gHistNP->GetCellError(iBinX,iBinY);
1835 gError += exy*exy;
1836 }
1837 if(nCounter != 0) {
1838 sum /= nCounter;
1839 gError = TMath::Sqrt(gError)/nCounter;
1840 }
1841 gHistNPprojection->SetBinContent(iBinX,sum);
1842 gHistNPprojection->SetBinError(iBinX,gError);
1843 }
1844 gHistNPprojection->GetXaxis()->SetRangeUser(-1.4,1.4);
41f2bc59 1845 if(kUseZYAM)
1846 gHistNPprojection->GetYaxis()->SetTitle("#frac{1}{N_{trig}}#frac{dN_{assoc}}{#Delta#eta} - b_{ZYAM}");
1847 else
1848 gHistNPprojection->GetYaxis()->SetTitle("C_{-+}(#Delta#eta)");
3b8a460b 1849 gHistNPprojection->GetXaxis()->SetTitle("#Delta#eta");
2b4abbc5 1850 }
3b8a460b 1851 else {
1852 gHistNPprojection = new TH1D("gHistNPprojection","",gHistNP->GetNbinsY(),gHistNP->GetYaxis()->GetXmin(),gHistNP->GetYaxis()->GetXmax());
1853 for(Int_t iBinY = 1; iBinY <= gHistNP->GetNbinsY(); iBinY++) {
1854 sum = 0.; gError = 0.0; nCounter = 0;
09b14834 1855 for(Int_t iBinX = binMin; iBinX <= binMax; iBinX++) {
1856 //for(Int_t iBinX = 1; iBinX <= gHistNP->GetNbinsX(); iBinX++) {
3b8a460b 1857 sum += gHistNP->GetBinContent(iBinX,iBinY);
1858 if(gHistNP->GetBinContent(iBinX,iBinY) != 0.) nCounter += 1;
1859 Double_t exy = gHistNP->GetCellError(iBinX,iBinY);
1860 gError += exy*exy;
1861 }
1862 if(nCounter != 0) {
1863 sum /= nCounter;
1864 gError = TMath::Sqrt(gError)/nCounter;
1865 }
1866 gHistNPprojection->SetBinContent(iBinY,sum);
1867 gHistNPprojection->SetBinError(iBinY,gError);
1868 }
41f2bc59 1869 if(kUseZYAM)
1870 gHistNPprojection->GetYaxis()->SetTitle("#frac{1}{N_{trig}}#frac{dN_{assoc}}{#Delta#varphi} - b_{ZYAM} (rad^{-1})");
1871 else
1872 gHistNPprojection->GetYaxis()->SetTitle("C_{-+}(#Delta#varphi)");
3b8a460b 1873 gHistNPprojection->GetXaxis()->SetTitle("#Delta#varphi (rad)");
2b4abbc5 1874 }
3b8a460b 1875 //ZYAM
41f2bc59 1876 if(kUseZYAM) {
1877 Double_t reference = gHistNPprojection->GetBinContent(gHistNPprojection->GetMinimumBin());
3b8a460b 1878 for(Int_t iBinX = 1; iBinX <= gHistNPprojection->GetNbinsX(); iBinX++)
1879 gHistNPprojection->SetBinContent(iBinX,gHistNPprojection->GetBinContent(iBinX) - reference);
41f2bc59 1880 }
3b8a460b 1881
2b4abbc5 1882 gHistNPprojection->GetYaxis()->SetTitleOffset(1.5);
3b8a460b 1883 gHistNPprojection->SetMarkerStyle(20);
1884 gHistNPprojection->SetStats(kFALSE);
1885 gHistNPprojection->DrawCopy("E");
2b4abbc5 1886 //================//
1887
1888 gPad->SetTheta(30); // default is 30
1889 gPad->SetPhi(-60); // default is 30
1890 gPad->Update();
1891
3b8a460b 1892 latexInfo1->DrawLatex(0.6,0.95,centralityLatex.Data());
1893 latexInfo1->DrawLatex(0.6,0.89,psiLatex.Data());
1894 latexInfo1->DrawLatex(0.6,0.83,pttLatex.Data());
1895 latexInfo1->DrawLatex(0.6,0.77,ptaLatex.Data());
2b4abbc5 1896
2340240c 1897 pngName = filename;
1898 if(kProjectInEta) pngName.ReplaceAll(".root","_DeltaEtaProjection");
1899 else pngName.ReplaceAll(".root","_DeltaPhiProjection");
2b4abbc5 1900 pngName += ".NegativePositive.png";
1901 cNP->SaveAs(pngName.Data());
1902
1903 //============================================================//
1904 //Get the ++ correlation function
1905 TH2D *gHistPP = dynamic_cast<TH2D *>(f->Get("gHistPPCorrelationFunctions"));
1906 gHistPP->SetStats(kFALSE);
1907 gHistPP->SetTitle("C_{++}(#Delta#eta,#Delta#varphi)");
3b8a460b 1908 gHistPP->GetXaxis()->SetRangeUser(-1.4,1.4);
2b4abbc5 1909 gHistPP->GetXaxis()->CenterTitle();
1910 gHistPP->GetXaxis()->SetTitleOffset(1.2);
1911 gHistPP->GetYaxis()->CenterTitle();
1912 gHistPP->GetYaxis()->SetTitleOffset(1.2);
1913 gHistPP->GetZaxis()->SetTitleOffset(1.5);
1914 TCanvas *cPP = new TCanvas("cPP","",100,100,600,500);
1915 cPP->SetFillColor(10); cPP->SetHighLightColor(10);
1916 cPP->SetLeftMargin(0.15);
3b8a460b 1917
2b4abbc5 1918 //=================//
1919 TH1D* gHistPPprojection = 0x0;
3b8a460b 1920 Double_t sum = 0.0;
1921 Double_t gError = 0.0;
1922 Int_t nCounter = 0;
1923 if(kProjectInEta) {
1924 gHistPPprojection = new TH1D("gHistPPprojection","",gHistPP->GetNbinsX(),gHistPP->GetXaxis()->GetXmin(),gHistPP->GetXaxis()->GetXmax());
1925 for(Int_t iBinX = 1; iBinX <= gHistPP->GetNbinsX(); iBinX++) {
1926 sum = 0.; gError = 0.0; nCounter = 0;
09b14834 1927 for(Int_t iBinY = binMin; iBinY <= binMax; iBinY++) {
1928 //for(Int_t iBinY = 1; iBinY <= gHistPP->GetNbinsY(); iBinY++) {
3b8a460b 1929 sum += gHistPP->GetBinContent(iBinX,iBinY);
1930 if(gHistPP->GetBinContent(iBinX,iBinY) != 0.) nCounter += 1;
1931 Double_t exy = gHistPP->GetCellError(iBinX,iBinY);
1932 gError += exy*exy;
1933 }
1934 if(nCounter != 0) {
1935 sum /= nCounter;
1936 gError = TMath::Sqrt(gError)/nCounter;
1937 }
1938 gHistPPprojection->SetBinContent(iBinX,sum);
1939 gHistPPprojection->SetBinError(iBinX,gError);
1940 }
1941 gHistPPprojection->GetXaxis()->SetRangeUser(-1.4,1.4);
41f2bc59 1942 if(kUseZYAM)
1943 gHistPPprojection->GetYaxis()->SetTitle("#frac{1}{N_{trig}}#frac{dN_{assoc}}{#Delta#eta} - b_{ZYAM}");
1944 else
1945 gHistPPprojection->GetYaxis()->SetTitle("C_{++}(#Delta#eta)");
3b8a460b 1946 gHistPPprojection->GetXaxis()->SetTitle("#Delta#eta");
2b4abbc5 1947 }
3b8a460b 1948 else {
1949 gHistPPprojection = new TH1D("gHistPPprojection","",gHistPP->GetNbinsY(),gHistPP->GetYaxis()->GetXmin(),gHistPP->GetYaxis()->GetXmax());
1950 for(Int_t iBinY = 1; iBinY <= gHistPP->GetNbinsY(); iBinY++) {
1951 sum = 0.; gError = 0.0; nCounter = 0;
09b14834 1952 for(Int_t iBinX = binMin; iBinX <= binMax; iBinX++) {
1953 //for(Int_t iBinX = 1; iBinX <= gHistPP->GetNbinsX(); iBinX++) {
3b8a460b 1954 sum += gHistPP->GetBinContent(iBinX,iBinY);
1955 if(gHistPP->GetBinContent(iBinX,iBinY) != 0.) nCounter += 1;
1956 Double_t exy = gHistPP->GetCellError(iBinX,iBinY);
1957 gError += exy*exy;
1958 }
1959 if(nCounter != 0) {
1960 sum /= nCounter;
1961 gError = TMath::Sqrt(gError)/nCounter;
1962 }
1963 gHistPPprojection->SetBinContent(iBinY,sum);
1964 gHistPPprojection->SetBinError(iBinY,gError);
1965 }
41f2bc59 1966 if(kUseZYAM)
1967 gHistPPprojection->GetYaxis()->SetTitle("#frac{1}{N_{trig}}#frac{dN_{assoc}}{#Delta#varphi} - b_{ZYAM} (rad^{-1})");
1968 else
1969 gHistPPprojection->GetYaxis()->SetTitle("C_{++}(#Delta#varphi)");
3b8a460b 1970 gHistPPprojection->GetXaxis()->SetTitle("#Delta#varphi (rad)");
2b4abbc5 1971 }
3b8a460b 1972 //ZYAM
41f2bc59 1973 if(kUseZYAM) {
1974 Double_t reference = gHistPPprojection->GetBinContent(gHistPPprojection->GetMinimumBin());
1975 for(Int_t iBinX = 1; iBinX <= gHistPPprojection->GetNbinsX(); iBinX++)
1976 gHistPPprojection->SetBinContent(iBinX,gHistPPprojection->GetBinContent(iBinX) - reference);
1977 }
1978
2b4abbc5 1979 gHistPPprojection->GetYaxis()->SetTitleOffset(1.5);
3b8a460b 1980 gHistPPprojection->SetMarkerStyle(20);
1981 gHistPPprojection->SetStats(kFALSE);
1982 gHistPPprojection->DrawCopy("E");
2b4abbc5 1983 //================//
1984
1985 gPad->SetTheta(30); // default is 30
1986 gPad->SetPhi(-60); // default is 30
1987 gPad->Update();
1988
3b8a460b 1989 latexInfo1->DrawLatex(0.6,0.95,centralityLatex.Data());
1990 latexInfo1->DrawLatex(0.6,0.89,psiLatex.Data());
1991 latexInfo1->DrawLatex(0.6,0.83,pttLatex.Data());
1992 latexInfo1->DrawLatex(0.6,0.77,ptaLatex.Data());
2b4abbc5 1993
2340240c 1994 pngName = filename;
1995 if(kProjectInEta) pngName.ReplaceAll(".root","_DeltaEtaProjection");
1996 else pngName.ReplaceAll(".root","_DeltaPhiProjection");
2b4abbc5 1997 pngName += ".PositivePositive.png";
1998 cPP->SaveAs(pngName.Data());
1999
2000 //============================================================//
2001 //Get the -- correlation function
2002 TH2D *gHistNN = dynamic_cast<TH2D *>(f->Get("gHistNNCorrelationFunctions"));
2003 gHistNN->SetStats(kFALSE);
2004 gHistNN->SetTitle("C_{--}(#Delta#eta,#Delta#varphi)");
3b8a460b 2005 gHistNN->GetXaxis()->SetRangeUser(-1.4,1.4);
2b4abbc5 2006 gHistNN->GetXaxis()->CenterTitle();
2007 gHistNN->GetXaxis()->SetTitleOffset(1.2);
2008 gHistNN->GetYaxis()->CenterTitle();
2009 gHistNN->GetYaxis()->SetTitleOffset(1.2);
2010 gHistNN->GetZaxis()->SetTitleOffset(1.5);
2011 TCanvas *cNN = new TCanvas("cNN","",150,150,600,500);
2012 cNN->SetFillColor(10); cNN->SetHighLightColor(10);
2013 cNN->SetLeftMargin(0.15);
2b4abbc5 2014
2015 //=================//
3b8a460b 2016 TH1D* gHistNNprojection = 0x0;
2017 Double_t sum = 0.0;
2018 Double_t gError = 0.0;
2019 Int_t nCounter = 0;
2020 if(kProjectInEta) {
2021 gHistNNprojection = new TH1D("gHistNNprojection","",gHistNN->GetNbinsX(),gHistNN->GetXaxis()->GetXmin(),gHistNN->GetXaxis()->GetXmax());
2022 for(Int_t iBinX = 1; iBinX <= gHistNN->GetNbinsX(); iBinX++) {
2023 sum = 0.; gError = 0.0; nCounter = 0;
09b14834 2024 for(Int_t iBinY = binMin; iBinY <= binMax; iBinY++) {
2025 //for(Int_t iBinY = 1; iBinY <= gHistNN->GetNbinsY(); iBinY++) {
3b8a460b 2026 sum += gHistNN->GetBinContent(iBinX,iBinY);
2027 if(gHistNN->GetBinContent(iBinX,iBinY) != 0.) nCounter += 1;
2028 Double_t exy = gHistNN->GetCellError(iBinX,iBinY);
2029 gError += exy*exy;
2030 }
2031 if(nCounter != 0) {
2032 sum /= nCounter;
2033 gError = TMath::Sqrt(gError)/nCounter;
2034 }
2035 gHistNNprojection->SetBinContent(iBinX,sum);
2036 gHistNNprojection->SetBinError(iBinX,gError);
2037 }
2038 gHistNNprojection->GetXaxis()->SetRangeUser(-1.4,1.4);
41f2bc59 2039 if(kUseZYAM)
2040 gHistNNprojection->GetYaxis()->SetTitle("#frac{1}{N_{trig}}#frac{dN_{assoc}}{#Delta#eta} - b_{ZYAM}");
2041 else
2042 gHistNNprojection->GetYaxis()->SetTitle("C_{--}(#Delta#eta)");
3b8a460b 2043 gHistNNprojection->GetXaxis()->SetTitle("#Delta#eta");
2b4abbc5 2044 }
3b8a460b 2045 else {
2046 gHistNNprojection = new TH1D("gHistNNprojection","",gHistNN->GetNbinsY(),gHistNN->GetYaxis()->GetXmin(),gHistNN->GetYaxis()->GetXmax());
2047 for(Int_t iBinY = 1; iBinY <= gHistNN->GetNbinsY(); iBinY++) {
2048 sum = 0.; gError = 0.0; nCounter = 0;
09b14834 2049 for(Int_t iBinX = binMin; iBinX <= binMax; iBinX++) {
2050 //for(Int_t iBinX = 1; iBinX <= gHistNN->GetNbinsX(); iBinX++) {
3b8a460b 2051 sum += gHistNN->GetBinContent(iBinX,iBinY);
2052 if(gHistNN->GetBinContent(iBinX,iBinY) != 0.) nCounter += 1;
2053 Double_t exy = gHistNN->GetCellError(iBinX,iBinY);
2054 gError += exy*exy;
2055 }
2056 if(nCounter != 0) {
2057 sum /= nCounter;
2058 gError = TMath::Sqrt(gError)/nCounter;
2059 }
2060 gHistNNprojection->SetBinContent(iBinY,sum);
2061 gHistNNprojection->SetBinError(iBinY,gError);
2062 }
41f2bc59 2063 if(kUseZYAM)
2064 gHistNNprojection->GetYaxis()->SetTitle("#frac{1}{N_{trig}}#frac{dN_{assoc}}{#Delta#varphi} - b_{ZYAM} (rad^{-1})");
2065 else
2066 gHistNNprojection->GetYaxis()->SetTitle("C_{--}(#Delta#varphi)");
3b8a460b 2067 gHistNNprojection->GetXaxis()->SetTitle("#Delta#varphi (rad)");
2b4abbc5 2068 }
3b8a460b 2069 //ZYAM
41f2bc59 2070 if(kUseZYAM) {
2071 Double_t reference = gHistNNprojection->GetBinContent(gHistNNprojection->GetMinimumBin());
2072 for(Int_t iBinX = 1; iBinX <= gHistNNprojection->GetNbinsX(); iBinX++)
2073 gHistNNprojection->SetBinContent(iBinX,gHistNNprojection->GetBinContent(iBinX) - reference);
2074 }
2b4abbc5 2075
2b4abbc5 2076 gHistNNprojection->GetYaxis()->SetTitleOffset(1.5);
3b8a460b 2077 gHistNNprojection->SetMarkerStyle(20);
2078 gHistNNprojection->SetStats(kFALSE);
2079 gHistNNprojection->DrawCopy("E");
2b4abbc5 2080 //=================//
2081
2082 gPad->SetTheta(30); // default is 30
2083 gPad->SetPhi(-60); // default is 30
2084 gPad->Update();
2085
3b8a460b 2086 latexInfo1->DrawLatex(0.6,0.95,centralityLatex.Data());
2087 latexInfo1->DrawLatex(0.6,0.89,psiLatex.Data());
2088 latexInfo1->DrawLatex(0.6,0.83,pttLatex.Data());
2089 latexInfo1->DrawLatex(0.6,0.77,ptaLatex.Data());
2b4abbc5 2090
2340240c 2091 pngName = filename;
2092 if(kProjectInEta) pngName.ReplaceAll(".root","_DeltaEtaProjection");
2093 else pngName.ReplaceAll(".root","_DeltaPhiProjection");
2b4abbc5 2094 pngName += ".NegativeNegative.png";
2095 cNN->SaveAs(pngName.Data());
2096
17c69b4e 2097 TString outFileName = filename;
2098 if(kProjectInEta) outFileName.ReplaceAll(".root","_DeltaEtaProjection.root");
2099 else outFileName.ReplaceAll(".root","_DeltaPhiProjection.root");
2100 TFile *fProjection = TFile::Open(outFileName.Data(),"recreate");
3b8a460b 2101 gHistNPprojection->Write();
2102 gHistPNprojection->Write();
2103 gHistNNprojection->Write();
2104 gHistPPprojection->Write();
2b4abbc5 2105 fProjection->Close();
2b4abbc5 2106}
742af4bd 2107
5c9b4733 2108//____________________________________________________________//
2109void fitCorrelationFunctions(Int_t gCentrality = 1,
2110 Double_t psiMin = -0.5, Double_t psiMax = 3.5,
20006629 2111 Double_t vertexZMin = -10.,Double_t vertexZMax = 10.,
5c9b4733 2112 Double_t ptTriggerMin = -1.,
2113 Double_t ptTriggerMax = -1.,
2114 Double_t ptAssociatedMin = -1.,
2115 Double_t ptAssociatedMax = -1.,
2116 TH2D *gHist) {
07d0a35c 2117
2118 cout<<"FITTING FUNCTION"<<endl;
2119
5c9b4733 2120 //near side peak: [1]*TMath::Exp(-TMath::Power((0.5*TMath::Power((x/[2]),2)+0.5*TMath::Power((y/[3]),2)),[4]))
2121 //away side ridge: [5]*TMath::Exp(-TMath::Power((0.5*TMath::Power(((y-TMath::Pi())/[6]),2)),[7]))
2122 //longitudinal ridge: [8]*TMath::Exp(-TMath::Power((0.5*TMath::Power((x/[9]),2)),[10]))
2123 //wing structures: [11]*TMath::Power(x,2)
2124 //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))
24ca71a6 2125 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.);
5c9b4733 2126 gFitFunction->SetName("gFitFunction");
742af4bd 2127
2128
5c9b4733 2129 //Normalization
742af4bd 2130 gFitFunction->SetParName(0,"N1");
5c9b4733 2131 //near side peak
742af4bd 2132 gFitFunction->SetParName(1,"N_{near side}");gFitFunction->FixParameter(1,0.0);
2133 gFitFunction->SetParName(2,"Sigma_{near side}(delta eta)"); gFitFunction->FixParameter(2,0.0);
2134 gFitFunction->SetParName(3,"Sigma_{near side}(delta phi)"); gFitFunction->FixParameter(3,0.0);
2135 gFitFunction->SetParName(4,"Exponent_{near side}"); gFitFunction->FixParameter(4,1.0);
5c9b4733 2136 //away side ridge
742af4bd 2137 gFitFunction->SetParName(5,"N_{away side}"); gFitFunction->FixParameter(5,0.0);
2138 gFitFunction->SetParName(6,"Sigma_{away side}(delta phi)"); gFitFunction->FixParameter(6,0.0);
2139 gFitFunction->SetParName(7,"Exponent_{away side}"); gFitFunction->FixParameter(7,1.0);
24ca71a6 2140 //longitudinal ridge
742af4bd 2141 gFitFunction->SetParName(8,"N_{long. ridge}"); gFitFunction->SetParameter(8,0.05);//
2142 gFitFunction->FixParameter(8,0.0);
2143 gFitFunction->SetParName(9,"Sigma_{long. ridge}(delta eta)"); gFitFunction->FixParameter(9,0.0);
2144 gFitFunction->SetParName(10,"Exponent_{long. ridge}"); gFitFunction->FixParameter(10,1.0);
5c9b4733 2145 //wing structures
742af4bd 2146 gFitFunction->SetParName(11,"N_{wing}"); gFitFunction->FixParameter(11,0.0);
5c9b4733 2147 //flow contribution
742af4bd 2148 gFitFunction->SetParName(12,"N_{flow}"); gFitFunction->SetParameter(12,0.2);gFitFunction->SetParLimits(12,0.0,10);
2149 gFitFunction->SetParName(13,"V1"); gFitFunction->SetParameter(13,0.005);gFitFunction->SetParLimits(13,0.0,10);
2150 gFitFunction->SetParName(14,"V2"); gFitFunction->SetParameter(14,0.1);gFitFunction->SetParLimits(14,0.0,10);
2151 gFitFunction->SetParName(15,"V3"); gFitFunction->SetParameter(15,0.05);gFitFunction->SetParLimits(15,0.0,10);
2152 gFitFunction->SetParName(16,"V4"); gFitFunction->SetParameter(16,0.005);gFitFunction->SetParLimits(16,0.0,10);
2153 gFitFunction->SetParName(17,"Offset"); gFitFunction->FixParameter(17,0.0);
2154
2155 // flow parameters
2156 Double_t fNV = 0.;
2157 Double_t fV1 = 0.;
2158 Double_t fV2 = 0.;
2159 Double_t fV3 = 0.;
2160 Double_t fV4 = 0.;
2161
2162 //Fitting the correlation function (first the edges to extract flow)
2163 gHist->Fit("gFitFunction","nm","",1.0,1.6);
2164
2165 fNV += gFitFunction->GetParameter(12);
2166 fV1 += gFitFunction->GetParameter(13);
2167 fV2 += gFitFunction->GetParameter(14);
2168 fV3 += gFitFunction->GetParameter(15);
2169 fV4 += gFitFunction->GetParameter(16);
2170
2171 gHist->Fit("gFitFunction","nm","",-1.6,-1.0);
2172
2173 fNV += gFitFunction->GetParameter(12);
2174 fV1 += gFitFunction->GetParameter(13);
2175 fV2 += gFitFunction->GetParameter(14);
2176 fV3 += gFitFunction->GetParameter(15);
2177 fV4 += gFitFunction->GetParameter(16);
2178
2179 // Now fit the whole with fixed flow
2180 gFitFunction->FixParameter(12,fNV/2.);
2181 gFitFunction->FixParameter(13,fV1/2.);
2182 gFitFunction->FixParameter(14,fV2/2.);
2183 gFitFunction->FixParameter(15,fV3/2.);
2184 gFitFunction->FixParameter(16,fV4/2.);
2185
2186 gFitFunction->ReleaseParameter(0);gFitFunction->SetParameter(0,1.0);
2187 gFitFunction->ReleaseParameter(1);gFitFunction->SetParameter(1,0.3);
2188 gFitFunction->ReleaseParameter(2);gFitFunction->SetParameter(2,0.3);gFitFunction->SetParLimits(2,0.05,0.7);
2189 gFitFunction->ReleaseParameter(3);gFitFunction->SetParameter(3,0.3);gFitFunction->SetParLimits(3,0.05,1.7);
2190 //gFitFunction->ReleaseParameter(4);gFitFunction->SetParameter(4,1.0);gFitFunction->SetParLimits(4,0.0,2.0);
2191 gFitFunction->ReleaseParameter(5);gFitFunction->SetParameter(5,1.0);//gFitFunction->SetParLimits(5,0.0,10);
2192 gFitFunction->ReleaseParameter(6);gFitFunction->SetParameter(6,0.5);//gFitFunction->SetParLimits(6,0.0,10);
2193 //gFitFunction->ReleaseParameter(7);gFitFunction->SetParameter(7,1.0);gFitFunction->SetParLimits(7,0.0,2.0);
478c95cf 2194 //gFitFunction->ReleaseParameter(8);gFitFunction->SetParameter(8,0.05);
2195 //gFitFunction->ReleaseParameter(9);gFitFunction->SetParameter(9,0.6);gFitFunction->SetParLimits(9,0.1,10.0);
742af4bd 2196 //gFitFunction->ReleaseParameter(10);gFitFunction->SetParameter(10,1.0);gFitFunction->SetParLimits(10,0.0,2.0);
2197 gFitFunction->ReleaseParameter(17);gFitFunction->SetParameter(17,0.7);gFitFunction->SetParLimits(17,0.0,0.9);
2198
5c9b4733 2199 gHist->Fit("gFitFunction","nm");
2200
742af4bd 2201
5c9b4733 2202 //Cloning the histogram
2203 TString histoName = gHist->GetName(); histoName += "Fit";
2204 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());
2205 TH2D *gHistResidual = dynamic_cast<TH2D *>(gHist->Clone());
2206 gHistResidual->SetName("gHistResidual");
2207 gHistResidual->Sumw2();
2208
2209 for(Int_t iBinDeltaEta = 1; iBinDeltaEta <= gHist->GetNbinsX(); iBinDeltaEta++) {
2210 for(Int_t iBinDeltaPhi = 1; iBinDeltaPhi <= gHist->GetNbinsY(); iBinDeltaPhi++) {
2211 gHistFit->SetBinContent(iBinDeltaEta,iBinDeltaPhi,gFitFunction->Eval(gHist->GetXaxis()->GetBinCenter(iBinDeltaEta),gHist->GetYaxis()->GetBinCenter(iBinDeltaPhi)));
2212 }
2213 }
2214 gHistResidual->Add(gHistFit,-1);
2215
2216 //Write to output file
2217 TString newFileName = "correlationFunctionFit";
2218 if(histoName.Contains("PN")) newFileName += "PN";
2219 else if(histoName.Contains("NP")) newFileName += "NP";
2220 else if(histoName.Contains("PP")) newFileName += "PP";
2221 else if(histoName.Contains("NN")) newFileName += "NN";
2222 newFileName += ".Centrality";
2223 newFileName += gCentrality; newFileName += ".Psi";
2224 if((psiMin == -0.5)&&(psiMax == 0.5)) newFileName += "InPlane.Ptt";
2225 else if((psiMin == 0.5)&&(psiMax == 1.5)) newFileName += "Intermediate.Ptt";
2226 else if((psiMin == 1.5)&&(psiMax == 2.5)) newFileName += "OutOfPlane.Ptt";
2227 else if((psiMin == 2.5)&&(psiMax == 3.5)) newFileName += "Rest.PttFrom";
2228 else newFileName += "All.PttFrom";
2229 newFileName += Form("%.1f",ptTriggerMin); newFileName += "To";
2230 newFileName += Form("%.1f",ptTriggerMax); newFileName += "PtaFrom";
2231 newFileName += Form("%.1f",ptAssociatedMin); newFileName += "To";
2232 newFileName += Form("%.1f",ptAssociatedMax);
09b14834 2233
2234 newFileName += "_";
2235 newFileName += Form("%.1f",psiMin);
2236 newFileName += "-";
2237 newFileName += Form("%.1f",psiMax);
2238
5c9b4733 2239 newFileName += ".root";
2240 TFile *newFile = TFile::Open(newFileName.Data(),"recreate");
5ea3b761 2241 gHist->Write();
5c9b4733 2242 gHistFit->Write();
2243 gHistResidual->Write();
2244 gFitFunction->Write();
2245 newFile->Close();
2246
2247
2248}
24ca71a6 2249
3b8a460b 2250//____________________________________________________________//
2251TH2D *convolute2D(TH2D* h1, TH2D* h2, TString hname) {
3ce45c12 2252 //
2253 // this function does the convolution of 2 eta or phi "efficiencies" in a deta or dphi "efficiency"
2254 // and returns a new histogram which is normalized to the number of bin combinations
2255
2256 // histogram definition
2257 TH2D *hConv = NULL;
1b5c6300 2258 hConv = new TH2D(hname.Data(),hname.Data(), h1->GetNbinsY(),-2,2,h1->GetNbinsX(),-TMath::Pi()/2.,3*TMath::Pi()/2.);
3ce45c12 2259
2260 Double_t x1 = 0.;
2261 Double_t x2 = 0.;
2262 Double_t y1 = 0.;
2263 Double_t y2 = 0.;
2264 Double_t z1 = 0.;
2265 Double_t z2 = 0.;
2266 Double_t dx = 0.;
2267 Double_t dy = 0.;
2268 Double_t dz = 0.;
2269
2270
2271 // convolution
2272 for(Int_t i = 0; i < h1->GetNbinsX(); i ++){
2273 cout<<"CONVOLUTING BIN = "<<i<<" of "<<h1->GetNbinsX()<<endl;
2274 for(Int_t j = 0; j < h1->GetNbinsY(); j ++){
2275 for(Int_t k = 0; k < h2->GetNbinsX(); k ++){
2276 for(Int_t l = 0; l < h2->GetNbinsY(); l ++){
2277
2278 x1 = (Double_t)h1->GetXaxis()->GetBinCenter(i+1);
2279 y1 = (Double_t)h1->GetYaxis()->GetBinCenter(j+1);
2280 x2 = (Double_t)h2->GetXaxis()->GetBinCenter(k+1);
2281 y2 = (Double_t)h2->GetYaxis()->GetBinCenter(l+1);
2282
2283 z1 = (Double_t)h1->GetBinContent(i+1,j+1);
2284 z2 = (Double_t)h2->GetBinContent(k+1,l+1);
2285
1b5c6300 2286 // need the gymnastics to keep the same binning
2287 dx = x1 - x2 - (h1->GetXaxis()->GetBinWidth(1)/2.);
2288 if(gRandom->Gaus() > 0) dy = y1 - y2 + (h1->GetYaxis()->GetBinWidth(1)/2.);
2289 else dy = y1 - y2 - (h1->GetYaxis()->GetBinWidth(1)/2.);
3ce45c12 2290
2291 if(dx>3./2.*TMath::Pi()) dx = dx - 2.*TMath::Pi();
2292 if(dx<-1./2.*TMath::Pi()) dx = 2*TMath::Pi() + dx;
2293
2294 dz = z1*z2;
2295
2296 hConv->Fill(dy,dx,dz);
2297
2298 }
2299 }
2300 }
2301 }
2302
2303 return hConv;
2304
2305}