Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 14:31:14

0001 //
0002 //
0003 // File: hitfit/Base_Constrainer.h
0004 // Purpose: Abstract base for the chisq fitter classes.
0005 //          This allows for different algorithms to be used.
0006 // Created: Jul, 2000, sss, based on run 1 mass analysis code.
0007 //
0008 // CMSSW File      : interface/Base_Constrainer.h
0009 // Original Author : Scott Stuart Snyder <snyder@bnl.gov> for D0
0010 // Imported to CMSSW by Haryo Sumowidagdo <Suharyo.Sumowidagdo@cern.ch>
0011 
0012 /**
0013 
0014     @file Base_Constrainer.h
0015 
0016     @brief Abstract base classes for the \f$\chi^{2}\f$ fitter classes.
0017 
0018     This file contains abstract base classes for the \f$\chi^{2}\f$ fitter
0019     classes.  By doing so, different fitting algorithm may be used.
0020 
0021     @par Creation date:
0022     July 2000.
0023 
0024     @author
0025     Scott Stuart Snyder <snyder@bnl.gov>.
0026 
0027     @par Modification History:
0028     Apr 2009: Haryo Sumowidagdo <Suharyo.Sumowidagdo@cern.ch>:
0029     Imported to CMSSW.<br>
0030     Oct 2009: Haryo Sumowidagdo <Suharyo.Sumowidagdo@cern.ch>:
0031     Added Doxygen tags for automatic generation of documentation.
0032 
0033     @par Terms of Usage:
0034     With consent from the original author (Scott Snyder).
0035 
0036  */
0037 #ifndef HITFIT_BASE_CONSTRAINER_H
0038 #define HITFIT_BASE_CONSTRAINER_H
0039 
0040 #include "TopQuarkAnalysis/TopHitFit/interface/matutil.h"
0041 #include <iosfwd>
0042 
0043 namespace hitfit {
0044 
0045   class Defaults;
0046 
0047   //*************************************************************************
0048 
0049   /**
0050     @class Base_Constrainer_Args
0051 
0052     @brief Hold on to parameters for the Base_Constrainer class.
0053  */
0054   class Base_Constrainer_Args
0055   //
0056   // Purpose: Hold on to parameters for the Base_Constrainer class.
0057   //
0058   // Parameters:
0059   //   bool test_gradient - If true, check the constraint gradient calculations
0060   //                        by also doing them numerically.
0061   //   float test_step    - When test_gradient is true, the step size to use
0062   //                        for numeric differentiation.
0063   //   float test_eps     - When test_gradient is true, the maximum relative
0064   //                        difference permitted between returned and
0065   //                        numerically calculated gradients.
0066   //
0067   {
0068   public:
0069     // Constructor.  Initialize from a Defaults object.
0070 
0071     /**
0072      Instantiate Base_Constrainer_Args from an instance of Defaults object.
0073 
0074      @param defs An instance of Defaults object. The instance must contain the
0075      variables of type and name:
0076      - bool <i>test_gradient</i>.
0077      - double <i>test_step</i>.
0078      - double <i>test_eps</i>.
0079    */
0080     Base_Constrainer_Args(const Defaults& defs);
0081 
0082     // Retrieve parameter values.
0083     /**
0084      Return the <i>_test_gradient</i> parameter.
0085    */
0086     bool test_gradient() const;
0087 
0088     /**
0089      Return the <i>_test_step</i> parameter.
0090    */
0091     double test_step() const;
0092 
0093     /**
0094      Return the <i>_test_eps</i> parameter.
0095    */
0096     double test_eps() const;
0097 
0098   private:
0099     // Hold on to parameter values.
0100 
0101     /**
0102      If true, check constraint gradient calculation by also doing
0103      them numerically.
0104    */
0105     bool _test_gradient;
0106 
0107     /**
0108      When <i>_test_gradient</i> is true, the step size use for numerical
0109      differentation.
0110    */
0111     double _test_step;
0112 
0113     /**
0114      When <i>_test_gradient</i> is true, the maximum relative difference
0115      permitted between returned and numerically calculated gradients.
0116    */
0117     double _test_eps;
0118   };
0119 
0120   //*************************************************************************
0121 
0122   /**
0123 
0124    @class Constraint_Calculator
0125 
0126    @brief Abstract base class for evaluating constraints.  Users derive
0127    from this class and implement the eval() method.
0128 
0129  */
0130   class Constraint_Calculator
0131   //
0132   // Purpose: Abstract base class for evaluating constraints.
0133   //          Derive from this and implement the eval() method.
0134   //
0135   {
0136   public:
0137     // Constructor, destructor.  Pass in the number of constraints.
0138 
0139     /**
0140      Constructor.
0141      @param nconstraints Number of constraint equations.
0142    */
0143     Constraint_Calculator(int nconstraints);
0144 
0145     /**
0146      Destructor.
0147    */
0148     virtual ~Constraint_Calculator() {}
0149 
0150     // Get back the number of constraints.
0151     /**
0152      Return the number of constraints.
0153    */
0154     int nconstraints() const;
0155 
0156     // Evaluate constraints at the point described by X and Y (well-measured
0157     // and poorly-measured variables, respectively).  The results should
0158     // be stored in F.  BX and BY should be set to the gradients of F with
0159     // respect to X and Y, respectively.
0160     //
0161     // Return true if the point X, Y is accepted.
0162     // Return false if it is rejected (i.e., in an unphysical region).
0163     // The constraints need not be evaluated in that case.
0164 
0165     /**
0166     @brief Evaluate constraints at the point described by <i>x</i> and
0167     <i>y</i>
0168     (well-measured and poorly-measured variables, respectively).  The results
0169     should be stored in <i>F</i>.  <i>Bx</i> and <i>By</i> should be set to
0170     the gradients of <i>F</i> with respect to <i>x</i> and <i>y</i>,
0171     respectively.
0172 
0173     @param x  Column_Vector of well-measured variables.
0174 
0175     @param y  Column_Vector of poorly-measured variables.
0176 
0177     @param F  Row_Vector contains the results of the constraint evaluation.
0178 
0179     @param Bx Gradients of <i>F</i> with respect to <i>x</i>
0180 
0181     @param By Gradients of <i>F</i> with respect to <i>y</i>
0182 
0183     @par Output:
0184     - <i>F</i>.
0185     - <i>Bx</i>.
0186     - <i>By</i>.
0187 
0188     @par Return:
0189     <b>true</b> if the point <i>(x,y)</i> is accepted.<br>
0190     <b>false</b> if the point <i>(x,y)</i> is rejected
0191     (i.e., in an unphysical region).  The constraints need not be
0192     evaluated in that case.
0193 
0194    */
0195     virtual bool eval(const Column_Vector& x, const Column_Vector& y, Row_Vector& F, Matrix& Bx, Matrix& By) = 0;
0196 
0197   private:
0198     // The number of constraint functions.
0199     /**
0200      Number of constraints functions.
0201    */
0202     int _nconstraints;
0203   };
0204 
0205   //*************************************************************************
0206 
0207   /**
0208     @class Base_Constrainer
0209     @brief Base class for \f$\chi^{2}\f$ constrained fitter.
0210  */
0211   class Base_Constrainer
0212   //
0213   // Purpose: Base class for chisq constrained fitter.
0214   //
0215   {
0216   public:
0217     // Constructor, destructor.
0218     // ARGS holds the parameter settings for this instance.
0219 
0220     /**
0221      Constructor.
0222      @param args Contains the parameter settings for this instance.
0223    */
0224     Base_Constrainer(const Base_Constrainer_Args& args);
0225 
0226     /**
0227      Destructor.
0228    */
0229     virtual ~Base_Constrainer() {}
0230 
0231     // Do the fit.
0232     // Call the number of well-measured variables Nw, the number of
0233     // poorly-measured variables Np, and the number of constraints Nc.
0234     // Inputs:
0235     //   CONSTRAINT_CALCULATOR is the object that will be used to evaluate
0236     //     the constraints.
0237     //   XM(Nw) and YM(Np) are the measured values of the well- and
0238     //     poorly-measured variables, respectively.
0239     //   X(Nw) and Y(Np) are the starting values for the fit.
0240     //   G_I(Nw,Nw) is the error matrix for the well-measured variables.
0241     //   Y(Np,Np) is the inverse error matrix for the poorly-measured variables.
0242     //
0243     // Outputs:
0244     //   X(Nw) and Y(Np) is the point at the minimum.
0245     //   PULLX(Nw) and PULLY(Np) are the pull quantities.
0246     //   Q(Nw,Nw), R(Np,Np), and S(Nw,Np) are the final error matrices
0247     //     between all the variables.
0248     //
0249     // The return value is the final chisq.  Returns a value < 0 if the
0250     // fit failed to converge.
0251 
0252     /**
0253      @brief Perform the \f$\chi^{2}\f$ constrained fit.
0254 
0255 
0256      @param constraint_calculator The object that will be used to evaluate
0257      the constraints.
0258 
0259      @param xm Measured values of the well-measured variables, has dimension
0260      <i>Nw</i>.
0261 
0262      @param x  Before the fit: starting value of the well-measured variables
0263      for the fit, has dimension <i>Nw</i>.
0264      After the fit: final values of the well-measured variables.
0265 
0266      @param ym Measured values of the poorly-measured variables, has
0267      dimension <i>Np</i>.
0268 
0269      @param y  Before the fit: starting value of the poorly-measured variables
0270      for the fit, has dimension <i>Np</i>.
0271      After the fit: final values of the poorly-measured variables.
0272 
0273      @param G_i Error matrix for the well-measured variables, has dimension
0274      <i>Nw,Nw</i>.
0275 
0276      @param Y Inverse error matrix for the poorly-measured variables, has
0277      dimension <i>Np,Np</i>.
0278 
0279      @param pullx Pull quantities for the well-measured variables, has
0280      dimension <i>Nw</i>.
0281 
0282      @param pully Pull quantities for the poorly-measured variables, has
0283      dimension <i>Np</i>.
0284 
0285      @param Q Error matrix for the well-measured variables, has dimension
0286      <i>Nw,Nw</i>.
0287 
0288      @param R Error matrix for the poorly-measured variables, has dimension
0289      <i>Np,Np</i>.
0290 
0291      @param S Error matrix for the correlation between well-measured variables
0292      and poorly-measured variables, has dimension <i>Nw,Np</i>.
0293 
0294      @par Input:
0295      <i>constraint_calculator, xm, x, ym, y, G_i, y</i>.
0296 
0297      @par Output:
0298      <i>x, y, pullx, pully, Q, R, S</i>.
0299      @par Return:
0300      \f$\chi^{2}\f$ of the fit.  Should returns a negative value if the fit
0301      does not converge.
0302    */
0303     virtual double fit(Constraint_Calculator& constraint_calculator,
0304                        const Column_Vector& xm,
0305                        Column_Vector& x,
0306                        const Column_Vector& ym,
0307                        Column_Vector& y,
0308                        const Matrix& G_i,
0309                        const Diagonal_Matrix& Y,
0310                        Column_Vector& pullx,
0311                        Column_Vector& pully,
0312                        Matrix& Q,
0313                        Matrix& R,
0314                        Matrix& S) = 0;
0315 
0316     // Print out any internal state to S.
0317     /**
0318     @brief Print out internal state to output stream.
0319 
0320     @param s Output stream to which to write.
0321 
0322     @par Return:
0323     The stream <i>s</i>;
0324 
0325    */
0326     virtual std::ostream& print(std::ostream& s) const;
0327 
0328     // Print out internal state to S.
0329     friend std::ostream& operator<<(std::ostream& s, const Base_Constrainer& f);
0330 
0331   private:
0332     // Parameter settings.
0333 
0334     /**
0335      Parameter settings for this instance of Base_Constrainer.
0336    */
0337     const Base_Constrainer_Args _args;
0338 
0339   protected:
0340     // Helper function to evaluate the constraints.
0341     // This takes care of checking what the user function returns against
0342     // numerical derivatives, if that was requested.
0343 
0344     /**
0345      @brief Helper function to evaluate constraints.  This takes care of
0346      checking what the user function returns againts numerical
0347      derivatives, it that was requested.
0348 
0349      @par Purpose
0350      Call the constraint function for the point <i>(x,y)</i>.  Return
0351      <i>F</i>, <i>Bx</i>, <i>By</i>, and a boolean flag saying if the point
0352      is acceptable.  If <i>test_gradient</i> is on, we verify the gradients
0353      returned by also computing them numerically.
0354 
0355      @param constraint_calculator the User-supplied object to evaluate the
0356      constraints.
0357 
0358      @param x Vector of well-measured quantities where we evaluate the
0359      constraints.  Has dimension of <i>Nw</i>.
0360 
0361      @param y Vector of poorly-measured quantities where we evaluate the
0362      the constraints.  Has dimension of <i>Np</i>.
0363 
0364      @param F Results of the constraint functions.
0365 
0366      @param Bx Gradients of <i>F</i> with respect to <i>x</i>, has dimension
0367      <i>(Nw,Nc)</i>.
0368 
0369      @param By Gradients of <i>F</i> with respect to <i>y</i>, has dimension
0370      <i>(Np,Nc)</i>.
0371 
0372      @par Input:
0373      <i>constraint_calculator, x, y.</i>
0374 
0375      @par Output:
0376      <i>F, Bx, By.</i>
0377 
0378      @par Return:
0379      <b>true</b> if the point <i>(x,y)</i> is accepted.<br>
0380      <b>false</b> if the point <i>(x,y)</i> is rejected
0381 
0382    */
0383     bool call_constraint_fcn(Constraint_Calculator& constraint_calculator,
0384                              const Column_Vector& x,
0385                              const Column_Vector& y,
0386                              Row_Vector& F,
0387                              Matrix& Bx,
0388                              Matrix& By) const;
0389   };
0390 
0391 }  // namespace hitfit
0392 
0393 #endif  // not HITFIT_BASE_CONSTRAINER_H