2 circular targets

Getting Started
72

Thank you @anna

To store the sum of energy

  1. In B1SteppingAction.cc
    a. I wrote

    if (volume == fScoringVolume) {
    G4double edepStep = step->GetTotalEnergyDeposit();
    fEventAction->AddEdep(edepStep);
    return (edepStep);
    }
    if (volume == fScoringVolume2) {
    G4double edepStep2 = step->GetTotalEnergyDeposit();
    fEventAction->AddEdep2(edepStep2);
    return (edepStep2);
    }
    if (volume == fScoringVolume || volume == fScoringVolume2) {
    G4double edepStep2tartot = edepStep+edepStep2;
    fEventAction->AddEdeptot(edepStep2tartot);
    }

B1SteppingAction.cc (6.8 KB)

  1. In B1EventAction.cc
    a. I added fEdep2TarSum = 0.;
    b. I added fRunAction->AddEdeptot(fEdep);
    c. I added analysisManager->FillNtupleDColumn(2,0, fEdep2TarSum);
    analysisManager->AddNtupleRow(2);

B1EventAction.cc (4.3 KB)

  1. In B1RunAction.cc
    a. I added
    analysisManager->CreateNtuple("B1VolumeTot", "Volume Tot");
    analysisManager->CreateNtupleDColumn("OutputEnergyTot");
    analysisManager->FinishNtuple();
    b. I added
    void B1RunAction::AddEdeptot(G4double edep)
    {
    fEdep += edep;
    }

B1RunAction.cc (8.8 KB)

  1. In B1EventAction.hh
    a. I added void void AddEdeptot(G4double edep) { fEdep2TarSum += edep; }
    b. I added G4double fEdep2TarSum;

B1EventAction.hh (3.0 KB)

  1. In B1RunAction.hh
    a. I added void AddEdeptot (G4double edep);

B1RunAction.hh (2.9 KB)

But compiling I get these errors

image
image1105×650 28.6 KB

I understand they are due to

return (edepStep); and return (edepStep2);…

To watch what is happening

In my macro file I had the /tracking/verbose 1 line yet…what should I do to see what is happening?

run1.mac.txt (559 Bytes)

73

I am not sure why you introduced those changes. Do you need a variable that returns the sum of first and second target? What I meant is that the ntuple you fill in event action could be declared as a third ntuple, with Id=2 - which means the columns " OutputEnergy" and " OutputEnergy2" will be created within analysisManager->CreateNtuple("EventVariables", "Event");.
And yes, the errors mean you try to return a value in a method that is not supposed to return anything. Also, return means exit from the method, so your third “if” will never happen as one of the first 2 will.

Regarding the tracking verbosity, you need to analyse the output that you obtain. It should start with something like that:

G4WT0 > *********************************************************************************************************
G4WT0 > * G4Track Information:   Particle = e+,   Track ID = 1,   Parent ID = 0
G4WT0 > *********************************************************************************************************
G4WT0 > 
G4WT0 > Step#    X(mm)    Y(mm)    Z(mm) KinE(MeV)  dE(MeV) StepLeng TrackLeng  NextVolume ProcName
G4WT0 >     0        0        0   -1e+03   4.5e+04        0        0         0       World initStep
G4WT0 >     1        0        0     -609   4.5e+04   0.0805      391       391       World eIoni
74

Thank you @anna. No I don’t think to need the sum of the released energies in the two targets (my supervisor didn’t ask me it).Reading yuor previous message

As I wrote above, each case = either target1, or target 2, or event information (sum of energy in targets, etc).

I understood that you were suggesting me to do it!

Well… the if now I understand what you mean, I’ve to create a third ntupla whitin I’ve a column for the released energy in target 1 and a second column for the released energy in target 2…true?
Then

  1. In B1RunAction.cc
    a. I added
    analysisManager->CreateNtuple("B1EventVariables", "Events");
    analysisManager->CreateNtupleDColumn("ReleasedEnergy");
    analysisManager->CreateNtupleDColumn("ReleasedEnergy2");
    analysisManager->FinishNtuple();

