How to detect Gamma photons through composite material By NaI Detector

Please fill out the following information to help in answering your question, and also see tips for posting code snippets. If you don’t provide this information it will take more time to help with your problem!

Geant4 Version: 11.01
Operating System: Ububntu
Compiler/Version: 9.4.0
CMake Version: 3.20

— As a new user in Geant4 please i need guidance in developing this. i wanted to shoot gamma passing my composite material to be detected by NaI detector, please guide me on which example is the approriate to work on in order to make my work more easier as a new user. because i want to see the relaibility of using the composite material in blocking the gamma, i want to calculate the attenuation coefficients, so i need to find the Io and I values by detecting how much photons are detected by the Detector with the material and without the material.

1 Like

Besides i modified exampleB1 in geant4 version 11.01, to my desired set-up. But i am confused on how to get the number of photons recorded with my material and with out. Attached are the modified codes, please guide me whether the codes are correctly modified in order to archive the above aim.
DetectorConstruction.cc (14.3 KB)
PrimaryGeneratorAction.cc (4.2 KB)
I also modified run2.mac file to this.

# Macro file for example B1
# 
# To be run preferably in batch, without graphics:
# % exampleB1 run2.mac
#
#/run/numberOfThreads 4
/run/initialize
#
/control/verbose 2
/run/verbose 2
#
# gamma 6 MeV to the direction (0.,0.,1.)
# 10000 events
#
/gun/particle gamma
/gun/energy 150 KeV
#
/run/printProgress 100
/run/beamOn 1000000
# 
# proton 210 MeV to the direction (0.,0.,1.)
# 1000 events
#
#/gun/particle proton
/#gun/energy 210 MeV
#
#/run/beamOn 1000

when i modified the run2.mac, i received this in the terminal, i could’nt figure out what does that mean.

G4WT1 > 
-------- WWWW ------- G4Exception-START -------- WWWW -------
*** G4Exception : GeomMgt1001
      issued by : G4LogicalVolumeStore::GetVolume()
Volume NOT found in store !
        Volume Envelope NOT found in store !
        Returning NULL pointer.
*** This is just a warning message. ***
-------- WWWW -------- G4Exception-END --------- WWWW -------

G4WT1 > 
-------- WWWW ------- G4Exception-START -------- WWWW -------
*** G4Exception : MyCode0002
      issued by : PrimaryGeneratorAction::GeneratePrimaries()
Envelope volume of box shape not found.
Perhaps you have changed geometry.
The gun will be place at the center.
*** This is just a warning message. ***
-------- WWWW -------- G4Exception-END --------- WWWW -------

G4WT0 > 
-------- WWWW ------- G4Exception-START -------- WWWW -------
*** G4Exception : GeomMgt1001
      issued by : G4LogicalVolumeStore::GetVolume()
Volume NOT found in store !
        Volume Envelope NOT found in store !
        Returning NULL pointer.
*** This is just a warning message. ***
-------- WWWW -------- G4Exception-END --------- WWWW -------

and these

