#include <GenEventBeamEffectBooster.h>

Inheritance diagram for Simulation::GenEventBeamEffectBooster:

Collaboration diagram for Simulation::GenEventBeamEffectBooster:

Public Member Functions
	GenEventBeamEffectBooster (const std::string &t, const std::string &n, const IInterface *p)
	Constructor with parameters. More...

StatusCode	initialize () override final
	Athena algtool's Hooks. More...

StatusCode	finalize () override final
	Athena algtool's Hooks. More...

StatusCode	initializeAthenaEvent ()

StatusCode	manipulate (HepMC::GenEvent &ge, const EventContext &ctx) const override final
	modifies the given GenEvent More...

Private Member Functions
StatusCode	initializeGenEvent (CLHEP::HepLorentzRotation &transform, const EventContext &ctx) const
	calculate the transformations that we want to apply to the particles in the current GenEvent More...

void	boostParticle (HepMC::GenParticlePtr &p, const CLHEP::HepLorentzRotation &transform) const
	apply boost to individual GenParticles More...

Private Attributes
ServiceHandle< IAthRNGSvc >	m_rndGenSvc {this, "RandomSvc", "AthRNGSvc"}

ATHRNG::RNGWrapper *m_randomEngine	ATLAS_THREAD_SAFE {}
	Slot-local RNG. More...

Gaudi::Property< std::string >	m_randomEngineName {this, "RandomStream", "BEAM"}
	Name of the random number stream. More...

Gaudi::Property< bool >	m_applyBoost {this, "ApplyBoost", true}

Gaudi::Property< bool >	m_applyDivergence {this, "ApplyDivergence", true}

Gaudi::Property< double >	m_sigma_px_b1 {this, "Sigma_px_b1", 20.e-6, "angular divergence in x of beam 1 [rad]"}

Gaudi::Property< double >	m_sigma_px_b2 {this, "Sigma_px_b2", 21.e-6, "angular divergence in x of beam 2 [rad]"}

Gaudi::Property< double >	m_sigma_py_b1 {this, "Sigma_py_b1", 22.e-6, "angular divergence in y of beam 1 [rad]"}

Gaudi::Property< double >	m_sigma_py_b2 {this, "Sigma_py_b2", 23.e-6, "angular divergence in y of beam 2 [rad]"}

Gaudi::Property< double >	m_xing_x_b1 {this, "HalfXing_x_b1", 19.e-6, "half-crossing in x for beam 1 at IP1, i.e. angle between beam 1 and z-axis"}

Gaudi::Property< double >	m_xing_x_b2 {this, "HalfXing_x_b2", 1.e-6, "half-crossing in x for beam 2 at IP1, i.e. angle between beam 2 and z-axis"}

Gaudi::Property< double >	m_xing_y_b1 {this, "HalfXing_y_b1", 150.e-6, "half-crossing in y for beam 1 at IP1, i.e. angle between beam 1 and z-axis"}

Gaudi::Property< double >	m_xing_y_b2 {this, "HalfXing_y_b2", -152.e-6, "half-crossing in y for beam 2 at IP1, i.e. angle between beam 2 and z-axis"}

Gaudi::Property< double >	m_dE {this, "deltaEOverE", 1.1E-4, "Delta_E/E theoretical value"}

Gaudi::Property< double >	m_pbeam1 {this, "pbeam1", 3500.0E3, "Beam 1 Momentum / MeV"}

Gaudi::Property< double >	m_pbeam2 {this, "pbeam2", 3500.0E3, "Beam 2 Momentum / MeV"}

Gaudi::Property< double >	m_beam1ParticleMass {this, "Beam1ParticleMass", CLHEP::proton_mass_c2, "Mass of particle used in beam 1"}

Gaudi::Property< double >	m_beam2ParticleMass {this, "Beam2ParticleMass", CLHEP::proton_mass_c2, "Mass of particle used in beam 2"}

Detailed Description

This tool takes a HepMC::GenEvent and applies boosts due to beam tilt, etc.

based on: https://svnweb.cern.ch/trac/atlasoff/browser/Simulation/G4Atlas/G4AtlasUtilities/trunk/src/BeamEffectTransformation.cxx

