Home » Groups » STAR Protected

How to Insert User Step Routine to GSTAR

Updated on Sun, 2010-10-17 21:52. Originally created by leun on 2010-10-16 13:02.

How to Insert User Step Routine to GSTAR

 

Geant manual states that a user can write a custom hit-calculation routine (GUSTEP) for a sensitive volume, which will be called by GTRACK at the end of every step, and when entering a new volume. The following is my understanding of how to do this in the context of GSTAR and starsim.

For the purpose of record keeping, I will describe the process here. When it comes to STAR software, I am very much structurally challenged, so please excuse my ignorance on that front.

 

The first part comes from an ancient document posted on Drupal plus my own investigation.

http://drupal.star.bnl.gov/STAR/pwg/other-groups/photon-reconstruction/software-and-simulation-documentation/advanced-geometry-interface-gstar

 

1. Write a subroutine that calculates hits with a specific naming convention. If the sensitive volume that you want to associate this routine with has the name YYYY, (always 4 character) then the name of the subroutine should be YYYYstep(pointer,hit). From the document, "Their first integer input argument pointer is the address of the hit descriptor array (10 words, real) in the GCBANK memory", and the second argument is the hit value to be returned.

2. Stick this subroutine with other geometry files. (geometry/fpdmgeo in my case). I'm not sure exactly where this file has to be, but emprically this works fine. There is one such file currently in my directory, /star/u/leun/sim/fpdsim6/export/pams/geometry/fpdmgeo/ffpdstep.g 

3. Modify the geometry file in use in two ways. (/star/u/leun/sim/fpdsim6/export/pams/geometry/fpdmgeo/fpdmgeo2.g) Let's called the sensitive volume of interest FPCT. First, declare the subroutine you wrote as external. (example: EXTERNAL  FPCTSTEP) Second, find the line that declares the method of hit calculation for the volume of interest (example: "HITS FPCT ELOS:0:(0,50)") and replace it with "HITS YYYY USER:<Nbits>:(<min>,<max>)" (example: "HITS FPCT USER:0:(0,100000)") Here, ELOS is one of the many standard methods of hit calculation. 0 for Nbit should work fine for the most part, and min and max are the minimum and maximum hit values. 

 

Compile by running cons, and now you have your own hit calculation that will be called by Geant. 

As an example, here is a subroutine used to count Cerenkov photons in the FPD photocathode. (Written by someone else, comments added by me)

  

      subroutine FPCTSTEP(JJ,HIT)
      * define hit element USER as # of cherenkov photon detected 
      +CDE,TYPING,GCBANK,GCONST,GCUNIT,GCTMED,GCTRAK,GCKINE,GCSETS,AGCSTEP.
 
      Integer JJ
      Real HIT
 
      * Cerenkov detection efficiency controls fraction of the stopped particles that dumps energy. It's either all or nothing.
      * deStep>0.0 allows counting of only those stopped particles that dumped energy. Otherwise, efficiency is no use.
 
      if(Ipart == 50 & Istop ==2 & deStep>0.0) then
        hit=1.0
        return
      endif
      hit = 0.0
      return
      end 
 
 
I was able to follow this recipe and add my own subroutine to whatever sensitive volume I needed to mess with. However, when I started looking into the tracking process of the optical photons in the Pb glass, I realized that my subroutine was called only when the photon was absorbed. (non-zero energy loss) This was in sharp contrast to the tracking of the charged particles in the same volume, where I could see that my subroutine was being called at every step regardless of what happened. It seemed that a piece of code somewhere was blocking my routine for photons when Eloss is zero, but I couldn't find where.
 
Eventually, I contacted Jason Webb regarding this problem, and he was extremely helpful by looking into the starsim codes and pointing me to the part that controls the calling of the user defined routines. It is located at  

$STAR/asps/Simulation/starsim/atgeant/aggstep.age 
 
In it, there is indeed a line that says if neutral has no energy loss, don't bother calculating the hit. So, I simply had to comment it out, and that fixed the problem. The same directory contains other interesting codes that has to do with hit calculations and such. 
 
 
LKE
 
 

 

 

Groups: