Message: Re: Broken Not Logged In (login)
 Next-in-Thread Next-in-Thread
 Next-in-Forum Next-in-Forum

None Re: Broken 

Forum: Event and Track Management
Re: Question Broken (Jaime Gomez)
Re: None Re: Broken (Marc Verderi)
Date: 19 Aug, 2014
From: Marc Verderi <Marc Verderi>

Hi Jaime,

I am not sure I understand your use case, but the hint I was suggesting 
you should work for your case: if the PreStepPoint GetStepStatus() is 
fGeomBoundary, then this track is incident to the quartz. The reason is 
that the (current) PreStepPoint was the PostStepPoint of the previous 
step : the fGeomBoundary flag then indicates that this (old) 
PostStepPoint corresponded to a step ending on the volume surface (the 
quartz). The track was hence incident to the quartz.

Hope this helps.

Cheers,
Marc

On 08/19/2014 04:27 PM, Jaime Arturo Gomez wrote:
> Hi Marc,
>
> Thanks for your reply. The reason I was using the PostStep volume was 
> because I wanted to gurntee that the particles are incident on the 
> quartz. So in my code I wanted to do something like
>
> if(PreMaterial->GetName()==Air && PostMaterial==Quartz)
>
> The quartz is the sensitive material that I have defined as such in my 
> CalorimeterSD.cc. But my worry is that if an electron (or other 
> particle) hits the face of the quartz it might decay to other things. 
> So I wanted to make sure I get only the hits that are on the absolute 
> face of the detector.
>
>
> On Tue, Aug 19, 2014 at 5:16 AM, Marc Verderi <verderi@in2p3.fr 
> <mailto:verderi@in2p3.fr>> wrote:
>
>
>     *** Discussion title: Event and Track Management
>
>     Hi Jaime,
>
>     Please find comments below.
>
>     Cheers,
>     Marc
>
>     On 08/18/2014 08:33 PM, Jaime Gomez wrote:
>     > *** Discussion title: Event and Track Management
>     >
>     > Hi all,
>     >
>     >   I am trying to make plots of particles that are incident on
>     the face of my detector and later will compare it to the read out
>     energy.
>     >
>     >   My detector is located at 60mm (z direction) in the code.
>     >
>     >
>     -----------------------------------CODE----------------------------------------
>     >
>     > const G4VTouchable* touchable =
>     > aStep->GetPostStepPoint()->GetTouchable(); G4VPhysicalVolume*
>     volume =
>     > touchable->GetVolume();
>     >
>     >    G4String name = volume->GetName();
>     Above line might certainly be responsible for the crash. A "usual"
>     trap
>     is that when a track reaches the boundary of the world volume, the
>     PostStepPoint has no (further) volume: the post step point has the
>     volume information of the next volume, but there is not next volume
>     after the world volume, hence the "volume" pointer is zero. Then
>     you are
>     calling a method on a null pointer here (I wonder if even
>     "touchable" is
>     null also).
>
>     Also, I understand that you use the PostStepPoint volume to see in
>     what
>     volume the step occurred: this is wrong, you must use the PreStepPoint
>     volume for this. As mentioned above, the PostStepPoint volume is the
>     "next" volume. Assume the step occurs in some volume A and finishes on
>     the boundary of A, and that after that there is some volume B. In this
>     case the volume in the PreStepPoint is A, while the one of the
>     PostStepPoint is B. Also, note that the energy deposit in such a
>     step is
>     in volume A (and not in volume B).
>     >    G4String det1 = "RPD1";
>     >    G4String det2 = "RPD2";
>     >    G4String det3 = "RPD3";
>     >    G4String det4 = "RPD4";
>     >    G4String det5 = "RPD5";
>     >    G4String det6 = "RPD6";
>     >    G4String det7 = "RPD7";
>     >    G4String det8 = "RPD8";
>     >    G4String det9 = "RPD9";
>     >    G4String det10 = "RPD10";
>     >    G4String det11 = "RPD11";
>     >    G4String det12 = "RPD12";
>     >    G4String det13 = "RPD13";
>     >    G4String det14 = "RPD14";
>     >    G4String det15 = "RPD15";
>     >    G4String det16 = "RPD16";
>     >
>     >    G4Material* PreMaterial =
>     aStep->GetPreStepPoint()->GetMaterial();
>     >    G4Material* PostMaterial =
>     aStep->GetPostStepPoint()->GetMaterial();
>     >
>     > G4ThreeVector PreHitPos = aStep->GetPreStepPoint()->GetPosition();
>     > G4ThreeVector ThisHitPos = aStep->GetPostStepPoint()->GetPosition();
>     >
>     >    if(ThisHitPos.z()<60.5 && PreHitPos.z()>58)
>     >       {
>     >        if(PreMaterial!=PostMaterial)
>     G4cout<<PreMaterial->GetName()<<"
>     "<<PostMaterial->GetName()<<G4endl; //check to see it goes from
>     Air to quartz
>     >
>     >            if( name==det1 || name==det2 || name==det3 ||
>     name==det4 ||
>     >                name==det5 || name==det6 || name==det7 ||
>     name==det8 ||
>     >                name==det9 || name==det10 || name==det11 ||
>     >                name==det12|| name==det13 || name==det14 ||
>     name==det15 ||
>     >                name==det16)
>     >               {
>     >                 ....I Fill histos here
>     >               }//end of detector hit filter
>     >        }//end of z-hit filter
>
>     Have you considered using a G4VSensitiveDetector ? It would make you
>     code much easier (and faster than the above one). Given you are
>     using a
>     stepping action here, you have to have all the logic to check if
>     you are
>     in one of the right volumes. Using a G4VSensitiveDetector, you will
>     attach it only to the volume(s) you want to collect information
>     in. You
>     will do this attachment(s) in your DetectorConstruction class. Your
>     G4VSensitiveDetector code will be guaranteed to be invoked only
>     when the
>     track is making a step in this(these) volume(s).
>     In this case, you may check the PreStepPoint "GetStepStatus()" status,
>     and if the track is just entering the volume, you will have this
>     status
>     to be fGeomBoundary. Your sensitive detector code may hence reduce to
>     something like:
>
>     G4bool MySensitiveDetector::ProcessHit(const G4Step* step,
>     G4TouchableHistory*)
>     {
>          if ( step->GetPreStepPoint()->GetStepStatus() == fGeomBoundary )
>          {
>                 ... fill histos here
>          }
>
>         return true;
>
>     }
>
>     ProcessHit, you guess, is the method called if the step occurs in the
>     volume to which the sensitive detector is attached.
>
>     >
>     >
>     _____________________________________________________________________
>     >
>     > It just seems that the 2nd, and 3rd if statements when used together
>     > cause a crash and I was wondering if anyone could help me figure out
>     > why.
>     >
>     > thanks for any advice you can give J
>     >
>     >     Attachment:
>     >
>     http://hypernews.slac.stanford.edu/HyperNews/geant4/get/AUX/2014/08/18/11.16-1005-SteppingAction.cc
>     >
>     > -------------------------------------------------------------
>     > Visit this GEANT4 at hypernews.slac.stanford.edu
>     <http://hypernews.slac.stanford.edu> message (to reply or
>     unsubscribe) at:
>     >
>     http://hypernews.slac.stanford.edu/HyperNews/geant4/get/eventtrackmanage/1221.html
>
>     -------------------------------------------------------------
>     Visit this GEANT4 at hypernews.slac.stanford.edu
>     <http://hypernews.slac.stanford.edu> message (to reply or
>     unsubscribe) at:
>     http://hypernews.slac.stanford.edu/HyperNews/geant4/get/eventtrackmanage/1221/1.html
>
>


 [ MIME part of type text/html without a name stripped ]

 Add Message Add Message
to: "Re: Broken"

 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 ]