Definition at line 38 of file GenEventBeamEffectBooster.h.

Constructor & Destructor Documentation

◆ GenEventBeamEffectBooster()

Simulation::GenEventBeamEffectBooster::GenEventBeamEffectBooster	(	const std::string &	t,
		const std::string &	n,
		const IInterface *	p
	)

Constructor with parameters.

Constructor.

Definition at line 23 of file GenEventBeamEffectBooster.cxx.

     : base_class(t,n,p)
   {
   }

Member Function Documentation

◆ boostParticle()

void Simulation::GenEventBeamEffectBooster::boostParticle	(	HepMC::GenParticlePtr &	p,
		const CLHEP::HepLorentzRotation &	transform
	)		const

private

apply boost to individual GenParticles

Definition at line 151 of file GenEventBeamEffectBooster.cxx.

   {
     // Apply the same transformation for EVERY HepMC::GenParticle
     const HepMC::FourVector &mom = p->momentum();
     CLHEP::HepLorentzVector hv(mom.px(), mom.py(), mom.pz(), mom.e()); //TODO check units
     ATH_MSG_VERBOSE("BEAMBOOST initial momentum " << hv );
     hv.transform(transform);
     ATH_MSG_VERBOSE("BEAMBOOST transformed momentum " << hv);
     p->set_momentum(HepMC::FourVector(hv.px(),hv.py(),hv.pz(),hv.e())); //TODO check units
   }

◆ finalize()

StatusCode Simulation::GenEventBeamEffectBooster::finalize ( )

finaloverride

Athena algtool's Hooks.

Definition at line 47 of file GenEventBeamEffectBooster.cxx.

   {
     ATH_MSG_VERBOSE("Finalizing ...");
     return StatusCode::SUCCESS;
   }

◆ initialize()

StatusCode Simulation::GenEventBeamEffectBooster::initialize ( )

finaloverride

Athena algtool's Hooks.

Definition at line 32 of file GenEventBeamEffectBooster.cxx.

   {
     ATH_MSG_VERBOSE("Initializing ...");
     // prepare the RandonNumber generation
     ATH_CHECK(m_rndGenSvc.retrieve());
     m_randomEngine = m_rndGenSvc->getEngine(this, m_randomEngineName); //FIXME
     if (!m_randomEngine) {
       ATH_MSG_ERROR("Could not get random number engine from RandomNumberService. Abort.");
       return StatusCode::FAILURE;
     }
     return StatusCode::SUCCESS;
   }

◆ initializeAthenaEvent()

StatusCode Simulation::GenEventBeamEffectBooster::initializeAthenaEvent ( )

Definition at line 53 of file GenEventBeamEffectBooster.cxx.

   {
     // Refresh properties retrieved from the BeamCondSvc here if
     // necessary.  It may be more thread-safe to retrieve them directly
     // from the BeamCondSvc during initializeGenEvent though?
  
     // m_sigma_px_b1 = 20.E-6; // angular divergence in x of beam 1 [rad]
     // m_sigma_px_b2 = 21.E-6; // angular divergence in x of beam 2 [rad]
     // m_sigma_py_b1 = 22.E-6; // angular divergence in y of beam 1 [rad]
     // m_sigma_py_b2 = 23.E-6; // angular divergence in y of beam 2 [rad]
     // // Crossing angles, note the convention here is taken from online https://atlas-dcs.cern.ch/index.php?page=LHC_INS::INS_BPM
     // m_xing_x_b1 = 19.E-6;   // half-crossing in x for beam 1 at IP1, i.e. angle between beam 1 and z-axis
     // m_xing_x_b2 = 1.E-6;    // half-crossing in x for beam 2 at IP1, i.e. angle between beam 2 and z-axis
     // m_xing_y_b1 =  150.E-6;  // half-crossing in y for beam 1 at IP1, i.e. angle between beam 1 and z-axis
     // m_xing_y_b2 = -152.E-6; // half-crossing in y for beam 2 at IP1, i.e. angle between beam 2 and z-axis
     // m_dE = 1.1E-4; // Delta_E/E theoretical value
     // m_pbeam1 = 3500.0E3; /// @todo Get from a service
     // m_pbeam2 = 3500.0E3; /// @todo Get from a service
     // m_beam1ParticleMass =  CLHEP::proton_mass_c2;
     // m_beam2ParticleMass =  CLHEP::proton_mass_c2;
  
     return StatusCode::SUCCESS;
   }

