include/root/RooMultiVarGaussian.h

0001 /*****************************************************************************
0002  * Project: RooFit                                                           *
0003  * Package: RooFitModels                                                     *
0004  *    File: $Id$
0005  * Authors:                                                                  *
0006  *   WV, Wouter Verkerke, UC Santa Barbara, verkerke@slac.stanford.edu       *
0007  *   DK, David Kirkby,    UC Irvine,         dkirkby@uci.edu                 *
0008  *                                                                           *
0009  * Copyright (c) 2000-2005, Regents of the University of California          *
0010  *                          and Stanford University. All rights reserved.    *
0011  *                                                                           *
0012  * Redistribution and use in source and binary forms,                        *
0013  * with or without modification, are permitted according to the terms        *
0014  * listed in LICENSE (http://roofit.sourceforge.net/license.txt)             *
0015  *****************************************************************************/
0016 #ifndef ROO_MULTI_VAR_GAUSSIAN
0017 #define ROO_MULTI_VAR_GAUSSIAN
0018
0019 #include "RooAbsPdf.h"
0020 #include "RooListProxy.h"
0021 #include "TMatrixDSym.h"
0022 #include "TMatrixD.h"
0023 #include "TVectorD.h"
0024
0025 class RooRealVar;
0026 class RooFitResult ;
0027
0028 #include <map>
0029 #include <vector>
0030
0031 class RooMultiVarGaussian : public RooAbsPdf {
0032 public:
0033
0034   RooMultiVarGaussian() {} ;
0035   RooMultiVarGaussian(const char *name, const char *title, const RooArgList& xvec, const RooArgList& mu, const TMatrixDBase& covMatrix) ;
0036   RooMultiVarGaussian(const char *name, const char *title, const RooArgList& xvec, const RooFitResult& fr, bool reduceToConditional=true) ;
0037   RooMultiVarGaussian(const char *name, const char *title, const RooArgList& xvec, const TVectorD& mu, const TMatrixDBase& covMatrix) ;
0038   RooMultiVarGaussian(const char *name, const char *title, const RooArgList& xvec,const TMatrixDBase& covMatrix) ;
0039   void setAnaIntZ(double z) { _z = z ; }
0040
0041   RooMultiVarGaussian(const RooMultiVarGaussian& other, const char* name=nullptr) ;
0042   TObject* clone(const char* newname) const override { return new RooMultiVarGaussian(*this,newname); }
0043
0044   Int_t getAnalyticalIntegral(RooArgSet& allVars, RooArgSet& analVars, const char* rangeName=nullptr) const override ;
0045   double analyticalIntegral(Int_t code, const char* rangeName=nullptr) const override ;
0046
0047   Int_t getGenerator(const RooArgSet& directVars, RooArgSet &generateVars, bool staticInitOK=true) const override;
0048   void initGenerator(Int_t code) override ;
0049   void generateEvent(Int_t code) override;
0050
0051   const TMatrixDSym& covarianceMatrix() const { return _cov ; }
0052   const RooArgList& xVec() const { return _x;}
0053   const RooArgList& muVec() const { return _mu; }
0054
0055   class AnaIntData {
0056   public:
0057     TMatrixD    S22bar ;
0058     double    S22det ;
0059     std::vector<int> pmap ;
0060     Int_t       nint ;
0061   } ;
0062
0063   class GenData {
0064   public:
0065     TMatrixD    UT ;
0066     std::vector<int> omap ;
0067     std::vector<int> pmap ;
0068     TVectorD    mu1 ;
0069     TVectorD    mu2 ;
0070     TMatrixD    S12S22I ;
0071   } ;
0072
0073   class BitBlock {
0074   public:
0075     void setBit(Int_t ibit) ;
0076     bool getBit(Int_t ibit) ;
0077     friend bool operator==(BitBlock const& lhs, BitBlock const& rhs);
0078
0079     Int_t b0 = 0;
0080     Int_t b1 = 0;
0081     Int_t b2 = 0;
0082     Int_t b3 = 0;
0083   } ;
0084
0085   static void blockDecompose(const TMatrixD& input, const std::vector<int>& map1, const std::vector<int>& map2, TMatrixDSym& S11, TMatrixD& S12, TMatrixD& S21, TMatrixDSym& S22) ;
0086
0087 protected:
0088
0089   void decodeCode(Int_t code, std::vector<int>& map1, std::vector<int>& map2) const;
0090   AnaIntData& anaIntData(Int_t code) const ;
0091   GenData& genData(Int_t code) const ;
0092
0093   mutable std::map<int,AnaIntData> _anaIntCache ; ///<!
0094   mutable std::map<int,GenData> _genCache ; ///<!
0095
0096   mutable std::vector<BitBlock> _aicMap ; ///<!
0097
0098   RooListProxy _x ;
0099   RooListProxy _mu ;
0100   TMatrixDSym _cov ;
0101   TMatrixDSym _covI ;
0102   double    _det ;
0103   double    _z ;
0104
0105   void syncMuVec() const ;
0106   mutable TVectorD _muVec ; //! Do not persist
0107
0108   double evaluate() const override ;
0109
0110 private:
0111
0112   ClassDefOverride(RooMultiVarGaussian,1) // Multivariate Gaussian PDF with correlations
0113 };
0114
0115 #endif