Message: Re: Storing neutron events using SD Not Logged In (login)
 Next-in-Thread Next-in-Thread
 Next-in-Forum Next-in-Forum

Warning Re: Storing neutron events using SD 

Forum: Event and Track Management
Re: Question Storing neutron events using SD (C Mos)
Date: 29 Oct, 2012
From: Michael H. Kelsey <Michael H. Kelsey>

On Mon, 29 Oct 2012 17:52:59 GMT, C Mos wrote:
> Hi I am using a sensitive detector to detect particles in a specific
> volume. This is what I have:
> 
> //=====================================
> 
>   G4bool APIdrdcSpatialSD::ProcessHits(G4Step* aStep,G4TouchableHistory*){
>   G4double edep = aStep->GetTotalEnergyDeposit();
>   G4int partPDGencoding=aTrack->GetDefinition()->GetPDGEncoding();
> 
>   //Next line for debugging purposes
>   if(partPDGencoding==G4Neutron::Definition()->GetPDGEncoding())  G4cout << "Neutron!..." << G4endl; 
> 
>   //Only collect data if there was energy deposit.
>   if(edep==0.) return false;
> 
>   //Only collect data for gamma and neutrons
>   if(partPDGencoding!=G4ParticleTable::GetParticleTable()->FindParticle("gamma")->GetPDGEncoding() &&
>   partPDGencoding!=G4Neutron::Definition()->GetPDGEncoding()) return false;  
>   //Data collection next...
>   }
> 
> //=====================================
> 
> I notice that I can identify neutron hits but they always seem to have
> 0. energy deposit. Is there a proper way to detect energy deposit by
> these wonderful particles?

In my role of writing the simulation code for CDMS backgrounds studies, I had exactly the same issue. The bottom line is that for the vast majority of cases, neutrons don't deposit energy directly. It is only the secondary particles from the interaction of a neutron with a nucleus which deposit energy.

1) They can scatter elastically off nuclei (which is most of what happens). If one of those scatters causes a visible nuclear recoil, the energy deposit comes from that nucleus stopping, not from the neutron.

2) They can "scatter" inelastically off nuclei. In that case, the original neutron (track) disappears, and one or more secondaries, often including a neutron :-), take its place. Those secondaries, or the recoiling nuclear fragment, may go on to deposit energy.

3) The neutron may be captured (especially if you have NeutronHP in your physics list). This is a special case of an "inelastic scatter," but where the secondaries are more constrained than in a typical hadronic interaction. You'll get one or more gammas and a nuclear fragment from de-excitation. That fragment could go on to decay radioactively if you have G4RadioactiveDecay in your physics list.

4) With some target materials, you can get neutron-induced fission.

In your code, you probably want to ignore the simple elastic scatters (by checking what process was responsible for the step), but you will have to get rid of the "edep != 0." requirement if you want to see any neutron hits at all.

By the way, you don't need to get the PDGencoding to identify neutrons. All of the particle definition pointers are singletons, so you can write:

  if (aTrack->GetDefinition() == G4Neutron::Definition())
    G4cout << "Neutron! ..." << G4endl;

and similarly for all the other particles.

    -- Michael Kelsey

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

1 More: Re: Storing neutron events using SD   (C Mos - 01 Nov, 2012)
(_ None: Re: Storing neutron events using SD   (Mike Kelsey - 01 Nov, 2012)
 Add Message Add Message
to: "Re: Storing neutron events using SD"

 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 ]