The target cannot record the decay spectrum correctly!

I am currently using the example of rdecay02. My geometry is set to the following style.


The yellow one is the LaBr3 detector and the white one is CsI detector. Now I am simulating the LaBr3’ internal radioaction. The radioactive source is Ac227. The following is my code to describe the uniform distribution of radioactive sources in the LaBr3.

PrimaryGeneratorAction::PrimaryGeneratorAction()
{
 // GPS
  // fParticleGun = new G4GeneralParticleSource();

   // particleGun
#if 1
  G4int n_particle = 1;
  fParticleGun  = new G4ParticleGun(n_particle);
  
  fParticleGun->SetParticleEnergy(0*eV);

  // uniform distribution inside the cylinder
    G4double H = 2.5 * cm;
    G4double R = 1.9 * cm;
    G4double x, y, z;
    G4double sigma1, sigma2, sigma3;
    G4double eta1, eta2;

    do {
        sigma1 = G4UniformRand();
        sigma2 = G4UniformRand();
        sigma3 = G4UniformRand();
        eta1 = 2 * sigma1 - 1;
        eta2 = 2 * sigma2 - 1;
    } while ((eta1 * eta1 + eta2 * eta2) > 1);

    x = R * eta1;
    y = R * eta2;
    z = H * (2 * sigma3 - 1);

    fParticleGun->SetParticlePosition(G4ThreeVector(x, y, z));


      // set isotropic direction,
    G4double px; G4double py; G4double pz;

    G4double aperture = 180 * G4UniformRand();// 180.;
    G4double phi = 360 * G4UniformRand();

    G4double costheta; G4double sintheta; G4double cosphi; G4double sinphi;
    costheta = cos(aperture *deg);  // isotropic in aperture
    sintheta = sin(aperture *deg);
    cosphi = cos(phi *deg);
    sinphi = sin(phi *deg);

    px = sintheta * cosphi;
    py = sintheta * sinphi;
    pz = costheta;

    fParticleGun->SetParticleMomentumDirection(G4ThreeVector(px,py, pz));
  
  #endif
  //fParticleGun->SetParticleMomentumDirection(G4ThreeVector(0.,0.,1.));
}

//-------------------------------------------------------------

#
# Macro file for "rdecay02.cc"
#
/control/verbose 1
###/control/cout/ignoreThreadsExcept 0
###/run/numberOfThreads 2
/run/verbose 1
#
/process/em/verbose 0
/process/had/verbose 0
#
/run/initialize
#
# Set a very high time threshold to allow all decays to happen
/process/had/rdm/thresholdForVeryLongDecayTime 1.0e+60 year
#
/gun/particle ion
/gun/ion 89 227
#
/analysis/setFileName rdecay02
/analysis/h1/set 0 100 0. 10. MeV
/analysis/h1/set 1 100 0. 10. MeV
/analysis/h1/set 2 100 0. 10. MeV
/analysis/h1/set 3 100 0. 10. MeV
/analysis/h1/set 4 100 0. 10. MeV
/analysis/h1/set 5 100 0. 10. MeV
/analysis/h1/set 6 100 0. 10. MeV
#
/run/printProgress 1000
/run/beamOn 10000

Below is my macro file:

histogram 0: The Pulse Height Spectrum (PHS) of the target.
histogram 1: The PHS of the detector.
histogram 2: The combined PHS of the target and detector.
histogram 3: The anti-coincidece PHS of the target.
histogram 4: The anti-coincidece PHS of the detector.
histogram 5: The coincidece PHS between the target and detector.
histogram 6: The emitted particle energy spectrum.

However,the recorded data is not correct in H10.The decay emission spectrum is right.

Hi Jifei

This will not fix your problem, but…

will not give you an isotropic angular distribution. An isotropic angular distribution is not uniform in theta (what you call aperture). The solid angle reduces towards the poles (0 and 180 deg) (think of a globe of Earth). It is uniform in cos(theta). In any case, always use

  auto r = G4RandomDirection();
  px = r.x();
  ...

See G4RandomDirection.hh for some explanation, including how to restrict the aperture.

Thanks for your reply,this is indeed an easy way to express isotropic. My method is based on the coordinate system. It’s no problem for me to check it visually. :grinning: