discretepdf.cpp
Go to the documentation of this file.
1 // $Id$
2 // Copyright (C) 2002 Klaas Gadeyne <first dot last at gmail dot com>
3 //
4 // This program is free software; you can redistribute it and/or modify
5 // it under the terms of the GNU Lesser General Public License as published by
6 // the Free Software Foundation; either version 2.1 of the License, or
7 // (at your option) any later version.
8 //
9 // This program is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12 // GNU Lesser General Public License for more details.
13 //
14 // You should have received a copy of the GNU Lesser General Public License
15 // along with this program; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
17 //
18 #include "discretepdf.h"
19 #include "../bfl_err.h"
20 #include "../wrappers/rng/rng.h"
21 #include <vector>
22 #include <iostream>
23 
24 
25 namespace BFL
26 {
27  using namespace std;
28  using namespace MatrixWrapper;
29 
30 
31  DiscretePdf::DiscretePdf(unsigned int num_states): Pdf<int>(1)
32  ,_num_states(num_states)
33  {
34  //discrete pdf has dimension 1
35  _Values_p = new vector<Probability>(num_states);
36  for (int i=0; i<NumStatesGet();i++)
37  {
38  (*_Values_p)[i] = (Probability)(1.0/NumStatesGet());
39  }
40  _CumPDF.insert(_CumPDF.begin(),num_states+1,0.0);
41  CumPDFUpdate();
42 #ifdef __CONSTRUCTOR__
43  cout << "DiscretePdf constructor\n";
44 #endif // __CONSTRUCTOR__
45  }
46 
47  DiscretePdf::DiscretePdf(const DiscretePdf & my_dpdf):Pdf<int>(my_dpdf)
48  ,_num_states(my_dpdf.NumStatesGet())
49  {
50  _Values_p = new vector<Probability>(this->NumStatesGet());
51  (*_Values_p) = my_dpdf.ProbabilitiesGet();
52  _CumPDF.insert(_CumPDF.begin(),NumStatesGet()+1,0.0);
53  CumPDFUpdate();
54 #ifdef __CONSTRUCTOR__
55  cout << "DiscretePdf copy constructor\n";
56 #endif // __CONSTRUCTOR__
57  }
58 
60  {
61 #ifdef __CONSTRUCTOR__
62  cout << "DiscretePdf destructor\n";
63 #endif
64  // Release memory!
65  delete _Values_p;
66  }
67 
68  //Clone function
70  {
71  return new DiscretePdf(*this);
72  }
73 
74  unsigned int DiscretePdf::NumStatesGet()const
75  {
76  return _num_states;
77  }
78 
79 
80  Probability DiscretePdf::ProbabilityGet(const int& state) const
81  {
82  assert((int)state >= 0 && state < NumStatesGet());
83  return (*_Values_p)[state];
84  }
85 
87  {
88  assert((int)state >= 0 && state < NumStatesGet());
89  assert(a<=1);
90 
91  // renormalize other probabilities such that sum of probabilities will be
92  // one after the probability of state is set to a
93  // This should keep the probabilities normalized
94  Probability old_prob_state = ProbabilityGet(state);
95  if (old_prob_state!=1) {
96  double normalization_factor = (1-a)/(1-old_prob_state);
97  for (int i=0; i<NumStatesGet();i++)
98  {
99  (*_Values_p)[i] = (Probability)( (double)( (*_Values_p)[i] )* normalization_factor);
100  }
101  }
102  else{
103  for (int i=0; i<NumStatesGet();i++)
104  {
105  (*_Values_p)[i] = (Probability)( (1-a)/(NumStatesGet()-1));
106  }
107  }
108  (*_Values_p)[state] = a;
109  return CumPDFUpdate();
110  }
111 
112  vector<Probability> DiscretePdf::ProbabilitiesGet() const
113  {
114  return *_Values_p;
115  }
116 
117  bool DiscretePdf::ProbabilitiesSet(vector<Probability> & v)
118  {
119  assert(v.size() == NumStatesGet());
120 
121  *_Values_p = v;
122  //normalize the probabilities and update the cumulative sum
123  return (NormalizeProbs() && CumPDFUpdate());
124  }
125 
126  // For optimal performance!
127  bool
128  DiscretePdf::SampleFrom (vector<Sample<int> >& list_samples,
129  const unsigned int num_samples,
130  int method,
131  void * args) const
132  {
133  switch(method)
134  {
135  case DEFAULT: // O(N log(N) efficiency)
136  return Pdf<int>::SampleFrom(list_samples, num_samples,method,args);
137 
138  case RIPLEY: // See mcpdf.cpp for more explanation
139  {
140  list_samples.resize(num_samples);
141  // GENERATE N ORDERED IID UNIFORM SAMPLES
142  std::vector<double> unif_samples(num_samples);
143  for ( unsigned int i = 0; i < num_samples ; i++)
144  unif_samples[i] = runif();
145 
146  /* take n-th racine of u_N */
147  unif_samples[num_samples-1] = pow(unif_samples[num_samples-1],
148  double (1.0/num_samples));
149  /* rescale samples */
150  for ( int i = num_samples-2; i >= 0 ; i--)
151  unif_samples[i] = pow(unif_samples[i], double (1.0/(i+1))) * unif_samples[i+1];
152 
153  // CHECK WHERE THESE SAMPLES ARE IN _CUMPDF
154  unsigned int index = 0;
155  unsigned int num_states = NumStatesGet();
156  vector<double>::const_iterator CumPDFit = _CumPDF.begin();
157  vector<Sample<int> >::iterator sit = list_samples.begin();
158 
159  for ( unsigned int i = 0; i < num_samples ; i++)
160  {
161  while ( unif_samples[i] > *CumPDFit )
162  {
163  // check for internal error
164  assert(index <= num_states);
165  index++; CumPDFit++;
166  }
167  int a = index - 1;
168  sit->ValueSet(a);
169  sit++;
170  }
171  return true;
172  }
173  default:
174  cerr << "DiscretePdf::Samplefrom(int, void *): No such sampling method" << endl;
175  return false;
176  }
177  }
178 
179 
180 
181  bool DiscretePdf::SampleFrom (Sample<int>& one_sample, int method, void * args) const
182  {
183  switch(method)
184  {
185  case DEFAULT:
186  {
187  // Sample from univariate uniform rng between 0 and 1;
188  double unif_sample; unif_sample = runif();
189  // Compare where we should be
190  unsigned int index = 0;
191  while ( unif_sample > _CumPDF[index] )
192  {
193  assert(index <= NumStatesGet());
194  index++;
195  }
196  int a = index - 1;
197  one_sample.ValueSet(a);
198  return true;
199  }
200  default:
201  cerr << "DiscretePdf::Samplefrom(int, void *): No such sampling method"
202  << endl;
203  return false;
204  }
205  }
206 
208  {
209  double SumOfProbs = 0.0;
210  for ( unsigned int i = 0; i < NumStatesGet() ; i++){
211  SumOfProbs += (*_Values_p)[i];
212  }
213  if (SumOfProbs > 0){
214  for ( unsigned int i = 0; i < NumStatesGet() ; i++){
215  (*_Values_p)[i] = (Probability)( (double) ( (*_Values_p)[i]) /SumOfProbs);
216  }
217  return true;
218  }
219  else{
220  cerr << "DiscretePdf::NormalizeProbs(): SumOfProbs = " << SumOfProbs << endl;
221  return false;
222  }
223  }
224 
226  {
227  // precondition: sum of probabilities should be 1
228  double CumSum=0.0;
229  static vector<double>::iterator CumPDFit;
230  CumPDFit = _CumPDF.begin();
231  *CumPDFit = 0.0;
232 
233  // Calculate the Cumulative PDF
234  for ( unsigned int i = 0; i < NumStatesGet(); i++)
235  {
236  CumPDFit++;
237  // Calculate the __normalised__ Cumulative sum!!!
238  CumSum += ( (*_Values_p)[i] );
239  *CumPDFit = CumSum;
240  }
241  // Check if last element of valuelist is +- 1
242  assert( (_CumPDF[NumStatesGet()] >= 1.0 - NUMERIC_PRECISION) &&
243  (_CumPDF[NumStatesGet()] <= 1.0 + NUMERIC_PRECISION) );
244 
245  _CumPDF[NumStatesGet()]=1;
246 
247  return true;
248  }
249 
251  {
252  int index_mostProbableState = -1;
253  Probability prob_mostProbableState = 0.0;
254  for (int state = 0 ; state < NumStatesGet() ; state++)
255  {
256  if ( (*_Values_p)[state] >= prob_mostProbableState)
257  {
258  index_mostProbableState = state;
259  prob_mostProbableState = (*_Values_p)[state];
260  }
261  }
262  return index_mostProbableState;
263  }
264 
265 
266 } // End namespace
virtual ~DiscretePdf()
Destructor.
#define DEFAULT
Class representing a probability (a double between 0 and 1)
bool ProbabilitiesSet(vector< Probability > &values)
Set all probabilities.
virtual DiscretePdf * Clone() const
Clone function.
bool ProbabilitySet(int state, Probability a)
Function to change/set the probability of a single state.
unsigned int _num_states
The number of discrete state.
Definition: discretepdf.h:38
void ValueSet(const T &value)
Set the value of the Sample.
bool CumPDFUpdate()
Updates the cumPDF.
vector< Probability > * _Values_p
Pointer to the discrete PDF-values, the sum of the elements = 1.
Definition: discretepdf.h:41
Class PDF: Virtual Base class representing Probability Density Functions.
int MostProbableStateGet()
Get the index of the most probable state.
bool SampleFrom(vector< Sample< int > > &list_samples, const unsigned int num_samples, int method=DEFAULT, void *args=NULL) const
vector< double > _CumPDF
STL-vector containing the Cumulative PDF (for efficient sampling)
Definition: discretepdf.h:47
#define RIPLEY
Class representing a PDF on a discrete variable.
virtual bool SampleFrom(vector< Sample< T > > &list_samples, const unsigned int num_samples, int method=DEFAULT, void *args=NULL) const
Draw multiple samples from the Pdf (overloaded)
Probability ProbabilityGet(const int &state) const
Implementation of virtual base class method.
unsigned int NumStatesGet() const
Get the number of discrete States.
double runif()
vector< Probability > ProbabilitiesGet() const
Get all probabilities.
#define NUMERIC_PRECISION
bool NormalizeProbs()
Normalize all the probabilities (eg. after setting a probability)
DiscretePdf(unsigned int num_states=0)
Constructor (dimension = number of classes) An equal probability is set for all classes.
vector< Probability > ProbabilitiesGet() const
Get all probabilities.


bfl
Author(s): Klaas Gadeyne, Wim Meeussen, Tinne Delaet and many others. See web page for a full contributor list. ROS package maintained by Wim Meeussen.
autogenerated on Mon Feb 28 2022 21:56:33