B1RunAction.cc (8.6 KB)

  1. In B1EventAction.cc
    a. I added
    analysisManager->FillNtupleDColumn(2,0, fEdep);
    analysisManager->FillNtupleDColumn(2,1, fEdep2Tar);
    analysisManager->AddNtupleRow(2);

B1EventAction.cc (4.0 KB)

Then I compiled without errors and I see there are 3 Ntuples, the first one for the Volume1, the second one for the Volume 2 and the last one that I just created to print both relased energy in first and second target. But I think there is a problem, because, you see in the first 2 Ntuples I get 12 events, instead in the third one I get just 5 events! So I’m not filling all the events ( it’s the same problem that I had before to create a second Ntupla!)

image
image715×596 18.8 KB

image
image714×599 18.3 KB

image
image718×292 8.04 KB

Regarding the tracking verbosity, yes when I run the simulation, I get the table that you showed me. You can see it in this file
log.txt (60.7 KB)

I read the are Events having Track ID >1 and Parent ID>0, but I’m seeing that in this output I’ve 68 tables, instead in the ROOT file I’ve just 12 events in the Volume 1 and 12 event in the Volume 2…so maybe I’m not storing all the events in the Root file?

Thank you

75

You save only particles that traverse your detector. The last but one column says which volume the step enters. Check this output to see how many particles (and which) create entry in your ntuple. The row with your volume will be exactly the step that you store in the ntuple (momentum direction, energy, PDG…)

76

ah ok @anna, then I’ve just 12 events, just because the others don’t pass go through the 2 targets, but out of them, true?

Instead, what should I do to fix the problem about the third ntupla?

77

If that’s how many you counted then yes.

It’s not a problem with the 3rd ntuple, but you may create it (in run action), and create OutputEnergy and OutputEnergy columns inside it, and fill this ntuple (with probably ID=2) in the event action. Just as you did the other two ntuples.

78

@anna, could be possible that I’m storing more events than how many events I’ve?

My two targets are defined

G4Tubs* solidEnv =
new G4Tubs("Envelope", //its name
pRMin, pRMax, pDz, pSPhi, pDPhi); //its size*/

    and

G4Tubs* solidEnv2 =
new G4Tubs("Envelope2", //its name
pRMin, pRMax, pDz, pSPhi, pDPhi); //its size*/

then, their names are " Envelope" and “Envelope2” in the log file I read the rows, for the target 1

  • G4Track Information: Particle = e+, Track ID = 1, Parent ID = 0

15 -0.0131 0.00543 -1.5 4.5e+04 0.0312 146 999 Envelope Transportation

  • G4Track Information: Particle = e+, Track ID = 1, Parent ID = 0

11 -0.00536 -0.00191 -1.5 4.5e+04 0.000176 0.639 999 Envelope Transportation

  • G4Track Information: Particle = e+, Track ID = 1, Parent ID = 0
8 -0.00152   -0.004     -1.5   4.5e+04   0.0152     73.6       999    Envelope Transportation
  • G4Track Information: Particle = e+, Track ID = 1, Parent ID = 0

11 -0.00126 0.00137 -1.5 4.5e+04 0.0264 140 999 Envelope Transportation
12 -0.00128 0.0014 1.28 4.5e+04 0.595 2.78 1e+03 Envelope eIoni

  • G4Track Information: Particle = e-, Track ID = 12, Parent ID = 1
0 -0.00128   0.0014     1.28      1.65        0        0         0    Envelope initStep
  • G4Track Information: Particle = e+, Track ID = 1, Parent ID = 0

12 -0.00435 0.00198 -1.5 4.5e+04 0.00186 15.2 999 Envelope Transportation

and, for the target 2

  • G4Track Information: Particle = e+, Track ID = 1, Parent ID = 0

17 -0.0139 0.00571 21.5 4.5e+04 0.00374 20 1.02e+03 Envelope2 Transportation

  • G4Track Information: Particle = e+, Track ID = 1, Parent ID = 0

