/////////////////////////////////////////////////////////////////////////// // // Copyright (c) 2002, Industrial Light & Magic, a division of Lucas // Digital Ltd. LLC // // All rights reserved. // // Redistribution and use in source and binary forms, with or without // modification, are permitted provided that the following conditions are // met: // * Redistributions of source code must retain the above copyright // notice, this list of conditions and the following disclaimer. // * Redistributions in binary form must reproduce the above // copyright notice, this list of conditions and the following disclaimer // in the documentation and/or other materials provided with the // distribution. // * Neither the name of Industrial Light & Magic nor the names of // its contributors may be used to endorse or promote products derived // from this software without specific prior written permission. // // THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS // "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT // LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR // A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT // OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, // SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT // LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, // DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY // THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT // (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE // OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. // /////////////////////////////////////////////////////////////////////////// #ifndef INCLUDED_IMATHRANDOM_H #define INCLUDED_IMATHRANDOM_H //----------------------------------------------------------------------------- // // Generators for uniformly distributed pseudo-random numbers and // functions that use those generators to generate numbers with // different distributions: // // class Rand32 // class Rand48 // solidSphereRand() // hollowSphereRand() // gaussRand() // gaussSphereRand() // //----------------------------------------------------------------------------- // // Here is the copyright for the *rand48() functions implemented for // PLATFORM_DARWIN_PPC and PLATFORM_WINDOWS // // // Copyright (c) 1993 Martin Birgmeier // All rights reserved. // // You may redistribute unmodified or modified versions of this source // code provided that the above copyright notice and this and the // following conditions are retained. // // This software is provided ``as is'', and comes with no warranties // of any kind. I shall in no event be liable for anything that happens // to anyone/anything when using this software. // #include #include namespace Imath { //----------------------------------------------- // Fast random-number generator that generates // a uniformly distributed sequence with a period // length of 2^32. //----------------------------------------------- class Rand32 { public: //------------ // Constructor //------------ Rand32 (unsigned long int seed = 0); //-------------------------------- // Re-initialize with a given seed //-------------------------------- void init (unsigned long int seed); //---------------------------------------------------------- // Get the next value in the sequence (range: [false, true]) //---------------------------------------------------------- bool nextb (); //--------------------------------------------------------------- // Get the next value in the sequence (range: [0 ... 0xffffffff]) //--------------------------------------------------------------- unsigned long int nexti (); //------------------------------------------------------ // Get the next value in the sequence (range: [0 ... 1[) //------------------------------------------------------ float nextf (); //------------------------------------------------------------------- // Get the next value in the sequence (range [rangeMin ... rangeMax[) //------------------------------------------------------------------- float nextf (float rangeMin, float rangeMax); private: void next (); unsigned long int _state; }; //-------------------------------------------------------- // Random-number generator based on the C Standard Library // functions drand48(), lrand48() & company; generates a // uniformly distributed sequence. //-------------------------------------------------------- class Rand48 { public: //------------ // Constructor //------------ Rand48 (unsigned long int seed = 0); //-------------------------------- // Re-initialize with a given seed //-------------------------------- void init (unsigned long int seed); //---------------------------------------------------------- // Get the next value in the sequence (range: [false, true]) //---------------------------------------------------------- bool nextb (); //--------------------------------------------------------------- // Get the next value in the sequence (range: [0 ... 0x7fffffff]) //--------------------------------------------------------------- long int nexti (); //------------------------------------------------------ // Get the next value in the sequence (range: [0 ... 1[) //------------------------------------------------------ double nextf (); //------------------------------------------------------------------- // Get the next value in the sequence (range [rangeMin ... rangeMax[) //------------------------------------------------------------------- double nextf (double rangeMin, double rangeMax); private: unsigned short int _state[3]; #if defined( PLATFORM_DARWIN_PPC ) || defined ( PLATFORM_WINDOWS ) void shiftState(); #endif }; //------------------------------------------------------------ // Return random points uniformly distributed in a sphere with // radius 1 around the origin (distance from origin <= 1). //------------------------------------------------------------ template Vec solidSphereRand (Rand &rand); //------------------------------------------------------------- // Return random points uniformly distributed on the surface of // a sphere with radius 1 around the origin. //------------------------------------------------------------- template Vec hollowSphereRand (Rand &rand); //----------------------------------------------- // Return random numbers with a normal (Gaussian) // distribution with zero mean and unit variance. //----------------------------------------------- template float gaussRand (Rand &rand); //---------------------------------------------------- // Return random points whose distance from the origin // has a normal (Gaussian) distribution with zero mean // and unit variance. //---------------------------------------------------- template Vec gaussSphereRand (Rand &rand); //--------------- // Implementation //--------------- inline void Rand32::init (unsigned long int seed) { _state = (seed * 0xa5a573a5L) ^ 0x5a5a5a5aL; } inline Rand32::Rand32 (unsigned long int seed) { init (seed); } inline void Rand32::next () { _state = 1664525L * _state + 1013904223L; } inline bool Rand32::nextb () { next (); // Return the 31st (most significant) bit, by and-ing with 2 ^ 31. return !!(_state & 2147483648UL); } inline unsigned long int Rand32::nexti () { next (); return _state & 0xffffffff; } inline float Rand32::nextf () { next (); return ((int) (_state & 0xffffff)) * ((float) (1.0F / 0x1000000)); } inline float Rand32::nextf (float rangeMin, float rangeMax) { return rangeMin + nextf() * (rangeMax - rangeMin); } inline void Rand48::init (unsigned long int seed) { seed = (seed * 0xa5a573a5L) ^ 0x5a5a5a5aL; _state[0] = (unsigned short int) (seed); _state[1] = (unsigned short int) (seed >> 16); _state[2] = (unsigned short int) (seed); } inline Rand48::Rand48 (unsigned long int seed) { init (seed); } #if defined( PLATFORM_DARWIN_PPC ) || defined ( PLATFORM_WINDOWS ) inline void Rand48::shiftState() { unsigned long accu; unsigned short temp[2]; accu = 0xe66dUL * ( unsigned long )_state[0] + 0x000bUL; temp[0] = ( unsigned short )accu; /* lower 16 bits */ accu >>= sizeof( unsigned short ) * 8; accu += 0xe66dUL * ( unsigned long )_state[1] + 0xdeecUL * ( unsigned long )_state[0]; temp[1] = ( unsigned short )accu; /* middle 16 bits */ accu >>= sizeof( unsigned short ) * 8; accu += 0xe66dUL * _state[2] + 0xdeecUL * _state[1] + 0x0005UL * _state[0]; _state[0] = temp[0]; _state[1] = temp[1]; _state[2] = ( unsigned short )accu; } #endif inline bool Rand48::nextb () { #if defined( PLATFORM_DARWIN_PPC ) || defined ( PLATFORM_WINDOWS ) shiftState(); return ( ( long( _state[2] ) << 15 ) + ( long( _state[1] ) >> 1 ) ) & 0x1; #else return nrand48 (_state) & 1; #endif } inline long int Rand48::nexti () { #if defined( PLATFORM_DARWIN_PPC ) || defined ( PLATFORM_WINDOWS ) shiftState(); return ( long( _state[2] ) << 15 ) + ( long( _state[1] ) >> 1 ); #else return nrand48 (_state); #endif } inline double Rand48::nextf () { #if defined( PLATFORM_DARWIN_PPC ) || defined ( PLATFORM_WINDOWS ) shiftState(); return ldexp( double( _state[0] ), -48 ) + ldexp( double( _state[1] ), -32 ) + ldexp( double( _state[2] ), -16 ); #else return erand48 (_state); #endif } inline double Rand48::nextf (double rangeMin, double rangeMax) { return rangeMin + nextf() * (rangeMax - rangeMin); } template Vec solidSphereRand (Rand &rand) { Vec v; do { for (unsigned int i = 0; i < Vec::dimensions(); i++) v[i] = (typename Vec::BaseType) rand.nextf (-1, 1); } while (v.length2() > 1); return v; } template Vec hollowSphereRand (Rand &rand) { Vec v; typename Vec::BaseType length; do { for (unsigned int i = 0; i < Vec::dimensions(); i++) v[i] = (typename Vec::BaseType) rand.nextf (-1, 1); length = v.length(); } while (length > 1 || length == 0); return v / length; } template float gaussRand (Rand &rand) { float x; // Note: to avoid numerical problems with very small float y; // numbers, we make these variables singe-precision float length2; // floats, but later we call the double-precision log() // and sqrt() functions instead of logf() and sqrtf(). do { x = float (rand.nextf (-1, 1)); y = float (rand.nextf (-1, 1)); length2 = x * x + y * y; } while (length2 >= 1 || length2 == 0); return x * sqrt (-2 * log (length2) / length2); } template Vec gaussSphereRand (Rand &rand) { return hollowSphereRand (rand) * gaussRand (rand); } double drand48(); long int lrand48(); } // namespace Imath #endif