include/root/RooFFTConvPdf.h

0001 /*****************************************************************************
0002  * Project: RooFit                                                           *
0003  *                                                                           *
0004  * Copyright (c) 2000-2005, Regents of the University of California          *
0005  *                          and Stanford University. All rights reserved.    *
0006  *                                                                           *
0007  * Redistribution and use in source and binary forms,                        *
0008  * with or without modification, are permitted according to the terms        *
0009  * listed in LICENSE (http://roofit.sourceforge.net/license.txt)             *
0010  *****************************************************************************/
0011
0012 #ifndef ROOFFTCONVPDF
0013 #define ROOFFTCONVPDF
0014
0015 #include "RooAbsCachedPdf.h"
0016 #include "RooRealProxy.h"
0017 #include "RooSetProxy.h"
0018 #include "RooAbsReal.h"
0019 #include "RooHistPdf.h"
0020 #include "TVirtualFFT.h"
0021
0022 class RooRealVar;
0023
0024 ///PDF for the numerical (FFT) convolution of two PDFs.
0025 class RooFFTConvPdf : public RooAbsCachedPdf {
0026 public:
0027
0028   RooFFTConvPdf() {
0029     // coverity[UNINIT_CTOR]
0030   } ;
0031   RooFFTConvPdf(const char *name, const char *title, RooRealVar& convVar, RooAbsPdf& pdf1, RooAbsPdf& pdf2, Int_t ipOrder=2);
0032   RooFFTConvPdf(const char *name, const char *title, RooAbsReal& pdfConvVar, RooRealVar& convVar, RooAbsPdf& pdf1, RooAbsPdf& pdf2, Int_t ipOrder=2);
0033   RooFFTConvPdf(const RooFFTConvPdf& other, const char* name=nullptr) ;
0034   TObject* clone(const char* newname) const override { return new RooFFTConvPdf(*this,newname); }
0035   ~RooFFTConvPdf() override ;
0036
0037   void setShift(double val1, double val2) { _shift1 = val1 ; _shift2 = val2 ; }
0038   void setCacheObservables(const RooArgSet& obs) { _cacheObs.removeAll() ; _cacheObs.add(obs) ; }
0039   const RooArgSet& cacheObservables() const { return _cacheObs ; }
0040
0041   /// Return value of buffer fraction applied in FFT calculation array beyond either
0042   /// end of the observable domain to reduce cyclical effects
0043   double bufferFraction() const {
0044     return _bufFrac ;
0045   }
0046
0047   enum BufStrat { Extend=0, Mirror=1, Flat=2 } ;
0048   /// Return the strategy currently used to fill the buffer:
0049   /// 'Extend' means is that the input p.d.f convolution observable range is widened to include the buffer range
0050   /// 'Flat' means that the buffer is filled with the p.d.f. value at the boundary of the observable range
0051   /// 'Mirror' means that the buffer is filled with a mirror image of the p.d.f. around the convolution observable boundary
0052   BufStrat bufferStrategy() const {
0053     return _bufStrat ;
0054   }
0055   void setBufferStrategy(BufStrat bs) ;
0056   void setBufferFraction(double frac) ;
0057
0058   void printMetaArgs(std::ostream& os) const override ;
0059
0060   // Propagate maximum value estimate of pdf1 as convolution can only result in lower max values
0061   Int_t getMaxVal(const RooArgSet& vars) const override { return _pdf1.arg().getMaxVal(vars) ; }
0062   double maxVal(Int_t code) const override { return _pdf1.arg().maxVal(code) ; }
0063
0064
0065 protected:
0066
0067   RooRealProxy _x ;      ///< Convolution observable
0068   RooRealProxy _xprime ; ///< Input function representing value of convolution observable
0069   RooRealProxy _pdf1 ;   ///< First input p.d.f
0070   RooRealProxy _pdf2 ;   ///< Second input p.d.f
0071   RooSetProxy _params ;  ///< Effective parameters of this p.d.f.
0072
0073   void calcParams() ;
0074
0075   std::vector<double>  scanPdf(RooRealVar& obs, RooAbsPdf& pdf, const RooDataHist& hist, const RooArgSet& slicePos, Int_t& N, Int_t& N2, Int_t& zeroBin, double shift) const ;
0076
0077   class FFTCacheElem : public PdfCacheElem {
0078   public:
0079     FFTCacheElem(const RooFFTConvPdf& self, const RooArgSet* nset) ;
0080
0081     RooArgList containedArgs(Action) override ;
0082
0083     std::unique_ptr<TVirtualFFT> fftr2c1;
0084     std::unique_ptr<TVirtualFFT> fftr2c2;
0085     std::unique_ptr<TVirtualFFT> fftc2r;
0086
0087     std::unique_ptr<RooAbsPdf> pdf1Clone;
0088     std::unique_ptr<RooAbsPdf> pdf2Clone;
0089
0090     std::unique_ptr<RooAbsBinning> histBinning;
0091     std::unique_ptr<RooAbsBinning> scanBinning;
0092   };
0093
0094   friend class FFTCacheElem ;
0095
0096   double evaluate() const override { RooArgSet dummy(_x.arg()) ; return getVal(&dummy) ; } ; // dummy
0097   const char* inputBaseName() const override ;
0098   RooFit::OwningPtr<RooArgSet> actualObservables(const RooArgSet& nset) const override ;
0099   RooFit::OwningPtr<RooArgSet> actualParameters(const RooArgSet& nset) const override ;
0100   RooAbsArg& pdfObservable(RooAbsArg& histObservable) const override ;
0101   void fillCacheObject(PdfCacheElem& cache) const override ;
0102   void fillCacheSlice(FFTCacheElem& cache, const RooArgSet& slicePosition) const ;
0103
0104   PdfCacheElem* createCache(const RooArgSet* nset) const override ;
0105   TString histNameSuffix() const override ;
0106
0107   // mutable std::map<const RooHistPdf*,CacheAuxInfo*> _cacheAuxInfo ; //! Auxiliary Cache information (do not persist)
0108   double _bufFrac ; // Sampling buffer size as fraction of domain size
0109   BufStrat _bufStrat ; // Strategy to fill the buffer
0110
0111   double  _shift1 ;
0112   double  _shift2 ;
0113
0114   RooAbsGenContext* genContext(const RooArgSet &vars, const RooDataSet *prototype=nullptr,
0115                                        const RooArgSet* auxProto=nullptr, bool verbose= false) const override ;
0116
0117   friend class RooConvGenContext ;
0118   RooSetProxy  _cacheObs ; ///< Non-convolution observables that are also cached
0119
0120 private:
0121
0122   void prepareFFTBinning(RooRealVar& convVar) const;
0123
0124   ClassDefOverride(RooFFTConvPdf,1) // Convolution operator p.d.f based on numeric Fourier transforms
0125 };
0126
0127 #endif