Adapt SimBeamSpotObjects to handle both Betafunc and Gauss EvtVtxGene…

…rators

CondCore/BeamSpotPlugins/interface/BeamSpotPayloadInspectorHelper.h

@@ -715,17 +715,22 @@ namespace beamSpotPI {
 namespace simBeamSpotPI {
 
   enum parameters {
-    X = 0,              // 0  - Positions
-    Y = 1,              // 1
-    Z = 2,              // 2
-    sigmaZ = 3,         // 3  - Widths
-    betaStar = 4,       // 4
-    emittance = 5,      // 5
-    expTransWidth = 6,  // 6  - from LPC-like calculation
-    phi = 7,            // 7  - Additional parameters
-    alpha = 8,          // 8
-    timeOffset = 9,     // 9
-    END_OF_TYPES = 10,
+    X = 0,               // 0  - Positions
+    Y = 1,               // 1
+    Z = 2,               // 2
+    meanX = 3,           // 3
+    meanY = 4,           // 4
+    meanZ = 5,           // 5
+    sigmaX = 6,          // 6  - Widths
+    sigmaY = 7,          // 7
+    sigmaZ = 8,          // 8
+    betaStar = 9,        // 9
+    emittance = 10,      // 10
+    expTransWidth = 11,  // 11  - from LPC-like calculation
+    phi = 12,            // 12  - Additional parameters
+    alpha = 13,          // 13
+    timeOffset = 14,     // 14
+    END_OF_TYPES = 15,
   };
 
   /************************************************/
@@ -737,6 +742,16 @@ namespace simBeamSpotPI {
         return (addUnits ? "Y [cm]" : "Y");
       case Z:
         return (addUnits ? "Z [cm]" : "Z");
+      case meanX:
+        return (addUnits ? "MeanX [cm]" : "meanX");
+      case meanY:
+        return (addUnits ? "MeanY [cm]" : "meanY");
+      case meanZ:
+        return (addUnits ? "MeanZ [cm]" : "meanZ");
+      case sigmaX:
+        return (addUnits ? "#sigma_{X} [cm]" : "sigmaX");
+      case sigmaY:
+        return (addUnits ? "#sigma_{Y} [cm]" : "sigmaY");
       case sigmaZ:
         return (addUnits ? "#sigma_{Z} [cm]" : "sigmaZ");
       case betaStar:
@@ -769,8 +784,11 @@ namespace simBeamSpotPI {
     SimBSParamsHelper(const std::shared_ptr<PayloadType>& bs) {
       // fill in the values
       m_values[parameters::X] = bs->x(), m_values[parameters::Y] = bs->y(), m_values[parameters::Z] = bs->z();
-      m_values[parameters::sigmaZ] = bs->sigmaZ(), m_values[parameters::betaStar] = bs->betaStar(),
-      m_values[parameters::emittance] = bs->emittance();
+      m_values[parameters::meanX] = bs->meanX(), m_values[parameters::meanY] = bs->meanY(),
+      m_values[parameters::meanZ] = bs->meanZ();
+      m_values[parameters::sigmaX] = bs->sigmaX(), m_values[parameters::sigmaY] = bs->sigmaY(),
+      m_values[parameters::sigmaZ] = bs->sigmaZ();
+      m_values[parameters::betaStar] = bs->betaStar(), m_values[parameters::emittance] = bs->emittance();
       m_values[parameters::expTransWidth] = (1 / std::sqrt(2)) * std::sqrt(bs->emittance() * bs->betaStar()) * 10000.f;
       m_values[parameters::phi] = bs->phi(), m_values[parameters::alpha] = bs->alpha(),
       m_values[parameters::timeOffset] = bs->timeOffset();
@@ -842,6 +860,16 @@ namespace simBeamSpotPI {
             return m_payload->y();
           case Z:
             return m_payload->z();
+          case meanX:
+            return m_payload->meanX();
+          case meanY:
+            return m_payload->meanY();
+          case meanZ:
+            return m_payload->meanZ();
+          case sigmaX:
+            return m_payload->sigmaX();
+          case sigmaY:
+            return m_payload->sigmaY();
           case sigmaZ:
             return m_payload->sigmaZ();
           case betaStar:

CondCore/BeamSpotPlugins/test/testBeamSpotPayloadInspector.cpp

@@ -84,7 +84,7 @@ int main(int argc, char** argv) {
   // SimBeamSpot
   std::string prepConnectionString("frontier://FrontierPrep/CMS_CONDITIONS");
 
-  tag = "SimBeamSpotObjects_Realistic25ns13p6TeVEOY2022Collision_v0_mc";
+  tag = "SimBeamSpot_Realistic25ns13p6TeVEOY2022CollisionVtxSmearingParameters_v1_mc";
   start = static_cast<unsigned long long>(1);
   end = static_cast<unsigned long long>(1);
 
@@ -98,8 +98,8 @@ int main(int argc, char** argv) {
   histoSimParametersDiff.process(prepConnectionString, PI::mk_input(tag, start, end));
   edm::LogPrint("testBeamSpotPayloadInspector") << histoSimParametersDiff.data() << std::endl;
 
-  tag1 = "SimBeamSpotObjects_Realistic25ns13p6TeVEOY2022Collision_v0_mc";
-  tag2 = "SimBeamSpotObjects_Realistic25ns13p6TeVEarly2023Collision_v0_mc";
+  tag1 = "SimBeamSpot_Realistic25ns13p6TeVEOY2022CollisionVtxSmearingParameters_v1_mc";
+  tag2 = "SimBeamSpot_Realistic25ns13p6TeVEarly2023CollisionVtxSmearingParameters_v1_mc";
   start = static_cast<unsigned long long>(1);
 
   SimBeamSpotParametersDiffTwoTags histoSimParametersDiffTwoTags;

CondFormats/BeamSpotObjects/interface/SimBeamSpotObjects.h

@@ -3,7 +3,21 @@
 
 /** \class SimBeamSpotObjects
  *
- * provide the vertex smearing parameters from DB
+ * Provide the vertex smearing parameters from DB
+ *
+ * This Object contains the parameters needed by the vtx smearing functions used up to Run 3:
+ *   - BetafuncEvtVtxGenerator (realistic Run1/Run2/Run3 conditions)
+ *     Parameters used:
+ *        - fX0, fY0, fZ0
+ *        - fSigmaZ
+ *        - fbetastar, femittance
+ *        - fPhi, fAlpha
+ *        - fTimeOffset
+ *   - GaussEvtVtxGenerator (design Run1/Run2/Run3 conditions)
+ *     Parameters used:
+ *        - fMeanX, fMeanY, fMeanZ
+ *        - fSigmaX, fSigmaY, fSigmaZ
+ *        - fTimeOffset
  *
  */
 
@@ -18,6 +32,11 @@ class SimBeamSpotObjects {
     fX0 = 0.0;
     fY0 = 0.0;
     fZ0 = 0.0;
+    fMeanX = 0.0;
+    fMeanY = 0.0;
+    fMeanZ = 0.0;
+    fSigmaX = 0.0;
+    fSigmaY = 0.0;
     fSigmaZ = 0.0;
     fbetastar = 0.0;
     femittance = 0.0;
@@ -32,7 +51,13 @@ class SimBeamSpotObjects {
   void setX(double val) { fX0 = val; }
   void setY(double val) { fY0 = val; }
   void setZ(double val) { fZ0 = val; }
-  /// set sigmaZ
+  /// set meanX, meanY, meanZ
+  void setMeanX(double val) { fMeanX = val; }
+  void setMeanY(double val) { fMeanY = val; }
+  void setMeanZ(double val) { fMeanZ = val; }
+  /// set sigmaX, sigmaY, sigmaZ
+  void setSigmaX(double val) { fSigmaX = val; }
+  void setSigmaY(double val) { fSigmaY = val; }
   void setSigmaZ(double val) { fSigmaZ = val; }
   /// set BetaStar and Emittance
   void setBetaStar(double val) { fbetastar = val; }
@@ -42,31 +67,33 @@ class SimBeamSpotObjects {
   void setAlpha(double val) { fAlpha = val; }
   void setTimeOffset(double val) { fTimeOffset = val; }
 
-  /// get X position
+  /// get X, Y, Z position
   double x() const { return fX0; }
-  /// get Y position
   double y() const { return fY0; }
-  /// get Z position
   double z() const { return fZ0; }
-  /// get sigmaZ
+  /// get meanX, meanY, meanZ position
+  double meanX() const { return fMeanX; }
+  double meanY() const { return fMeanY; }
+  double meanZ() const { return fMeanZ; }
+  /// get sigmaX, sigmaY, sigmaZ
+  double sigmaX() const { return fSigmaX; }
+  double sigmaY() const { return fSigmaY; }
   double sigmaZ() const { return fSigmaZ; }
-  /// get BetaStar
+  /// get BetaStar and Emittance
   double betaStar() const { return fbetastar; }
-  /// get Emittance
   double emittance() const { return femittance; }
-  /// get Phi
+  /// get Phi, Alpha and TimeOffset
   double phi() const { return fPhi; }
-  /// get Alpha
   double alpha() const { return fAlpha; }
-  /// get TimeOffset
   double timeOffset() const { return fTimeOffset; }
 
   /// print sim beam spot parameters
   void print(std::stringstream& ss) const;
 
 private:
   double fX0, fY0, fZ0;
-  double fSigmaZ;
+  double fMeanX, fMeanY, fMeanZ;
+  double fSigmaX, fSigmaY, fSigmaZ;
   double fbetastar, femittance;
   double fPhi, fAlpha;
   double fTimeOffset;

CondFormats/BeamSpotObjects/src/SimBeamSpotObjects.cc

@@ -8,6 +8,11 @@ void SimBeamSpotObjects::print(std::stringstream& ss) const {
      << "       X0     = " << x() << " [cm]\n"
      << "       Y0     = " << y() << " [cm]\n"
      << "       Z0     = " << z() << " [cm]\n"
+     << "    MeanX     = " << meanX() << " [cm]\n"
+     << "    MeanY     = " << meanY() << " [cm]\n"
+     << "    MeanZ     = " << meanZ() << " [cm]\n"
+     << " Sigma X0     = " << sigmaX() << " [cm]\n"
+     << " Sigma Y0     = " << sigmaY() << " [cm]\n"
      << " Sigma Z0     = " << sigmaZ() << " [cm]\n"
      << " Beta star    = " << betaStar() << " [cm]\n"
      << " Emittance X  = " << emittance() << " [cm]\n"

CondTools/BeamSpot/plugins/BeamProfile2DBReader.cc

@@ -56,7 +56,8 @@ class BeamProfile2DBReader : public edm::one::EDAnalyzer<edm::one::SharedResourc
     int run;
     int ls;
     double fX0, fY0, fZ0;
-    double fSigmaZ;
+    double fMeanX, fMeanY, fMeanZ;
+    double fSigmaX, fSigmaY, fSigmaZ;
     double fbetastar, femittance;
     double fPhi, fAlpha;
     double fTimeOffset;
@@ -102,6 +103,11 @@ void BeamProfile2DBReader::TheBSfromDB::init() {
   fX0 = dummy_double;
   fY0 = dummy_double;
   fZ0 = dummy_double;
+  fMeanX = dummy_double;
+  fMeanY = dummy_double;
+  fMeanZ = dummy_double;
+  fSigmaX = dummy_double;
+  fSigmaY = dummy_double;
   fSigmaZ = dummy_double;
   fbetastar = dummy_double;
   femittance = dummy_double;
@@ -130,6 +136,11 @@ void BeamProfile2DBReader::analyze(const edm::Event& iEvent, const edm::EventSet
     theBSfromDB_.fX0 = mybeamspot->x();
     theBSfromDB_.fY0 = mybeamspot->y();
     theBSfromDB_.fZ0 = mybeamspot->z();
+    theBSfromDB_.fMeanX = mybeamspot->meanX();
+    theBSfromDB_.fMeanY = mybeamspot->meanY();
+    theBSfromDB_.fMeanZ = mybeamspot->meanZ();
+    theBSfromDB_.fSigmaX = mybeamspot->sigmaX();
+    theBSfromDB_.fSigmaY = mybeamspot->sigmaY();
     theBSfromDB_.fSigmaZ = mybeamspot->sigmaZ();
     theBSfromDB_.fbetastar = mybeamspot->betaStar();
     theBSfromDB_.femittance = mybeamspot->emittance();
@@ -157,6 +168,11 @@ void BeamProfile2DBReader::beginJob() {
   bstree_->Branch("BSx0", &theBSfromDB_.fX0, "BSx0/F");
   bstree_->Branch("BSy0", &theBSfromDB_.fY0, "BSy0/F");
   bstree_->Branch("BSz0", &theBSfromDB_.fZ0, "BSz0/F");
+  bstree_->Branch("BSmeanX", &theBSfromDB_.fMeanX, "BSmeanX/F");
+  bstree_->Branch("BSmeanY", &theBSfromDB_.fMeanY, "BSmeanY/F");
+  bstree_->Branch("BSmeanZ", &theBSfromDB_.fMeanZ, "BSmeanZ/F");
+  bstree_->Branch("Beamsigmax", &theBSfromDB_.fSigmaX, "Beamsigmax/F");
+  bstree_->Branch("Beamsigmay", &theBSfromDB_.fSigmaY, "Beamsigmay/F");
   bstree_->Branch("Beamsigmaz", &theBSfromDB_.fSigmaZ, "Beamsigmaz/F");
   bstree_->Branch("BetaStar", &theBSfromDB_.fbetastar, "BetaStar/F");
   bstree_->Branch("Emittance", &theBSfromDB_.femittance, "Emittance/F");

CondTools/BeamSpot/plugins/BeamProfile2DBWriter.cc

@@ -60,6 +60,11 @@ BeamProfile2DBWriter::BeamProfile2DBWriter(const edm::ParameterSet& iConfig)
   beamSpot_.setX(iConfig.getParameter<double>("X0"));
   beamSpot_.setY(iConfig.getParameter<double>("Y0"));
   beamSpot_.setZ(iConfig.getParameter<double>("Z0"));
+  beamSpot_.setMeanX(iConfig.getParameter<double>("MeanX"));
+  beamSpot_.setMeanY(iConfig.getParameter<double>("MeanY"));
+  beamSpot_.setMeanZ(iConfig.getParameter<double>("MeanZ"));
+  beamSpot_.setSigmaX(iConfig.getParameter<double>("SigmaX"));
+  beamSpot_.setSigmaY(iConfig.getParameter<double>("SigmaY"));
   beamSpot_.setSigmaZ(iConfig.getParameter<double>("SigmaZ"));
   beamSpot_.setAlpha(iConfig.getParameter<double>("Alpha"));
   beamSpot_.setPhi(iConfig.getParameter<double>("Phi"));
@@ -88,12 +93,17 @@ void BeamProfile2DBWriter::fillDescriptions(edm::ConfigurationDescriptions& desc
   edm::ParameterSetDescription desc;
   desc.add<std::string>("recordName", "SimBeamSpotObjectsRcd")
       ->setComment("name of the record to use for the PoolDBOutputService");
-  desc.add<double>("X0")->setComment("in cm");
-  desc.add<double>("Y0")->setComment("in cm");
-  desc.add<double>("Z0")->setComment("in cm");
-  desc.add<double>("SigmaZ")->setComment("in cm");
-  desc.add<double>("BetaStar")->setComment("in cm");
-  desc.add<double>("Emittance")->setComment("in cm");
+  desc.add<double>("X0", 0.0)->setComment("in cm");
+  desc.add<double>("Y0", 0.0)->setComment("in cm");
+  desc.add<double>("Z0", 0.0)->setComment("in cm");
+  desc.add<double>("MeanX", 0.0)->setComment("in cm");
+  desc.add<double>("MeanY", 0.0)->setComment("in cm");
+  desc.add<double>("MeanZ", 0.0)->setComment("in cm");
+  desc.add<double>("SigmaX", 0.0)->setComment("in cm");
+  desc.add<double>("SigmaY", 0.0)->setComment("in cm");
+  desc.add<double>("SigmaZ", 0.0)->setComment("in cm");
+  desc.add<double>("BetaStar", 0.0)->setComment("in cm");
+  desc.add<double>("Emittance", 0.0)->setComment("in cm");
   desc.add<double>("Alpha", 0.0)->setComment("in radians");
   desc.add<double>("Phi", 0.0)->setComment("in radians");
   desc.add<double>("TimeOffset", 0.0)->setComment("in ns");

CondTools/BeamSpot/test/BeamProfile2DBWriterAll_cfg.py

@@ -31,20 +31,20 @@
 from CondTools.BeamSpot.beamProfile2DBWriter_cfi import beamProfile2DBWriter
 process.load("IOMC.EventVertexGenerators.VtxSmearedParameters_cfi")
 
-# Get the PSets that have BetaStar as a parameter
-psets_with_BetaStar = []
+# Get the PSets that have BetaStar (BetafuncEvtVtxGenerator) or SigmaX (GaussEvtVtxGenerator) as a parameter
+psets_SimBeamSpotObjects = []
 psets_names = []
 for psetName in process.__dict__:
     pset = getattr(process, psetName)
-    if isinstance(pset, cms.PSet) and "BetaStar" in pset.parameterNames_():
+    if isinstance(pset, cms.PSet) and ("BetaStar" in pset.parameterNames_() or "SigmaX" in pset.parameterNames_()):
         print(psetName)
         psets_names.append(psetName)
-        psets_with_BetaStar.append(pset)
+        psets_SimBeamSpotObjects.append(pset)
 
 # Create a VPSet to store the parameter sets
 myVPSet = cms.VPSet()
 
-for i, pset in enumerate(psets_with_BetaStar):
+for i, pset in enumerate(psets_SimBeamSpotObjects):
     cloneName = 'BeamProfile2DBWriter_' + str(i+1)  # Unique clone name
     setattr(process, cloneName, beamProfile2DBWriter.clone(pset,
                                                            recordName = cms.string(psets_names[i])))
@@ -74,7 +74,7 @@
 process.end = cms.EndPath()
 
 process.schedule = cms.Schedule()
-for i, pset in enumerate(psets_with_BetaStar):
+for i, pset in enumerate(psets_SimBeamSpotObjects):
     pathName = 'Path_' + str(i+1)  # Unique path name
     process.schedule.append(getattr(process, pathName))
 process.schedule.append(process.end)

CondTools/BeamSpot/test/extractAndUploadAll.py

@@ -31,7 +31,7 @@
     data = {
         "destinationDatabase": "oracle://cms_orcoff_prep/CMS_CONDITIONS",
         "destinationTags": {
-            "SimBeamSpot_" + value + "_v0_mc": {}
+            "SimBeamSpot_" + value + "_v1_mc": {}
         },
         "inputTag": value,
         "since": None,

IOMC/EventVertexGenerators/interface/BetafuncEvtVtxGenerator.h

@@ -35,7 +35,7 @@ class BetafuncEvtVtxGenerator : public BaseEvtVtxGenerator {
   BetafuncEvtVtxGenerator(const BetafuncEvtVtxGenerator& p) = delete;
   /** Copy assignment operator */
   BetafuncEvtVtxGenerator& operator=(const BetafuncEvtVtxGenerator& rhs) = delete;
-  ~BetafuncEvtVtxGenerator() override;
+  ~BetafuncEvtVtxGenerator() override = default;
 
   void beginLuminosityBlock(edm::LuminosityBlock const&, edm::EventSetup const&) override;
 

IOMC/EventVertexGenerators/interface/GaussEvtVtxGenerator.h

@@ -8,6 +8,10 @@
  */
 
 #include "IOMC/EventVertexGenerators/interface/BaseEvtVtxGenerator.h"
+#include "FWCore/Framework/interface/ESWatcher.h"
+#include "FWCore/Utilities/interface/ESGetToken.h"
+#include "CondFormats/DataRecord/interface/SimBeamSpotObjectsRcd.h"
+#include "CondFormats/BeamSpotObjects/interface/SimBeamSpotObjects.h"
 
 namespace CLHEP {
   class HepRandomEngine;
@@ -20,7 +24,9 @@ class GaussEvtVtxGenerator : public BaseEvtVtxGenerator {
   GaussEvtVtxGenerator(const GaussEvtVtxGenerator& p) = delete;
   /** Copy assignment operator */
   GaussEvtVtxGenerator& operator=(const GaussEvtVtxGenerator& rhs) = delete;
-  ~GaussEvtVtxGenerator() override;
+  ~GaussEvtVtxGenerator() override = default;
+
+  void beginLuminosityBlock(edm::LuminosityBlock const&, edm::EventSetup const&) override;
 
   /// return a new event vertex
   //virtual CLHEP::Hep3Vector* newVertex();
@@ -43,9 +49,15 @@ class GaussEvtVtxGenerator : public BaseEvtVtxGenerator {
   void meanZ(double m = 0) { fMeanZ = m; }
 
 private:
+  bool readDB_;
+
   double fSigmaX, fSigmaY, fSigmaZ;
   double fMeanX, fMeanY, fMeanZ;
   double fTimeOffset;
+
+  void update(const edm::EventSetup& iEventSetup);
+  edm::ESWatcher<SimBeamSpotObjectsRcd> parameterWatcher_;
+  edm::ESGetToken<SimBeamSpotObjects, SimBeamSpotObjectsRcd> beamToken_;
 };
 
 #endif

IOMC/EventVertexGenerators/src/BetafuncEvtVtxGenerator.cc

@@ -38,7 +38,7 @@ BetafuncEvtVtxGenerator::BetafuncEvtVtxGenerator(const edm::ParameterSet& p) : B
     fSigmaZ = p.getParameter<double>("SigmaZ") * cm;
     fbetastar = p.getParameter<double>("BetaStar") * cm;
     femittance = p.getParameter<double>("Emittance") * cm;              // this is not the normalized emittance
-    fTimeOffset = p.getParameter<double>("TimeOffset") * ns * c_light;  // HepMC distance units are mm
+    fTimeOffset = p.getParameter<double>("TimeOffset") * ns * c_light;  // HepMC distance units are in mm
 
     setBoost(p.getParameter<double>("Alpha") * radian, p.getParameter<double>("Phi") * radian);
     if (fSigmaZ <= 0) {
@@ -53,8 +53,6 @@ BetafuncEvtVtxGenerator::BetafuncEvtVtxGenerator(const edm::ParameterSet& p) : B
   }
 }
 
-BetafuncEvtVtxGenerator::~BetafuncEvtVtxGenerator() {}
-
 void BetafuncEvtVtxGenerator::beginLuminosityBlock(edm::LuminosityBlock const&, edm::EventSetup const& iEventSetup) {
   update(iEventSetup);
 }