Refactoring of DB Monitoring code

CalibCalorimetry/HcalAlgos/interface/AbsElectronicODERHS.h

@@ -0,0 +1,144 @@
+#ifndef CalibCalorimetry_HcalAlgos_AbsElectronicODERHS_h_
+#define CalibCalorimetry_HcalAlgos_AbsElectronicODERHS_h_
+
+#include <vector>
+#include <climits>
+#include <algorithm>
+
+#include "CalibCalorimetry/HcalAlgos/interface/AbsODERHS.h"
+#include "CondFormats/HcalObjects/interface/HcalInterpolatedPulse.h"
+
+//
+// Modeling of electronic circuits always involves an input pulse that
+// determines the circuit output. This class adds an input pulse to
+// AbsODERHS and establishes a uniform interface to circuit parameters.
+//
+class AbsElectronicODERHS : public AbsODERHS
+{
+public:
+    static const unsigned invalidNode = UINT_MAX - 1U;
+
+    inline AbsElectronicODERHS() : initialized_(false), allSet_(false) {}
+
+    inline explicit AbsElectronicODERHS(const HcalInterpolatedPulse& pulse)
+        : inputPulse_(pulse) {}
+
+    inline virtual ~AbsElectronicODERHS() {}
+
+    inline const HcalInterpolatedPulse& inputPulse() const {return inputPulse_;}
+
+    inline HcalInterpolatedPulse& inputPulse() {return inputPulse_;}
+
+    template<class Pulse>
+    inline void setInputPulse(const Pulse& pulse) {inputPulse_ = pulse;}
+
+    // The following methods must be overriden by derived classes.
+    // Total number of nodes included in the simulation:
+    virtual unsigned numberOfNodes() const = 0;
+
+    // The node which counts as "output" (preamp output in case
+    // of QIE8, which is not necessarily the node which accumulates
+    // the charge):
+    virtual unsigned outputNode() const = 0;
+
+    // The node which counts as "control". If this method returns
+    // "invalidNode" then there is no such node in the circuit.
+    virtual unsigned controlNode() const {return invalidNode;}
+
+    // The number of simulation parameters:
+    virtual unsigned nParameters() const = 0;
+
+    // Check if all parameters have been set
+    inline bool allParametersSet() const
+    {
+        // Raise "allSet_" flag if all parameters have been set
+        if (!allSet_)
+        {
+            const unsigned nExpected = this->nParameters();
+            if (nExpected)
+            {
+                if (paramMask_.size() != nExpected)
+                    return false;
+                unsigned count = 0;
+                const unsigned char* mask = &paramMask_[0];
+                for (unsigned i=0; i<nExpected; ++i)
+                    count += mask[i];
+                allSet_ = count == nExpected;
+            }
+            else
+                allSet_ = true;
+        }
+        return allSet_;
+    }
+
+    inline void setParameter(const unsigned which, const double value)
+    {
+        if (!initialized_)
+            initialize();
+        paramMask_.at(which) = 1;
+        params_[which] = value;
+    }
+
+    inline double getParameter(const unsigned which) const
+    {
+        if (!paramMask_.at(which)) throw cms::Exception(
+            "In AbsElectronicODERHS::getParameter: no such parameter or "
+            "parameter value is not established yet");
+        return params_[which];
+    }
+
+    inline const std::vector<double>& getAllParameters() const
+    {
+        if (!allParametersSet()) throw cms::Exception(
+            "In AbsElectronicODERHS::getAllParameters: "
+            "some parameter values were not established yet");
+        return params_;
+    }
+
+    inline void setLeadingParameters(const double* values, const unsigned len)
+    {
+        if (len)
+        {
+            assert(values);
+            if (!initialized_)
+                initialize();
+            const unsigned sz = params_.size();
+            const unsigned imax = std::min(sz, len);
+            for (unsigned i=0; i<imax; ++i)
+            {
+                params_[i] = values[i];
+                paramMask_[i] = 1;
+            }
+        }
+    }
+
+    inline void setLeadingParameters(const std::vector<double>& values)
+    {
+        if (!values.empty())
+            setLeadingParameters(&values[0], values.size());
+    }
+
+protected:
+    HcalInterpolatedPulse inputPulse_;
+    std::vector<double> params_;
+
+private:
+    std::vector<unsigned char> paramMask_;
+    bool initialized_;
+    mutable bool allSet_;
+
+    inline void initialize()
+    {
+        const unsigned nExpected = this->nParameters();
+        if (nExpected)
+        {
+            params_.resize(nExpected);
+            paramMask_.resize(nExpected);
+            for (unsigned i=0; i<nExpected; ++i)
+                paramMask_[i] = 0;
+        }
+        initialized_ = true;
+    }
+};
+
+#endif // CalibCalorimetry_HcalAlgos_AbsElectronicODERHS_h_

