Clean-up
[u/mrichter/AliRoot.git] / EVGEN / AliGenReaderEMD.cxx
CommitLineData
f1d1affa 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
16
17// Class to read events from external (TNtupla) file
18// Events -> neutron removal by EM dissociation of Pb nuclei
19// Data from RELDIS code (by I. Pshenichov)
20
21#include <TFile.h>
22#include <TParticle.h>
23#include <TTree.h>
24#include <TVirtualMC.h>
25#include <TDatabasePDG.h>
26#include <TPDGCode.h>
f1d1affa 27#include "AliGenReaderEMD.h"
2078b0f7 28#include "AliStack.h"
29
f1d1affa 30
706938e6 31ClassImp(AliGenReaderEMD)
f1d1affa 32
1c56e311 33AliGenReaderEMD::AliGenReaderEMD():
34 fStartEvent(0),
35 fNcurrent(0),
36 fNparticle(0),
37 fTreeNtuple(0),
1c56e311 38 fPcToTrack(0),
a16b785f 39 fOffset(0),
2078b0f7 40 fNnAside(0),
41 fEnAside(0),
a16b785f 42 fnPDGCode(0),
2078b0f7 43 fNnCside(0),
44 fEnCside(0),
45 fNpAside(0),
46 fEtapAside(0),
a16b785f 47 fpPDGCode(0),
2078b0f7 48 fNpCside(0),
a16b785f 49 fEtapCside(0),
50 fNppAside(0),
51 fEtappAside(0),
52 fppPDGCode(0),
53 fNppCside(0),
54 fEtappCside(0),
55 fNpmAside(0),
56 fEtapmAside(0),
57 fpmPDGCode(0),
58 fNpmCside(0),
59 fEtapmCside(0),
60 fNp0Aside(0),
61 fEtap0Aside(0),
62 fp0PDGCode(0),
63 fNp0Cside(0),
64 fEtap0Cside(0),
65 fNetaAside(0),
66 fEtaetaAside(0),
67 fetaPDGCode(0),
68 fNetaCside(0),
69 fEtaetaCside(0),
70 fNomegaAside(0),
71 fEtaomegaAside(0),
72 fomegaPDGCode(0),
73 fNomegaCside(0),
74 fEtaomegaCside(0)
f1d1affa 75{
a16b785f 76// Std constructor
77 for(int i=0; i<70; i++){
78 fPxnAside[i] = fPynAside[i] = fPznAside[i] = 0.;
79 fPxnCside[i] = fPynCside[i] = fPznCside[i] = 0.;
80 if(i<50){
81 fPxpAside[i] = fPypAside[i] = fPzpAside[i] = 0.;
82 fPxpCside[i] = fPypCside[i] = fPzpCside[i] = 0.;
83 if(i<30){
84 fPxppAside[i] = fPyppAside[i] = fPzppAside[i] = 0.;
85 fPxppCside[i] = fPyppCside[i] = fPzppCside[i] = 0.;
86 fPxpmAside[i] = fPypmAside[i] = fPzpmAside[i] = 0.;
87 fPxpmCside[i] = fPypmCside[i] = fPzpmCside[i] = 0.;
88 fPxp0Aside[i] = fPyp0Aside[i] = fPzp0Aside[i] = 0.;
89 fPxp0Cside[i] = fPyp0Cside[i] = fPzp0Cside[i] = 0.;
90 if(i<15){
91 fPxetaAside[i] = fPyetaAside[i] = fPzetaAside[i] = 0.;
92 fPxetaCside[i] = fPyetaCside[i] = fPzetaCside[i] = 0.;
93 fPxomegaAside[i] = fPyomegaAside[i] = fPzomegaAside[i] = 0.;
94 fPxomegaCside[i] = fPyomegaCside[i] = fPzomegaCside[i] = 0.;
95 }
96 }
97 }
98 }
99 if(fPcToTrack==kAll) printf("\n\t *** AliGenReaderEMD will track all produced particles \n\n");
100 else if(fPcToTrack==kNotNucleons) printf("\n\t *** AliGenReaderEMD will track all produced particles except nucleons\n\n");
101 else if(fPcToTrack==kOnlyNucleons) printf("\n\t *** AliGenReaderEMD will track only nucleons\n\n");
f1d1affa 102}
103
a16b785f 104
1c56e311 105AliGenReaderEMD::AliGenReaderEMD(const AliGenReaderEMD &reader):
106 AliGenReader(reader),
107 fStartEvent(0),
108 fNcurrent(0),
109 fNparticle(0),
110 fTreeNtuple(0),
1c56e311 111 fPcToTrack(0),
a16b785f 112 fOffset(0),
2078b0f7 113 fNnAside(0),
114 fEnAside(0),
a16b785f 115 fnPDGCode(0),
2078b0f7 116 fNnCside(0),
117 fEnCside(0),
118 fNpAside(0),
119 fEtapAside(0),
a16b785f 120 fpPDGCode(0),
2078b0f7 121 fNpCside(0),
a16b785f 122 fEtapCside(0),
123 fNppAside(0),
124 fEtappAside(0),
125 fppPDGCode(0),
126 fNppCside(0),
127 fEtappCside(0),
128 fNpmAside(0),
129 fEtapmAside(0),
130 fpmPDGCode(0),
131 fNpmCside(0),
132 fEtapmCside(0),
133 fNp0Aside(0),
134 fEtap0Aside(0),
135 fp0PDGCode(0),
136 fNp0Cside(0),
137 fEtap0Cside(0),
138 fNetaAside(0),
139 fEtaetaAside(0),
140 fetaPDGCode(0),
141 fNetaCside(0),
142 fEtaetaCside(0),
143 fNomegaAside(0),
144 fEtaomegaAside(0),
145 fomegaPDGCode(0),
146 fNomegaCside(0),
147 fEtaomegaCside(0)
1c56e311 148{
149 // Copy Constructor
150 reader.Copy(*this);
151}
f1d1affa 152 // -----------------------------------------------------------------------------------
153AliGenReaderEMD::~AliGenReaderEMD()
154{
155 delete fTreeNtuple;
156}
157
158// -----------------------------------------------------------------------------------
159AliGenReaderEMD& AliGenReaderEMD::operator=(const AliGenReaderEMD& rhs)
160{
161// Assignment operator
162 rhs.Copy(*this);
163 return *this;
164}
165
166// -----------------------------------------------------------------------------------
167void AliGenReaderEMD::Copy(TObject&) const
168{
169 //
170 // Copy
171 //
172 Fatal("Copy","Not implemented!\n");
173}
174
175// -----------------------------------------------------------------------------------
176void AliGenReaderEMD::Init()
177{
178//
179// Reset the existing file environment and open a new root file
180
181 TFile *pFile=0;
182 if (!pFile) {
183 pFile = new TFile(fFileName);
184 pFile->cd();
a16b785f 185 printf("\n %s file opened to read RELDIS EMD events\n\n", fFileName);
f1d1affa 186 }
187 fTreeNtuple = (TTree*)gDirectory->Get("h2032");
188 fNcurrent = fStartEvent;
189
190 TTree *Ntu=fTreeNtuple;
191 //
192 // Set branch addresses
193 // **** neutrons
2078b0f7 194 Ntu->SetBranchAddress("Nleft",&fNnAside);
195 Ntu->SetBranchAddress("Eleft",&fEnAside);
a16b785f 196 Ntu->SetBranchAddress("Ipdg_l_n",&fnPDGCode);
2078b0f7 197 Ntu->SetBranchAddress("Pxl", fPxnAside);
198 Ntu->SetBranchAddress("Pyl", fPynAside);
199 Ntu->SetBranchAddress("Pzl", fPznAside);
200 Ntu->SetBranchAddress("Nright",&fNnCside);
201 Ntu->SetBranchAddress("Eright",&fEnCside);
202 Ntu->SetBranchAddress("Pxr", fPxnCside);
203 Ntu->SetBranchAddress("Pyr", fPynCside);
204 Ntu->SetBranchAddress("Pzr", fPznCside);
f1d1affa 205 // **** protons
2078b0f7 206 Ntu->SetBranchAddress("Nleft_p",&fNpAside);
207 Ntu->SetBranchAddress("Etaleft_p",&fEtapAside);
a16b785f 208 Ntu->SetBranchAddress("Ipdg_l_p",&fpPDGCode);
2078b0f7 209 Ntu->SetBranchAddress("Pxl_p", fPxpAside);
210 Ntu->SetBranchAddress("Pyl_p", fPypAside);
211 Ntu->SetBranchAddress("Pzl_p", fPzpAside);
212 Ntu->SetBranchAddress("Nright_p",&fNpCside);
213 Ntu->SetBranchAddress("Etaright_p",&fEtapCside);
214 Ntu->SetBranchAddress("Pxr_p", fPxpCside);
215 Ntu->SetBranchAddress("Pyr_p", fPypCside);
216 Ntu->SetBranchAddress("Pzr_p", fPzpCside);
a16b785f 217 // **** pi+
218 Ntu->SetBranchAddress("Nleft_pp",&fNppAside);
219 Ntu->SetBranchAddress("Etaleft_pp",&fEtappAside);
220 Ntu->SetBranchAddress("Ipdg_l_pp",&fppPDGCode);
221 Ntu->SetBranchAddress("Pxl_pp", fPxppAside);
222 Ntu->SetBranchAddress("Pyl_pp", fPyppAside);
223 Ntu->SetBranchAddress("Pzl_pp", fPzppAside);
224 Ntu->SetBranchAddress("Nright_pp",&fNppCside);
225 Ntu->SetBranchAddress("Etaright_pp",&fEtappCside);
226 Ntu->SetBranchAddress("Pxr_pp", fPxppCside);
227 Ntu->SetBranchAddress("Pyr_pp", fPyppCside);
228 Ntu->SetBranchAddress("Pzr_pp", fPzppCside);
229 // **** pi-
230 Ntu->SetBranchAddress("Nleft_pm",&fNpmAside);
231 Ntu->SetBranchAddress("Etaleft_pm",&fEtapmAside);
232 Ntu->SetBranchAddress("Ipdg_l_pm",&fpmPDGCode);
233 Ntu->SetBranchAddress("Pxl_pm", fPxpmAside);
234 Ntu->SetBranchAddress("Pyl_pm", fPypmAside);
235 Ntu->SetBranchAddress("Pzl_pm", fPzpmAside);
236 Ntu->SetBranchAddress("Nright_pm",&fNpmCside);
237 Ntu->SetBranchAddress("Etaright_pm",&fEtapmCside);
238 Ntu->SetBranchAddress("Pxr_pm", fPxpmCside);
239 Ntu->SetBranchAddress("Pyr_pm", fPypmCside);
240 Ntu->SetBranchAddress("Pzr_pm", fPzpmCside);
241 // **** pi0
242 Ntu->SetBranchAddress("Nleft_p0",&fNp0Aside);
243 Ntu->SetBranchAddress("Etaleft_p0",&fEtap0Aside);
244 Ntu->SetBranchAddress("Ipdg_l_p0",&fp0PDGCode);
245 Ntu->SetBranchAddress("Pxl_p0", fPxp0Aside);
246 Ntu->SetBranchAddress("Pyl_p0", fPyp0Aside);
247 Ntu->SetBranchAddress("Pzl_p0", fPzp0Aside);
248 Ntu->SetBranchAddress("Nright_p0",&fNp0Cside);
249 Ntu->SetBranchAddress("Etaright_p0",&fEtap0Cside);
250 Ntu->SetBranchAddress("Pxr_p0", fPxp0Cside);
251 Ntu->SetBranchAddress("Pyr_p0", fPyp0Cside);
252 Ntu->SetBranchAddress("Pzr_p0", fPzp0Cside);
253 // **** eta
254 Ntu->SetBranchAddress("Nleft_et",&fNetaAside);
255 Ntu->SetBranchAddress("Etaleft_et",&fEtaetaAside);
256 Ntu->SetBranchAddress("Ipdg_l_et",&fetaPDGCode);
257 Ntu->SetBranchAddress("Pxl_et", fPxetaAside);
258 Ntu->SetBranchAddress("Pyl_et", fPyetaAside);
259 Ntu->SetBranchAddress("Pzl_et", fPzetaAside);
260 Ntu->SetBranchAddress("Nright_et",&fNetaCside);
261 Ntu->SetBranchAddress("Etaright_et",&fEtaetaCside);
262 Ntu->SetBranchAddress("Pxr_et", fPxetaCside);
263 Ntu->SetBranchAddress("Pyr_et", fPyetaCside);
264 Ntu->SetBranchAddress("Pzr_et", fPzetaCside);
265 // **** omega
266 Ntu->SetBranchAddress("Nleft_om",&fNomegaAside);
267 Ntu->SetBranchAddress("Etaleft_om",&fEtaomegaAside);
268 Ntu->SetBranchAddress("Ipdg_l_om",&fomegaPDGCode);
269 Ntu->SetBranchAddress("Pxl_om", fPxomegaAside);
270 Ntu->SetBranchAddress("Pyl_om", fPyomegaAside);
271 Ntu->SetBranchAddress("Pzl_om", fPzomegaAside);
272 Ntu->SetBranchAddress("Nright_om",&fNomegaCside);
273 Ntu->SetBranchAddress("Etaright_om",&fEtaomegaCside);
274 Ntu->SetBranchAddress("Pxr_om", fPxomegaCside);
275 Ntu->SetBranchAddress("Pyr_om", fPyomegaCside);
276 Ntu->SetBranchAddress("Pzr_om", fPzomegaCside);
f1d1affa 277}
278
279// -----------------------------------------------------------------------------------
280Int_t AliGenReaderEMD::NextEvent()
281{
282 // Read the next event
283 Int_t nTracks=0;
a16b785f 284 fNparticle = 0; fOffset=0;
f1d1affa 285
286 TFile* pFile = fTreeNtuple->GetCurrentFile();
287 pFile->cd();
288
289
290 Int_t nentries = (Int_t) fTreeNtuple->GetEntries();
291 if(fNcurrent < nentries) {
292 fTreeNtuple->GetEvent(fNcurrent);
a16b785f 293 if(fNcurrent%100 == 0) printf("\n *** Reading event %d ***\n",fNcurrent);
f1d1affa 294 //
a16b785f 295 if(fPcToTrack==kAll || fPcToTrack==kOnlyNucleons){ // nucleons
296 nTracks = fNnCside+fNnAside+fNpCside+fNpAside;
f1d1affa 297 }
a16b785f 298 if(fPcToTrack==kAll || fPcToTrack==kNotNucleons){ //pions,eta,omega
299 nTracks += fNppCside+fNpmCside+fNppAside+fNpmAside+fNp0Aside+fNp0Cside+
300 fNetaAside+fNetaCside+fNomegaAside+fNomegaCside;
f1d1affa 301 }
2078b0f7 302 fNcurrent++;
a16b785f 303 printf("\t #### Putting %d particles in the stack\n", nTracks);
067ae68e 304 /*if(fPcToTrack==kAll || fPcToTrack==kOnlyNucleons) printf("\t\t %d+%d neutrons, %d+%d protons\n",
a16b785f 305 fNnAside,fNnCside, fNpAside,fNpCside);
306 if(fPcToTrack==kAll || fPcToTrack==kNotNucleons) printf("\t %d+%d pi+, %d+%d pi-, %d+%d pi0, %d+%d eta, %d+%d omega\n",
307 fNppAside,fNppCside,fNpmAside,fNpmCside,
067ae68e 308 fNp0Aside,fNp0Cside,fNetaAside,fNetaCside, fNomegaAside,fNomegaCside);*/
f1d1affa 309 return nTracks;
310 }
311
312 return 0;
313}
314
315// -----------------------------------------------------------------------------------
316TParticle* AliGenReaderEMD::NextParticle()
317{
318 // Read the next particle
2078b0f7 319 Float_t p[4]={0.,0.,0.,0.};
a16b785f 320 int pdgCode=0;
2078b0f7 321
a16b785f 322 if(fPcToTrack==kAll || fPcToTrack==kOnlyNucleons){//***********************************************
323 if(fNparticle<fNnAside){
324 p[0] = fPxnAside[fNparticle];
325 p[1] = fPynAside[fNparticle];
326 p[2] = fPznAside[fNparticle];
327 pdgCode = fnPDGCode;
328// printf(" pc%d n sideA: PDG code %d, momentum (%f, %f, %f) \n", fNparticle, pdgCode, p[0],p[1],p[2]);
2078b0f7 329 }
a16b785f 330 else if(fNparticle>=fNnAside && fNparticle<(fNnAside+fNnCside)){
2078b0f7 331 p[0] = fPxnCside[fNparticle];
332 p[1] = fPynCside[fNparticle];
333 p[2] = fPznCside[fNparticle];
a16b785f 334 pdgCode = fnPDGCode;
335// printf(" pc%d n sideC: PDG code %d, momentum (%f, %f, %f) \n", fNparticle, pdgCode, p[0],p[1],p[2]);
2078b0f7 336 }
a16b785f 337 else if(fNparticle>=fNnAside+fNnCside && fNparticle<(fNnAside+fNnCside+fNpAside)){
338 p[0] = fPxpAside[fNparticle];
339 p[1] = fPypAside[fNparticle];
340 p[2] = fPzpAside[fNparticle];
341 pdgCode = fpPDGCode;
342// printf(" pc%d p sideA: PDG code %d, momentum (%f, %f, %f) \n", fNparticle, pdgCode, p[0],p[1],p[2]);
2078b0f7 343 }
a16b785f 344 else if(fNparticle>=fNnAside+fNnCside+fNpAside && fNparticle<(fNnAside+fNnCside+fNpCside+fNpAside)){
2078b0f7 345 p[0] = fPxpCside[fNparticle];
346 p[1] = fPypCside[fNparticle];
347 p[2] = fPzpCside[fNparticle];
a16b785f 348 pdgCode = fpPDGCode;
349// printf(" pc%d p sideC: PDG code %d, momentum (%f, %f, %f) \n", fNparticle, pdgCode, p[0],p[1],p[2]);
2078b0f7 350 }
a16b785f 351 fOffset = fNnAside+fNnCside+fNpCside+fNpAside;
352 } //**********************************************************************************************
353 if(fPcToTrack==kAll || fPcToTrack==kNotNucleons){
354 if(fNparticle>=fOffset && fNparticle<fOffset+fNppAside){ // *** pi +
355 p[0] = fPxppAside[fNparticle];
356 p[1] = fPyppAside[fNparticle];
357 p[2] = fPzppAside[fNparticle];
358 pdgCode = fppPDGCode;
359// printf(" pc%d pi+ sideA: PDG code %d, momentum (%f, %f, %f) \n", fNparticle, pdgCode, p[0],p[1],p[2]);
360 }
361 if(fNparticle>=fOffset+fNppAside && fNparticle<fOffset+fNppAside+fNppCside){
362 p[0] = fPxppCside[fNparticle];
363 p[1] = fPyppCside[fNparticle];
364 p[2] = fPzppCside[fNparticle];
365 pdgCode = fppPDGCode;
366// printf(" pc%d pi+ sideC: PDG code %d, momentum (%f, %f, %f) \n", fNparticle, pdgCode, p[0],p[1],p[2]);
367 }
368 if(fNparticle>=fOffset+fNppAside+fNppCside && fNparticle<fOffset+fNppAside+fNppCside+fNpmAside){ // *** pi -
369 p[0] = fPxpmAside[fNparticle];
370 p[1] = fPypmAside[fNparticle];
371 p[2] = fPzpmAside[fNparticle];
372 pdgCode = fpmPDGCode;
373// printf(" pc%d pi- sideA: PDG code %d, momentum (%f, %f, %f) \n", fNparticle, pdgCode, p[0],p[1],p[2]);
374 }
375 if(fNparticle>=fOffset+fNppAside+fNppCside+fNpmAside && fNparticle<fOffset+fNppAside+fNppCside+fNpmAside+fNpmCside){
376 p[0] = fPxpmCside[fNparticle];
377 p[1] = fPypmCside[fNparticle];
378 p[2] = fPzpmCside[fNparticle];
379 pdgCode = fpmPDGCode;
380// printf(" pc%d pi- sideC: PDG code %d, momentum (%f, %f, %f) \n", fNparticle, pdgCode, p[0],p[1],p[2]);
381 }
382 if(fNparticle>=fOffset+fNppAside+fNppCside+fNpmAside+fNpmCside &&
383 fNparticle<fOffset+fNppAside+fNppCside+fNpmAside+fNpmCside+fNp0Aside){ // *** pi 0
384 p[0] = fPxp0Aside[fNparticle];
385 p[1] = fPyp0Aside[fNparticle];
386 p[2] = fPzp0Aside[fNparticle];
387 pdgCode = fp0PDGCode;
388// printf(" pc%d pi0 sideA: PDG code %d, momentum (%f, %f, %f) \n", fNparticle, pdgCode, p[0],p[1],p[2]);
389 }
390 if(fNparticle>=fOffset+fNppAside+fNppCside+fNpmAside+fNpmCside+fNp0Aside &&
391 fNparticle<fOffset+fNppAside+fNppCside+fNpmAside+fNpmCside+fNp0Aside+fNp0Cside){
392 p[0] = fPxp0Cside[fNparticle];
393 p[1] = fPyp0Cside[fNparticle];
394 p[2] = fPzp0Cside[fNparticle];
395 pdgCode = fp0PDGCode;
396// printf(" pc%d pi0 sideC: PDG code %d, momentum (%f, %f, %f) \n", fNparticle, pdgCode, p[0],p[1],p[2]);
2078b0f7 397 }
a16b785f 398 if(fNparticle>=fOffset+fNppAside+fNppCside+fNpmAside+fNpmCside+fNp0Aside+fNp0Cside &&
399 fNparticle<fOffset+fNppAside+fNppCside+fNpmAside+fNpmCside+fNp0Aside+fNp0Cside+fNetaAside){ // *** eta
400 p[0] = fPxetaAside[fNparticle];
401 p[1] = fPyetaAside[fNparticle];
402 p[2] = fPzetaAside[fNparticle];
403 pdgCode = fetaPDGCode;
404// printf(" pc%d eta sideA: PDG code %d, momentum (%f, %f, %f) \n", fNparticle, pdgCode, p[0],p[1],p[2]);
405 }
406 if(fNparticle>=fOffset+fNppAside+fNppCside+fNpmAside+fNpmCside+fNp0Aside+fNp0Cside+fNetaAside &&
407 fNparticle<fOffset+fNppAside+fNppCside+fNpmAside+fNpmCside+fNp0Aside+fNp0Cside+fNetaAside+fNetaCside){
408 p[0] = fPxetaCside[fNparticle];
409 p[1] = fPyetaCside[fNparticle];
410 p[2] = fPzetaCside[fNparticle];
411 pdgCode = fetaPDGCode;
412// printf(" pc%d eta sideC: PDG code %d, momentum (%f, %f, %f) \n", fNparticle, pdgCode, p[0],p[1],p[2]);
413 }
414 if(fNparticle>=fOffset+fNppAside+fNppCside+fNpmAside+fNpmCside+fNp0Aside+fNp0Cside+fNetaAside+fNetaCside &&
415 fNparticle<fOffset+fNppAside+fNppCside+fNpmAside+fNpmCside+fNp0Aside+fNp0Cside+fNetaAside+fNetaCside+fNomegaAside){ // *** omega
416 p[0] = fPxomegaAside[fNparticle];
417 p[1] = fPyomegaAside[fNparticle];
418 p[2] = fPzomegaAside[fNparticle];
419 pdgCode = fomegaPDGCode;
420// printf(" pc%d omega sideA: PDG code %d, momentum (%f, %f, %f) \n", fNparticle, pdgCode, p[0],p[1],p[2]);
421 }
422 if(fNparticle>=fOffset+fNppAside+fNppCside+fNpmAside+fNpmCside+fNp0Aside+fNp0Cside+fNetaAside+fNetaCside+fNomegaAside
423 && fNparticle<fOffset+fNppAside+fNppCside+fNpmAside+fNpmCside+fNp0Aside+fNp0Cside+fNetaAside+fNetaCside+fNomegaAside+fNomegaCside){
424 p[0] = fPxomegaCside[fNparticle];
425 p[1] = fPyomegaCside[fNparticle];
426 p[2] = fPzomegaCside[fNparticle];
427 pdgCode = fomegaPDGCode;
428// printf(" pc%d omega sideC: PDG code %d, momentum (%f, %f, %f) \n", fNparticle, pdgCode, p[0],p[1],p[2]);
429 }
430
2078b0f7 431 }
a16b785f 432
f1d1affa 433 Float_t ptot = TMath::Sqrt(p[0]*p[0]+p[1]*p[1]+p[2]*p[2]);
a16b785f 434 Double_t amass = TDatabasePDG::Instance()->GetParticle(pdgCode)->Mass();
f1d1affa 435 p[3] = TMath::Sqrt(ptot*ptot+amass*amass);
f1d1affa 436
067ae68e 437 if(p[3]<=amass){
438 Warning("Generate","Particle %d E = %f GeV mass = %f GeV ",pdgCode,p[3],amass);
439 }
440
441 //printf(" Pc %d: PDGcode %d p(%1.2f, %1.2f, %1.2f, %1.3f)\n",
442 // fNparticle,pdgCode,p[0], p[1], p[2], p[3]);
a16b785f 443
444 TParticle* particle = new TParticle(pdgCode, 0, -1, -1, -1, -1,
f1d1affa 445 p[0], p[1], p[2], p[3], 0., 0., 0., 0.);
067ae68e 446 if((p[0]*p[0]+p[1]*p[1]+p[2]*p[2])>1e-5) particle->SetBit(kTransportBit);
f1d1affa 447 fNparticle++;
448 return particle;
449}