1 | <HTML> |
---|
2 | <TITLE> |
---|
3 | </TITLE> |
---|
4 | <!-- Changed by: Katsuya Amako, 21-Sep-1998 --> |
---|
5 | <!-- Changed by: Dennis Wright, 25-Jun-2002 --> |
---|
6 | <!-- Proof read by: Joe Chuma, 28-Jun-1999 --> |
---|
7 | |
---|
8 | <BODY> |
---|
9 | <TABLE WIDTH="100%"><TR> |
---|
10 | <TD> |
---|
11 | |
---|
12 | |
---|
13 | <A HREF="index.html"> |
---|
14 | <IMG SRC="../../../../resources/html/IconsGIF/Contents.gif" ALT="Contents"></A> |
---|
15 | <A HREF="event.html"> |
---|
16 | <IMG SRC="../../../../resources/html/IconsGIF/Previous.gif" ALT="Previous"></A> |
---|
17 | <A HREF="biasing.html"> |
---|
18 | <IMG SRC="../../../../resources/html/IconsGIF/Next.gif" ALT="Next"></A> |
---|
19 | </TD> |
---|
20 | <TD ALIGN="Right"> |
---|
21 | <FONT SIZE="-1" COLOR="#238E23"> |
---|
22 | <B>Geant4 User's Guide</B> |
---|
23 | <BR> |
---|
24 | <B>For Application Developers</B> |
---|
25 | <BR> |
---|
26 | <B>Toolkit Fundamentals</B> |
---|
27 | </FONT> |
---|
28 | </TD> |
---|
29 | </TR></TABLE> |
---|
30 | <BR><BR> |
---|
31 | |
---|
32 | <P ALIGN="Center"> |
---|
33 | <FONT SIZE="+3" COLOR="#238E23"> |
---|
34 | <B>3.6 Event Generator Interface</B> |
---|
35 | </FONT> |
---|
36 | <BR><BR> |
---|
37 | |
---|
38 | <HR ALIGN="Center" SIZE="7%"> |
---|
39 | <p> |
---|
40 | |
---|
41 | <a name="3.6.1"> |
---|
42 | <H2>3.6.1 Structure of a primary event</H2></a> |
---|
43 | |
---|
44 | <b>Primary vertex and primary particle</b> |
---|
45 | <p> |
---|
46 | The <i>G4Event</i> class object should have a set of primary particles |
---|
47 | when it is sent to <i>G4EventManager</i> via <tt>processOneEvent()</tt> |
---|
48 | method. It is the mandate of your <i>G4VUserPrimaryGeneratorAction</i> |
---|
49 | concrete class to send primary particles to the <i>G4Event</i> |
---|
50 | object. |
---|
51 | <p> |
---|
52 | The <i>G4PrimaryParticle</i> class represents a primary particle |
---|
53 | with which Geant4 starts simulating an event. This class object |
---|
54 | has information on particle type and its three momenta. The |
---|
55 | positional and time information of primary particle(s) are stored |
---|
56 | in the <i>G4PrimaryVertex</i> class object and, thus, this class object |
---|
57 | can have one or more <i>G4PrimaryParticle</i> class objects which |
---|
58 | share the same vertex. As shown in Fig.?.?, primary vertexes and |
---|
59 | primary particles are associated with the <i>G4Event</i> object by |
---|
60 | a form of linked list. |
---|
61 | <p> |
---|
62 | A concrete class of <i>G4VPrimaryGenerator</i>, the |
---|
63 | <i>G4PrimaryParticle</i> object is constructed with either |
---|
64 | a pointer to <i>G4ParticleDefinition</i> or an integer number |
---|
65 | which represents P.D.G. particle code. For the case of some |
---|
66 | artificial particles, e.g., geantino, optical photon, etc., |
---|
67 | or exotic nuclear fragments, which the P.D.G. particle code does |
---|
68 | not cover, the <i>G4PrimaryParticle</i> should be constructed |
---|
69 | by <i>G4ParticleDefinition</i> pointer. On the other hand, |
---|
70 | elementary particles with very short life time, e.g., weak |
---|
71 | bosons, or quarks/gluons, can be instantiated as <i>G4PrimaryParticle</i> |
---|
72 | objects using the P.D.G. particle code. It should be noted that, |
---|
73 | even though primary particles with such a very short life time |
---|
74 | are defined, Geant4 will simulate only the particles which |
---|
75 | are defined as <i>G4ParticleDefinition</i> class objects. Other |
---|
76 | primary particles will be simply ignored by <i>G4EventManager</i>. |
---|
77 | But it may still be useful to construct such "intermediate" |
---|
78 | particles for recording the origin of the primary event. |
---|
79 | <p> |
---|
80 | <b>Forced decay channel</b> |
---|
81 | <p> |
---|
82 | The <i>G4PrimaryParticle</i> class object can have a list of its |
---|
83 | daughter particles. If the parent particle is an "intermediate" |
---|
84 | particle, which Geant4 does not have a corresponding <i>G4ParticleDefinition</i>, |
---|
85 | this parent particle is ignored and daughters |
---|
86 | are assumed to start from the vertex with which their parent is associated. |
---|
87 | For example, a Z boson is associated with a vertex and it has |
---|
88 | positive and negative muons as its daughters, these muons will |
---|
89 | start from that vertex. |
---|
90 | <p> |
---|
91 | There are some kinds of particles which should fly some reasonable |
---|
92 | distances and, thus, should be simulated by Geant4, but you still |
---|
93 | want to follow the decay channel generated by an event generator. |
---|
94 | A typical case of these particles is B meson. Even for the case of |
---|
95 | a primary particle which has a corresponding <i>G4ParticleDefinition</i>, |
---|
96 | it can have daughter primary particles. Geant4 will trace the |
---|
97 | parent particle until it comes to decay, obeying multiple |
---|
98 | scattering, ionization loss, rotation with the magnetic field, etc. |
---|
99 | according to its particle type. When the parent comes to decay, |
---|
100 | instead of randomly choosing its decay channel, it follows the |
---|
101 | "pre-assigned" decay channel. To conserve the energy and the momentum |
---|
102 | of the parent, daughters will be Lorentz transformed according to |
---|
103 | their parent's frame. |
---|
104 | <p> |
---|
105 | |
---|
106 | <HR> |
---|
107 | <a name="3.6.2"> |
---|
108 | <H2>3.6.2 Interface to a primary generator</H2></a> |
---|
109 | |
---|
110 | <b><i>G4HEPEvtInterface</i></b> |
---|
111 | <p> |
---|
112 | Unfortunately, almost all event generators presently in use, |
---|
113 | commonly are written in FORTRAN. For Geant4, it was decided to not link with any |
---|
114 | FORTRAN program or library, even though the C++ language syntax itself |
---|
115 | allows such a link. Linking to a FORTRAN package might be convenient in |
---|
116 | some cases, but we will lose many advantages of object-oriented |
---|
117 | features of C++, such as robustness. Instead, Geant4 provides |
---|
118 | an ASCII file interface for such event generators. |
---|
119 | <p> |
---|
120 | <i>G4HEPEvtInterface</i> is one of <i>G4VPrimaryGenerator</i> |
---|
121 | concrete class and thus it can be used in your |
---|
122 | <i>G4VUserPrimaryGeneratorAction</i> concrete class. |
---|
123 | <i>G4HEPEvtInterface</i> reads an ASCII file produced by |
---|
124 | an event generator and reproduces <i>G4PrimaryParticle</i> objects |
---|
125 | associated with a <i>G4PrimaryVertex</i> object. It reproduces |
---|
126 | a full production chain of the event generator, starting with |
---|
127 | primary quarks, etc. In other words, <i>G4HEPEvtInterface</i> |
---|
128 | converts information stored in the <tt>/HEPEVT/</tt> common block to |
---|
129 | an object-oriented data structure. Because the <tt>/HEPEVT/</tt> common block |
---|
130 | is commonly used by almost all event generators written in |
---|
131 | FORTRAN, <i>G4HEPEvtInterface</i> can interface to almost |
---|
132 | all event generators currently used in the HEP community. |
---|
133 | The constructor of <i>G4HEPEvtInterface</i> takes the file name. |
---|
134 | Source listing 3.6.1 shows an example how to use <i>G4HEPEvtInterface</i>. |
---|
135 | Note that an event generator is not assumed to give a place of the |
---|
136 | primary particles, the interaction point must be set before |
---|
137 | invoking <i>GeneratePrimaryVertex()</i> method. |
---|
138 | <p> |
---|
139 | <center> |
---|
140 | <table border=2 cellpadding=10> |
---|
141 | <tr> |
---|
142 | <td> |
---|
143 | <PRE> |
---|
144 | <FONT FACE="Courier New"> |
---|
145 | #ifndef ExN04PrimaryGeneratorAction_h |
---|
146 | #define ExN04PrimaryGeneratorAction_h 1 |
---|
147 | |
---|
148 | #include "G4VUserPrimaryGeneratorAction.hh" |
---|
149 | #include "globals.hh" |
---|
150 | |
---|
151 | class G4VPrimaryGenerator; |
---|
152 | class G4Event; |
---|
153 | |
---|
154 | class ExN04PrimaryGeneratorAction : public G4VUserPrimaryGeneratorAction |
---|
155 | { |
---|
156 | public: |
---|
157 | ExN04PrimaryGeneratorAction(); |
---|
158 | ~ExN04PrimaryGeneratorAction(); |
---|
159 | |
---|
160 | public: |
---|
161 | void GeneratePrimaries(G4Event* anEvent); |
---|
162 | |
---|
163 | private: |
---|
164 | G4VPrimaryGenerator* HEPEvt; |
---|
165 | }; |
---|
166 | |
---|
167 | #endif |
---|
168 | |
---|
169 | |
---|
170 | #include "ExN04PrimaryGeneratorAction.hh" |
---|
171 | |
---|
172 | #include "G4Event.hh" |
---|
173 | #include "G4HEPEvtInterface.hh" |
---|
174 | |
---|
175 | ExN04PrimaryGeneratorAction::ExN04PrimaryGeneratorAction() |
---|
176 | { |
---|
177 | HEPEvt = new G4HEPEvtInterface("pythia_event.data"); |
---|
178 | } |
---|
179 | |
---|
180 | ExN04PrimaryGeneratorAction::~ExN04PrimaryGeneratorAction() |
---|
181 | { |
---|
182 | delete HEPEvt; |
---|
183 | } |
---|
184 | |
---|
185 | void ExN04PrimaryGeneratorAction::GeneratePrimaries(G4Event* anEvent) |
---|
186 | { |
---|
187 | HEPEvt->SetParticlePosition(G4ThreeVector(0.*cm,0.*cm,0.*cm)); |
---|
188 | HEPEvt->GeneratePrimaryVertex(anEvent); |
---|
189 | } |
---|
190 | </font> |
---|
191 | </pre> |
---|
192 | </td> |
---|
193 | </tr> |
---|
194 | <tr> |
---|
195 | <td align=center> |
---|
196 | Source listing 3.6.1<BR> |
---|
197 | An example code for <i>G4HEPEvtInterface</i>. |
---|
198 | </td> |
---|
199 | </tr> |
---|
200 | </table></center> |
---|
201 | |
---|
202 | |
---|
203 | <p> |
---|
204 | <b>Format of the ASCII file</b> |
---|
205 | <p> |
---|
206 | An ASCII file, which will be fed by <i>G4HEPEvtInterface</i> |
---|
207 | should have the following format. |
---|
208 | <ul> |
---|
209 | <li>The first line of each primary event should be an integer |
---|
210 | which represents the number of the following lines of primary |
---|
211 | particles. |
---|
212 | <li>Each line in an event corresponds to a particle in the |
---|
213 | <tt>/HEPEVT/</tt> common. Each line has <tt>ISTHEP, IDHEP, JDAHEP(1), JDAHEP(2), |
---|
214 | PHEP(1), PHEP(2), PHEP(3), PHEP(5).</tt> Refer to the <tt>/HEPEVT/</tt> manual for |
---|
215 | the meanings of these variables. |
---|
216 | </ul> |
---|
217 | Source listing 3.6.2 shows an example FORTRAN code to generate an ASCII file. |
---|
218 | <p> |
---|
219 | <center> |
---|
220 | <table border=2 cellpadding=10> |
---|
221 | <tr> |
---|
222 | <td> |
---|
223 | <pre> |
---|
224 | *********************************************************** |
---|
225 | SUBROUTINE HEP2G4 |
---|
226 | * |
---|
227 | * Convert /HEPEVT/ event structure to an ASCII file |
---|
228 | * to be fed by G4HEPEvtInterface |
---|
229 | * |
---|
230 | *********************************************************** |
---|
231 | PARAMETER (NMXHEP=2000) |
---|
232 | COMMON/HEPEVT/NEVHEP,NHEP,ISTHEP(NMXHEP),IDHEP(NMXHEP), |
---|
233 | >JMOHEP(2,NMXHEP),JDAHEP(2,NMXHEP),PHEP(5,NMXHEP),VHEP(4,NMXHEP) |
---|
234 | DOUBLE PRECISION PHEP,VHEP |
---|
235 | * |
---|
236 | WRITE(6,*) NHEP |
---|
237 | DO IHEP=1,NHEP |
---|
238 | WRITE(6,10) |
---|
239 | > ISTHEP(IHEP),IDHEP(IHEP),JDAHEP(1,IHEP),JDAHEP(2,IHEP), |
---|
240 | > PHEP(1,IHEP),PHEP(2,IHEP),PHEP(3,IHEP),PHEP(5,IHEP) |
---|
241 | 10 FORMAT(4I10,4(1X,D15.8)) |
---|
242 | ENDDO |
---|
243 | * |
---|
244 | RETURN |
---|
245 | END |
---|
246 | </pre> |
---|
247 | </td> |
---|
248 | </tr> |
---|
249 | <tr> |
---|
250 | <td align=center> |
---|
251 | Source listing 3.6.2<BR> |
---|
252 | A FORTRAN example using the <tt>/HEPEVT/</tt> common. |
---|
253 | </td> |
---|
254 | </tr> |
---|
255 | </table></center> |
---|
256 | <p> |
---|
257 | <b>Future interface to the new generation generators</b> |
---|
258 | <p> |
---|
259 | Several activities have already been started for developing |
---|
260 | object-oriented event generators. Such new generators can |
---|
261 | be easily linked and used with a Geant4 based simulation. |
---|
262 | Furthermore, we need not distinguish a primary generator from |
---|
263 | the physics processes used in Geant4. Future generators can be |
---|
264 | a kind of physics process plugged-in by inheriting |
---|
265 | <i>G4VProcess</i>. |
---|
266 | <p> |
---|
267 | |
---|
268 | <HR> |
---|
269 | <a name="3.6.3"> |
---|
270 | <H2>3.6.3 Event overlap using multiple generators</H2></a> |
---|
271 | |
---|
272 | Your <i>G4VUserPrimaryGeneratorAction</i> concrete class can |
---|
273 | have more than one <i>G4VPrimaryGenerator</i> concrete class. |
---|
274 | Each <i>G4VPrimaryGenerator</i> concrete class can be accessed |
---|
275 | more than once per event. Using these class objects, one |
---|
276 | event can have more than one primary event. |
---|
277 | <p> |
---|
278 | One possible use is the following. Within an event, |
---|
279 | a <i>G4HEPEvtInterface</i> class object instantiated with a |
---|
280 | minimum bias event file is accessed 20 times and another |
---|
281 | <i>G4HEPEvtInterface</i> class object instantiated with a |
---|
282 | signal event file is accessed once. Thus, this event represents |
---|
283 | a typical signal event of LHC overlapping 20 minimum bias |
---|
284 | events. It should be noted that a simulation of event |
---|
285 | overlapping can be done by merging hits and/or digits |
---|
286 | associated with several events, and these events can be |
---|
287 | simulated independently. Digitization over multiple events |
---|
288 | will be mentioned in <a href="../Detector/digitization.html">Section 4.5</a>. |
---|
289 | |
---|
290 | <BR><BR> |
---|
291 | <HR> |
---|
292 | <A HREF="../../../../Authors/html/subjectsToAuthors.html"> |
---|
293 | <I>About the authors</I></A> |
---|
294 | |
---|
295 | </BODY> |
---|
296 | </HTML> |
---|