File indexing completed on 2025-01-11 03:38:36
0001
0002
0003
0004
0005
0006
0007 #include <iostream>
0008 #include <algorithm> // because we use std::swap
0009
0010 #include <HepMC/GenEvent.h>
0011 #include <HepMC/SimpleVector.h>
0012 #include "SimDataFormats/GeneratorProducts/interface/HepMCProduct.h"
0013
0014 using namespace edm;
0015 using namespace std;
0016
0017 HepMCProduct::HepMCProduct(HepMC::GenEvent* evt)
0018 : evt_(evt), isVtxGenApplied_(false), isVtxBoostApplied_(false), isPBoostApplied_(false) {}
0019
0020 HepMCProduct::~HepMCProduct() {
0021 delete evt_;
0022 evt_ = nullptr;
0023 isVtxGenApplied_ = false;
0024 isVtxBoostApplied_ = false;
0025 isPBoostApplied_ = false;
0026 }
0027
0028 void HepMCProduct::addHepMCData(HepMC::GenEvent* evt) {
0029
0030
0031 evt_ = evt;
0032
0033
0034
0035 }
0036
0037 void HepMCProduct::applyVtxGen(HepMC::FourVector const& vtxShift) {
0038
0039
0040 if (isVtxGenApplied())
0041 return;
0042
0043 for (HepMC::GenEvent::vertex_iterator vt = evt_->vertices_begin(); vt != evt_->vertices_end(); ++vt) {
0044 double x = (*vt)->position().x() + vtxShift.x();
0045 double y = (*vt)->position().y() + vtxShift.y();
0046 double z = (*vt)->position().z() + vtxShift.z();
0047 double t = (*vt)->position().t() + vtxShift.t();
0048
0049 (*vt)->set_position(HepMC::FourVector(x, y, z, t));
0050 }
0051
0052 isVtxGenApplied_ = true;
0053
0054 return;
0055 }
0056
0057 void HepMCProduct::boostToLab(TMatrixD const* lorentz, std::string const& type) {
0058
0059
0060 if (lorentz == nullptr) {
0061
0062 return;
0063 }
0064
0065
0066
0067 TMatrixD tmplorentz(*lorentz);
0068
0069
0070 if (type == "vertex") {
0071 if (isVtxBoostApplied()) {
0072
0073 return;
0074 }
0075
0076 for (HepMC::GenEvent::vertex_iterator vt = evt_->vertices_begin(); vt != evt_->vertices_end(); ++vt) {
0077
0078 TMatrixD p4(4, 1);
0079 p4(0, 0) = (*vt)->position().t();
0080 p4(1, 0) = (*vt)->position().x();
0081 p4(2, 0) = (*vt)->position().z();
0082 p4(3, 0) = (*vt)->position().y();
0083
0084 TMatrixD p4lab(4, 1);
0085 p4lab = tmplorentz * p4;
0086
0087 (*vt)->set_position(HepMC::FourVector(p4lab(1, 0), p4lab(3, 0), p4lab(2, 0), p4lab(0, 0)));
0088 }
0089
0090 isVtxBoostApplied_ = true;
0091 } else if (type == "momentum") {
0092 if (isPBoostApplied()) {
0093
0094 return;
0095 }
0096
0097 for (HepMC::GenEvent::particle_iterator part = evt_->particles_begin(); part != evt_->particles_end(); ++part) {
0098
0099 TMatrixD p4(4, 1);
0100 p4(0, 0) = (*part)->momentum().e();
0101 p4(1, 0) = (*part)->momentum().x();
0102 p4(2, 0) = (*part)->momentum().z();
0103 p4(3, 0) = (*part)->momentum().y();
0104
0105 TMatrixD p4lab(4, 1);
0106 p4lab = tmplorentz * p4;
0107
0108 (*part)->set_momentum(HepMC::FourVector(p4lab(1, 0), p4lab(3, 0), p4lab(2, 0), p4lab(0, 0)));
0109 }
0110
0111 isPBoostApplied_ = true;
0112 } else {
0113 std::cout << " no type found for boostToLab(std::string), options are vertex or momentum" << std::endl;
0114 }
0115
0116 return;
0117 }
0118
0119 const HepMC::GenEvent& HepMCProduct::getHepMCData() const { return *evt_; }
0120
0121
0122 HepMCProduct::HepMCProduct(HepMCProduct const& other) : evt_(nullptr) {
0123 if (other.evt_)
0124 evt_ = new HepMC::GenEvent(*other.evt_);
0125 isVtxGenApplied_ = other.isVtxGenApplied_;
0126 isVtxBoostApplied_ = other.isVtxBoostApplied_;
0127 isPBoostApplied_ = other.isPBoostApplied_;
0128 }
0129
0130
0131 void HepMCProduct::swap(HepMCProduct& other) {
0132 std::swap(evt_, other.evt_);
0133 std::swap(isVtxGenApplied_, other.isVtxGenApplied_);
0134 std::swap(isVtxBoostApplied_, other.isVtxBoostApplied_);
0135 std::swap(isPBoostApplied_, other.isPBoostApplied_);
0136 }
0137
0138
0139 HepMCProduct& HepMCProduct::operator=(HepMCProduct const& other) {
0140 HepMCProduct temp(other);
0141 swap(temp);
0142 return *this;
0143 }
0144
0145
0146 HepMCProduct::HepMCProduct(HepMCProduct&& other)
0147 : evt_(nullptr), isVtxGenApplied_(false), isVtxBoostApplied_(false), isPBoostApplied_(false) {
0148 swap(other);
0149 }
0150 HepMCProduct& HepMCProduct::operator=(HepMCProduct&& other) {
0151 swap(other);
0152 return *this;
0153 }