Notification/src/TmpSimEvent.cc

0001 #include "SimG4Core/Notification/interface/TmpSimEvent.h"
0002 #include "SimDataFormats/EncodedEventId/interface/EncodedEventId.h"
0003
0004 #include <CLHEP/Units/SystemOfUnits.h>
0005
0006 class IdSort {
0007 public:
0008   bool operator()(const SimTrack& a, const SimTrack& b) { return a.trackId() < b.trackId(); }
0009 };
0010
0011 TmpSimEvent::TmpSimEvent() {
0012   g4vertices_.reserve(2000);
0013   g4tracks_.reserve(4000);
0014 }
0015
0016 TmpSimEvent::~TmpSimEvent() { clear(); }
0017
0018 void TmpSimEvent::clear() {
0019   for (auto& ptr : g4tracks_) {
0020     delete ptr;
0021   }
0022   g4tracks_.clear();
0023   for (auto& ptr : g4vertices_) {
0024     delete ptr;
0025   }
0026   g4vertices_.clear();
0027 }
0028
0029 void TmpSimEvent::load(edm::SimTrackContainer& c) const {
0030   const double invgev = 1.0 / CLHEP::GeV;
0031   for (auto& trk : g4tracks_) {
0032     int ip = trk->part();
0033     const math::XYZVectorD& mom = trk->momentum();
0034     math::XYZTLorentzVectorD p(mom.x() * invgev, mom.y() * invgev, mom.z() * invgev, trk->energy() * invgev);
0035     int iv = trk->ivert();
0036     int ig = trk->igenpart();
0037     int id = trk->id();
0038     // ip = particle ID as PDG
0039     // pp = 4-momentum in GeV
0040     // iv = corresponding TmpSimVertex index
0041     // ig = corresponding GenParticle index
0042     SimTrack t = SimTrack(ip, p, iv, ig, trk->trackerSurfacePosition(), trk->trackerSurfaceMomentum());
0043     t.setTrackId(id);
0044     t.setEventId(EncodedEventId(0));
0045     t.setCrossedBoundaryVars(
0046         trk->crossedBoundary(), trk->getIDAtBoundary(), trk->getPositionAtBoundary(), trk->getMomentumAtBoundary());
0047     c.push_back(t);
0048   }
0049   std::stable_sort(c.begin(), c.end(), IdSort());
0050 }
0051
0052 void TmpSimEvent::load(edm::SimVertexContainer& c) const {
0053   const double invcm = 1.0 / CLHEP::cm;
0054   // index of the vertex is needed to make SimVertex object
0055   for (unsigned int i = 0; i < g4vertices_.size(); ++i) {
0056     TmpSimVertex* vtx = g4vertices_[i];
0057     auto pos = vtx->vertexPosition();
0058     math::XYZVectorD v3(pos.x() * invcm, pos.y() * invcm, pos.z() * invcm);
0059     float t = (float)(vtx->vertexGlobalTime() / CLHEP::second);
0060     int iv = vtx->parentIndex();
0061     // v3 = position in cm
0062     // t  = global time in second
0063     // iv = index of the parent in the SimEvent SimTrack container (-1 if no parent)
0064     SimVertex v = SimVertex(v3, t, iv, i);
0065     v.setProcessType((unsigned int)vtx->processType());
0066     v.setEventId(EncodedEventId(0));
0067     c.push_back(v);
0068   }
0069 }