]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGCF/EBYE/BalanceFunctions/AliBalancePsi.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
fix for AddTaskChaoticity function calls
[u/mrichter/AliRoot.git] / PWGCF / EBYE / BalanceFunctions / AliBalancePsi.cxx
CommitLineData
0879e280 1/**************************************************************************
2 * Author: Panos Christakoglou. *
3 * Contributors are mentioned in the code where appropriate. *
4 * *
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 **************************************************************************/
13
14/* $Id: AliBalancePsi.cxx 54125 2012-01-24 21:07:41Z miweber $ */
15
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//-----------------------------------------------------------------
21
22
23//ROOT
24#include <Riostream.h>
621329de 25#include <TCanvas.h>
0879e280 26#include <TMath.h>
27#include <TAxis.h>
28#include <TH2D.h>
29#include <TH3D.h>
30#include <TLorentzVector.h>
31#include <TObjArray.h>
32#include <TGraphErrors.h>
33#include <TString.h>
34
35#include "AliVParticle.h"
36#include "AliMCParticle.h"
37#include "AliESDtrack.h"
38#include "AliAODTrack.h"
9afe3098 39#include "AliTHn.h"
0879e280 40
41#include "AliBalancePsi.h"
42
43ClassImp(AliBalancePsi)
44
45//____________________________________________________________________//
46AliBalancePsi::AliBalancePsi() :
47 TObject(),
9fd4b54e 48 fShuffle(kFALSE),
0879e280 49 fAnalysisLevel("ESD"),
50 fAnalyzedEvents(0) ,
51 fCentralityId(0) ,
52 fCentStart(0.),
53 fCentStop(0.),
9afe3098 54 fHistP(0),
55 fHistN(0),
56 fHistPN(0),
57 fHistNP(0),
58 fHistPP(0),
59 fHistNN(0),
73609a48 60 fHistHBTbefore(0),
61 fHistHBTafter(0),
62 fHistConversionbefore(0),
63 fHistConversionafter(0),
c8ad9635 64 fHistPsiMinusPhi(0),
73609a48 65 fPsiInterval(15.),
66 fHBTCut(kFALSE),
67 fConversionCut(kFALSE) {
0879e280 68 // Default constructor
0879e280 69}
70
0879e280 71//____________________________________________________________________//
72AliBalancePsi::AliBalancePsi(const AliBalancePsi& balance):
9fd4b54e 73 TObject(balance), fShuffle(balance.fShuffle),
0879e280 74 fAnalysisLevel(balance.fAnalysisLevel),
75 fAnalyzedEvents(balance.fAnalyzedEvents),
76 fCentralityId(balance.fCentralityId),
77 fCentStart(balance.fCentStart),
78 fCentStop(balance.fCentStop),
9afe3098 79 fHistP(balance.fHistP),
80 fHistN(balance.fHistN),
81 fHistPN(balance.fHistPN),
82 fHistNP(balance.fHistNP),
83 fHistPP(balance.fHistPP),
84 fHistNN(balance.fHistNN),
73609a48 85 fHistHBTbefore(balance.fHistHBTbefore),
86 fHistHBTafter(balance.fHistHBTafter),
87 fHistConversionbefore(balance.fHistConversionbefore),
88 fHistConversionafter(balance.fHistConversionafter),
c8ad9635 89 fHistPsiMinusPhi(balance.fHistPsiMinusPhi),
73609a48 90 fPsiInterval(balance.fPsiInterval),
91 fHBTCut(balance.fHBTCut),
92 fConversionCut(balance.fConversionCut) {
0879e280 93 //copy constructor
0879e280 94}
95
96//____________________________________________________________________//
97AliBalancePsi::~AliBalancePsi() {
98 // Destructor
9afe3098 99 delete fHistP;
100 delete fHistN;
101 delete fHistPN;
102 delete fHistNP;
103 delete fHistPP;
104 delete fHistNN;
73609a48 105
106 delete fHistHBTbefore;
107 delete fHistHBTafter;
108 delete fHistConversionbefore;
109 delete fHistConversionafter;
c8ad9635 110 delete fHistPsiMinusPhi;
0879e280 111}
112
113//____________________________________________________________________//
9afe3098 114void AliBalancePsi::InitHistograms() {
115 // single particle histograms
211b716d 116
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);
121
9afe3098 122 Int_t anaSteps = 1; // analysis steps
9fd4b54e 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
9afe3098 126
127 // two particle histograms
9fd4b54e 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
9afe3098 131
132 //centrality
621329de 133 /*const Int_t kNCentralityBins = 9;
9afe3098 134 Double_t centralityBins[kNCentralityBins+1] = {0.,5.,10.,20.,30.,40.,50.,60.,70.,80.};
9afe3098 135 iBinSingle[0] = kNCentralityBins;
136 dBinsSingle[0] = centralityBins;
137 axisTitleSingle[0] = "Centrality percentile [%]";
138 iBinPair[0] = kNCentralityBins;
139 dBinsPair[0] = centralityBins;
621329de 140 axisTitlePair[0] = "Centrality percentile [%]"; */
9afe3098 141
c8ad9635 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)
6acdbcb2 143 const Int_t kNPsi2Bins = 4;
144 Double_t psi2Bins[kNPsi2Bins+1] = {-0.5,0.5,1.5,2.5,3.5};
621329de 145 iBinSingle[0] = kNPsi2Bins;
146 dBinsSingle[0] = psi2Bins;
c8ad9635 147 axisTitleSingle[0] = "#varphi - #Psi_{2} (a.u.)";
621329de 148 iBinPair[0] = kNPsi2Bins;
149 dBinsPair[0] = psi2Bins;
c8ad9635 150 axisTitlePair[0] = "#varphi - #Psi_{2} (a.u.)";
9afe3098 151
621329de 152 // delta eta
153 const Int_t kNDeltaEtaBins = 80;
9afe3098 154 Double_t deltaEtaBins[kNDeltaEtaBins+1];
155 for(Int_t i = 0; i < kNDeltaEtaBins+1; i++)
621329de 156 deltaEtaBins[i] = -2.0 + i * 0.05;
157 iBinPair[1] = kNDeltaEtaBins;
158 dBinsPair[1] = deltaEtaBins;
c8ad9635 159 axisTitlePair[1] = "#Delta#eta";
621329de 160
161 // delta phi
9afe3098 162 const Int_t kNDeltaPhiBins = 72;
163 Double_t deltaPhiBins[kNDeltaPhiBins+1];
164 for(Int_t i = 0; i < kNDeltaPhiBins+1; i++){
c8ad9635 165 //deltaPhiBins[i] = -180.0 + i * 5.;
166 deltaPhiBins[i] = -TMath::Pi()/2. + i * 5.*TMath::Pi()/180.;
9afe3098 167 }
621329de 168 iBinPair[2] = kNDeltaPhiBins;
169 dBinsPair[2] = deltaPhiBins;
c8ad9635 170 axisTitlePair[2] = "#Delta#varphi (rad)";
9afe3098 171
a38fd7f3 172 // pt(trigger-associated)
8795e993 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;
177 //}
621329de 178 iBinSingle[1] = kNPtBins;
179 dBinsSingle[1] = ptBins;
c8ad9635 180 axisTitleSingle[1] = "p_{T,trig.} (GeV/c)";
621329de 181
182 iBinPair[3] = kNPtBins;
183 dBinsPair[3] = ptBins;
c8ad9635 184 axisTitlePair[3] = "p_{T,trig.} (GeV/c)";
a38fd7f3 185
186 iBinPair[4] = kNPtBins;
187 dBinsPair[4] = ptBins;
c8ad9635 188 axisTitlePair[4] = "p_{T,assoc.} (GeV/c)";
9afe3098 189
190 TString histName;
191 //+ triggered particles
192 histName = "fHistP";
9fd4b54e 193 if(fShuffle) histName.Append("_shuffle");
9afe3098 194 if(fCentralityId) histName += fCentralityId.Data();
9fd4b54e 195 fHistP = new AliTHn(histName.Data(),histName.Data(),anaSteps,kTrackVariablesSingle,iBinSingle);
196 for (Int_t j=0; j<kTrackVariablesSingle; j++) {
9afe3098 197 fHistP->SetBinLimits(j, dBinsSingle[j]);
198 fHistP->SetVarTitle(j, axisTitleSingle[j]);
0879e280 199 }
9afe3098 200
201 //- triggered particles
202 histName = "fHistN";
9fd4b54e 203 if(fShuffle) histName.Append("_shuffle");
9afe3098 204 if(fCentralityId) histName += fCentralityId.Data();
9fd4b54e 205 fHistN = new AliTHn(histName.Data(),histName.Data(),anaSteps,kTrackVariablesSingle,iBinSingle);
206 for (Int_t j=0; j<kTrackVariablesSingle; j++) {
9afe3098 207 fHistN->SetBinLimits(j, dBinsSingle[j]);
208 fHistN->SetVarTitle(j, axisTitleSingle[j]);
0879e280 209 }
9afe3098 210
211 //+- pairs
212 histName = "fHistPN";
9fd4b54e 213 if(fShuffle) histName.Append("_shuffle");
9afe3098 214 if(fCentralityId) histName += fCentralityId.Data();
9fd4b54e 215 fHistPN = new AliTHn(histName.Data(),histName.Data(),anaSteps, kTrackVariablesPair, iBinPair);
216 for (Int_t j=0; j<kTrackVariablesPair; j++) {
9afe3098 217 fHistPN->SetBinLimits(j, dBinsPair[j]);
218 fHistPN->SetVarTitle(j, axisTitlePair[j]);
0879e280 219 }
0879e280 220
9afe3098 221 //-+ pairs
222 histName = "fHistNP";
9fd4b54e 223 if(fShuffle) histName.Append("_shuffle");
9afe3098 224 if(fCentralityId) histName += fCentralityId.Data();
9fd4b54e 225 fHistNP = new AliTHn(histName.Data(),histName.Data(),anaSteps, kTrackVariablesPair, iBinPair);
226 for (Int_t j=0; j<kTrackVariablesPair; j++) {
9afe3098 227 fHistNP->SetBinLimits(j, dBinsPair[j]);
228 fHistNP->SetVarTitle(j, axisTitlePair[j]);
0879e280 229 }
0879e280 230
9afe3098 231 //++ pairs
232 histName = "fHistPP";
9fd4b54e 233 if(fShuffle) histName.Append("_shuffle");
9afe3098 234 if(fCentralityId) histName += fCentralityId.Data();
9fd4b54e 235 fHistPP = new AliTHn(histName.Data(),histName.Data(),anaSteps, kTrackVariablesPair, iBinPair);
236 for (Int_t j=0; j<kTrackVariablesPair; j++) {
9afe3098 237 fHistPP->SetBinLimits(j, dBinsPair[j]);
238 fHistPP->SetVarTitle(j, axisTitlePair[j]);
239 }
240
241 //-- pairs
242 histName = "fHistNN";
9fd4b54e 243 if(fShuffle) histName.Append("_shuffle");
9afe3098 244 if(fCentralityId) histName += fCentralityId.Data();
9fd4b54e 245 fHistNN = new AliTHn(histName.Data(),histName.Data(),anaSteps, kTrackVariablesPair, iBinPair);
246 for (Int_t j=0; j<kTrackVariablesPair; j++) {
9afe3098 247 fHistNN->SetBinLimits(j, dBinsPair[j]);
248 fHistNN->SetVarTitle(j, axisTitlePair[j]);
0879e280 249 }
f06d59b3 250 AliInfo("Finished setting up the AliTHn");
73609a48 251
252 // QA histograms
c8ad9635 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());
211b716d 258
259 TH1::AddDirectory(oldStatus);
260
0879e280 261}
262
263//____________________________________________________________________//
f06d59b3 264void AliBalancePsi::CalculateBalance(Double_t gReactionPlane,
265 TObjArray *particles,
73609a48 266 TObjArray *particlesMixed,
267 Float_t bSign) {
0879e280 268 // Calculates the balance function
269 fAnalyzedEvents++;
0879e280 270
271 // Initialize histograms if not done yet
9afe3098 272 if(!fHistPN){
0879e280 273 AliWarning("Histograms not yet initialized! --> Will be done now");
274 AliWarning("This works only in local mode --> Add 'gBalance->InitHistograms()' in your configBalanceFunction");
275 InitHistograms();
276 }
277
9fd4b54e 278 Double_t trackVariablesSingle[kTrackVariablesSingle];
279 Double_t trackVariablesPair[kTrackVariablesPair];
f06d59b3 280
281 if (!particles){
282 AliWarning("particles TObjArray is NULL pointer --> return");
283 return;
284 }
9afe3098 285
f06d59b3 286 // define end of particle loops
287 Int_t iMax = particles->GetEntriesFast();
288 Int_t jMax = iMax;
289 if (particlesMixed)
290 jMax = particlesMixed->GetEntriesFast();
291
292 // Eta() is extremely time consuming, therefore cache it for the inner loop here:
293 TObjArray* particlesSecond = (particlesMixed) ? particlesMixed : particles;
294
295 TArrayF secondEta(jMax);
296 TArrayF secondPhi(jMax);
297 TArrayF secondPt(jMax);
298
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();
303 }
304
305 // 1st particle loop
73609a48 306 for (Int_t i = 0; i < iMax; i++) {
6acdbcb2 307 AliVParticle* firstParticle = (AliVParticle*) particles->At(i);
308
309 // some optimization
310 Float_t firstEta = firstParticle->Eta();
311 Float_t firstPhi = firstParticle->Phi();
312 Float_t firstPt = firstParticle->Pt();
313
314 // Event plane (determine psi bin)
315 Double_t gPsiMinusPhi = 0.;
316 Double_t gPsiMinusPhiBin = -10.;
317 gPsiMinusPhi = TMath::Abs(firstPhi - gReactionPlane);
318 //in-plane
c8ad9635 319 if((gPsiMinusPhi <= 7.5*TMath::DegToRad())||
320 ((172.5*TMath::DegToRad() <= gPsiMinusPhi)&&(gPsiMinusPhi <= 187.5*TMath::DegToRad())))
6acdbcb2 321 gPsiMinusPhiBin = 0.0;
322 //intermediate
c8ad9635 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())))
6acdbcb2 327 gPsiMinusPhiBin = 1.0;
328 //out of plane
c8ad9635 329 else if(((82.5*TMath::DegToRad() <= gPsiMinusPhi)&&(gPsiMinusPhi <= 97.5*TMath::DegToRad()))||
330 ((262.5*TMath::DegToRad() <= gPsiMinusPhi)&&(gPsiMinusPhi <= 277.5*TMath::DegToRad())))
6acdbcb2 331 gPsiMinusPhiBin = 2.0;
332 //everything else
333 else
334 gPsiMinusPhiBin = 3.0;
335
c8ad9635 336 fHistPsiMinusPhi->Fill(gPsiMinusPhiBin,gPsiMinusPhi);
73609a48 337
338 Short_t charge1 = (Short_t) firstParticle->Charge();
6acdbcb2 339
340 trackVariablesSingle[0] = gPsiMinusPhiBin;
341 trackVariablesSingle[1] = firstPt;
342
343 //fill single particle histograms
73609a48 344 if(charge1 > 0) fHistP->Fill(trackVariablesSingle,0,1.);
345 else if(charge1 < 0) fHistN->Fill(trackVariablesSingle,0,1.);
6acdbcb2 346
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 ?!
a86235ae 351 for(Int_t j = 0; j < iMax; j++) {
a86235ae 352 if (particlesMixed && j == jMax-1 ) // if the mixed track number is smaller than the main event one
6acdbcb2 353 break;
a86235ae 354
355 if(j == i) continue; // no auto correlations
9afe3098 356
6acdbcb2 357 AliVParticle* secondParticle = (AliVParticle*) particlesSecond->At(j);
6acdbcb2 358 Short_t charge2 = (Short_t) secondParticle->Charge();
9afe3098 359
6acdbcb2 360 trackVariablesPair[0] = gPsiMinusPhiBin;
361 trackVariablesPair[1] = firstEta - secondEta[j]; // delta eta
362 trackVariablesPair[2] = firstPhi - secondPhi[j]; // delta phi
c8ad9635 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();
9afe3098 373
6acdbcb2 374 trackVariablesPair[3] = firstPt; // pt trigger
375 trackVariablesPair[4] = secondPt[j]; // pt
376 // trackVariablesPair[5] = fCentrality; // centrality
73609a48 377
378 // HBT like cut
379 if(fHBTCut && charge1 * charge2 > 0){
380 //if( dphi < 3 || deta < 0.01 ){ // VERSION 1
381 // continue;
382
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);
389
390 // optimization
391 if (TMath::Abs(deta) < 0.02 * 2.5 * 3) //twoTrackEfficiencyCutValue = 0.02 [default for dphicorrelations]
392 {
393 // phi in rad
c8ad9635 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];
73609a48 398
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);
402
403 const Float_t kLimit = 0.02 * 3;
404
405 Float_t dphistarminabs = 1e5;
406 Float_t dphistarmin = 1e5;
407
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);
412
413 if (dphistarabs < dphistarminabs) {
414 dphistarmin = dphistar;
415 dphistarminabs = dphistarabs;
416 }
417 }
418
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));
421 continue;
422 }
423 }
424 }
425 fHistHBTafter->Fill(deta,dphi);
426 }//HBT cut
427
428 // conversions
429 if(fConversionCut) {
430 if (charge1 * charge2 < 0) {
431 Double_t deta = firstEta - secondEta[j];
432 Double_t dphi = firstPhi - secondPhi[j];
433 fHistConversionbefore->Fill(deta,dphi);
434
435 Float_t m0 = 0.510e-3;
436 Float_t tantheta1 = 1e10;
437
438 // phi in rad
c8ad9635 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];
73609a48 443
444 if (firstEta < -1e-10 || firstEta > 1e-10)
445 tantheta1 = 2 * TMath::Exp(-firstEta) / ( 1 - TMath::Exp(-2*firstEta));
446
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]));
450
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);
453
454 Float_t masssqu = 2 * m0 * m0 + 2 * ( TMath::Sqrt(e1squ * e2squ) - ( firstPt * secondPt[j] * ( TMath::Cos(phi1rad - phi2rad) + 1.0 / tantheta1 / tantheta2 ) ) );
455
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));
458 continue;
459 }
460 fHistConversionafter->Fill(deta,dphi);
461 }
462 }//conversion cut
9afe3098 463
73609a48 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.);
8795e993 468 else {
469 //AliWarning(Form("Wrong charge combination: charge1 = %d and charge2 = %d",charge,charge2));
470 continue;
471 }
6acdbcb2 472 }//end of 2nd particle loop
473 }//end of 1st particle loop
0879e280 474}
475
0879e280 476//____________________________________________________________________//
9afe3098 477TH1D *AliBalancePsi::GetBalanceFunctionHistogram(Int_t iVariableSingle,
478 Int_t iVariablePair,
9afe3098 479 Double_t psiMin,
6acdbcb2 480 Double_t psiMax,
481 Double_t ptTriggerMin,
482 Double_t ptTriggerMax,
483 Double_t ptAssociatedMin,
484 Double_t ptAssociatedMax) {
9afe3098 485 //Returns the BF histogram, extracted from the 6 AliTHn objects
0879e280 486 //(private members) of the AliBalancePsi class.
621329de 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
9afe3098 489 // Psi_2
621329de 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);
9afe3098 496
6acdbcb2 497 // pt trigger
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);
505 }
506
507 // pt associated
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);
513 }
514
a38fd7f3 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));
516
9afe3098 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);
524
525 TH1D *gHistBalanceFunctionHistogram = 0x0;
526 if((hTemp1)&&(hTemp2)&&(hTemp3)&&(hTemp4)&&(hTemp5)&&(hTemp6)) {
527 gHistBalanceFunctionHistogram = (TH1D*)hTemp1->Clone();
528 gHistBalanceFunctionHistogram->Reset();
529
530 switch(iVariablePair) {
621329de 531 case 1:
c8ad9635 532 gHistBalanceFunctionHistogram->GetXaxis()->SetTitle("#Delta#eta");
533 gHistBalanceFunctionHistogram->GetYaxis()->SetTitle("B(#Delta#eta)");
9afe3098 534 break;
621329de 535 case 2:
c8ad9635 536 gHistBalanceFunctionHistogram->GetXaxis()->SetTitle("#Delta#varphi (rad)");
537 gHistBalanceFunctionHistogram->GetYaxis()->SetTitle("B(#Delta#varphi)");
9afe3098 538 break;
9afe3098 539 default:
540 break;
541 }
0879e280 542
0879e280 543 hTemp1->Sumw2();
544 hTemp2->Sumw2();
545 hTemp3->Sumw2();
546 hTemp4->Sumw2();
a38fd7f3 547 hTemp1->Add(hTemp3,-1.);
0879e280 548 hTemp1->Scale(1./hTemp5->GetEntries());
a38fd7f3 549 hTemp2->Add(hTemp4,-1.);
0879e280 550 hTemp2->Scale(1./hTemp6->GetEntries());
551 gHistBalanceFunctionHistogram->Add(hTemp1,hTemp2,1.,1.);
9afe3098 552 gHistBalanceFunctionHistogram->Scale(0.5);
0879e280 553 }
554
0879e280 555 return gHistBalanceFunctionHistogram;
556}
a38fd7f3 557
f0c5040c 558//____________________________________________________________________//
559TH1D *AliBalancePsi::GetBalanceFunctionHistogram2pMethod(Int_t iVariableSingle,
560 Int_t iVariablePair,
561 Double_t psiMin,
562 Double_t psiMax,
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
573 // Psi_2
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);
580
581 // pt trigger
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);
589 }
590
591 // pt associated
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);
597 }
598
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();
607
608 // Psi_2
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);
615
616 // pt trigger
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);
624 }
625
626 // pt associated
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);
632 }
633 // ============================================================================================
634
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));
636
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);
644
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 // ============================================================================================
654
655 TH1D *gHistBalanceFunctionHistogram = 0x0;
656 if((hTemp1)&&(hTemp2)&&(hTemp3)&&(hTemp4)&&(hTemp5)&&(hTemp6)) {
657 gHistBalanceFunctionHistogram = (TH1D*)hTemp1->Clone();
658 gHistBalanceFunctionHistogram->Reset();
659
660 switch(iVariablePair) {
661 case 1:
662 gHistBalanceFunctionHistogram->GetXaxis()->SetTitle("#Delta#eta");
663 gHistBalanceFunctionHistogram->GetYaxis()->SetTitle("B(#Delta#eta)");
664 break;
665 case 2:
666 gHistBalanceFunctionHistogram->GetXaxis()->SetTitle("#Delta#varphi (rad)");
667 gHistBalanceFunctionHistogram->GetYaxis()->SetTitle("B(#Delta#varphi)");
668 break;
669 default:
670 break;
671 }
672
673 hTemp1->Sumw2();
674 hTemp2->Sumw2();
675 hTemp3->Sumw2();
676 hTemp4->Sumw2();
677 hTemp1Mix->Sumw2();
678 hTemp2Mix->Sumw2();
679 hTemp3Mix->Sumw2();
680 hTemp4Mix->Sumw2();
681
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());
690
691 hTemp1->Divide(hTemp1Mix);
692 hTemp2->Divide(hTemp2Mix);
693 hTemp3->Divide(hTemp3Mix);
694 hTemp4->Divide(hTemp4Mix);
695
696 hTemp1->Add(hTemp3,-1.);
697 hTemp2->Add(hTemp4,-1.);
698
699 gHistBalanceFunctionHistogram->Add(hTemp1,hTemp2,1.,1.);
700 gHistBalanceFunctionHistogram->Scale(0.5);
701 }
702
703 return gHistBalanceFunctionHistogram;
704}
705
a38fd7f3 706//____________________________________________________________________//
621329de 707TH2D *AliBalancePsi::GetBalanceFunctionDeltaEtaDeltaPhi(Double_t psiMin,
6acdbcb2 708 Double_t psiMax,
709 Double_t ptTriggerMin,
710 Double_t ptTriggerMax,
711 Double_t ptAssociatedMin,
712 Double_t ptAssociatedMax) {
621329de 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";
719
720 // Psi_2
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);
727
6acdbcb2 728 // pt trigger
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);
736 }
737
738 // pt associated
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);
744 }
745
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)));
621329de 747
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);
755
756 TH2D *gHistBalanceFunctionHistogram = 0x0;
757 if((hTemp1)&&(hTemp2)&&(hTemp3)&&(hTemp4)&&(hTemp5)&&(hTemp6)) {
758 gHistBalanceFunctionHistogram = (TH2D*)hTemp1->Clone();
759 gHistBalanceFunctionHistogram->Reset();
c8ad9635 760 gHistBalanceFunctionHistogram->GetXaxis()->SetTitle("#Delta#eta");
761 gHistBalanceFunctionHistogram->GetYaxis()->SetTitle("#Delta#varphi (rad)");
762 gHistBalanceFunctionHistogram->GetZaxis()->SetTitle("B(#Delta#eta,#Delta#varphi)");
621329de 763
764 hTemp1->Sumw2();
765 hTemp2->Sumw2();
766 hTemp3->Sumw2();
767 hTemp4->Sumw2();
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);
774 }
775
776 return gHistBalanceFunctionHistogram;
777}
778
f0c5040c 779//____________________________________________________________________//
780TH2D *AliBalancePsi::GetBalanceFunctionDeltaEtaDeltaPhi2pMethod(Double_t psiMin,
781 Double_t psiMax,
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";
794
795 if(!bfMix){
796 AliError("balance function object for event mixing not available");
797 return NULL;
798 }
799
800 // Psi_2
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);
807
808 // pt trigger
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);
816 }
817
818 // pt associated
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);
824 }
825
826
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();
835
836 // Psi_2
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);
843
844 // pt trigger
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);
852 }
853
854 // pt associated
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);
860 }
861 // ============================================================================================
862
863
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)));
865
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);
873
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 // ============================================================================================
883
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)");
891
892 hTemp1->Sumw2();
893 hTemp2->Sumw2();
894 hTemp3->Sumw2();
895 hTemp4->Sumw2();
896 hTemp1Mix->Sumw2();
897 hTemp2Mix->Sumw2();
898 hTemp3Mix->Sumw2();
899 hTemp4Mix->Sumw2();
900
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());
909
910 hTemp1->Divide(hTemp1Mix);
911 hTemp2->Divide(hTemp2Mix);
912 hTemp3->Divide(hTemp3Mix);
913 hTemp4->Divide(hTemp4Mix);
914
915 hTemp1->Add(hTemp3,-1.);
916 hTemp2->Add(hTemp4,-1.);
917
918 gHistBalanceFunctionHistogram->Add(hTemp1,hTemp2,1.,1.);
919 gHistBalanceFunctionHistogram->Scale(0.5);
920 }
921
922 return gHistBalanceFunctionHistogram;
923}
924
925
621329de 926//____________________________________________________________________//
927TH2D *AliBalancePsi::GetCorrelationFunctionPN(Double_t psiMin,
6acdbcb2 928 Double_t psiMax,
929 Double_t ptTriggerMin,
930 Double_t ptTriggerMax,
931 Double_t ptAssociatedMin,
73dbfa54 932 Double_t ptAssociatedMax,
933 Bool_t normToTrig) {
a38fd7f3 934 //Returns the 2D correlation function for (+-) pairs
621329de 935 // Psi_2: axis 0
936 fHistP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
937 fHistPN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
621329de 938
6acdbcb2 939 // pt trigger
940 if((ptTriggerMin != -1.)&&(ptTriggerMax != -1.)) {
941 fHistP->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
942 fHistPN->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
5ef8eaa7 943 }
621329de 944
6acdbcb2 945 // pt associated
946 if((ptAssociatedMin != -1.)&&(ptAssociatedMax != -1.))
947 fHistPN->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
948
949 //fHistP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(-0.5,2.5);
950 //fHistPN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(-0.5,2.5);
951
952 //TH2D *gHistTest = dynamic_cast<TH2D *>(fHistP->Project(0,0,1));
953 //TCanvas *c1 = new TCanvas("c1","");
954 //c1->cd();
955 //if(!gHistTest){
956 //AliError("Projection of fHistP = NULL");
957 //return gHistTest;
958 //}
959 //else{
960 //gHistTest->DrawCopy("colz");
961 //}
962
963 //0:step, 1: Delta eta, 2: Delta phi
621329de 964 TH2D *gHist = dynamic_cast<TH2D *>(fHistPN->Project(0,1,2));
89c00c43 965 if(!gHist){
966 AliError("Projection of fHistPN = NULL");
967 return gHist;
968 }
969
6acdbcb2 970 //AliInfo(Form("Entries (test): %lf",(Double_t)(gHistTest->GetEntries())));
971 //AliInfo(Form("Entries (1D): %lf",(Double_t)(fHistP->Project(0,1)->GetEntries())));
972 //AliInfo(Form("Entries (2D): %lf",(Double_t)(fHistPN->Project(0,1,2)->GetEntries())));
621329de 973
6acdbcb2 974 //TCanvas *c2 = new TCanvas("c2","");
975 //c2->cd();
976 //fHistPN->Project(0,1,2)->DrawCopy("colz");
621329de 977
73dbfa54 978 if(normToTrig && (Double_t)(fHistP->Project(0,1)->GetEntries())!=0)
621329de 979 gHist->Scale(1./(Double_t)(fHistP->Project(0,1)->GetEntries()));
a38fd7f3 980
981 return gHist;
982}
983
984//____________________________________________________________________//
621329de 985TH2D *AliBalancePsi::GetCorrelationFunctionNP(Double_t psiMin,
6acdbcb2 986 Double_t psiMax,
987 Double_t ptTriggerMin,
988 Double_t ptTriggerMax,
989 Double_t ptAssociatedMin,
73dbfa54 990 Double_t ptAssociatedMax,
991 Bool_t normToTrig) {
a38fd7f3 992 //Returns the 2D correlation function for (+-) pairs
621329de 993 // Psi_2: axis 0
994 fHistN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
995 fHistNP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
a38fd7f3 996
6acdbcb2 997 // pt trigger
998 if((ptTriggerMin != -1.)&&(ptTriggerMax != -1.)) {
999 fHistN->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1000 fHistNP->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1001 }
1002
1003 // pt associated
1004 if((ptAssociatedMin != -1.)&&(ptAssociatedMax != -1.))
1005 fHistNP->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
1006
1007 //0:step, 1: Delta eta, 2: Delta phi
621329de 1008 TH2D *gHist = dynamic_cast<TH2D *>(fHistNP->Project(0,1,2));
89c00c43 1009 if(!gHist){
1010 AliError("Projection of fHistPN = NULL");
1011 return gHist;
1012 }
1013
a38fd7f3 1014 //Printf("Entries (1D): %lf",(Double_t)(fHistN->Project(0,2)->GetEntries()));
1015 //Printf("Entries (2D): %lf",(Double_t)(fHistNP->Project(0,2,3)->GetEntries()));
73dbfa54 1016 if(normToTrig && (Double_t)(fHistN->Project(0,1)->GetEntries())!=0)
621329de 1017 gHist->Scale(1./(Double_t)(fHistN->Project(0,1)->GetEntries()));
a38fd7f3 1018
1019 return gHist;
1020}
1021
1022//____________________________________________________________________//
621329de 1023TH2D *AliBalancePsi::GetCorrelationFunctionPP(Double_t psiMin,
6acdbcb2 1024 Double_t psiMax,
1025 Double_t ptTriggerMin,
1026 Double_t ptTriggerMax,
1027 Double_t ptAssociatedMin,
73dbfa54 1028 Double_t ptAssociatedMax,
1029 Bool_t normToTrig) {
a38fd7f3 1030 //Returns the 2D correlation function for (+-) pairs
621329de 1031 // Psi_2: axis 0
1032 fHistP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
1033 fHistPP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
6acdbcb2 1034
1035 // pt trigger
1036 if((ptTriggerMin != -1.)&&(ptTriggerMax != -1.)) {
1037 fHistP->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1038 fHistPP->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1039 }
1040
1041 // pt associated
1042 if((ptAssociatedMin != -1.)&&(ptAssociatedMax != -1.))
1043 fHistPP->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
a38fd7f3 1044
6acdbcb2 1045 //0:step, 1: Delta eta, 2: Delta phi
621329de 1046 TH2D *gHist = dynamic_cast<TH2D *>(fHistPP->Project(0,1,2));
89c00c43 1047 if(!gHist){
1048 AliError("Projection of fHistPN = NULL");
1049 return gHist;
1050 }
1051
a38fd7f3 1052 //Printf("Entries (1D): %lf",(Double_t)(fHistP->Project(0,2)->GetEntries()));
1053 //Printf("Entries (2D): %lf",(Double_t)(fHistPP->Project(0,2,3)->GetEntries()));
73dbfa54 1054 if(normToTrig && (Double_t)(fHistP->Project(0,1)->GetEntries())!=0)
621329de 1055 gHist->Scale(1./(Double_t)(fHistP->Project(0,1)->GetEntries()));
a38fd7f3 1056
1057 return gHist;
1058}
1059
1060//____________________________________________________________________//
621329de 1061TH2D *AliBalancePsi::GetCorrelationFunctionNN(Double_t psiMin,
6acdbcb2 1062 Double_t psiMax,
1063 Double_t ptTriggerMin,
1064 Double_t ptTriggerMax,
1065 Double_t ptAssociatedMin,
73dbfa54 1066 Double_t ptAssociatedMax,
1067 Bool_t normToTrig) {
a38fd7f3 1068 //Returns the 2D correlation function for (+-) pairs
621329de 1069 // Psi_2: axis 0
1070 fHistN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
1071 fHistNN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
6acdbcb2 1072
1073 // pt trigger
1074 if((ptTriggerMin != -1.)&&(ptTriggerMax != -1.)) {
1075 fHistN->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1076 fHistNN->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1077 }
1078
1079 // pt associated
1080 if((ptAssociatedMin != -1.)&&(ptAssociatedMax != -1.))
1081 fHistNN->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
a38fd7f3 1082
6acdbcb2 1083 //0:step, 1: Delta eta, 2: Delta phi
621329de 1084 TH2D *gHist = dynamic_cast<TH2D *>(fHistNN->Project(0,1,2));
89c00c43 1085 if(!gHist){
1086 AliError("Projection of fHistPN = NULL");
1087 return gHist;
1088 }
1089
a38fd7f3 1090 //Printf("Entries (1D): %lf",(Double_t)(fHistN->Project(0,2)->GetEntries()));
1091 //Printf("Entries (2D): %lf",(Double_t)(fHistNN->Project(0,2,3)->GetEntries()));
73dbfa54 1092 if(normToTrig && (Double_t)(fHistN->Project(0,1)->GetEntries())!=0)
621329de 1093 gHist->Scale(1./(Double_t)(fHistN->Project(0,1)->GetEntries()));
a38fd7f3 1094
1095 return gHist;
1096}
621329de 1097
7b610f27 1098//____________________________________________________________________//
1099TH2D *AliBalancePsi::GetCorrelationFunctionChargeIndependent(Double_t psiMin,
1100 Double_t psiMax,
1101 Double_t ptTriggerMin,
1102 Double_t ptTriggerMax,
1103 Double_t ptAssociatedMin,
1104 Double_t ptAssociatedMax,
1105 Bool_t normToTrig) {
1106 //Returns the 2D correlation function for the sum of all charge combination pairs
1107 // Psi_2: axis 0
1108 fHistN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
1109 fHistP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
1110 fHistNN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
1111 fHistPP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
1112 fHistNP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
1113 fHistPN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
1114
1115 // pt trigger
1116 if((ptTriggerMin != -1.)&&(ptTriggerMax != -1.)) {
1117 fHistN->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1118 fHistP->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1119 fHistNN->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1120 fHistPP->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1121 fHistNP->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1122 fHistPN->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1123 }
1124
1125 // pt associated
1126 if((ptAssociatedMin != -1.)&&(ptAssociatedMax != -1.))
1127 fHistNN->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
1128 fHistPP->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
1129 fHistNP->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
1130 fHistPN->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
1131
1132 //0:step, 1: Delta eta, 2: Delta phi
1133 TH2D *gHistNN = dynamic_cast<TH2D *>(fHistNN->Project(0,1,2));
1134 if(!gHistNN){
1135 AliError("Projection of fHistNN = NULL");
1136 return gHistNN;
1137 }
1138 TH2D *gHistPP = dynamic_cast<TH2D *>(fHistPP->Project(0,1,2));
1139 if(!gHistPP){
1140 AliError("Projection of fHistPP = NULL");
1141 return gHistPP;
1142 }
1143 TH2D *gHistNP = dynamic_cast<TH2D *>(fHistNP->Project(0,1,2));
1144 if(!gHistNP){
1145 AliError("Projection of fHistNP = NULL");
1146 return gHistNP;
1147 }
1148 TH2D *gHistPN = dynamic_cast<TH2D *>(fHistPN->Project(0,1,2));
1149 if(!gHistPN){
1150 AliError("Projection of fHistPN = NULL");
1151 return gHistPN;
1152 }
1153
1154 // sum all 2 particle histograms
1155 gHistNN->Add(gHistPP);
1156 gHistNN->Add(gHistNP);
1157 gHistNN->Add(gHistPN);
1158
1159 // if normalization to trigger particle, divide by sum of + and - triggers
1160 if(normToTrig && (Double_t)(fHistN->Project(0,1)->GetEntries())!=0 && (Double_t)(fHistP->Project(0,1)->GetEntries())!=0)
1161 gHistNN->Scale(1./(Double_t)(fHistN->Project(0,1)->GetEntries() + fHistN->Project(0,1)->GetEntries()));
1162
1163 return gHistNN;
1164}
1165
73609a48 1166//____________________________________________________________________//
1167Float_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) {
1168 //
1169 // calculates dphistar
1170 //
1171 Float_t dphistar = phi1 - phi2 - charge1 * bSign * TMath::ASin(0.075 * radius / pt1) + charge2 * bSign * TMath::ASin(0.075 * radius / pt2);
1172
1173 static const Double_t kPi = TMath::Pi();
1174
1175 // circularity
1176// if (dphistar > 2 * kPi)
1177// dphistar -= 2 * kPi;
1178// if (dphistar < -2 * kPi)
1179// dphistar += 2 * kPi;
1180
1181 if (dphistar > kPi)
1182 dphistar = kPi * 2 - dphistar;
1183 if (dphistar < -kPi)
1184 dphistar = -kPi * 2 - dphistar;
1185 if (dphistar > kPi) // might look funny but is needed
1186 dphistar = kPi * 2 - dphistar;
1187
1188 return dphistar;
1189}
1190