◆ initializeGenEvent()

StatusCode Simulation::GenEventBeamEffectBooster::initializeGenEvent	(	CLHEP::HepLorentzRotation &	transform,
		const EventContext &	ctx
	)		const

private

calculate the transformations that we want to apply to the particles in the current GenEvent

Definition at line 77 of file GenEventBeamEffectBooster.cxx.

   {
     // Reset the transformation
     transform = CLHEP::HepLorentzRotation(); //TODO drop this
  
     // Prepare the random engine
     m_randomEngine->setSeed( name(), ctx );
     CLHEP::HepRandomEngine* randomEngine(m_randomEngine->getEngine(ctx));
  
     // Create beam 1 and 2 momenta, including divergence and crossing-angle effects
     const double px1 = m_pbeam1 * CLHEP::RandGaussZiggurat::shoot(randomEngine, m_xing_x_b2, m_sigma_px_b2);
     const double py1 = m_pbeam1 * CLHEP::RandGaussZiggurat::shoot(randomEngine, m_xing_y_b2, m_sigma_py_b2);
     const double pz1 = CLHEP::RandGaussZiggurat::shoot(randomEngine, sqrt(m_pbeam1*m_pbeam1 - px1*px1 - py1*py1), m_dE);
     const double e1 = sqrt(px1*px1 + py1*py1 + pz1*pz1 + m_beam1ParticleMass*m_beam1ParticleMass);
     CLHEP::HepLorentzVector pp1(px1, py1, pz1, e1);
  
     const double px2 = -m_pbeam2 * CLHEP::RandGaussZiggurat::shoot(randomEngine, m_xing_x_b1, m_sigma_px_b1); // crossing angle & divergence
     const double py2 = -m_pbeam2 * CLHEP::RandGaussZiggurat::shoot(randomEngine, m_xing_y_b1, m_sigma_py_b1);
     const double pz2 = CLHEP::RandGaussZiggurat::shoot(randomEngine, -sqrt(m_pbeam2*m_pbeam2 - px2*px2 - py2*py2), m_dE); // longitudinal component in + direction energy smeared
     const double e2 = sqrt(px2*px2 + py2*py2 + pz2*pz2 + m_beam2ParticleMass*m_beam2ParticleMass);
     CLHEP::HepLorentzVector pp2(px2, py2, pz2, e2);
  
     // Now set-up rotation-boost matrix
     const CLHEP::HepLorentzVector psum = pp1 + pp2;
     const CLHEP::HepLorentzVector& dir = pp1;
     // Boost psum back on the direction of dir, adapted from bstback & fromCMframe PYTHIA8, credit to T.Sjostrand
     const double betaX = -psum.x() / psum.t();
     const double betaY = -psum.y() / psum.t();
     const double betaZ = -psum.z() / psum.t();
     const double beta2 = betaX*betaX + betaY*betaY + betaZ*betaZ;
     const double gamma = 1. / sqrt(1. - beta2);
     const double prod1 = betaX * dir.x() + betaY * dir.y() + betaZ * dir.z();
     const double prod2 = gamma * (gamma * prod1 / (1. + gamma) + dir.t());
     // Get the angle of rotation
     const CLHEP::HepLorentzVector back(dir.x() + prod2 * betaX,
                                        dir.y() + prod2 * betaY,
                                        dir.z() + prod2 * betaZ,
                                        gamma * (dir.t() + prod1));
     const double thback = back.theta();
     const double phback = back.phi();
     // Setting up two rotation matrices, via 3x3 rep. rather than even more messy EulerAngles
     const double sph = sin(phback), cph = cos(phback), sth = sin(thback), cth = cos(thback);
     const CLHEP::HepRotation rot1(CLHEP::HepRep3x3(cph, sph, 0.0, -sph, cph, 0.0, 0.0, 0.0, 1.0));
     const CLHEP::HepRotation rot2(CLHEP::HepRep3x3(cth*cph, -sph, sth*cph, cth*sph, cph, sth*sph, -sth, 0.0, cth));
     const CLHEP::HepBoost bst(psum.px()/psum.e(), psum.py()/psum.e(), psum.pz()/psum.e());
     // Combine the two rotations and the boost; matrix multiplication order matters
     const CLHEP::HepRotation rot = rot2*rot1;
     if(m_applyBoost && m_applyDivergence) {
       transform.set(rot,bst);
     }
     else {
       if(m_applyBoost) {
         transform.set(bst);
       }
       if(m_applyDivergence) {
         transform.set(rot);
       }
     }
     return StatusCode::SUCCESS;
   }

