Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2025-01-28 23:24:34

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 
0009 using CLHEP::cm;
0010 using CLHEP::ns;
0011 
0012 GaussEvtVtxGenerator::GaussEvtVtxGenerator(const edm::ParameterSet& p) : BaseEvtVtxGenerator(p) {
0013   readDB_ = p.getParameter<bool>("readDB");
0014   if (!readDB_) {
0015     fMeanX = p.getParameter<double>("MeanX") * cm;
0016     fMeanY = p.getParameter<double>("MeanY") * cm;
0017     fMeanZ = p.getParameter<double>("MeanZ") * cm;
0018     fSigmaX = p.getParameter<double>("SigmaX") * cm;
0019     fSigmaY = p.getParameter<double>("SigmaY") * cm;
0020     fSigmaZ = p.getParameter<double>("SigmaZ") * cm;
0021     fTimeOffset = p.getParameter<double>("TimeOffset") * ns * c_light;  // HepMC distance units are in mm
0022 
0023     if (fSigmaX < 0) {
0024       throw cms::Exception("Configuration") << "Error in GaussEvtVtxGenerator: "
0025                                             << "Illegal resolution in X (SigmaX is negative)";
0026     }
0027     if (fSigmaY < 0) {
0028       throw cms::Exception("Configuration") << "Error in GaussEvtVtxGenerator: "
0029                                             << "Illegal resolution in Y (SigmaY is negative)";
0030     }
0031     if (fSigmaZ < 0) {
0032       throw cms::Exception("Configuration") << "Error in GaussEvtVtxGenerator: "
0033                                             << "Illegal resolution in Z (SigmaZ is negative)";
0034     }
0035   }
0036   if (readDB_) {
0037     // NOTE: this is currently watching LS transitions, while it should watch Run transitions,
0038     // even though in reality there is no Run Dependent MC (yet) in CMS
0039     beamToken_ = esConsumes<SimBeamSpotObjects, SimBeamSpotObjectsRcd, edm::Transition::BeginLuminosityBlock>();
0040   }
0041 }
0042 
0043 void GaussEvtVtxGenerator::beginLuminosityBlock(edm::LuminosityBlock const&, edm::EventSetup const& iEventSetup) {
0044   update(iEventSetup);
0045 }
0046 
0047 void GaussEvtVtxGenerator::update(const edm::EventSetup& iEventSetup) {
0048   if (readDB_ && parameterWatcher_.check(iEventSetup)) {
0049     edm::ESHandle<SimBeamSpotObjects> beamhandle = iEventSetup.getHandle(beamToken_);
0050     if (beamhandle->isGaussian()) {
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     } else {
0059       throw cms::Exception("Configuration")
0060           << "Error in GaussEvtVtxGenerator::update: The provided SimBeamSpotObjects is not Gaussian.\n"
0061           << "Please check the configuration and ensure that the beam spot parameters are appropriate for a Gaussian "
0062              "distribution.";
0063     }
0064   }
0065 }
0066 
0067 ROOT::Math::XYZTVector GaussEvtVtxGenerator::vertexShift(CLHEP::HepRandomEngine* engine) const {
0068   double X, Y, Z, T;
0069   X = CLHEP::RandGaussQ::shoot(engine, fMeanX, fSigmaX);
0070   Y = CLHEP::RandGaussQ::shoot(engine, fMeanY, fSigmaY);
0071   Z = CLHEP::RandGaussQ::shoot(engine, fMeanZ, fSigmaZ);
0072   T = CLHEP::RandGaussQ::shoot(engine, fTimeOffset, fSigmaZ);
0073 
0074   return ROOT::Math::XYZTVector(X, Y, Z, T);
0075 }
0076 
0077 void GaussEvtVtxGenerator::sigmaX(double s) {
0078   if (s >= 0) {
0079     fSigmaX = s;
0080   } else {
0081     throw cms::Exception("LogicError") << "Error in GaussEvtVtxGenerator::sigmaX: "
0082                                        << "Illegal resolution in X (negative)";
0083   }
0084 }
0085 
0086 void GaussEvtVtxGenerator::sigmaY(double s) {
0087   if (s >= 0) {
0088     fSigmaY = s;
0089   } else {
0090     throw cms::Exception("LogicError") << "Error in GaussEvtVtxGenerator::sigmaY: "
0091                                        << "Illegal resolution in Y (negative)";
0092   }
0093 }
0094 
0095 void GaussEvtVtxGenerator::sigmaZ(double s) {
0096   if (s >= 0) {
0097     fSigmaZ = s;
0098   } else {
0099     throw cms::Exception("LogicError") << "Error in GaussEvtVtxGenerator::sigmaZ: "
0100                                        << "Illegal resolution in Z (negative)";
0101   }
0102 }
0103 
0104 void GaussEvtVtxGenerator::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0105   edm::ParameterSetDescription desc;
0106   desc.add<double>("MeanX", 0.0)->setComment("in cm");
0107   desc.add<double>("MeanY", 0.0)->setComment("in cm");
0108   desc.add<double>("MeanZ", 0.0)->setComment("in cm");
0109   desc.add<double>("SigmaX", 0.0)->setComment("in cm");
0110   desc.add<double>("SigmaY", 0.0)->setComment("in cm");
0111   desc.add<double>("SigmaZ", 0.0)->setComment("in cm");
0112   desc.add<double>("TimeOffset", 0.0)->setComment("in ns");
0113   desc.add<edm::InputTag>("src");
0114   desc.add<bool>("readDB");
0115   descriptions.add("GaussEvtVtxGenerator", desc);
0116 }