G4WT2 > [thread 2] Thread-local run terminated.
G4WT2 > [thread 2] Run Summary
G4WT1 > [thread 1] Thread-local run terminated.
G4WT1 > [thread 1] Run Summary
G4WT1 > [thread 1]   Number of events processed : 248
G4WT3 > [thread 3] Thread-local run terminated.
G4WT0 > [thread 0] Thread-local run terminated.
G4WT3 > [thread 3] Run Summary
G4WT0 > [thread 0] Run Summary
G4WT0 > [thread 0]   Number of events processed : 217
G4WT3 > [thread 3]   Number of events processed : 256
G4WT1 > [thread 1]   User=11.150000s Real=3.075644s Sys=0.080000s [Cpu=365.1%]
G4WT1 > 
G4WT2 > [thread 2]   Number of events processed : 279
G4WT3 > [thread 3]   User=0.430000s Real=0.185412s Sys=0.000000s [Cpu=231.9%]
G4WT0 > [thread 0]   User=11.170000s Real=3.083001s Sys=0.080000s [Cpu=364.9%]
G4WT3 > 
G4WT2 > [thread 2]   User=11.170000s Real=3.082571s Sys=0.080000s [Cpu=365.0%]
G4WT0 > 
G4WT3 > --------------------End of Local Run------------------------
G4WT3 >  The run consists of 256 proton of 210 MeV
G4WT3 >  Cumulated dose per run, in scoring volume : 492.042 picoGy  rms = 6.54833 picoGy 
G4WT3 > ------------------------------------------------------------
G4WT3 > 
G4WT2 > 
G4WT0 > --------------------End of Local Run------------------------
G4WT1 > --------------------End of Local Run------------------------
G4WT0 >  The run consists of 217 proton of 210 MeV
G4WT2 > --------------------End of Local Run------------------------
G4WT1 >  The run consists of 248 proton of 210 MeV
G4WT0 >  Cumulated dose per run, in scoring volume : 414.417 picoGy  rms = 6.84741 picoGy 
G4WT2 >  The run consists of 279 proton of 210 MeV
G4WT0 > ------------------------------------------------------------
G4WT0 > 
G4WT2 >  Cumulated dose per run, in scoring volume : 538.28 picoGy  rms = 6.82819 picoGy 
G4WT1 >  Cumulated dose per run, in scoring volume : 475.73 picoGy  rms = 6.44467 picoGy 
G4WT2 > ------------------------------------------------------------
G4WT2 > 
G4WT1 > ------------------------------------------------------------
G4WT1 > 
 Run terminated.
Run Summary
  Number of events processed : 1000
  User=11.190000s Real=3.092077s Sys=0.080000s [Cpu=364.5%]

--------------------End of Global Run-----------------------
 The run consists of 1000 
 Cumulated dose per run, in scoring volume : 1.92047 nanoGy  rms = 13.3407 picoGy 
------------------------------------------------------------

please is there any body to Help me on this? iam still waiting hope to get your feedback soon thank you.

Today, i also able to modify example B4. to my desired set-up as follows:
DetectorConstruction.cc (14.8 KB)
PrimaryGeneratorAction.cc (4.0 KB)
DetectorConstruction.hh (2.4 KB)
in the Example B4a i only modify the above codes, below is the result of the setup.


and the result from the terminal are.

G4WT3 > 
-------- WWWW ------- G4Exception-START -------- WWWW -------
*** G4Exception : MyCode0002
      issued by : PrimaryGeneratorAction::GeneratePrimaries()
World volume of box shape not found.
Perhaps you have changed geometry.
The gun will be place in the center.
*** This is just a warning message. ***
-------- WWWW -------- G4Exception-END --------- WWWW -------

