LDMX Software
ldmx::PhotonuclearInteraction Class Reference

Stores detailed information about a photonuclear interaction. More...

#include <PhotonuclearInteraction.h>

Classes

struct  ParticleInfo
 Stores kinematic and identity information for a particle. More...
 

Public Member Functions

 PhotonuclearInteraction ()=default
 Constructor.
 
virtual ~PhotonuclearInteraction ()=default
 Destructor.
 
void clear ()
 Clear the interaction data.
 
void setIncidentPhoton (int track_id, int pdg_id, double energy, double px, double py, double pz, double x, double y, double z, double time)
 Set the incident photon information.
 
ParticleInfo getIncidentPhoton () const
 Get the incident photon information.
 
void setTarget (int Z, int A, const std::string &material)
 Set the target nucleus information.
 
int getTargetZ () const
 Get the target atomic number.
 
int getTargetA () const
 Get the target mass number.
 
std::string getTargetMaterial () const
 Get the target material name.
 
void addImmediateSecondary (const ParticleInfo &particle)
 Add an immediate secondary particle from the PN interaction.
 
std::vector< ParticleInfogetImmediateSecondaries () const
 Get all immediate secondary particles.
 
int getNumImmediateSecondaries () const
 Get the number of immediate secondary particles.
 
void addDescendant (int immediate_secondary_id, int final_state_track_id)
 Add a final state descendant track for an immediate secondary.
 
std::map< int, std::vector< int > > getDescendantMap () const
 Get the descendant map.
 
std::vector< int > getDescendants (int immediate_secondary_id) const
 Get descendants for a specific immediate secondary.
 
void setInteractionVolume (const std::string &volume)
 Set the volume where the interaction occurred.
 
std::string getInteractionVolume () const
 Get the interaction volume name.
 
void setProcessName (const std::string &process)
 Set the process name.
 
std::string getProcessName () const
 Get the process name.
 
void print () const
 Print information about this interaction.
 

Private Member Functions

 ClassDef (PhotonuclearInteraction, 1)
 ROOT dictionary generation.
 

Private Attributes

ParticleInfo incident_photon_
 Incident photon information.
 
int target_z_ {0}
 Target atomic number.
 
int target_a_ {0}
 Target mass number.
 
std::string target_material_ {""}
 Target material name.
 
std::vector< ParticleInfoimmediate_secondaries_
 Immediate secondary particles from the cascade.
 
std::map< int, std::vector< int > > descendant_map_
 Maps immediate secondary track ID -> final state track IDs.
 
std::string interaction_volume_ {""}
 Volume where interaction occurred.
 
std::string process_name_ {""}
 Process name (e.g., "photonNuclear")
 

Friends

std::ostream & operator<< (std::ostream &o, const PhotonuclearInteraction &pn)
 Stream operator for printing.
 

Detailed Description

Stores detailed information about a photonuclear interaction.

This class records comprehensive details about photonuclear (PN) interactions to enable validation of cascade mechanisms, identification of edge cases, and analysis of final state particle production.

For each PN interaction, it captures:

  • Incident photon information (energy, momentum, position, time)
  • Target nucleus (Z, A, material)
  • All immediate secondary particles from the cascade
  • Mapping from cascade products to their final state descendants

Definition at line 32 of file PhotonuclearInteraction.h.

Member Function Documentation

◆ addDescendant()

void ldmx::PhotonuclearInteraction::addDescendant ( int immediate_secondary_id,
int final_state_track_id )

Add a final state descendant track for an immediate secondary.

This builds the map from immediate cascade products to their eventual final state particles.

Parameters
immediate_secondary_idTrack ID of the immediate secondary
final_state_track_idTrack ID of the final state descendant

Definition at line 50 of file PhotonuclearInteraction.cxx.

51 {
52 descendant_map_[immediate_secondary_id].push_back(final_state_track_id);
53}
std::map< int, std::vector< int > > descendant_map_
Maps immediate secondary track ID -> final state track IDs.

References descendant_map_.

◆ addImmediateSecondary()

void ldmx::PhotonuclearInteraction::addImmediateSecondary ( const ParticleInfo & particle)

Add an immediate secondary particle from the PN interaction.

Parameters
particleParticleInfo for the secondary

Definition at line 45 of file PhotonuclearInteraction.cxx.

46 {
47 immediate_secondaries_.push_back(particle);
48}
std::vector< ParticleInfo > immediate_secondaries_
Immediate secondary particles from the cascade.

References immediate_secondaries_.

Referenced by simcore::PhotonuclearTracker::stepping().

◆ clear()

void ldmx::PhotonuclearInteraction::clear ( )

Clear the interaction data.

Definition at line 10 of file PhotonuclearInteraction.cxx.

10 {
11 incident_photon_ = ParticleInfo();
12 target_z_ = 0;
13 target_a_ = 0;
16 descendant_map_.clear();
18 process_name_ = "";
19}
std::string target_material_
Target material name.
ParticleInfo incident_photon_
Incident photon information.
std::string process_name_
Process name (e.g., "photonNuclear")
std::string interaction_volume_
Volume where interaction occurred.

