EventVertexGenerators/src/GaussEvtVtxGenerator.cc

0001 #include "IOMC/EventVertexGenerators/interface/GaussEvtVtxGenerator.h"
0002 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0003 #include "FWCore/Utilities/interface/Exception.h"
0004
0005 #include <CLHEP/Random/RandGaussQ.h>
0006 #include <CLHEP/Units/SystemOfUnits.h>
0007 #include <CLHEP/Units/GlobalPhysicalConstants.h>
0008 #include "HepMC/SimpleVector.h"
0009
0010 using CLHEP::cm;
0011 using CLHEP::ns;
0012
0013 GaussEvtVtxGenerator::GaussEvtVtxGenerator(const edm::ParameterSet& p) : BaseEvtVtxGenerator(p) {
0014   readDB_ = p.getParameter<bool>("readDB");
0015   if (!readDB_) {
0016     fMeanX = p.getParameter<double>("MeanX") * cm;
0017     fMeanY = p.getParameter<double>("MeanY") * cm;
0018     fMeanZ = p.getParameter<double>("MeanZ") * cm;
0019     fSigmaX = p.getParameter<double>("SigmaX") * cm;
0020     fSigmaY = p.getParameter<double>("SigmaY") * cm;
0021     fSigmaZ = p.getParameter<double>("SigmaZ") * cm;
0022     fTimeOffset = p.getParameter<double>("TimeOffset") * ns * c_light;  // HepMC distance units are in mm
0023
0024     if (fSigmaX < 0) {
0025       throw cms::Exception("Configuration") << "Error in GaussEvtVtxGenerator: "
0026                                             << "Illegal resolution in X (SigmaX is negative)";
0027     }
0028     if (fSigmaY < 0) {
0029       throw cms::Exception("Configuration") << "Error in GaussEvtVtxGenerator: "
0030                                             << "Illegal resolution in Y (SigmaY is negative)";
0031     }
0032     if (fSigmaZ < 0) {
0033       throw cms::Exception("Configuration") << "Error in GaussEvtVtxGenerator: "
0034                                             << "Illegal resolution in Z (SigmaZ is negative)";
0035     }
0036   }
0037   if (readDB_) {
0038     // NOTE: this is currently watching LS transitions, while it should watch Run transitions,
0039     // even though in reality there is no Run Dependent MC (yet) in CMS
0040     beamToken_ = esConsumes<SimBeamSpotObjects, SimBeamSpotObjectsRcd, edm::Transition::BeginLuminosityBlock>();
0041   }
0042 }
0043
0044 void GaussEvtVtxGenerator::beginLuminosityBlock(edm::LuminosityBlock const&, edm::EventSetup const& iEventSetup) {
0045   update(iEventSetup);
0046 }
0047
0048 void GaussEvtVtxGenerator::update(const edm::EventSetup& iEventSetup) {
0049   if (readDB_ && parameterWatcher_.check(iEventSetup)) {
0050     edm::ESHandle<SimBeamSpotObjects> beamhandle = iEventSetup.getHandle(beamToken_);
0051     fMeanX = beamhandle->meanX() * cm;
0052     fMeanY = beamhandle->meanY() * cm;
0053     fMeanZ = beamhandle->meanZ() * cm;
0054     fSigmaX = beamhandle->sigmaX() * cm;
0055     fSigmaY = beamhandle->sigmaY() * cm;
0056     fSigmaZ = beamhandle->sigmaZ() * cm;
0057     fTimeOffset = beamhandle->timeOffset() * ns * c_light;  // HepMC distance units are in mm
0058   }
0059 }
0060
0061 HepMC::FourVector GaussEvtVtxGenerator::newVertex(CLHEP::HepRandomEngine* engine) const {
0062   double X, Y, Z, T;
0063   X = CLHEP::RandGaussQ::shoot(engine, fMeanX, fSigmaX);
0064   Y = CLHEP::RandGaussQ::shoot(engine, fMeanY, fSigmaY);
0065   Z = CLHEP::RandGaussQ::shoot(engine, fMeanZ, fSigmaZ);
0066   T = CLHEP::RandGaussQ::shoot(engine, fTimeOffset, fSigmaZ);
0067
0068   return HepMC::FourVector(X, Y, Z, T);
0069 }
0070
0071 void GaussEvtVtxGenerator::sigmaX(double s) {
0072   if (s >= 0) {
0073     fSigmaX = s;
0074   } else {
0075     throw cms::Exception("LogicError") << "Error in GaussEvtVtxGenerator::sigmaX: "
0076                                        << "Illegal resolution in X (negative)";
0077   }
0078 }
0079
0080 void GaussEvtVtxGenerator::sigmaY(double s) {
0081   if (s >= 0) {
0082     fSigmaY = s;
0083   } else {
0084     throw cms::Exception("LogicError") << "Error in GaussEvtVtxGenerator::sigmaY: "
0085                                        << "Illegal resolution in Y (negative)";
0086   }
0087 }
0088
0089 void GaussEvtVtxGenerator::sigmaZ(double s) {
0090   if (s >= 0) {
0091     fSigmaZ = s;
0092   } else {
0093     throw cms::Exception("LogicError") << "Error in GaussEvtVtxGenerator::sigmaZ: "
0094                                        << "Illegal resolution in Z (negative)";
0095   }
0096 }
0097
0098 void GaussEvtVtxGenerator::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0099   edm::ParameterSetDescription desc;
0100   desc.add<double>("MeanX", 0.0)->setComment("in cm");
0101   desc.add<double>("MeanY", 0.0)->setComment("in cm");
0102   desc.add<double>("MeanZ", 0.0)->setComment("in cm");
0103   desc.add<double>("SigmaX", 0.0)->setComment("in cm");
0104   desc.add<double>("SigmaY", 0.0)->setComment("in cm");
0105   desc.add<double>("SigmaZ", 0.0)->setComment("in cm");
0106   desc.add<double>("TimeOffset", 0.0)->setComment("in ns");
0107   desc.add<edm::InputTag>("src");
0108   desc.add<bool>("readDB");
0109   descriptions.add("GaussEvtVtxGenerator", desc);
0110 }