Line data Source code
1 : #include "TauolaHEPEVTParticle.h"
2 :
3 : #include "TauolaLog.h"
4 :
5 : namespace Tauolapp
6 : {
7 :
8 : TauolaHEPEVTParticle::~TauolaHEPEVTParticle()
9 0 : {
10 : // Cleanup particles that do not belong to event
11 0 : for(unsigned int i=0;i<cache.size();i++)
12 0 : if(cache[i]->m_barcode<0)
13 0 : delete cache[i];
14 0 : }
15 :
16 0 : TauolaHEPEVTParticle::TauolaHEPEVTParticle(int pdgid, int status, double px, double py, double pz, double e, double m, int ms, int me, int ds, int de){
17 0 : m_px = px;
18 0 : m_py = py;
19 0 : m_pz = pz;
20 0 : m_e = e;
21 0 : m_generated_mass = m;
22 :
23 0 : m_pdgid = pdgid;
24 0 : m_status = status;
25 :
26 0 : m_first_mother = ms;
27 0 : m_second_mother = me;
28 0 : m_daughter_start = ds;
29 0 : m_daughter_end = de;
30 :
31 0 : m_barcode = -1;
32 0 : m_event = NULL;
33 0 : }
34 :
35 : void TauolaHEPEVTParticle::undecay(){
36 :
37 0 : Log::Info()<<"TauolaHEPEVTParticle::undecay not implemented for HEPEVT"<<endl;
38 0 : }
39 :
40 : void TauolaHEPEVTParticle::setMothers(vector<TauolaParticle*> mothers){
41 :
42 : // If this particle has not yet been added to the event record
43 : // then add it to the mothers' event record
44 0 : if(m_barcode<0 && mothers.size()>0)
45 : {
46 0 : TauolaHEPEVTEvent *evt = ((TauolaHEPEVTParticle*)mothers[0])->m_event;
47 0 : evt->addParticle(this);
48 0 : }
49 :
50 0 : if(mothers.size()>2) Log::Fatal("TauolaHEPEVTParticle::setMothers: HEPEVT does not allow more than two mothers!");
51 :
52 0 : if(mothers.size()>0) m_first_mother = mothers[0]->getBarcode();
53 0 : if(mothers.size()>1) m_second_mother = mothers[1]->getBarcode();
54 0 : }
55 :
56 : void TauolaHEPEVTParticle::setDaughters(vector<TauolaParticle*> daughters){
57 :
58 : // This particle must be inside some event record to be able to add daughters
59 0 : if(m_event==NULL) Log::Fatal("TauolaHEPEVTParticle::setDaughters: particle not inside event record.");
60 :
61 : int beg = 65535, end = -1;
62 :
63 0 : for(unsigned int i=0;i<daughters.size();i++)
64 : {
65 0 : int bc = daughters[i]->getBarcode();
66 0 : if(bc<0) Log::Fatal("TauolaHEPEVTParticle::setDaughters: all daughters has to be in event record first");
67 0 : if(bc<beg) beg = bc;
68 0 : if(bc>end) end = bc;
69 : }
70 0 : if(end == -1) beg = -1;
71 :
72 0 : m_daughter_start = beg;
73 0 : m_daughter_end = end;
74 0 : }
75 :
76 : std::vector<TauolaParticle*> TauolaHEPEVTParticle::getMothers(){
77 :
78 0 : std::vector<TauolaParticle*> mothers;
79 :
80 0 : TauolaParticle *p1 = NULL;
81 0 : TauolaParticle *p2 = NULL;
82 :
83 0 : if(m_first_mother>=0) p1 = m_event->getParticle(m_first_mother);
84 0 : if(m_second_mother>=0) p2 = m_event->getParticle(m_second_mother);
85 :
86 0 : if(p1) mothers.push_back(p1);
87 0 : if(p2) mothers.push_back(p2);
88 :
89 : return mothers;
90 0 : }
91 :
92 : // WARNING: this method also corrects daughter indices
93 : // if such were not defined
94 : std::vector<TauolaParticle*> TauolaHEPEVTParticle::getDaughters(){
95 :
96 0 : std::vector<TauolaParticle*> daughters;
97 :
98 0 : if(!m_event) return daughters;
99 :
100 : // Check if m_daughter_start and m_daughter_end are set
101 : // If not - try to get list of daughters from event
102 0 : if(m_daughter_end<0)
103 : {
104 : int min_d=65535, max_d=-1;
105 0 : for(int i=0;i<m_event->getParticleCount();i++)
106 : {
107 0 : if(m_event->getParticle(i)->isDaughterOf(this))
108 : {
109 0 : if(i<min_d) min_d = i;
110 0 : if(i>max_d) max_d = i;
111 : }
112 : }
113 0 : if(max_d>=0)
114 : {
115 0 : m_daughter_start = min_d;
116 0 : m_daughter_end = max_d;
117 0 : m_status = 2;
118 0 : }
119 0 : }
120 :
121 : // If m_daughter_end is still not set - there are no daughters
122 : // Otherwsie - get daughters
123 0 : if(m_daughter_end>=0)
124 : {
125 0 : for(int i=m_daughter_start;i<=m_daughter_end;i++)
126 : {
127 0 : TauolaParticle *p = m_event->getParticle(i);
128 0 : if(p==NULL)
129 : {
130 0 : Log::Warning()<<"TauolaHEPEVTParticle::getDaughters(): No particle with index "<<i<<endl;
131 0 : return daughters;
132 : }
133 :
134 0 : daughters.push_back(p);
135 0 : }
136 : }
137 :
138 0 : return daughters;
139 0 : }
140 :
141 : void TauolaHEPEVTParticle::checkMomentumConservation(){
142 :
143 0 : if(!m_event) return;
144 0 : if(m_daughter_end < 0) return;
145 :
146 0 : TauolaHEPEVTParticle *buf = m_event->getParticle(m_daughter_start);
147 :
148 0 : int first_mother_idx = buf->getFirstMotherIndex();
149 0 : int second_mother_idx = buf->getSecondMotherIndex();
150 :
151 : double px =0.0, py =0.0, pz =0.0, e =0.0;
152 : double px2=0.0, py2=0.0, pz2=0.0, e2=0.0;
153 :
154 0 : for(int i=m_daughter_start;i<=m_daughter_end;i++)
155 : {
156 0 : buf = m_event->getParticle(i);
157 0 : px += buf->getPx();
158 0 : py += buf->getPy();
159 0 : pz += buf->getPz();
160 0 : e += buf->getE ();
161 : }
162 :
163 0 : if(first_mother_idx>=0)
164 : {
165 0 : buf = m_event->getParticle(first_mother_idx);
166 0 : px2 += buf->getPx();
167 0 : py2 += buf->getPy();
168 0 : pz2 += buf->getPz();
169 0 : e2 += buf->getE();
170 0 : }
171 :
172 0 : if(second_mother_idx>=0)
173 : {
174 0 : buf = m_event->getParticle(second_mother_idx);
175 0 : px2 += buf->getPx();
176 0 : py2 += buf->getPy();
177 0 : pz2 += buf->getPz();
178 0 : e2 += buf->getE();
179 0 : }
180 : // 3-momentum // test HepMC style
181 0 : double dp = sqrt( (px-px2)*(px-px2) + (py-py2)*(py-py2) + (pz-pz2)*(pz-pz2) );
182 : // virtuality test as well.
183 0 : double m1 = sqrt( fabs( e*e - px*px - py*py - pz*pz ) );
184 0 : double m2 = sqrt( fabs( e2*e2 - px2*px2 - py2*py2 - pz2*pz2 ) );
185 :
186 0 : if( fabs(m1-m2) > 0.0001 || dp > 0.0001*(e+e2))
187 : {
188 0 : Log::RedirectOutput( Log::Warning()<<"Momentum not conserved in vertex: " );
189 0 : if(first_mother_idx >=0) m_event->getParticle(first_mother_idx) ->print();
190 0 : if(second_mother_idx>=0) m_event->getParticle(second_mother_idx)->print();
191 0 : for(int i=m_daughter_start;i<=m_daughter_end;i++) m_event->getParticle(i)->print();
192 0 : Log::RevertOutput();
193 0 : }
194 0 : }
195 :
196 : TauolaHEPEVTParticle * TauolaHEPEVTParticle::createNewParticle(
197 : int pdg_id, int status, double mass,
198 : double px, double py, double pz, double e){
199 :
200 : // New particles created using this method are added to cache
201 : // They will be deleted when this particle will be deleted
202 :
203 0 : cache.push_back(new TauolaHEPEVTParticle(pdg_id,status,px,py,pz,e,mass,-1,-1,-1,-1));
204 0 : return cache.back();
205 : }
206 :
207 : bool TauolaHEPEVTParticle::isDaughterOf(TauolaHEPEVTParticle *p)
208 : {
209 0 : int bc = p->getBarcode();
210 0 : if(bc==m_first_mother || bc==m_second_mother) return true;
211 :
212 0 : return false;
213 0 : }
214 :
215 : bool TauolaHEPEVTParticle::isMotherOf (TauolaHEPEVTParticle *p)
216 : {
217 0 : int bc = p->getBarcode();
218 0 : if(bc>=m_daughter_start && bc<=m_daughter_end) return true;
219 :
220 0 : return false;
221 0 : }
222 :
223 : void TauolaHEPEVTParticle::print(){
224 0 : char buf[256];
225 0 : sprintf(buf,"P: (%2i) %6i %2i | %11.4e %11.4e %11.4e %11.4e | %11.4e | M: %2i %2i | D: %2i %2i\n",
226 0 : m_barcode, m_pdgid, m_status, m_px, m_py, m_pz, m_e, m_generated_mass,
227 0 : m_first_mother, m_second_mother, m_daughter_start, m_daughter_end);
228 :
229 0 : cout<<buf;
230 0 : }
231 :
232 : /******** Getter and Setter methods: ***********************/
233 :
234 : void TauolaHEPEVTParticle::setPdgID(int pdg_id){
235 0 : m_pdgid = pdg_id;
236 0 : }
237 :
238 : void TauolaHEPEVTParticle::setStatus(int status){
239 0 : m_status = status;
240 0 : }
241 :
242 : void TauolaHEPEVTParticle::setMass(double mass){
243 0 : m_generated_mass = mass;
244 0 : }
245 :
246 : int TauolaHEPEVTParticle::getPdgID(){
247 0 : return m_pdgid;
248 : }
249 :
250 : int TauolaHEPEVTParticle::getStatus(){
251 0 : return m_status;
252 : }
253 :
254 : double TauolaHEPEVTParticle::getMass(){
255 0 : return m_generated_mass;
256 : }
257 :
258 : inline double TauolaHEPEVTParticle::getPx(){
259 0 : return m_px;
260 : }
261 :
262 : inline double TauolaHEPEVTParticle::getPy(){
263 0 : return m_py;
264 : }
265 :
266 : double TauolaHEPEVTParticle::getPz(){
267 0 : return m_pz;
268 : }
269 :
270 : double TauolaHEPEVTParticle::getE(){
271 0 : return m_e;
272 : }
273 :
274 : void TauolaHEPEVTParticle::setPx(double px){
275 0 : m_px = px;
276 0 : }
277 :
278 : void TauolaHEPEVTParticle::setPy(double py){
279 0 : m_py = py;
280 0 : }
281 :
282 :
283 : void TauolaHEPEVTParticle::setPz(double pz){
284 0 : m_pz = pz;
285 0 : }
286 :
287 : void TauolaHEPEVTParticle::setE(double e){
288 0 : m_e = e;
289 0 : }
290 :
291 : int TauolaHEPEVTParticle::getBarcode(){
292 0 : return m_barcode;
293 : }
294 :
295 : void TauolaHEPEVTParticle::setBarcode(int barcode){
296 0 : m_barcode = barcode;
297 0 : }
298 :
299 : void TauolaHEPEVTParticle::setEvent(TauolaHEPEVTEvent *event){
300 0 : m_event = event;
301 0 : }
302 :
303 : int TauolaHEPEVTParticle::getFirstMotherIndex(){
304 0 : return m_first_mother;
305 : }
306 :
307 : int TauolaHEPEVTParticle::getSecondMotherIndex(){
308 0 : return m_second_mother;
309 : }
310 :
311 : int TauolaHEPEVTParticle::getDaughterRangeStart(){
312 0 : return m_daughter_start;
313 : }
314 :
315 : int TauolaHEPEVTParticle::getDaughterRangeEnd(){
316 0 : return m_daughter_end;
317 : }
318 :
319 : } // namespace Tauolapp
|