References descendant_map_, immediate_secondaries_, incident_photon_, interaction_volume_, process_name_, target_a_, target_material_, and target_z_.

◆ getDescendantMap()

std::map< int, std::vector< int > > ldmx::PhotonuclearInteraction::getDescendantMap ( ) const
inline

Get the descendant map.

Maps immediate secondary track ID -> list of final state track IDs

Returns
Map of secondary IDs to their descendant IDs

Definition at line 177 of file PhotonuclearInteraction.h.

177 {
178 return descendant_map_;
179 }

References descendant_map_.

◆ getDescendants()

std::vector< int > ldmx::PhotonuclearInteraction::getDescendants ( int immediate_secondary_id) const

Get descendants for a specific immediate secondary.

Parameters
immediate_secondary_idTrack ID of the immediate secondary
Returns
Vector of track IDs for descendants (empty if none)

Definition at line 55 of file PhotonuclearInteraction.cxx.

56 {
57 auto it = descendant_map_.find(immediate_secondary_id);
58 if (it != descendant_map_.end()) {
59 return it->second;
60 }
61 return std::vector<int>();
62}

References descendant_map_.

◆ getImmediateSecondaries()

std::vector< ParticleInfo > ldmx::PhotonuclearInteraction::getImmediateSecondaries ( ) const
inline

Get all immediate secondary particles.

Returns
Vector of ParticleInfo for all immediate secondaries

Definition at line 146 of file PhotonuclearInteraction.h.

146 {
148 }

References immediate_secondaries_.

◆ getIncidentPhoton()

ParticleInfo ldmx::PhotonuclearInteraction::getIncidentPhoton ( ) const
inline

Get the incident photon information.

Returns
ParticleInfo for the incident photon

Definition at line 102 of file PhotonuclearInteraction.h.

102{ return incident_photon_; }

References incident_photon_.

◆ getInteractionVolume()

std::string ldmx::PhotonuclearInteraction::getInteractionVolume ( ) const
inline

Get the interaction volume name.

Returns
Volume name

Definition at line 203 of file PhotonuclearInteraction.h.

203{ return interaction_volume_; }

References interaction_volume_.

◆ getNumImmediateSecondaries()

int ldmx::PhotonuclearInteraction::getNumImmediateSecondaries ( ) const
inline

Get the number of immediate secondary particles.

Returns
Number of immediate secondaries

Definition at line 155 of file PhotonuclearInteraction.h.

155 {
156 return immediate_secondaries_.size();
157 }

References immediate_secondaries_.

◆ getProcessName()

std::string ldmx::PhotonuclearInteraction::getProcessName ( ) const
inline

Get the process name.

Returns
Process name

Definition at line 217 of file PhotonuclearInteraction.h.

217{ return process_name_; }

References process_name_.

◆ getTargetA()

int ldmx::PhotonuclearInteraction::getTargetA ( ) const
inline

Get the target mass number.

Returns
Mass number A

Definition at line 125 of file PhotonuclearInteraction.h.

125{ return target_a_; }

References target_a_.

◆ getTargetMaterial()

std::string ldmx::PhotonuclearInteraction::getTargetMaterial ( ) const
inline

Get the target material name.

Returns
Material name

Definition at line 132 of file PhotonuclearInteraction.h.

132{ return target_material_; }

References target_material_.

◆ getTargetZ()

int ldmx::PhotonuclearInteraction::getTargetZ ( ) const
inline

Get the target atomic number.

Returns
Atomic number Z

Definition at line 118 of file PhotonuclearInteraction.h.

118{ return target_z_; }

References target_z_.

◆ print()

void ldmx::PhotonuclearInteraction::print ( ) const

Print information about this interaction.

Definition at line 86 of file PhotonuclearInteraction.cxx.

86{ std::cout << *this; }

◆ setIncidentPhoton()

void ldmx::PhotonuclearInteraction::setIncidentPhoton ( int track_id,
int pdg_id,
double energy,
double px,
double py,
double pz,
double x,
double y,
double z,
double time )

Set the incident photon information.

Parameters
track_idGeant4 track ID
energyTotal energy [MeV]
pxMomentum x-component [MeV]
pyMomentum y-component [MeV]
pzMomentum z-component [MeV]
xPosition x-component [mm]
yPosition y-component [mm]
zPosition z-component [mm]
timeGlobal time [ns]

Definition at line 21 of file PhotonuclearInteraction.cxx.

References ldmx::PhotonuclearInteraction::ParticleInfo::energy_, incident_photon_, ldmx::PhotonuclearInteraction::ParticleInfo::pdg_id_, ldmx::PhotonuclearInteraction::ParticleInfo::px_, ldmx::PhotonuclearInteraction::ParticleInfo::time_, ldmx::PhotonuclearInteraction::ParticleInfo::track_id_, and ldmx::PhotonuclearInteraction::ParticleInfo::x_.

Referenced by simcore::PhotonuclearTracker::stepping().

