Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

calibration target upstream of everything in ldmx-sw #5

Closed
2 of 3 tasks
tomeichlersmith opened this issue Jan 10, 2024 · 9 comments
Closed
2 of 3 tasks

calibration target upstream of everything in ldmx-sw #5

tomeichlersmith opened this issue Jan 10, 2024 · 9 comments

Comments

@tomeichlersmith
Copy link
Owner

tomeichlersmith commented Jan 10, 2024

Another calibration possibility is to just put the thicker tungsten target in front of everything that is currently apart of the LDMX detector apparatus. This would resolve a few issues we are seeing

  • Allows longer baseline for the muons to spread out to the edges of the ECal flower
  • Gives opportunity to align LDMX tracker systems with (hopefully cleaner) muon pairs
  • Doing the simulation with ldmx-sw allows us to see the effect of the magnetic field and the partitioning of the ECal modules into their cells

To Dos

If neither of these require C++ changes, we could do an initial study of this proposal without having to edit ldmx-sw.

  • update detector geometry: hopefully, I can just drop a hunk into detector.gdml and have the rest stay the same, but if the world volume is not big enough, I may need to enlarge that as well.
  • biasing/filtering ecosystem: need to make sure I can select this upstream hunk for the biasing/filtering ecosystem. I suspect this to not be too difficult as long as I name the geometry in a good way. Hopefully, I won't need to change the C++ either in this case.
  • turning of the magnetic field: we may want to be able to run calibration without the magnetic field. Need to do this in the simulation.
@tomeichlersmith
Copy link
Owner Author

So, I've been able to write a detector.gdml that includes a hunk of calibration tungsten material upstream of everything, throw electrons at it, and generate muons within it using the biasing/filtering system already in existence.

The issue is that when I do this procedure, all events, including ones with actual muons produced, are aborted. I caved and did a debug build of ldmx-sw and was able to observe the following (normal running printout omitted).

gdb setup

ldmx gdb
file ../ldmx-sw/install/bin/fire
run calib-hunk.py 1 # run once quick to load the libraries
break G4RunManager::AbortEvent
break G4Event::SetEventAborted
run calib-hunk.py 1

normal bt output from abort due to lack of muon energy

[ MidShowerDiMuonBkgdFilter ]: (5) Not enough energy went to the input process. Aborting event.

Breakpoint 3, G4RunManager::AbortEvent (this=0x562568d44240) at /src/source/run/src/G4RunManager.cc:625
625     /src/source/run/src/G4RunManager.cc: No such file or directory.
(gdb) bt
#0  G4RunManager::AbortEvent (this=0x562568d44240) at /src/source/run/src/G4RunManager.cc:625
#1  0x00007fd48f843ee7 in biasing::MidShowerDiMuonBkgdFilter::AbortEvent (this=0x56256d4ce950, reason="Not enough energy went to the input process.") at /home/tom/code/ldmx/ldmx-sw/Biasing/src/Biasing/MidShowerDiMuonBkgdFilter.cxx:114
#2  0x00007fd48f843d12 in biasing::MidShowerDiMuonBkgdFilter::NewStage (this=0x56256d4ce950) at /home/tom/code/ldmx/ldmx-sw/Biasing/src/Biasing/MidShowerDiMuonBkgdFilter.cxx:85

abnormal output

[ PartialEnergySorter ] : Classifying track 7 with energy 6740.01 MeV.

Breakpoint 3, G4RunManager::AbortEvent (this=0x562568d44240) at /src/source/run/src/G4RunManager.cc:625
625     /src/source/run/src/G4RunManager.cc: No such file or directory.
(gdb) bt
#0  G4RunManager::AbortEvent (this=0x562568d44240) at /src/source/run/src/G4RunManager.cc:625
#1  0x00007fd49290e214 in G4ExceptionHandler::Notify (this=<optimized out>, originOfException=<optimized out>, 
    exceptionCode=<optimized out>, severity=<optimized out>, description=<optimized out>)
    at /src/source/run/src/G4ExceptionHandler.cc:117
#2  0x00007fd49073c748 in G4Exception (originOfException=0x7fd490bca3c8 "G4VParticleChange::CheckSecondary", 
    exceptionCode=0x7fd490bca50f "TRACK001", severity=EventMustBeAborted, 
    description=0x7fd490bca3a0 "Secondary with illegal energy/momentum ")
    at /src/source/global/management/src/G4Exception.cc:51
