3

I am using GATE(which uses Geant4) to do MC studies on dosimetric output. I am using a cylindrical cobalt source at 80 cm SAD to measure the PDD in a water phantom and dose at depth of 10 cm.

I now want to simulate a smaller source (say, r/2 and h/2) and compare the dosimetric output at a depth of 10 cm. Besides the geometry, I see that I am able to control the number of particles and time of the simulation. What would be the best way to change these two parameters to mimic the lower output from a smaller source? Or is there any other parameter that can be changed to mimic a smaller source? I am trying to calculate the output factor of the smaller source w.r.t. to the original source.

Severin Pappadeux
  • 18,636
  • 3
  • 38
  • 64
MadPhysicist
  • 175
  • 2
  • 11
  • What is the problem changing (r,h)? It is unclear how you're going to calibrate it. At the end of the day you should normalize your PDD per Bq, single decay/second (or fixed activity per some Bq and/or Ci of the source) – Severin Pappadeux Feb 17 '22 at 14:38
  • Thanks for responding. I am able to change r,h. I am using a fixed activity. However, I am getting an output factor of 1(w.r.t. ref dose at d=10 cm of the previous case). I was thinking it is because I am using the same number of particles(?). – MadPhysicist Feb 17 '22 at 15:39
  • You have to normalize to the same number of decays. Just in case I put in answer my code for Co60 cylindrical source in Geant4, 10 events at once – Severin Pappadeux Feb 17 '22 at 16:32

1 Answers1

1

Not sure if it helps, this is cylindrical source with Co60

Source::Source():
    _particleGun{nullptr},
    _sourceMessenger{nullptr},
    _radius{-1.0},
    _halfz{-1.0},
    _nof_particles{10}
{
    _particleGun = new G4ParticleGun( 1 );

    G4ParticleTable* particleTable = G4ParticleTable::GetParticleTable();
    G4String particleName = "gamma"; // "geantino"
    _particleGun->SetParticleDefinition(particleTable->FindParticle(particleName));

    _particleGun->SetParticlePosition(G4ThreeVector(0., 0., 0.));
    _particleGun->SetParticleMomentumDirection(G4ThreeVector(0., 0., 1.));
    _particleGun->SetParticleEnergy(1000.0*MeV);

    _sourceMessenger = new SourceMessenger(this);
}


Source::~Source()
{
    delete _particleGun;
    delete _sourceMessenger;
}


troika Source::sample_direction()
{
    double phi   = 2.0 * M_PI * G4UniformRand();

    double cos_z = 2.0 * G4UniformRand() - 1.0;
    double sin_z = sqrt( (1.0 - cos_z) * (1.0 + cos_z) );

    return troika{ sin_z * cos(phi), sin_z * sin(phi), cos_z };
}


double Source::sample_energy()
{
    return (G4UniformRand() < P_lo) ? E_lo : E_hi;
}


void Source::GeneratePrimaries(G4Event* anEvent)
{
    for(int k = 0; k != _nof_particles; ++k) // we generate _nof_particles at once
    {
        // here we sample spatial decay vertex uniformly in the cylinder
        double z   = _halfz * ( 2.0*G4UniformRand() - 1.0 );
        double phi = 2.0 * M_PI * G4UniformRand();
        double r   = _radius * sqrt(G4UniformRand());

        auto x = r * cos(phi);
        auto y = r * sin(phi);
        _particleGun->SetParticlePosition(G4ThreeVector(x, y, z));

        // now uniform-on-the-sphere direction
        auto dir = sample_direction();
        _particleGun->SetParticleMomentumDirection(G4ThreeVector(dir._wx, dir._wy, dir._wz));

        // energy 50/50 1.17 or 1.33
        auto e = sample_energy();
        _particleGun->SetParticleEnergy(e);

        // all together in a vertex
        _particleGun->GeneratePrimaryVertex(anEvent);
    }
}
Severin Pappadeux
  • 18,636
  • 3
  • 38
  • 64
  • Hey! Thanks again! Could you please explain how the normalization of the number of particles takes place? I am using GATE where the commands/statements are different and I need to translate from C++. I am using GPS inside the physical volume of the source and forcing the particle(2 Energies) to emit within a solid angle pointed towards the phantom. What I understood from your code is that you are defining a source, and emitting a particle, giving the particle defined energy and initial momentum. How particles were normalized is not clear. Could you please explain more? – MadPhysicist Feb 18 '22 at 03:53
  • Here is my **[code](https://drive.google.com/file/d/1yfjk5vI4-5CHHcSGSn6bIs8rg9-ycD1w/view?usp=sharing)** for reference. – MadPhysicist Feb 18 '22 at 03:55
  • Sure. This is part (source term) of the Co60 collimation system. I have two collimators, one large and one smaller. I sample 10bln decays (this is, btw, I need 10 photons per event, counter is 32bit integer) for both. I check yield (writing phase space file) after the collimator but before phantom. Naturally, for smaller collimator I have less particles in the file, smaller yield. Then I calculate dose at given point in the phantom. Dose is per one Bq (one decay), ca 5*10^-15. For larger collimator I set OF=1 (multiplier of 1/5 10^-15), for smaller collimator I get lower dose and OF=0.93 or so – Severin Pappadeux Feb 18 '22 at 14:28
  • Thanks again! So my question is that for the smaller collimator, you say you use a smaller number of particles. How do you decide on this number? Is that decided by a ratio of the collimator sizes? So in my case, should it be a ratio of volumes? I am not using a phase space file but directly setting my number of primary particles and I think setting the right number of particles should be the solution. – MadPhysicist Feb 18 '22 at 21:30
  • No, no, no, I sample the same number of decays (10bln) in the source. But after collimation there are smaller number of particles going toward a phantom (I checked it collecting PhSF after collimator before phantom), and smaller dose in the area of interest. Thus, naturally I got smaller OF for the small collimator, .93 vs1. – Severin Pappadeux Feb 19 '22 at 00:10
  • 'How do you decide on this number? Is that decided by a ratio of the collimator sizes?' It happens naturally as smaller number of particles are leaving the source in the collimator opening and be able to reach phantom. Look, I simulate all decays. I waste CPU time a lot because decay photons are mostly absorbed in the source and surrounding shielding. – Severin Pappadeux Feb 19 '22 at 00:12
  • And difference with you, I venture, is that you photons are somehow already sampled in the proper theta angle, not like me where I'm always doing proper 4pi directions. Frankly, I'm curious how GATE is doing inside-the-source proper absorption and scattering, if angle is limited to 12 degrees or so. – Severin Pappadeux Feb 19 '22 at 00:15
  • If your source is similar to mine, namely, made from pellets packed into ss shell, and pellets have the same Co60 content (atoms per volume), then looks like number decays would be proportional to the volume, and therefore you have to normalize by volumes ratio. Still curious how sampling in the fixed theta was dose. – Severin Pappadeux Feb 19 '22 at 17:18
  • "was dose" -> "was done" – Severin Pappadeux Feb 21 '22 at 13:00