1 /**************************************************************************
2 * Author: Panos Christakoglou. *
3 * Contributors are mentioned in the code where appropriate. *
5 * Permission to use, copy, modify and distribute this software and its *
6 * documentation strictly for non-commercial purposes is hereby granted *
7 * without fee, provided that the above copyright notice appears in all *
8 * copies and that both the copyright notice and this permission notice *
9 * appear in the supporting documentation. The authors make no claims *
10 * about the suitability of this software for any purpose. It is *
11 * provided "as is" without express or implied warranty. *
12 **************************************************************************/
14 /* $Id: AliBalancePsi.cxx 54125 2012-01-24 21:07:41Z miweber $ */
16 //-----------------------------------------------------------------
17 // Balance Function class
18 // This is the class to deal with the Balance Function wrt Psi analysis
19 // Origin: Panos Christakoglou, Nikhef, Panos.Christakoglou@cern.ch
20 //-----------------------------------------------------------------
24 #include <Riostream.h>
30 #include <TLorentzVector.h>
31 #include <TObjArray.h>
32 #include <TGraphErrors.h>
35 #include "AliVParticle.h"
36 #include "AliMCParticle.h"
37 #include "AliESDtrack.h"
38 #include "AliAODTrack.h"
41 #include "AliBalancePsi.h"
43 ClassImp(AliBalancePsi)
45 //____________________________________________________________________//
46 AliBalancePsi::AliBalancePsi() :
49 fAnalysisLevel("ESD"),
62 fHistConversionbefore(0),
63 fHistConversionafter(0),
67 fConversionCut(kFALSE) {
68 // Default constructor
71 //____________________________________________________________________//
72 AliBalancePsi::AliBalancePsi(const AliBalancePsi& balance):
73 TObject(balance), fShuffle(balance.fShuffle),
74 fAnalysisLevel(balance.fAnalysisLevel),
75 fAnalyzedEvents(balance.fAnalyzedEvents),
76 fCentralityId(balance.fCentralityId),
77 fCentStart(balance.fCentStart),
78 fCentStop(balance.fCentStop),
79 fHistP(balance.fHistP),
80 fHistN(balance.fHistN),
81 fHistPN(balance.fHistPN),
82 fHistNP(balance.fHistNP),
83 fHistPP(balance.fHistPP),
84 fHistNN(balance.fHistNN),
85 fHistHBTbefore(balance.fHistHBTbefore),
86 fHistHBTafter(balance.fHistHBTafter),
87 fHistConversionbefore(balance.fHistConversionbefore),
88 fHistConversionafter(balance.fHistConversionafter),
89 fHistPsiMinusPhi(balance.fHistPsiMinusPhi),
90 fPsiInterval(balance.fPsiInterval),
91 fHBTCut(balance.fHBTCut),
92 fConversionCut(balance.fConversionCut) {
96 //____________________________________________________________________//
97 AliBalancePsi::~AliBalancePsi() {
106 delete fHistHBTbefore;
107 delete fHistHBTafter;
108 delete fHistConversionbefore;
109 delete fHistConversionafter;
110 delete fHistPsiMinusPhi;
113 //____________________________________________________________________//
114 void AliBalancePsi::InitHistograms() {
115 // single particle histograms
117 // global switch disabling the reference
118 // (to avoid "Replacing existing TH1" if several wagons are created in train)
119 Bool_t oldStatus = TH1::AddDirectoryStatus();
120 TH1::AddDirectory(kFALSE);
122 Int_t anaSteps = 1; // analysis steps
123 Int_t iBinSingle[kTrackVariablesSingle]; // binning for track variables
124 Double_t* dBinsSingle[kTrackVariablesSingle]; // bins for track variables
125 TString axisTitleSingle[kTrackVariablesSingle]; // axis titles for track variables
127 // two particle histograms
128 Int_t iBinPair[kTrackVariablesPair]; // binning for track variables
129 Double_t* dBinsPair[kTrackVariablesPair]; // bins for track variables
130 TString axisTitlePair[kTrackVariablesPair]; // axis titles for track variables
133 /*const Int_t kNCentralityBins = 9;
134 Double_t centralityBins[kNCentralityBins+1] = {0.,5.,10.,20.,30.,40.,50.,60.,70.,80.};
135 iBinSingle[0] = kNCentralityBins;
136 dBinsSingle[0] = centralityBins;
137 axisTitleSingle[0] = "Centrality percentile [%]";
138 iBinPair[0] = kNCentralityBins;
139 dBinsPair[0] = centralityBins;
140 axisTitlePair[0] = "Centrality percentile [%]"; */
142 //Psi_2: -0.5->0.5 (in plane), 0.5->1.5 (intermediate), 1.5->2.5 (out of plane), 2.5->3.5 (rest)
143 const Int_t kNPsi2Bins = 4;
144 Double_t psi2Bins[kNPsi2Bins+1] = {-0.5,0.5,1.5,2.5,3.5};
145 iBinSingle[0] = kNPsi2Bins;
146 dBinsSingle[0] = psi2Bins;
147 axisTitleSingle[0] = "#varphi - #Psi_{2} (a.u.)";
148 iBinPair[0] = kNPsi2Bins;
149 dBinsPair[0] = psi2Bins;
150 axisTitlePair[0] = "#varphi - #Psi_{2} (a.u.)";
153 const Int_t kNDeltaEtaBins = 80;
154 Double_t deltaEtaBins[kNDeltaEtaBins+1];
155 for(Int_t i = 0; i < kNDeltaEtaBins+1; i++)
156 deltaEtaBins[i] = -2.0 + i * 0.05;
157 iBinPair[1] = kNDeltaEtaBins;
158 dBinsPair[1] = deltaEtaBins;
159 axisTitlePair[1] = "#Delta#eta";
162 const Int_t kNDeltaPhiBins = 72;
163 Double_t deltaPhiBins[kNDeltaPhiBins+1];
164 for(Int_t i = 0; i < kNDeltaPhiBins+1; i++){
165 //deltaPhiBins[i] = -180.0 + i * 5.;
166 deltaPhiBins[i] = -TMath::Pi()/2. + i * 5.*TMath::Pi()/180.;
168 iBinPair[2] = kNDeltaPhiBins;
169 dBinsPair[2] = deltaPhiBins;
170 axisTitlePair[2] = "#Delta#varphi (rad)";
172 // pt(trigger-associated)
173 const Int_t kNPtBins = 16;
174 Double_t ptBins[kNPtBins+1] = {0.2,0.6,1.0,1.5,2.0,2.5,3.0,3.5,4.0,5.0,6.0,7.0,8.0,10.,12.,15.,20.};
175 //for(Int_t i = 0; i < kNPtBins+1; i++){
176 //ptBins[i] = 0.2 + i * 0.5;
178 iBinSingle[1] = kNPtBins;
179 dBinsSingle[1] = ptBins;
180 axisTitleSingle[1] = "p_{T,trig.} (GeV/c)";
182 iBinPair[3] = kNPtBins;
183 dBinsPair[3] = ptBins;
184 axisTitlePair[3] = "p_{T,trig.} (GeV/c)";
186 iBinPair[4] = kNPtBins;
187 dBinsPair[4] = ptBins;
188 axisTitlePair[4] = "p_{T,assoc.} (GeV/c)";
191 //+ triggered particles
193 if(fShuffle) histName.Append("_shuffle");
194 if(fCentralityId) histName += fCentralityId.Data();
195 fHistP = new AliTHn(histName.Data(),histName.Data(),anaSteps,kTrackVariablesSingle,iBinSingle);
196 for (Int_t j=0; j<kTrackVariablesSingle; j++) {
197 fHistP->SetBinLimits(j, dBinsSingle[j]);
198 fHistP->SetVarTitle(j, axisTitleSingle[j]);
201 //- triggered particles
203 if(fShuffle) histName.Append("_shuffle");
204 if(fCentralityId) histName += fCentralityId.Data();
205 fHistN = new AliTHn(histName.Data(),histName.Data(),anaSteps,kTrackVariablesSingle,iBinSingle);
206 for (Int_t j=0; j<kTrackVariablesSingle; j++) {
207 fHistN->SetBinLimits(j, dBinsSingle[j]);
208 fHistN->SetVarTitle(j, axisTitleSingle[j]);
212 histName = "fHistPN";
213 if(fShuffle) histName.Append("_shuffle");
214 if(fCentralityId) histName += fCentralityId.Data();
215 fHistPN = new AliTHn(histName.Data(),histName.Data(),anaSteps, kTrackVariablesPair, iBinPair);
216 for (Int_t j=0; j<kTrackVariablesPair; j++) {
217 fHistPN->SetBinLimits(j, dBinsPair[j]);
218 fHistPN->SetVarTitle(j, axisTitlePair[j]);
222 histName = "fHistNP";
223 if(fShuffle) histName.Append("_shuffle");
224 if(fCentralityId) histName += fCentralityId.Data();
225 fHistNP = new AliTHn(histName.Data(),histName.Data(),anaSteps, kTrackVariablesPair, iBinPair);
226 for (Int_t j=0; j<kTrackVariablesPair; j++) {
227 fHistNP->SetBinLimits(j, dBinsPair[j]);
228 fHistNP->SetVarTitle(j, axisTitlePair[j]);
232 histName = "fHistPP";
233 if(fShuffle) histName.Append("_shuffle");
234 if(fCentralityId) histName += fCentralityId.Data();
235 fHistPP = new AliTHn(histName.Data(),histName.Data(),anaSteps, kTrackVariablesPair, iBinPair);
236 for (Int_t j=0; j<kTrackVariablesPair; j++) {
237 fHistPP->SetBinLimits(j, dBinsPair[j]);
238 fHistPP->SetVarTitle(j, axisTitlePair[j]);
242 histName = "fHistNN";
243 if(fShuffle) histName.Append("_shuffle");
244 if(fCentralityId) histName += fCentralityId.Data();
245 fHistNN = new AliTHn(histName.Data(),histName.Data(),anaSteps, kTrackVariablesPair, iBinPair);
246 for (Int_t j=0; j<kTrackVariablesPair; j++) {
247 fHistNN->SetBinLimits(j, dBinsPair[j]);
248 fHistNN->SetVarTitle(j, axisTitlePair[j]);
250 AliInfo("Finished setting up the AliTHn");
253 fHistHBTbefore = new TH2D("fHistHBTbefore","before HBT cut",200,0,2,200,0,2.*TMath::Pi());
254 fHistHBTafter = new TH2D("fHistHBTafter","after HBT cut",200,0,2,200,0,2.*TMath::Pi());
255 fHistConversionbefore = new TH2D("fHistConversionbefore","before Conversion cut",200,0,2,200,0,2.*TMath::Pi());
256 fHistConversionafter = new TH2D("fHistConversionafter","after Conversion cut",200,0,2,200,0,2.*TMath::Pi());
257 fHistPsiMinusPhi = new TH2D("fHistPsiMinusPhi","",4,-0.5,3.5,100,0,2.*TMath::Pi());
259 TH1::AddDirectory(oldStatus);
263 //____________________________________________________________________//
264 void AliBalancePsi::CalculateBalance(Double_t gReactionPlane,
265 TObjArray *particles,
266 TObjArray *particlesMixed,
268 // Calculates the balance function
271 // Initialize histograms if not done yet
273 AliWarning("Histograms not yet initialized! --> Will be done now");
274 AliWarning("This works only in local mode --> Add 'gBalance->InitHistograms()' in your configBalanceFunction");
278 Double_t trackVariablesSingle[kTrackVariablesSingle];
279 Double_t trackVariablesPair[kTrackVariablesPair];
282 AliWarning("particles TObjArray is NULL pointer --> return");
286 // define end of particle loops
287 Int_t iMax = particles->GetEntriesFast();
290 jMax = particlesMixed->GetEntriesFast();
292 // Eta() is extremely time consuming, therefore cache it for the inner loop here:
293 TObjArray* particlesSecond = (particlesMixed) ? particlesMixed : particles;
295 TArrayF secondEta(jMax);
296 TArrayF secondPhi(jMax);
297 TArrayF secondPt(jMax);
299 for (Int_t i=0; i<jMax; i++){
300 secondEta[i] = ((AliVParticle*) particlesSecond->At(i))->Eta();
301 secondPhi[i] = ((AliVParticle*) particlesSecond->At(i))->Phi();
302 secondPt[i] = ((AliVParticle*) particlesSecond->At(i))->Pt();
306 for (Int_t i = 0; i < iMax; i++) {
307 AliVParticle* firstParticle = (AliVParticle*) particles->At(i);
310 Float_t firstEta = firstParticle->Eta();
311 Float_t firstPhi = firstParticle->Phi();
312 Float_t firstPt = firstParticle->Pt();
314 // Event plane (determine psi bin)
315 Double_t gPsiMinusPhi = 0.;
316 Double_t gPsiMinusPhiBin = -10.;
317 gPsiMinusPhi = TMath::Abs(firstPhi - gReactionPlane);
319 if((gPsiMinusPhi <= 7.5*TMath::DegToRad())||
320 ((172.5*TMath::DegToRad() <= gPsiMinusPhi)&&(gPsiMinusPhi <= 187.5*TMath::DegToRad())))
321 gPsiMinusPhiBin = 0.0;
323 else if(((37.5*TMath::DegToRad() <= gPsiMinusPhi)&&(gPsiMinusPhi <= 52.5*TMath::DegToRad()))||
324 ((127.5*TMath::DegToRad() <= gPsiMinusPhi)&&(gPsiMinusPhi <= 142.5*TMath::DegToRad()))||
325 ((217.5*TMath::DegToRad() <= gPsiMinusPhi)&&(gPsiMinusPhi <= 232.5*TMath::DegToRad()))||
326 ((307.5*TMath::DegToRad() <= gPsiMinusPhi)&&(gPsiMinusPhi <= 322.5*TMath::DegToRad())))
327 gPsiMinusPhiBin = 1.0;
329 else if(((82.5*TMath::DegToRad() <= gPsiMinusPhi)&&(gPsiMinusPhi <= 97.5*TMath::DegToRad()))||
330 ((262.5*TMath::DegToRad() <= gPsiMinusPhi)&&(gPsiMinusPhi <= 277.5*TMath::DegToRad())))
331 gPsiMinusPhiBin = 2.0;
334 gPsiMinusPhiBin = 3.0;
336 fHistPsiMinusPhi->Fill(gPsiMinusPhiBin,gPsiMinusPhi);
338 Short_t charge1 = (Short_t) firstParticle->Charge();
340 trackVariablesSingle[0] = gPsiMinusPhiBin;
341 trackVariablesSingle[1] = firstPt;
343 //fill single particle histograms
344 if(charge1 > 0) fHistP->Fill(trackVariablesSingle,0,1.);
345 else if(charge1 < 0) fHistN->Fill(trackVariablesSingle,0,1.);
347 // 2nd particle loop (only for j < i for non double counting in the same pT region)
348 // --> SAME pT region for trigger and assoc: NO double counting with this
349 // --> DIFF pT region for trigger and assoc: Missing assoc. particles with j > i to a trigger i
350 // --> can be handled afterwards by using assoc. as trigger as well ?!
351 for(Int_t j = 0; j < iMax; j++) {
352 if (particlesMixed && j == jMax-1 ) // if the mixed track number is smaller than the main event one
355 if(j == i) continue; // no auto correlations
357 AliVParticle* secondParticle = (AliVParticle*) particlesSecond->At(j);
358 Short_t charge2 = (Short_t) secondParticle->Charge();
360 trackVariablesPair[0] = gPsiMinusPhiBin;
361 trackVariablesPair[1] = firstEta - secondEta[j]; // delta eta
362 trackVariablesPair[2] = firstPhi - secondPhi[j]; // delta phi
363 //if (trackVariablesPair[2] > 180.) // delta phi between -180 and 180
364 //trackVariablesPair[2] -= 360.;
365 //if (trackVariablesPair[2] < - 180.)
366 //trackVariablesPair[2] += 360.;
367 if (trackVariablesPair[2] > TMath::Pi()) // delta phi between -pi and pi
368 trackVariablesPair[2] -= 2.*TMath::Pi();
369 if (trackVariablesPair[2] < - TMath::Pi())
370 trackVariablesPair[2] += 2.*TMath::Pi();
371 if (trackVariablesPair[2] < - TMath::Pi()/2.)
372 trackVariablesPair[2] += 2.*TMath::Pi();
374 trackVariablesPair[3] = firstPt; // pt trigger
375 trackVariablesPair[4] = secondPt[j]; // pt
376 // trackVariablesPair[5] = fCentrality; // centrality
379 if(fHBTCut && charge1 * charge2 > 0){
380 //if( dphi < 3 || deta < 0.01 ){ // VERSION 1
383 Double_t deta = firstEta - secondEta[j];
384 Double_t dphi = firstPhi - secondPhi[j];
385 // VERSION 2 (Taken from DPhiCorrelations)
386 // the variables & cuthave been developed by the HBT group
387 // see e.g. https://indico.cern.ch/materialDisplay.py?contribId=36&sessionId=6&materialId=slides&confId=142700
388 fHistHBTbefore->Fill(deta,dphi);
391 if (TMath::Abs(deta) < 0.02 * 2.5 * 3) //twoTrackEfficiencyCutValue = 0.02 [default for dphicorrelations]
394 //Float_t phi1rad = firstPhi*TMath::DegToRad();
395 //Float_t phi2rad = secondPhi[j]*TMath::DegToRad();
396 Float_t phi1rad = firstPhi;
397 Float_t phi2rad = secondPhi[j];
399 // check first boundaries to see if is worth to loop and find the minimum
400 Float_t dphistar1 = GetDPhiStar(phi1rad, firstPt, charge1, phi2rad, secondPt[j], charge2, 0.8, bSign);
401 Float_t dphistar2 = GetDPhiStar(phi1rad, firstPt, charge1, phi2rad, secondPt[j], charge2, 2.5, bSign);
403 const Float_t kLimit = 0.02 * 3;
405 Float_t dphistarminabs = 1e5;
406 Float_t dphistarmin = 1e5;
408 if (TMath::Abs(dphistar1) < kLimit || TMath::Abs(dphistar2) < kLimit || dphistar1 * dphistar2 < 0 ) {
409 for (Double_t rad=0.8; rad<2.51; rad+=0.01) {
410 Float_t dphistar = GetDPhiStar(phi1rad, firstPt, charge1, phi2rad, secondPt[j], charge2, rad, bSign);
411 Float_t dphistarabs = TMath::Abs(dphistar);
413 if (dphistarabs < dphistarminabs) {
414 dphistarmin = dphistar;
415 dphistarminabs = dphistarabs;
419 if (dphistarminabs < 0.02 && TMath::Abs(deta) < 0.02) {
420 //AliInfo(Form("HBT: Removed track pair %d %d with [[%f %f]] %f %f %f | %f %f %d %f %f %d %f", i, j, deta, dphi, dphistarminabs, dphistar1, dphistar2, phi1rad, pt1, charge1, phi2rad, pt2, charge2, bSign));
425 fHistHBTafter->Fill(deta,dphi);
430 if (charge1 * charge2 < 0) {
431 Double_t deta = firstEta - secondEta[j];
432 Double_t dphi = firstPhi - secondPhi[j];
433 fHistConversionbefore->Fill(deta,dphi);
435 Float_t m0 = 0.510e-3;
436 Float_t tantheta1 = 1e10;
439 //Float_t phi1rad = firstPhi*TMath::DegToRad();
440 //Float_t phi2rad = secondPhi[j]*TMath::DegToRad();
441 Float_t phi1rad = firstPhi;
442 Float_t phi2rad = secondPhi[j];
444 if (firstEta < -1e-10 || firstEta > 1e-10)
445 tantheta1 = 2 * TMath::Exp(-firstEta) / ( 1 - TMath::Exp(-2*firstEta));
447 Float_t tantheta2 = 1e10;
448 if (secondEta[j] < -1e-10 || secondEta[j] > 1e-10)
449 tantheta2 = 2 * TMath::Exp(-secondEta[j]) / ( 1 - TMath::Exp(-2*secondEta[j]));
451 Float_t e1squ = m0 * m0 + firstPt * firstPt * (1.0 + 1.0 / tantheta1 / tantheta1);
452 Float_t e2squ = m0 * m0 + secondPt[j] * secondPt[j] * (1.0 + 1.0 / tantheta2 / tantheta2);
454 Float_t masssqu = 2 * m0 * m0 + 2 * ( TMath::Sqrt(e1squ * e2squ) - ( firstPt * secondPt[j] * ( TMath::Cos(phi1rad - phi2rad) + 1.0 / tantheta1 / tantheta2 ) ) );
456 if (masssqu < 0.04*0.04){
457 //AliInfo(Form("Conversion: Removed track pair %d %d with [[%f %f] %f %f] %d %d <- %f %f %f %f %f %f ", i, j, deta, dphi, masssqu, charge1, charge2,eta1,eta2,phi1,phi2,pt1,pt2));
460 fHistConversionafter->Fill(deta,dphi);
464 if( charge1 > 0 && charge2 < 0) fHistPN->Fill(trackVariablesPair,0,1.);
465 else if( charge1 < 0 && charge2 > 0) fHistNP->Fill(trackVariablesPair,0,1.);
466 else if( charge1 > 0 && charge2 > 0) fHistPP->Fill(trackVariablesPair,0,1.);
467 else if( charge1 < 0 && charge2 < 0) fHistNN->Fill(trackVariablesPair,0,1.);
469 //AliWarning(Form("Wrong charge combination: charge1 = %d and charge2 = %d",charge,charge2));
472 }//end of 2nd particle loop
473 }//end of 1st particle loop
476 //____________________________________________________________________//
477 TH1D *AliBalancePsi::GetBalanceFunctionHistogram(Int_t iVariableSingle,
481 Double_t ptTriggerMin,
482 Double_t ptTriggerMax,
483 Double_t ptAssociatedMin,
484 Double_t ptAssociatedMax) {
485 //Returns the BF histogram, extracted from the 6 AliTHn objects
486 //(private members) of the AliBalancePsi class.
487 //iVariableSingle: 0(phi-Psi), 1(pt-trigger)
488 //iVariablePair: 0(phi-Psi) 1(Delta eta), 2(Delta phi), 3(pt-trigger), 4(pt-associated
490 fHistP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
491 fHistN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
492 fHistPN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
493 fHistNP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
494 fHistPP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
495 fHistNN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
498 if((ptTriggerMin != -1.)&&(ptTriggerMax != -1.)) {
499 fHistP->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
500 fHistN->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
501 fHistPN->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
502 fHistNP->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
503 fHistPP->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
504 fHistNN->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
508 if((ptAssociatedMin != -1.)&&(ptAssociatedMax != -1.)) {
509 fHistPN->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
510 fHistNP->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
511 fHistPP->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
512 fHistNN->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
515 //Printf("P:%lf - N:%lf - PN:%lf - NP:%lf - PP:%lf - NN:%lf",fHistP->GetEntries(0),fHistN->GetEntries(0),fHistPN->GetEntries(0),fHistNP->GetEntries(0),fHistPP->GetEntries(0),fHistNN->GetEntries(0));
517 // Project into the wanted space (1st: analysis step, 2nd: axis)
518 TH1D* hTemp1 = (TH1D*)fHistPN->Project(0,iVariablePair);
519 TH1D* hTemp2 = (TH1D*)fHistNP->Project(0,iVariablePair);
520 TH1D* hTemp3 = (TH1D*)fHistPP->Project(0,iVariablePair);
521 TH1D* hTemp4 = (TH1D*)fHistNN->Project(0,iVariablePair);
522 TH1D* hTemp5 = (TH1D*)fHistP->Project(0,iVariableSingle);
523 TH1D* hTemp6 = (TH1D*)fHistN->Project(0,iVariableSingle);
525 TH1D *gHistBalanceFunctionHistogram = 0x0;
526 if((hTemp1)&&(hTemp2)&&(hTemp3)&&(hTemp4)&&(hTemp5)&&(hTemp6)) {
527 gHistBalanceFunctionHistogram = (TH1D*)hTemp1->Clone();
528 gHistBalanceFunctionHistogram->Reset();
530 switch(iVariablePair) {
532 gHistBalanceFunctionHistogram->GetXaxis()->SetTitle("#Delta#eta");
533 gHistBalanceFunctionHistogram->GetYaxis()->SetTitle("B(#Delta#eta)");
536 gHistBalanceFunctionHistogram->GetXaxis()->SetTitle("#Delta#varphi (rad)");
537 gHistBalanceFunctionHistogram->GetYaxis()->SetTitle("B(#Delta#varphi)");
547 hTemp1->Add(hTemp3,-1.);
548 hTemp1->Scale(1./hTemp5->GetEntries());
549 hTemp2->Add(hTemp4,-1.);
550 hTemp2->Scale(1./hTemp6->GetEntries());
551 gHistBalanceFunctionHistogram->Add(hTemp1,hTemp2,1.,1.);
552 gHistBalanceFunctionHistogram->Scale(0.5);
555 return gHistBalanceFunctionHistogram;
558 //____________________________________________________________________//
559 TH1D *AliBalancePsi::GetBalanceFunctionHistogram2pMethod(Int_t iVariableSingle,
563 Double_t ptTriggerMin,
564 Double_t ptTriggerMax,
565 Double_t ptAssociatedMin,
566 Double_t ptAssociatedMax,
567 AliBalancePsi *bfMix) {
568 //Returns the BF histogram, extracted from the 6 AliTHn objects
569 //after dividing each correlation function by the Event Mixing one
570 //(private members) of the AliBalancePsi class.
571 //iVariableSingle: 0(phi-Psi), 1(pt-trigger)
572 //iVariablePair: 0(phi-Psi) 1(Delta eta), 2(Delta phi), 3(pt-trigger), 4(pt-associated
574 fHistP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
575 fHistN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
576 fHistPN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
577 fHistNP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
578 fHistPP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
579 fHistNN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
582 if((ptTriggerMin != -1.)&&(ptTriggerMax != -1.)) {
583 fHistP->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
584 fHistN->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
585 fHistPN->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
586 fHistNP->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
587 fHistPP->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
588 fHistNN->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
592 if((ptAssociatedMin != -1.)&&(ptAssociatedMax != -1.)) {
593 fHistPN->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
594 fHistNP->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
595 fHistPP->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
596 fHistNN->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
599 // ============================================================================================
600 // the same for event mixing
601 AliTHn *fHistPMix = bfMix->GetHistNp();
602 AliTHn *fHistNMix = bfMix->GetHistNn();
603 AliTHn *fHistPNMix = bfMix->GetHistNpn();
604 AliTHn *fHistNPMix = bfMix->GetHistNnp();
605 AliTHn *fHistPPMix = bfMix->GetHistNpp();
606 AliTHn *fHistNNMix = bfMix->GetHistNnn();
609 fHistPMix->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
610 fHistNMix->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
611 fHistPNMix->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
612 fHistNPMix->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
613 fHistPPMix->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
614 fHistNNMix->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
617 if((ptTriggerMin != -1.)&&(ptTriggerMax != -1.)) {
618 fHistPMix->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
619 fHistNMix->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
620 fHistPNMix->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
621 fHistNPMix->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
622 fHistPPMix->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
623 fHistNNMix->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
627 if((ptAssociatedMin != -1.)&&(ptAssociatedMax != -1.)) {
628 fHistPNMix->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
629 fHistNPMix->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
630 fHistPPMix->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
631 fHistNNMix->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
633 // ============================================================================================
635 //Printf("P:%lf - N:%lf - PN:%lf - NP:%lf - PP:%lf - NN:%lf",fHistP->GetEntries(0),fHistN->GetEntries(0),fHistPN->GetEntries(0),fHistNP->GetEntries(0),fHistPP->GetEntries(0),fHistNN->GetEntries(0));
637 // Project into the wanted space (1st: analysis step, 2nd: axis)
638 TH1D* hTemp1 = (TH1D*)fHistPN->Project(0,iVariablePair);
639 TH1D* hTemp2 = (TH1D*)fHistNP->Project(0,iVariablePair);
640 TH1D* hTemp3 = (TH1D*)fHistPP->Project(0,iVariablePair);
641 TH1D* hTemp4 = (TH1D*)fHistNN->Project(0,iVariablePair);
642 TH1D* hTemp5 = (TH1D*)fHistP->Project(0,iVariableSingle);
643 TH1D* hTemp6 = (TH1D*)fHistN->Project(0,iVariableSingle);
645 // ============================================================================================
646 // the same for event mixing
647 TH2D* hTemp1Mix = (TH2D*)fHistPNMix->Project(0,iVariablePair);
648 TH2D* hTemp2Mix = (TH2D*)fHistNPMix->Project(0,iVariablePair);
649 TH2D* hTemp3Mix = (TH2D*)fHistPPMix->Project(0,iVariablePair);
650 TH2D* hTemp4Mix = (TH2D*)fHistNNMix->Project(0,iVariablePair);
651 TH1D* hTemp5Mix = (TH1D*)fHistPMix->Project(0,iVariableSingle);
652 TH1D* hTemp6Mix = (TH1D*)fHistNMix->Project(0,iVariableSingle);
653 // ============================================================================================
655 TH1D *gHistBalanceFunctionHistogram = 0x0;
656 if((hTemp1)&&(hTemp2)&&(hTemp3)&&(hTemp4)&&(hTemp5)&&(hTemp6)) {
657 gHistBalanceFunctionHistogram = (TH1D*)hTemp1->Clone();
658 gHistBalanceFunctionHistogram->Reset();
660 switch(iVariablePair) {
662 gHistBalanceFunctionHistogram->GetXaxis()->SetTitle("#Delta#eta");
663 gHistBalanceFunctionHistogram->GetYaxis()->SetTitle("B(#Delta#eta)");
666 gHistBalanceFunctionHistogram->GetXaxis()->SetTitle("#Delta#varphi (rad)");
667 gHistBalanceFunctionHistogram->GetYaxis()->SetTitle("B(#Delta#varphi)");
682 hTemp1->Scale(1./hTemp5->GetEntries());
683 hTemp3->Scale(1./hTemp5->GetEntries());
684 hTemp2->Scale(1./hTemp6->GetEntries());
685 hTemp4->Scale(1./hTemp6->GetEntries());
686 hTemp1Mix->Scale(1./hTemp5Mix->GetEntries());
687 hTemp3Mix->Scale(1./hTemp5Mix->GetEntries());
688 hTemp2Mix->Scale(1./hTemp6Mix->GetEntries());
689 hTemp4Mix->Scale(1./hTemp6Mix->GetEntries());
691 hTemp1->Divide(hTemp1Mix);
692 hTemp2->Divide(hTemp2Mix);
693 hTemp3->Divide(hTemp3Mix);
694 hTemp4->Divide(hTemp4Mix);
696 hTemp1->Add(hTemp3,-1.);
697 hTemp2->Add(hTemp4,-1.);
699 gHistBalanceFunctionHistogram->Add(hTemp1,hTemp2,1.,1.);
700 gHistBalanceFunctionHistogram->Scale(0.5);
703 return gHistBalanceFunctionHistogram;
706 //____________________________________________________________________//
707 TH2D *AliBalancePsi::GetBalanceFunctionDeltaEtaDeltaPhi(Double_t psiMin,
709 Double_t ptTriggerMin,
710 Double_t ptTriggerMax,
711 Double_t ptAssociatedMin,
712 Double_t ptAssociatedMax) {
713 //Returns the BF histogram in Delta eta vs Delta phi,
714 //extracted from the 6 AliTHn objects
715 //(private members) of the AliBalancePsi class.
716 //iVariableSingle: 0(phi-Psi), 1(pt-trigger)
717 //iVariablePair: 0(phi-Psi) 1(Delta eta), 2(Delta phi), 3(pt-trigger), 4(pt-associated
718 TString histName = "gHistBalanceFunctionHistogram2D";
721 fHistP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
722 fHistN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
723 fHistPN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
724 fHistNP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
725 fHistPP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
726 fHistNN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
729 if((ptTriggerMin != -1.)&&(ptTriggerMax != -1.)) {
730 fHistP->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
731 fHistN->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
732 fHistPN->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
733 fHistNP->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
734 fHistPP->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
735 fHistNN->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
739 if((ptAssociatedMin != -1.)&&(ptAssociatedMax != -1.)) {
740 fHistPN->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
741 fHistNP->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
742 fHistPP->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
743 fHistNN->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
746 //AliInfo(Form("P:%lf - N:%lf - PN:%lf - NP:%lf - PP:%lf - NN:%lf",fHistP->GetEntries(0),fHistN->GetEntries(0),fHistPN->GetEntries(0),fHistNP->GetEntries(0),fHistPP->GetEntries(0),fHistNN->GetEntries(0)));
748 // Project into the wanted space (1st: analysis step, 2nd: axis)
749 TH2D* hTemp1 = (TH2D*)fHistPN->Project(0,1,2);
750 TH2D* hTemp2 = (TH2D*)fHistNP->Project(0,1,2);
751 TH2D* hTemp3 = (TH2D*)fHistPP->Project(0,1,2);
752 TH2D* hTemp4 = (TH2D*)fHistNN->Project(0,1,2);
753 TH1D* hTemp5 = (TH1D*)fHistP->Project(0,1);
754 TH1D* hTemp6 = (TH1D*)fHistN->Project(0,1);
756 TH2D *gHistBalanceFunctionHistogram = 0x0;
757 if((hTemp1)&&(hTemp2)&&(hTemp3)&&(hTemp4)&&(hTemp5)&&(hTemp6)) {
758 gHistBalanceFunctionHistogram = (TH2D*)hTemp1->Clone();
759 gHistBalanceFunctionHistogram->Reset();
760 gHistBalanceFunctionHistogram->GetXaxis()->SetTitle("#Delta#eta");
761 gHistBalanceFunctionHistogram->GetYaxis()->SetTitle("#Delta#varphi (rad)");
762 gHistBalanceFunctionHistogram->GetZaxis()->SetTitle("B(#Delta#eta,#Delta#varphi)");
768 hTemp1->Add(hTemp3,-1.);
769 hTemp1->Scale(1./hTemp5->GetEntries());
770 hTemp2->Add(hTemp4,-1.);
771 hTemp2->Scale(1./hTemp6->GetEntries());
772 gHistBalanceFunctionHistogram->Add(hTemp1,hTemp2,1.,1.);
773 gHistBalanceFunctionHistogram->Scale(0.5);
776 return gHistBalanceFunctionHistogram;
779 //____________________________________________________________________//
780 TH2D *AliBalancePsi::GetBalanceFunctionDeltaEtaDeltaPhi2pMethod(Double_t psiMin,
782 Double_t ptTriggerMin,
783 Double_t ptTriggerMax,
784 Double_t ptAssociatedMin,
785 Double_t ptAssociatedMax,
786 AliBalancePsi *bfMix) {
787 //Returns the BF histogram in Delta eta vs Delta phi,
788 //after dividing each correlation function by the Event Mixing one
789 //extracted from the 6 AliTHn objects
790 //(private members) of the AliBalancePsi class.
791 //iVariableSingle: 0(phi-Psi), 1(pt-trigger)
792 //iVariablePair: 0(phi-Psi) 1(Delta eta), 2(Delta phi), 3(pt-trigger), 4(pt-associated
793 TString histName = "gHistBalanceFunctionHistogram2D";
796 AliError("balance function object for event mixing not available");
801 fHistP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
802 fHistN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
803 fHistPN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
804 fHistNP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
805 fHistPP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
806 fHistNN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
809 if((ptTriggerMin != -1.)&&(ptTriggerMax != -1.)) {
810 fHistP->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
811 fHistN->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
812 fHistPN->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
813 fHistNP->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
814 fHistPP->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
815 fHistNN->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
819 if((ptAssociatedMin != -1.)&&(ptAssociatedMax != -1.)) {
820 fHistPN->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
821 fHistNP->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
822 fHistPP->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
823 fHistNN->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
827 // ============================================================================================
828 // the same for event mixing
829 AliTHn *fHistPMix = bfMix->GetHistNp();
830 AliTHn *fHistNMix = bfMix->GetHistNn();
831 AliTHn *fHistPNMix = bfMix->GetHistNpn();
832 AliTHn *fHistNPMix = bfMix->GetHistNnp();
833 AliTHn *fHistPPMix = bfMix->GetHistNpp();
834 AliTHn *fHistNNMix = bfMix->GetHistNnn();
837 fHistPMix->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
838 fHistNMix->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
839 fHistPNMix->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
840 fHistNPMix->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
841 fHistPPMix->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
842 fHistNNMix->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
845 if((ptTriggerMin != -1.)&&(ptTriggerMax != -1.)) {
846 fHistPMix->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
847 fHistNMix->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
848 fHistPNMix->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
849 fHistNPMix->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
850 fHistPPMix->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
851 fHistNNMix->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
855 if((ptAssociatedMin != -1.)&&(ptAssociatedMax != -1.)) {
856 fHistPNMix->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
857 fHistNPMix->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
858 fHistPPMix->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
859 fHistNNMix->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
861 // ============================================================================================
864 //AliInfo(Form("P:%lf - N:%lf - PN:%lf - NP:%lf - PP:%lf - NN:%lf",fHistP->GetEntries(0),fHistN->GetEntries(0),fHistPN->GetEntries(0),fHistNP->GetEntries(0),fHistPP->GetEntries(0),fHistNN->GetEntries(0)));
866 // Project into the wanted space (1st: analysis step, 2nd: axis)
867 TH2D* hTemp1 = (TH2D*)fHistPN->Project(0,1,2);
868 TH2D* hTemp2 = (TH2D*)fHistNP->Project(0,1,2);
869 TH2D* hTemp3 = (TH2D*)fHistPP->Project(0,1,2);
870 TH2D* hTemp4 = (TH2D*)fHistNN->Project(0,1,2);
871 TH1D* hTemp5 = (TH1D*)fHistP->Project(0,1);
872 TH1D* hTemp6 = (TH1D*)fHistN->Project(0,1);
874 // ============================================================================================
875 // the same for event mixing
876 TH2D* hTemp1Mix = (TH2D*)fHistPNMix->Project(0,1,2);
877 TH2D* hTemp2Mix = (TH2D*)fHistNPMix->Project(0,1,2);
878 TH2D* hTemp3Mix = (TH2D*)fHistPPMix->Project(0,1,2);
879 TH2D* hTemp4Mix = (TH2D*)fHistNNMix->Project(0,1,2);
880 TH1D* hTemp5Mix = (TH1D*)fHistPMix->Project(0,1);
881 TH1D* hTemp6Mix = (TH1D*)fHistNMix->Project(0,1);
882 // ============================================================================================
884 TH2D *gHistBalanceFunctionHistogram = 0x0;
885 if((hTemp1)&&(hTemp2)&&(hTemp3)&&(hTemp4)&&(hTemp5)&&(hTemp6)) {
886 gHistBalanceFunctionHistogram = (TH2D*)hTemp1->Clone();
887 gHistBalanceFunctionHistogram->Reset();
888 gHistBalanceFunctionHistogram->GetXaxis()->SetTitle("#Delta#eta");
889 gHistBalanceFunctionHistogram->GetYaxis()->SetTitle("#Delta#varphi (rad)");
890 gHistBalanceFunctionHistogram->GetZaxis()->SetTitle("B(#Delta#eta,#Delta#varphi)");
901 hTemp1->Scale(1./hTemp5->GetEntries());
902 hTemp3->Scale(1./hTemp5->GetEntries());
903 hTemp2->Scale(1./hTemp6->GetEntries());
904 hTemp4->Scale(1./hTemp6->GetEntries());
905 hTemp1Mix->Scale(1./hTemp5Mix->GetEntries());
906 hTemp3Mix->Scale(1./hTemp5Mix->GetEntries());
907 hTemp2Mix->Scale(1./hTemp6Mix->GetEntries());
908 hTemp4Mix->Scale(1./hTemp6Mix->GetEntries());
910 hTemp1->Divide(hTemp1Mix);
911 hTemp2->Divide(hTemp2Mix);
912 hTemp3->Divide(hTemp3Mix);
913 hTemp4->Divide(hTemp4Mix);
915 hTemp1->Add(hTemp3,-1.);
916 hTemp2->Add(hTemp4,-1.);
918 gHistBalanceFunctionHistogram->Add(hTemp1,hTemp2,1.,1.);
919 gHistBalanceFunctionHistogram->Scale(0.5);
922 return gHistBalanceFunctionHistogram;
926 //____________________________________________________________________//
927 TH2D *AliBalancePsi::GetCorrelationFunctionPN(Double_t psiMin,
929 Double_t ptTriggerMin,
930 Double_t ptTriggerMax,
931 Double_t ptAssociatedMin,
932 Double_t ptAssociatedMax) {
933 //Returns the 2D correlation function for (+-) pairs
935 fHistP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
936 fHistPN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
939 if((ptTriggerMin != -1.)&&(ptTriggerMax != -1.)) {
940 fHistP->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
941 fHistPN->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
945 if((ptAssociatedMin != -1.)&&(ptAssociatedMax != -1.))
946 fHistPN->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
948 //fHistP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(-0.5,2.5);
949 //fHistPN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(-0.5,2.5);
951 //TH2D *gHistTest = dynamic_cast<TH2D *>(fHistP->Project(0,0,1));
952 //TCanvas *c1 = new TCanvas("c1","");
955 //AliError("Projection of fHistP = NULL");
959 //gHistTest->DrawCopy("colz");
962 //0:step, 1: Delta eta, 2: Delta phi
963 TH2D *gHist = dynamic_cast<TH2D *>(fHistPN->Project(0,1,2));
965 AliError("Projection of fHistPN = NULL");
969 //AliInfo(Form("Entries (test): %lf",(Double_t)(gHistTest->GetEntries())));
970 //AliInfo(Form("Entries (1D): %lf",(Double_t)(fHistP->Project(0,1)->GetEntries())));
971 //AliInfo(Form("Entries (2D): %lf",(Double_t)(fHistPN->Project(0,1,2)->GetEntries())));
973 //TCanvas *c2 = new TCanvas("c2","");
975 //fHistPN->Project(0,1,2)->DrawCopy("colz");
977 if((Double_t)(fHistP->Project(0,1)->GetEntries())!=0)
978 gHist->Scale(1./(Double_t)(fHistP->Project(0,1)->GetEntries()));
983 //____________________________________________________________________//
984 TH2D *AliBalancePsi::GetCorrelationFunctionNP(Double_t psiMin,
986 Double_t ptTriggerMin,
987 Double_t ptTriggerMax,
988 Double_t ptAssociatedMin,
989 Double_t ptAssociatedMax) {
990 //Returns the 2D correlation function for (+-) pairs
992 fHistN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
993 fHistNP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
996 if((ptTriggerMin != -1.)&&(ptTriggerMax != -1.)) {
997 fHistN->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
998 fHistNP->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1002 if((ptAssociatedMin != -1.)&&(ptAssociatedMax != -1.))
1003 fHistNP->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
1005 //0:step, 1: Delta eta, 2: Delta phi
1006 TH2D *gHist = dynamic_cast<TH2D *>(fHistNP->Project(0,1,2));
1008 AliError("Projection of fHistPN = NULL");
1012 //Printf("Entries (1D): %lf",(Double_t)(fHistN->Project(0,2)->GetEntries()));
1013 //Printf("Entries (2D): %lf",(Double_t)(fHistNP->Project(0,2,3)->GetEntries()));
1014 if((Double_t)(fHistN->Project(0,1)->GetEntries())!=0)
1015 gHist->Scale(1./(Double_t)(fHistN->Project(0,1)->GetEntries()));
1020 //____________________________________________________________________//
1021 TH2D *AliBalancePsi::GetCorrelationFunctionPP(Double_t psiMin,
1023 Double_t ptTriggerMin,
1024 Double_t ptTriggerMax,
1025 Double_t ptAssociatedMin,
1026 Double_t ptAssociatedMax) {
1027 //Returns the 2D correlation function for (+-) pairs
1029 fHistP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
1030 fHistPP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
1033 if((ptTriggerMin != -1.)&&(ptTriggerMax != -1.)) {
1034 fHistP->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1035 fHistPP->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1039 if((ptAssociatedMin != -1.)&&(ptAssociatedMax != -1.))
1040 fHistPP->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
1042 //0:step, 1: Delta eta, 2: Delta phi
1043 TH2D *gHist = dynamic_cast<TH2D *>(fHistPP->Project(0,1,2));
1045 AliError("Projection of fHistPN = NULL");
1049 //Printf("Entries (1D): %lf",(Double_t)(fHistP->Project(0,2)->GetEntries()));
1050 //Printf("Entries (2D): %lf",(Double_t)(fHistPP->Project(0,2,3)->GetEntries()));
1051 if((Double_t)(fHistP->Project(0,1)->GetEntries())!=0)
1052 gHist->Scale(1./(Double_t)(fHistP->Project(0,1)->GetEntries()));
1057 //____________________________________________________________________//
1058 TH2D *AliBalancePsi::GetCorrelationFunctionNN(Double_t psiMin,
1060 Double_t ptTriggerMin,
1061 Double_t ptTriggerMax,
1062 Double_t ptAssociatedMin,
1063 Double_t ptAssociatedMax) {
1064 //Returns the 2D correlation function for (+-) pairs
1066 fHistN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
1067 fHistNN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
1070 if((ptTriggerMin != -1.)&&(ptTriggerMax != -1.)) {
1071 fHistN->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1072 fHistNN->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1076 if((ptAssociatedMin != -1.)&&(ptAssociatedMax != -1.))
1077 fHistNN->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
1079 //0:step, 1: Delta eta, 2: Delta phi
1080 TH2D *gHist = dynamic_cast<TH2D *>(fHistNN->Project(0,1,2));
1082 AliError("Projection of fHistPN = NULL");
1086 //Printf("Entries (1D): %lf",(Double_t)(fHistN->Project(0,2)->GetEntries()));
1087 //Printf("Entries (2D): %lf",(Double_t)(fHistNN->Project(0,2,3)->GetEntries()));
1088 if((Double_t)(fHistN->Project(0,1)->GetEntries())!=0)
1089 gHist->Scale(1./(Double_t)(fHistN->Project(0,1)->GetEntries()));
1094 //____________________________________________________________________//
1095 TH2D *AliBalancePsi::GetCorrelationFunctionChargeIndependent(Double_t psiMin,
1097 Double_t ptTriggerMin,
1098 Double_t ptTriggerMax,
1099 Double_t ptAssociatedMin,
1100 Double_t ptAssociatedMax) {
1101 //Returns the 2D correlation function for the sum of all charge combination pairs
1103 fHistN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
1104 fHistP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
1105 fHistNN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
1106 fHistPP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
1107 fHistNP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
1108 fHistPN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
1111 if((ptTriggerMin != -1.)&&(ptTriggerMax != -1.)) {
1112 fHistN->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1113 fHistP->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1114 fHistNN->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1115 fHistPP->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1116 fHistNP->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1117 fHistPN->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1121 if((ptAssociatedMin != -1.)&&(ptAssociatedMax != -1.))
1122 fHistNN->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
1123 fHistPP->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
1124 fHistNP->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
1125 fHistPN->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
1127 //0:step, 1: Delta eta, 2: Delta phi
1128 TH2D *gHistNN = dynamic_cast<TH2D *>(fHistNN->Project(0,1,2));
1130 AliError("Projection of fHistNN = NULL");
1133 TH2D *gHistPP = dynamic_cast<TH2D *>(fHistPP->Project(0,1,2));
1135 AliError("Projection of fHistPP = NULL");
1138 TH2D *gHistNP = dynamic_cast<TH2D *>(fHistNP->Project(0,1,2));
1140 AliError("Projection of fHistNP = NULL");
1143 TH2D *gHistPN = dynamic_cast<TH2D *>(fHistPN->Project(0,1,2));
1145 AliError("Projection of fHistPN = NULL");
1149 // sum all 2 particle histograms
1150 gHistNN->Add(gHistPP);
1151 gHistNN->Add(gHistNP);
1152 gHistNN->Add(gHistPN);
1154 // divide by sum of + and - triggers
1155 if((Double_t)(fHistN->Project(0,1)->GetEntries())!=0 && (Double_t)(fHistP->Project(0,1)->GetEntries())!=0)
1156 gHistNN->Scale(1./(Double_t)(fHistN->Project(0,1)->GetEntries() + fHistN->Project(0,1)->GetEntries()));
1161 //____________________________________________________________________//
1162 Float_t AliBalancePsi::GetDPhiStar(Float_t phi1, Float_t pt1, Float_t charge1, Float_t phi2, Float_t pt2, Float_t charge2, Float_t radius, Float_t bSign) {
1164 // calculates dphistar
1166 Float_t dphistar = phi1 - phi2 - charge1 * bSign * TMath::ASin(0.075 * radius / pt1) + charge2 * bSign * TMath::ASin(0.075 * radius / pt2);
1168 static const Double_t kPi = TMath::Pi();
1171 // if (dphistar > 2 * kPi)
1172 // dphistar -= 2 * kPi;
1173 // if (dphistar < -2 * kPi)
1174 // dphistar += 2 * kPi;
1177 dphistar = kPi * 2 - dphistar;
1178 if (dphistar < -kPi)
1179 dphistar = -kPi * 2 - dphistar;
1180 if (dphistar > kPi) // might look funny but is needed
1181 dphistar = kPi * 2 - dphistar;