Message: Re: Precision of geometry definition and tracking in replicas Not Logged In (login)
 Next-in-Thread Next-in-Thread
 Next-in-Forum Next-in-Forum

None Re: Precision of geometry definition and tracking in replicas 

Forum: Event and Track Management
Re: Question Precision of geometry definition and tracking in replicas (Matthias Böhm)
Re: None Re: Precision of geometry definition and tracking in replicas (Vladimir IVANTCHENKO )
Date: 04 Feb, 2005
From: Gumplinger Peter <Gumplinger Peter>

Please, see /examples/extended/optical/LXe/src/LXeScintSD::ProcessHits as a partial example of what you want to do:

G4double edep = aStep->GetTotalEnergyDeposit();

if (edep==0.) return false; // No edep so dont count as hit

and:

  G4StepPoint* thePrePoint = aStep->GetPreStepPoint();
  G4StepPoint* thePostPoint = aStep->GetPostStepPoint();

  //Get the average position of the hit
  G4ThreeVector pos = thePrePoint->GetPosition() + thePostPoint->GetPosition();
  pos/=2.;

For your information, I post here a summary of G4's kinematical/geometrical design by Katsuya Amako

===============================================================

You have to remember two very basic design decisions in Geant4 particle transportation to understand the design:

  1. Geant4 transports "particle" in stepwise
  ============================================
     From this decision:
       a) we need an object to describe 'step', which we call 'G4Step',
       b) because 'step' is a line segment (i.e. geometrically extend 
          object), it has two endpoints, which we call 'PreStepPoint' 
          and 'PostStepPoint'.   These two endpoints keep 
          kinematical/geometrical information of the "particle"  
          before and during/after one step.

  2. Geant4 transports 'track' instead of so-called "particle"
  ============================================================
     From this decision:
       a) we need an object to describe 'track', which we call 'G4Track',
       b) 'track' is transported through geometry taking into account
          physics interactions,
       c) 'track' keeps kinematical/geometrical information of the 
          "particle",
       d) 'track' knows the most recent step taken before coming to the 
          current position, i.e. 'track' has 'step' information. 

The design diagram is as follows:

   +---------+           +--------+         +-------------+
   | G4Track |o----------| G4Step |o--------| G4StepPoint |
   |         |----------o|        |       2 |             |
   +---------+           +--------+         +-------------+

The design described above is explained in the Section '5.1 Tracking -- 5.1.1 Basic Concepts' of 'Geant4 User's Guide'

Based on the above design Geant4 provides two possible ways for users to access to kinematical/geometrical information of the "particle":

  - 'track' oriented way:  To retrieve information basically through 'track'.
  - 'step'  oriented way:  To retrieve information basically through 'step'.

In addition to above two design decisions, you have to know another important design decision described below on 'track' and 'step'. This is related to the volume relocation when the "particle" is transported on the volume boundary. You need to understand this when you want to access geometry related information like "current" volume, touchable, material, etc, in UserAction/Hits.

  1. Geometrical information in 'step' at geometrical boundary
  ============================================================
     Because 'step' is a line segment it is natural to design so that the 
     head of the endpoint of 'step' (i.e. 'PostStepPoint') belongs to the 
     new volume when the "particle" is transported to a geometrical 
     boundary. This means the 'step' (i.e. line segment) itself belongs to 
     the old volume but only 'PostStepPoint' (i.e. point) belongs to the 
     new volume. 

  2. Geometrical information in 'track' at geometrical boundary
  ==============================================================
     Because 'track' is conceptually an object which keeps "particle"
     information at a particular space point, there are two
     choices to assign its volume when it is transported to the geometrical
     boundary. The natural choice here in the design was that it belongs 
     to old volume because the 'track' has a 'step' which belongs to the 
     old volume. You can memorize this by saying that the "particle" (i.e.
     'track') always moves in the same volume during and after one step. 

     Except for this volume related information (volume, touchable, material,
     ...) the information kept by 'track' is same as one kept by 'PostStepPoint'.   

'track' oriented approach =========================

 [Note] 
   - An example of this approach can be found in Novice/ExN03.
   - The statements below which starts with ">>" are questions.

>>>> 1) Being in SteppingAction (or ProcessHits), what is the *simplest* and
>>>> safe way to get the current physical volume ?

   In this question you have to be very careful what you mean by "current" 
   physical volume (also true for 'touchable', 'material' etc). 

   Case 1:
   -------
   If "current" volume means that you always want to access the volume 
   where its 'step' belongs (this is the most case you want to know), 
   you should use the method G4Track::GetVolume(). For 'Touchable' and 
   'Material' G4Track::GetTouchable(), G4Track::GetMaterial().

   Case 2:
   -------
   If "current volume" means that you want the volume automatically points 
   to the new volume when the relocation occurred at the geometry  boundary 
   (i.e. the new volume is one 'particle' gets into when the next step starts), 
   then you should always use G4Track::GetNextVolume(). For 'Touchable' 
   and 'Material' G4TracK::GetNextTouchable(), G4Track::GetNextMaterial().

