effc++ corrections (alberto)
[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>
27
28#include "AliGenReaderEMD.h"
29
af12eeb3 30ClassImp(AliGenReaderEMD)
f1d1affa 31
32
33 // -----------------------------------------------------------------------------------
34AliGenReaderEMD::AliGenReaderEMD()
35{
36// Default constructor
37 fStartEvent = 0;
38 fTreeNtuple = 0;
39 fIPSide = 0;
40 fPcToTrack = 0;
41}
42
43 // -----------------------------------------------------------------------------------
44AliGenReaderEMD::~AliGenReaderEMD()
45{
46 delete fTreeNtuple;
47}
48
49// -----------------------------------------------------------------------------------
50AliGenReaderEMD& AliGenReaderEMD::operator=(const AliGenReaderEMD& rhs)
51{
52// Assignment operator
53 rhs.Copy(*this);
54 return *this;
55}
56
57// -----------------------------------------------------------------------------------
58void AliGenReaderEMD::Copy(TObject&) const
59{
60 //
61 // Copy
62 //
63 Fatal("Copy","Not implemented!\n");
64}
65
66// -----------------------------------------------------------------------------------
67void AliGenReaderEMD::Init()
68{
69//
70// Reset the existing file environment and open a new root file
71
72 TFile *pFile=0;
73 if (!pFile) {
74 pFile = new TFile(fFileName);
75 pFile->cd();
76 printf("\n %s file opened\n\n", fFileName);
77 }
78 fTreeNtuple = (TTree*)gDirectory->Get("h2032");
79 fNcurrent = fStartEvent;
80
81 TTree *Ntu=fTreeNtuple;
82 //
83 // Set branch addresses
84 // **** neutrons
85 Ntu->SetBranchAddress("Nleft",&fNnLeft);
86 Ntu->SetBranchAddress("Eleft",&fEnLeft);
2524c56f 87 Ntu->SetBranchAddress("Pxl", fPxnLeft);
88 Ntu->SetBranchAddress("Pyl", fPynLeft);
89 Ntu->SetBranchAddress("Pzl", fPznLeft);
f1d1affa 90 Ntu->SetBranchAddress("Nright",&fNnRight);
91 Ntu->SetBranchAddress("Eright",&fEnRight);
2524c56f 92 Ntu->SetBranchAddress("Pxr", fPxnRight);
93 Ntu->SetBranchAddress("Pyr", fPynRight);
94 Ntu->SetBranchAddress("Pzr", fPznRight);
f1d1affa 95 // **** protons
96 Ntu->SetBranchAddress("Nleft_p",&fNpLeft);
97 Ntu->SetBranchAddress("Etaleft_p",&fEtapLeft);
2524c56f 98 Ntu->SetBranchAddress("Pxl_p", fPxpLeft);
99 Ntu->SetBranchAddress("Pyl_p", fPypLeft);
100 Ntu->SetBranchAddress("Pzl_p", fPzpLeft);
f1d1affa 101 Ntu->SetBranchAddress("Nright_p",&fNpRight);
102 Ntu->SetBranchAddress("Etaright_p",&fEtapRight);
2524c56f 103 Ntu->SetBranchAddress("Pxr_p", fPxpRight);
104 Ntu->SetBranchAddress("Pyr_p", fPypRight);
105 Ntu->SetBranchAddress("Pzr_p", fPzpRight);
f1d1affa 106}
107
108// -----------------------------------------------------------------------------------
109Int_t AliGenReaderEMD::NextEvent()
110{
111 // Read the next event
112 Int_t nTracks=0;
113
114 TFile* pFile = fTreeNtuple->GetCurrentFile();
115 pFile->cd();
116
117
118 Int_t nentries = (Int_t) fTreeNtuple->GetEntries();
119 if(fNcurrent < nentries) {
120 fTreeNtuple->GetEvent(fNcurrent);
121 printf("\n *** Event %d read ***\n\n",fNcurrent);
122 fNcurrent++;
123 //printf("\n \t \t LEFT-> %d neutrons and %d protons emitted\n", fNnLeft, fNpLeft);
124 //printf("\t \t RIGHT-> %d neutrons and %d protons emitted\n\n", fNnRight, fNpRight);
125 //
126 // **** IPSIde =0->RIGHT side, =1->LEFT side
127 // #### fPcToTrack =0->neutrons, =1->protons
128 if(fIPSide==0){
129 if(fPcToTrack==0){
130 printf("\n \t \t Tracking %d neutrons emitted on RIGHT side\n\n", fNnRight);
131 nTracks = fNnRight;
132 }
133 else if(fPcToTrack==1){
134 printf("\n \t \t Tracking %d protons emitted on RIGHT side\n\n", fNpRight);
135 nTracks = fNpRight;
136 }
137 }
138 else if(fIPSide==1){
139 if(fPcToTrack==0){
140 printf("\n \t \t Tracking %d neutrons emitted on LEFT side\n", fNnLeft);
141 nTracks = fNnLeft;
142 }
143 else if(fPcToTrack==1){
144 printf("\n \t \t Tracking %d protons emitted on LEFT side\n", fNpLeft);
145 nTracks = fNpLeft;
146 }
147 }
148 fNparticle = 0;
149 return nTracks;
150 }
151
152 return 0;
153}
154
155// -----------------------------------------------------------------------------------
156TParticle* AliGenReaderEMD::NextParticle()
157{
158 // Read the next particle
159 Float_t p[4];
160 Int_t ipart = kNeutron;
161 Double_t amass = TDatabasePDG::Instance()->GetParticle(kNeutron)->Mass();
162 p[0] = fPxnRight[fNparticle];
163 p[1] = fPynRight[fNparticle];
164 p[2] = fPznRight[fNparticle];
165 Float_t ptot = TMath::Sqrt(p[0]*p[0]+p[1]*p[1]+p[2]*p[2]);
166 p[3] = TMath::Sqrt(ptot*ptot+amass*amass);
167 //printf("\t Nucleon momentum: px = %f, py = %f, pz = %f\n", p[0],p[1],p[2]);
168
169 if(p[3]<=amass)
170 Warning("Generate","Particle %d E = %f mass = %f \n",ipart,p[3],amass);
171
172 TParticle* particle = new TParticle(ipart, 0, -1, -1, -1, -1,
173 p[0], p[1], p[2], p[3], 0., 0., 0., 0.);
174 //fNcurrent++;
175 fNparticle++;
176 return particle;
177}