Back to home page

sPhenix code displayed by LXR

 
 

    


File indexing completed on 2025-08-05 08:16:00

0001 // $Id: $
0002 
0003 /*!
0004  * \file PHHepMCGenHelper.h
0005  * \brief
0006  * \author Jin Huang <jhuang@bnl.gov>
0007  * \version $Revision:   $
0008  * \date $Date: $
0009  */
0010 
0011 #ifndef PHHEPMC_PHHEPMCGENHELPER_H
0012 #define PHHEPMC_PHHEPMCGENHELPER_H
0013 
0014 #include <CLHEP/Vector/ThreeVector.h>
0015 
0016 #include <gsl/gsl_rng.h>
0017 
0018 #include <cmath>
0019 #include <string>
0020 #include <utility>
0021 #include <vector>
0022 
0023 class PHCompositeNode;
0024 class PHHepMCGenEvent;
0025 class PHHepMCGenEventMap;
0026 
0027 namespace HepMC
0028 {
0029   class GenEvent;
0030 }
0031 
0032 /*!
0033  * \brief PHHepMCGenHelper provides service of DST upload of HepMC subevent, vertex assignment and random generator
0034  */
0035 class PHHepMCGenHelper
0036 {
0037  public:
0038   PHHepMCGenHelper();
0039   virtual ~PHHepMCGenHelper();
0040 
0041   //! supported function distributions
0042   enum VTXFUNC
0043   {
0044     //! uniform distribution with half width set via set_vertex_distribution_width()
0045     Uniform,
0046     //! normal distribution with sigma width set via set_vertex_distribution_width()
0047     Gaus
0048   };
0049 
0050   //! toss a new vertex according to a Uniform or Gaus distribution
0051   void set_vertex_distribution_function(VTXFUNC x, VTXFUNC y, VTXFUNC z, VTXFUNC t);
0052 
0053   //! set the mean value of the vertex distribution, use PHENIX units of cm, ns
0054   void set_vertex_distribution_mean(const double x, const double y, const double z, const double t);
0055 
0056   //! set the width of the vertex distribution function about the mean, use PHENIX units of cm, ns
0057   void set_vertex_distribution_width(const double x, const double y, const double z, const double t);
0058 
0059   //! embedding ID for the event
0060   //! positive ID is the embedded event of interest, e.g. jetty event from pythia
0061   //! negative IDs are backgrounds, .e.g out of time pile up collisions
0062   //! Usually, ID = 0 means the primary Au+Au collision background
0063   int get_embedding_id() const { return _embedding_id; }
0064   //
0065   //! embedding ID for the event
0066   //! positive ID is the embedded event of interest, e.g. jetty event from pythia
0067   //! negative IDs are backgrounds, .e.g out of time pile up collisions
0068   //! Usually, ID = 0 means the primary Au+Au collision background
0069   void set_embedding_id(int id) { _embedding_id = id; }
0070   //
0071   //! reuse vertex from another PHHepMCGenEvent with embedding_id = src_embedding_id Additional smearing and shift possible with set_vertex_distribution_*()
0072   void set_reuse_vertex(int src_embedding_id)
0073   {
0074     _reuse_vertex = true;
0075     _reuse_vertex_embedding_id = src_embedding_id;
0076   }
0077 
0078   //
0079   //! init interface nodes
0080   virtual int create_node_tree(PHCompositeNode *topNode);
0081 
0082   //! choice of reference version of the PHHepMCGenEvent
0083   static const PHHepMCGenEvent *get_PHHepMCGenEvent_template();
0084 
0085   //! send HepMC::GenEvent to DST tree. This function takes ownership of evt
0086   PHHepMCGenEvent *insert_event(HepMC::GenEvent *evt);
0087 
0088   const PHHepMCGenEventMap *get_geneventmap() const
0089   {
0090     return _geneventmap;
0091   }
0092 
0093   PHHepMCGenEventMap *get_geneventmap()
0094   {
0095     return _geneventmap;
0096   }
0097 
0098   gsl_rng *get_random_generator()
0099   {
0100     return RandomGenerator;
0101   }
0102 
0103   void set_geneventmap(PHHepMCGenEventMap *geneventmap)
0104   {
0105     _geneventmap = geneventmap;
0106   }
0107 
0108   //! Beam angle in lab polar coordinate.
0109   //! @param[in] beamA_theta beamA, in pair of Theta-Phi. BeamA is aimed to +z direction in the HepMC event generator's coordinate
0110   //! @param[in] beamA_phi beamA, in pair of Theta-Phi. BeamA is aimed to +z direction in the HepMC event generator's coordinate
0111   //! @param[in] beamB_theta  beamB, in pair of Theta-Phi. BeamA is aimed to -z direction in the HepMC event generator's coordinate
0112   //! @param[in] beamB_phi  beamB, in pair of Theta-Phi. BeamA is aimed to -z direction in the HepMC event generator's coordinate
0113   void set_beam_direction_theta_phi(
0114       const double beamA_theta,
0115       const double beamA_phi,
0116       const double beamB_theta,
0117       const double beamB_phi)
0118   {
0119     m_beam_direction_theta_phi = {{beamA_theta, beamA_phi}, {beamB_theta, beamB_phi}};
0120   }
0121 
0122   //! Beam angle divergence in accelerator beam coordinate.
0123   //! @param[in] beamA_divergence_h beamA, horizontal divergence Gaussian Sigma. BeamA is aimed to +z direction in the HepMC event generator's coordinate
0124   //! @param[in] beamA_divergence_v beamA, vertical divergence Gaussian Sigma. BeamA is aimed to +z direction in the HepMC event generator's coordinate
0125   //! @param[in] beamB_divergence_h beamB, horizontal divergence Gaussian Sigma. BeamA is aimed to -z direction in the HepMC event generator's coordinate
0126   //! @param[in] beamB_divergence_v beamB, vertical divergence Gaussian Sigma. BeamA is aimed to -z direction in the HepMC event generator's coordinate
0127   void set_beam_angular_divergence_hv(
0128       const double beamA_divergence_h,
0129       const double beamA_divergence_v,
0130       const double beamB_divergence_h,
0131       const double beamB_divergence_v)
0132   {
0133     m_beam_angular_divergence_hv = {{beamA_divergence_h, beamA_divergence_v}, {beamB_divergence_h, beamB_divergence_v}};
0134   }
0135 
0136   //! Central beam angle shift as linear function of longitudinal vertex position, d_shift/dz,
0137   //! which is used to represent leading order effect of crab cavity momentum kick on the beam bunch
0138   //! @param[in] beamA_h beamA, horizontal angle dh/dz. BeamA is aimed to +z direction in the HepMC event generator's coordinate
0139   //! @param[in] beamA_v beamA, vertical angle dv/dz. BeamA is aimed to +z direction in the HepMC event generator's coordinate
0140   //! @param[in] beamB_h beamB, horizontal angle dh/dz. BeamA is aimed to -z direction in the HepMC event generator's coordinate
0141   //! @param[in] beamB_v beamB, vertical angle dv/dz. BeamA is aimed to -z direction in the HepMC event generator's coordinate
0142   void set_beam_angular_z_coefficient_hv(
0143       const double beamA_h,
0144       const double beamA_v,
0145       const double beamB_h,
0146       const double beamB_v)
0147   {
0148     m_beam_angular_z_coefficient_hv = {{beamA_h, beamA_v}, {beamB_h, beamB_v}};
0149   }
0150 
0151   //! simulate bunch interaction instead of applying vertex distributions
0152   void use_beam_bunch_sim(bool b) { m_use_beam_bunch_sim = b; }
0153 
0154   //! Beam bunch geometry as 3D Gauss width
0155   //! First element is beamA, in vector of Gaussian Sigma H,V,Longitudinal
0156   //! Second element is beamB, in vector of Gaussian Sigma H,V,Longitudinal
0157   void set_beam_bunch_width(const std::vector<double> &beamA, const std::vector<double> &beamB);
0158 
0159   void CopySettings(PHHepMCGenHelper &helper);
0160 
0161   //! copy setting to helper_dest
0162   void CopySettings(PHHepMCGenHelper *helper_dest);
0163 
0164   //! copy setting from helper_src
0165   void CopyHelperSettings(PHHepMCGenHelper *helper_src);
0166 
0167   void Print(const std::string &what = "ALL") const;
0168 
0169   void PHHepMCGenHelper_Verbosity(int v) { m_verbosity = v; }
0170 
0171   int PHHepMCGenHelper_Verbosity() { return m_verbosity; }
0172 
0173  protected:
0174   //! Record the translation,boost,rotation for HepMC frame to lab frame according to collision settings
0175   void HepMC2Lab_boost_rotation_translation(PHHepMCGenEvent *genevent);
0176 
0177   //! move vertex in translation according to vertex settings
0178   void move_vertex(PHHepMCGenEvent *genevent);
0179 
0180   //! generate vertx with bunch interaction according to
0181   //! https://github.com/eic/documents/blob/d06b5597a0a89dcad215bab50fe3eefa17a097a5/reports/general/Note-Simulations-BeamEffects.pdf
0182   //! \return pair of bunch local z position for beam A and beam B
0183   std::pair<double, double> generate_vertx_with_bunch_interaction(PHHepMCGenEvent *genevent);
0184 
0185  private:
0186   gsl_rng *RandomGenerator{nullptr};
0187 
0188   double smear(const double position, const double width, VTXFUNC dist) const;
0189 
0190   //! function to convert spherical coordinate to Hep3Vector in x-y-z
0191   static CLHEP::Hep3Vector pair2Hep3Vector(const std::pair<double, double> &theta_phi);
0192 
0193   VTXFUNC _vertex_func_x{Gaus};
0194   VTXFUNC _vertex_func_y{Gaus};
0195   VTXFUNC _vertex_func_z{Gaus};
0196   VTXFUNC _vertex_func_t{Gaus};
0197 
0198   double _vertex_x{0.};
0199   double _vertex_y{0.};
0200   double _vertex_z{0.};
0201   double _vertex_t{0.};
0202 
0203   double _vertex_width_x{0.};
0204   double _vertex_width_y{0.};
0205   double _vertex_width_z{0.};
0206   double _vertex_width_t{0.};
0207 
0208   //! Beam angle in lab polar coordinate.
0209   //! First element is beamA, in pair of Theta-Phi. BeamA is aimed to +z direction in the HepMC event generator's coordinate
0210   //! Second element is beamB, in pair of Theta-Phi. BeamA is aimed to -z direction in the HepMC event generator's coordinate
0211   std::pair<std::pair<double, double>, std::pair<double, double>> m_beam_direction_theta_phi = {
0212       {0, 0},    //+z beam
0213       {M_PI, 0}  //-z beam
0214   };
0215 
0216   //! Beam angle divergence in accelerator beam coordinate.
0217   //! First element is beamA, in pair of Gaussian Sigma_H Sigma_V. BeamA is aimed to +z direction in the HepMC event generator's coordinate
0218   //! Second element is beamB, in pair of Gaussian Sigma_H Sigma_V. BeamA is aimed to -z direction in the HepMC event generator's coordinate
0219   std::pair<std::pair<double, double>, std::pair<double, double>> m_beam_angular_divergence_hv = {
0220       {0, 0},  //+z beam
0221       {0, 0}   //-z beam
0222   };
0223 
0224   //! Central beam angle shift as linear function of longitudinal vertex position, d_shift/dz,
0225   //! which is used to represent leading order effect of crab cavity momentum kick on the beam bunch
0226   //! First element is beamA, in pair of dh/dz, dv/dz. BeamA is aimed to +z direction in the HepMC event generator's coordinate
0227   //! Second element is beamB, in pair of dh/dz, dv/dz. BeamA is aimed to -z direction in the HepMC event generator's coordinate
0228   std::pair<std::pair<double, double>, std::pair<double, double>> m_beam_angular_z_coefficient_hv = {
0229       {0, 0},  //+z beam
0230       {0, 0}   //-z beam
0231   };
0232 
0233   //! positive ID is the embedded event of interest, e.g. jetty event from pythia
0234   //! negative IDs are backgrounds, .e.g out of time pile up collisions
0235   //! Usually, ID = 0 means the primary Au+Au collision background
0236   int _embedding_id{0};
0237 
0238   //! whether reuse vertex from another PHHepMCGenEvent
0239   bool _reuse_vertex{false};
0240 
0241   //! if _reuse_vertex, which embedding_id provide the source vertex. Additional smearing and shift possible with set_vertex_distribution_*()
0242   int _reuse_vertex_embedding_id{std::numeric_limits<int>::min()};
0243 
0244   //! pointer to the output container
0245   PHHepMCGenEventMap *_geneventmap{nullptr};
0246 
0247   //! verbosity
0248   int m_verbosity{0};
0249 
0250   //! simulate bunch interaction instead of applying vertex distributions
0251   bool m_use_beam_bunch_sim{false};
0252 
0253   //! Beam bunch geometry as 3D Gauss width
0254   //! First element is beamA, in vector of Gaussian Sigma H,V,Longitudinal
0255   //! Second element is beamB, in vector of Gaussian Sigma H,V,Longitudinal
0256   std::pair<std::vector<double>, std::vector<double>> m_beam_bunch_width = {
0257       {0, 0, 0},  //+z beam
0258       {0, 0, 0}   //-z beam
0259   };
0260 
0261   //! use m_beam_bunch_width to calculate horizontal and vertical collision width
0262   //! \param[in] hv_index 0: horizontal. 1: vertical
0263   double get_collision_width(unsigned int hv_index);
0264 };
0265 
0266 #endif /* PHHEPMC_PHHEPMCGENHELPER_H */