◆ manipulate()

StatusCode Simulation::GenEventBeamEffectBooster::manipulate	(	HepMC::GenEvent &	ge,
		const EventContext &	ctx
	)		const

finaloverride

modifies the given GenEvent

modifies (displaces) the given GenEvent

Definition at line 139 of file GenEventBeamEffectBooster.cxx.

   {
     // Obtain the transformation
     CLHEP::HepLorentzRotation transform = CLHEP::HepLorentzRotation();
     ATH_CHECK( initializeGenEvent(transform, ctx) );
  
     for( auto particleIter: ge) {
       this->boostParticle(particleIter, transform);
     }
     return StatusCode::SUCCESS;
   }

Member Data Documentation

◆ ATLAS_THREAD_SAFE

ATHRNG::RNGWrapper* m_randomEngine Simulation::GenEventBeamEffectBooster::ATLAS_THREAD_SAFE {}

private

Slot-local RNG.

Definition at line 60 of file GenEventBeamEffectBooster.h.

◆ m_applyBoost

Gaudi::Property<bool> Simulation::GenEventBeamEffectBooster::m_applyBoost {this, "ApplyBoost", true}

private

Definition at line 62 of file GenEventBeamEffectBooster.h.

◆ m_applyDivergence

Gaudi::Property<bool> Simulation::GenEventBeamEffectBooster::m_applyDivergence {this, "ApplyDivergence", true}

private

Definition at line 63 of file GenEventBeamEffectBooster.h.

◆ m_beam1ParticleMass

Gaudi::Property<double> Simulation::GenEventBeamEffectBooster::m_beam1ParticleMass {this, "Beam1ParticleMass", CLHEP::proton_mass_c2, "Mass of particle used in beam 1"}

private

Todo:: Get from a service

Definition at line 78 of file GenEventBeamEffectBooster.h.

◆ m_beam2ParticleMass

Gaudi::Property<double> Simulation::GenEventBeamEffectBooster::m_beam2ParticleMass {this, "Beam2ParticleMass", CLHEP::proton_mass_c2, "Mass of particle used in beam 2"}

private

Definition at line 79 of file GenEventBeamEffectBooster.h.

◆ m_dE

Gaudi::Property<double> Simulation::GenEventBeamEffectBooster::m_dE {this, "deltaEOverE", 1.1E-4, "Delta_E/E theoretical value"}

private

Definition at line 75 of file GenEventBeamEffectBooster.h.

◆ m_pbeam1

Gaudi::Property<double> Simulation::GenEventBeamEffectBooster::m_pbeam1 {this, "pbeam1", 3500.0E3, "Beam 1 Momentum / MeV"}

private

Definition at line 76 of file GenEventBeamEffectBooster.h.

◆ m_pbeam2

Gaudi::Property<double> Simulation::GenEventBeamEffectBooster::m_pbeam2 {this, "pbeam2", 3500.0E3, "Beam 2 Momentum / MeV"}

private

Todo:: Get from a service

Definition at line 77 of file GenEventBeamEffectBooster.h.

◆ m_randomEngineName

Gaudi::Property<std::string> Simulation::GenEventBeamEffectBooster::m_randomEngineName {this, "RandomStream", "BEAM"}