G4WT0 > [thread 0] Thread-local run terminated.
G4WT1 > [thread 1] Thread-local run terminated.
G4WT3 > [thread 3] Thread-local run terminated.
G4WT2 > [thread 2] Thread-local run terminated.
G4WT0 > [thread 0] Run Summary
G4WT0 > [thread 0]   Number of events processed : 10
G4WT0 > [thread 0]   User=0.040000s Real=0.040465s Sys=0.010000s [Cpu=123.6%]
G4WT2 > [thread 2] Run Summary
G4WT2 > [thread 2]   Number of events processed : 50
G4WT2 > [thread 2]   User=0.670000s Real=0.823868s Sys=0.250000s [Cpu=111.7%]
G4WT2 > 
G4WT2 >  ----> print histograms statistic for the local thread 
G4WT2 > 
G4WT2 >  EAbs : mean = 0 eV  rms = 0 eV 
G4WT1 > [thread 1] Run Summary
G4WT3 > [thread 3] Run Summary
G4WT3 > [thread 3]   Number of events processed : 20
G4WT3 > [thread 3]   User=0.140000s Real=0.146578s Sys=0.060000s [Cpu=136.4%]
G4WT3 > 
G4WT3 >  ----> print histograms statistic for the local thread 
G4WT3 > 
G4WT3 >  EAbs : mean = 0 eV  rms = 0 eV 
G4WT0 > 
G4WT2 >  EGap : mean = 0 eV  rms = 0 eV 
G4WT1 > [thread 1]   Number of events processed : 20
G4WT1 > [thread 1]   User=0.270000s Real=0.331259s Sys=0.130000s [Cpu=120.8%]
G4WT1 > 
G4WT1 >  ----> print histograms statistic for the local thread 
G4WT1 > 
G4WT1 >  EAbs : mean = 0 eV  rms = 0 eV 
G4WT1 >  EGap : mean = 0 eV  rms = 0 eV 
G4WT1 >  LAbs : mean = 0 fm  rms = 0 fm 
G4WT1 >  LGap : mean = 0 fm  rms = 0 fm 
G4WT1 > ... merge all h1 - done
G4WT3 >  EGap : mean = 0 eV  rms = 0 eV 
G4WT3 >  LAbs : mean = 0 fm  rms = 0 fm 
G4WT3 >  LGap : mean = 0 fm  rms = 0 fm 
G4WT0 >  ----> print histograms statistic for the local thread 
G4WT0 > 
G4WT0 >  EAbs : mean = 0 eV  rms = 0 eV 
G4WT0 >  EGap : mean = 0 eV  rms = 0 eV 
G4WT0 >  LAbs : mean = 0 fm  rms = 0 fm 
G4WT0 >  LGap : mean = 0 fm  rms = 0 fm 
G4WT2 >  LAbs : mean = 0 fm  rms = 0 fm 
G4WT2 >  LGap : mean = 0 fm  rms = 0 fm 
G4WT1 > ... merge all h2 - done
G4WT1 > ... merge all h3 - done
G4WT1 > ... merge all p1 - done
G4WT1 > ... merge all p2 - done
G4WT0 > ... merge all h1 - done
G4WT1 > ... merge slave ntuples - done
G4WT0 > ... merge all h2 - done
G4WT0 > ... merge all h3 - done
G4WT0 > ... merge all p1 - done
G4WT0 > ... merge all p2 - done
G4WT2 > ... merge all h1 - done
G4WT0 > ... merge slave ntuples - done
G4WT2 > ... merge all h2 - done
G4WT2 > ... merge all h3 - done
G4WT2 > ... merge all p1 - done
G4WT2 > ... merge all p2 - done
G4WT3 > ... merge all h1 - done
G4WT2 > ... merge slave ntuples - done
G4WT3 > ... merge all h2 - done
G4WT3 > ... merge all h3 - done
G4WT3 > ... merge all p1 - done
G4WT3 > ... merge all p2 - done
G4WT3 > ... merge slave ntuples - done
 Run terminated.
Run Summary
  Number of events processed : 100
  User=0.880000s Real=1.057696s Sys=0.320000s [Cpu=113.5%]

 ----> print histograms statistic for the entire run 

 EAbs : mean = 0 eV  rms = 0 eV 
 EGap : mean = 0 eV  rms = 0 eV 
 LAbs : mean = 0 fm  rms = 0 fm 
 LGap : mean = 0 fm  rms = 0 fm 
... merge main ntuples - done
... write file : B4.root - done
... close file : B4.root - done
100 events have been kept for refreshing and/or reviewing.
  "/vis/reviewKeptEvents" to review one by one.
  To see accumulated, "/vis/enable", then "/vis/viewer/flush" or "/vis/viewer/rebuild".
There are 4 h1 histograms
  0 with      0 entries: Edep in absorber
  1 with      0 entries: Edep in gap
  2 with      0 entries: trackL in absorber
  3 with      0 entries: trackL in gap
List them with "/analysis/list".
View them with "/vis/plot" or "/vis/reviewPlots".
WARNING: Viewpoint direction is very close to the up vector direction.
  Change the up vector or "/vis/viewer/set/rotationStyle freeRotation".



that was the result when i run 10 events, my confusion here is how can i make the set up to record the photons that pass the material and detected by my detector. i will be very grateful if i can get any help.

