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=nullptr) 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 TMatrixDSym& covarianceMatrixInverse() const { return _covI ; }
0053   const RooArgList& xVec() const { return _x;}
0054   const RooArgList& muVec() const { return _mu; }
0055
0056   class AnaIntData {
0057   public:
0058     TMatrixD    S22bar ;
0059     double    S22det ;
0060     std::vector<int> pmap ;
0061     Int_t       nint ;
0062   } ;
0063
0064   class GenData {
0065   public:
0066     TMatrixD    UT ;
0067     std::vector<int> omap ;
0068     std::vector<int> pmap ;
0069     TVectorD    mu1 ;
0070     TVectorD    mu2 ;
0071     TMatrixD    S12S22I ;
0072   } ;
0073
0074   class BitBlock {
0075   public:
0076     void setBit(Int_t ibit) ;
0077     bool getBit(Int_t ibit) ;
0078     friend bool operator==(BitBlock const& lhs, BitBlock const& rhs);
0079
0080     Int_t b0 = 0;
0081     Int_t b1 = 0;
0082     Int_t b2 = 0;
0083     Int_t b3 = 0;
0084   } ;
0085
0086   static void blockDecompose(const TMatrixD& input, const std::vector<int>& map1, const std::vector<int>& map2, TMatrixDSym& S11, TMatrixD& S12, TMatrixD& S21, TMatrixDSym& S22) ;
0087
0088 protected:
0089
0090   void decodeCode(Int_t code, std::vector<int>& map1, std::vector<int>& map2) const;
0091   AnaIntData& anaIntData(Int_t code) const ;
0092   GenData& genData(Int_t code) const ;
0093
0094   mutable std::map<int,AnaIntData> _anaIntCache ; ///<!
0095   mutable std::map<int,GenData> _genCache ; ///<!
0096
0097   mutable std::vector<BitBlock> _aicMap ; ///<!
0098
0099   RooListProxy _x ;
0100   RooListProxy _mu ;
0101   TMatrixDSym _cov ;
0102   TMatrixDSym _covI ;
0103   double    _det ;
0104   double    _z ;
0105
0106   void syncMuVec() const ;
0107   mutable TVectorD _muVec ; //! Do not persist
0108
0109   double evaluate() const override ;
0110
0111 private:
0112
0113   ClassDefOverride(RooMultiVarGaussian,1) // Multivariate Gaussian PDF with correlations
0114 };
0115
0116 #endif