#3  0x00007fd490bc698d in G4VParticleChange::CheckSecondary (this=0x56256d476d80, aTrack=...)
    at /src/source/track/src/G4VParticleChange.cc:482
#4  0x00007fd490bc72f5 in G4VParticleChange::AddSecondary (this=0x56256d476d80, aTrack=0x562571c188d0)
    at /src/source/track/src/G4VParticleChange.cc:173
#5  0x00007fd4910c64ed in G4ParticleChangeForOccurenceBiasing::StealSecondaries (this=0x56256d476d80)
    at /src/source/processes/biasing/generic/src/G4ParticleChangeForOccurenceBiasing.cc:51
#6  0x00007fd4910b94d8 in G4BiasingProcessInterface::PostStepDoIt (this=0x56256d4dbec0, track=..., step=...)
    at /src/source/processes/biasing/generic/src/G4BiasingProcessInterface.cc:497

Conclusion

It looks like, concerningly, that the gamma conversion to muons process is producing illegal energy/momenta.

#2  0x00007fd49073c748 in G4Exception (originOfException=0x7fd490bca3c8 "G4VParticleChange::CheckSecondary", 
    exceptionCode=0x7fd490bca50f "TRACK001", severity=EventMustBeAborted, 
    description=0x7fd490bca3a0 "Secondary with illegal energy/momentum ")

The G4VParticleChange::CheckSecondary function throws this exception if

  1. The momentum direction is not a unit vector
  2. The KE is negative
  3. The track time is going backwards relative to parent time

@tomeichlersmith
Copy link
Owner Author

I doubt its (1) since the momentum directions are computed using cosine and sines, but it could be (2). I doubt (3) since this process does nothing to alter the default time behavior.

https://github.com/LDMX-Software/geant4/blob/428860e83aad9b562dd835fe5b924bdb9df6d448/source/processes/electromagnetic/highenergy/src/G4GammaConversionToMuons.cc#L471-L493

@tomeichlersmith
Copy link
Owner Author

I am very puzzled on how this error comes to be.

The G4VParticleChange::CheckSecondary is only called if debugFlag is true. debugFlag is constructed with a value of false and is only set to true if

  • G4VParticleChange::SetDebugFlag is called
  • Geant4 was compiled with G4VERBOSE defined at compile time
  • We are copying the data from a G4VParticleChange which has debugFlag set to true

I cannot find any place within the geant4 source code that changes the debug flag via the public member function

tom@zuko:~/code/ldmx/geant4$ rg --no-ignore SetDebugFlag
source/track/include/G4VParticleChange.hh
303:    void   SetDebugFlag();

source/track/include/G4VParticleChange.icc
296: void G4VParticleChange::SetDebugFlag()

nor any other place that directly modifies the (protected) member variabled

tom@zuko:~/code/ldmx/geant4$ rg --no-ignore debugFlag
source/geometry/magneticfield/include/G4MagIntegratorDriver.hh
167:                              G4int     debugFlag);

source/track/include/G4VParticleChange.hh
314:    G4bool   debugFlag;

source/track/src/G4ParticleChangeForGamma.cc
53:  debugFlag = false;

source/track/src/G4ParticleChangeForTransport.cc
202:  if (debugFlag) CheckIt(*aTrack);

source/track/include/G4VParticleChange.icc
292:  debugFlag = false;
298:  debugFlag = true; #SetDebugFlag function definition
304:  return debugFlag;

source/track/src/G4ParticleChange.cc
341:  if (debugFlag) CheckIt(*aTrack);
394:  if (debugFlag) CheckIt(*aTrack);
434:  if (debugFlag) CheckIt(*aTrack);

source/track/src/G4ParticleChangeForDecay.cc
177:  if (debugFlag) CheckIt(*aTrack);

source/track/src/G4VParticleChange.cc
64:   debugFlag(false)
68:  debugFlag = true; # constructor wrapped with #ifdef G4VERBOSE
105:   debugFlag(right.debugFlag)
154:    debugFlag = right.debugFlag;
173:  if (debugFlag) CheckSecondary(*aTrack);

source/track/src/G4ParticleChangeForLoss.cc
57:  debugFlag = false;

source/processes/parameterisation/History
213:3) Call to CheckIt controled by the G4VParticleChange::debugFlag attribute.
419:      debugFlag
426:  G4bool debugFlag;

source/processes/parameterisation/src/G4FastStep.cc
392:  if (debugFlag) CheckIt(*aTrack);
425:  if (debugFlag) CheckIt(*aTrack);

There is no reference to the G4VERBOSE compile-time flag in any of the CMakeLists

