]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGCF/EBYE/BalanceFunctions/AliBalancePsi.cxx
updated
[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>
25#include <TMath.h>
26#include <TAxis.h>
27#include <TH2D.h>
28#include <TH3D.h>
29#include <TLorentzVector.h>
30#include <TObjArray.h>
31#include <TGraphErrors.h>
32#include <TString.h>
33
34#include "AliVParticle.h"
35#include "AliMCParticle.h"
36#include "AliESDtrack.h"
37#include "AliAODTrack.h"
9afe3098 38#include "AliTHn.h"
0879e280 39
40#include "AliBalancePsi.h"
41
42ClassImp(AliBalancePsi)
43
44//____________________________________________________________________//
45AliBalancePsi::AliBalancePsi() :
46 TObject(),
47 bShuffle(kFALSE),
48 fAnalysisLevel("ESD"),
49 fAnalyzedEvents(0) ,
50 fCentralityId(0) ,
51 fCentStart(0.),
52 fCentStop(0.),
9afe3098 53 fHistP(0),
54 fHistN(0),
55 fHistPN(0),
56 fHistNP(0),
57 fHistPP(0),
58 fHistNN(0),
59 fPsiInterval(15.) {
0879e280 60 // Default constructor
0879e280 61}
62
0879e280 63//____________________________________________________________________//
64AliBalancePsi::AliBalancePsi(const AliBalancePsi& balance):
65 TObject(balance), bShuffle(balance.bShuffle),
66 fAnalysisLevel(balance.fAnalysisLevel),
67 fAnalyzedEvents(balance.fAnalyzedEvents),
68 fCentralityId(balance.fCentralityId),
69 fCentStart(balance.fCentStart),
70 fCentStop(balance.fCentStop),
9afe3098 71 fHistP(balance.fHistP),
72 fHistN(balance.fHistN),
73 fHistPN(balance.fHistPN),
74 fHistNP(balance.fHistNP),
75 fHistPP(balance.fHistPP),
76 fHistNN(balance.fHistNN),
77 fPsiInterval(balance.fPsiInterval) {
0879e280 78 //copy constructor
0879e280 79}
80
81//____________________________________________________________________//
82AliBalancePsi::~AliBalancePsi() {
83 // Destructor
9afe3098 84 delete fHistP;
85 delete fHistN;
86 delete fHistPN;
87 delete fHistNP;
88 delete fHistPP;
89 delete fHistNN;
0879e280 90}
91
92//____________________________________________________________________//
9afe3098 93void AliBalancePsi::InitHistograms() {
94 // single particle histograms
95 Int_t anaSteps = 1; // analysis steps
96 Int_t iBinSingle[nTrackVariablesSingle]; // binning for track variables
97 Double_t* dBinsSingle[nTrackVariablesSingle]; // bins for track variables
98 TString axisTitleSingle[nTrackVariablesSingle]; // axis titles for track variables
99
100 // two particle histograms
101 Int_t iBinPair[nTrackVariablesPair]; // binning for track variables
102 Double_t* dBinsPair[nTrackVariablesPair]; // bins for track variables
103 TString axisTitlePair[nTrackVariablesPair]; // axis titles for track variables
104
105 //centrality
106 const Int_t kNCentralityBins = 9;
107 Double_t centralityBins[kNCentralityBins+1] = {0.,5.,10.,20.,30.,40.,50.,60.,70.,80.};
108 //const Int_t kNCentralityBins = 200;
109 //Double_t centralityBins[kNCentralityBins+1];
110 //for(Int_t i = 0; i < kNCentralityBins+1; i++)
111 //centralityBins[i] = 0.0 + i * 0.5;
112 iBinSingle[0] = kNCentralityBins;
113 dBinsSingle[0] = centralityBins;
114 axisTitleSingle[0] = "Centrality percentile [%]";
115 iBinPair[0] = kNCentralityBins;
116 dBinsPair[0] = centralityBins;
117 axisTitlePair[0] = "Centrality percentile [%]";
118
119 //Psi_2
5439ee18 120 const Int_t kNPsi2Bins = 24;
9afe3098 121 Double_t psi2Bins[kNPsi2Bins+1];
122 for(Int_t i = 0; i < kNPsi2Bins+1; i++)
5439ee18 123 psi2Bins[i] = 0.0 + i * 7.5;
9afe3098 124 iBinSingle[1] = kNPsi2Bins;
125 dBinsSingle[1] = psi2Bins;
126 axisTitleSingle[1] = "#phi - #Psi_{2} (#circ)";
127 iBinPair[1] = kNPsi2Bins;
128 dBinsPair[1] = psi2Bins;
129 axisTitlePair[1] = "#phi - #Psi_{2} (#circ)";
130
131 // eta
532ff49a 132 /*const Int_t kNEtaBins = 10;
9afe3098 133 Double_t etaBins[kNEtaBins+1];
134 for(Int_t i = 0; i < kNEtaBins+1; i++)
a38fd7f3 135 etaBins[i] = -1.0 + i * 0.2;
9afe3098 136 iBinSingle[2] = kNEtaBins;
137 dBinsSingle[2] = etaBins;
532ff49a 138 axisTitleSingle[2] = "#eta"; */
9afe3098 139
140 //#eta of triggered particle
141 /*iBinPair[2] = kNEtaBins;
142 dBinsPair[2] = etaBins;
143 axisTitlePair[2] = "#eta"; */
144
145 // delta eta
532ff49a 146 const Int_t kNDeltaEtaBins = 40;
9afe3098 147 Double_t deltaEtaBins[kNDeltaEtaBins+1];
148 for(Int_t i = 0; i < kNDeltaEtaBins+1; i++)
532ff49a 149 deltaEtaBins[i] = -2.0 + i * 0.1;
9afe3098 150 iBinPair[2] = kNDeltaEtaBins;
151 dBinsPair[2] = deltaEtaBins;
152 axisTitlePair[2] = "#Delta #eta";
153
154 // y
155 /*const Int_t kNYBins = 40;
156 Double_t yBins[kNYBins+1];
157 for(Int_t i = 0; i < kNYBins+1; i++)
158 yBins[i] = -1.0 + i * 0.05;
159 iBinSingle[3] = kNYBins;
160 dBinsSingle[3] = yBins;
161 axisTitleSingle[3] = "y";
162
163 //y of triggered particle
164 iBinPair[3] = kNYBins;
165 dBinsPair[3] = yBins;
166 axisTitlePair[3] = "y"; */
167
168 // delta y
169 /*const Int_t kNDeltaYBins = 40;
170 Double_t deltaYBins[kNDeltaYBins+1];
171 for(Int_t i = 0; i < kNDeltaYBins+1; i++)
172 deltaYBins[i] = -2.0 + i * 0.1;
173 iBinPair[8] = kNDeltaYBins;
174 dBinsPair[8] = deltaYBins;
175 axisTitlePair[8] = "#Delta y"; */
176
177 // phi
532ff49a 178 /*const Int_t kNPhiBins = 18;
9afe3098 179 Double_t phiBins[kNPhiBins+1];
180 for(Int_t i = 0; i < kNPhiBins+1; i++){
a38fd7f3 181 phiBins[i] = 0.0 + i * 20.;
9afe3098 182 }
183 iBinSingle[3] = kNPhiBins;
184 dBinsSingle[3] = phiBins;
532ff49a 185 axisTitleSingle[3] = "#phi (#circ)"; */
9afe3098 186
9afe3098 187 // delta phi
188 const Int_t kNDeltaPhiBins = 72;
189 Double_t deltaPhiBins[kNDeltaPhiBins+1];
190 for(Int_t i = 0; i < kNDeltaPhiBins+1; i++){
191 deltaPhiBins[i] = -180.0 + i * 5.;
192 }
193 iBinPair[3] = kNDeltaPhiBins;
194 dBinsPair[3] = deltaPhiBins;
195 axisTitlePair[3] = "#Delta #phi (#circ)";
196
a38fd7f3 197 // pt(trigger-associated)
532ff49a 198 const Int_t kNPtBins = 10;
a38fd7f3 199 Double_t ptBins[kNPtBins+1];
200 for(Int_t i = 0; i < kNPtBins+1; i++){
532ff49a 201 ptBins[i] = 0.0 + i * 2.0;
a38fd7f3 202 }
532ff49a 203 iBinSingle[2] = kNPtBins;
204 dBinsSingle[2] = ptBins;
205 axisTitleSingle[2] = "p_{t}^{trig.} (GeV/c)";
a38fd7f3 206
207 iBinPair[4] = kNPtBins;
208 dBinsPair[4] = ptBins;
209 axisTitlePair[4] = "p_{t}^{trig.} (GeV/c)";
210
211 iBinPair[5] = kNPtBins;
212 dBinsPair[5] = ptBins;
213 axisTitlePair[5] = "p_{t}^{assoc.} (GeV/c)";
214
9afe3098 215 // qside
216 /*const Int_t kNQSideBins = 200;
217 Double_t qSideBins[kNQSideBins+1];
218 for(Int_t i = 0; i < kNQSideBins+1; i++)
219 qSideBins[i] = 0.0 + i * 0.02;
220 iBinPair[10] = kNQSideBins;
221 dBinsPair[10] = qSideBins;
222 axisTitlePair[10] = "q_{side} (GeV/c)";
223
224 // qout
225 const Int_t kNQoutBins = 200;
226 Double_t qoutBins[kNQoutBins+1];
227 for(Int_t i = 0; i < kNQoutBins+1; i++)
228 qoutBins[i] = 0.0 + i * 0.02;
229 iBinPair[11] = kNQoutBins;
230 dBinsPair[11] = qoutBins;
231 axisTitlePair[11] = "q_{out} (GeV/c)";
232
233 // qlong
234 const Int_t kNQlongBins = 200;
235 Double_t qlongBins[kNQlongBins+1];
236 for(Int_t i = 0; i < kNQlongBins+1; i++)
237 qlongBins[i] = 0.0 + i * 0.02;
238 iBinPair[12] = kNQlongBins;
239 dBinsPair[12] = qlongBins;
240 axisTitlePair[12] = "q_{long} (GeV/c)";
241
242 // qinv
243 const Int_t kNQinvBins = 200;
244 Double_t qinvBins[kNQinvBins+1];
245 for(Int_t i = 0; i < kNQinvBins+1; i++)
246 qinvBins[i] = 0.0 + i * 0.02;
247 iBinPair[13] = kNQinvBins;
248 dBinsPair[13] = qinvBins;
249 axisTitlePair[13] = "q_{inv} (GeV/c)";*/
250
251 TString histName;
252 //+ triggered particles
253 histName = "fHistP";
254 if(bShuffle) histName.Append("_shuffle");
255 if(fCentralityId) histName += fCentralityId.Data();
256 fHistP = new AliTHn(histName.Data(),histName.Data(),anaSteps,nTrackVariablesSingle,iBinSingle);
257 for (Int_t j=0; j<nTrackVariablesSingle; j++) {
258 fHistP->SetBinLimits(j, dBinsSingle[j]);
259 fHistP->SetVarTitle(j, axisTitleSingle[j]);
0879e280 260 }
9afe3098 261
262 //- triggered particles
263 histName = "fHistN";
264 if(bShuffle) histName.Append("_shuffle");
265 if(fCentralityId) histName += fCentralityId.Data();
266 fHistN = new AliTHn(histName.Data(),histName.Data(),anaSteps,nTrackVariablesSingle,iBinSingle);
267 for (Int_t j=0; j<nTrackVariablesSingle; j++) {
268 fHistN->SetBinLimits(j, dBinsSingle[j]);
269 fHistN->SetVarTitle(j, axisTitleSingle[j]);
0879e280 270 }
9afe3098 271
272 //+- pairs
273 histName = "fHistPN";
274 if(bShuffle) histName.Append("_shuffle");
275 if(fCentralityId) histName += fCentralityId.Data();
276 fHistPN = new AliTHn(histName.Data(),histName.Data(),anaSteps, nTrackVariablesPair, iBinPair);
277 for (Int_t j=0; j<nTrackVariablesPair; j++) {
278 fHistPN->SetBinLimits(j, dBinsPair[j]);
279 fHistPN->SetVarTitle(j, axisTitlePair[j]);
0879e280 280 }
0879e280 281
9afe3098 282 //-+ pairs
283 histName = "fHistNP";
284 if(bShuffle) histName.Append("_shuffle");
285 if(fCentralityId) histName += fCentralityId.Data();
286 fHistNP = new AliTHn(histName.Data(),histName.Data(),anaSteps, nTrackVariablesPair, iBinPair);
287 for (Int_t j=0; j<nTrackVariablesPair; j++) {
288 fHistNP->SetBinLimits(j, dBinsPair[j]);
289 fHistNP->SetVarTitle(j, axisTitlePair[j]);
0879e280 290 }
0879e280 291
9afe3098 292 //++ pairs
293 histName = "fHistPP";
294 if(bShuffle) histName.Append("_shuffle");
295 if(fCentralityId) histName += fCentralityId.Data();
296 fHistPP = new AliTHn(histName.Data(),histName.Data(),anaSteps, nTrackVariablesPair, iBinPair);
297 for (Int_t j=0; j<nTrackVariablesPair; j++) {
298 fHistPP->SetBinLimits(j, dBinsPair[j]);
299 fHistPP->SetVarTitle(j, axisTitlePair[j]);
300 }
301
302 //-- pairs
303 histName = "fHistNN";
304 if(bShuffle) histName.Append("_shuffle");
305 if(fCentralityId) histName += fCentralityId.Data();
306 fHistNN = new AliTHn(histName.Data(),histName.Data(),anaSteps, nTrackVariablesPair, iBinPair);
307 for (Int_t j=0; j<nTrackVariablesPair; j++) {
308 fHistNN->SetBinLimits(j, dBinsPair[j]);
309 fHistNN->SetVarTitle(j, axisTitlePair[j]);
0879e280 310 }
9afe3098 311 Printf("Finished setting up the AliTHn");
0879e280 312}
313
314//____________________________________________________________________//
315void AliBalancePsi::CalculateBalance(Float_t fCentrality,
316 Double_t gReactionPlane,
317 vector<Double_t> **chargeVector) {
318 // Calculates the balance function
319 fAnalyzedEvents++;
320 Int_t i = 0 , j = 0;
0879e280 321
322 // Initialize histograms if not done yet
9afe3098 323 if(!fHistPN){
0879e280 324 AliWarning("Histograms not yet initialized! --> Will be done now");
325 AliWarning("This works only in local mode --> Add 'gBalance->InitHistograms()' in your configBalanceFunction");
326 InitHistograms();
327 }
328
9afe3098 329 Double_t trackVarsSingle[nTrackVariablesSingle];
330 Double_t trackVarsPair[nTrackVariablesPair];
331
0879e280 332 Int_t gNtrack = chargeVector[0]->size();
333 //Printf("(AliBalancePsi) Number of tracks: %d",gNtrack);
35cd3e21 334 Double_t gPsiMinusPhi = 0.;
0879e280 335 Double_t dy = 0., deta = 0.;
336 Double_t qLong = 0., qOut = 0., qSide = 0., qInv = 0.;
337 Double_t dphi = 0.;
338
339 Double_t charge1 = 0;
340 Double_t eta1 = 0., rap1 = 0.;
341 Double_t px1 = 0., py1 = 0., pz1 = 0.;
342 Double_t pt1 = 0.;
343 Double_t energy1 = 0.;
344 Double_t phi1 = 0.;
345
346 Double_t charge2 = 0;
347 Double_t eta2 = 0., rap2 = 0.;
348 Double_t px2 = 0., py2 = 0., pz2 = 0.;
349 Double_t pt2 = 0.;
350 Double_t energy2 = 0.;
351 Double_t phi2 = 0.;
9afe3098 352
353 for(i = 0; i < gNtrack; i++) {
354 charge1 = chargeVector[0]->at(i);
355 rap1 = chargeVector[1]->at(i);
356 eta1 = chargeVector[2]->at(i);
357 phi1 = chargeVector[3]->at(i);
358 px1 = chargeVector[4]->at(i);
359 py1 = chargeVector[5]->at(i);
360 pz1 = chargeVector[6]->at(i);
361 pt1 = chargeVector[7]->at(i);
362 energy1 = chargeVector[8]->at(i);
5439ee18 363 gPsiMinusPhi = TMath::Abs(phi1 - gReactionPlane);
364 if(gPsiMinusPhi > 180.) gPsiMinusPhi = 360. - gPsiMinusPhi;
365 //if(gPsiMinusPhi < -fPsiInterval/2) gPsiMinusPhi = 360. + gPsiMinusPhi;
9afe3098 366
367 trackVarsSingle[0] = fCentrality;
368 trackVarsSingle[1] = gPsiMinusPhi;
532ff49a 369 //trackVarsSingle[2] = eta1;
9afe3098 370 //trackVarsSingle[3] = rap1;
532ff49a 371 //trackVarsSingle[3] = phi1;
372 trackVarsSingle[2] = pt1;
9afe3098 373
374 //Printf("Track(a) %d - phi-Psi: %lf",i+1,trackVarsSingle[1]);
375 //fill single particle histograms
376 if(charge1 > 0)
377 fHistP->Fill(trackVarsSingle,0,1.);
378 else if(charge1 < 0)
379 fHistN->Fill(trackVarsSingle,0,1.);
380
381 for(j = 0; j < i; j++) {
382 charge2 = chargeVector[0]->at(j);
383 rap2 = chargeVector[1]->at(j);
384 eta2 = chargeVector[2]->at(j);
385 phi2 = chargeVector[3]->at(j);
386 px2 = chargeVector[4]->at(j);
387 py2 = chargeVector[5]->at(j);
388 pz2 = chargeVector[6]->at(j);
389 pt2 = chargeVector[7]->at(j);
390 energy2 = chargeVector[8]->at(j);
391 //Printf("Track(b) %d - pt: %lf",j+1,pt2);
392
393 // filling the arrays
394 // RAPIDITY
395 dy = rap1 - rap2;
396
397 // Eta
398 deta = eta1 - eta2;
399
400 //qlong
401 Double_t eTot = energy1 + energy2;
402 Double_t pxTot = px1 + px2;
403 Double_t pyTot = py1 + py2;
404 Double_t pzTot = pz1 + pz2;
405 Double_t q0Tot = energy1 - energy2;
406 Double_t qxTot = px1 - px2;
407 Double_t qyTot = py1 - py2;
408 Double_t qzTot = pz1 - pz2;
409
410 Double_t eTot2 = eTot*eTot;
411 Double_t pTot2 = pxTot*pxTot + pyTot*pyTot + pzTot*pzTot;
412 Double_t pzTot2 = pzTot*pzTot;
413
414 Double_t q0Tot2 = q0Tot*q0Tot;
415 Double_t qTot2 = qxTot*qxTot + qyTot*qyTot + qzTot*qzTot;
416
417 Double_t snn = eTot2 - pTot2;
418 Double_t ptTot2 = pTot2 - pzTot2 ;
419 Double_t ptTot = TMath::Sqrt( ptTot2 );
420
421 qLong = TMath::Abs(eTot*qzTot - pzTot*q0Tot)/TMath::Sqrt(snn + ptTot2);
422
423 //qout
424 qOut = TMath::Sqrt(snn/(snn + ptTot2)) * TMath::Abs(pxTot*qxTot + pyTot*qyTot)/ptTot;
425
426 //qside
427 qSide = TMath::Abs(pxTot*qyTot - pyTot*qxTot)/ptTot;
428
429 //qinv
430 qInv = TMath::Sqrt(TMath::Abs(-q0Tot2 + qTot2 ));
431
432 //phi
433 dphi = phi2 - phi1;
434 if(dphi < -180.) dphi += 360.; //dphi should be between -180 and 180!
435 else if(dphi > 180.) dphi -= 360.; //dphi should be between -180 and 180!
436
437 trackVarsPair[0] = fCentrality;
438 trackVarsPair[1] = gPsiMinusPhi;
439 //trackVarsPair[2] = eta1;
440 //trackVarsPair[3] = rap1;
441 //trackVarsPair[4] = phi1;
9afe3098 442 trackVarsPair[2] = deta;
443 //trackVarsPair[8] = dy;
444 trackVarsPair[3] = dphi;
a38fd7f3 445 trackVarsPair[4] = pt1;
446 trackVarsPair[5] = pt2;
9afe3098 447 //trackVarsPair[10] = qSide;
448 //trackVarsPair[11] = qOut;
449 //trackVarsPair[12] = qLong;
450 //trackVarsPair[13] = qInv;
451
452 if( charge1 > 0 && charge2 < 0) fHistPN->Fill(trackVarsPair,0,1.);
453 else if( charge1 < 0 && charge2 > 0) fHistNP->Fill(trackVarsPair,0,1.);
454 else if( charge1 > 0 && charge2 > 0) fHistPP->Fill(trackVarsPair,0,1.);
455 else if( charge1 < 0 && charge2 < 0) fHistNN->Fill(trackVarsPair,0,1.);
456 else AliWarning("Wrong charge combination!");
457
0879e280 458 }//end of 2nd particle loop
459 }//end of 1st particle loop
0879e280 460}
461
0879e280 462//____________________________________________________________________//
9afe3098 463TH1D *AliBalancePsi::GetBalanceFunctionHistogram(Int_t iVariableSingle,
464 Int_t iVariablePair,
465 Double_t centrMin,
466 Double_t centrMax,
467 Double_t psiMin,
468 Double_t psiMax) {
469 //Returns the BF histogram, extracted from the 6 AliTHn objects
0879e280 470 //(private members) of the AliBalancePsi class.
9afe3098 471 //iVariableSingle: 2(eta), 3(y), 4(phi)
472 //iVariablePair: 2(Delta eta), 3(Delta y), 4(Delta phi), 5(qside), 6(qout), 7(qlong) 8(qinv)
473 TString gAnalysisType[ANALYSIS_TYPES] = {"y","eta","phi","qlong","qout","qside","qinv"};
0879e280 474 TString histName = "gHistBalanceFunctionHistogram";
9afe3098 475 histName += gAnalysisType[iVariablePair];
476
477 // centrality
478 fHistP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(centrMin,centrMax);
479 fHistN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(centrMin,centrMax);
480 fHistPN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(centrMin,centrMax);
481 fHistNP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(centrMin,centrMax);
482 fHistPP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(centrMin,centrMax);
483 fHistNN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(centrMin,centrMax);
484
485 // Psi_2
486 fHistP->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(psiMin,psiMax);
487 fHistN->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(psiMin,psiMax);
488 fHistPN->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(psiMin,psiMax);
489 fHistNP->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(psiMin,psiMax);
490 fHistPP->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(psiMin,psiMax);
491 fHistNN->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(psiMin,psiMax);
492
a38fd7f3 493 //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));
494
9afe3098 495 // Project into the wanted space (1st: analysis step, 2nd: axis)
496 TH1D* hTemp1 = (TH1D*)fHistPN->Project(0,iVariablePair);
497 TH1D* hTemp2 = (TH1D*)fHistNP->Project(0,iVariablePair);
498 TH1D* hTemp3 = (TH1D*)fHistPP->Project(0,iVariablePair);
499 TH1D* hTemp4 = (TH1D*)fHistNN->Project(0,iVariablePair);
500 TH1D* hTemp5 = (TH1D*)fHistP->Project(0,iVariableSingle);
501 TH1D* hTemp6 = (TH1D*)fHistN->Project(0,iVariableSingle);
502
503 TH1D *gHistBalanceFunctionHistogram = 0x0;
504 if((hTemp1)&&(hTemp2)&&(hTemp3)&&(hTemp4)&&(hTemp5)&&(hTemp6)) {
505 gHistBalanceFunctionHistogram = (TH1D*)hTemp1->Clone();
506 gHistBalanceFunctionHistogram->Reset();
507
508 switch(iVariablePair) {
509 case 2:
510 gHistBalanceFunctionHistogram->GetXaxis()->SetTitle("#Delta #eta");
511 gHistBalanceFunctionHistogram->GetYaxis()->SetTitle("B(#Delta #eta)");
512 break;
513 case 3:
514 gHistBalanceFunctionHistogram->GetXaxis()->SetTitle("#Delta y");
515 gHistBalanceFunctionHistogram->GetYaxis()->SetTitle("B(#Delta y)");
516 break;
517 case 4:
518 gHistBalanceFunctionHistogram->GetXaxis()->SetTitle("#Delta #phi (deg.)");
519 gHistBalanceFunctionHistogram->GetYaxis()->SetTitle("B(#Delta #phi)");
520 break;
521 case 5:
522 gHistBalanceFunctionHistogram->GetXaxis()->SetTitle("q_{side} (GeV/c)");
523 gHistBalanceFunctionHistogram->GetYaxis()->SetTitle("B(q_{side})");
524 break;
525 case 6:
526 gHistBalanceFunctionHistogram->GetXaxis()->SetTitle("q_{out} (GeV/c)");
527 gHistBalanceFunctionHistogram->GetYaxis()->SetTitle("B(q_{out})");
528 break;
529 case 7:
530 gHistBalanceFunctionHistogram->GetXaxis()->SetTitle("q_{long} (GeV/c)");
531 gHistBalanceFunctionHistogram->GetYaxis()->SetTitle("B(q_{long})");
532 break;
533 case 8:
534 gHistBalanceFunctionHistogram->GetXaxis()->SetTitle("q_{inv} (GeV/c)");
535 gHistBalanceFunctionHistogram->GetYaxis()->SetTitle("B(q_{inv})");
536 break;
537 default:
538 break;
539 }
0879e280 540
0879e280 541 hTemp1->Sumw2();
542 hTemp2->Sumw2();
543 hTemp3->Sumw2();
544 hTemp4->Sumw2();
a38fd7f3 545 hTemp1->Add(hTemp3,-1.);
0879e280 546 hTemp1->Scale(1./hTemp5->GetEntries());
a38fd7f3 547 hTemp2->Add(hTemp4,-1.);
0879e280 548 hTemp2->Scale(1./hTemp6->GetEntries());
549 gHistBalanceFunctionHistogram->Add(hTemp1,hTemp2,1.,1.);
9afe3098 550 gHistBalanceFunctionHistogram->Scale(0.5);
0879e280 551 }
552
0879e280 553 return gHistBalanceFunctionHistogram;
554}
a38fd7f3 555
556//____________________________________________________________________//
557TH2D *AliBalancePsi::GetCorrelationFunctionPN(Double_t centrMin,
558 Double_t centrMax,
559 Double_t psiMin,
560 Double_t psiMax) {
561 //Returns the 2D correlation function for (+-) pairs
562 // centrality: axis 0
563 fHistP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(centrMin,centrMax);
564 fHistPN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(centrMin,centrMax);
565
566 // Psi_2: axis 1
567 fHistP->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(psiMin,psiMax);
568 fHistPN->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(psiMin,psiMax);
569
570 //0:step, 2: Delta eta, 3: Delta phi
571 TH2D *gHist = dynamic_cast<TH2D *>(fHistPN->Project(0,2,3));
89c00c43 572 if(!gHist){
573 AliError("Projection of fHistPN = NULL");
574 return gHist;
575 }
576
a38fd7f3 577 //Printf("Entries (1D): %lf",(Double_t)(fHistP->Project(0,2)->GetEntries()));
578 //Printf("Entries (2D): %lf",(Double_t)(fHistPN->Project(0,2,3)->GetEntries()));
579 if((Double_t)(fHistP->Project(0,2)->GetEntries())!=0)
580 gHist->Scale(1./(Double_t)(fHistP->Project(0,2)->GetEntries()));
581
582 return gHist;
583}
584
585//____________________________________________________________________//
586TH2D *AliBalancePsi::GetCorrelationFunctionNP(Double_t centrMin,
587 Double_t centrMax,
588 Double_t psiMin,
589 Double_t psiMax) {
590 //Returns the 2D correlation function for (+-) pairs
591 // centrality: axis 0
592 fHistN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(centrMin,centrMax);
593 fHistNP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(centrMin,centrMax);
594
595 // Psi_2: axis 1
596 fHistN->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(psiMin,psiMax);
597 fHistNP->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(psiMin,psiMax);
598
599 //0:step, 2: Delta eta, 3: Delta phi
600 TH2D *gHist = dynamic_cast<TH2D *>(fHistNP->Project(0,2,3));
89c00c43 601 if(!gHist){
602 AliError("Projection of fHistPN = NULL");
603 return gHist;
604 }
605
a38fd7f3 606 //Printf("Entries (1D): %lf",(Double_t)(fHistN->Project(0,2)->GetEntries()));
607 //Printf("Entries (2D): %lf",(Double_t)(fHistNP->Project(0,2,3)->GetEntries()));
608 if((Double_t)(fHistN->Project(0,2)->GetEntries())!=0)
609 gHist->Scale(1./(Double_t)(fHistN->Project(0,2)->GetEntries()));
610
611 return gHist;
612}
613
614//____________________________________________________________________//
615TH2D *AliBalancePsi::GetCorrelationFunctionPP(Double_t centrMin,
616 Double_t centrMax,
617 Double_t psiMin,
618 Double_t psiMax) {
619 //Returns the 2D correlation function for (+-) pairs
620 // centrality: axis 0
621 fHistP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(centrMin,centrMax);
622 fHistPP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(centrMin,centrMax);
623
624 // Psi_2: axis 1
625 fHistP->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(psiMin,psiMax);
626 fHistPP->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(psiMin,psiMax);
627
628 //0:step, 2: Delta eta, 3: Delta phi
629 TH2D *gHist = dynamic_cast<TH2D *>(fHistPP->Project(0,2,3));
89c00c43 630 if(!gHist){
631 AliError("Projection of fHistPN = NULL");
632 return gHist;
633 }
634
a38fd7f3 635 //Printf("Entries (1D): %lf",(Double_t)(fHistP->Project(0,2)->GetEntries()));
636 //Printf("Entries (2D): %lf",(Double_t)(fHistPP->Project(0,2,3)->GetEntries()));
637 if((Double_t)(fHistP->Project(0,2)->GetEntries())!=0)
638 gHist->Scale(1./(Double_t)(fHistP->Project(0,2)->GetEntries()));
639
640 return gHist;
641}
642
643//____________________________________________________________________//
644TH2D *AliBalancePsi::GetCorrelationFunctionNN(Double_t centrMin,
645 Double_t centrMax,
646 Double_t psiMin,
647 Double_t psiMax) {
648 //Returns the 2D correlation function for (+-) pairs
649 // centrality: axis 0
650 fHistN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(centrMin,centrMax);
651 fHistNN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(centrMin,centrMax);
652
653 // Psi_2: axis 1
654 fHistN->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(psiMin,psiMax);
655 fHistNN->GetGrid(0)->GetGrid()->GetAxis(1)->SetRangeUser(psiMin,psiMax);
656
657 //0:step, 2: Delta eta, 3: Delta phi
658 TH2D *gHist = dynamic_cast<TH2D *>(fHistNN->Project(0,2,3));
89c00c43 659 if(!gHist){
660 AliError("Projection of fHistPN = NULL");
661 return gHist;
662 }
663
a38fd7f3 664 //Printf("Entries (1D): %lf",(Double_t)(fHistN->Project(0,2)->GetEntries()));
665 //Printf("Entries (2D): %lf",(Double_t)(fHistNN->Project(0,2,3)->GetEntries()));
666 if((Double_t)(fHistN->Project(0,2)->GetEntries())!=0)
667 gHist->Scale(1./(Double_t)(fHistN->Project(0,2)->GetEntries()));
668
669 return gHist;
670}