]> 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
Update class
[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),
f0fb4ac1 67 fConversionCut(kFALSE),
68 fEventClass("EventPlane"){
0879e280 69 // Default constructor
0879e280 70}
71
0879e280 72//____________________________________________________________________//
73AliBalancePsi::AliBalancePsi(const AliBalancePsi& balance):
9fd4b54e 74 TObject(balance), fShuffle(balance.fShuffle),
0879e280 75 fAnalysisLevel(balance.fAnalysisLevel),
76 fAnalyzedEvents(balance.fAnalyzedEvents),
77 fCentralityId(balance.fCentralityId),
78 fCentStart(balance.fCentStart),
79 fCentStop(balance.fCentStop),
9afe3098 80 fHistP(balance.fHistP),
81 fHistN(balance.fHistN),
82 fHistPN(balance.fHistPN),
83 fHistNP(balance.fHistNP),
84 fHistPP(balance.fHistPP),
85 fHistNN(balance.fHistNN),
73609a48 86 fHistHBTbefore(balance.fHistHBTbefore),
87 fHistHBTafter(balance.fHistHBTafter),
88 fHistConversionbefore(balance.fHistConversionbefore),
89 fHistConversionafter(balance.fHistConversionafter),
c8ad9635 90 fHistPsiMinusPhi(balance.fHistPsiMinusPhi),
73609a48 91 fPsiInterval(balance.fPsiInterval),
92 fHBTCut(balance.fHBTCut),
f0fb4ac1 93 fConversionCut(balance.fConversionCut),
94 fEventClass("EventPlane"){
0879e280 95 //copy constructor
0879e280 96}
97
98//____________________________________________________________________//
99AliBalancePsi::~AliBalancePsi() {
100 // Destructor
9afe3098 101 delete fHistP;
102 delete fHistN;
103 delete fHistPN;
104 delete fHistNP;
105 delete fHistPP;
106 delete fHistNN;
73609a48 107
108 delete fHistHBTbefore;
109 delete fHistHBTafter;
110 delete fHistConversionbefore;
111 delete fHistConversionafter;
c8ad9635 112 delete fHistPsiMinusPhi;
f0fb4ac1 113
0879e280 114}
115
116//____________________________________________________________________//
9afe3098 117void AliBalancePsi::InitHistograms() {
118 // single particle histograms
211b716d 119
120 // global switch disabling the reference
121 // (to avoid "Replacing existing TH1" if several wagons are created in train)
122 Bool_t oldStatus = TH1::AddDirectoryStatus();
123 TH1::AddDirectory(kFALSE);
124
9afe3098 125 Int_t anaSteps = 1; // analysis steps
9fd4b54e 126 Int_t iBinSingle[kTrackVariablesSingle]; // binning for track variables
127 Double_t* dBinsSingle[kTrackVariablesSingle]; // bins for track variables
128 TString axisTitleSingle[kTrackVariablesSingle]; // axis titles for track variables
9afe3098 129
130 // two particle histograms
9fd4b54e 131 Int_t iBinPair[kTrackVariablesPair]; // binning for track variables
132 Double_t* dBinsPair[kTrackVariablesPair]; // bins for track variables
133 TString axisTitlePair[kTrackVariablesPair]; // axis titles for track variables
f0fb4ac1 134 /**********************************************************
135
136 ======> Modification: Change Event Classification Scheme
137
138 ---> fEventClass == "EventPlane"
139
140 Default operation with Event Plane
141
142 ---> fEventClass == "Multiplicity"
143
144 Work with kTPCITStracklet multiplicity (from GetReferenceMultiplicity)
145
146 ---> fEventClass == "Centrality"
147
148 Work with Centrality Bins
149
150 ***********************************************************/
151
152 //--- Multiplicity Bins ------------------------------------
153 const Int_t kMultBins = 8;
154 //A first rough attempt at four bins
155 Double_t kMultBinLimits[kMultBins+1]={0,10,20,30,40,50,60,70,80};
156 //----------------------------------------------------------
157
158 //--- Centrality Bins --------------------------------------
159 const Int_t kNCentralityBins = 9;
160 Double_t centralityBins[kNCentralityBins+1] = {0.,5.,10.,20.,30.,40.,50.,60.,70.,80.};
161 //----------------------------------------------------------
162
163 //--- Event Plane Bins -------------------------------------
164 //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)
165 const Int_t kNPsi2Bins = 4;
166 Double_t psi2Bins[kNPsi2Bins+1] = {-0.5,0.5,1.5,2.5,3.5};
167 //----------------------------------------------------------
168
169 //Depending on fEventClass Variable, do one thing or the other...
170 if(fEventClass == "Multiplicity"){
171 iBinSingle[0] = kMultBins;
172 dBinsSingle[0] = kMultBinLimits;
173 axisTitleSingle[0] = "kTPCITStracklet multiplicity";
174 iBinPair[0] = kMultBins;
175 dBinsPair[0] = kMultBinLimits;
176 axisTitlePair[0] = "kTPCITStracklet multiplicity";
177 }
178 if(fEventClass == "Centrality"){
179 iBinSingle[0] = kNCentralityBins;
180 dBinsSingle[0] = centralityBins;
181 axisTitleSingle[0] = "Centrality percentile [%]";
182 iBinPair[0] = kNCentralityBins;
183 dBinsPair[0] = centralityBins;
184 axisTitlePair[0] = "Centrality percentile [%]";
185 }
186 if(fEventClass == "EventPlane"){
187 iBinSingle[0] = kNPsi2Bins;
188 dBinsSingle[0] = psi2Bins;
189 axisTitleSingle[0] = "#varphi - #Psi_{2} (a.u.)";
190 iBinPair[0] = kNPsi2Bins;
191 dBinsPair[0] = psi2Bins;
192 axisTitlePair[0] = "#varphi - #Psi_{2} (a.u.)";
193 }
9afe3098 194
621329de 195 // delta eta
196 const Int_t kNDeltaEtaBins = 80;
9afe3098 197 Double_t deltaEtaBins[kNDeltaEtaBins+1];
198 for(Int_t i = 0; i < kNDeltaEtaBins+1; i++)
621329de 199 deltaEtaBins[i] = -2.0 + i * 0.05;
200 iBinPair[1] = kNDeltaEtaBins;
201 dBinsPair[1] = deltaEtaBins;
c8ad9635 202 axisTitlePair[1] = "#Delta#eta";
621329de 203
204 // delta phi
9afe3098 205 const Int_t kNDeltaPhiBins = 72;
206 Double_t deltaPhiBins[kNDeltaPhiBins+1];
207 for(Int_t i = 0; i < kNDeltaPhiBins+1; i++){
c8ad9635 208 //deltaPhiBins[i] = -180.0 + i * 5.;
209 deltaPhiBins[i] = -TMath::Pi()/2. + i * 5.*TMath::Pi()/180.;
9afe3098 210 }
621329de 211 iBinPair[2] = kNDeltaPhiBins;
212 dBinsPair[2] = deltaPhiBins;
c8ad9635 213 axisTitlePair[2] = "#Delta#varphi (rad)";
9afe3098 214
a38fd7f3 215 // pt(trigger-associated)
8795e993 216 const Int_t kNPtBins = 16;
217 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.};
218 //for(Int_t i = 0; i < kNPtBins+1; i++){
219 //ptBins[i] = 0.2 + i * 0.5;
220 //}
621329de 221 iBinSingle[1] = kNPtBins;
222 dBinsSingle[1] = ptBins;
c8ad9635 223 axisTitleSingle[1] = "p_{T,trig.} (GeV/c)";
621329de 224
225 iBinPair[3] = kNPtBins;
226 dBinsPair[3] = ptBins;
c8ad9635 227 axisTitlePair[3] = "p_{T,trig.} (GeV/c)";
a38fd7f3 228
229 iBinPair[4] = kNPtBins;
230 dBinsPair[4] = ptBins;
c8ad9635 231 axisTitlePair[4] = "p_{T,assoc.} (GeV/c)";
9afe3098 232
233 TString histName;
234 //+ triggered particles
235 histName = "fHistP";
9fd4b54e 236 if(fShuffle) histName.Append("_shuffle");
9afe3098 237 if(fCentralityId) histName += fCentralityId.Data();
9fd4b54e 238 fHistP = new AliTHn(histName.Data(),histName.Data(),anaSteps,kTrackVariablesSingle,iBinSingle);
239 for (Int_t j=0; j<kTrackVariablesSingle; j++) {
9afe3098 240 fHistP->SetBinLimits(j, dBinsSingle[j]);
241 fHistP->SetVarTitle(j, axisTitleSingle[j]);
0879e280 242 }
9afe3098 243
244 //- triggered particles
245 histName = "fHistN";
9fd4b54e 246 if(fShuffle) histName.Append("_shuffle");
9afe3098 247 if(fCentralityId) histName += fCentralityId.Data();
9fd4b54e 248 fHistN = new AliTHn(histName.Data(),histName.Data(),anaSteps,kTrackVariablesSingle,iBinSingle);
249 for (Int_t j=0; j<kTrackVariablesSingle; j++) {
9afe3098 250 fHistN->SetBinLimits(j, dBinsSingle[j]);
251 fHistN->SetVarTitle(j, axisTitleSingle[j]);
0879e280 252 }
9afe3098 253
254 //+- pairs
255 histName = "fHistPN";
9fd4b54e 256 if(fShuffle) histName.Append("_shuffle");
9afe3098 257 if(fCentralityId) histName += fCentralityId.Data();
9fd4b54e 258 fHistPN = new AliTHn(histName.Data(),histName.Data(),anaSteps, kTrackVariablesPair, iBinPair);
259 for (Int_t j=0; j<kTrackVariablesPair; j++) {
9afe3098 260 fHistPN->SetBinLimits(j, dBinsPair[j]);
261 fHistPN->SetVarTitle(j, axisTitlePair[j]);
0879e280 262 }
0879e280 263
9afe3098 264 //-+ pairs
265 histName = "fHistNP";
9fd4b54e 266 if(fShuffle) histName.Append("_shuffle");
9afe3098 267 if(fCentralityId) histName += fCentralityId.Data();
9fd4b54e 268 fHistNP = new AliTHn(histName.Data(),histName.Data(),anaSteps, kTrackVariablesPair, iBinPair);
269 for (Int_t j=0; j<kTrackVariablesPair; j++) {
9afe3098 270 fHistNP->SetBinLimits(j, dBinsPair[j]);
271 fHistNP->SetVarTitle(j, axisTitlePair[j]);
0879e280 272 }
0879e280 273
9afe3098 274 //++ pairs
275 histName = "fHistPP";
9fd4b54e 276 if(fShuffle) histName.Append("_shuffle");
9afe3098 277 if(fCentralityId) histName += fCentralityId.Data();
9fd4b54e 278 fHistPP = new AliTHn(histName.Data(),histName.Data(),anaSteps, kTrackVariablesPair, iBinPair);
279 for (Int_t j=0; j<kTrackVariablesPair; j++) {
9afe3098 280 fHistPP->SetBinLimits(j, dBinsPair[j]);
281 fHistPP->SetVarTitle(j, axisTitlePair[j]);
282 }
283
284 //-- pairs
285 histName = "fHistNN";
9fd4b54e 286 if(fShuffle) histName.Append("_shuffle");
9afe3098 287 if(fCentralityId) histName += fCentralityId.Data();
9fd4b54e 288 fHistNN = new AliTHn(histName.Data(),histName.Data(),anaSteps, kTrackVariablesPair, iBinPair);
289 for (Int_t j=0; j<kTrackVariablesPair; j++) {
9afe3098 290 fHistNN->SetBinLimits(j, dBinsPair[j]);
291 fHistNN->SetVarTitle(j, axisTitlePair[j]);
0879e280 292 }
f06d59b3 293 AliInfo("Finished setting up the AliTHn");
73609a48 294
295 // QA histograms
c8ad9635 296 fHistHBTbefore = new TH2D("fHistHBTbefore","before HBT cut",200,0,2,200,0,2.*TMath::Pi());
297 fHistHBTafter = new TH2D("fHistHBTafter","after HBT cut",200,0,2,200,0,2.*TMath::Pi());
298 fHistConversionbefore = new TH2D("fHistConversionbefore","before Conversion cut",200,0,2,200,0,2.*TMath::Pi());
299 fHistConversionafter = new TH2D("fHistConversionafter","after Conversion cut",200,0,2,200,0,2.*TMath::Pi());
300 fHistPsiMinusPhi = new TH2D("fHistPsiMinusPhi","",4,-0.5,3.5,100,0,2.*TMath::Pi());
211b716d 301
302 TH1::AddDirectory(oldStatus);
303
0879e280 304}
305
306//____________________________________________________________________//
f06d59b3 307void AliBalancePsi::CalculateBalance(Double_t gReactionPlane,
308 TObjArray *particles,
73609a48 309 TObjArray *particlesMixed,
f0fb4ac1 310 Float_t bSign,
311 Double_t kMultorCent) {
0879e280 312 // Calculates the balance function
313 fAnalyzedEvents++;
f0fb4ac1 314
0879e280 315 // Initialize histograms if not done yet
9afe3098 316 if(!fHistPN){
0879e280 317 AliWarning("Histograms not yet initialized! --> Will be done now");
318 AliWarning("This works only in local mode --> Add 'gBalance->InitHistograms()' in your configBalanceFunction");
319 InitHistograms();
320 }
321
9fd4b54e 322 Double_t trackVariablesSingle[kTrackVariablesSingle];
323 Double_t trackVariablesPair[kTrackVariablesPair];
f06d59b3 324
325 if (!particles){
326 AliWarning("particles TObjArray is NULL pointer --> return");
327 return;
328 }
9afe3098 329
f06d59b3 330 // define end of particle loops
331 Int_t iMax = particles->GetEntriesFast();
332 Int_t jMax = iMax;
333 if (particlesMixed)
334 jMax = particlesMixed->GetEntriesFast();
335
336 // Eta() is extremely time consuming, therefore cache it for the inner loop here:
337 TObjArray* particlesSecond = (particlesMixed) ? particlesMixed : particles;
338
339 TArrayF secondEta(jMax);
340 TArrayF secondPhi(jMax);
341 TArrayF secondPt(jMax);
c1847400 342 TArrayS secondCharge(jMax);
f06d59b3 343
344 for (Int_t i=0; i<jMax; i++){
345 secondEta[i] = ((AliVParticle*) particlesSecond->At(i))->Eta();
346 secondPhi[i] = ((AliVParticle*) particlesSecond->At(i))->Phi();
347 secondPt[i] = ((AliVParticle*) particlesSecond->At(i))->Pt();
c1847400 348 secondCharge[i] = (Short_t)((AliVParticle*) particlesSecond->At(i))->Charge();
f06d59b3 349 }
350
351 // 1st particle loop
73609a48 352 for (Int_t i = 0; i < iMax; i++) {
6acdbcb2 353 AliVParticle* firstParticle = (AliVParticle*) particles->At(i);
354
355 // some optimization
356 Float_t firstEta = firstParticle->Eta();
357 Float_t firstPhi = firstParticle->Phi();
358 Float_t firstPt = firstParticle->Pt();
359
360 // Event plane (determine psi bin)
361 Double_t gPsiMinusPhi = 0.;
362 Double_t gPsiMinusPhiBin = -10.;
363 gPsiMinusPhi = TMath::Abs(firstPhi - gReactionPlane);
364 //in-plane
c8ad9635 365 if((gPsiMinusPhi <= 7.5*TMath::DegToRad())||
366 ((172.5*TMath::DegToRad() <= gPsiMinusPhi)&&(gPsiMinusPhi <= 187.5*TMath::DegToRad())))
6acdbcb2 367 gPsiMinusPhiBin = 0.0;
368 //intermediate
c8ad9635 369 else if(((37.5*TMath::DegToRad() <= gPsiMinusPhi)&&(gPsiMinusPhi <= 52.5*TMath::DegToRad()))||
370 ((127.5*TMath::DegToRad() <= gPsiMinusPhi)&&(gPsiMinusPhi <= 142.5*TMath::DegToRad()))||
371 ((217.5*TMath::DegToRad() <= gPsiMinusPhi)&&(gPsiMinusPhi <= 232.5*TMath::DegToRad()))||
372 ((307.5*TMath::DegToRad() <= gPsiMinusPhi)&&(gPsiMinusPhi <= 322.5*TMath::DegToRad())))
6acdbcb2 373 gPsiMinusPhiBin = 1.0;
374 //out of plane
c8ad9635 375 else if(((82.5*TMath::DegToRad() <= gPsiMinusPhi)&&(gPsiMinusPhi <= 97.5*TMath::DegToRad()))||
376 ((262.5*TMath::DegToRad() <= gPsiMinusPhi)&&(gPsiMinusPhi <= 277.5*TMath::DegToRad())))
6acdbcb2 377 gPsiMinusPhiBin = 2.0;
378 //everything else
379 else
380 gPsiMinusPhiBin = 3.0;
381
c8ad9635 382 fHistPsiMinusPhi->Fill(gPsiMinusPhiBin,gPsiMinusPhi);
73609a48 383
384 Short_t charge1 = (Short_t) firstParticle->Charge();
6acdbcb2 385
386 trackVariablesSingle[0] = gPsiMinusPhiBin;
f0fb4ac1 387 trackVariablesSingle[1] = firstPt;
bb473a33 388 if(fEventClass=="Multiplicity" || fEventClass == "Centrality" ) trackVariablesSingle[0] = kMultorCent;
6acdbcb2 389
390 //fill single particle histograms
73609a48 391 if(charge1 > 0) fHistP->Fill(trackVariablesSingle,0,1.);
392 else if(charge1 < 0) fHistN->Fill(trackVariablesSingle,0,1.);
6acdbcb2 393
c1847400 394 // 2nd particle loop
395 for(Int_t j = 0; j < jMax; j++) {
396
397 if(!particlesMixed && j == i) continue; // no auto correlations (only for non mixing)
398
399 Short_t charge2 = secondCharge[j];
9afe3098 400
f0fb4ac1 401 trackVariablesPair[0] = trackVariablesSingle[0];
6acdbcb2 402 trackVariablesPair[1] = firstEta - secondEta[j]; // delta eta
403 trackVariablesPair[2] = firstPhi - secondPhi[j]; // delta phi
c8ad9635 404 //if (trackVariablesPair[2] > 180.) // delta phi between -180 and 180
405 //trackVariablesPair[2] -= 360.;
406 //if (trackVariablesPair[2] < - 180.)
407 //trackVariablesPair[2] += 360.;
408 if (trackVariablesPair[2] > TMath::Pi()) // delta phi between -pi and pi
409 trackVariablesPair[2] -= 2.*TMath::Pi();
410 if (trackVariablesPair[2] < - TMath::Pi())
411 trackVariablesPair[2] += 2.*TMath::Pi();
412 if (trackVariablesPair[2] < - TMath::Pi()/2.)
413 trackVariablesPair[2] += 2.*TMath::Pi();
9afe3098 414
6acdbcb2 415 trackVariablesPair[3] = firstPt; // pt trigger
416 trackVariablesPair[4] = secondPt[j]; // pt
417 // trackVariablesPair[5] = fCentrality; // centrality
73609a48 418
419 // HBT like cut
420 if(fHBTCut && charge1 * charge2 > 0){
421 //if( dphi < 3 || deta < 0.01 ){ // VERSION 1
422 // continue;
423
424 Double_t deta = firstEta - secondEta[j];
425 Double_t dphi = firstPhi - secondPhi[j];
426 // VERSION 2 (Taken from DPhiCorrelations)
427 // the variables & cuthave been developed by the HBT group
428 // see e.g. https://indico.cern.ch/materialDisplay.py?contribId=36&sessionId=6&materialId=slides&confId=142700
429 fHistHBTbefore->Fill(deta,dphi);
430
431 // optimization
432 if (TMath::Abs(deta) < 0.02 * 2.5 * 3) //twoTrackEfficiencyCutValue = 0.02 [default for dphicorrelations]
433 {
434 // phi in rad
c8ad9635 435 //Float_t phi1rad = firstPhi*TMath::DegToRad();
436 //Float_t phi2rad = secondPhi[j]*TMath::DegToRad();
437 Float_t phi1rad = firstPhi;
438 Float_t phi2rad = secondPhi[j];
73609a48 439
440 // check first boundaries to see if is worth to loop and find the minimum
441 Float_t dphistar1 = GetDPhiStar(phi1rad, firstPt, charge1, phi2rad, secondPt[j], charge2, 0.8, bSign);
442 Float_t dphistar2 = GetDPhiStar(phi1rad, firstPt, charge1, phi2rad, secondPt[j], charge2, 2.5, bSign);
443
444 const Float_t kLimit = 0.02 * 3;
445
446 Float_t dphistarminabs = 1e5;
447 Float_t dphistarmin = 1e5;
448
449 if (TMath::Abs(dphistar1) < kLimit || TMath::Abs(dphistar2) < kLimit || dphistar1 * dphistar2 < 0 ) {
450 for (Double_t rad=0.8; rad<2.51; rad+=0.01) {
451 Float_t dphistar = GetDPhiStar(phi1rad, firstPt, charge1, phi2rad, secondPt[j], charge2, rad, bSign);
452 Float_t dphistarabs = TMath::Abs(dphistar);
453
454 if (dphistarabs < dphistarminabs) {
455 dphistarmin = dphistar;
456 dphistarminabs = dphistarabs;
457 }
458 }
459
460 if (dphistarminabs < 0.02 && TMath::Abs(deta) < 0.02) {
461 //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));
462 continue;
463 }
464 }
465 }
466 fHistHBTafter->Fill(deta,dphi);
467 }//HBT cut
468
469 // conversions
470 if(fConversionCut) {
471 if (charge1 * charge2 < 0) {
472 Double_t deta = firstEta - secondEta[j];
473 Double_t dphi = firstPhi - secondPhi[j];
474 fHistConversionbefore->Fill(deta,dphi);
475
476 Float_t m0 = 0.510e-3;
477 Float_t tantheta1 = 1e10;
478
479 // phi in rad
c8ad9635 480 //Float_t phi1rad = firstPhi*TMath::DegToRad();
481 //Float_t phi2rad = secondPhi[j]*TMath::DegToRad();
482 Float_t phi1rad = firstPhi;
483 Float_t phi2rad = secondPhi[j];
73609a48 484
485 if (firstEta < -1e-10 || firstEta > 1e-10)
486 tantheta1 = 2 * TMath::Exp(-firstEta) / ( 1 - TMath::Exp(-2*firstEta));
487
488 Float_t tantheta2 = 1e10;
489 if (secondEta[j] < -1e-10 || secondEta[j] > 1e-10)
490 tantheta2 = 2 * TMath::Exp(-secondEta[j]) / ( 1 - TMath::Exp(-2*secondEta[j]));
491
492 Float_t e1squ = m0 * m0 + firstPt * firstPt * (1.0 + 1.0 / tantheta1 / tantheta1);
493 Float_t e2squ = m0 * m0 + secondPt[j] * secondPt[j] * (1.0 + 1.0 / tantheta2 / tantheta2);
494
495 Float_t masssqu = 2 * m0 * m0 + 2 * ( TMath::Sqrt(e1squ * e2squ) - ( firstPt * secondPt[j] * ( TMath::Cos(phi1rad - phi2rad) + 1.0 / tantheta1 / tantheta2 ) ) );
496
497 if (masssqu < 0.04*0.04){
498 //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));
499 continue;
500 }
501 fHistConversionafter->Fill(deta,dphi);
502 }
503 }//conversion cut
9afe3098 504
73609a48 505 if( charge1 > 0 && charge2 < 0) fHistPN->Fill(trackVariablesPair,0,1.);
506 else if( charge1 < 0 && charge2 > 0) fHistNP->Fill(trackVariablesPair,0,1.);
507 else if( charge1 > 0 && charge2 > 0) fHistPP->Fill(trackVariablesPair,0,1.);
508 else if( charge1 < 0 && charge2 < 0) fHistNN->Fill(trackVariablesPair,0,1.);
8795e993 509 else {
510 //AliWarning(Form("Wrong charge combination: charge1 = %d and charge2 = %d",charge,charge2));
511 continue;
512 }
6acdbcb2 513 }//end of 2nd particle loop
514 }//end of 1st particle loop
0879e280 515}
516
0879e280 517//____________________________________________________________________//
9afe3098 518TH1D *AliBalancePsi::GetBalanceFunctionHistogram(Int_t iVariableSingle,
519 Int_t iVariablePair,
9afe3098 520 Double_t psiMin,
6acdbcb2 521 Double_t psiMax,
522 Double_t ptTriggerMin,
523 Double_t ptTriggerMax,
524 Double_t ptAssociatedMin,
525 Double_t ptAssociatedMax) {
9afe3098 526 //Returns the BF histogram, extracted from the 6 AliTHn objects
0879e280 527 //(private members) of the AliBalancePsi class.
621329de 528 //iVariableSingle: 0(phi-Psi), 1(pt-trigger)
529 //iVariablePair: 0(phi-Psi) 1(Delta eta), 2(Delta phi), 3(pt-trigger), 4(pt-associated
9afe3098 530 // Psi_2
621329de 531 fHistP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
532 fHistN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
533 fHistPN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
534 fHistNP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
535 fHistPP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
536 fHistNN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
9afe3098 537
6acdbcb2 538 // pt trigger
539 if((ptTriggerMin != -1.)&&(ptTriggerMax != -1.)) {
540 fHistP->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
541 fHistN->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
542 fHistPN->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
543 fHistNP->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
544 fHistPP->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
545 fHistNN->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
546 }
547
548 // pt associated
549 if((ptAssociatedMin != -1.)&&(ptAssociatedMax != -1.)) {
550 fHistPN->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
551 fHistNP->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
552 fHistPP->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
553 fHistNN->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
554 }
555
a38fd7f3 556 //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));
557
9afe3098 558 // Project into the wanted space (1st: analysis step, 2nd: axis)
559 TH1D* hTemp1 = (TH1D*)fHistPN->Project(0,iVariablePair);
560 TH1D* hTemp2 = (TH1D*)fHistNP->Project(0,iVariablePair);
561 TH1D* hTemp3 = (TH1D*)fHistPP->Project(0,iVariablePair);
562 TH1D* hTemp4 = (TH1D*)fHistNN->Project(0,iVariablePair);
563 TH1D* hTemp5 = (TH1D*)fHistP->Project(0,iVariableSingle);
564 TH1D* hTemp6 = (TH1D*)fHistN->Project(0,iVariableSingle);
565
566 TH1D *gHistBalanceFunctionHistogram = 0x0;
567 if((hTemp1)&&(hTemp2)&&(hTemp3)&&(hTemp4)&&(hTemp5)&&(hTemp6)) {
568 gHistBalanceFunctionHistogram = (TH1D*)hTemp1->Clone();
569 gHistBalanceFunctionHistogram->Reset();
570
571 switch(iVariablePair) {
621329de 572 case 1:
c8ad9635 573 gHistBalanceFunctionHistogram->GetXaxis()->SetTitle("#Delta#eta");
574 gHistBalanceFunctionHistogram->GetYaxis()->SetTitle("B(#Delta#eta)");
9afe3098 575 break;
621329de 576 case 2:
c8ad9635 577 gHistBalanceFunctionHistogram->GetXaxis()->SetTitle("#Delta#varphi (rad)");
578 gHistBalanceFunctionHistogram->GetYaxis()->SetTitle("B(#Delta#varphi)");
9afe3098 579 break;
9afe3098 580 default:
581 break;
582 }
0879e280 583
0879e280 584 hTemp1->Sumw2();
585 hTemp2->Sumw2();
586 hTemp3->Sumw2();
587 hTemp4->Sumw2();
a38fd7f3 588 hTemp1->Add(hTemp3,-1.);
0879e280 589 hTemp1->Scale(1./hTemp5->GetEntries());
a38fd7f3 590 hTemp2->Add(hTemp4,-1.);
0879e280 591 hTemp2->Scale(1./hTemp6->GetEntries());
592 gHistBalanceFunctionHistogram->Add(hTemp1,hTemp2,1.,1.);
9afe3098 593 gHistBalanceFunctionHistogram->Scale(0.5);
0879e280 594 }
595
0879e280 596 return gHistBalanceFunctionHistogram;
597}
a38fd7f3 598
f0c5040c 599//____________________________________________________________________//
600TH1D *AliBalancePsi::GetBalanceFunctionHistogram2pMethod(Int_t iVariableSingle,
601 Int_t iVariablePair,
602 Double_t psiMin,
603 Double_t psiMax,
604 Double_t ptTriggerMin,
605 Double_t ptTriggerMax,
606 Double_t ptAssociatedMin,
607 Double_t ptAssociatedMax,
608 AliBalancePsi *bfMix) {
609 //Returns the BF histogram, extracted from the 6 AliTHn objects
610 //after dividing each correlation function by the Event Mixing one
611 //(private members) of the AliBalancePsi class.
612 //iVariableSingle: 0(phi-Psi), 1(pt-trigger)
613 //iVariablePair: 0(phi-Psi) 1(Delta eta), 2(Delta phi), 3(pt-trigger), 4(pt-associated
614 // Psi_2
615 fHistP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
616 fHistN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
617 fHistPN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
618 fHistNP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
619 fHistPP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
620 fHistNN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
621
622 // pt trigger
623 if((ptTriggerMin != -1.)&&(ptTriggerMax != -1.)) {
624 fHistP->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
625 fHistN->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
626 fHistPN->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
627 fHistNP->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
628 fHistPP->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
629 fHistNN->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
630 }
631
632 // pt associated
633 if((ptAssociatedMin != -1.)&&(ptAssociatedMax != -1.)) {
634 fHistPN->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
635 fHistNP->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
636 fHistPP->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
637 fHistNN->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
638 }
639
640 // ============================================================================================
641 // the same for event mixing
642 AliTHn *fHistPMix = bfMix->GetHistNp();
643 AliTHn *fHistNMix = bfMix->GetHistNn();
644 AliTHn *fHistPNMix = bfMix->GetHistNpn();
645 AliTHn *fHistNPMix = bfMix->GetHistNnp();
646 AliTHn *fHistPPMix = bfMix->GetHistNpp();
647 AliTHn *fHistNNMix = bfMix->GetHistNnn();
648
649 // Psi_2
650 fHistPMix->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
651 fHistNMix->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
652 fHistPNMix->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
653 fHistNPMix->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
654 fHistPPMix->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
655 fHistNNMix->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
656
657 // pt trigger
658 if((ptTriggerMin != -1.)&&(ptTriggerMax != -1.)) {
659 fHistPMix->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
660 fHistNMix->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
661 fHistPNMix->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
662 fHistNPMix->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
663 fHistPPMix->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
664 fHistNNMix->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
665 }
666
667 // pt associated
668 if((ptAssociatedMin != -1.)&&(ptAssociatedMax != -1.)) {
669 fHistPNMix->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
670 fHistNPMix->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
671 fHistPPMix->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
672 fHistNNMix->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
673 }
674 // ============================================================================================
675
676 //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));
677
678 // Project into the wanted space (1st: analysis step, 2nd: axis)
679 TH1D* hTemp1 = (TH1D*)fHistPN->Project(0,iVariablePair);
680 TH1D* hTemp2 = (TH1D*)fHistNP->Project(0,iVariablePair);
681 TH1D* hTemp3 = (TH1D*)fHistPP->Project(0,iVariablePair);
682 TH1D* hTemp4 = (TH1D*)fHistNN->Project(0,iVariablePair);
683 TH1D* hTemp5 = (TH1D*)fHistP->Project(0,iVariableSingle);
684 TH1D* hTemp6 = (TH1D*)fHistN->Project(0,iVariableSingle);
685
686 // ============================================================================================
687 // the same for event mixing
f0fb4ac1 688 TH2D* hTemp1Mix = (TH2D*)fHistPNMix->Project(0,iVariablePair);
689 TH2D* hTemp2Mix = (TH2D*)fHistNPMix->Project(0,iVariablePair);
690 TH2D* hTemp3Mix = (TH2D*)fHistPPMix->Project(0,iVariablePair);
691 TH2D* hTemp4Mix = (TH2D*)fHistNNMix->Project(0,iVariablePair);
692 TH1D* hTemp5Mix = (TH1D*)fHistPMix->Project(0,iVariableSingle);
693 TH1D* hTemp6Mix = (TH1D*)fHistNMix->Project(0,iVariableSingle);
f0c5040c 694 // ============================================================================================
695
696 TH1D *gHistBalanceFunctionHistogram = 0x0;
697 if((hTemp1)&&(hTemp2)&&(hTemp3)&&(hTemp4)&&(hTemp5)&&(hTemp6)) {
698 gHistBalanceFunctionHistogram = (TH1D*)hTemp1->Clone();
699 gHistBalanceFunctionHistogram->Reset();
700
701 switch(iVariablePair) {
702 case 1:
703 gHistBalanceFunctionHistogram->GetXaxis()->SetTitle("#Delta#eta");
704 gHistBalanceFunctionHistogram->GetYaxis()->SetTitle("B(#Delta#eta)");
705 break;
706 case 2:
707 gHistBalanceFunctionHistogram->GetXaxis()->SetTitle("#Delta#varphi (rad)");
708 gHistBalanceFunctionHistogram->GetYaxis()->SetTitle("B(#Delta#varphi)");
709 break;
710 default:
711 break;
712 }
713
714 hTemp1->Sumw2();
715 hTemp2->Sumw2();
716 hTemp3->Sumw2();
717 hTemp4->Sumw2();
718 hTemp1Mix->Sumw2();
719 hTemp2Mix->Sumw2();
720 hTemp3Mix->Sumw2();
721 hTemp4Mix->Sumw2();
722
723 hTemp1->Scale(1./hTemp5->GetEntries());
724 hTemp3->Scale(1./hTemp5->GetEntries());
725 hTemp2->Scale(1./hTemp6->GetEntries());
726 hTemp4->Scale(1./hTemp6->GetEntries());
f0fb4ac1 727 hTemp1Mix->Scale(1./hTemp5Mix->GetEntries());
728 hTemp3Mix->Scale(1./hTemp5Mix->GetEntries());
729 hTemp2Mix->Scale(1./hTemp6Mix->GetEntries());
730 hTemp4Mix->Scale(1./hTemp6Mix->GetEntries());
f0c5040c 731
732 hTemp1->Divide(hTemp1Mix);
733 hTemp2->Divide(hTemp2Mix);
734 hTemp3->Divide(hTemp3Mix);
735 hTemp4->Divide(hTemp4Mix);
736
737 hTemp1->Add(hTemp3,-1.);
738 hTemp2->Add(hTemp4,-1.);
739
740 gHistBalanceFunctionHistogram->Add(hTemp1,hTemp2,1.,1.);
741 gHistBalanceFunctionHistogram->Scale(0.5);
742 }
743
744 return gHistBalanceFunctionHistogram;
745}
746
a38fd7f3 747//____________________________________________________________________//
621329de 748TH2D *AliBalancePsi::GetBalanceFunctionDeltaEtaDeltaPhi(Double_t psiMin,
6acdbcb2 749 Double_t psiMax,
750 Double_t ptTriggerMin,
751 Double_t ptTriggerMax,
752 Double_t ptAssociatedMin,
753 Double_t ptAssociatedMax) {
621329de 754 //Returns the BF histogram in Delta eta vs Delta phi,
755 //extracted from the 6 AliTHn objects
756 //(private members) of the AliBalancePsi class.
757 //iVariableSingle: 0(phi-Psi), 1(pt-trigger)
758 //iVariablePair: 0(phi-Psi) 1(Delta eta), 2(Delta phi), 3(pt-trigger), 4(pt-associated
759 TString histName = "gHistBalanceFunctionHistogram2D";
760
761 // Psi_2
762 fHistP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
763 fHistN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
764 fHistPN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
765 fHistNP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
766 fHistPP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
767 fHistNN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
768
6acdbcb2 769 // pt trigger
770 if((ptTriggerMin != -1.)&&(ptTriggerMax != -1.)) {
771 fHistP->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
772 fHistN->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
773 fHistPN->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
774 fHistNP->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
775 fHistPP->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
776 fHistNN->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
777 }
778
779 // pt associated
780 if((ptAssociatedMin != -1.)&&(ptAssociatedMax != -1.)) {
781 fHistPN->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
782 fHistNP->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
783 fHistPP->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
784 fHistNN->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
785 }
786
787 //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 788
789 // Project into the wanted space (1st: analysis step, 2nd: axis)
790 TH2D* hTemp1 = (TH2D*)fHistPN->Project(0,1,2);
791 TH2D* hTemp2 = (TH2D*)fHistNP->Project(0,1,2);
792 TH2D* hTemp3 = (TH2D*)fHistPP->Project(0,1,2);
793 TH2D* hTemp4 = (TH2D*)fHistNN->Project(0,1,2);
794 TH1D* hTemp5 = (TH1D*)fHistP->Project(0,1);
795 TH1D* hTemp6 = (TH1D*)fHistN->Project(0,1);
796
797 TH2D *gHistBalanceFunctionHistogram = 0x0;
798 if((hTemp1)&&(hTemp2)&&(hTemp3)&&(hTemp4)&&(hTemp5)&&(hTemp6)) {
799 gHistBalanceFunctionHistogram = (TH2D*)hTemp1->Clone();
800 gHistBalanceFunctionHistogram->Reset();
c8ad9635 801 gHistBalanceFunctionHistogram->GetXaxis()->SetTitle("#Delta#eta");
802 gHistBalanceFunctionHistogram->GetYaxis()->SetTitle("#Delta#varphi (rad)");
803 gHistBalanceFunctionHistogram->GetZaxis()->SetTitle("B(#Delta#eta,#Delta#varphi)");
621329de 804
805 hTemp1->Sumw2();
806 hTemp2->Sumw2();
807 hTemp3->Sumw2();
808 hTemp4->Sumw2();
809 hTemp1->Add(hTemp3,-1.);
810 hTemp1->Scale(1./hTemp5->GetEntries());
811 hTemp2->Add(hTemp4,-1.);
812 hTemp2->Scale(1./hTemp6->GetEntries());
813 gHistBalanceFunctionHistogram->Add(hTemp1,hTemp2,1.,1.);
814 gHistBalanceFunctionHistogram->Scale(0.5);
815 }
816
817 return gHistBalanceFunctionHistogram;
818}
819
f0c5040c 820//____________________________________________________________________//
821TH2D *AliBalancePsi::GetBalanceFunctionDeltaEtaDeltaPhi2pMethod(Double_t psiMin,
822 Double_t psiMax,
823 Double_t ptTriggerMin,
824 Double_t ptTriggerMax,
825 Double_t ptAssociatedMin,
826 Double_t ptAssociatedMax,
827 AliBalancePsi *bfMix) {
828 //Returns the BF histogram in Delta eta vs Delta phi,
829 //after dividing each correlation function by the Event Mixing one
830 //extracted from the 6 AliTHn objects
831 //(private members) of the AliBalancePsi class.
832 //iVariableSingle: 0(phi-Psi), 1(pt-trigger)
833 //iVariablePair: 0(phi-Psi) 1(Delta eta), 2(Delta phi), 3(pt-trigger), 4(pt-associated
834 TString histName = "gHistBalanceFunctionHistogram2D";
835
836 if(!bfMix){
837 AliError("balance function object for event mixing not available");
838 return NULL;
839 }
840
841 // Psi_2
842 fHistP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
843 fHistN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
844 fHistPN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
845 fHistNP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
846 fHistPP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
847 fHistNN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
848
849 // pt trigger
850 if((ptTriggerMin != -1.)&&(ptTriggerMax != -1.)) {
851 fHistP->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
852 fHistN->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
853 fHistPN->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
854 fHistNP->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
855 fHistPP->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
856 fHistNN->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
857 }
858
859 // pt associated
860 if((ptAssociatedMin != -1.)&&(ptAssociatedMax != -1.)) {
861 fHistPN->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
862 fHistNP->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
863 fHistPP->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
864 fHistNN->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
865 }
866
867
868 // ============================================================================================
869 // the same for event mixing
870 AliTHn *fHistPMix = bfMix->GetHistNp();
871 AliTHn *fHistNMix = bfMix->GetHistNn();
872 AliTHn *fHistPNMix = bfMix->GetHistNpn();
873 AliTHn *fHistNPMix = bfMix->GetHistNnp();
874 AliTHn *fHistPPMix = bfMix->GetHistNpp();
875 AliTHn *fHistNNMix = bfMix->GetHistNnn();
876
877 // Psi_2
878 fHistPMix->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
879 fHistNMix->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
880 fHistPNMix->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
881 fHistNPMix->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
882 fHistPPMix->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
883 fHistNNMix->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
884
885 // pt trigger
886 if((ptTriggerMin != -1.)&&(ptTriggerMax != -1.)) {
887 fHistPMix->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
888 fHistNMix->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
889 fHistPNMix->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
890 fHistNPMix->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
891 fHistPPMix->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
892 fHistNNMix->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
893 }
894
895 // pt associated
896 if((ptAssociatedMin != -1.)&&(ptAssociatedMax != -1.)) {
897 fHistPNMix->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
898 fHistNPMix->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
899 fHistPPMix->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
900 fHistNNMix->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
901 }
902 // ============================================================================================
903
904
905 //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)));
906
907 // Project into the wanted space (1st: analysis step, 2nd: axis)
908 TH2D* hTemp1 = (TH2D*)fHistPN->Project(0,1,2);
909 TH2D* hTemp2 = (TH2D*)fHistNP->Project(0,1,2);
910 TH2D* hTemp3 = (TH2D*)fHistPP->Project(0,1,2);
911 TH2D* hTemp4 = (TH2D*)fHistNN->Project(0,1,2);
912 TH1D* hTemp5 = (TH1D*)fHistP->Project(0,1);
913 TH1D* hTemp6 = (TH1D*)fHistN->Project(0,1);
914
915 // ============================================================================================
916 // the same for event mixing
917 TH2D* hTemp1Mix = (TH2D*)fHistPNMix->Project(0,1,2);
918 TH2D* hTemp2Mix = (TH2D*)fHistNPMix->Project(0,1,2);
919 TH2D* hTemp3Mix = (TH2D*)fHistPPMix->Project(0,1,2);
920 TH2D* hTemp4Mix = (TH2D*)fHistNNMix->Project(0,1,2);
f0fb4ac1 921 TH1D* hTemp5Mix = (TH1D*)fHistPMix->Project(0,1);
922 TH1D* hTemp6Mix = (TH1D*)fHistNMix->Project(0,1);
f0c5040c 923 // ============================================================================================
924
925 TH2D *gHistBalanceFunctionHistogram = 0x0;
926 if((hTemp1)&&(hTemp2)&&(hTemp3)&&(hTemp4)&&(hTemp5)&&(hTemp6)) {
927 gHistBalanceFunctionHistogram = (TH2D*)hTemp1->Clone();
928 gHistBalanceFunctionHistogram->Reset();
929 gHistBalanceFunctionHistogram->GetXaxis()->SetTitle("#Delta#eta");
930 gHistBalanceFunctionHistogram->GetYaxis()->SetTitle("#Delta#varphi (rad)");
931 gHistBalanceFunctionHistogram->GetZaxis()->SetTitle("B(#Delta#eta,#Delta#varphi)");
932
933 hTemp1->Sumw2();
934 hTemp2->Sumw2();
935 hTemp3->Sumw2();
936 hTemp4->Sumw2();
937 hTemp1Mix->Sumw2();
938 hTemp2Mix->Sumw2();
939 hTemp3Mix->Sumw2();
940 hTemp4Mix->Sumw2();
941
942 hTemp1->Scale(1./hTemp5->GetEntries());
943 hTemp3->Scale(1./hTemp5->GetEntries());
944 hTemp2->Scale(1./hTemp6->GetEntries());
945 hTemp4->Scale(1./hTemp6->GetEntries());
f0fb4ac1 946 hTemp1Mix->Scale(1./hTemp5Mix->GetEntries());
947 hTemp3Mix->Scale(1./hTemp5Mix->GetEntries());
948 hTemp2Mix->Scale(1./hTemp6Mix->GetEntries());
949 hTemp4Mix->Scale(1./hTemp6Mix->GetEntries());
f0c5040c 950
951 hTemp1->Divide(hTemp1Mix);
952 hTemp2->Divide(hTemp2Mix);
953 hTemp3->Divide(hTemp3Mix);
954 hTemp4->Divide(hTemp4Mix);
955
956 hTemp1->Add(hTemp3,-1.);
957 hTemp2->Add(hTemp4,-1.);
958
959 gHistBalanceFunctionHistogram->Add(hTemp1,hTemp2,1.,1.);
960 gHistBalanceFunctionHistogram->Scale(0.5);
961 }
962
963 return gHistBalanceFunctionHistogram;
964}
965
966
621329de 967//____________________________________________________________________//
968TH2D *AliBalancePsi::GetCorrelationFunctionPN(Double_t psiMin,
6acdbcb2 969 Double_t psiMax,
970 Double_t ptTriggerMin,
971 Double_t ptTriggerMax,
972 Double_t ptAssociatedMin,
47e040b5 973 Double_t ptAssociatedMax) {
a38fd7f3 974 //Returns the 2D correlation function for (+-) pairs
621329de 975 // Psi_2: axis 0
976 fHistP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
977 fHistPN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
621329de 978
6acdbcb2 979 // pt trigger
980 if((ptTriggerMin != -1.)&&(ptTriggerMax != -1.)) {
981 fHistP->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
982 fHistPN->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
5ef8eaa7 983 }
621329de 984
6acdbcb2 985 // pt associated
986 if((ptAssociatedMin != -1.)&&(ptAssociatedMax != -1.))
987 fHistPN->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
988
989 //fHistP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(-0.5,2.5);
990 //fHistPN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(-0.5,2.5);
991
992 //TH2D *gHistTest = dynamic_cast<TH2D *>(fHistP->Project(0,0,1));
993 //TCanvas *c1 = new TCanvas("c1","");
994 //c1->cd();
995 //if(!gHistTest){
996 //AliError("Projection of fHistP = NULL");
997 //return gHistTest;
998 //}
999 //else{
1000 //gHistTest->DrawCopy("colz");
1001 //}
1002
1003 //0:step, 1: Delta eta, 2: Delta phi
621329de 1004 TH2D *gHist = dynamic_cast<TH2D *>(fHistPN->Project(0,1,2));
89c00c43 1005 if(!gHist){
1006 AliError("Projection of fHistPN = NULL");
1007 return gHist;
1008 }
1009
6acdbcb2 1010 //AliInfo(Form("Entries (test): %lf",(Double_t)(gHistTest->GetEntries())));
1011 //AliInfo(Form("Entries (1D): %lf",(Double_t)(fHistP->Project(0,1)->GetEntries())));
1012 //AliInfo(Form("Entries (2D): %lf",(Double_t)(fHistPN->Project(0,1,2)->GetEntries())));
621329de 1013
6acdbcb2 1014 //TCanvas *c2 = new TCanvas("c2","");
1015 //c2->cd();
1016 //fHistPN->Project(0,1,2)->DrawCopy("colz");
621329de 1017
47e040b5 1018 if((Double_t)(fHistP->Project(0,1)->GetEntries())!=0)
621329de 1019 gHist->Scale(1./(Double_t)(fHistP->Project(0,1)->GetEntries()));
a38fd7f3 1020
1021 return gHist;
1022}
1023
1024//____________________________________________________________________//
621329de 1025TH2D *AliBalancePsi::GetCorrelationFunctionNP(Double_t psiMin,
6acdbcb2 1026 Double_t psiMax,
1027 Double_t ptTriggerMin,
1028 Double_t ptTriggerMax,
1029 Double_t ptAssociatedMin,
47e040b5 1030 Double_t ptAssociatedMax) {
a38fd7f3 1031 //Returns the 2D correlation function for (+-) pairs
621329de 1032 // Psi_2: axis 0
1033 fHistN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
1034 fHistNP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
a38fd7f3 1035
6acdbcb2 1036 // pt trigger
1037 if((ptTriggerMin != -1.)&&(ptTriggerMax != -1.)) {
1038 fHistN->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1039 fHistNP->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1040 }
1041
1042 // pt associated
1043 if((ptAssociatedMin != -1.)&&(ptAssociatedMax != -1.))
1044 fHistNP->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
1045
1046 //0:step, 1: Delta eta, 2: Delta phi
621329de 1047 TH2D *gHist = dynamic_cast<TH2D *>(fHistNP->Project(0,1,2));
89c00c43 1048 if(!gHist){
1049 AliError("Projection of fHistPN = NULL");
1050 return gHist;
1051 }
1052
a38fd7f3 1053 //Printf("Entries (1D): %lf",(Double_t)(fHistN->Project(0,2)->GetEntries()));
1054 //Printf("Entries (2D): %lf",(Double_t)(fHistNP->Project(0,2,3)->GetEntries()));
47e040b5 1055 if((Double_t)(fHistN->Project(0,1)->GetEntries())!=0)
621329de 1056 gHist->Scale(1./(Double_t)(fHistN->Project(0,1)->GetEntries()));
a38fd7f3 1057
1058 return gHist;
1059}
1060
1061//____________________________________________________________________//
621329de 1062TH2D *AliBalancePsi::GetCorrelationFunctionPP(Double_t psiMin,
6acdbcb2 1063 Double_t psiMax,
1064 Double_t ptTriggerMin,
1065 Double_t ptTriggerMax,
1066 Double_t ptAssociatedMin,
47e040b5 1067 Double_t ptAssociatedMax) {
a38fd7f3 1068 //Returns the 2D correlation function for (+-) pairs
621329de 1069 // Psi_2: axis 0
1070 fHistP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
1071 fHistPP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
6acdbcb2 1072
1073 // pt trigger
1074 if((ptTriggerMin != -1.)&&(ptTriggerMax != -1.)) {
1075 fHistP->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1076 fHistPP->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1077 }
1078
1079 // pt associated
1080 if((ptAssociatedMin != -1.)&&(ptAssociatedMax != -1.))
1081 fHistPP->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 *>(fHistPP->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)(fHistP->Project(0,2)->GetEntries()));
1091 //Printf("Entries (2D): %lf",(Double_t)(fHistPP->Project(0,2,3)->GetEntries()));
47e040b5 1092 if((Double_t)(fHistP->Project(0,1)->GetEntries())!=0)
621329de 1093 gHist->Scale(1./(Double_t)(fHistP->Project(0,1)->GetEntries()));
a38fd7f3 1094
1095 return gHist;
1096}
1097
1098//____________________________________________________________________//
621329de 1099TH2D *AliBalancePsi::GetCorrelationFunctionNN(Double_t psiMin,
6acdbcb2 1100 Double_t psiMax,
1101 Double_t ptTriggerMin,
1102 Double_t ptTriggerMax,
1103 Double_t ptAssociatedMin,
47e040b5 1104 Double_t ptAssociatedMax) {
a38fd7f3 1105 //Returns the 2D correlation function for (+-) pairs
621329de 1106 // Psi_2: axis 0
1107 fHistN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
1108 fHistNN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
6acdbcb2 1109
1110 // pt trigger
1111 if((ptTriggerMin != -1.)&&(ptTriggerMax != -1.)) {
1112 fHistN->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1113 fHistNN->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1114 }
1115
1116 // pt associated
1117 if((ptAssociatedMin != -1.)&&(ptAssociatedMax != -1.))
1118 fHistNN->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
a38fd7f3 1119
6acdbcb2 1120 //0:step, 1: Delta eta, 2: Delta phi
621329de 1121 TH2D *gHist = dynamic_cast<TH2D *>(fHistNN->Project(0,1,2));
89c00c43 1122 if(!gHist){
1123 AliError("Projection of fHistPN = NULL");
1124 return gHist;
1125 }
1126
a38fd7f3 1127 //Printf("Entries (1D): %lf",(Double_t)(fHistN->Project(0,2)->GetEntries()));
1128 //Printf("Entries (2D): %lf",(Double_t)(fHistNN->Project(0,2,3)->GetEntries()));
47e040b5 1129 if((Double_t)(fHistN->Project(0,1)->GetEntries())!=0)
621329de 1130 gHist->Scale(1./(Double_t)(fHistN->Project(0,1)->GetEntries()));
a38fd7f3 1131
1132 return gHist;
1133}
621329de 1134
7b610f27 1135//____________________________________________________________________//
1136TH2D *AliBalancePsi::GetCorrelationFunctionChargeIndependent(Double_t psiMin,
1137 Double_t psiMax,
1138 Double_t ptTriggerMin,
1139 Double_t ptTriggerMax,
1140 Double_t ptAssociatedMin,
47e040b5 1141 Double_t ptAssociatedMax) {
7b610f27 1142 //Returns the 2D correlation function for the sum of all charge combination pairs
1143 // Psi_2: axis 0
1144 fHistN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
1145 fHistP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
1146 fHistNN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
1147 fHistPP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
1148 fHistNP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
1149 fHistPN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax);
1150
1151 // pt trigger
1152 if((ptTriggerMin != -1.)&&(ptTriggerMax != -1.)) {
1153 fHistN->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1154 fHistP->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1155 fHistNN->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1156 fHistPP->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1157 fHistNP->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1158 fHistPN->GetGrid(0)->GetGrid()->GetAxis(3)->SetRangeUser(ptTriggerMin,ptTriggerMax);
1159 }
1160
1161 // pt associated
1162 if((ptAssociatedMin != -1.)&&(ptAssociatedMax != -1.))
1163 fHistNN->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
1164 fHistPP->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
1165 fHistNP->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
1166 fHistPN->GetGrid(0)->GetGrid()->GetAxis(4)->SetRangeUser(ptAssociatedMin,ptAssociatedMax);
1167
1168 //0:step, 1: Delta eta, 2: Delta phi
1169 TH2D *gHistNN = dynamic_cast<TH2D *>(fHistNN->Project(0,1,2));
1170 if(!gHistNN){
1171 AliError("Projection of fHistNN = NULL");
1172 return gHistNN;
1173 }
1174 TH2D *gHistPP = dynamic_cast<TH2D *>(fHistPP->Project(0,1,2));
1175 if(!gHistPP){
1176 AliError("Projection of fHistPP = NULL");
1177 return gHistPP;
1178 }
1179 TH2D *gHistNP = dynamic_cast<TH2D *>(fHistNP->Project(0,1,2));
1180 if(!gHistNP){
1181 AliError("Projection of fHistNP = NULL");
1182 return gHistNP;
1183 }
1184 TH2D *gHistPN = dynamic_cast<TH2D *>(fHistPN->Project(0,1,2));
1185 if(!gHistPN){
1186 AliError("Projection of fHistPN = NULL");
1187 return gHistPN;
1188 }
1189
1190 // sum all 2 particle histograms
1191 gHistNN->Add(gHistPP);
1192 gHistNN->Add(gHistNP);
1193 gHistNN->Add(gHistPN);
1194
47e040b5 1195 // divide by sum of + and - triggers
1196 if((Double_t)(fHistN->Project(0,1)->GetEntries())!=0 && (Double_t)(fHistP->Project(0,1)->GetEntries())!=0)
7b610f27 1197 gHistNN->Scale(1./(Double_t)(fHistN->Project(0,1)->GetEntries() + fHistN->Project(0,1)->GetEntries()));
1198
1199 return gHistNN;
1200}
1201
73609a48 1202//____________________________________________________________________//
1203Float_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) {
1204 //
1205 // calculates dphistar
1206 //
1207 Float_t dphistar = phi1 - phi2 - charge1 * bSign * TMath::ASin(0.075 * radius / pt1) + charge2 * bSign * TMath::ASin(0.075 * radius / pt2);
1208
1209 static const Double_t kPi = TMath::Pi();
1210
1211 // circularity
1212// if (dphistar > 2 * kPi)
1213// dphistar -= 2 * kPi;
1214// if (dphistar < -2 * kPi)
1215// dphistar += 2 * kPi;
1216
1217 if (dphistar > kPi)
1218 dphistar = kPi * 2 - dphistar;
1219 if (dphistar < -kPi)
1220 dphistar = -kPi * 2 - dphistar;
1221 if (dphistar > kPi) // might look funny but is needed
1222 dphistar = kPi * 2 - dphistar;
1223
1224 return dphistar;
1225}
1226