Logo ROOT   6.10/09
Reference Guide
SequentialProposal.cxx
Go to the documentation of this file.
1 // @(#)root/roostats:$Id$
2 // Authors: Giovanni Petrucciani 4/21/2011
3 /*************************************************************************
4  * Copyright (C) 1995-2008, Rene Brun and Fons Rademakers. *
5  * All rights reserved. *
6  * *
7  * For the licensing terms see $ROOTSYS/LICENSE. *
8  * For the list of contributors see $ROOTSYS/README/CREDITS. *
9  *************************************************************************/
10 
11 /** \class RooStats::SequentialProposal
12  \ingroup Roostats
13 
14 Class implementing a proposal function that samples the parameter space
15 by moving only in one coordinate (chosen randomly) at each step
16 */
17 
19 #include <RooArgSet.h>
20 #include <iostream>
21 #include <memory>
22 #include <TIterator.h>
23 #include <RooRandom.h>
24 #include <RooStats/RooStatsUtils.h>
25 
26 using namespace std;
27 
29 
30 namespace RooStats {
31 
32 ////////////////////////////////////////////////////////////////////////////////
33 
34 SequentialProposal::SequentialProposal(double divisor) :
36  fDivisor(1./divisor)
37 {
38 }
39 
40 ////////////////////////////////////////////////////////////////////////////////
41 /// Populate xPrime with a new proposed point
42 
44 {
45  RooStats::SetParameters(&x, &xPrime);
46  RooLinkedListIter it(xPrime.iterator());
47  RooRealVar* var;
48  int n = xPrime.getSize();
49  int j = int( floor(RooRandom::uniform()*n) );
50  for (int i = 0; (var = (RooRealVar*)it.Next()) != NULL; ++i) {
51  if (i == j) {
52  double val = var->getVal(), max = var->getMax(), min = var->getMin(), len = max - min;
53  val += RooRandom::gaussian() * len * fDivisor;
54  while (val > max) val -= len;
55  while (val < min) val += len;
56  var->setVal(val);
57  //std::cout << "Proposing a step along " << var->GetName() << std::endl;
58  }
59  }
60 }
61 
62 ////////////////////////////////////////////////////////////////////////////////
63 /// Return the probability of proposing the point x1 given the starting
64 /// point x2
65 
67  return true;
68 }
69 
71  RooArgSet& )
72 {
73  return 1.0; // should not be needed
74 }
75 
76 }
ProposalFunction is an interface for all proposal functions that would be used with a Markov Chain Mo...
static Double_t gaussian(TRandom *generator=randomGenerator())
Return a Gaussian random variable with mean 0 and variance 1.
Definition: RooRandom.cxx:111
virtual Double_t GetProposalDensity(RooArgSet &x1, RooArgSet &x2)
Return the probability of proposing the point x1 given the starting point x2.
void SetParameters(const RooArgSet *desiredVals, RooArgSet *paramsToChange)
Definition: RooStatsUtils.h:58
bool Bool_t
Definition: RtypesCore.h:59
STL namespace.
#define NULL
Definition: RtypesCore.h:88
Double_t x[n]
Definition: legend1.C:17
RooRealVar represents a fundamental (non-derived) real valued object.
Definition: RooRealVar.h:36
Int_t getSize() const
Class implementing a proposal function that samples the parameter space by moving only in one coordin...
double floor(double)
Namespace for the RooStats classes.
Definition: Asimov.h:20
#define ClassImp(name)
Definition: Rtypes.h:336
static Double_t uniform(TRandom *generator=randomGenerator())
Return a number uniformly distributed from (0,1)
Definition: RooRandom.cxx:84
double Double_t
Definition: RtypesCore.h:55
virtual void Propose(RooArgSet &xPrime, RooArgSet &x)
Populate xPrime with a new proposed point.
RooLinkedListIter iterator(Bool_t dir=kIterForward) const
const Int_t n
Definition: legend1.C:16
RooLinkedListIter is the TIterator implementation for RooLinkedList.
virtual Bool_t IsSymmetric(RooArgSet &x1, RooArgSet &x2)
Determine whether or not the proposal density is symmetric for points x1 and x2 - that is...