1 | <HTML> |
---|
2 | <TITLE> |
---|
3 | </TITLE> |
---|
4 | <!-- Changed by: Dennis Wright, 25-Jun-2002 --> |
---|
5 | |
---|
6 | <BODY> |
---|
7 | <TABLE WIDTH="100%"><TR> |
---|
8 | <TD> |
---|
9 | |
---|
10 | |
---|
11 | <A HREF="index.html"> |
---|
12 | <IMG SRC="../../../../resources/html/IconsGIF/Contents.gif" ALT="Contents"></A> |
---|
13 | <A HREF="eventGenerator.html"> |
---|
14 | <IMG SRC="../../../../resources/html/IconsGIF/Previous.gif" ALT="Previous"></A> |
---|
15 | <IMG SRC="../../../../resources/html/IconsGIF/NextGR.gif" ALT="Next"></A> |
---|
16 | </TD> |
---|
17 | <TD ALIGN="Right"> |
---|
18 | <FONT SIZE="-1" COLOR="#238E23"> |
---|
19 | <B>Geant4 User's Guide</B> |
---|
20 | <BR> |
---|
21 | <B>For Application Developers</B> |
---|
22 | <BR> |
---|
23 | <B>Toolkit Fundamentals</B> |
---|
24 | </FONT> |
---|
25 | </TD> |
---|
26 | </TR></TABLE> |
---|
27 | <BR><BR> |
---|
28 | |
---|
29 | <P ALIGN="Center"> |
---|
30 | <FONT SIZE="+3" COLOR="#238E23"> |
---|
31 | <B>3.7 Event Biasing Techniques</B> |
---|
32 | </FONT> |
---|
33 | <BR><BR> |
---|
34 | |
---|
35 | <HR ALIGN="Center" SIZE="7%"> |
---|
36 | |
---|
37 | <p> |
---|
38 | <a name="3.7.1"> |
---|
39 | <H2>3.7.1 Scoring, Geometrical Importance Sampling and Weight Roulette</H2></a> |
---|
40 | |
---|
41 | Geant4 provides event biasing techniques which may be used to save computing time |
---|
42 | in such applications as the simulation of radiation shielding. These are |
---|
43 | <I>geometrical splitting</I> and <I>Russian roulette</I> (also called geometrical |
---|
44 | importance sampling), and <I>weight roulette</I>. A basic scoring system is also |
---|
45 | provided in order to monitor the sampling. In this chapter, it is assumed that |
---|
46 | the reader is familiar with both the usage of Geant4 and the concepts of importance |
---|
47 | sampling. More detailed documentation may be found at |
---|
48 | <a href="#imp-biasing1">[1]</a>. |
---|
49 | A detailed description of different use-cases which employ the sampling and scoring |
---|
50 | techniques can be found at <a href="#imp-biasing2">[2]</a>. |
---|
51 | <BR> |
---|
52 | The goals of the scoring technique are to: |
---|
53 | <UL> |
---|
54 | <LI> appraise particle quantities related to special regions or surfaces, |
---|
55 | <LI> be applicable to all "cells" (physical volumes or replicas) of a given geometry, |
---|
56 | <LI> be customizable. |
---|
57 | </UL> |
---|
58 | Standard scoring must be provided for quantities such as tracks entering a |
---|
59 | cell, average weight of entering tracks, energy of entering tracks, and |
---|
60 | collisions inside the cell. |
---|
61 | <BR> |
---|
62 | The purpose of importance sampling is to save computing time by sampling less |
---|
63 | often the particle histories entering "less important" geometry regions, and |
---|
64 | more often in more "important" regions. Given the same amount of computing |
---|
65 | time, an importance-sampled and an analogue-sampled simulation must show |
---|
66 | equal mean values, while the importance-sampled simulation will have a |
---|
67 | decreased variance. |
---|
68 | <BR> |
---|
69 | The implementation of scoring is independent of the implementation of |
---|
70 | importance sampling. However both share common concepts. <I>Scoring and |
---|
71 | importance sampling apply to particle types chosen by the user</I>. |
---|
72 | <BR> |
---|
73 | Examples on how to use scoring and importance sampling may be found |
---|
74 | in <tt>examples/extended/biasing</tt> and <tt>examples/advanced/Tiara</tt>. |
---|
75 | |
---|
76 | <h3><a name="geometries"></a>3.7.1.1 Geometries</h3> |
---|
77 | |
---|
78 | The kind of scoring referred to in this note and the importance sampling |
---|
79 | apply to spatial cells provided by the user. |
---|
80 | <P> |
---|
81 | A <B>cell</B> is a physical volume (further specified by it's replica number, |
---|
82 | if the volume is a replica). Cells may be defined in two kinds of |
---|
83 | geometries: |
---|
84 | <OL> |
---|
85 | <LI><B>mass geometry</B>: the geometry setup of the experiment to be |
---|
86 | simulated. Physics processes apply to this geometry.</LI> |
---|
87 | <LI><B>parallel-geometry</B>: a geometry constructed to define the physical |
---|
88 | volumes according to which scoring and/or importance sampling is applied.</LI> |
---|
89 | </OL> |
---|
90 | The user has the choice to score and/or sample by importance the particles |
---|
91 | of the chosen type, according to mass geometry or to parallel geometry. |
---|
92 | It is possible to utilize several parallel geometries in addition to the mass |
---|
93 | geometry. This provides the user with a lot of flexibility to define separate |
---|
94 | geometries for different particle types in order to apply scoring or/and |
---|
95 | importance sampling. |
---|
96 | <P> |
---|
97 | <B><a name="remarkparallel"></a>Limitations of "parallel" geometries</B> |
---|
98 | <UL> |
---|
99 | <LI>A special kind of "transportation" inside the "parallel" geometry has been |
---|
100 | developed. In the current implementation this does not take fields into |
---|
101 | account. Therefore "parallel" geometries can only be used in field-free |
---|
102 | simulations.<BR> |
---|
103 | The world volume of the parallel geometry must overlap the world volume |
---|
104 | of the mass geometry.<BR> |
---|
105 | Particles crossing a boundary in the parallel geometry where there is |
---|
106 | vacuum in the mass geometry are also biased. This may be optimized in |
---|
107 | later versions.</LI> |
---|
108 | <LI>Scoring cells must not share boundaries with the world volume.</LI> |
---|
109 | </UL> |
---|
110 | |
---|
111 | <h3><a name="sampler"></a>3.7.1.2 Changing the Sampling</h3> |
---|
112 | <p> |
---|
113 | Samplers are higher level tools which perform the necessary changes |
---|
114 | of the Geant4 sampling in order to apply importance sampling and weight |
---|
115 | roulette. The samplers may also be used to apply scoring.</p> |
---|
116 | |
---|
117 | <p>Scoring and variance reduction may be combined arbitrarily |
---|
118 | for chosen particle types and may be applied to the mass or to parallel |
---|
119 | geometries. The samplers support all combinations.</p> |
---|
120 | |
---|
121 | <p>Different samplers are implemented for the mass and parallel |
---|
122 | geometries. Both implement the interface <tt>G4VSampler</tt>. |
---|
123 | <tt>G4VSampler</tt> allows the chosen combinations of scoring and variance |
---|
124 | reduction to be prepared and configures the sampling appropriately. To this |
---|
125 | end <tt>G4VSampler</tt> provides <tt>Prepare...</tt> methods and a |
---|
126 | <tt>Configure</tt> method:</p> |
---|
127 | |
---|
128 | <a name="G4VSampler"> |
---|
129 | <PRE> |
---|
130 | class G4VSampler |
---|
131 | { |
---|
132 | public: |
---|
133 | G4VSampler(); |
---|
134 | virtual ~G4VSampler(); |
---|
135 | virtual void PrepareScoring(G4VScorer *Scorer) = 0; |
---|
136 | virtual void PrepareImportanceSampling(G4VIStore *istore, |
---|
137 | const G4VImportanceAlgorithm |
---|
138 | *ialg = 0) = 0; |
---|
139 | virtual void PrepareWeightRoulett(G4double wsurvive = 0.5, |
---|
140 | G4double wlimit = 0.25, |
---|
141 | G4double isource = 1) = 0; |
---|
142 | virtual void PrepareWeightWindow(G4VWeightWindowStore *wwstore, |
---|
143 | G4VWeightWindowAlgorithm *wwAlg = 0, |
---|
144 | G4PlaceOfAction placeOfAction = |
---|
145 | onBoundary) = 0; |
---|
146 | virtual void Configure() = 0; |
---|
147 | virtual void ClearSampling() = 0; |
---|
148 | virtual G4bool IsConfigured() const = 0; |
---|
149 | }; |
---|
150 | </PRE> |
---|
151 | |
---|
152 | <p> |
---|
153 | The methods for setting up the desired combination need specific information: |
---|
154 | </p> |
---|
155 | <div class="itemizedlist"> |
---|
156 | <ul type="disc"> |
---|
157 | <li>Scoring: message <tt>PrepareScoring</tt> with a |
---|
158 | <a href="#G4VScorer"><tt>G4VScorer</tt></a></li> |
---|
159 | <li>Importance sampling: message <tt>PrepareImportanceSampling</tt> with a |
---|
160 | <a href="#G4VIStore"><tt>G4VIStore</tt></a> and optionally a |
---|
161 | <tt>G4VImportanceAlgorithm</tt></li> |
---|
162 | <li><p>Weight window: message <tt>PrepareWeightWindow</tt> with the |
---|
163 | arguments:</p> |
---|
164 | <div class="itemizedlist"> |
---|
165 | <ul type="round"> |
---|
166 | <li><p><span class="emphasis"><em>*wwstore</em></span>: a |
---|
167 | <tt>G4VWeightWindowStore</tt> for retrieving the lower weight |
---|
168 | bound s for the energy-space cells</p></li> |
---|
169 | <li><p><span class="emphasis"><em>*wwAlg</em></span>: a |
---|
170 | <tt>G4VWeightWindowAlgorithm</tt> |
---|
171 | if a customized algorithm should be used</p></li> |
---|
172 | <li><p><span class="emphasis"><em>placeOfAction</em></span>: a |
---|
173 | <tt>G4PlaceOfAction</tt> specifying where to perform the |
---|
174 | biasing</p></li> |
---|
175 | </ul></div></li> |
---|
176 | <li>Weight roulette: message <tt>PrepareWeightRoulett</tt> with the |
---|
177 | optional parameters: |
---|
178 | <div class="itemizedlist"> |
---|
179 | <ul type="round"> |
---|
180 | <li><span class="emphasis"><em>wsurvive</em></span>: survival weight</li> |
---|
181 | <li><span class="emphasis"><em>wlimit</em></span>: minimal allowed value of |
---|
182 | weight * source importance / cell importance</li> |
---|
183 | <li><span class="emphasis"><em>isource</em></span>: importance of the |
---|
184 | source cell</li> |
---|
185 | </ul> |
---|
186 | </div></li> |
---|
187 | </ul> |
---|
188 | </div> |
---|
189 | <p> |
---|
190 | Each object of a sampler class is responsible for one particle type. |
---|
191 | The particle type is given to the constructor of the sampler classes via the |
---|
192 | particle type name, e.g. "neutron". Depending on the specific purpose, |
---|
193 | the <tt>Configure()</tt> of a sampler will set up specialized processes |
---|
194 | (derived from <tt>G4VProcess</tt>) for transportation in the parallel geometry, |
---|
195 | scoring, importance sampling and weight roulette for the given particle type. |
---|
196 | When <tt>Configure()</tt> is invoked the sampler places the processes in the |
---|
197 | correct order independent of the order in which user invoked the |
---|
198 | <tt>Prepare...</tt> methods.</p> |
---|
199 | <div class="note" style="margin-left: 0.5in; margin-right: 0.5in;"> |
---|
200 | <h3 class="title">Notes</h3> |
---|
201 | <ul> |
---|
202 | <li>The <tt>Prepare...()</tt> functions may each only be invoked once.</li> |
---|
203 | <li>To configure the sampling the function <tt>Configure()</tt> must be |
---|
204 | called <span class="emphasis"><em>after</em></span> the |
---|
205 | <tt>G4RunManager</tt> has been initialized.</li> |
---|
206 | </ul> |
---|
207 | </div></div> |
---|
208 | <p>Two classes implementing the interface |
---|
209 | <a href="#G4VSampler"><tt>G4VSampler</tt></a> are provided for the mass and |
---|
210 | a parallel geometry, respectively:</p> |
---|
211 | <div class="itemizedlist"> |
---|
212 | <ul type="disc"> |
---|
213 | <li><tt>G4MassGeometrySampler</tt></li> |
---|
214 | <li><tt>G4ParallelGeometrySampler</tt></li> |
---|
215 | </ul> |
---|
216 | </div> |
---|
217 | <p> |
---|
218 | The constructors of both classes take the particle type name (e.g. "neutron") |
---|
219 | as an argument. In addition, the constructor of |
---|
220 | <tt>G4ParallelGeometrySampler</tt> needs a reference to the physical world |
---|
221 | volume of the parallel geometry.</p> |
---|
222 | |
---|
223 | <h3><a name="scoring"></a>3.7.1.3 Scoring</h3> |
---|
224 | |
---|
225 | <p> |
---|
226 | Scoring is provided by a framework and is done according to particle type. |
---|
227 | Nevertheless it is also possible to score particles of different |
---|
228 | types into the same score. The framework may also be easily used for |
---|
229 | customized scoring.</p> |
---|
230 | <p> Scoring may be applied to a mass or a parallel geometry. It is done with |
---|
231 | an object genericly called a scorer using a sampler described above. The scorer |
---|
232 | receives the information about every step taken by a particle of chosen type. |
---|
233 | This information consists of an object of the Geant4 kernel class |
---|
234 | <tt>G4Step</tt> and an object of the class <tt>G4GeometryCellStep</tt> provided |
---|
235 | specifically for the purpose of scoring and importance sampling. |
---|
236 | <tt>G4GeometryCellStep</tt> provides information about the previous and current |
---|
237 | "cell" of the particle track. |
---|
238 | </p> |
---|
239 | |
---|
240 | <p> |
---|
241 | A "scorer" class derives from the interface <tt>G4VScorer</tt>. Users may |
---|
242 | create customized "scorers" or use the standard scoring. |
---|
243 | </p> |
---|
244 | <p> |
---|
245 | Classes involved in the framework:</p> |
---|
246 | <div class="itemizedlist"> |
---|
247 | <ul type="disc"> |
---|
248 | <li><a name="G4VScorer"></a><tt>G4VScorer</tt><br> |
---|
249 | Classes derived from <tt>G4VScorer</tt> may use the Geant4 framework for |
---|
250 | scoring. This is the minimum dependence a customized scoring has to have |
---|
251 | to use the scoring frame. The frame will provide a process messaging a |
---|
252 | scorer derived from <tt>G4VScorer</tt> and classes for the setup.<br> |
---|
253 | In particular an object of this type may be given to a <tt>G4VSampler</tt>. |
---|
254 | In order to score different particle types in one object of <tt>G4VScorer</tt> |
---|
255 | the object may be given to different sampler objects.<br> |
---|
256 | The interface: |
---|
257 | <pre class="programlisting"> |
---|
258 | class G4VScorer |
---|
259 | { |
---|
260 | public: |
---|
261 | G4VScorer(); |
---|
262 | virtual ~G4VScorer(); |
---|
263 | virtual void Score(const G4Step &step, const G4GeometryCellStep &gstep) = 0; |
---|
264 | }; |
---|
265 | </pre> |
---|
266 | The member function |
---|
267 | <tt>Score(const G4Step &step, const G4GeometryCellStep &gstep)</tt> |
---|
268 | is invoked for every "PostStep" of a chosen particle before the |
---|
269 | "<tt>PostStepDoIt()</tt>"-functions of the physical processes are invoked.</li> |
---|
270 | <li><tt>G4GeometryCellStep</tt>: |
---|
271 | <pre class="programlisting"> |
---|
272 | class G4GeometryCellStep |
---|
273 | { |
---|
274 | public: |
---|
275 | G4GeometryCellStep(const G4GeometryCell &preCell, |
---|
276 | const G4GeometryCell &postCell); |
---|
277 | ~G4GeometryCellStep(); |
---|
278 | const G4GeometryCell &GetPreGeometryCell() const; |
---|
279 | const G4GeometryCell &GetPostGeometryCell() const; |
---|
280 | G4bool GetCrossBoundary() const; |
---|
281 | void SetPreGeometryCell(const G4GeometryCell &preCell); |
---|
282 | void SetPostGeometryCell(const G4GeometryCell &postCell); |
---|
283 | void SetCrossBoundary(G4bool b); |
---|
284 | private: .... |
---|
285 | }; |
---|
286 | </pre> |
---|
287 | |
---|
288 | The classes and fuctions used in <tt>G4GeometryCellStep</tt>: |
---|
289 | <div class="itemizedlist"> |
---|
290 | <ul type="round"> |
---|
291 | <li><tt>G4GeometryCell()</tt> identifies a "cell" as a physical volume with |
---|
292 | a replica number. The "cell" is in the "parallel" geometry if a sampler |
---|
293 | for "parallel" geometries is used, otherwise it is in the mass geometry.</li> |
---|
294 | <li><p><tt>GetPreGeometryCell()</tt> returns the previous "cell" touched by the |
---|
295 | particle, or it is equal to <tt>GetPostGeometryCell()</tt> if the track |
---|
296 | has not crossed a boundary.</li> |
---|
297 | <li><tt>GetPostGeometryCell()</tt> refers to the current "cell".</li> |
---|
298 | <li><tt>GetCrossBoundary()</tt> returns <tt>true</tt> if the particle |
---|
299 | crosses a boundary in the parallel or mass geometry with respect to the |
---|
300 | sampler used.</li> |
---|
301 | </ul> |
---|
302 | </div></li> |
---|
303 | </ul> |
---|
304 | </div> |
---|
305 | |
---|
306 | <div class="note" style="margin-left: 0.5in; margin-right: 0.5in;"> |
---|
307 | <h3 class="title">Note</h3> |
---|
308 | <ul> |
---|
309 | <li>When scoring is done in a "parallel" geometry, special action must be taken |
---|
310 | to prevent counting of "collisions" with boundaries of the mass geometry as |
---|
311 | interactions. This is handled differently when scoring is done in the mass |
---|
312 | geometry.</li> |
---|
313 | </ul> |
---|
314 | </div> |
---|
315 | |
---|
316 | <h3><a name="importancesampling"></a>3.7.1.4 Importance Sampling</h3> |
---|
317 | |
---|
318 | |
---|
319 | <p>Importance sampling acts on particles crossing boundaries between |
---|
320 | "importance cells". The action taken depends on the importance values |
---|
321 | assigned to the cells. In general a particle history is either split |
---|
322 | or Russian roulette is played if the importance increases or decreases, |
---|
323 | respectively. A weight assigned to the history is changed according to |
---|
324 | the action taken.</p> |
---|
325 | |
---|
326 | <p>The tools provided for importance sampling require the user |
---|
327 | to have a good understanding of the physics in the problem. This is |
---|
328 | because the user has to decide which particle types require importance |
---|
329 | sampled, define the cells, and assign importance values to the cells. If |
---|
330 | this is not done properly the results cannot be expected to describe a |
---|
331 | real experiment.</p> |
---|
332 | |
---|
333 | <p>The assignment of importance values to a cell is done using an |
---|
334 | importance store described below.</p> |
---|
335 | |
---|
336 | <p>An "importance store" with the interface <tt>G4VIStore</tt> is used |
---|
337 | to store importance values related to cells. In order to do importance |
---|
338 | sampling the user has to create an object (e.g. of class <tt>G4IStore</tt>) |
---|
339 | of type <tt>G4VIStore</tt>. The samplers may be given a <tt>G4VIStore</tt>. |
---|
340 | The user fills the store with cells and their importance values.</p> |
---|
341 | |
---|
342 | <p>An importance store has to be constructed with a reference to |
---|
343 | the world volume of the geometry used for importance sampling. This may |
---|
344 | be the world volume of the mass or of a parallel geometry. Importance |
---|
345 | stores derive from the interface <tt>G4VIStore</tt>:</p> |
---|
346 | <a name="G4VIStore"></a> |
---|
347 | <pre class="programlisting"> |
---|
348 | class G4VIStore |
---|
349 | { |
---|
350 | public: |
---|
351 | G4VIStore(); |
---|
352 | virtual ~G4VIStore(); |
---|
353 | virtual G4double GetImportance(const G4GeometryCell &gCell) const = 0; |
---|
354 | virtual G4bool IsKnown(const G4GeometryCell &gCell) const = 0; |
---|
355 | virtual const G4VPhysicalVolume &GetWorldVolume() const = 0; |
---|
356 | }; |
---|
357 | </pre> |
---|
358 | |
---|
359 | <p> |
---|
360 | A concrete implementation of an importance store is provided by the class |
---|
361 | <tt>G4VStore</tt>. The <span class="emphasis"><em>public</em></span> |
---|
362 | part of the class is:</p> |
---|
363 | <pre class="programlisting"> |
---|
364 | class G4IStore : public G4VIStore |
---|
365 | { |
---|
366 | public: |
---|
367 | explicit G4IStore(const G4VPhysicalVolume &worldvolume); |
---|
368 | virtual ~G4IStore(); |
---|
369 | virtual G4double GetImportance(const G4GeometryCell &gCell) const; |
---|
370 | virtual G4bool IsKnown(const G4GeometryCell &gCell) const; |
---|
371 | virtual const G4VPhysicalVolume &GetWorldVolume() const; |
---|
372 | void AddImportanceGeometryCell(G4double importance, |
---|
373 | const G4GeometryCell &gCell); |
---|
374 | void AddImportanceGeometryCell(G4double importance, |
---|
375 | const G4VPhysicalVolume &, |
---|
376 | G4int aRepNum = 0); |
---|
377 | void ChangeImportance(G4double importance, |
---|
378 | const G4GeometryCell &gCell); |
---|
379 | void ChangeImportance(G4double importance, |
---|
380 | const G4VPhysicalVolume &, |
---|
381 | G4int aRepNum = 0); |
---|
382 | G4double GetImportance(const G4VPhysicalVolume &, |
---|
383 | G4int aRepNum = 0) const ; |
---|
384 | private: ..... |
---|
385 | }; |
---|
386 | </pre> |
---|
387 | <p> |
---|
388 | The member function <tt>AddImportanceGeometryCell()</tt> enters a cell and |
---|
389 | an importance value into the importance store. The importance values may be |
---|
390 | returned either according to a physical volume and a replica number or |
---|
391 | according to a <tt>G4GeometryCell</tt>. The user mus be aware of the |
---|
392 | interpretation of assigning importance values to a cell.</p> |
---|
393 | |
---|
394 | <div class="note" style="margin-left: 0.5in; margin-right: 0.5in;"> |
---|
395 | <h3 class="title">Notes</h3> |
---|
396 | <ul> |
---|
397 | <li>An importance value must be assigned to every cell.</li> |
---|
398 | </ul> |
---|
399 | </div> |
---|
400 | |
---|
401 | <p>The different cases:</p> |
---|
402 | <div class="itemizedlist"> |
---|
403 | <ul type="disc"> |
---|
404 | <li><span class="emphasis"><em>Cell is not in store</em></span><br> |
---|
405 | Not filling a certain cell in the store will cause an exception.</li> |
---|
406 | <li><span class="emphasis"><em>Importance value = zero</em></span><br> |
---|
407 | Tracks of the chosen particle type will be killed.</li> |
---|
408 | <li><span class="emphasis"><em>importance values > 0</em></span><br> |
---|
409 | Normal allowed values</li> |
---|
410 | <li><span class="emphasis"><em>Importance value smaller zero</em></span><br> |
---|
411 | Not allowed!</li> |
---|
412 | </ul> |
---|
413 | </div></div> |
---|
414 | |
---|
415 | <h3><a name="importancealgorithm"></a>3.7.1.5 The Importance Sampling Algorithm</h3> |
---|
416 | |
---|
417 | <p> |
---|
418 | Importance sampling supports using a customized importance sampling algorithm. |
---|
419 | To this end, the sampler interface |
---|
420 | <a href="#G4VSampler"><tt>G4VSampler</tt></a> may be given a pointer to |
---|
421 | the interface <tt>G4VImportanceAlgorithm</tt>:</p> |
---|
422 | <a name="G4VImportanceAlgorithm"></a> |
---|
423 | <pre class="programlisting"> |
---|
424 | class G4VImportanceAlgorithm |
---|
425 | { |
---|
426 | public: |
---|
427 | G4VImportanceAlgorithm(); |
---|
428 | virtual ~G4VImportanceAlgorithm(); |
---|
429 | virtual G4Nsplit_Weight Calculate(G4double ipre, |
---|
430 | G4double ipost, |
---|
431 | G4double init_w) const = 0; |
---|
432 | }; |
---|
433 | </pre> |
---|
434 | |
---|
435 | <p>The method <tt>Calculate()</tt> takes the arguments:</p> |
---|
436 | <div class="itemizedlist"> |
---|
437 | <ul type="disc"> |
---|
438 | <li><span class="emphasis"><em>ipre, ipost</em></span>: importance |
---|
439 | of the previous cell and the importance of the current cell, respectively. |
---|
440 | </li> |
---|
441 | <li><span class="emphasis"><em>init_w</em></span>: the particles weight</li> |
---|
442 | </ul> |
---|
443 | </div> |
---|
444 | |
---|
445 | <p>It returns the struct:</p> |
---|
446 | <pre class="programlisting"> |
---|
447 | class G4Nsplit_Weight |
---|
448 | { |
---|
449 | public: |
---|
450 | |
---|
451 | G4int fN; |
---|
452 | G4double fW; |
---|
453 | }; |
---|
454 | </pre> |
---|
455 | |
---|
456 | <div class="itemizedlist"> |
---|
457 | <ul type="disc"> |
---|
458 | <li><span class="emphasis"><em>fN</em></span>: the calculated number of |
---|
459 | particles to exit the importance sampling</li> |
---|
460 | <li><span class="emphasis"><em>fW</em></span>: the weight of the particles</li> |
---|
461 | </ul> |
---|
462 | </div> |
---|
463 | |
---|
464 | <p> |
---|
465 | The user may have a customized algorithm used by providing a class |
---|
466 | inheriting from <tt>G4VImportanceAlgorithm</tt>. <br> |
---|
467 | If no customized algorithm is given to the sampler the default importance sampling |
---|
468 | algorithm is used. This algorithm is implemented in <tt>G4ImportanceAlgorithm</tt>.</p> |
---|
469 | |
---|
470 | <h3><a name="weighwindow"></a>3.7.1.6 The Weight Window Technique</h3> |
---|
471 | |
---|
472 | <p> The weight window technique is a weight-based alternative to |
---|
473 | importance sampling:</p> |
---|
474 | <ul type="disc"> |
---|
475 | <li><p>applies splitting and Russian roulette depending on space (cells) and |
---|
476 | energy</p></li> |
---|
477 | <li><p>user defines weight windows in contrast to defining |
---|
478 | importance values as in importance sampling</p></li></ul> |
---|
479 | <p> |
---|
480 | In contrast to importance sampling this technique is not weight blind. Instead |
---|
481 | the technique is applied according to the particle weight with respect to the |
---|
482 | current energy-space cell.</p> |
---|
483 | <p> |
---|
484 | Therefore the technique is convenient to apply in combination with other |
---|
485 | variance reduction techniques such as cross-section biasing and implicit capture.</p> |
---|
486 | |
---|
487 | <p> |
---|
488 | A weight window may be specified for every cell and for several energy regions: |
---|
489 | <span class="emphasis"><em>space-energy cell</em> |
---|
490 | </span>.</p> |
---|
491 | <div class="informalfigure"><div class="mediaobject"><table border="0" summary="manufactured viewport for HTML img" cellspacing="0" cellpadding="0"><tr><td><img src="wwconcept.gif" alt="Weight window concept
"></td></tr></table><div class="caption"><p>Weight window concept</p></div></div></div> |
---|
492 | <p>The user specifies a <span class="emphasis"><em>lower weight bound W_L |
---|
493 | </em></span> for every space-energy cell.</p> |
---|
494 | <ul type="disc"><li><p>The upper weight bound W_U and the survival weight |
---|
495 | W_S are calculated as: </p><p>W_U = C_U <span class="emphasis"><em>W_L |
---|
496 | </em></span></p><p>and</p><p>W_S = C_S <span class="emphasis"> |
---|
497 | <em>W_L</em></span>.</p></li> |
---|
498 | <li><p>The user specifies C_S and C_U once for the whole problem.</p></li> |
---|
499 | <li><p>The user may give different sets of energy bounds for |
---|
500 | every cell or one set for all geometrical cells </p></li> |
---|
501 | <li><p>Special case: if C_S = C_U = 1 for all energies then |
---|
502 | weight window is equivalent to importance sampling</p></li> |
---|
503 | <li><p>The user can choose to apply the technique: at boundaries, |
---|
504 | on collisions or on boundaries and collisions</p></li></ul> |
---|
505 | |
---|
506 | |
---|
507 | <p> |
---|
508 | The energy-space cells are realized by <tt>G4GeometryCell</tt> as in importance |
---|
509 | sampling. The cells are stored in a weight window store defined by |
---|
510 | <tt>G4VWeightWindowStore</tt>:</p> |
---|
511 | |
---|
512 | <a name="G4VWeightWindowStore"></a><pre class="programlisting"> |
---|
513 | class G4VWeightWindowStore { |
---|
514 | public: |
---|
515 | G4VWeightWindowStore(); |
---|
516 | virtual ~G4VWeightWindowStore(); |
---|
517 | virtual G4double GetLowerWeitgh(const G4GeometryCell &gCell, |
---|
518 | G4double partEnergy) const = 0; |
---|
519 | virtual G4bool IsKnown(const G4GeometryCell &gCell) const = 0; |
---|
520 | virtual const G4VPhysicalVolume &GetWorldVolume() const = 0; |
---|
521 | }; |
---|
522 | </pre> |
---|
523 | <p>A concrete implementation is provided:</p><a name="G4WeightWindowStore"></a><pre class="programlisting"> |
---|
524 | class G4WeightWindowStore: public G4VWeightWindowStore { |
---|
525 | public: |
---|
526 | explicit G4WeightWindowStore(const G4VPhysicalVolume &worldvolume); |
---|
527 | virtual ~G4WeightWindowStore(); |
---|
528 | virtual G4double GetLowerWeitgh(const G4GeometryCell &gCell, |
---|
529 | G4double partEnergy) const; |
---|
530 | virtual G4bool IsKnown(const G4GeometryCell &gCell) const; |
---|
531 | virtual const G4VPhysicalVolume &GetWorldVolume() const; |
---|
532 | void AddLowerWeights(const G4GeometryCell &gCell, |
---|
533 | const std::vector<G4double> &lowerWeights); |
---|
534 | void AddUpperEboundLowerWeightPairs(const G4GeometryCell &gCell, |
---|
535 | const G4UpperEnergyToLowerWeightMap& |
---|
536 | enWeMap); |
---|
537 | void SetGeneralUpperEnergyBounds(const |
---|
538 | std::set<G4double, std::less<G4double> > & enBounds); |
---|
539 | |
---|
540 | private:: |
---|
541 | ... |
---|
542 | }; |
---|
543 | </pre> |
---|
544 | <p> |
---|
545 | The user may choose equal energy bounds for all cells. In this case a set of |
---|
546 | upper energy bounds must be given to the store using the method |
---|
547 | <tt>SetGeneralUpperEnergyBounds</tt>. If a general set of energy bounds have |
---|
548 | been set <tt>AddLowerWeights</tt> can be used to add the cells.</p> |
---|
549 | <p> |
---|
550 | Alternatively, the user may chose different energy regions for different cells. |
---|
551 | In this case the user must provide a mapping of upper energy bounds to lower |
---|
552 | weight bounds for every cell using the method |
---|
553 | <tt>AddUpperEboundLowerWeightPairs</tt>.</p> |
---|
554 | <p> |
---|
555 | Weight window algorithms implementing the interface class |
---|
556 | <tt>G4VWeightWindowAlgorithm</tt> can be used to define a customized algorithm:</p> |
---|
557 | <a name="G4VWeightWindowAlgorithm"></a><pre class="programlisting"> |
---|
558 | class G4VWeightWindowAlgorithm { |
---|
559 | public: |
---|
560 | G4VWeightWindowAlgorithm(); |
---|
561 | virtual ~G4VWeightWindowAlgorithm(); |
---|
562 | virtual G4Nsplit_Weight Calculate(G4double init_w, |
---|
563 | G4double lowerWeightBound) const = 0; |
---|
564 | }; |
---|
565 | </pre> |
---|
566 | <p> |
---|
567 | A concrete implementation is provided and used as a default:</p> |
---|
568 | <a name="G4WeightWindowAlgorithm"></a><pre class="programlisting"> |
---|
569 | class G4WeightWindowAlgorithm : public G4VWeightWindowAlgorithm { |
---|
570 | public: |
---|
571 | G4WeightWindowAlgorithm(G4double upperLimitFaktor = 5, |
---|
572 | G4double survivalFaktor = 3, |
---|
573 | G4int maxNumberOfSplits = 5); |
---|
574 | virtual ~G4WeightWindowAlgorithm(); |
---|
575 | virtual G4Nsplit_Weight Calculate(G4double init_w, |
---|
576 | G4double lowerWeightBound) const; |
---|
577 | private: |
---|
578 | ... |
---|
579 | }; |
---|
580 | </pre> |
---|
581 | <p> |
---|
582 | The constructor takes three parameters which are used to: |
---|
583 | calculate the upper weight bound (upperLimitFaktor), calculate the survival weight |
---|
584 | (survivalFaktor), and introduce a maximal number (maxNumberOfSplits) of copies to |
---|
585 | be created in one go.</p> |
---|
586 | <p> |
---|
587 | In addition, the inverse of the maxNumberOfSplits is used to specify the minimum |
---|
588 | survival probability in case of Russian roulette.</p></div></div> |
---|
589 | |
---|
590 | <h3><a name="weightroulette"></a>3.7.1.7 The Weight Roulette Technique</h3> |
---|
591 | |
---|
592 | <p>Weight roulette (also called weight cutoff) is usually applied if importance |
---|
593 | sampling and implicit capture are used together. Implicit capture is not |
---|
594 | described here but it is useful to note that this procedure reduces a particle |
---|
595 | weight in every collision instead of killing the particle with some probability.</p> |
---|
596 | |
---|
597 | <p> |
---|
598 | Together with importance sampling the weight of a particle may become so low |
---|
599 | that it does not change any result significantly. Hence tracking a very low |
---|
600 | weight particle is a waste of computing time. Weight roulette is applied in |
---|
601 | order to solve this problem.</p> |
---|
602 | |
---|
603 | <b>The weight roulette concept</b> |
---|
604 | <p>Weight roulette takes into account the importance "Ic" of the current cell |
---|
605 | and the importance "Is" of the cell in which the source is located, by using |
---|
606 | the ratio "R=Is/Ic".</p> |
---|
607 | |
---|
608 | <p> |
---|
609 | Weight roulette uses a relative minimal weight limit and a relative survival |
---|
610 | weight. When a particle falls below the weight limit Russian roulette is applied. |
---|
611 | If the particle survives, tracking will be continued and the particle weight will |
---|
612 | be set to the survival weight.</p> |
---|
613 | |
---|
614 | <p> |
---|
615 | The weight roulette uses the following parameters with their default values:</p> |
---|
616 | <div class="itemizedlist"> |
---|
617 | <ul type="disc"> |
---|
618 | <li><p><span class="emphasis"><em>wsurvival</em></span>: 0.5</p></li> |
---|
619 | <li><p><span class="emphasis"><em>wlimit</em></span>: 0.25</p></li> |
---|
620 | <li><p><span class="emphasis"><em>isource</em></span>: 1</p></li> |
---|
621 | </ul> |
---|
622 | </div> |
---|
623 | <p>The following algorithm is applied:</p> |
---|
624 | <p>If a particle weight "w" is lower than R*wlimit:</p> |
---|
625 | <div class="itemizedlist"> |
---|
626 | <ul type="disc"> |
---|
627 | <li> |
---|
628 | <p>the weight of the particle will be changed to "ws = wsurvival*R"</p> |
---|
629 | </li> |
---|
630 | <li> |
---|
631 | <p>the probability for the particle to survive is "p = w/ws"</p> |
---|
632 | </li> |
---|
633 | </ul> |
---|
634 | </div> |
---|
635 | |
---|
636 | <p> |
---|
637 | <a name="3.7.2"> |
---|
638 | <H2>3.7.2 Physics Based Biasing</H2></a> |
---|
639 | Geant4 supports physics based biasing through a number of general use, built in biasing |
---|
640 | techniques. A utility class, G4WrapperProcess, is also available to support user defined biasing. |
---|
641 | |
---|
642 | <h3><a name="3.7.2.1"></a>3.7.2.1 Built in Biasing Options</h3> |
---|
643 | <p> |
---|
644 | <h3><a name="3.7.2.1.1"></a>3.7.2.1.1 Primary Particle Biasing</h3> |
---|
645 | <p> |
---|
646 | Primary particle biasing can be used to increase the number of primary particles generated in a particular phase space region of interest. The weight of the primary particle is modified as appropriate. A general implementation is provided through the G4GeneralParticleSource class. It is possible to bias position, angular and energy distributions. |
---|
647 | |
---|
648 | <p> |
---|
649 | G4GeneralParticleSource is a concrete implementation of G4VPrimaryGenerator. To use, instantiate G4GeneralParticleSource in the G4VUserPrimaryGeneratorAction class, as demonstrated below. |
---|
650 | |
---|
651 | <PRE> |
---|
652 | MyPrimaryGeneratorAction::MyPrimaryGeneratorAction() { |
---|
653 | generator = new G4GeneralParticleSource; |
---|
654 | } |
---|
655 | |
---|
656 | void |
---|
657 | MyPrimaryGeneratorAction::GeneratePrimaries(G4Event*anEvent){ |
---|
658 | generator->GeneratePrimaryVertex(anEvent); |
---|
659 | } |
---|
660 | |
---|
661 | </PRE> |
---|
662 | The biasing can be configured through interactive commands. Extensive documentation can be found in <a href="#primary-particle-biasing">[3]</a>. Examples are also distributed with the Geant4 distribution in <B>examples/extended/eventgenerator/exgps</B>. |
---|
663 | |
---|
664 | <h3><a name="3.7.2.1.2"></a>3.7.2.1.2 Radioactive Decay Biasing</h3> |
---|
665 | <p> |
---|
666 | The G4RadioactiveDecay class simulates the decay of radioactive nuclei and implements the following biasing options: |
---|
667 | <UL> |
---|
668 | <LI> Increase the sampling rate of radionuclides within observation times through a user defined probability distribution function |
---|
669 | <LI> Nuclear splitting, where the parent nuclide is split into a user defined number of nuclides |
---|
670 | <LI> Branching ratio biasing where branching ratios are sampled with equal probability |
---|
671 | </UL> |
---|
672 | G4RadioactiveDecay is a process which must be registered with a process manager, as demonstrated below. |
---|
673 | |
---|
674 | <PRE> |
---|
675 | void MyPhysicsList::ConstructProcess() |
---|
676 | { |
---|
677 | ... |
---|
678 | G4RadioactiveDecay* theRadioactiveDecay = |
---|
679 | new G4RadioactiveDecay(); |
---|
680 | |
---|
681 | G4ProcessManager* pmanager = ... |
---|
682 | pmanager ->AddProcess(theRadioactiveDecay); |
---|
683 | ... |
---|
684 | } |
---|
685 | </PRE> |
---|
686 | <p> |
---|
687 | The biasing can be controlled either in compiled code or through interactive commands. Extensive documentation can be found in <a href="#exrdm">[4]</a> and <a href="#rdm">[5]</a>. Radioactive decay biasing examples are also distributed with the Geant4 distribution in <B>examples/extended/radioactivedecay/exrdm</B>. |
---|
688 | |
---|
689 | <h3><a name="3.7.2.1.3"></a>3.7.2.1.3 Hadronic Leading Particle Biasing</h3> |
---|
690 | <p>Two built in hadronic leading particle techniques are available. One is through G4Mars5GeV. This is an inclusive event generator for hadron interactions with nuclei. It is translated from the Mars23(98) version of the MARS code system <a href="#mars">[6]</a>. A fixed number of particles are generated at each vertex, with weights assigned as appropriate. G4Mars5GeV is valid with energies less than 5 GeV with these particle types: pi-, pi+, K+, K-, K0L, K0S, proton, anti-proton and gamma. To use this generator, create and register a G4Mars5GeV object with an appropriate inelastic process, as demonstrated below. |
---|
691 | |
---|
692 | <a name="G4Mars5GeV"> |
---|
693 | <PRE> |
---|
694 | MyPhysicsList::ConstructProcess() |
---|
695 | { |
---|
696 | ... |
---|
697 | G4Mars5Gev* leadModel = new G4Mars5GeV(); |
---|
698 | |
---|
699 | G4ProtonInelasticProcess* inelProcess = |
---|
700 | new G4ProtonInelasticProcess(); |
---|
701 | inelProcess->RegisterMe(leadModel); |
---|
702 | |
---|
703 | processManager->AdddiscreteProcess(inelProcess); |
---|
704 | } |
---|
705 | |
---|
706 | </PRE> |
---|
707 | |
---|
708 | <p> |
---|
709 | More examples can be found in the LHEP_LEAD, LHEP_LEAD_HP, QGSC_LEAD and QGSC_LEAD_HP physics |
---|
710 | lists. Additional documentation can be found in <a href="#g4mars5gev">[7]</a>. |
---|
711 | |
---|
712 | <P> |
---|
713 | Another hadronic leading particle biasing technique is implemented in the G4HadLeadBias utility. This method keeps only the most important part of the event, as well as representative tracks of each given particle type. So the track with the highest energy as well as one of each of Baryon, pi0, mesons and leptons. As usual, appropriate weights are assigned to the particles. Setting the <B>SwitchLeadBiasOn</B> environmental variable will activate this utility. |
---|
714 | |
---|
715 | <h3><a name="3.7.2.1.4"></a>3.7.2.1.4 Hadronic Cross Section Biasing</h3> |
---|
716 | Cross section biasing artificially enhances/reduces the cross section of a process. This may be useful for studying thin layer interactions or thick layer shielding. The built in hadronic cross section biasing applies to photon inelastic, electron nuclear and positron nuclear processes. |
---|
717 | </P> |
---|
718 | <P> |
---|
719 | The biasing is controlled through the <B>BiasCrossSectionByFactor</B> method in G4HadronicProcess, as demonstrated below. |
---|
720 | |
---|
721 | <a name="G4Mars5GeV"> |
---|
722 | <PRE> |
---|
723 | void MyPhysicsList::ConstructProcess() |
---|
724 | { |
---|
725 | ... |
---|
726 | G4ElectroNuclearReaction * theElectroReaction = |
---|
727 | new G4ElectroNuclearReaction; |
---|
728 | |
---|
729 | G4ElectronNuclearProcess theElectronNuclearProcess; |
---|
730 | theElectronNuclearProcess.RegisterMe(theElectroReaction); |
---|
731 | theElectronNuclearProcess.BiasCrossSectionByFactor(100); |
---|
732 | |
---|
733 | pManager->AddDiscreteProcess(&theElectronNuclearProcess); |
---|
734 | ... |
---|
735 | } |
---|
736 | </PRE> |
---|
737 | |
---|
738 | <p>More details can be found in <a href="#hadronicscrosssection">[8]</a>. |
---|
739 | |
---|
740 | <h3><a name="3.7.2.1"></a>3.7.2.2 G4WrapperProcess</h3> |
---|
741 | <p> |
---|
742 | G4WrapperProcess can be used to implement user defined event biasing. G4WrapperProcess, which is a process itself, wraps an existing process. By default, all function calls are forwared to the wrapped process. It is a non-invasive way to modify the behaviour of an existing process. |
---|
743 | |
---|
744 | <p> To use this utility, first create a derived class inheriting from G4WrapperProcess. Override the methods whose behaviour you would like to modify, for example, PostStepDoIt, and register the derived class in place of the process to be wrapped. Finally, register the wrapped process with G4WrapperProcess. The code snippets below demonstrate its use. |
---|
745 | |
---|
746 | |
---|
747 | <PRE> |
---|
748 | class MyWrapperProcess : public G4WrapperProcess { |
---|
749 | ... |
---|
750 | G4VParticleChange* PostStepDoIt(const G4Track& track, |
---|
751 | const G4Step& step) { |
---|
752 | // Do something interesting |
---|
753 | } |
---|
754 | }; |
---|
755 | |
---|
756 | </PRE> |
---|
757 | |
---|
758 | <PRE> |
---|
759 | void MyPhysicsList::ConstructProcess() |
---|
760 | { |
---|
761 | ... |
---|
762 | G4LowEnergyBremsstrahlung* bremProcess = |
---|
763 | new G4LowEnergyBremsstrahlung(); |
---|
764 | |
---|
765 | MyWrapperProcess* wrapper = new MyWrapperProcess(); |
---|
766 | wrapper->RegisterProcess(bremProcess); |
---|
767 | |
---|
768 | processManager->AddProcess(wrapper, -1, -1, 3); |
---|
769 | } |
---|
770 | |
---|
771 | </PRE> |
---|
772 | <P> |
---|
773 | <table> |
---|
774 | <tr><td valign=top><a name="imp-biasing1">[1]</a> |
---|
775 | <td><a href="http://cern.ch/geant4/collaboration/working_groups/geometry/biasing/Sampling.html">Scoring, |
---|
776 | geometrical importance sampling and weight roulette in Geant4</a> |
---|
777 | <tr><td valign=top><a name="imp-biasing2">[2]</a> |
---|
778 | <td><a href="http://cern.ch/geant4/collaboration/working_groups/geometry/biasing/BiasScoreUseCases.html">Use cases for |
---|
779 | importance biasing and scoring technique</a> |
---|
780 | <tr><td valign=top><a name="primary-particle-biasing">[3]</a> |
---|
781 | <td><a href="http://reat.space.qinetiq.com/gps/">Primary particle biasing</a> |
---|
782 | <tr><td valign=top><a name="exrdm">[4]</a> |
---|
783 | <td><a href="http://reat.space.qinetiq.com/septimess/exrdm/">Radioactive decay biasing example</a> |
---|
784 | <tr><td valign=top><a name="rdm">[5]</a> |
---|
785 | <td><a href="http://www.space.qinetiq.com/geant4/rdm.html">Radioactive decay biasing</a> |
---|
786 | |
---|
787 | <tr><td valign=top><a name="mars">[6]</a> |
---|
788 | <td><a href="http://www-ap.fnal.gov/MARS">MARS code system</a> |
---|
789 | |
---|
790 | <tr><td valign=top><a name="g4mars5gev">[7]</a> |
---|
791 | <td><a href="http://geant4.web.cern.ch/geant4/support/proc_mod_catalog/models/hadronic/LeadParticleBias.html">G4Mars5GeV</a> |
---|
792 | |
---|
793 | <tr><td valign=top><a name="hadronicscrosssection">[8]</a> |
---|
794 | <td><a href="http://www.triumf.ca/geant4-03/talks/03-Wednesday-AM-1/03-J.Wellisch/biasing.hadronics.pdf">Hadronic cross section biasing</a> |
---|
795 | |
---|
796 | </table> |
---|
797 | <P> |
---|
798 | |
---|
799 | |
---|
800 | </BODY> |
---|
801 | </HTML> |
---|