[u/mrichter/AliRoot.git] / TEvtGen / Tauola / TauolaHepMCParticle.cxx

#include "TauolaHepMCParticle.h"
#include "TauolaLog.h"

namespace Tauolapp
{

TauolaHepMCParticle::TauolaHepMCParticle(){
  m_particle = new HepMC::GenParticle();
}

TauolaHepMCParticle::~TauolaHepMCParticle(){
  
  //delete the mother and daughter pointers
  while(m_mothers.size()!=0){
    TauolaParticle * temp = m_mothers.back();
    m_mothers.pop_back();
    delete temp;
  }
  while(m_daughters.size()!=0){
    TauolaParticle * temp = m_daughters.back();
    m_daughters.pop_back();
    delete temp;
  }

  while(m_created_particles.size()!=0){
    TauolaHepMCParticle * temp = (TauolaHepMCParticle*) m_created_particles.back();
    m_created_particles.pop_back();
    if(temp->getHepMC()->barcode()==0) delete temp->getHepMC();
    delete temp;
  } 

}

// NOTE: Not executed by release examples
TauolaHepMCParticle::TauolaHepMCParticle(int pdg_id, int status, double mass){
  m_particle = new HepMC::GenParticle();
  m_particle->set_pdg_id(pdg_id);
  m_particle->set_status(status);
  m_particle->set_generated_mass(mass);
}

TauolaHepMCParticle::TauolaHepMCParticle(HepMC::GenParticle * particle){
  m_particle = particle;
}

HepMC::GenParticle * TauolaHepMCParticle::getHepMC(){
  return m_particle;
}

void TauolaHepMCParticle::undecay(){
  std::vector<TauolaParticle*> daughters = getDaughters();
  std::vector<TauolaParticle*>::iterator dIter = daughters.begin();

  for(; dIter != daughters.end(); dIter++)
    (*dIter)->undecay();

  if(m_particle->end_vertex())
  {
  while(m_particle->end_vertex()->particles_out_size())
  {
    HepMC::GenParticle *p = m_particle->end_vertex()->remove_particle(*(m_particle->end_vertex()->particles_out_const_begin()));
    delete p;
  }
  delete m_particle->end_vertex();
  }

  m_daughters.clear();
  m_particle->set_status(TauolaParticle::STABLE);

  for(unsigned int i=0;i<daughters.size();i++)
    delete daughters[i];
}

void TauolaHepMCParticle::setMothers(vector<TauolaParticle*> mothers){

  /******** Deal with mothers ***********/

  //If there are mothers
  if(mothers.size()>0){

    HepMC::GenParticle * part;
    part=dynamic_cast<TauolaHepMCParticle*>(mothers.at(0))->getHepMC();

    //Use end vertex of first mother as production vertex for particle
    HepMC::GenVertex * production_vertex = part->end_vertex();
    HepMC::GenVertex * orig_production_vertex = production_vertex;

    //If production_vertex does not exist - create it
    //If it's tau decay - set the time and position including the tau lifetime correction
    //otherwise - copy the time and position of decaying particle
    if(!production_vertex){
      production_vertex = new HepMC::GenVertex();
      HepMC::FourVector point = part->production_vertex()->position();
      production_vertex->set_position(point);
      part->parent_event()->add_vertex(production_vertex);
    }

    //Loop over all mothers to check that the end points to the right place
    vector<TauolaParticle*>::iterator mother_itr;
    for(mother_itr = mothers.begin(); mother_itr != mothers.end();
        mother_itr++){

      HepMC::GenParticle * moth;
      moth = dynamic_cast<TauolaHepMCParticle*>(*mother_itr)->getHepMC();

      if(moth->end_vertex()!=orig_production_vertex)
        Log::Fatal("Mother production_vertices point to difference places. Can not override. Please delete vertices first.",1);
      else
        production_vertex->add_particle_in(moth);

      //update status info
      if(moth->status()==TauolaParticle::STABLE)
        moth->set_status(TauolaParticle::DECAYED);
    }
    production_vertex->add_particle_out(m_particle);
  }
}

void TauolaHepMCParticle::setDaughters(vector<TauolaParticle*> daughters){

  if(!m_particle->parent_event())
    Log::Fatal("New particle needs the event set before it's daughters can be added",2);

  //If there are daughters
  if(daughters.size()>0){
    // NOTE: Not executed by release examples
    //       because daughters.size() is always 0

    //Use production vertex of first daughter as end vertex for particle
    HepMC::GenParticle * first_daughter;
    first_daughter = (dynamic_cast<TauolaHepMCParticle*>(daughters.at(0)))->getHepMC();

    HepMC::GenVertex * end_vertex;
    end_vertex=first_daughter->production_vertex();
    HepMC::GenVertex * orig_end_vertex = end_vertex;

    if(!end_vertex){ //if it does not exist create it
      end_vertex = new HepMC::GenVertex();
      m_particle->parent_event()->add_vertex(end_vertex);
    }

    //Loop over all daughters to check that the end points to the right place
    vector<TauolaParticle*>::iterator daughter_itr;
    for(daughter_itr = daughters.begin(); daughter_itr != daughters.end();
        daughter_itr++){

      HepMC::GenParticle * daug;
      daug = dynamic_cast<TauolaHepMCParticle*>(*daughter_itr)->getHepMC();


      if(daug->production_vertex()!=orig_end_vertex)
        Log::Fatal("Daughter production_vertices point to difference places. Can not override. Please delete vertices first.",3);
      else
        end_vertex->add_particle_out(daug);
    }
    end_vertex->add_particle_in(m_particle);
  }
}

std::vector<TauolaParticle*> TauolaHepMCParticle::getMothers(){

  if(m_mothers.size()==0&&m_particle->production_vertex()){
    HepMC::GenVertex::particles_in_const_iterator pcle_itr;
    pcle_itr=m_particle->production_vertex()->particles_in_const_begin();
    
    HepMC::GenVertex::particles_in_const_iterator pcle_itr_end;
    pcle_itr_end=m_particle->production_vertex()->particles_in_const_end();
    
    for(;pcle_itr != pcle_itr_end; pcle_itr++){
      m_mothers.push_back(new TauolaHepMCParticle(*pcle_itr));
    }
  }
  return m_mothers;
}

std::vector<TauolaParticle*> TauolaHepMCParticle::getDaughters(){
  
  if(m_daughters.size()==0&&m_particle->end_vertex()){
    HepMC::GenVertex::particles_out_const_iterator pcle_itr;
    pcle_itr=m_particle->end_vertex()->particles_out_const_begin();
    
    HepMC::GenVertex::particles_out_const_iterator pcle_itr_end;
    pcle_itr_end=m_particle->end_vertex()->particles_out_const_end();
    
    for(;pcle_itr != pcle_itr_end; pcle_itr++){
      m_daughters.push_back(new TauolaHepMCParticle(*pcle_itr));
    }
  }
  return m_daughters;
}

void TauolaHepMCParticle::checkMomentumConservation(){

  if(!m_particle->end_vertex()) return;
  
  // HepMC version of check_momentum_conservation
  // with added energy check

  double sumpx = 0, sumpy = 0, sumpz = 0, sume = 0;
  for( HepMC::GenVertex::particles_in_const_iterator part1 = m_particle->end_vertex()->particles_in_const_begin();
       part1 != m_particle->end_vertex()->particles_in_const_end(); part1++ ){

    sumpx += (*part1)->momentum().px();
    sumpy += (*part1)->momentum().py();
    sumpz += (*part1)->momentum().pz();
    sume  += (*part1)->momentum().e();
  }
  
  for( HepMC::GenVertex::particles_out_const_iterator part2 = m_particle->end_vertex()->particles_out_const_begin();
       part2 != m_particle->end_vertex()->particles_out_const_end(); part2++ ){

    sumpx -= (*part2)->momentum().px();
    sumpy -= (*part2)->momentum().py();
    sumpz -= (*part2)->momentum().pz();
    sume  -= (*part2)->momentum().e();
  }

  if( sqrt( sumpx*sumpx + sumpy*sumpy + sumpz*sumpz + sume*sume) > Tauola::momentum_conservation_threshold ) {
    Log::Warning()<<"Momentum not conserved in the vertex:"<<endl;
    Log::RedirectOutput(Log::Warning(false));
    m_particle->end_vertex()->print();
    Log::RevertOutput();
    return;
  }
  
  return;
}

// NOTE: Not executed by release examples
void TauolaHepMCParticle::setPdgID(int pdg_id){
  m_particle->set_pdg_id(pdg_id);
}

void TauolaHepMCParticle::setMass(double mass){
  m_particle->set_generated_mass(mass);
}

// NOTE: Not executed by release examples
void TauolaHepMCParticle::setStatus(int status){
  m_particle->set_status(status);
}

int TauolaHepMCParticle::getPdgID(){
  return m_particle->pdg_id();
}

int TauolaHepMCParticle::getStatus(){
  return m_particle->status();
}

int TauolaHepMCParticle::getBarcode(){
  return m_particle->barcode();
}

// Set (X,T) Position of tau decay trees
void TauolaHepMCParticle::decayEndgame(){

  double lifetime = Tauola::tau_lifetime * (-log( Tauola::randomDouble() ));
  HepMC::FourVector tau_momentum = m_particle->momentum();

  double mass     = sqrt(abs(  tau_momentum.e()*tau_momentum.e()
                             - tau_momentum.px()*tau_momentum.px()
                             - tau_momentum.py()*tau_momentum.py()
                             - tau_momentum.pz()*tau_momentum.pz()
                            ) );

  // Get previous position
  HepMC::FourVector previous_position = m_particle->production_vertex()->position();

  // Calculate new position
  HepMC::FourVector new_position(previous_position.x()+tau_momentum.px()/mass*lifetime,
                                 previous_position.y()+tau_momentum.py()/mass*lifetime,
                                 previous_position.z()+tau_momentum.pz()/mass*lifetime,
                                 previous_position.t()+tau_momentum.e() /mass*lifetime);

  // Set new position
  m_particle->end_vertex()->set_position(new_position);
  recursiveSetPosition(m_particle,new_position);
}

void TauolaHepMCParticle::recursiveSetPosition(HepMC::GenParticle *p, HepMC::FourVector pos){

  if(!p->end_vertex()) return;

  // Iterate over all outgoing particles
  for(HepMC::GenVertex::particles_out_const_iterator pp = p->end_vertex()->particles_out_const_begin();
      pp != p->end_vertex()->particles_out_const_end();
      ++pp){
    if( !(*pp)->end_vertex() ) continue;

    // Set position
    (*pp)->end_vertex()->set_position(pos);
    recursiveSetPosition(*pp,pos);
  }
}

TauolaHepMCParticle * TauolaHepMCParticle::createNewParticle(
                        int pdg_id, int status, double mass,
                        double px, double py, double pz, double e){

  TauolaHepMCParticle * new_particle = new TauolaHepMCParticle();
  new_particle->getHepMC()->set_pdg_id(pdg_id);
  new_particle->getHepMC()->set_status(status);
  new_particle->getHepMC()->set_generated_mass(mass);

  HepMC::FourVector momentum(px,py,pz,e);
  new_particle->getHepMC()->set_momentum(momentum);

  m_created_particles.push_back(new_particle);
  
  return new_particle;
}

void TauolaHepMCParticle::print(){
  m_particle->print();
}


/******** Getter and Setter methods: ***********************/

inline double TauolaHepMCParticle::getPx(){
  return m_particle->momentum().px();
}

inline double TauolaHepMCParticle::getPy(){
  return m_particle->momentum().py();
}

double TauolaHepMCParticle::getPz(){
  return m_particle->momentum().pz();
}

double TauolaHepMCParticle::getE(){
  return m_particle->momentum().e();
}

void TauolaHepMCParticle::setPx(double px){
  //make new momentum as something is wrong with
  //the HepMC momentum setters

  HepMC::FourVector momentum(m_particle->momentum());
  momentum.setPx(px);
  m_particle->set_momentum(momentum);
}

void TauolaHepMCParticle::setPy(double py){
  HepMC::FourVector momentum(m_particle->momentum());
  momentum.setPy(py);
  m_particle->set_momentum(momentum);
}


void TauolaHepMCParticle::setPz(double pz){
  HepMC::FourVector momentum(m_particle->momentum());
  momentum.setPz(pz);
  m_particle->set_momentum(momentum);
}

void TauolaHepMCParticle::setE(double e){
  HepMC::FourVector momentum(m_particle->momentum());
  momentum.setE(e);
  m_particle->set_momentum(momentum);
}

} // namespace Tauolapp
Commit	Line	Data
0ca57c2f	1	#include "TauolaHepMCParticle.h"
e1938fed	2	#include "TauolaLog.h"
0ca57c2f	3
	4	namespace Tauolapp
	5	{
	6
	7	TauolaHepMCParticle::TauolaHepMCParticle(){
	8	m_particle = new HepMC::GenParticle();
	9	}
	10
	11	TauolaHepMCParticle::~TauolaHepMCParticle(){
	12
	13	//delete the mother and daughter pointers
	14	while(m_mothers.size()!=0){
	15	TauolaParticle * temp = m_mothers.back();
	16	m_mothers.pop_back();
	17	delete temp;
	18	}
	19	while(m_daughters.size()!=0){
	20	TauolaParticle * temp = m_daughters.back();
	21	m_daughters.pop_back();
	22	delete temp;
	23	}
	24
	25	while(m_created_particles.size()!=0){
	26	TauolaHepMCParticle * temp = (TauolaHepMCParticle*) m_created_particles.back();
	27	m_created_particles.pop_back();
	28	if(temp->getHepMC()->barcode()==0) delete temp->getHepMC();
	29	delete temp;
	30	}
	31
	32	}
	33
	34	// NOTE: Not executed by release examples
	35	TauolaHepMCParticle::TauolaHepMCParticle(int pdg_id, int status, double mass){
	36	m_particle = new HepMC::GenParticle();
	37	m_particle->set_pdg_id(pdg_id);
	38	m_particle->set_status(status);
	39	m_particle->set_generated_mass(mass);
	40	}
	41
	42	TauolaHepMCParticle::TauolaHepMCParticle(HepMC::GenParticle * particle){
	43	m_particle = particle;
	44	}
	45
	46	HepMC::GenParticle * TauolaHepMCParticle::getHepMC(){
	47	return m_particle;
	48	}
	49
	50	void TauolaHepMCParticle::undecay(){
	51	std::vector<TauolaParticle*> daughters = getDaughters();
	52	std::vector<TauolaParticle*>::iterator dIter = daughters.begin();
	53
	54	for(; dIter != daughters.end(); dIter++)
	55	(*dIter)->undecay();
	56
	57	if(m_particle->end_vertex())
	58	{
	59	while(m_particle->end_vertex()->particles_out_size())
	60	{
	61	HepMC::GenParticle p = m_particle->end_vertex()->remove_particle((m_particle->end_vertex()->particles_out_const_begin()));
	62	delete p;
	63	}
	64	delete m_particle->end_vertex();
	65	}
	66
67	m_daughters.clear();
68	m_particle->set_status(TauolaParticle::STABLE);
69
70	for(unsigned int i=0;i<daughters.size();i++)
71	delete daughters[i];
72	}
73
74	void TauolaHepMCParticle::setMothers(vector<TauolaParticle*> mothers){
75
76	/****** Deal with mothers *********/
77
78	//If there are mothers
79	if(mothers.size()>0){
80
81	HepMC::GenParticle * part;
82	part=dynamic_cast<TauolaHepMCParticle*>(mothers.at(0))->getHepMC();
83
84	//Use end vertex of first mother as production vertex for particle
85	HepMC::GenVertex * production_vertex = part->end_vertex();
86	HepMC::GenVertex * orig_production_vertex = production_vertex;
87
88	//If production_vertex does not exist - create it
89	//If it's tau decay - set the time and position including the tau lifetime correction
90	//otherwise - copy the time and position of decaying particle
91	if(!production_vertex){
92	production_vertex = new HepMC::GenVertex();
93	HepMC::FourVector point = part->production_vertex()->position();
94	production_vertex->set_position(point);
95	part->parent_event()->add_vertex(production_vertex);
96	}
97
98	//Loop over all mothers to check that the end points to the right place
99	vector<TauolaParticle*>::iterator mother_itr;
100	for(mother_itr = mothers.begin(); mother_itr != mothers.end();
101	mother_itr++){
102
103	HepMC::GenParticle * moth;
104	moth = dynamic_cast<TauolaHepMCParticle>(mother_itr)->getHepMC();
105
106	if(moth->end_vertex()!=orig_production_vertex)
107	Log::Fatal("Mother production_vertices point to difference places. Can not override. Please delete vertices first.",1);
108	else
109	production_vertex->add_particle_in(moth);
110
111	//update status info
112	if(moth->status()==TauolaParticle::STABLE)
113	moth->set_status(TauolaParticle::DECAYED);
114	}
115	production_vertex->add_particle_out(m_particle);
116	}
117	}
118
119	void TauolaHepMCParticle::setDaughters(vector<TauolaParticle*> daughters){
120
121	if(!m_particle->parent_event())
122	Log::Fatal("New particle needs the event set before it's daughters can be added",2);
123
124	//If there are daughters
125	if(daughters.size()>0){
126	// NOTE: Not executed by release examples
127	// because daughters.size() is always 0
128
129	//Use production vertex of first daughter as end vertex for particle
130	HepMC::GenParticle * first_daughter;
131	first_daughter = (dynamic_cast<TauolaHepMCParticle*>(daughters.at(0)))->getHepMC();
132
133	HepMC::GenVertex * end_vertex;
134	end_vertex=first_daughter->production_vertex();
135	HepMC::GenVertex * orig_end_vertex = end_vertex;
136
137	if(!end_vertex){ //if it does not exist create it
138	end_vertex = new HepMC::GenVertex();
139	m_particle->parent_event()->add_vertex(end_vertex);
140	}
141
142	//Loop over all daughters to check that the end points to the right place
143	vector<TauolaParticle*>::iterator daughter_itr;
144	for(daughter_itr = daughters.begin(); daughter_itr != daughters.end();
145	daughter_itr++){
146
147	HepMC::GenParticle * daug;
148	daug = dynamic_cast<TauolaHepMCParticle>(daughter_itr)->getHepMC();
149
150
151	if(daug->production_vertex()!=orig_end_vertex)
152	Log::Fatal("Daughter production_vertices point to difference places. Can not override. Please delete vertices first.",3);
153	else
154	end_vertex->add_particle_out(daug);
155	}
156	end_vertex->add_particle_in(m_particle);
157	}
158	}
159
160	std::vector<TauolaParticle*> TauolaHepMCParticle::getMothers(){
161
162	if(m_mothers.size()==0&&m_particle->production_vertex()){
163	HepMC::GenVertex::particles_in_const_iterator pcle_itr;
164	pcle_itr=m_particle->production_vertex()->particles_in_const_begin();
165
166	HepMC::GenVertex::particles_in_const_iterator pcle_itr_end;
167	pcle_itr_end=m_particle->production_vertex()->particles_in_const_end();
168
169	for(;pcle_itr != pcle_itr_end; pcle_itr++){
170	m_mothers.push_back(new TauolaHepMCParticle(*pcle_itr));
171	}
172	}
173	return m_mothers;
174	}
175
176	std::vector<TauolaParticle*> TauolaHepMCParticle::getDaughters(){
177
178	if(m_daughters.size()==0&&m_particle->end_vertex()){
179	HepMC::GenVertex::particles_out_const_iterator pcle_itr;
180	pcle_itr=m_particle->end_vertex()->particles_out_const_begin();
181
182	HepMC::GenVertex::particles_out_const_iterator pcle_itr_end;
183	pcle_itr_end=m_particle->end_vertex()->particles_out_const_end();
184
185	for(;pcle_itr != pcle_itr_end; pcle_itr++){
186	m_daughters.push_back(new TauolaHepMCParticle(*pcle_itr));
187	}
188	}
189	return m_daughters;
190	}
191
192	void TauolaHepMCParticle::checkMomentumConservation(){
193
194	if(!m_particle->end_vertex()) return;
195
196	// HepMC version of check_momentum_conservation
197	// with added energy check
198
199	double sumpx = 0, sumpy = 0, sumpz = 0, sume = 0;
200	for( HepMC::GenVertex::particles_in_const_iterator part1 = m_particle->end_vertex()->particles_in_const_begin();
201	part1 != m_particle->end_vertex()->particles_in_const_end(); part1++ ){
202
203	sumpx += (*part1)->momentum().px();
204	sumpy += (*part1)->momentum().py();
205	sumpz += (*part1)->momentum().pz();
206	sume += (*part1)->momentum().e();
207	}
208
209	for( HepMC::GenVertex::particles_out_const_iterator part2 = m_particle->end_vertex()->particles_out_const_begin();
210	part2 != m_particle->end_vertex()->particles_out_const_end(); part2++ ){
211
212	sumpx -= (*part2)->momentum().px();
213	sumpy -= (*part2)->momentum().py();
214	sumpz -= (*part2)->momentum().pz();
215	sume -= (*part2)->momentum().e();
216	}
217
218	if( sqrt( sumpxsumpx + sumpysumpy + sumpzsumpz + sumesume) > Tauola::momentum_conservation_threshold ) {
219	Log::Warning()<<"Momentum not conserved in the vertex:"<<endl;
220	Log::RedirectOutput(Log::Warning(false));
221	m_particle->end_vertex()->print();
222	Log::RevertOutput();
223	return;
224	}
225
226	return;
227	}
228
229	// NOTE: Not executed by release examples
230	void TauolaHepMCParticle::setPdgID(int pdg_id){
231	m_particle->set_pdg_id(pdg_id);
232	}
233
234	void TauolaHepMCParticle::setMass(double mass){
235	m_particle->set_generated_mass(mass);
236	}
237
238	// NOTE: Not executed by release examples
239	void TauolaHepMCParticle::setStatus(int status){
240	m_particle->set_status(status);
241	}
242
243	int TauolaHepMCParticle::getPdgID(){
244	return m_particle->pdg_id();
245	}
246
247	int TauolaHepMCParticle::getStatus(){
248	return m_particle->status();
249	}
250
251	int TauolaHepMCParticle::getBarcode(){
252	return m_particle->barcode();
253	}
254
255	// Set (X,T) Position of tau decay trees
256	void TauolaHepMCParticle::decayEndgame(){
257
258	double lifetime = Tauola::tau_lifetime * (-log( Tauola::randomDouble() ));
259	HepMC::FourVector tau_momentum = m_particle->momentum();
260
261	double mass = sqrt(abs( tau_momentum.e()*tau_momentum.e()
262	- tau_momentum.px()*tau_momentum.px()
263	- tau_momentum.py()*tau_momentum.py()
264	- tau_momentum.pz()*tau_momentum.pz()
265	) );
266
267	// Get previous position
268	HepMC::FourVector previous_position = m_particle->production_vertex()->position();
269
270	// Calculate new position
271	HepMC::FourVector new_position(previous_position.x()+tau_momentum.px()/mass*lifetime,
272	previous_position.y()+tau_momentum.py()/mass*lifetime,
273	previous_position.z()+tau_momentum.pz()/mass*lifetime,
274	previous_position.t()+tau_momentum.e() /mass*lifetime);
275
276	// Set new position
277	m_particle->end_vertex()->set_position(new_position);
278	recursiveSetPosition(m_particle,new_position);
279	}
280
281	void TauolaHepMCParticle::recursiveSetPosition(HepMC::GenParticle *p, HepMC::FourVector pos){
282
283	if(!p->end_vertex()) return;
284
285	// Iterate over all outgoing particles
286	for(HepMC::GenVertex::particles_out_const_iterator pp = p->end_vertex()->particles_out_const_begin();
287	pp != p->end_vertex()->particles_out_const_end();
288	++pp){
289	if( !(*pp)->end_vertex() ) continue;
290
291	// Set position
292	(*pp)->end_vertex()->set_position(pos);
293	recursiveSetPosition(*pp,pos);
294	}
295	}
296
297	TauolaHepMCParticle * TauolaHepMCParticle::createNewParticle(
298	int pdg_id, int status, double mass,
299	double px, double py, double pz, double e){
300
301	TauolaHepMCParticle * new_particle = new TauolaHepMCParticle();
302	new_particle->getHepMC()->set_pdg_id(pdg_id);
303	new_particle->getHepMC()->set_status(status);
304	new_particle->getHepMC()->set_generated_mass(mass);
305
306	HepMC::FourVector momentum(px,py,pz,e);
307	new_particle->getHepMC()->set_momentum(momentum);
308
309	m_created_particles.push_back(new_particle);
310
311	return new_particle;
312	}
313
314	void TauolaHepMCParticle::print(){
315	m_particle->print();
316	}
317
318
319	/****** Getter and Setter methods: *********************/
320
321	inline double TauolaHepMCParticle::getPx(){
322	return m_particle->momentum().px();
323	}
324
325	inline double TauolaHepMCParticle::getPy(){
326	return m_particle->momentum().py();
327	}
328
329	double TauolaHepMCParticle::getPz(){
330	return m_particle->momentum().pz();
331	}
332
333	double TauolaHepMCParticle::getE(){
334	return m_particle->momentum().e();
335	}
336
337	void TauolaHepMCParticle::setPx(double px){
338	//make new momentum as something is wrong with
339	//the HepMC momentum setters
340
341	HepMC::FourVector momentum(m_particle->momentum());
342	momentum.setPx(px);
343	m_particle->set_momentum(momentum);
344	}
345
346	void TauolaHepMCParticle::setPy(double py){
347	HepMC::FourVector momentum(m_particle->momentum());
348	momentum.setPy(py);
349	m_particle->set_momentum(momentum);
350	}
351
352
353	void TauolaHepMCParticle::setPz(double pz){
354	HepMC::FourVector momentum(m_particle->momentum());
355	momentum.setPz(pz);
356	m_particle->set_momentum(momentum);
357	}
358
359	void TauolaHepMCParticle::setE(double e){
360	HepMC::FourVector momentum(m_particle->momentum());
361	momentum.setE(e);
362	m_particle->set_momentum(momentum);
363	}
364
365	} // namespace Tauolapp