Logo ROOT   6.12/07
Reference Guide
ProposalHelper.cxx
Go to the documentation of this file.
1 // @(#)root/roostats:$Id$
2 // Authors: Kevin Belasco 7/22/2009
3 // Authors: Kyle Cranmer 7/22/2009
4 /*************************************************************************
5  * Copyright (C) 1995-2008, Rene Brun and Fons Rademakers. *
6  * All rights reserved. *
7  * *
8  * For the licensing terms see $ROOTSYS/LICENSE. *
9  * For the list of contributors see $ROOTSYS/README/CREDITS. *
10  *************************************************************************/
11 
12 /** \class RooStats::ProposalHelper
13  \ingroup Roostats
14 */
15 
16 #include "Rtypes.h"
18 #include "RooStats/PdfProposal.h"
19 #include "RooStats/RooStatsUtils.h"
20 #include "RooArgSet.h"
21 #include "RooDataSet.h"
22 #include "RooAbsPdf.h"
23 #include "RooAddPdf.h"
24 #include "RooNDKeysPdf.h"
25 #include "RooUniform.h"
26 #include "RooMsgService.h"
27 #include "RooRealVar.h"
28 #include "TIterator.h"
29 #include "RooMultiVarGaussian.h"
30 #include "RooConstVar.h"
31 #include "TString.h"
32 
33 #include <map>
34 
35 namespace RooStats {
36  class ProposalFunction;
37 }
38 
40 
41 using namespace RooFit;
42 using namespace RooStats;
43 using namespace std;
44 
45 //static const Double_t DEFAULT_UNI_FRAC = 0.10;
46 static const Double_t DEFAULT_CLUES_FRAC = 0.20;
47 //static const Double_t SIGMA_RANGE_DIVISOR = 6;
48 static const Double_t SIGMA_RANGE_DIVISOR = 5;
49 //static const Int_t DEFAULT_CACHE_SIZE = 100;
50 //static const Option_t* CLUES_OPTIONS = "a";
51 
52 ////////////////////////////////////////////////////////////////////////////////
53 
55 {
56  fPdfProp = new PdfProposal();
57  fVars = NULL;
58  fOwnsPdfProp = kTRUE;
59  fOwnsPdf = kFALSE;
60  fOwnsCluesPdf = kFALSE;
61  fOwnsVars = kFALSE;
62  fUseUpdates = kFALSE;
63  fPdf = NULL;
64  fSigmaRangeDivisor = SIGMA_RANGE_DIVISOR;
65  fCluesPdf = NULL;
66  fUniformPdf = NULL;
67  fClues = NULL;
68  fCovMatrix = NULL;
69  fCluesFrac = -1;
70  fUniFrac = -1;
71  fCacheSize = -1;
72  fCluesOptions = NULL;
73 }
74 
75 ////////////////////////////////////////////////////////////////////////////////
76 
78 {
79  if (fPdf == NULL)
80  CreatePdf();
81  // kbelasco: check here for memory leaks: does RooAddPdf make copies or
82  // take ownership of components, coeffs
83  RooArgList* components = new RooArgList();
84  RooArgList* coeffs = new RooArgList();
85  if (fCluesPdf == NULL)
86  CreateCluesPdf();
87  if (fCluesPdf != NULL) {
88  if (fCluesFrac < 0)
89  fCluesFrac = DEFAULT_CLUES_FRAC;
90  printf("added clues from dataset %s with fraction %g\n",
91  fClues->GetName(), fCluesFrac);
92  components->add(*fCluesPdf);
93  coeffs->add(RooConst(fCluesFrac));
94  }
95  if (fUniFrac > 0.) {
96  CreateUniformPdf();
97  components->add(*fUniformPdf);
98  coeffs->add(RooConst(fUniFrac));
99  }
100  components->add(*fPdf);
101  RooAddPdf* addPdf = new RooAddPdf("proposalFunction", "Proposal Density",
102  *components, *coeffs);
103  fPdfProp->SetPdf(*addPdf);
104  fPdfProp->SetOwnsPdf(kTRUE);
105  if (fCacheSize > 0)
106  fPdfProp->SetCacheSize(fCacheSize);
107  fOwnsPdfProp = kFALSE;
108  return fPdfProp;
109 }
110 
111 ////////////////////////////////////////////////////////////////////////////////
112 
114 {
115  // kbelasco: check here for memory leaks:
116  // does RooMultiVarGaussian make copies of xVec and muVec?
117  // or should we delete them?
118  if (fVars == NULL) {
119  coutE(InputArguments) << "ProposalHelper::CreatePdf(): " <<
120  "Variables to create proposal function for are not set." << endl;
121  return;
122  }
123  RooArgList* xVec = new RooArgList();
124  RooArgList* muVec = new RooArgList();
125  TIterator* it = fVars->createIterator();
126  RooRealVar* r;
127  RooRealVar* clone;
128  while ((r = (RooRealVar*)it->Next()) != NULL) {
129  xVec->add(*r);
130  TString cloneName = TString::Format("%s%s", "mu__", r->GetName());
131  clone = (RooRealVar*)r->clone(cloneName.Data());
132  muVec->add(*clone);
133  if (fUseUpdates)
134  fPdfProp->AddMapping(*clone, *r);
135  }
136  if (fCovMatrix == NULL)
137  CreateCovMatrix(*xVec);
138  fPdf = new RooMultiVarGaussian("mvg", "MVG Proposal", *xVec, *muVec,
139  *fCovMatrix);
140  delete xVec;
141  delete muVec;
142  delete it;
143 }
144 
145 ////////////////////////////////////////////////////////////////////////////////
146 
148 {
149  Int_t size = xVec.getSize();
150  fCovMatrix = new TMatrixDSym(size);
151  RooRealVar* r;
152  for (Int_t i = 0; i < size; i++) {
153  r = (RooRealVar*)xVec.at(i);
154  Double_t range = r->getMax() - r->getMin();
155  (*fCovMatrix)(i,i) = range / fSigmaRangeDivisor;
156  }
157 }
158 
159 ////////////////////////////////////////////////////////////////////////////////
160 
162 {
163  if (fClues != NULL) {
164  if (fCluesOptions == NULL)
165  fCluesPdf = new RooNDKeysPdf("cluesPdf", "Clues PDF", *fVars, *fClues);
166  else
167  fCluesPdf = new RooNDKeysPdf("cluesPdf", "Clues PDF", *fVars, *fClues,
168  fCluesOptions);
169  }
170 }
171 
172 ////////////////////////////////////////////////////////////////////////////////
173 
175 {
176  fUniformPdf = new RooUniform("uniform", "Uniform Proposal PDF",
177  RooArgSet(*fVars));
178 }
virtual Double_t getMin(const char *name=0) const
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
RooAddPdf is an efficient implementation of a sum of PDFs of the form.
Definition: RooAddPdf.h:29
ProposalFunction is an interface for all proposal functions that would be used with a Markov Chain Mo...
#define coutE(a)
Definition: RooMsgService.h:34
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Add the specified argument to list.
virtual Double_t getMax(const char *name=0) const
virtual TObject * clone(const char *newname) const
Definition: RooRealVar.h:47
static const Double_t DEFAULT_CLUES_FRAC
virtual ProposalFunction * GetProposalFunction()
int Int_t
Definition: RtypesCore.h:41
STL namespace.
Iterator abstract base class.
Definition: TIterator.h:30
static TString Format(const char *fmt,...)
Static method which formats a string using a printf style format descriptor and return a TString...
Definition: TString.cxx:2365
RooRealVar represents a fundamental (non-derived) real valued object.
Definition: RooRealVar.h:36
Int_t getSize() const
ROOT::R::TRInterface & r
Definition: Object.C:4
RooAbsArg * at(Int_t idx) const
Definition: RooArgList.h:84
static const Double_t SIGMA_RANGE_DIVISOR
Flat p.d.f.
Definition: RooUniform.h:24
Generic N-dimensional implementation of a kernel estimation p.d.f.
Definition: RooNDKeysPdf.h:48
const Bool_t kFALSE
Definition: RtypesCore.h:88
Multivariate Gaussian p.d.f.
Namespace for the RooStats classes.
Definition: Asimov.h:20
#define ClassImp(name)
Definition: Rtypes.h:359
double Double_t
Definition: RtypesCore.h:55
TMatrixTSym< Double_t > TMatrixDSym
PdfProposal is a concrete implementation of the ProposalFunction interface.
Definition: PdfProposal.h:30
RooConstVar & RooConst(Double_t val)
virtual TObject * Next()=0
void CreateCovMatrix(RooArgList &xVec)
const Bool_t kTRUE
Definition: RtypesCore.h:87