CalibCalorimetry/HcalAlgos/interface/AbsODERHS.h

@@ -0,0 +1,18 @@
+#ifndef CalibCalorimetry_HcalAlgos_AbsODERHS_h_
+#define CalibCalorimetry_HcalAlgos_AbsODERHS_h_
+
+//
+// Base class for the ODE right hand sides
+//
+class AbsODERHS
+{
+public:
+    inline virtual ~AbsODERHS() {}
+
+    virtual AbsODERHS* clone() const = 0;
+
+    virtual void calc(double t, const double* x, unsigned lenX,
+                      double* derivative) = 0;
+};
+
+#endif // CalibCalorimetry_HcalAlgos_AbsODERHS_h_

CalibCalorimetry/HcalAlgos/interface/ConstantStepOdeSolver.h

@@ -0,0 +1,184 @@
+#ifndef CalibCalorimetry_HcalAlgos_ConstantStepOdeSolver_h_
+#define CalibCalorimetry_HcalAlgos_ConstantStepOdeSolver_h_
+
+#include <vector>
+#include <iostream>
+#include "FWCore/Utilities/interface/Exception.h"
+
+#include "CalibCalorimetry/HcalAlgos/interface/AbsODERHS.h"
+
+//
+// ODE solver with a constant time step. The derived classes are supposed
+// to implement concrete ODE solving algorithms (Runge-Kutta, etc).
+//
+class ConstantStepOdeSolver
+{
+public:
+    inline ConstantStepOdeSolver()
+        : rhs_(0), dt_(0.0), dim_(0), runLen_(0), lastIntegrated_(0) {}
+
+    inline ConstantStepOdeSolver(const AbsODERHS& rhs) :
+        rhs_(0), dt_(0.0), dim_(0), runLen_(0), lastIntegrated_(0)
+    {
+        rhs_ = rhs.clone();
+    }
+
+    // The copy constructor and the assignment operator are explicitly provided
+    ConstantStepOdeSolver(const ConstantStepOdeSolver& r);
+    ConstantStepOdeSolver& operator=(const ConstantStepOdeSolver& r);
+
+    inline virtual ~ConstantStepOdeSolver() {delete rhs_;}
+
+    // Access the equation right hand side
+    inline void setRHS(const AbsODERHS& rhs)
+    {
+        delete rhs_;
+        rhs_ = rhs.clone();
+    }
+    inline const AbsODERHS* getRHS() const {return rhs_;}
+    inline AbsODERHS* getRHS() {return rhs_;}
+
+    // Inspectors (will be valid after at least one "run" call)
+    inline unsigned lastDim() const {return dim_;}
+    inline unsigned lastRunLength() const {return runLen_;}
+    inline double lastDeltaT() const {return dt_;}
+    inline double lastMaxT() const {return runLen_ ? dt_*(runLen_-1U) : 0.0;}
+
+    inline double getTime(const unsigned idx) const
+    {
+        if (idx >= runLen_) throw cms::Exception(
+            "In ConstantStepOdeSolver::getTime: index out of range");
+        return idx*dt_;
+    }
+
+    inline double getCoordinate(const unsigned which, const unsigned idx) const
+    {
+        if (which >= dim_ || idx >= runLen_) throw cms::Exception(
+            "In ConstantStepOdeSolver::getCoordinate: index out of range");
+        return historyBuffer_[dim_*idx + which];
+    }
+
+    // Integrate the node with the given number and get
+    // the value of the integral at the given history point
+    double getIntegrated(unsigned which, unsigned idx) const;
+
+    // Truncate some coordinate
+    void truncateCoordinate(unsigned which, double minValue, double maxValue);
+
+    // Linear interpolation methods will be used in case the "cubic"
+    // argument is false, and cubic in case it is true
+    double interpolateCoordinate(unsigned which, double t,
+                                 bool cubic = false) const;
+
+    // Interpolate the integrated node
+    double interpolateIntegrated(unsigned which, double t,
+                                 bool cubic = false) const;
+
+    // Get the time of the peak position
+    double getPeakTime(unsigned which) const;
+
+    // Solve the ODE and remember the history
+    void run(const double* initialConditions, unsigned lenConditions,
+             double dt, unsigned nSteps);
+
+    // Set the history from some external source. The size
+    // of the data array should be at least dim*runLen.
+    void setHistory(double dt, const double* data,
+                    unsigned dim, unsigned runLen);
+
+    // Write out the history
+    void writeHistory(std::ostream& os, double dt, bool cubic = false) const;
+
+    // Write out the integrated node
+    void writeIntegrated(std::ostream& os, unsigned which,
+                         double dt, bool cubic = false) const;
+
+    // The following method must be overriden by derived classes
+    virtual const char* methodName() const = 0;
+
+protected:
+    AbsODERHS* rhs_;
+
+private:
+    // The following method must be overriden by derived classes
+    virtual void step(double t, double dt,
+                      const double* x, unsigned lenX,
+                      double* coordIncrement) const = 0;
+
+    // The following integration corresponds to the cubic
+    // interpolation of the coordinate
+    void integrateCoordinate(const unsigned which);
+
+    double dt_;
+    unsigned dim_;
+    unsigned runLen_;
+    unsigned lastIntegrated_;
+
+    std::vector<double> historyBuffer_;
+    std::vector<double> chargeBuffer_;
+};
+
+
+// Dump the coordinate history as it was collected
+inline std::ostream& operator<<(std::ostream& os,
+                                const ConstantStepOdeSolver& s)
+{
+    s.writeHistory(os, s.lastDeltaT());
+    return os;
+}
+
+
+// A few concrete ODE solvers
+class EulerOdeSolver : public ConstantStepOdeSolver
+{
+public:
+    inline EulerOdeSolver() : ConstantStepOdeSolver() {}
+
+    inline explicit EulerOdeSolver(const AbsODERHS& rhs)
+        : ConstantStepOdeSolver(rhs) {}
+
+    inline const char* methodName() const {return "Euler";}
+
+private:
+    void step(double t, double dt,
+              const double* x, unsigned lenX,
+              double* coordIncrement) const;
+};
+
+
+class RK2 : public ConstantStepOdeSolver
+{
+public:
+    inline RK2() : ConstantStepOdeSolver() {}
+
+    inline explicit RK2(const AbsODERHS& rhs) : ConstantStepOdeSolver(rhs) {}
+
+    inline const char* methodName() const {return "2nd order Runge-Kutta";}
+
+private:
+    void step(double t, double dt,
+              const double* x, unsigned lenX,
+              double* coordIncrement) const;
+
+    mutable std::vector<double> buf_;
+};
+
+
+class RK4 : public ConstantStepOdeSolver
+{
+public:
+    inline RK4() : ConstantStepOdeSolver() {}
+
+    inline explicit RK4(const AbsODERHS& rhs) : ConstantStepOdeSolver(rhs) {}
+
+    inline const char* methodName() const {return "4th order Runge-Kutta";}
+
+private:
+    void step(double t, double dt,
+              const double* x, unsigned lenX,
+              double* coordIncrement) const;
+
+    mutable std::vector<double> buf_;
+};
+
+#endif // CalibCalorimetry_HcalAlgos_ConstantStepOdeSolver_h_