Line data Source code
1 : //////////////////////////////////////////////////////////////////////////
2 : // Matt.Dobbs@Cern.CH June 30, 2000
3 : // Generic Wrapper for the fortran HEPEVT common block
4 : //
5 : // The static data member's initializations must be separate from .h file.
6 : //
7 : //////////////////////////////////////////////////////////////////////////
8 :
9 : #include "HepMC/HEPEVT_Wrapper.h"
10 :
11 : namespace HepMC {
12 :
13 : ////////////////////////////////////////
14 : // static data member initializations //
15 : ////////////////////////////////////////
16 :
17 : unsigned int HEPEVT_Wrapper::s_sizeof_int = 4;
18 :
19 : unsigned int HEPEVT_Wrapper::s_sizeof_real = sizeof(double);
20 :
21 : unsigned int HEPEVT_Wrapper::s_max_number_entries = 4000;
22 :
23 : ///////////////////
24 : // Print Methods //
25 : ///////////////////
26 :
27 : void HEPEVT_Wrapper::print_hepevt( std::ostream& ostr )
28 : {
29 : /// dumps the content of this HEPEVT event to ostr (Width is 80)
30 0 : ostr << "________________________________________"
31 0 : << "________________________________________" << std::endl;
32 0 : ostr << "***** HEPEVT Common Event#: "
33 0 : << event_number()
34 0 : << ", " << number_entries() << " particles (max "
35 0 : << max_number_entries() << ") *****";
36 0 : if ( is_double_precision() ) {
37 0 : ostr << " Double Precision" << std::endl;
38 0 : } else {
39 0 : ostr << " Single Precision" << std::endl;
40 : }
41 0 : ostr << sizeof_int() << "-byte integers, "
42 0 : << sizeof_real() << "-byte floating point numbers, "
43 0 : << max_number_entries() << "-allocated entries."
44 0 : << std::endl;
45 0 : print_legend(ostr);
46 0 : ostr << "________________________________________"
47 0 : << "________________________________________" << std::endl;
48 0 : for ( int i=1; i <= number_entries(); ++i ) {
49 0 : print_hepevt_particle( i, ostr );
50 : }
51 0 : ostr << "________________________________________"
52 0 : << "________________________________________" << std::endl;
53 0 : }
54 :
55 : void HEPEVT_Wrapper::print_legend( std::ostream& ostr )
56 : {
57 0 : char outline[81];
58 0 : sprintf( outline,"%4s %4s %4s %5s %10s, %9s, %9s, %9s, %10s",
59 : "Indx","Stat","Par-","chil-",
60 : "( P_x","P_y","P_z","Energy","M ) ");
61 0 : ostr << outline << std::endl;
62 0 : sprintf( outline,"%9s %4s %4s %10s, %9s, %9s, %9s) %9s",
63 : "ID ","ents","dren",
64 : "Prod ( X","Y","Z","cT", "[mm]");
65 0 : ostr << outline << std::endl;
66 0 : }
67 :
68 : void HEPEVT_Wrapper::print_hepevt_particle( int i, std::ostream& ostr )
69 : {
70 : /// dumps the content HEPEVT particle entry i (Width is 120)
71 : /// here i is the C array index (i.e. it starts at 0 ... whereas the
72 : /// fortran array index starts at 1) So if there's 100 particles, the
73 : /// last valid index is 100-1=99
74 0 : char outline[81];
75 0 : sprintf( outline,
76 : "%4d %+4d %4d %4d (%9.3g, %9.3g, %9.3g, %9.3g, %9.3g)"
77 0 : ,i, status(i), first_parent(i), first_child(i),
78 0 : px(i), py(i), pz(i), e(i), m(i) );
79 0 : ostr << outline << "\n";
80 0 : sprintf( outline,"%+9d %4d %4d (%9.3g, %9.3g, %9.3g, %9.3g)",
81 : // old version was:" (%+9.2e, %+9.2e, %+9.2e, %+9.2e)"
82 0 : id(i), last_parent(i), last_child(i),
83 0 : x(i), y(i), z(i), t(i) );
84 0 : ostr << outline << std::endl;
85 0 : }
86 :
87 :
88 : bool HEPEVT_Wrapper::check_hepevt_consistency( std::ostream& os )
89 : {
90 : /// This method inspects the HEPEVT common block and looks for
91 : /// inconsistencies in the mother/daughter pointers
92 : bool isConsistent=true;
93 0 : char header[81];
94 0 : sprintf( header,
95 : "\n\n\t**** WARNINGInconsistent HEPEVT input, Event %10d ****"
96 0 : , HEPEVT_Wrapper::event_number() );
97 :
98 0 : for ( int i = 1; i <= HEPEVT_Wrapper::number_entries(); ++i ) {
99 : // 1. check its mothers
100 0 : int moth1 = HEPEVT_Wrapper::first_parent( i );
101 0 : int moth2 = HEPEVT_Wrapper::last_parent( i );
102 0 : if ( moth2<moth1 ) {
103 0 : if ( isConsistent ) {
104 0 : os << header << std::endl;
105 : isConsistent = false;
106 0 : print_legend(os);
107 0 : }
108 0 : os << "Inconsistent entry " << i
109 0 : << " first parent > last parent " << std::endl;
110 0 : HEPEVT_Wrapper::print_hepevt_particle( i, os );
111 0 : }
112 0 : for ( int m = moth1; m<=moth2 && m!=0; ++m ) {
113 0 : if ( m>HEPEVT_Wrapper::number_entries() || m < 0 ) {
114 0 : if ( isConsistent ) {
115 0 : os << header << std::endl;
116 : isConsistent = false;
117 0 : print_legend(os);
118 0 : }
119 0 : os << "Inconsistent entry " << i
120 0 : << " mother points out of range " << std::endl;
121 0 : HEPEVT_Wrapper::print_hepevt_particle( i, os );
122 0 : }
123 0 : int mChild1 = HEPEVT_Wrapper::first_child(m);
124 0 : int mChild2 = HEPEVT_Wrapper::last_child(m);
125 : // we don't consider null pointers as inconsistent
126 0 : if ( mChild1==0 && mChild2==0 ) continue;
127 0 : if ( i<mChild1 || i>mChild2 ) {
128 0 : if ( isConsistent ) {
129 0 : os << header << std::endl;
130 : isConsistent = false;
131 0 : print_legend(os);
132 0 : }
133 0 : os << "Inconsistent mother-daughter relationship between "
134 0 : << i << " & " << m
135 0 : << " (try !trust_mother)" << std::endl;
136 0 : HEPEVT_Wrapper::print_hepevt_particle( i, os );
137 0 : HEPEVT_Wrapper::print_hepevt_particle( m, os );
138 0 : }
139 0 : }
140 : // 2. check its daughters
141 0 : int dau1 = HEPEVT_Wrapper::first_child( i );
142 0 : int dau2 = HEPEVT_Wrapper::last_child( i );
143 0 : if ( dau2<dau1 ) {
144 0 : if ( isConsistent ) {
145 0 : os << header << std::endl;
146 : isConsistent = false;
147 0 : print_legend(os);
148 0 : }
149 0 : os << "Inconsistent entry " << i
150 0 : << " first child > last child " << std::endl;
151 0 : HEPEVT_Wrapper::print_hepevt_particle( i, os );
152 0 : }
153 0 : for ( int d = dau1; d<=dau2 && d!=0; ++d ) {
154 0 : if ( d>HEPEVT_Wrapper::number_entries() || d < 0 ) {
155 0 : if ( isConsistent ) {
156 0 : os << header << std::endl;
157 : isConsistent = false;
158 0 : print_legend(os);
159 0 : }
160 0 : os << "Inconsistent entry " << i
161 0 : << " child points out of range " << std::endl;
162 0 : HEPEVT_Wrapper::print_hepevt_particle( i, os );
163 0 : }
164 0 : int d_moth1 = HEPEVT_Wrapper::first_parent(d);
165 0 : int d_moth2 = HEPEVT_Wrapper::last_parent(d);
166 : // we don't consider null pointers as inconsistent
167 0 : if ( d_moth1==0 && d_moth2==0 ) continue;
168 0 : if ( i<d_moth1 || i>d_moth2 ) {
169 0 : if ( isConsistent ) {
170 0 : os << header << std::endl;
171 : isConsistent = false;
172 0 : print_legend(os);
173 0 : }
174 0 : os << "Inconsistent mother-daughter relationship between "
175 0 : << i << " & " << d
176 0 : << " (try trust_mothers)"<< std::endl;
177 0 : HEPEVT_Wrapper::print_hepevt_particle( i, os );
178 0 : HEPEVT_Wrapper::print_hepevt_particle( d, os );
179 0 : }
180 0 : }
181 : }
182 0 : if (!isConsistent) {
183 0 : os << "Above lists all the inconsistencies in the HEPEVT common "
184 0 : << "\n block which has been provided as input to HepMC. "
185 0 : << "\n HepMC WILL have trouble interpreting the mother-daughter"
186 0 : << "\n relationships ... but all other information "
187 0 : << "\n (4-vectors etc) will be correctly transferred."
188 0 : << "\n In order for HepMC to be able to interpret the mother/"
189 0 : << "\n daughter hierachy, it MUST be given consistent input."
190 0 : << "\n This is one of the design criteria of HepMC: "
191 0 : << "\n consistency is enforced by the code.";
192 0 : os << "\nThere is a switch in IO_HEPEVT, set-able using "
193 0 : << "\n IO_HEPEVT::set_trust_mothers_before_daughters( bool )"
194 0 : << "\n which you may want to try.";
195 0 : os << "\nNote: if HEPEVT common block has been filled by pythia"
196 0 : << "\n pyhepc, then the switch MSTP(128)=2 should be used in"
197 0 : << "\n pythia, which instructs pythia not to put multiple "
198 0 : << "\n copies of resonances in the event record.\n";
199 0 : os << "To obtain a file summarizing the inconsistency, you should:"
200 0 : << "\n\t ofstream myFile(\"myInconsistentEvent.txt\"); "
201 0 : << "\n\t HEPEVT_Wrapper::check_hepevt_consistency(myFile); "
202 0 : << "\n\t HEPEVT_Wrapper::print_hepevt(myFile); "
203 0 : << "\n[now write the event to HepMC using something like"
204 0 : << "\n\t\t myIO_HEPEVT->write_event(myEvent); ]"
205 0 : << "\n\t myEvent->print( myFile ); "
206 0 : << " // print event as HepMC sees it"
207 0 : << "\n ------------------------- Thank-you. \n\n" << std::endl;
208 0 : }
209 0 : return isConsistent;
210 0 : }
211 :
212 : void HEPEVT_Wrapper::zero_everything()
213 : {
214 0 : set_event_number( 0 );
215 0 : set_number_entries( 0 );
216 0 : for ( int i = 1; i<=max_number_entries(); ++i ) {
217 0 : set_status( i, 0 );
218 0 : set_id( i, 0 );
219 0 : set_parents( i, 0, 0 );
220 0 : set_children( i, 0, 0 );
221 0 : set_momentum( i, 0, 0, 0, 0 );
222 0 : set_mass( i, 0 );
223 0 : set_position( i, 0, 0, 0, 0 );
224 : }
225 0 : }
226 :
227 : } // HepMC
228 :
|