I think the best way to do this would be through the use of G4MultiFunctionalDetector attached to the detector volume together with the appropriate G4PS... primitive scorer classes registered with the multifunctional detector. The scored quantities can be extracted/persisted as needed in an EndOfEventAction.

Use of these is demonstrated in the B3 and B4d basic examples, with more advanced use in the RE02, RE03 and RE06 examples in extended/runAndEvent.

You meant i should go through B3, B4d and extended examples RE02, RE03 and RE06 i can understand how to develop it? But iam still new to geant4 i don’t know if you can help me pointed out in the examples what i should do or what i should modify or combine from the examples. Thank you so much for the response.

Yes, though to simplify things, let’s just focus on the B3a example for the first step of adding the scoring. Here, look at:

  1. B3a/src/DetectorConstruction.cc and the DetectorConstruction::ConstructSDandField() member function. This is where the G4MultiFunctionalDetector is constructed together with the primitive scorers. It is then attached to the logical volume where scoring is needed using the SetSensitiveDetector function, passing the name of the logical volume (so whatever you’ve called the LV for your NaI crystal) and the pointer to the MFD.
  2. B3a/src/EventAction.cc and the EventAction::EndOfEventAction member function. This is where the scored quantities for the event are extracted. In particular, see the first few lines where the hits collections for the primitive scorers are extracted and how these relate to the names these were given in Step 1. The remainder of the function deals with iterating over those collections and extracting wanted results. The best first step here is simply to print out numbers to G4cout, then decide how best to process/persist/accumulate these as needed for your use case.

The other examples noted follow on from these two key areas, demonstrating other uses/persistency cases. All of these should be looked at together with the Hits section of the Application Developer Guide, in particular the subsection on Multifunctional detectors and primitive scorers and the proceeding sections.

thanks so much for the response i will look at it and give the feedback.

i was able to follow what you have highlighted. i looked in to B3a/src/DetectorConstruction.cc and the DetectorConstruction::ConstructSDandField() member function. so also the B3a/src/EventAction.cc and the EventAction::EndOfEventAction member function. attached are the modified codes. please can you check and give me the further comments. Thank you so much i appreciate your kind guidance.
DetectorConstruction.cc (14.5 KB)

EventAction.cc (3.7 KB)
EventAction.hh (2.5 KB)

@bmorgan please come to my aid i appreciate your guidance

It looks reasonable from a quick scan, but have you compiled/run it?

yes i did that, this it the visualisation when i run 100 evennts.



from the above it looks like the particle is going in different direction. while from (case) my setup i want the particle photon to start from the source geometry then pass through the material and be detected.
looking at the PrimaryGeneratorAction.cc file the particle is chargegeantino. while for my case i want to shoot gamma photons and neutron. How to modify this to searve my purpose? furthermore how to easily make a macro file to run a high number of events?
PrimaryGeneratorAction.cc


//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

PrimaryGeneratorAction::PrimaryGeneratorAction()
{
  G4int n_particle = 1;
  fParticleGun  = new G4ParticleGun(n_particle);

  // default particle kinematic

  G4ParticleTable* particleTable = G4ParticleTable::GetParticleTable();
  G4ParticleDefinition* particle
                    = particleTable->FindParticle("chargedgeantino");
  fParticleGun->SetParticleDefinition(particle);
  fParticleGun->SetParticlePosition(G4ThreeVector(0.,0.,0.));
  fParticleGun->SetParticleEnergy(1*eV);
  fParticleGun->SetParticleMomentumDirection(G4ThreeVector(1.,0.,0.));
}

//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

PrimaryGeneratorAction::~PrimaryGeneratorAction()
{
  delete fParticleGun;
}

//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

