01-dev-2004 NvE Installation scripts modified to make use of environment variable...
[u/mrichter/AliRoot.git] / RALICE / AliInvmass.cxx
CommitLineData
4c039060 1/**************************************************************************
2 * Copyright(c) 1998-1999, 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
f531a546 16// $Id$
4c039060 17
959fbac5 18////////////////////////////////////////////////////////////////////////////////
19// Class AliInvmass
20// Construction of invariant mass and combinatorial background.
21//
22// Example :
23// ---------
24//
25// TObjArray* photons=new TObjArray(); // Array with photon tracks for pi0 rec.
26//
27// // Code to create some photon tracks from pi0 decays
28// Int_t ntracks=200;
29// for (Int_t i=0; i<ntracks; i++)
30// {
31// photons->Add(new Alitrack);
32// ...
33// ...
34// ...
35// }
36//
37// // Perform the invariant mass and comb. bkg. reconstruction
38//
39// TObjArray* allm=q.Invmass(photons,photons); // All reconstructed invariant masses
40//
41// TH1F* hall=new TH1F("hall","hall",200,0,2); // Histo with M_inv of all combinations
42//
43// Int_t nall=0;
44// if (allm) nall=allm->GetEntries();
45//
46// AliTrack* t;
47// Float_t minv;
48// for (Int_t j=0; j<nall; j++)
49// {
50// t=(AliTrack*)allm->At(j);
51// if (t)
52// {
53// minv=t->GetMass();
54// hall->Fill(minv);
55// }
56// }
57//
58// TObjArray* bkgm=q.CombBkg(photons,photons); // Reconstructed comb. background
59//
60// TH1F* hbkg=new TH1F("hbkg","hbkg",200,0,2); // Histo with M_inv. of comb. background
61//
62// Int_t nbkg=0;
63// if (bkgm) nbkg=bkgm->GetEntries();
64//
65// for (Int_t j=0; j<nbkg; j++)
66// {
67// t=(AliTrack*)bkgm->At(j);
68// if (t)
69// {
70// minv=t->GetMass();
71// hbkg->Fill(minv);
72// }
73// }
74//
75// TH1F* hsig=new TH1F("sig","sig",200,0,2); // Histo with the bkg. subtracted signal
76// hsig->Sumw2();
77// hsig->Add(hall,hbkg,1,-1);
78//
79//
80// Note : By default the storage of the reconstructed information is performed
81// in separate TObjArrays for the signal and comb. background resp.
82// In order to limit the memory usage, AliInvmass::SetStorageMode(1) may be
83// used to activate only a single TObjArray to store the reconstructed information.
84// Consequently, the following statements
85//
86// TObjArray* allm=q.Invmass(photons,photons);
87// TObjArray* bkgm=q.CombBkg(photons,photons);
88//
89// will result in the fact that after he invokation of CombBkg
90// the information of "allm" is lost due to the fact that the storage is
91// is re-used for "bkgm" in case the "single storage" option has been selected.
92// Usage of the, in that case invalid, pointer "allm" may cause your
93// program to crash.
94//
95// * Thus : In case of single storage usage, all invokations of the returned
96// array pointer have to be completed before invoking any memberfunction
97// of the same AliInvmass object again.
98//
99//
100//
101//--- Author: Nick van Eijndhoven 12-apr-1999 UU-SAP Utrecht
f531a546 102//- Modified: NvE $Date$ UU-SAP Utrecht
959fbac5 103////////////////////////////////////////////////////////////////////////////////
104
d88f97cc 105#include "AliInvmass.h"
c72198f1 106#include "Riostream.h"
d88f97cc 107
108ClassImp(AliInvmass) // Class implementation to enable ROOT I/O
109
110AliInvmass::AliInvmass()
111{
112// Creation of an AliInvmass object and initialisation of parameters
113 fPi=acos(-1.);
114 fMode=2;
115 fBkg=0;
116 fNewtheta=1;
117 fNewphi=1;
118 fMinv=0;
119 fMbkg=0;
120}
121////////////////////////////////////////////////////////////////////////////////
122AliInvmass::~AliInvmass()
123{
124// Destructor to delete dynamically allocated memory
125 if (fMinv)
126 {
d88f97cc 127 delete fMinv;
128 fMinv=0;
129 }
130
131 if (fMbkg)
132 {
d88f97cc 133 delete fMbkg;
134 fMbkg=0;
135 }
136}
137////////////////////////////////////////////////////////////////////////////////
138void AliInvmass::SetStorageMode(Int_t m)
139{
140// Set storage mode for the result arrays for inv. mass and comb. background
141 fMode=2;
142 if (m==1) fMode=1;
143}
144////////////////////////////////////////////////////////////////////////////////
145void AliInvmass::SetThetaSwitch(Int_t i)
146{
147// Enable/Disable (1/0) switching of theta angle in comb. bkg. reconstruction.
148// Default : Switching of theta is enabled.
149 fNewtheta=1;
150 if (i==0) fNewtheta=0;
151}
152////////////////////////////////////////////////////////////////////////////////
153void AliInvmass::SetPhiSwitch(Int_t i)
154{
155// Enable/Disable (1/0) switching of phi angle in comb. bkg. reconstruction.
156// Default : Switching of phi is enabled.
157 fNewphi=1;
158 if (i==0) fNewphi=0;
159}
160////////////////////////////////////////////////////////////////////////////////
261c0caf 161Int_t AliInvmass::GetStorageMode() const
d88f97cc 162{
163// Provide mode of storage for the result arrays for inv. mass and comb. background
164 return fMode;
165}
166////////////////////////////////////////////////////////////////////////////////
261c0caf 167Int_t AliInvmass::GetThetaSwitch() const
d88f97cc 168{
169// Provide the theta switching flag
170 return fNewtheta;
171}
172////////////////////////////////////////////////////////////////////////////////
261c0caf 173Int_t AliInvmass::GetPhiSwitch() const
d88f97cc 174{
175// Provide the phi switching flag
176 return fNewphi;
177}
178////////////////////////////////////////////////////////////////////////////////
179TObjArray* AliInvmass::Invmass(TObjArray* a1,TObjArray* a2)
180{
181// Perform two-particle invariant mass reconstruction
182 fBkg=0;
183 Combine(a1,a2);
184 return fMinv;
185}
186////////////////////////////////////////////////////////////////////////////////
187TObjArray* AliInvmass::CombBkg(TObjArray* a1,TObjArray* a2)
188{
189// Perform two-particle combinatorial background reconstruction
190 fBkg=1;
191 Combine(a1,a2);
192 if (fMode != 1)
193 {
194 return fMbkg;
195 }
196 else
197 {
198 return fMinv;
199 }
200}
201////////////////////////////////////////////////////////////////////////////////
202void AliInvmass::Combine(TObjArray* a1,TObjArray* a2)
203{
204// Perform two-particle invariant mass reconstruction
205
206 if ((!fBkg || fMode==1) && fMinv)
207 {
d88f97cc 208 delete fMinv;
209 fMinv=0;
210 }
211
212 if (fBkg && (fMode !=1) && fMbkg)
213 {
d88f97cc 214 delete fMbkg;
215 fMbkg=0;
216 }
217
218 Int_t isame; // Indicates whether both lists are identical
219 isame=0;
220 if (a1==a2) isame=1;
221
222 // Index i must loop over the shortest of a1 and a2
223 TObjArray* listi=a1;
224 TObjArray* listj=a2;
225 Int_t ni=a1->GetEntries();
226 Int_t nj=a2->GetEntries();
227 if (nj < ni)
228 {
229 ni=a2->GetEntries();
230 nj=a1->GetEntries();
231 listi=a2;
232 listj=a1;
233 }
234
235 AliTrack* p1=0;
236 AliTrack* p2=0;
237 AliTrack* px=0;
238 Ali4Vector ptot;
239 AliTrack* t=0;
240 Double_t v2[4],vx[4];
241 Float_t q1,q2;
242
243 Int_t jmin; // Start index for list j
244 Int_t jx; // Index for randomly picked particle for comb. bkg. reconstruction
245
246 for (Int_t i=0; i<ni; i++) // Select first a particle from list i
247 {
248 p1=(AliTrack*)listi->At(i);
249 p2=0;
250
251 if (!p1) continue;
252
253 jmin=0;
254 if (isame) jmin=i+1;
255 for (Int_t j=jmin; j<nj; j++) // Select also a particle from list j
256 {
257 p2=(AliTrack*)listj->At(j);
258 if (p1==p2) p2=0; // Don't combine particle with itself
259
260 if (!p2) continue;
261
262 p2->GetVector(v2,"sph");
263
264 // Take theta and phi from randomly chosen other list j particle for bkg. reconstr.
265 if (fBkg)
266 {
267 px=0;
268 if ((!isame && nj>1) || (isame && nj>2))
269 {
270 jx=int(fRndm.Uniform(0,float(nj)));
271 px=(AliTrack*)listj->At(jx);
272
273 while (!px || px==p2 || px==p1)
274 {
275 jx++;
276 if (jx >= nj) jx=0;
277 px=(AliTrack*)listj->At(jx);
278 }
279
280 px->GetVector(vx,"sph");
281 if (fNewtheta) v2[2]=vx[2]; // Replace the theta angle in the v2 vector
282 if (fNewphi) v2[3]=vx[3]; // Replace the phi angle in the v2 vector
283 }
284 }
285
286 if ((!fBkg && p2) || (fBkg && px))
287 {
288 // Store the data of this two-particle combination
289 ptot.SetVector(v2,"sph");
290 ptot=(Ali4Vector)(ptot+(*p1));
291 q1=p1->GetCharge();
292 q2=p2->GetCharge();
293 t=new AliTrack;
294 t->Set4Momentum(ptot);
295 t->SetCharge(q1+q2);
296 if (!fBkg || fMode==1)
297 {
6516b62d 298 if (!fMinv)
299 {
300 fMinv=new TObjArray();
301 fMinv->SetOwner();
302 }
d88f97cc 303 fMinv->Add(t);
304 }
305 else
306 {
6516b62d 307 if (!fMbkg)
308 {
309 fMbkg=new TObjArray();
310 fMbkg->SetOwner();
311 }
d88f97cc 312 fMbkg->Add(t);
313 }
314 }
315
316 } // End of second particle loop
317
318 } // End of first particle loop
319}
320////////////////////////////////////////////////////////////////////////////////