private

Name of the random number stream.

Definition at line 61 of file GenEventBeamEffectBooster.h.

◆ m_rndGenSvc

ServiceHandle<IAthRNGSvc> Simulation::GenEventBeamEffectBooster::m_rndGenSvc {this, "RandomSvc", "AthRNGSvc"}

private

Definition at line 59 of file GenEventBeamEffectBooster.h.

◆ m_sigma_px_b1

Gaudi::Property<double> Simulation::GenEventBeamEffectBooster::m_sigma_px_b1 {this, "Sigma_px_b1", 20.e-6, "angular divergence in x of beam 1 [rad]"}

private

Definition at line 66 of file GenEventBeamEffectBooster.h.

◆ m_sigma_px_b2

Gaudi::Property<double> Simulation::GenEventBeamEffectBooster::m_sigma_px_b2 {this, "Sigma_px_b2", 21.e-6, "angular divergence in x of beam 2 [rad]"}

private

Definition at line 67 of file GenEventBeamEffectBooster.h.

◆ m_sigma_py_b1

Gaudi::Property<double> Simulation::GenEventBeamEffectBooster::m_sigma_py_b1 {this, "Sigma_py_b1", 22.e-6, "angular divergence in y of beam 1 [rad]"}

private

Definition at line 68 of file GenEventBeamEffectBooster.h.

◆ m_sigma_py_b2

Gaudi::Property<double> Simulation::GenEventBeamEffectBooster::m_sigma_py_b2 {this, "Sigma_py_b2", 23.e-6, "angular divergence in y of beam 2 [rad]"}

private

Definition at line 69 of file GenEventBeamEffectBooster.h.

◆ m_xing_x_b1

Gaudi::Property<double> Simulation::GenEventBeamEffectBooster::m_xing_x_b1 {this, "HalfXing_x_b1", 19.e-6, "half-crossing in x for beam 1 at IP1, i.e. angle between beam 1 and z-axis"}

private

Definition at line 71 of file GenEventBeamEffectBooster.h.

◆ m_xing_x_b2

Gaudi::Property<double> Simulation::GenEventBeamEffectBooster::m_xing_x_b2 {this, "HalfXing_x_b2", 1.e-6, "half-crossing in x for beam 2 at IP1, i.e. angle between beam 2 and z-axis"}

private

Definition at line 72 of file GenEventBeamEffectBooster.h.

◆ m_xing_y_b1

Gaudi::Property<double> Simulation::GenEventBeamEffectBooster::m_xing_y_b1 {this, "HalfXing_y_b1", 150.e-6, "half-crossing in y for beam 1 at IP1, i.e. angle between beam 1 and z-axis"}

private

Definition at line 73 of file GenEventBeamEffectBooster.h.

◆ m_xing_y_b2

Gaudi::Property<double> Simulation::GenEventBeamEffectBooster::m_xing_y_b2 {this, "HalfXing_y_b2", -152.e-6, "half-crossing in y for beam 2 at IP1, i.e. angle between beam 2 and z-axis"}

private

Definition at line 74 of file GenEventBeamEffectBooster.h.

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

Public Member Functions

Private Member Functions

Private Attributes

Detailed Description

Constructor & Destructor Documentation

◆ GenEventBeamEffectBooster()

Member Function Documentation

◆ boostParticle()

◆ finalize()

◆ initialize()

◆ initializeAthenaEvent()

◆ initializeGenEvent()

◆ manipulate()

Member Data Documentation

◆ ATLAS_THREAD_SAFE

◆ m_applyBoost

◆ m_applyDivergence

◆ m_beam1ParticleMass

◆ m_beam2ParticleMass

◆ m_dE

◆ m_pbeam1

◆ m_pbeam2

◆ m_randomEngineName

◆ m_rndGenSvc

◆ m_sigma_px_b1

◆ m_sigma_px_b2

◆ m_sigma_py_b1

◆ m_sigma_py_b2

◆ m_xing_x_b1

◆ m_xing_x_b2

◆ m_xing_y_b1

◆ m_xing_y_b2