00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027 #include "CeylanRandomGeneratorFromPDF.h"
00028
00029
00030 #include "CeylanOperators.h"
00031 #include "CeylanWhiteNoiseRandomGenerator.h"
00032 #include "CeylanStringUtils.h"
00033 #include "CeylanLogPlug.h"
00034
00035 #ifdef CEYLAN_USES_CONFIG_H
00036 #include "CeylanConfig.h"
00037 #endif // CEYLAN_USES_CONFIG_H
00038
00039
00040 #include <cstdlib>
00041 #include <list>
00042
00043
00044
00045 using std::string ;
00046 using std::list ;
00047
00048 using namespace Ceylan ;
00049 using namespace Ceylan::Log ;
00050
00051 using namespace Ceylan::Maths ;
00052 using namespace Ceylan::Maths::Random ;
00053
00054
00055
00056
00057
00058 RandomGeneratorFromPDF::RandomGeneratorFromPDF(
00059 const ProbabilityFunction & pdf,
00060 Sample lowerLimit, Sample upperLimit, Seed aSeed ) :
00061 RandomGenerator( lowerLimit, upperLimit, aSeed ),
00062 _pdf( & pdf ),
00063 _whiteNoiseGenerator( 0 ),
00064 _probabilitiesTable( 0 )
00065 {
00066
00067
00068
00069 preCompute() ;
00070
00071 }
00072
00073
00074
00075 RandomGeneratorFromPDF::~RandomGeneratorFromPDF() throw()
00076 {
00077
00078
00079
00080 if ( _whiteNoiseGenerator != 0 )
00081 delete _whiteNoiseGenerator ;
00082
00083 if ( _probabilitiesTable != 0 )
00084 delete [] _probabilitiesTable ;
00085
00086 if ( _sampleRangesTable != 0 )
00087 delete [] _sampleRangesTable ;
00088
00089 }
00090
00091
00092
00093 RandomValue RandomGeneratorFromPDF::getNewValue()
00094 {
00095
00096 #if CEYLAN_DEBUG
00097
00098 if ( _whiteNoiseGenerator == 0 )
00099 throw MathsException( "RandomGeneratorFromPDF::getNewValue: "
00100 "no white noise generator available." ) ;
00101
00102 #endif // CEYLAN_DEBUG
00103
00104 Sample uniformRandomValue = _whiteNoiseGenerator->getNewValue() ;
00105
00106
00107
00108
00109
00110
00111
00112 Ceylan::Uint32 count = 0 ;
00113
00114
00115
00116
00117
00118
00119 while ( uniformRandomValue > _sampleRangesTable[ count ] )
00120 count++ ;
00121
00122 RandomValue res = count + _lowerLimit ;
00123
00124 #if CEYLAN_DEBUG
00125
00126 if ( static_cast<Sample>( res ) >= _upperLimit )
00127 {
00128
00129 LogPlug::error( "RandomGeneratorFromPDF::getNewValue: value "
00130 + Ceylan::toString( res )
00131 + " was drawn, whereas upper excluded limit is "
00132 + Ceylan::toString( _upperLimit ) ) ;
00133
00134 try
00135 {
00136
00137 displayProbabilities() ;
00138
00139 }
00140 catch ( const Ceylan::Exception & e )
00141 {
00142 LogPlug::error( "RandomGeneratorFromPDF::getNewValue: "
00143 "additional error: " + e.toString() ) ;
00144 }
00145
00146
00147 throw MathsException( "RandomGeneratorFromPDF::getNewValue: "
00148 "incorrect random value returned: "
00149 + Ceylan::toString( res ) + "." ) ;
00150
00151 }
00152 #endif // CEYLAN_DEBUG
00153
00154 return res ;
00155
00156 }
00157
00158
00159
00160 void RandomGeneratorFromPDF::reset( Seed newSeed )
00161 {
00162
00163 if ( _whiteNoiseGenerator != 0 )
00164 _whiteNoiseGenerator->reset( newSeed ) ;
00165 else
00166 LogPlug::warning( "RandomGeneratorFromPDF::reset: "
00167 "no white noise generator available, none reset." ) ;
00168
00169 }
00170
00171
00172
00173 const string RandomGeneratorFromPDF::displayProbabilities() const
00174 {
00175
00176
00177 if ( _probabilitiesTable == 0 )
00178 throw MathsException(
00179 "RandomGeneratorFromPDF::displayProbabilities: "
00180 "no probability table available." ) ;
00181
00182 list<string> l ;
00183
00184 Probability sum = 0 ;
00185 Probability currentProbability ;
00186
00187
00188
00189 Probability minProbability = 1000 ;
00190 Sample sampleMin = 0 ;
00191
00192 Probability maxProbability = 0 ;
00193 Sample sampleMax = 0 ;
00194
00195 string result = "Probability table is: " ;
00196
00197 for ( Sample sample = _lowerLimit; sample < _upperLimit; sample++ )
00198 {
00199
00200 currentProbability = _probabilitiesTable[ sample - _lowerLimit ] ;
00201
00202 l.push_back( "Probability to realize sample "
00203 + Ceylan::toString( sample )
00204 + ": " + Ceylan::toString( 100 * currentProbability )
00205 + " %." ) ;
00206
00207 sum += _probabilitiesTable[ sample - _lowerLimit ] ;
00208
00209 if ( currentProbability > maxProbability )
00210 {
00211 maxProbability = currentProbability ;
00212 sampleMax = sample ;
00213 }
00214
00215 if ( currentProbability < minProbability )
00216 {
00217 minProbability = currentProbability ;
00218 sampleMin = sample ;
00219 }
00220
00221
00222 }
00223
00224
00225
00226 #if CEYLAN_DEBUG
00227
00228 LogPlug::debug(
00229 "Displaying corresponding look-up sample range table:" ) ;
00230
00231 for ( Sample sample = _lowerLimit; sample < _upperLimit; sample++ )
00232 {
00233 LogPlug::debug( "Sample " + Ceylan::toString( sample )
00234 + " stops at " + _sampleRangesTable[ sample - _lowerLimit ] ) ;
00235 }
00236
00237 LogPlug::debug( "Upper bound to White noise random values is "
00238 + Ceylan::toString( RAND_MAX ) ) ;
00239
00240 #endif // CEYLAN_DEBUG
00241
00242
00243 return result + Ceylan::formatStringList( l )
00244 + "Probabilities sum up to "
00245 + Ceylan::toString( 100 * sum )
00246 + "%. Maximum probability ("
00247 + Ceylan::toString( 100 * maxProbability )
00248 + "%) reached for sample "
00249 + Ceylan::toString( sampleMax )
00250 + ", minimum probability ("
00251 + Ceylan::toString( 100 * minProbability )
00252 + "%) reached for sample "
00253 + Ceylan::toString( sampleMin ) + "." ;
00254
00255 }
00256
00257
00258
00259 const string RandomGeneratorFromPDF::toString( VerbosityLevels level ) const
00260 {
00261
00262 string base = "Probability Density Function-based random generator, " ;
00263
00264 if ( _pdf != 0 )
00265 {
00266 return base + "whose PDF is " + _pdf->toString() + ". "
00267 + RandomGenerator::toString( level) ;
00268 }
00269 else
00270 {
00271 return base + "with no PDF registered. "
00272 + RandomGenerator::toString( level ) ;
00273 }
00274
00275
00276 return base ;
00277
00278 }
00279
00280
00281
00282 void RandomGeneratorFromPDF::preCompute()
00283 {
00284
00285
00286
00287
00288
00289
00290
00291
00292 _whiteNoiseGenerator = new WhiteNoiseGenerator( 0, RAND_MAX ) ;
00293
00294 if ( _whiteNoiseGenerator == 0 )
00295 throw MathsException( "RandomGeneratorFromPDF::preCompute: "
00296 "allocation of internal white noise generator failed." ) ;
00297
00298 if ( _pdf == 0 )
00299 throw MathsException( "RandomGeneratorFromPDF::preCompute: "
00300 "no probability density function available" ) ;
00301
00302
00303
00304
00305
00306
00307
00308
00309
00310 _probabilitiesTable = new Probability[ _upperLimit - _lowerLimit ] ;
00311
00312 if ( _probabilitiesTable == 0 )
00313 throw MathsException( "RandomGeneratorFromPDF::preCompute: "
00314 "not enough memory to allocate probability table." ) ;
00315
00316 Probability probaTemp ;
00317 Probability toNormalize = 0 ;
00318
00319 for ( Sample sample = _lowerLimit; sample < _upperLimit; sample++ )
00320 {
00321 probaTemp = (*_pdf)( sample ) ;
00322 _probabilitiesTable[ sample - _lowerLimit ] = probaTemp ;
00323 toNormalize += probaTemp ;
00324 }
00325
00326
00327
00328
00329
00330
00331
00332 for ( Sample sample = _lowerLimit; sample < _upperLimit; sample ++ )
00333 _probabilitiesTable[ sample - _lowerLimit ] /= toNormalize ;
00334
00335
00336
00337
00338 _sampleRangesTable = new Sample[ _upperLimit - _lowerLimit ] ;
00339
00340 if ( _sampleRangesTable == 0 )
00341 throw MathsException( "RandomGeneratorFromPDF::preCompute: "
00342 "not enough memory to allocate probability table." ) ;
00343
00344 Sample sampleIndex = 0 ;
00345
00346 for ( Sample sample = _lowerLimit; sample < _upperLimit; sample++ )
00347 {
00348 sampleIndex += static_cast<Sample>(
00349 _probabilitiesTable[ sample - _lowerLimit ] * RAND_MAX ) ;
00350
00351 _sampleRangesTable[ sample - _lowerLimit ] = sampleIndex ;
00352 }
00353
00354
00355
00356
00357
00358
00359
00360 _sampleRangesTable[ _upperLimit - _lowerLimit - 1 ] = RAND_MAX ;
00361
00362
00363
00364
00365
00366
00367
00368
00369 }
00370