| |
| /////////////////////////////////////////////////////////////////////////// |
| // |
| // 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. |
| // |
| /////////////////////////////////////////////////////////////////////////// |
| |
| //----------------------------------------------------------------------------- |
| // |
| // Routines that generate pseudo-random numbers compatible |
| // with the standard erand48(), nrand48(), etc. functions. |
| // |
| //----------------------------------------------------------------------------- |
| |
| #include "ImathRandom.h" |
| #include "ImathInt64.h" |
| |
| namespace Imath { |
| namespace { |
| |
| // |
| // Static state used by Imath::drand48(), Imath::lrand48() and Imath::srand48() |
| // |
| |
| unsigned short staticState[3] = {0, 0, 0}; |
| |
| |
| void |
| rand48Next (unsigned short state[3]) |
| { |
| // |
| // drand48() and friends are all based on a linear congruential |
| // sequence, |
| // |
| // x[n+1] = (a * x[n] + c) % m, |
| // |
| // where a and c are as specified below, and m == (1 << 48) |
| // |
| |
| static const Int64 a = Int64 (0x5deece66dLL); |
| static const Int64 c = Int64 (0xbLL); |
| |
| // |
| // Assemble the 48-bit value x[n] from the |
| // three 16-bit values stored in state. |
| // |
| |
| Int64 x = (Int64 (state[2]) << 32) | |
| (Int64 (state[1]) << 16) | |
| Int64 (state[0]); |
| |
| // |
| // Compute x[n+1], except for the "modulo m" part. |
| // |
| |
| x = a * x + c; |
| |
| // |
| // Disassemble the 48 least significant bits of x[n+1] into |
| // three 16-bit values. Discard the 16 most significant bits; |
| // this takes care of the "modulo m" operation. |
| // |
| // We assume that sizeof (unsigned short) == 2. |
| // |
| |
| state[2] = x >> 32; |
| state[1] = x >> 16; |
| state[0] = x; |
| } |
| |
| } // namespace |
| |
| |
| double |
| erand48 (unsigned short state[3]) |
| { |
| // |
| // Generate double-precision floating-point values between 0.0 and 1.0: |
| // |
| // The exponent is set to 0x3ff, which indicates a value greater |
| // than or equal to 1.0, and less than 2.0. The 48 most significant |
| // bits of the significand (mantissa) are filled with pseudo-random |
| // bits generated by rand48Next(). The remaining 4 bits are a copy |
| // of the 4 most significant bits of the significand. This results |
| // in bit patterns between 0x3ff0000000000000 and 0x3fffffffffffffff, |
| // which correspond to uniformly distributed floating-point values |
| // between 1.0 and 1.99999999999999978. Subtracting 1.0 from those |
| // values produces numbers between 0.0 and 0.99999999999999978, that |
| // is, between 0.0 and 1.0-DBL_EPSILON. |
| // |
| |
| rand48Next (state); |
| |
| union {double d; Int64 i;} u; |
| |
| u.i = (Int64 (0x3ff) << 52) | // sign and exponent |
| (Int64 (state[2]) << 36) | // significand |
| (Int64 (state[1]) << 20) | |
| (Int64 (state[0]) << 4) | |
| (Int64 (state[2]) >> 12); |
| |
| return u.d - 1; |
| } |
| |
| |
| double |
| drand48 () |
| { |
| return Imath::erand48 (staticState); |
| } |
| |
| |
| long int |
| nrand48 (unsigned short state[3]) |
| { |
| // |
| // Generate uniformly distributed integers between 0 and 0x7fffffff. |
| // |
| |
| rand48Next (state); |
| |
| return ((long int) (state[2]) << 15) | |
| ((long int) (state[1]) >> 1); |
| } |
| |
| |
| long int |
| lrand48 () |
| { |
| return Imath::nrand48 (staticState); |
| } |
| |
| |
| void |
| srand48 (long int seed) |
| { |
| staticState[2] = seed >> 16; |
| staticState[1] = seed; |
| staticState[0] = 0x330e; |
| } |
| |
| |
| float |
| Rand32::nextf () |
| { |
| // |
| // Generate single-precision floating-point values between 0.0 and 1.0: |
| // |
| // The exponent is set to 0x7f, which indicates a value greater than |
| // or equal to 1.0, and less than 2.0. The 23 bits of the significand |
| // (mantissa) are filled with pseudo-random bits generated by |
| // Rand32::next(). This results in in bit patterns between 0x3f800000 |
| // and 0x3fffffff, which correspond to uniformly distributed floating- |
| // point values between 1.0 and 1.99999988. Subtracting 1.0 from |
| // those values produces numbers between 0.0 and 0.99999988, that is, |
| // between 0.0 and 1.0-FLT_EPSILON. |
| // |
| |
| next (); |
| |
| union {float f; unsigned int i;} u; |
| |
| u.i = 0x3f800000 | (_state & 0x7fffff); |
| return u.f - 1; |
| } |
| |
| } // namespace Imath |