source: trunk/documents/UserDoc/UsersGuides/ForApplicationDeveloper/html/Fundamentals/biasing.html @ 1358

Last change on this file since 1358 was 1208, checked in by garnier, 15 years ago

CVS update
File size: 33.9 KB
Line 
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
41Geant4 provides event biasing techniques which may be used to save computing time
42in such applications as the simulation of radiation shielding.  These are
43<I>geometrical splitting</I> and <I>Russian roulette</I> (also called geometrical
44importance sampling), and <I>weight roulette</I>.  A basic scoring system is also
45provided in order to monitor the sampling.  In this chapter, it is assumed that 
46the reader is familiar with both the usage of Geant4 and the concepts of importance
47sampling.  More detailed documentation may be found at
48<a href="#imp-biasing1">[1]</a>.
49A detailed description of different use-cases which employ the sampling and scoring
50techniques can be found at <a href="#imp-biasing2">[2]</a>.
51<BR>
52The 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>
58Standard scoring must be provided for quantities such as tracks entering a
59cell, average weight of entering tracks, energy of entering tracks, and
60collisions inside the cell.
61<BR>
62The purpose of importance sampling is to save computing time by sampling less
63often the particle histories entering "less important" geometry regions, and
64more often in more "important" regions. Given the same amount of computing
65time, an importance-sampled and an analogue-sampled simulation must show
66equal mean values, while the importance-sampled simulation will have a
67decreased variance.
68<BR>
69The implementation of scoring is independent of the implementation of
70importance sampling. However both share common concepts. <I>Scoring and
71importance sampling apply to particle types chosen by the user</I>.
72<BR>
73Examples on how to use scoring and importance sampling may be found
74in <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
78The kind of scoring referred to in this note and the importance sampling
79apply to spatial cells provided by the user.
80<P>
81A <B>cell</B> is a physical volume (further specified by it's replica number,
82if the volume is a replica). Cells may be defined in two kinds of
83geometries:
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>
90The user has the choice to score and/or sample by importance the particles
91of the chosen type, according to mass geometry or to parallel geometry.
92It is possible to utilize several parallel geometries in addition to the mass
93geometry. This provides the user with a lot of flexibility to define separate
94geometries for different particle types in order to apply scoring or/and
95importance 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>
113Samplers are higher level tools which perform the necessary changes
114of the Geant4 sampling in order to apply importance sampling and weight
115roulette. The samplers may also be used to apply scoring.</p>
116
117<p>Scoring and variance reduction may be combined arbitrarily
118for chosen particle types and may be applied to the mass or to parallel
119geometries. The samplers support all combinations.</p>
120
121<p>Different samplers are implemented for the mass and parallel
122geometries. Both implement the interface <tt>G4VSampler</tt>.
123<tt>G4VSampler</tt> allows the chosen combinations of scoring and variance
124reduction to be prepared and configures the sampling appropriately.  To this
125end <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>
153The 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
163arguments:</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>
190Each object of a sampler class is responsible for one particle type.
191The particle type is given to the constructor of the sampler classes via the
192particle type name, e.g. "neutron". Depending on the specific purpose,
193the <tt>Configure()</tt> of a sampler will set up specialized processes
194(derived from <tt>G4VProcess</tt>) for transportation in the parallel geometry,
195scoring, importance sampling and weight roulette for the given particle type.
196When <tt>Configure()</tt> is invoked the sampler places the processes in the
197correct 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
210a 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>
218The constructors of both classes take the particle type name (e.g. "neutron")
219as an argument. In addition, the constructor of
220<tt>G4ParallelGeometrySampler</tt> needs a reference to the physical world
221volume of the parallel geometry.</p>
222
223<h3><a name="scoring"></a>3.7.1.3 Scoring</h3>
224
225<p>
226Scoring is provided by a framework and is done according to particle type.
227Nevertheless it is also possible to score particles of different
228types into the same score. The framework may also be easily used for
229customized scoring.</p>
230<p> Scoring may be applied to a mass or a parallel geometry. It is done with
231an object genericly called a scorer using a sampler described above. The scorer
232receives the information about every step taken by a particle of chosen type.
233This information consists of an object of the Geant4 kernel class
234<tt>G4Step</tt> and an object of the class <tt>G4GeometryCellStep</tt> provided
235specifically 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>
241A "scorer" class derives from the interface <tt>G4VScorer</tt>. Users may
242create customized "scorers" or use the standard scoring.
243</p>
244<p>
245Classes 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 &amp;step, const G4GeometryCellStep &amp;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
321assigned to the cells. In general a particle history is either split
322or Russian roulette is played if the importance increases or decreases,
323respectively. A weight assigned to the history is changed according to
324the action taken.</p>
325
326<p>The tools provided for importance sampling require the user
327to have a good understanding of the physics in the problem. This is
328because the user has to decide which particle types require importance
329sampled, define the cells, and assign importance values to the cells. If
330this is not done properly the results cannot be expected to describe a
331real experiment.</p>
332
333<p>The assignment of importance values to a cell is done using an
334importance store described below.</p>
335
336<p>An "importance store" with the interface <tt>G4VIStore</tt> is used
337to store importance values related to cells. In order to do importance
338sampling the user has to create an object (e.g. of class <tt>G4IStore</tt>)
339of type <tt>G4VIStore</tt>. The samplers may be given a <tt>G4VIStore</tt>.
340The 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
343the world volume of the geometry used for importance sampling. This may
344be the world volume of the mass or of a parallel geometry. Importance
345stores 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>
360A concrete implementation of an importance store is provided by the class
361<tt>G4VStore</tt>. The <span class="emphasis"><em>public</em></span>
362part 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>
388The member function <tt>AddImportanceGeometryCell()</tt> enters a cell and
389an importance value into the importance store. The importance values may be
390returned either according to a physical volume and a replica number or
391according to a <tt>G4GeometryCell</tt>. The user mus be aware of the
392interpretation 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 &gt; 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>
418Importance sampling supports using a customized importance sampling algorithm.
419To this end, the sampler interface
420<a href="#G4VSampler"><tt>G4VSampler</tt></a> may be given a pointer to
421the 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>
465The user may have a customized algorithm used by providing a class
466inheriting from <tt>G4VImportanceAlgorithm</tt>. <br>
467If no customized algorithm is given to the sampler the default importance sampling
468algorithm 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
473importance 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 
478importance values as in importance sampling</p></li></ul>
479<p> 
480In contrast to importance sampling this technique is not weight blind. Instead
481the technique is applied according to the particle weight with respect to the
482current energy-space cell.</p>
483<p>
484Therefore the technique is convenient to apply in combination with other
485variance reduction techniques such as cross-section biasing and implicit capture.</p>
486
487<p>
488A 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&#xA;"></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>
508The energy-space cells are realized by <tt>G4GeometryCell</tt> as in importance
509sampling. 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 &amp;gCell,
518                                  G4double partEnergy) const = 0;
519    virtual G4bool IsKnown(const G4GeometryCell &amp;gCell) const = 0;
520    virtual const G4VPhysicalVolume &amp;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 &amp;worldvolume);
527     virtual ~G4WeightWindowStore();
528     virtual G4double GetLowerWeitgh(const G4GeometryCell &amp;gCell,
529                                     G4double partEnergy) const;
530     virtual G4bool IsKnown(const G4GeometryCell &amp;gCell) const;
531     virtual const G4VPhysicalVolume &amp;GetWorldVolume() const;
532     void AddLowerWeights(const G4GeometryCell &amp;gCell,
533                          const std::vector&lt;G4double&gt; &amp;lowerWeights);
534     void AddUpperEboundLowerWeightPairs(const G4GeometryCell &amp;gCell,
535                                         const G4UpperEnergyToLowerWeightMap&amp;
536                                         enWeMap);
537     void SetGeneralUpperEnergyBounds(const
538         std::set&lt;G4double, std::less&lt;G4double&gt; &gt; &amp; enBounds);
539 
540   private::
541   ... 
542  };
543</pre>
544<p>
545The user may choose equal energy bounds for all cells. In this case a set of
546upper energy bounds must be given to the store using the method
547<tt>SetGeneralUpperEnergyBounds</tt>.  If a general set of energy bounds have
548been set <tt>AddLowerWeights</tt> can be used to add the cells.</p>
549<p>
550Alternatively, the user may chose different energy regions for different cells.
551In this case the user must provide a mapping of upper energy bounds to lower
552weight bounds for every cell using the method
553<tt>AddUpperEboundLowerWeightPairs</tt>.</p>
554<p>
555Weight 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>
567A 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>
582The constructor takes three parameters which are used to:
583calculate the upper weight bound (upperLimitFaktor), calculate the survival weight
584(survivalFaktor), and introduce a maximal number (maxNumberOfSplits) of copies to
585be created in one go.</p>
586<p>
587In addition, the inverse of the maxNumberOfSplits is used to specify the minimum
588survival 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
593sampling and implicit capture are used together.  Implicit capture is not
594described here but it is useful to note that this procedure reduces a particle
595weight in every collision instead of killing the particle with some probability.</p>
596
597<p>
598Together with importance sampling the weight of a particle may become so low
599that it does not change any result significantly. Hence tracking a very low
600weight particle is a waste of computing time. Weight roulette is applied in
601order 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
605and the importance "Is" of the cell in which the source is located, by using
606the ratio "R=Is/Ic".</p>
607
608<p>
609Weight roulette uses a relative minimal weight limit and a relative survival
610weight. When a particle falls below the weight limit Russian roulette is applied.
611If the particle survives, tracking will be continued and the particle weight will
612be set to the survival weight.</p>
613
614<p>
615The 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>
639Geant4 supports physics based biasing through a number of general use, built in biasing
640techniques. 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>
646Primary 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>
649G4GeneralParticleSource 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>
662The 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>
666The 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>
672G4RadioactiveDecay 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>
709More examples can be found in the LHEP_LEAD, LHEP_LEAD_HP, QGSC_LEAD and QGSC_LEAD_HP physics
710lists. Additional documentation can be found in <a href="#g4mars5gev">[7]</a>.
711
712<P>
713Another 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>
716Cross 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>
719The 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>
742G4WrapperProcess 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>
Note: See TracBrowser for help on using the repository browser.