>>>> -its pointer

    const G4Track* track = aStep->GetTrack();    // aStep comes from method argument
    G4VPhysicalVolume* volumeOld = track->GetVolume();                // Case 1
    G4VPhysicalVolume* volumeNew = track->GetNextVolume();            // Case 2

>>>> -its name

   G4String name = volumeOld->GetName();                              // Case 1
   G4String name = volumeNew->GetName();                              // Case 2

>>>> -its copy number

   G4int copynumber = track->GetTouchable()->GetReplicaNumber();      // Case 1
   G4int copynumber = track->GetNextTouchable()->GetReplicaNumber();  // Case 2

>>>> -its material

   G4Material* material = track->GetMaterial();                       // Case 1
   G4Material* material = track->GetNextMaterial();                   // Case 2

>>>> -its logical volume

   G4LogicalVolume* logVolOld = volumeOld->GetLogicalVolume();        // Case 1
   G4LogicalVolume* logVolNew = volumeNew->GetLogicalVolume();        // Case 2

>>>> -its region

   G4Region* logVolOld->GetRegion();                                  // Case 1
   G4Region* logVolNew->GetRegion();                                  // Case 2

>>>> 2) Same questions for its mother.

   G4VPhysicalVolume* mphysVolOld = track->GetTouchable()->GetVolume(1);     // Case 1
   G4VPhysicalVolume* mphysVolNew = track->GetNextTouchable()->GetVolume(1); // Case 2

   G4String mname = mphysVolOld->GetName();                           // Case 1
   G4String mname = mphysVolNew->GetName();                           // Case 2

>>>> 3) How to check that we reach a bondary ?

   Here you have to access 'step' to get this information.

   G4StepPoint* pPostStepP = aStep->GetPostStepPoint(); // aStep comes from method argument
   if (pPostStepP->GetStepStatus() = fGeomBoundary) {......}

>>>> 4) In that case, same questions as above for the next volume.

   See codes in 'Case 2' above.

'step' oriented approach =========================

 [Note] 
   - An example of this approach can be found in Novice/ExN04.
   - The statements below which starts with ">>" are original questions.

>>>> 1) Being in SteppingAction (or ProcessHits), what is the *simplest* and
>>>> safe way to get the current physical volume ?

   Again you have to be careful what you mean by "current".

   Case 1:
   -------
   If "current" volume means that you always want to access to the volume 
   where its 'step' belongs (this is the most case user want to know), 
   you should use PreStepPoint.

   Case 2:
   -------
   If "current volume" means that you want the volume automatically points 
   to the new volume when the relocation occurred at the geometry  boundary 
   (i.e. the new volume is one 'particle' gets into after reaching a boundary), 
   then you should use PostStepPoint.

>>>> -its pointer

   G4StepPoint* pPreStepP  = aStep->GetPreStepPoint();  // aStep comes from method argument
   G4StepPoint* pPostStepP = aStep->GetPostStepPoint(); // aStep comes from method argument

   G4VPhysicalVolume* volumeOld = pPreStepP->GetPhysicalVolume();     // Case 1
   G4VPhysicalVolume* volumeNew = pPostStepP->GetPhysicalVolume();    // Case 2

>>>> -its name

   G4String name = volumeOld->GetName();                              // Case 1
   G4String name = volumeNew->GetName();                              // Case 2

>>>> -its copy number

   G4int copynumber = pPreStepP->GetTouchable()->GetReplicaNumber();  // Case 1
   G4int copynumber = pPostStepP->GetTouchable()->GetReplicaNumber(); // Case 2

>>>> -its material

   G4Material* material = pPreStepP->GetMaterial();                   // Case 1
   G4Material* material = pPostStepP->GetMaterial();                  // Case 2

>>>> -its logical volume

   G4LogicalVolume* logVolOld = volumeOld->GetLogicalVolume();        // Case 1
   G4LogicalVolume* logVolNew = volumeNew->GetLogicalVolume();        // Case 2

>>>> -its region

   G4Region* logVolOld->GetRegion();                                  // Case 1
   G4Region* logVolNew->GetRegion();                                  // Case 2

>>>> 2) Same questions for its mother.

   G4VPhysicalVolume* mphysVolOld = pPreStepP->GetTouchable()->GetVolume(1);  // Case 1
   G4VPhysicalVolume* mphysVolNew = pPostStepP->GetTouchable()->GetVolume(1); // Case 2

   G4String mname = mphysVolOld->GetName();                           // Case 1
   G4String mname = mphysVolNew->GetName();                           // Case 2

>>>> 3) How to check that we reach a bondary ?

   G4StepPoint* pPostStepP = aStep->GetPostStepPoint(); // aStep comes from method argument
   if (pPostStepPoint->GetStepStatus() = fGeomBoundary) {......}

>>>> 4) In that case, same questions as above for the next volume.

   See codes in 'Case 2' above.

Inline Depth:
 1 1
 All All
Outline Depth:
 1 1
 2 2
 All All
Add message: (add)

 Add Message Add Message
to: "Re: Precision of geometry definition and tracking in replicas"

 Subscribe Subscribe

This site runs SLAC HyperNews version 1.11-slac-98, derived from the original HyperNews


[ Geant 4 Home | Geant 4 HyperNews | Search | Request New Forum | Feedback ]