Skip to content

GitLab

Commit 5dd89f63 authored by Lynn Garren's avatar Lynn Garren
Browse files

more changes from HepMC 1.24

parent 8f1948a6
Hide whitespace changes
Inline Side-by-side
HepMC/HepMC/CBhepevt.h
View file @ 5dd89f63
// $Id: CBhepevt.h,v 1.1.1.1 2003/07/15 20:15:05 garren Exp $
// $Id: CBhepevt.h,v 1.2 2003/10/08 19:36:47 garren Exp $
// ----------------------------------------------------------------------
// CBhepevt.h
// ----------------------------------------------------------------------
......@@ -37,7 +37,8 @@ public:
hptr3( & hepev3_ ),
hptr4( & hepev4_ ),
hptr5( & hepev5_ ),
itsTrustMothers( true )
itsTrustMothers( true ),
itsTrustMothersAndDaughters( false )
{ ; }
hepevt_t * hepevt() { return hptr; }
......@@ -47,17 +48,32 @@ public:
hepev5_t * hepev5() { return hptr5; }
static int max_number_entries() { return NMXHEP; }
static int max_multiple_interactions() { return NMXMLT; }
int event_number() { return hptr->nevhep; }
int number_entries() { return hptr->nhep; }
int first_parent( int index ); // index of 1st mother
int last_parent( int index ); // index of last mother
int first_child( int index ); // index of 1st daughter
int last_child( int index ); // index of last daughter
int number_children( int index );
int number_parents( int index );
void clean( );
void zero_everything() { clean(); }
bool check_hepevt_consistency( std::ostream & os = std::cout );
bool toGenEvent( GenEvent *, bool printErrors );
bool fromGenEvent( const GenEvent * );
bool addtoHEPEVT( const GenEvent * ); // not yet implemented
void print( std::ostream & os = std::cout ) const;
void print_hepevt( std::ostream & os = std::cout );
void print_legend( std::ostream & os = std::cout );
void print_hepevt_particle( int index, std::ostream & os = std::cout );
// decide how to deal with HEPEVT
bool trustMothers() const { return itsTrustMothers; }
bool trustMothersAndDaughters() const { return itsTrustMothersAndDaughters; }
void setTrustMothers( bool b ) { itsTrustMothers = b; }
void setTrustMothersAndDaughters( bool b ) { itsTrustMothersAndDaughters = b; }
private:
......@@ -67,6 +83,7 @@ private:
hepev4_t * hptr4;
hepev5_t * hptr5;
bool itsTrustMothers;
bool itsTrustMothersAndDaughters;
// internal functions
GenParticle* createParticle( int index );
......
HepMC/HepMC/CBhepevt.icc
View file @ 5dd89f63
// $Id: CBhepevt.icc,v 1.1.1.1 2003/07/15 20:15:05 garren Exp $
// $Id: CBhepevt.icc,v 1.2 2003/10/08 19:36:47 garren Exp $
// ----------------------------------------------------------------------
// CBhepevt.icc
// ----------------------------------------------------------------------
......@@ -22,6 +22,234 @@ int CBhepevt::number_parents( int index ) {
return secondparent - firstparent + 1;
}
inline int CBhepevt::first_parent( int index )
{
int parent = hptr->jmohep[index][0];
return ( parent > 0 && parent <= number_entries() ) ?
parent : 0;
}
inline int CBhepevt::last_parent( int index )
{
// Returns the Index of the LAST parent in the HEPEVT record
// for particle with Index index.
// If there is only one parent, the last parent is forced to
// be the same as the first parent.
// If there are no parents for this particle, both the first_parent
// and the last_parent with return 0.
// Error checking is done to ensure the parent is always
// within range ( 0 <= parent <= nhep )
//
int firstparent = first_parent(index);
int parent = hptr->jmohep[index][1];
return ( parent > firstparent && parent <= number_entries() )
? parent : firstparent;
}
inline int CBhepevt::first_child( int index )
{
int child = hptr->jdahep[index][0];
return ( child > 0 && child <= number_entries() ) ?
child : 0;
}
inline int CBhepevt::last_child( int index )
{
// Returns the Index of the LAST child in the HEPEVT record
// for particle with Index index.
// If there is only one child, the last child is forced to
// be the same as the first child.
// If there are no children for this particle, both the first_child
// and the last_child with return 0.
// Error checking is done to ensure the child is always
// within range ( 0 <= parent <= nhep )
//
int firstchild = first_child(index);
int child = hptr->jdahep[index][1];
return ( child > firstchild && child <= number_entries() )
? child : firstchild;
}
void CBhepevt::print( std::ostream & os ) const
{
os << "CBhepevt (gets an event from HEPEVT common block)" << std::endl;
os << " trustMothers: " << itsTrustMothers << std::endl;
os << " trustMothersAndDaughters: " << itsTrustMothersAndDaughters << std::endl;
}
void CBhepevt::print_hepevt( std::ostream & os )
{
// dumps the content of this HEPEVT event to os (Width is 80)
os << "________________________________________"
<< "________________________________________" << std::endl;
os << "***** HEPEVT Common Event#: "
<< event_number()
<< ", " << number_entries() << " particles (max "
<< max_number_entries() << ") *****";
os << " Double Precision" << std::endl;
print_legend(os);
os << "________________________________________"
<< "________________________________________" << std::endl;
for ( int i=1; i <= number_entries(); ++i ) {
print_hepevt_particle( i, os );
}
os << "________________________________________"
<< "________________________________________" << std::endl;
}
void CBhepevt::print_legend( std::ostream & os )
{
os << "Indx Stat Par- chil- "
<< "( P_x, P_y, P_z, Energy M ) [GeV]" << std::endl;
os << " ID ents dren "
<< "Prod ( X, Y, Z, cT ) [mm]" << std::endl;
}
void CBhepevt::print_hepevt_particle( int i, std::ostream & os )
{
// dumps the content HEPEVT particle entry i (Width is 120)
// here i is the C array index (i.e. it starts at 0 ... whereas the
// fortran array index starts at 1) So if there's 100 particles, the
// last valid index is 100-1=99
char outline[81];
sprintf( outline,
"%4d %+4d %4d %4d (%9.3g, %9.3g, %9.3g, %9.3g, %9.3g)"
,i, hptr->isthep[i], first_parent(i), first_child(i),
hptr->phep[i][0], hptr->phep[i][1], hptr->phep[i][2],
hptr->phep[i][3], hptr->phep[i][4] );
os << outline << "\n";
sprintf( outline,"%+9d %4d %4d (%9.3g, %9.3g, %9.3g, %9.3g)",
// old version was:" (%+9.2e, %+9.2e, %+9.2e, %+9.2e)"
hptr->idhep[i], last_parent(i), last_child(i),
hptr->vhep[i][0], hptr->vhep[i][1], hptr->vhep[i][2],
hptr->vhep[i][3] );
os << outline << std::endl;
}
bool CBhepevt::check_hepevt_consistency( std::ostream & os )
{
// This method inspects the HEPEVT common block and looks for
// inconsistencies in the mother/daughter pointers
bool isConsistent=true;
char header[81];
sprintf( header,
"\n\n\t**** WARNINGInconsistent HEPEVT input, Event %10d ****"
, CBhepevt::event_number() );
for ( int i = 1; i <= CBhepevt::number_entries(); ++i ) {
// 1. check its mothers
int moth1 = CBhepevt::first_parent( i );
int moth2 = CBhepevt::last_parent( i );
if ( moth2<moth1 ) {
if ( isConsistent ) {
os << header << std::endl;
isConsistent = false;
print_legend(os);
}
os << "Inconsistent entry " << i
<< " first parent > last parent " << std::endl;
CBhepevt::print_hepevt_particle( i, os );
}
for ( int m = moth1; m<=moth2 && m!=0; ++m ) {
if ( m>CBhepevt::number_entries() || m < 0 ) {
if ( isConsistent ) {
os << header << std::endl;
isConsistent = false;
print_legend(os);
}
os << "Inconsistent entry " << i
<< " mother points out of range " << std::endl;
CBhepevt::print_hepevt_particle( i, os );
}
int mChild1 = CBhepevt::first_child(m);
int mChild2 = CBhepevt::last_child(m);
// we don't consider null pointers as inconsistent
if ( mChild1==0 && mChild2==0 ) continue;
if ( i<mChild1 || i>mChild2 ) {
if ( isConsistent ) {
os << header << std::endl;
isConsistent = false;
print_legend(os);
}
os << "Inconsistent mother-daughter relationship between "
<< i << " & " << m
<< " (try !trustMother)" << std::endl;
CBhepevt::print_hepevt_particle( i, os );
CBhepevt::print_hepevt_particle( m, os );
}
}
// 2. check its daughters
int dau1 = CBhepevt::first_child( i );
int dau2 = CBhepevt::last_child( i );
if ( dau2<dau1 ) {
if ( isConsistent ) {
os << header << std::endl;
isConsistent = false;
print_legend(os);
}
os << "Inconsistent entry " << i
<< " first child > last child " << std::endl;
CBhepevt::print_hepevt_particle( i, os );
}
for ( int d = dau1; d<=dau2 && d!=0; ++d ) {
if ( d>CBhepevt::number_entries() || d < 0 ) {
if ( isConsistent ) {
os << header << std::endl;
isConsistent = false;
print_legend(os);
}
os << "Inconsistent entry " << i
<< " child points out of range " << std::endl;
CBhepevt::print_hepevt_particle( i, os );
}
int d_moth1 = CBhepevt::first_parent(d);
int d_moth2 = CBhepevt::last_parent(d);
// we don't consider null pointers as inconsistent
if ( d_moth1==0 && d_moth2==0 ) continue;
if ( i<d_moth1 || i>d_moth2 ) {
if ( isConsistent ) {
os << header << std::endl;
isConsistent = false;
print_legend(os);
}
os << "Inconsistent mother-daughter relationship between "
<< i << " & " << d
<< " (try trustMothers)"<< std::endl;
CBhepevt::print_hepevt_particle( i, os );
CBhepevt::print_hepevt_particle( d, os );
}
}
}
if (!isConsistent) {
os << "Above lists all the inconsistencies in the HEPEVT common "
<< "\n block which has been provided as input to HepMC. "
<< "\n HepMC WILL have trouble interpreting the mother-daughter"
<< "\n relationships ... but all other information "
<< "\n (4-vectors etc) will be correctly transferred."
<< "\n In order for HepMC to be able to interpret the mother/"
<< "\n daughter hierachy, it MUST be given consistent input."
<< "\n This is one of the design criteria of HepMC: "
<< "\n consistency is enforced by the code.";
os << "\nThere is a switch in CBhepevt, set-able using "
<< "\n CBhepevt::setTrustMothers( bool )"
<< "\n which you may want to try.";
os << "\nNote: if HEPEVT common block has been filled by pythia"
<< "\n pyhepc, then the switch MSTP(128)=2 should be used in"
<< "\n pythia, which instructs pythia not to put multiple "
<< "\n copies of resonances in the event record.\n";
os << "To obtain a file summarizing the inconsistency, you should:"
<< "\n\t ofstream myFile(\"myInconsistentEvent.txt\"); "
<< "\n\t CBhepevt::check_hepevt_consistency(myFile); "
<< "\n\t CBhepevt::print_hepevt(myFile); "
<< "\n[now write the event to HepMC using something like"
<< "\n\t\t myIO_HEPEVT->write_event(myEvent); ]"
<< "\n\t myEvent->print( myFile ); "
<< " // print event as HepMC sees it"
<< "\n ------------------------- Thank-you. \n\n" << std::endl;
}
return isConsistent;
}
void CBhepevt::clean( )
{
int i, k;
......
HepMC/HepMC/GenEventConvert.h
View file @ 5dd89f63
// $Id: GenEventConvert.h,v 1.1.1.1 2003/07/15 20:15:05 garren Exp $
// $Id: GenEventConvert.h,v 1.2 2003/10/08 19:36:47 garren Exp $
// ----------------------------------------------------------------------
//
// GenEventConvert.h
......@@ -41,6 +41,7 @@ public:
bool toGenEvent( GenEvent * );
bool fromGenEvent( const GenEvent * );
bool addtoHEPEVT( const GenEvent * );
bool fill_next_event( GenEvent * evt ) { toGenEvent(evt); }
// decide how to deal with HEPEVT
bool printInconsistencyErrors() const { return itsInconsitencyErrors; }
......
HepMC/HepMC/GenEventConvert.icc
View file @ 5dd89f63
// $Id: GenEventConvert.icc,v 1.1.1.1 2003/07/15 20:15:05 garren Exp $
// $Id: GenEventConvert.icc,v 1.2 2003/10/08 19:36:47 garren Exp $
// ----------------------------------------------------------------------
//
// GenEventConvert.icc
......@@ -37,7 +37,7 @@ bool GenEventConvert< S >::toGenEvent( GenEvent* evt ) {
// read one event from the HEPEVT common block and fill GenEvent
// return T/F =success/failure
//
// 1. test that evt pointer is not null and set event number
// 1. test that evt pointer is not null
if ( !evt ) {
std::cerr
<< "GenEventConvert:toGenEvent: error - passed null event."
......@@ -46,6 +46,7 @@ bool GenEventConvert< S >::toGenEvent( GenEvent* evt ) {
}
CBInterface< S > * hep = CBInterface< S >::instance();
bool perr = printInconsistencyErrors();
// now fill GenEvent
return (*hep)->toGenEvent( evt, perr );
}
......@@ -61,12 +62,14 @@ bool GenEventConvert< S >::fromGenEvent( const GenEvent* evt ) {
// get the common block and clean it
CBInterface< S > * hep = CBInterface< S >::instance();
(*hep)->clean( );
// now fill HEPEVT
return (*hep)->fromGenEvent( evt );
}
template<class S>
bool GenEventConvert< S >::addtoHEPEVT( const GenEvent* evt ) {
CBInterface< S > * hep = CBInterface< S >::instance();
// add this event to existing GenEvent
return (*hep)->addtoHEPEVT( evt );
}
......
HepMC/HepMC/HEPEVTtoGenEvent.icc
View file @ 5dd89f63
// $Id: HEPEVTtoGenEvent.icc,v 1.1.1.1 2003/07/15 20:15:05 garren Exp $
// $Id: HEPEVTtoGenEvent.icc,v 1.2 2003/10/08 19:36:47 garren Exp $
// ----------------------------------------------------------------------
//
// toGenEvent.icc
......@@ -21,13 +21,14 @@ bool CBhepevt::toGenEvent( GenEvent* evt, bool printErrors ) {
// the children pointers are not always correct (i.e. there is
// oftentimes an internal inconsistency between the parents and
// children pointers). The parent pointers always seem to be correct.
// Thus the switch trust_mothers_before_daughters=1 is appropriate for
// Thus the switch trustMothers=1 is appropriate for
// pythia. NOTE: you should also set the switch MSTP(128) = 2 in
// pythia (not the default!), so that pythia doesn't
// store two copies of resonances in the event record.
// The situation is opposite for the HEPEVT which comes from Isajet
// via stdhep, so then use the switch trust_mothers_before_daughters=0
// via stdhep, so then use the switch trustMothers=0
//
// 1. set event number
evt->set_event_number( hptr->nevhep );
// these do not have the correct defaults if hepev4 is not filled
//evt->set_event_scale( hepev4()->scalelh[1] );
......@@ -54,13 +55,13 @@ bool CBhepevt::toGenEvent( GenEvent* evt, bool printErrors ) {
/// sufficient to do one or the other.
//
// 3. Build the production_vertex (if necessary)
if ( trustMothers() ) {
if ( trustMothers() || trustMothersAndDaughters() ) {
buildProductionVertex( i, hepevt_particle, evt, printErrors );
}
//
// 4. Build the end_vertex (if necessary)
// Identical steps as for production vertex
if ( !trustMothers() ) {
if ( !trustMothers() || trustMothersAndDaughters() ) {
buildEndVertex( i, hepevt_particle, evt, printErrors );
}
}
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment