TargetDarkBremFilter.cxx

/*
 * @file TargetDarkBremFilter.cxx
 * @class TargetDarkBremFilter
 * @brief Class defining a UserActionPlugin that allows a user to filter out
 *        events that don't result in a dark brem within a given volume
 * @author Michael Revering, University of Minnesota
 * @author Tom Eichlersmith, University of Minnesota
 */
 
#include "Biasing/TargetDarkBremFilter.h"
 
#include "G4DarkBreM/G4APrime.h"  //checking if particles match A'
#include "G4Electron.hh"          //to check if track is electron
 
namespace biasing {
 
TargetDarkBremFilter::TargetDarkBremFilter(
    const std::string& name, framework::config::Parameters& parameters)
    : simcore::UserAction(name, parameters) {
  threshold_ = parameters.get<double>("threshold");
}
 
void TargetDarkBremFilter::beginOfEventAction(const G4Event*) {
  found_aprime_ = false;
}
 
void TargetDarkBremFilter::stepping(const G4Step* step) {
  // don't process if event is already aborted
  if (G4EventManager::GetEventManager()->GetConstCurrentEvent()->IsAborted())
    return;
 
  // Get the track associated with this step.
  auto track{step->GetTrack()};
 
  // Leave if track is not an electron
  auto particle_def{track->GetParticleDefinition()};
  if (particle_def != G4Electron::Electron()) {
    if (particle_def == G4APrime::APrime() and
        track->GetCurrentStepNumber() == 1) {
      if (isOutsideTargetRegion(track->GetLogicalVolumeAtVertex()))
        abortEvent("A' was not created within target.");
      // A' was found and originated in correct region
      found_aprime_ = true;
    }  // first step of A'
    return;
  }
 
  // Leave if track is not primary
  if (track->GetParentID() != 0) return;
 
  // Leave if track is not in target region
  if (isOutsideTargetRegion(track->GetVolume())) return;
 
  if (isOutsideTargetRegion(track->GetNextVolume())  // leaving target region
      or
      track->GetTrackStatus() == fStopAndKill  // stopping within target region
      or track->GetKineticEnergy() == 0.       // stopping within target region
  ) {
    // Get the electron secondries
    const std::vector<G4Track*>* secondaries{step->GetSecondary()};
    if (not secondaries or secondaries->size() == 0) {
      abortEvent("Primary electron did not create secondaries.");
      return;
    } else {
      // check secondaries to see if we made a dark brem
      for (auto& secondary_track : *secondaries) {
        if (secondary_track->GetParticleDefinition() == G4APrime::APrime()) {
          // we found an A', woo-hoo!
 
          if (secondary_track->GetTotalEnergy() < threshold_) {
            abortEvent(
                "A' was not created with total energy above input threshold.");
            return;
          }
 
          if (isOutsideTargetRegion(secondary_track->GetVolume())) {
            abortEvent("A' was not created within target.");
            return;
          }
 
          /* debug message
          std::cout << "[ TargetDarkBremFilter ] : "
              << "(" <<
          G4EventManager::GetEventManager()->GetConstCurrentEvent()->GetEventID()
          << ") "
              << "Found an A' weighted " << secondary_track->GetWeight()
              << " originating in " <<
          secondary_track->GetVolume()->GetLogicalVolume()->GetName()
              << " with energy " << secondary_track->GetTotalEnergy() << " MeV."
              << std::endl;
           */
 
          // we found a good A', so we can leave
          return;
        }  // this secondary is an A'
      }  // loop through secondaries
    }  // are there secondaries to loop through
 
    // got here without finding A'
    abortEvent("Primary electron did not create A'.");
    return;
  }  // should check secondaries of primary
 
  // get here if we are the primary in the target region,
  //  but shouldn't check the secondaries yet
  return;
}
 
void TargetDarkBremFilter::endOfEventAction(const G4Event* event) {
  if (not event->IsAborted() and not found_aprime_) {
    abortEvent("Did not find an A' after entire event was simulated.");
  }
}
 
}  // namespace biasing
 
DECLARE_ACTION(biasing::TargetDarkBremFilter)