]>
Commit | Line | Data |
---|---|---|
621688e4 | 1 | /************************************************************************* |
2 | * Copyright(c) 1998-2048, ALICE Experiment at CERN, All rights reserved. * | |
3 | * * | |
4 | * Author: The ALICE Off-line Project. * | |
5 | * Contributors are mentioned in the code where appropriate. * | |
6 | * * | |
7 | * Permission to use, copy, modify and distribute this software and its * | |
8 | * documentation strictly for non-commercial purposes is hereby granted * | |
9 | * without fee, provided that the above copyright notice appears in all * | |
10 | * copies and that both the copyright notice and this permission notice * | |
11 | * appear in the supporting documentation. The authors make no claims * | |
12 | * about the suitability of this software for any purpose. It is * | |
13 | * provided "as is" without express or implied warranty. * | |
14 | **************************************************************************/ | |
15 | ||
16 | // Author: Dariusz Miskowiec <mailto:d.miskowiec@gsi.de> 2005 | |
17 | ||
18 | //============================================================================= | |
19 | // two-particle correlation analyzer | |
20 | // Loop over pairs and fill pair histograms. The first particle is always | |
21 | // taken from ev0 and the second from ev1 so for ev0=ev1 one is getting true | |
22 | // pairs, otherwise mixed ones. | |
23 | //============================================================================= | |
24 | ||
25 | #include <TROOT.h> | |
26 | #include <TMath.h> | |
27 | #include <TRandom2.h> | |
28 | #include "AliUnicorEvent.h" | |
29 | #include "AliUnicorHN.h" | |
30 | #include "AliUnicorAnalCorrel.h" | |
31 | ||
32 | ClassImp(AliUnicorAnalCorrel) | |
33 | ||
34 | //============================================================================= | |
61e4657c | 35 | AliUnicorAnalCorrel::AliUnicorAnalCorrel(const char *nam, Double_t emi, Double_t ema, |
28eee19b | 36 | Int_t pid0, Int_t pid1, AnalysisFrame frame): |
37 | AliUnicorAnal(nam), fPid0(pid0), fPid1(pid1), fMass0(0), fMass1(0), fZ0(0), fZ1(0), | |
38 | fFrame(frame), fPa() { | |
621688e4 | 39 | // constructor |
40 | // emi and ema define the rapidity range for histogram | |
41 | ||
42 | TParticlePDG *part0 = AliUnicorAnal::fgPDG.GetParticle(fPid0); | |
43 | TParticlePDG *part1 = AliUnicorAnal::fgPDG.GetParticle(fPid1); | |
44 | fMass0 = part0? part0->Mass() : 0; | |
45 | fMass1 = part1? part1->Mass() : 0; | |
28eee19b | 46 | fZ0 = part0? part0->Charge()/3.0 : 0; |
47 | fZ1 = part1? part1->Charge()/3.0 : 0; | |
621688e4 | 48 | double pi = TMath::Pi(); |
c6fc7f72 | 49 | |
50 | // correlation function | |
51 | ||
28eee19b | 52 | int nce = 6; double cebins[]={0,0.05,0.1,0.2,0.4,0.6,1.0}; // centrality bins |
53 | ||
76d78859 | 54 | //int npt = 7; double ptbins[]={0,0.1,0.2,0.3,0.4,0.5,0.7,1.0}; |
55 | //int npt = 6; double ptbins[]={0,0.1,0.25,0.35,0.55,1.0,2.0}; // like Adam, except last bin split | |
61e4657c | 56 | //int npt = 7; double ptbins[]={0,0.1,0.25,0.40,0.55,0.7,1.0,2.0}; // like Adam in Mar-2010, + first+last |
28eee19b | 57 | int npt = 10; double ptbins[]={0,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,1.0,2.0}; // for Pb |
76d78859 | 58 | |
61e4657c | 59 | double qbins[200] = {0}; |
28eee19b | 60 | // for (int i=0;i<60;i++) qbins[i]=i*0.005; |
61e4657c | 61 | // for (int i=0;i<45;i++) qbins[20+i]=0.1+i*0.02; |
28eee19b | 62 | // for (int i=0;i<30;i++) qbins[i]=i*0.010; |
63 | for (int i=0;i<60;i++) qbins[i]=i*0.005; | |
64 | for (int i=0;i<35;i++) qbins[60+i]=0.3+i*0.02; | |
65 | for (int i=0;i<20;i++) qbins[95+i]=1.0+i*0.05; | |
66 | for (int i=0;i<11;i++) qbins[115+i]=2.0+i*0.20; | |
c6fc7f72 | 67 | |
621688e4 | 68 | TAxis *ax[8]; |
28eee19b | 69 | ax[0] = new TAxis(4,-0.5,3.5);ax[0]->SetTitle("trumixrot"); |
70 | // ax[1] = new TAxis(5,0,1.0); ax[1]->SetTitle("centrality"); | |
71 | ax[1] = new TAxis(nce,cebins);ax[1]->SetTitle("centrality"); | |
76d78859 | 72 | ax[2] = new TAxis(3,emi,ema); ax[2]->SetTitle("y"); // pair y |
73 | // ax[3] = new TAxis(8,-pi,pi); ax[3]->SetTitle("phi"); // wrt event plane | |
74 | ax[3] = new TAxis(1,-pi,pi); ax[3]->SetTitle("phi"); // wrt event plane | |
75 | ax[4] = new TAxis(npt,ptbins);ax[4]->SetTitle("kt (GeV/c)"); // pair pt/2 | |
621688e4 | 76 | ax[5] = new TAxis(8,0,pi); ax[5]->SetTitle("q-theta"); |
77 | ax[6] = new TAxis(16,-pi,pi); ax[6]->SetTitle("q-phi"); | |
28eee19b | 78 | ax[7] = new TAxis(125,qbins); ax[7]->SetTitle("q (GeV/c)"); |
61e4657c | 79 | // ax[7] = new TAxis(700,0,3.5); ax[7]->SetTitle("q (GeV/c)"); |
621688e4 | 80 | AliUnicorHN *pair = new AliUnicorHN("pair",8,ax); |
621688e4 | 81 | fHistos.Add(pair); |
c6fc7f72 | 82 | |
76d78859 | 83 | // correlation function bin monitor (needed to get <kt> etc.) |
84 | ||
28eee19b | 85 | TAxis *bx[3]={0}; |
86 | //bx[0] = new TAxis(*(pair->GetAxis(1))); // wait until root bug (516-00..527-06) fixed. For now, do: | |
87 | bx[0] = new TAxis(ax[1]->GetNbins(), ax[1]->GetXmin(), ax[1]->GetXmax()); bx[0]->SetTitle("centrality"); | |
88 | bx[1] = new TAxis(10*ax[2]->GetNbins(),emi,ema); bx[1]->SetTitle("y"); // pair y | |
89 | bx[2] = new TAxis(100,0,2); bx[2]->SetTitle("kt (GeV/c)"); // pair pt/2 | |
90 | AliUnicorHN *bimo = new AliUnicorHN("bimo",3,bx); | |
91 | for (int i=0; i<8; i++) delete ax[i]; | |
92 | for (int i=0; i<3; i++) delete bx[i]; | |
76d78859 | 93 | fHistos.Add(bimo); |
94 | ||
c6fc7f72 | 95 | // two-track resolution monitoring histogram |
96 | ||
97 | ax[0] = new TAxis(3,-0.5,2.5); ax[0]->SetTitle("trumixrot"); | |
98 | ax[1] = new TAxis(2,-0.5,1.5); ax[1]->SetTitle("cut applied"); | |
76d78859 | 99 | ax[2] = new TAxis(npt,ptbins); ax[2]->SetTitle("(pair pt)/2 (GeV)"); |
28eee19b | 100 | ax[3] = new TAxis(80,-0.08,0.08); ax[3]->SetTitle("dtheta"); |
101 | ax[4] = new TAxis(80,-0.20,0.20); ax[4]->SetTitle("dphi"); | |
c6fc7f72 | 102 | AliUnicorHN *twot = new AliUnicorHN("twot",5,ax); |
103 | for (int i=0; i<5; i++) delete ax[i]; | |
104 | fHistos.Add(twot); | |
105 | ||
621688e4 | 106 | gROOT->cd(); |
107 | printf("%s object named %s created\n",ClassName(),GetName()); | |
108 | } | |
109 | //============================================================================= | |
c6fc7f72 | 110 | void AliUnicorAnalCorrel::Process(Int_t tmr, const AliUnicorEvent * const ev0, const AliUnicorEvent * const ev1, Double_t phirot) |
621688e4 | 111 | { |
112 | // process pairs from one or two (if mixing) events | |
113 | // tmr tells which histogram (bins) to fill: tru,mix,rot | |
114 | ||
115 | // Could be possibly accelerated by checking the "good particle" only once | |
116 | // and caching the result. (Maybe the optimizer does it already.) | |
117 | ||
118 | static TRandom2 ran; | |
119 | AliUnicorHN *pair = (AliUnicorHN*) fHistos.At(0); | |
76d78859 | 120 | AliUnicorHN *bimo = (AliUnicorHN*) fHistos.At(1); |
121 | AliUnicorHN *twot = (AliUnicorHN*) fHistos.At(2); | |
621688e4 | 122 | |
123 | // mixing-and-rotating-proof centrality and reaction plane angle | |
124 | // (but not rotation-proof for rotation angles much different from 0 and 180) | |
125 | // true and rotated pairs are within the triangle (j<i), mixed - all | |
126 | // thus, proper rotation is either by 180, or by 170 AND 190, etc. | |
127 | ||
128 | double cent = (ev0->Centrality()+ev1->Centrality())/2.0; | |
61e4657c | 129 | double q0x=0,q0y=0,q1x=0,q1y=0; |
621688e4 | 130 | ev0->RP(q0x,q0y); |
131 | ev1->RP(q1x,q1y); | |
132 | double rpphi = atan2(q0y+q1y,q0x+q1x); | |
133 | ||
134 | // loop over pairs | |
135 | ||
136 | for (int i=0; i<ev0->NParticles(); i++) { | |
137 | if (!ev0->ParticleGood(i,fPid0)) continue; | |
138 | for (int j=0; j<ev1->NParticles(); j++) { | |
139 | if (ev0==ev1 && j<i && fPid0==fPid1 ) continue; | |
140 | if (ev0==ev1 && j==i) continue; // beware, not even when rotated or non-identical | |
141 | if (!ev1->ParticleGood(j,fPid1)) continue; | |
621688e4 | 142 | fPa.Set0(fMass0,ev0->ParticleP(i),ev0->ParticleTheta(i),ev0->ParticlePhi(i)); |
143 | fPa.Set1(fMass1,ev1->ParticleP(j),ev1->ParticleTheta(j),ev1->ParticlePhi(j)+phirot); | |
144 | if (ev0==ev1 && fPid0==fPid1 && ran.Rndm()>=0.5) fPa.Swap(); | |
c6fc7f72 | 145 | twot->Fill((double) tmr, 0.0, fPa.Pt()/2.0, fPa.DTheta(), fPa.DPhi(),1.0); |
28eee19b | 146 | if (!ev0->PairGood(ev0->ParticleP(i),ev0->ParticleTheta(i),ev0->ParticlePhi(i),fZ0, |
147 | ev1->ParticleP(j),ev1->ParticleTheta(j),ev1->ParticlePhi(j)+phirot,fZ1)) continue; | |
c6fc7f72 | 148 | twot->Fill((double) tmr, 1.0, fPa.Pt()/2.0, fPa.DTheta(), fPa.DPhi(),1.0); |
28eee19b | 149 | fPa.CalcLAB(); // this could be organized better. AliUnicorPair should give k*? |
150 | fPa.CalcPairCM(); | |
151 | double qcm = fPa.QCM(); // momdif in pair cm - argument for Coulomb correction | |
152 | ||
621688e4 | 153 | fPa.CalcLAB(); |
61e4657c | 154 | if (fFrame == kPairFrame) fPa.CalcPairCM(); |
155 | if (fFrame == kLCMS) fPa.CalcLcmsCM(); | |
156 | if (fPa.QCM()==0) {printf("AliUnicorAnalCorrel: Q=0\n"); return;} // should not be too frequent | |
621688e4 | 157 | double phi = TVector2::Phi_mpi_pi(fPa.Phi()-rpphi); |
040da09c | 158 | double weigth = 1.0; |
28eee19b | 159 | /* |
160 | static TH2D *coul = 0; | |
161 | if (!coul) { | |
162 | TFile::Open("coulomb.root","read"); | |
163 | coul = (TH2D*) gFile->Get("co"); | |
164 | coul->SetDirectory(gROOT); | |
165 | gFile->Close(); | |
166 | } | |
167 | if (tmr==0 && fPid0==fPid1) { | |
168 | double co = 0; | |
169 | if (qcm>0.999) co = 1; | |
170 | else if (qcm>0.001) co = coul->Interpolate(7,qcm); | |
171 | weigth = 1.0-0.5+0.5*co*(1+exp(-pow(fPa.QCM()*7/0.197,2))); | |
172 | } | |
173 | */ | |
c6fc7f72 | 174 | pair->Fill((double) tmr, // 0 for tru, 1 for mix, 2 for rot |
621688e4 | 175 | cent, // centrality |
176 | fPa.Rapidity(), // pair rapidity | |
177 | phi, // pair phi wrt reaction plane | |
178 | fPa.Pt()/2.0, // half of pair pt | |
179 | fPa.QCMTheta(), // polar angle of Q | |
180 | fPa.QCMPhiOut(), // azimuthal angle of Q w.r.t. out | |
181 | fPa.QCM(), // |p2-p1| in c.m.s. | |
040da09c | 182 | weigth); // weigth |
28eee19b | 183 | if (tmr==0 && fPa.QCM()<0.2) bimo->Fill(cent, fPa.Rapidity(), fPa.Pt()/2.0, weigth); |
184 | if (tmr==0) pair->Fill((double) 3, // this is for Coulomb correction, maybe not necessary | |
185 | cent, // centrality | |
186 | fPa.Rapidity(), // pair rapidity | |
187 | phi, // pair phi wrt reaction plane | |
188 | fPa.Pt()/2.0, // half of pair pt | |
189 | fPa.QCMTheta(), // polar angle of Q | |
190 | fPa.QCMPhiOut(), // azimuthal angle of Q w.r.t. out | |
191 | fPa.QCM(), // |p2-p1| in c.m.s. | |
192 | weigth*qcm); // weigth*Q | |
621688e4 | 193 | } |
194 | } | |
195 | } | |
196 | //============================================================================= |