◆ setInteractionVolume()

void ldmx::PhotonuclearInteraction::setInteractionVolume ( const std::string & volume)
inline

Set the volume where the interaction occurred.

Parameters
volumeVolume name

Definition at line 194 of file PhotonuclearInteraction.h.

194 {
195 interaction_volume_ = volume;
196 }

References interaction_volume_.

Referenced by simcore::PhotonuclearTracker::stepping().

◆ setProcessName()

void ldmx::PhotonuclearInteraction::setProcessName ( const std::string & process)
inline

Set the process name.

Parameters
processProcess name (e.g., "photonNuclear")

Definition at line 210 of file PhotonuclearInteraction.h.

210{ process_name_ = process; }

References process_name_.

Referenced by simcore::PhotonuclearTracker::stepping().

◆ setTarget()

void ldmx::PhotonuclearInteraction::setTarget ( int Z,
int A,
const std::string & material )

Set the target nucleus information.

Parameters
ZAtomic number
AMass number
materialMaterial name

Definition at line 38 of file PhotonuclearInteraction.cxx.

39 {
40 target_z_ = Z;
41 target_a_ = A;
42 target_material_ = material;
43}

References target_a_, target_material_, and target_z_.

Referenced by simcore::PhotonuclearTracker::stepping().

Friends And Related Symbol Documentation

◆ operator<<

std::ostream & operator<< ( std::ostream & o,
const PhotonuclearInteraction & pn )
friend

Stream operator for printing.

Definition at line 64 of file PhotonuclearInteraction.cxx.

64 {
65 o << "PhotonuclearInteraction:" << " Process: " << pn.process_name_
66 << " Volume: " << pn.interaction_volume_
67 << " Incident photon: track_id=" << pn.incident_photon_.track_id_
68 << " pdg=" << pn.incident_photon_.pdg_id_
69 << " E=" << pn.incident_photon_.energy_ << " MeV"
70 << " Target: Z=" << pn.target_z_ << " A=" << pn.target_a_
71 << " material=" << pn.target_material_
72 << " Immediate secondaries: " << pn.immediate_secondaries_.size()
73 << std::endl;
74 for (const auto& sec : pn.immediate_secondaries_) {
75 o << " track_id=" << sec.track_id_ << " pdg=" << sec.pdg_id_
76 << " E=" << sec.energy_ << " MeV" << std::endl;
77 }
78 o << " Descendant map entries: " << pn.descendant_map_.size() << std::endl;
79 for (const auto& [sec_id, descendants] : pn.descendant_map_) {
80 o << " Secondary " << sec_id << " -> " << descendants.size()
81 << " descendants" << std::endl;
82 }
83 return o;
84}

Member Data Documentation

◆ descendant_map_

std::map<int, std::vector<int> > ldmx::PhotonuclearInteraction::descendant_map_
private

Maps immediate secondary track ID -> final state track IDs.

Definition at line 247 of file PhotonuclearInteraction.h.

Referenced by addDescendant(), clear(), getDescendantMap(), and getDescendants().

◆ immediate_secondaries_

std::vector<ParticleInfo> ldmx::PhotonuclearInteraction::immediate_secondaries_
private

Immediate secondary particles from the cascade.

Definition at line 244 of file PhotonuclearInteraction.h.

Referenced by addImmediateSecondary(), clear(), getImmediateSecondaries(), and getNumImmediateSecondaries().

◆ incident_photon_

ParticleInfo ldmx::PhotonuclearInteraction::incident_photon_
private

Incident photon information.

Definition at line 232 of file PhotonuclearInteraction.h.

Referenced by clear(), getIncidentPhoton(), and setIncidentPhoton().

◆ interaction_volume_

std::string ldmx::PhotonuclearInteraction::interaction_volume_ {""}
private

Volume where interaction occurred.

Definition at line 250 of file PhotonuclearInteraction.h.

250{""};

Referenced by clear(), getInteractionVolume(), and setInteractionVolume().

◆ process_name_

std::string ldmx::PhotonuclearInteraction::process_name_ {""}
private

Process name (e.g., "photonNuclear")

Definition at line 253 of file PhotonuclearInteraction.h.

253{""};

Referenced by clear(), getProcessName(), and setProcessName().

◆ target_a_

int ldmx::PhotonuclearInteraction::target_a_ {0}
private

Target mass number.

Definition at line 238 of file PhotonuclearInteraction.h.

238{0};

Referenced by clear(), getTargetA(), and setTarget().

◆ target_material_

std::string ldmx::PhotonuclearInteraction::target_material_ {""}
private

Target material name.

Definition at line 241 of file PhotonuclearInteraction.h.

241{""};

Referenced by clear(), getTargetMaterial(), and setTarget().

◆ target_z_

int ldmx::PhotonuclearInteraction::target_z_ {0}
private

Target atomic number.

Definition at line 235 of file PhotonuclearInteraction.h.

235{0};

Referenced by clear(), getTargetZ(), and setTarget().


The documentation for this class was generated from the following files: