generators/phhepmc/PHHepMCGenEventv1.cc

0001 #include "PHHepMCGenEventv1.h"
0002
0003 #include <HepMC/GenEvent.h>
0004 #include <HepMC/SimpleVector.h>  // for FourVector
0005
0006 #include <CLHEP/Vector/Boost.h>
0007 #include <CLHEP/Vector/LorentzRotation.h>
0008 #include <CLHEP/Vector/LorentzVector.h>
0009 #include <CLHEP/Vector/Rotation.h>
0010
0011 #include <iostream>  // for cout
0012 #include <map>       // for map
0013 #include <sstream>
0014 #include <utility>  // for swap
0015
0016 PHHepMCGenEventv1::PHHepMCGenEventv1()
0017   : m_boost_beta_vector(0, 0, 0)
0018   , m_rotation_vector(0, 0, 1)
0019   , m_rotation_angle(0)
0020
0021 {
0022 }
0023
0024 PHHepMCGenEventv1::PHHepMCGenEventv1(const PHHepMCGenEventv1& event)
0025   : PHHepMCGenEvent(event)
0026   , m_boost_beta_vector(event.get_boost_beta_vector())
0027   , m_rotation_vector(event.get_rotation_vector())
0028   , m_rotation_angle(event.get_rotation_angle())
0029   , m_psi_n(event.get_flow_psi_map())
0030 {
0031   return;
0032 }
0033
0034 PHHepMCGenEventv1& PHHepMCGenEventv1::operator=(const PHHepMCGenEventv1& event)
0035 {
0036   if (&event == this)
0037   {
0038     return *this;
0039   }
0040
0041   Reset();
0042
0043   _embedding_id = event.get_embedding_id();
0044   _isSimulated = event.is_simulated();
0045   _theEvt = new HepMC::GenEvent(*event.getEvent());
0046
0047   return *this;
0048 }
0049
0050 void PHHepMCGenEventv1::Reset()
0051 {
0052   PHHepMCGenEvent::Reset();
0053
0054   m_boost_beta_vector.set(0, 0, 0);
0055   m_rotation_vector.set(0, 0, 1);
0056   m_rotation_angle = 0;
0057   m_psi_n.clear();
0058 }
0059
0060 //_____________________________________________________________________________
0061 void PHHepMCGenEventv1::identify(std::ostream& os) const
0062 {
0063   PHHepMCGenEvent::identify(os);
0064
0065   os << " m_boost_beta_vector = (" << m_boost_beta_vector.x() << "," << m_boost_beta_vector.y() << "," << m_boost_beta_vector.z() << ") " << std::endl;
0066   os << " m_rotation_vector = (" << m_rotation_vector.x() << "," << m_rotation_vector.y() << "," << m_rotation_vector.z() << ") by " << m_rotation_angle << " rad" << std::endl;
0067   if (!m_psi_n.empty())
0068   {
0069     os << " Reaction plane angles psi_n from hijing flowAfterburner: ";
0070     for (const auto& it : m_psi_n)
0071     {
0072       os << " n=" << it.first << " psi_n=" << it.second << "; ";
0073     }
0074     os << std::endl;
0075   }
0076
0077   static const CLHEP::HepLorentzVector zp_lightcone(0, 0, 1, 1);
0078   static const CLHEP::HepLorentzVector zm_lightcone(0, 0, -1, 1);
0079
0080   os << " HepMC Frame unit light cone vector along +z axis " << zp_lightcone << " translate to lab at : "
0081      << (get_LorentzRotation_EvtGen2Lab() * zp_lightcone)
0082      << std::endl;
0083   os << " HepMC Frame unit light cone vector along -z axis  " << zm_lightcone << " translate to lab at : "
0084      << (get_LorentzRotation_EvtGen2Lab() * zm_lightcone)
0085      << std::endl;
0086
0087   return;
0088 }
0089
0090 CLHEP::HepLorentzRotation PHHepMCGenEventv1::get_LorentzRotation_EvtGen2Lab() const
0091 {
0092   CLHEP::HepBoost boost(m_boost_beta_vector.x(), m_boost_beta_vector.y(), m_boost_beta_vector.z());
0093   CLHEP::Hep3Vector axis(m_rotation_vector.x(), m_rotation_vector.y(), m_rotation_vector.z());
0094   CLHEP::HepRotation rotation(axis, m_rotation_angle);
0095
0096   return CLHEP::HepLorentzRotation(boost, rotation);
0097 }
0098
0099 CLHEP::HepLorentzRotation PHHepMCGenEventv1::get_LorentzRotation_Lab2EvtGen() const
0100 {
0101   return get_LorentzRotation_EvtGen2Lab().inverse();
0102 }
0103
0104 float PHHepMCGenEventv1::get_flow_psi(unsigned int n) const
0105 {
0106   auto it = m_psi_n.find(n);
0107   if (it != m_psi_n.end())
0108   {
0109     return it->second;
0110   }
0111
0112   std::cout << "PHHepMCGenEventv1::get_flow_psi - Warning - requested reaction plane angle psi_n for n=" << n << " does not exist. Returning 0.0" << std::endl;
0113   return 0.0F;
0114 }