12 #include "HepMC/GenEvent.h"
13 #include "HepMC/IO_HEPEVT.h"
14 #include "HepMC/HEPEVT_Wrapper.h"
93 for (
int ip=0; ip<4000; ip++ )
102 for (
int iv=1; iv<=evt->vertices_size(); iv++ )
105 bool legalVtx =
false;
109 HepMC::GenVertex* vtx = evt->barcode_to_vertex( -iv ) ;
111 if ( vtx->particles_in_size() != 1 )
continue;
112 if ( vtx->particles_out_size() <= 1 )
continue;
114 if ( (*(vtx->particles_in_const_begin()))->pdg_id() == 111 )
continue;
116 if (
fOnlyPDG != 1 && (*(vtx->particles_in_const_begin()))->pdg_id() !=
fOnlyPDG )
127 for ( HepMC::GenVertex::particle_iterator pitr=vtx->particles_begin(HepMC::children);
128 pitr != vtx->particles_end(HepMC::children); ++pitr)
130 if ( (*pitr)->pdg_id() ==
fOnlyPDG )
136 if ( same )
continue;
157 for ( HepMC::GenVertex::particle_iterator pitr=vtx->particles_begin(HepMC::children);
158 pitr != vtx->particles_end(HepMC::children); ++pitr)
162 if (
abs((*pitr)->pdg_id()) >=1 &&
abs((*pitr)->pdg_id()) <=8 )
break;
163 if (
abs((*pitr)->pdg_id()) == 21 )
break;
165 if ( (*pitr)->status() == 1 || (*pitr)->end_vertex() )
173 if ( !legalVtx )
continue;
180 HepMC::HEPEVT_Wrapper::set_event_number( evt->event_number() );
189 HepMC::GenVertex* vtx = evt->barcode_to_vertex( vtxbcode );
199 HepMC::HEPEVT_Wrapper::zero_everything();
205 HepMC::HEPEVT_Wrapper::set_id( index, (*(vtx->particles_in_const_begin()))->pdg_id() );
206 HepMC::FourVector
vec4;
207 vec4 = (*(vtx->particles_in_const_begin()))->momentum();
208 HepMC::HEPEVT_Wrapper::set_momentum( index, vec4.x(), vec4.y(), vec4.z(), vec4.e() );
209 HepMC::HEPEVT_Wrapper::set_mass( index, (*(vtx->particles_in_const_begin()))->generated_mass() );
210 HepMC::HEPEVT_Wrapper::set_position( index, vtx->position().x(), vtx->position().y(),
211 vtx->position().z(), vtx->position().t() );
212 HepMC::HEPEVT_Wrapper::set_status( index, (*(vtx->particles_in_const_begin()))->status() );
213 HepMC::HEPEVT_Wrapper::set_parents( index, 0, 0 );
214 fBarcodes.push_back( (*(vtx->particles_in_const_begin()))->barcode() );
219 for ( HepMC::GenVertex::particle_iterator pitr=vtx->particles_begin(HepMC::children);
220 pitr != vtx->particles_end(HepMC::children); ++pitr)
223 if ( (*pitr)->status() == 1 || (*pitr)->end_vertex() )
226 vec4 = (*pitr)->momentum();
227 HepMC::HEPEVT_Wrapper::set_id( index, (*pitr)->pdg_id() );
228 HepMC::HEPEVT_Wrapper::set_momentum( index, vec4.x(), vec4.y(), vec4.z(), vec4.e() );
229 HepMC::HEPEVT_Wrapper::set_mass( index, (*pitr)->generated_mass() );
230 vec4 = (*pitr)->production_vertex()->position();
231 HepMC::HEPEVT_Wrapper::set_position( index, vec4.x(), vec4.y(), vec4.z(), vec4.t() );
232 HepMC::HEPEVT_Wrapper::set_status( index, (*pitr)->status() );
233 HepMC::HEPEVT_Wrapper::set_parents( index, 1, 1 );
234 fBarcodes.push_back( (*pitr)->barcode() );
237 if ( (*pitr)->end_vertex() )
239 fSecVtxStore.push_back( (*pitr)->end_vertex()->barcode() );
245 int nentries =
index;
249 HepMC::HEPEVT_Wrapper::set_children ( index, 2, lastDau );
254 HepMC::HEPEVT_Wrapper::set_number_entries( nentries );
292 unsigned int vcounter = 0;
309 if ( HepMC::HEPEVT_Wrapper::number_entries() > nentries )
320 int largestBarcode = -1;
323 for (
int ip=1; ip<Nbcodes; ip++ )
328 if ( bcode > largestBarcode ) largestBarcode = bcode;
329 double px = HepMC::HEPEVT_Wrapper::px(ip+1);
330 double py = HepMC::HEPEVT_Wrapper::py(ip+1);
331 double pz = HepMC::HEPEVT_Wrapper::pz(ip+1);
335 if ( prt->end_vertex() )
338 HepMC::GenVertex* endVtx = prt->end_vertex();
340 std::vector<int> secVtxStorage;
341 secVtxStorage.clear();
343 secVtxStorage.push_back( endVtx->barcode() );
345 HepMC::FourVector mom4 = prt->momentum();
348 double bet1[3], bet2[3], gam1, gam2, pb;
349 double mass = mom4.m();
350 bet1[0] = -(mom4.px()/mass);
351 bet1[1] = -(mom4.py()/mass);
352 bet1[2] = -(mom4.pz()/mass);
356 gam1 = mom4.e()/mass;
359 unsigned int vcounter = 0;
361 while ( vcounter < secVtxStorage.size() )
364 HepMC::GenVertex* theVtx = evt->barcode_to_vertex( secVtxStorage[vcounter] );
366 for ( HepMC::GenVertex::particle_iterator pitr=theVtx->particles_begin(HepMC::children);
367 pitr != theVtx->particles_end(HepMC::children); ++pitr)
370 if ( (*pitr)->end_vertex() )
372 secVtxStorage.push_back( (*pitr)->end_vertex()->barcode() );
375 if ( theVtx->particles_out_size() == 1 && (*pitr)->pdg_id() == prt->pdg_id() )
378 (*pitr)->set_momentum( HepMC::FourVector(px,py,pz,e) );
382 HepMC::FourVector dmom4 = (*pitr)->momentum();
385 pb = bet1[0]*dmom4.px() + bet1[1]*dmom4.py() + bet1[2]*dmom4.pz();
386 double dpx = dmom4.px() + bet1[0] * (dmom4.e() + pb/(gam1+1.) );
387 double dpy = dmom4.py() + bet1[1] * (dmom4.e() + pb/(gam1+1.) );
388 double dpz = dmom4.pz() + bet1[2] * (dmom4.e() + pb/(gam1+1.) );
389 double de = gam1*dmom4.e() + pb;
391 pb = bet2[0]*dpx + bet2[1]*dpy + bet2[2]*dpz;
392 dpx += bet2[0] * ( de + pb/(gam2+1.) );
393 dpy += bet2[1] * ( de + pb/(gam2+1.) );
394 dpz += bet2[2] * ( de + pb/(gam2+1.) );
398 (*pitr)->set_momentum( HepMC::FourVector(dpx,dpy,dpz,de) );
404 secVtxStorage.clear();
407 prt->set_momentum( HepMC::FourVector(px,py,pz,e) );
411 int newlyGen = HepMC::HEPEVT_Wrapper::number_entries() - nentries;
413 if ( largestBarcode < evt->particles_size() )
419 for (
int ipp=evt->particles_size(); ipp>largestBarcode; ipp-- )
421 (evt->barcode_to_particle(ipp))->suggest_barcode( ipp+newlyGen );
427 for (
int ipnw=1; ipnw<=newlyGen; ipnw++ )
429 int nbcode = largestBarcode+ipnw;
430 int pdg_id = HepMC::HEPEVT_Wrapper::id( nentries+ipnw );
432 double px = HepMC::HEPEVT_Wrapper::px( nentries+ipnw );
433 double py = HepMC::HEPEVT_Wrapper::py( nentries+ipnw );
434 double pz = HepMC::HEPEVT_Wrapper::pz( nentries+ipnw );
440 NewPart->set_generated_mass( m );
441 NewPart->suggest_barcode( nbcode );
442 vtx->add_particle_out( NewPart ) ;
456 if (
abs((*(vtx->particles_in_const_begin()))->pdg_id()) != 15 )
return false;
458 for ( HepMC::GenVertex::particle_iterator pitr=vtx->particles_begin(HepMC::children);
459 pitr != vtx->particles_end(HepMC::children); ++pitr)
461 if (
abs((*pitr)->pdg_id()) == 11 ||
abs((*pitr)->pdg_id()) == 13 )
void attachParticles(HepMC::GenEvent *, HepMC::GenVertex *, int)
double phoran_(int *idummy)
std::vector< int > fBarcodes
bool fAvoidTauLeptonicDecays
void applyToBranch(HepMC::GenEvent *, int)
std::vector< std::string > fSpecialSettings
std::vector< int > fSecVtxStore
void applyToVertex(HepMC::GenEvent *, int)
CLHEP::HepRandomEngine * decayRandomEngine
bool isTauLeptonicDecay(HepMC::GenVertex *)
HepMC::GenEvent * apply(HepMC::GenEvent *)
void configureOnlyFor(int)