tom@zuko:~/code/ldmx/geant4$ fd CMakeLists -X rg G4VERBOSE
tom@zuko:~/code/ldmx/geant4$ 

and so the only way to enable it would be to define it at the cmake time which we don't do.

So I'm truly perplexed.

@tomeichlersmith
Copy link
Owner Author

I added a breakpoint on the G4VParticleChange::CheckSecondary function and it is only being triggered in these situations. Somehow, the debugFlag is being set to true and the track's time is negative globally.

[ PartialEnergySorter ] : Classifying track 17 with energy 4056.53 MeV.

Breakpoint 3, 0x00007fd75164fad0 in G4VParticleChange::CheckSecondary(G4Track&)@plt () from /usr/local/lib/libG4track.so
(gdb) n
Single stepping until exit from function _ZN17G4VParticleChange14CheckSecondaryER7G4Track@plt,
which has no line number information.

Breakpoint 3, G4VParticleChange::CheckSecondary (this=0x55dbc858dc30, aTrack=...) at /src/source/track/src/G4VParticleChange.cc:399
399     /src/source/track/src/G4VParticleChange.cc: No such file or directory.
(gdb) aTrack.GetParentID()
Undefined command: "aTrack.GetParentID".  Try "help".
(gdb) print aTrack.GetParentID()
$3 = 0
(gdb) print aTrack.GetTrackID()
$4 = 0
(gdb) print aTrack.GetKineticEnergy()
$5 = 2473.7504098196055
(gdb) print aTrack.GetGlobalTime()
$6 = -5.0283618986054899
(gdb) print aTrack.GetMomentumDirection()
$7 = (const G4ThreeVector &) @0x55dbcc20a8c8: {dx = 0.0086461913125100986, dy = 0.013992621393402247, dz = 0.99986471581025815}

@tomeichlersmith
Copy link
Owner Author

I can't figure out how the debugFlag ends up being true, but I did find the issue.

We time-shift primaries so that global time t=0 corresponds to the beam arriving at the target. In this calibration running mode, we end up calling a process of importance on a track with a negative global time (since its upstream of the target) which fails the checking mechanism and aborts the event.

@tomeichlersmith
Copy link
Owner Author

It seems like setting the number of secondaries modifies the debug flag member as well.

Hardware watchpoint 3: ((G4GammaConversionToMuons * const) 0x55cdfa52b050)->aParticleChange.debugFlag

Old value = false
New value = 2
G4VParticleChange::SetNumberOfSecondaries (this=0x55cdfa52b0b0, totSecondaries=2) at /home/tom/code/ldmx/dimuon/geant4/source/track/include/G4VParticleChange.icc:274
274       theListOfSecondaries->Initialize(totSecondaries);
(gdb) bt
#0  G4VParticleChange::SetNumberOfSecondaries (this=0x55cdfa52b0b0, totSecondaries=2)
    at /home/tom/code/ldmx/dimuon/geant4/source/track/include/G4VParticleChange.icc:274
#1  0x00007f9f1f2743ef in G4GammaConversionToMuons::PostStepDoIt (this=0x55cdfa52b098, aTrack=..., aStep=...)
    at /home/tom/code/ldmx/dimuon/geant4/source/processes/electromagnetic/highenergy/src/G4GammaConversionToMuons.cc:485
#2  0x00007f9f1f0ba7e7 in G4BiasingProcessInterface::PostStepDoIt (this=0x55cdfec56b30, track=..., step=...)
    at /home/tom/code/ldmx/dimuon/geant4/source/processes/biasing/generic/src/G4BiasingProcessInterface.cc:430

Not sure why this would be, but the good news is that this appears to be GammaConversion specific

@tomeichlersmith
Copy link
Owner Author

I've fleshed out the reasoning for this more thoroughly and opened an issue in our copy of Geant4.

LDMX-Software/geant4#13

My initial analysis of the code gives no evidence that this bug (or its patch) could affect physics distributions and so I hope that it can be passed into the development image with ease.

@tomeichlersmith
Copy link
Owner Author

The updated version of Geant4 is now in 4.2.1 of the development image. This enables this repository to run with a build of ldmx-sw:trunk (or recent version) within this new image without any modifications to the config. As such, I am going to close this issue and transition to generating larger samples for further analysis.

As a side note, I believe we can simply remove the magnetic field from the detector.gdml we are using and that is enough to disable it within the simulation. I have yet to test it, but I'm pretty confident that is the case.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

1 participant