14 -0.00571 -0.00213 21.5 4.5e+04 0.00198 18.7 1.02e+03 Envelope2 Transportation
15 -0.00573 -0.00214 22.3 4.5e+04 0.159 0.774 1.02e+03 Envelope2 eIoni

  • G4Track Information: Particle = e-, Track ID = 13, Parent ID = 1
0 -0.00573 -0.00214     22.3       5.9        0        0         0   Envelope2 initStep
  • G4Track Information: Particle = e+, Track ID = 1, Parent ID = 0

10 -0.00128 -0.00456 21.5 4.5e+04 0.00537 20 1.02e+03 Envelope2 Transportation

  • G4Track Information: Particle = e+, Track ID = 1, Parent ID = 0

14 -0.00204 0.00215 21.5 4.5e+04 0.00683 20 1.02e+03 Envelope2 Transportation

  • G4Track Information: Particle = e+, Track ID = 1, Parent ID = 0

15 -0.0039 0.00171 21.5 4.5e+04 0.00222 17.9 1.02e+03 Envelope2 Transportation

Then, I read 7 steps in the target 1 and 7 steps in the target 2…but in the Volume 1 and Volume 2 of the ROOT file I’ve 12 events.…

Regarding of the third Ntupla, I create it in RunAction by the code

analysisManager->CreateNtuple("B1EventVariables", "Events");
analysisManager->CreateNtupleDColumn("ReleasedEnergy");
analysisManager->CreateNtupleDColumn("ReleasedEnergy2");
`analysisManager->FinishNtuple();

and I fill it in EventAction by

analysisManager->FillNtupleDColumn(2,0, fEdep);
analysisManager->FillNtupleDColumn(2,1, fEdep2Tar);
alysisManager->AddNtupleRow(2);

but it doesn’t fill 12 events (as in Volume1 and Volume2), but just 5 events…

79

5 events is probably correct. Isn’t that what you wrote in your macro (after /run/beamOn) ?
What you store in the other ntuples is the number of steps that particles make inside the (sensitive) volumes.

80

Yes @anna in the run1.mac I wrote /run/beamOn 5 , so do you think that it’s correct to get 12 events in the first and second ntupla and to get 5 events in the third one?

And what about the difference between the 7 events counting the ling of the tracking verbosity and the 12 events stored in the ROOT file? Why do I have 7 events in the tracking verbosity and 12 in the Volume 1 and Volume 2 ntupla of the ROOT file?
Thank you

81

5 events in the third ntuple is coming from the initial positrons (which we call 5 events).
If you say you see 7 particles within volume, it fits the scenario where you call AddNtupleRow in the stepping action (7 times = once per particle step within volume) AND in the event action (5 times = once per event). 5+7=12. You need to fill those two first ntuples only in the stepping action, not also in the event action.

82

thank you @anna then

  1. In B1Eventaction.cc
    a. I deleted the lines
    analysisManager->FillNtupleDColumn(0,5, fEdep);
    analysisManager->AddNtupleRow();
    analysisManager->FillNtupleDColumn(1,5, fEdep2Tar);
    analysisManager->AddNtupleRow(1);

B1EventAction.hh (2.7 KB)

  1. In B1SteppingAction.cc
    a. In the if (volume == fScoringVolume) I added the line
    analysisManager->FillNtupleDColumn(0,5, fEdep);
    b. In the if (volume == fScoringVolume2) I added the line
    analysisManager->FillNtupleDColumn(1,5, fEdep2Tar);

B1SteppingAction.hh (2.6 KB)

  1. In B1SteppingAction.hh
    a. I added the lines
    G4double fEdep;
    G4double fEdep2Tar;

the I compiled and runned without errors.

log.txt (97.9 KB)
In the tracking verbosity I read there are 7 steps in the Volume 1 and 7 steps in the Volume 2. Then I opened the ROOT file but you see…

image
image713×989 29.5 KB

I still get 12 events in the first two ntuples, instead I’ve 7 steps in the tracking verbosity…

Then in B1EventAction.cc

  1. I also deleted the lines
    a. fRunAction->AddEdep(fEdep);
    fRunAction->AddEdep2(fEdep);

B1EventAction.hh (2.7 KB)

I compiled and runned again without erros

log.txt (97.9 KB)

and, in the tracking verbosity, I still get 7 steps in Volume 1 and 7 steps in Volume 2. Then I opened the ROOT file and now I get just 7 events in the ROOT file too!

image
image759×849 41.1 KB

But….given that to get 7 events in the ROOT file I had to delete the lines

fRunAction->AddEdep(fEdep);
fRunAction->AddEdep2(fEdep);

is it a problem for the simulation? Will I miss something when I will start the simulations? Near this line there was this comment

accumulate statistics in run action

Thank you

83

You didn’t use the run action, it was a leftover from the example.

84

Thank you @anna. Then if I undestand,

  1. Currently in the ROOT file I’m storing each step in which there is an interaction between the target (1 or 2)
  2. I store these interactions (7 interactions because, currently I’ve a beam just constitued by only 5 e+) and I store
    a. the ID of interacting particle
    b. the kinetic energy of interacting particle
    c. The momentum (x,y,z) of interacting particles
    d. The released energy in the targets in each step by particles.

Right?

Now, to end the program, I should:

  1. Set the beam dimension, so that I can study the differences of relased energy in the targets in function of the beam dimension
  2. Store in the ROOT file the information about secondary particles after the second target; I mean
    a. The PDG ID of secondary particles outgoing from the second target
    b. The kinetic energy of the secondary particles outgoing the secondary particles

Can you help me for this too, please?

85

To modify the beam you can read in the user guide on the general particle source.
To check for particles exiting the volume, you need to introduce another check into the stepping action, I believe it is touched in another thread on the forum.

88

Thank you @anna, but before to speak about the beam size and the outgoing particle…I just notice a problem!
Moving the

analysisManager->FillNtupleDColumn(0,5, fEdep);
analysisManager->FillNtupleDColumn(1,5, fEdep2Tar);

from the B1EventAction to the B1SteppingAction to get just 7 events as I read by the tracking verbosity, I have a problem in the released energy, i.e. I get a 0GeV released energy in the two targets!
I just noticed it

image
image856×617 30.6 KB

image
image856×609 30.6 KB

So, what should I do to have 7 events, but no zero energy?

90

Does anyone know how to fix it please?
I’m sorry to update my help request, but this is my first Gean4 simulation and I need it to be admitted to next PhD year…so if someone can help me to finish it, please do it. Currently I just need

  1. to fix this problem of 0 GeV released energy;
  2. to get information about outgoing particles
  3. to set the beam size
    and the program will be finished.
91

For the answer to 2) and 3) please see my previous post.
Regarding 1), I cannot reproduce your problem. Attaching header files does not help much with the code, it’s the source files that you edit (the ones ending with .cc)
However, let me show explicitly how it can be done. Energy accumulated in the target should not be filled within step, but rather after all the particles are simulated (it’s a sum of energy of all particles within an event). Therefore, they should be filled from event action. To do so you need to create a separate (third) ntuple, as I wrote in many previous messages:

  analysisManager->CreateNtuple("B1Event", "Event");
  analysisManager->CreateNtupleDColumn("OutputEnergy");
  analysisManager->CreateNtupleDColumn("OutputEnergy2");
  analysisManager->FinishNtuple();

(removing “OutputEnergy” and “OutputEnergy2” from other ntuples)
And fill it in in ::EndOfEventAction(...):

   auto analysisManager = G4AnalysisManager::Instance();
  analysisManager->FillNtupleDColumn(2,0, fEdep);
  analysisManager->FillNtupleDColumn(2,1, fEdep2Tar);
  analysisManager->AddNtupleRow(2);

Of course you need to keep accumulating the energy in the stepping action:

  G4double edepStep = step->GetTotalEnergyDeposit();
  if (volume == fScoringVolume) {
    fEventAction->AddEdep(edepStep);
  } else if (volume == fScoringVolume2) {
    fEventAction->AddEdep2(edepStep);
  }

That way you have an ntuple that should have as many entries as events that you run (5).

93

Thank you @anna ,I know you replied me about the outgoing particles and the beam size, but I first want to solve the 0 energy problem.

I created the third ntupla to store just the released energy in the two targets (I showed you it). But, as you wrote, in this way I just store the released energy by primary particles. Instead, as I wrote previous, my supervisor wants to store the released energy not just by primary particles, but also by secondary particle…this is my problem and this is the reason because you wrote me to see what happen in the verbosity tracker and to move the storing energy from the event action to the stepping action.
So, do you know how can I store both primary and secondary particles released energy?

Ps. Here there is all the code https://we.tl/t-iYl3DPpNFg

94

Have you looked at any of the basic examples? I am fairly sure that accumulating results for a full event (not just one step), as well as accumulating totals across many events in a run, are part of those examples.

98

Hi @mkelsey , thanks. I’m modifying the B1 example and @anna was helping me to do it.

My supervisror wants that the program

  1. Store the released energy in my two targets both by primary and secondary particles.
  2. Store the PDG particle ID of the particle releasing energy
  3. Store the kinetic energy of the particle releasing energy
  4. Store the directions of the particle releasing energy
  5. Store PDG particle ID of the particle outgoing the second target
  6. Store kinetic energy of the particle outgoing the second target

Let’s go step by step and let’s think first about what happens inside the two targets for now.

I think that in the B1 example the accumulating energy is for the event, because in the B1SteppingAction.cc I have the code:

  1. For the first target:
    if (volume == fScoringVolume) { G4double edepStep = step->GetTotalEnergyDeposit(); fEventAction->AddEdep(edepStep); }

  2. For the second target:

     else if (volume == fScoringVolume2) {
   G4double edepStep2 = step->GetTotalEnergyDeposit();
   fEventAction->AddEdep2(edepStep2);
 	}

then I think that I’m accumulating energy in the full event. Am I wrong?

Given that I need both the released energy by primary and secondary particles in the target, @anna said me to create a first Ntupla to store in the first one the released energy in the first target and in the second one the released energy in the second target.

Here the plots about the PDG, the kinetic energy, momentum and the released energy in the 2 targets:

TARGET 1
Particle ID

image
image873×620 28.4 KB

Kinetic Energy
image
image854×621 31.8 KB

Momentum x
image
image856×616 31.4 KB

Momentum y
image
image850×613 32.3 KB

Momentum z
image
image859×621 32.9 KB

Released Energy
image
image857×606 33.9 KB

TARGET 2

Particle ID

image
image852×616 28.2 KB

Kinetic Energy
image
image859×621 32.7 KB

Momentum x
image
image857×625 35.7 KB

Momentum y
image
image863×619 36.8 KB

Momentum z
image
image859×615 34.3 KB

Released energy
image
image855×621 34.2 KB

Then in the ROOT plots I read that I have 12 particles (i.e 6 e+, 5 photons and 1 e-) in each target. But it looks like to me that it doesn’t fit with the verbosity tracker

verbosity.txt (60.6 KB)

Or not? @anna and @mkelsey what do you think? Are the result fitting with the verbosity tracker or should I fix something?

@anna, last week, wrote me

image
image1061×353 26.7 KB

but
if I move the lines
analysisManager->FillNtupleDColumn(0,5, fEdep);
analysisManager->FillNtupleDColumn(1,5, fEdep2Tar);
In the stepping action, I get 0 GeV released energy!

She also wrote me

image
image764×219 12.3 KB

but in this way I just store the released energy by prirmary particles (not the released enrergy by primary and secondary)…or not? @anna did I understand wrotng your message?

What I mean is:
If my beam is constituted by five 45GeV e+ and hitting the first target, one of them produce a couple e+e-, I must accumulate both the released energy by the 5 primary e+, and the secondary e+ and e-. Does my code is doing it correctly?