Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-04-20 23:57:15

0001 /*************************
0002  *
0003  *  Date: 2005/08 $
0004  *  \author J Weng - F. Moortgat'
0005  */
0006 
0007 #include <iostream>
0008 #include <algorithm>  // because we use std::swap
0009 
0010 #include <HepMC3/GenEvent.h>
0011 #include <HepMC3/GenVertex.h>
0012 #include <HepMC3/GenParticle.h>
0013 #include "SimDataFormats/GeneratorProducts/interface/HepMC3Product.h"
0014 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0015 
0016 using namespace edm;
0017 using namespace std;
0018 
0019 HepMC3Product::HepMC3Product(HepMC3::GenEvent* evt)
0020     : evt_(evt), isVtxGenApplied_(false), isVtxBoostApplied_(false), isPBoostApplied_(false) {}
0021 
0022 HepMC3Product::~HepMC3Product() {
0023   delete evt_;
0024   evt_ = nullptr;
0025   isVtxGenApplied_ = false;
0026   isVtxBoostApplied_ = false;
0027   isPBoostApplied_ = false;
0028 }
0029 
0030 void HepMC3Product::addHepMCData(HepMC3::GenEvent* evt) {
0031   //  evt->print();
0032   // cout <<sizeof (evt)  <<"  " << sizeof ( HepMC::GenEvent)   << endl;
0033   evt_ = evt;
0034 
0035   // same story about vertex smearing - GenEvent won't know it...
0036   // in fact, would be better to implement CmsGenEvent...
0037 }
0038 
0039 void HepMC3Product::applyVtxGen(HepMC3::FourVector const& vtxShift) {
0040   //std::cout<< " applyVtxGen called " << isVtxGenApplied_ << endl;
0041   //fTimeOffset = 0;
0042   if (isVtxGenApplied())
0043     return;
0044   evt_->shift_position_by(vtxShift);
0045   isVtxGenApplied_ = true;
0046   return;
0047 }
0048 
0049 void HepMC3Product::boostToLab(TMatrixD const* lorentz, std::string const& type) {
0050   //std::cout << "from boostToLab:" << std::endl;
0051 
0052   if (lorentz == nullptr) {
0053     //std::cout << " lorentz = 0 " << std::endl;
0054     return;
0055   }
0056 
0057   //lorentz->Print();
0058 
0059   TMatrixD tmplorentz(*lorentz);
0060   //tmplorentz.Print();
0061 
0062   if (type == "vertex") {
0063     if (isVtxBoostApplied()) {
0064       //std::cout << " isVtxBoostApplied true " << std::endl;
0065       return;
0066     }
0067 
0068     for (const HepMC3::GenVertexPtr& vt : evt_->vertices()) {
0069       // change basis to lorentz boost definition: (t,x,z,y)
0070       TMatrixD p4(4, 1);
0071       p4(0, 0) = vt->position().t();
0072       p4(1, 0) = vt->position().x();
0073       p4(2, 0) = vt->position().z();
0074       p4(3, 0) = vt->position().y();
0075 
0076       TMatrixD p4lab(4, 1);
0077       p4lab = tmplorentz * p4;
0078       //std::cout << " vertex lorentz: " << p4lab(1,0) << " " << p4lab(3,0) << " " << p4lab(2,0) << std::endl;
0079       vt->set_position(HepMC3::FourVector(p4lab(1, 0), p4lab(3, 0), p4lab(2, 0), p4lab(0, 0)));
0080     }
0081 
0082     isVtxBoostApplied_ = true;
0083   } else if (type == "momentum") {
0084     if (isPBoostApplied()) {
0085       //std::cout << " isPBoostApplied true " << std::endl;
0086       return;
0087     }
0088 
0089     for (const HepMC3::GenParticlePtr& part : evt_->particles()) {
0090       // change basis to lorentz boost definition: (E,Px,Pz,Py)
0091       TMatrixD p4(4, 1);
0092       p4(0, 0) = part->momentum().e();
0093       p4(1, 0) = part->momentum().x();
0094       p4(2, 0) = part->momentum().z();
0095       p4(3, 0) = part->momentum().y();
0096 
0097       TMatrixD p4lab(4, 1);
0098       p4lab = tmplorentz * p4;
0099       //std::cout << " momentum lorentz: " << p4lab(1,0) << " " << p4lab(3,0) << " " << p4lab(2,0) << std::endl;
0100       part->set_momentum(HepMC3::FourVector(p4lab(1, 0), p4lab(3, 0), p4lab(2, 0), p4lab(0, 0)));
0101     }
0102 
0103     isPBoostApplied_ = true;
0104   } else {
0105     edm::LogWarning("GeneratorProducts") << "no type found for boostToLab(std::string), options are vertex or momentum";
0106     //throw edm::Exception(edm::errors::Configuration, "GeneratorProducts")
0107     //  << "no type found for boostToLab(std::string), options are vertex or momentum \n";
0108   }
0109 
0110   return;
0111 }
0112 
0113 const HepMC3::GenEvent& HepMC3Product::getHepMCData() const { return *evt_; }
0114 
0115 // copy constructor
0116 HepMC3Product::HepMC3Product(HepMC3Product const& other) : evt_(nullptr) {
0117   if (other.evt_)
0118     evt_ = new HepMC3::GenEvent(*other.evt_);
0119   isVtxGenApplied_ = other.isVtxGenApplied_;
0120   isVtxBoostApplied_ = other.isVtxBoostApplied_;
0121   isPBoostApplied_ = other.isPBoostApplied_;
0122   //fTimeOffset = other.fTimeOffset;
0123 }
0124 
0125 // swap
0126 void HepMC3Product::swap(HepMC3Product& other) {
0127   std::swap(evt_, other.evt_);
0128   std::swap(isVtxGenApplied_, other.isVtxGenApplied_);
0129   std::swap(isVtxBoostApplied_, other.isVtxBoostApplied_);
0130   std::swap(isPBoostApplied_, other.isPBoostApplied_);
0131   //std::swap(fTimeOffset, other.fTimeOffset);
0132 }
0133 
0134 // assignment: use copy/swap idiom for exception safety.
0135 HepMC3Product& HepMC3Product::operator=(HepMC3Product const& other) {
0136   HepMC3Product temp(other);
0137   swap(temp);
0138   return *this;
0139 }
0140 
0141 // move, needed explicitly as we have raw pointer...
0142 HepMC3Product::HepMC3Product(HepMC3Product&& other) : evt_(nullptr) { swap(other); }
0143 HepMC3Product& HepMC3Product::operator=(HepMC3Product&& other) {
0144   swap(other);
0145   return *this;
0146 }