void PrimaryGeneratorAction::GeneratePrimaries(G4Event* anEvent)
{
  G4ParticleDefinition* particle = fParticleGun->GetParticleDefinition();
  if (particle == G4ChargedGeantino::ChargedGeantino()) {
    //fluorine
    G4int Z = 9, A = 18;
    G4double ionCharge   = 0.*eplus;
    G4double excitEnergy = 0.*keV;

    G4ParticleDefinition* ion
       = G4IonTable::GetIonTable()->GetIon(Z,A,excitEnergy);
    fParticleGun->SetParticleDefinition(ion);
    fParticleGun->SetParticleCharge(ionCharge);
  }

  // randomized position
  //
  ///G4double x0  = 0*cm, y0  = 0*cm, z0  = 0*cm;
  ///G4double dx0 = 0*cm, dy0 = 0*cm, dz0 = 0*cm;
  G4double x0  = 4*cm, y0  = 4*cm, z0  = 4*cm;
  G4double dx0 = 1*cm, dy0 = 1*cm, dz0 = 1*cm;
  x0 += dx0*(G4UniformRand()-0.5);
  y0 += dy0*(G4UniformRand()-0.5);
  z0 += dz0*(G4UniformRand()-0.5);
  fParticleGun->SetParticlePosition(G4ThreeVector(x0,y0,z0));

  //create vertex
  //
  fParticleGun->GeneratePrimaryVertex(anEvent);
}

//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

}

thank you so much for your kind guidance i hope to hear from you soon.

@bmorgan still waiting for the comments thanks.

Hi @bmorgan i was able to set all the setup to my requirement. but i want to know how can i get the energy deposited/photon on the detector.
attached are the relevant important files of my codes and the visualization of the setup.
RunAction.cc (4.8 KB)
PrimaryGeneratorAction.cc (5.3 KB)
EventAction.cc (3.7 KB)
DetectorConstruction.cc (14.5 KB)
ActionInitialization.cc (2.5 KB)
RunAction.hh (2.4 KB)
PrimaryGeneratorAction.hh (2.9 KB)
EventAction.hh (2.5 KB)
DetectorConstruction.hh (2.5 KB)
ActionInitialization.hh (2.3 KB)
Also when you look at the the terminal below, why am i seeing dose on the patient, although i changed everything to my requirements. i only wanted the dose/photon deposited on the detector.


thank you so much for the guidance. Hope to hear from you soon.

I think you need to take a few steps back first, as that code is rather complex and is copy-pasting from several different examples (I’d also suggest deleting code rather than commenting it out when posting it here as this will make it easier to read). Only a few quick points:

  1. Modify the PrimaryGeneratorAction::GeneratePrimaries member function so that it simply reads

    void PrimaryGeneratorAction::GeneratePrimaries(G4Event* anEvent)
    {
      fParticleGun->GeneratePrimaryVertex(anEvent);
    }
    

    That will allow you to control particle type/position/direction/energy etc through the /gun/<cmd> UI commands at runtime through the UI or scripts.

  2. The output from “patient” is purely because you’ve copy-pasted and not removed code that’s not relevant to your use case/.

  3. To get the energy deposition in the NaI, the scoring looks set up o.k. in the detector construction, and the energy deposit, if any, should be being extracted in EventAction::EndOfEventAction(const G4Event* evt ) at the lines:

     for (itr = evtMap->GetMap()->begin(); itr != evtMap->GetMap()->end(); itr++) {
        ///G4int copyNb  = (itr->first);
        G4double edep = *(itr->second); // EDEP EXTRACTED HERE
        if (edep > eThreshold) nbOfFired++;
        ///G4cout << "\n  cryst" << copyNb << ": " << edep/keV << " keV ";
      }
    

    I would guess that for your use case, then simply accumulating the energy should be enough as there’s only one crystal e.g.

    G4double edep = 0.0;
    for (itr = evtMap->GetMap()->begin(); itr != evtMap->GetMap()->end(); itr++) {
        edep += *(itr->second); // EDEP EXTRACTED HERE
    }
    // Do whatever you need with edep here, e.g. use G4cout to output it to check

ok i will surely do that. thank you so much.

This topic was automatically closed 30 days after the last reply. New replies are no longer allowed.