ITK 6.0.0
Insight Toolkit
 
Loading...
Searching...
No Matches
itkMersenneTwisterRandomVariateGenerator.h
Go to the documentation of this file.
1/*=========================================================================
2 *
3 * Copyright NumFOCUS
4 *
5 * Licensed under the Apache License, Version 2.0 (the "License");
6 * you may not use this file except in compliance with the License.
7 * You may obtain a copy of the License at
8 *
9 * https://www.apache.org/licenses/LICENSE-2.0.txt
10 *
11 * Unless required by applicable law or agreed to in writing, software
12 * distributed under the License is distributed on an "AS IS" BASIS,
13 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14 * See the License for the specific language governing permissions and
15 * limitations under the License.
16 *
17 *=========================================================================*/
18#ifndef itkMersenneTwisterRandomVariateGenerator_h
19#define itkMersenneTwisterRandomVariateGenerator_h
20
21#include "itkMacro.h"
22#include "itkObjectFactory.h"
24#include "itkIntTypes.h"
25#include "itkMath.h"
26#include "itkSingletonMacro.h"
27
28#include <atomic>
29#include <mutex>
30#include <climits>
31#include <ctime>
32#include <limits>
33
35{
120
121struct MersenneTwisterGlobals;
122
124{
125public:
131
132 using IntegerType = uint32_t;
133
135 itkOverrideGetNameOfClassMacro(MersenneTwisterRandomVariateGenerator);
136
145 static Pointer
147
156 static Pointer
158
161 static void
163
171 static constexpr IntegerType DefaultSeed = 121212;
172
174 static constexpr IntegerType StateVectorLength = 624;
175
176#if !defined(ITK_FUTURE_LEGACY_REMOVE) && !defined(ITK_WRAPPING_PARSER)
180 ITK_FUTURE_DEPRECATED("ITK 6 discourages using this member function. Please use `SetSeed` instead!")
181 void
182 Initialize(const IntegerType seed = Self::CreateRandomSeed())
183 {
184 this->SetSeed(seed);
185 }
186#endif
187
189 double
191 {
192 return static_cast<double>(GetIntegerVariate()) * (1.0 / double{ std::numeric_limits<uint32_t>::max() });
193 }
194
196 double
198 {
199 return GetVariateWithClosedRange() * n;
200 }
201
203 double
205 {
206 return static_cast<double>(GetIntegerVariate()) / double{ 1ULL << 32ULL };
207 }
208
210 double
212 {
213 return GetVariateWithOpenUpperRange() * n;
214 }
215
217 double
219 {
220 return (static_cast<double>(GetIntegerVariate()) + 0.5) / double{ 1ULL << 32ULL };
221 }
222
224 double
226 {
227 return GetVariateWithOpenRange() * n;
228 }
229
231 IntegerType
233
237
240 double
242
243 /* Access to a normal random number distribution
244 * TODO: Compare with vnl_sample_normal */
245 double
246 GetNormalVariate(const double mean = 0.0, const double variance = 1.0);
247
248 /* Access to a uniform random number distribution in the range [a, b)
249 * TODO: Compare with vnl_sample_uniform */
250 double
251 GetUniformVariate(const double a, const double b);
252
258 double
259 GetVariate() override
260 {
262 }
263
265 double
267 {
268 return GetVariate();
269 }
270
275 void
277
283 GetSeed() const
284 {
285 return m_Seed;
286 }
287
293 static IntegerType
295
296protected:
299 void
300 PrintSelf(std::ostream & os, Indent indent) const override;
301
303 void
305
306 static IntegerType
307 hash(time_t t, clock_t c);
308
309private:
311 itkGetGlobalDeclarationMacro(MersenneTwisterGlobals, PimplGlobals);
312
314 static Pointer
316
318 static IntegerType
320 {
321 return hash(time(nullptr), clock());
322 }
323
325 void
327
328 // Internal state
330
331 // Next value to get from state
333
334 // Number of values left before reload is needed
335 int m_Left{};
336
337 // Seed value
338 std::atomic<IntegerType> m_Seed{};
339
340 // Local lock to enable concurrent access to singleton
341 std::mutex m_InstanceMutex{};
342
343 // Static/Global Variable need to be thread-safely accessed
344
345 static MersenneTwisterGlobals * m_PimplGlobals;
346
347}; // end of class
348
349// Declare inlined functions.... (must be declared in the header)
350
353inline double
355{
356 const IntegerType a = GetIntegerVariate() >> 5;
357 const IntegerType b = GetIntegerVariate() >> 6;
358
359 return (a * 67108864.0 + b) * (1.0 / 9007199254740992.0); // by Isaku
360 // Wada
361}
362
364// TODO: Compare with vnl_sample_normal
365inline double
366MersenneTwisterRandomVariateGenerator::GetNormalVariate(const double mean, const double variance)
367{
368 // Return a real number from a normal (Gaussian) distribution with given
369 // mean and variance by Box-Muller method
370 const double r = std::sqrt(-2.0 * std::log(1.0 - GetVariateWithOpenRange()) * variance);
371 const double phi = 2.0 * itk::Math::pi * GetVariateWithOpenUpperRange();
372
373 return mean + r * std::cos(phi);
374}
375
377// TODO: Compare with vnl_sample_uniform
378inline double
380{
381 const double u = GetVariateWithOpenUpperRange();
382
383 return ((1.0 - u) * a + u * b);
384}
385
386/* Change log from MTRand.h */
387// Change log:
388//
389// v0.1 - First release on 15 May 2000
390// - Based on code by Makoto Matsumoto, Takuji Nishimura, and Shawn Cokus
391// - Translated from C to C++
392// - Made completely ANSI compliant
393// - Designed convenient interface for initialization, seeding, and
394// obtaining numbers in default or user-defined ranges
395// - Added automatic seeding from /dev/urandom or time() and clock()
396// - Provided functions for saving and loading generator state
397//
398// v0.2 - Fixed bug which reloaded generator one step too late
399//
400// v0.3 - Switched to clearer, faster reload() code from Matthew Bellew
401//
402// v0.4 - Removed trailing newline in saved generator format to be consistent
403// with output format of built-in types
404//
405// v0.5 - Improved portability by replacing static const int's with enum's and
406// clarifying return values in seed(); suggested by Eric Heimburg
407// - Removed MAXINT constant; use 0xffffffffUL instead
408//
409// v0.6 - Eliminated seed overflow when uint32 is larger than 32 bits
410// - Changed integer [0,n] generator to give better uniformity
411//
412// v0.7 - Fixed operator precedence ambiguity in reload()
413// - Added access for real numbers in (0,1) and (0,n)
414//
415// v0.8 - Included time.h header to properly support time_t and clock_t
416//
417// v1.0 - Revised seeding to match 26 Jan 2002 update of Nishimura and Matsumoto
418// - Allowed for seeding with arrays of any length
419// - Added access for real numbers in [0,1) with 53-bit resolution
420// - Added access for real numbers from normal (Gaussian) distributions
421// - Increased overall speed by optimizing twist()
422// - Doubled speed of integer [0,n] generation
423// - Fixed out-of-range number generation on 64-bit machines
424// - Improved portability by substituting literal constants for long enum's
425// - Changed license from GNU LGPL to BSD
426} // namespace itk::Statistics
427
428#endif
Control indentation during Print() invocation.
Definition itkIndent.h:50
Implements transparent reference counting.
IntegerType GetIntegerVariate(const IntegerType n)
static Pointer New()
Method for creation through the object factory.
void PrintSelf(std::ostream &os, Indent indent) const override
static IntegerType hash(time_t t, clock_t c)
itkGetGlobalDeclarationMacro(MersenneTwisterGlobals, PimplGlobals)
void SetSeed(const IntegerType seed=Self::CreateRandomSeed())
void InitializeWithoutMutexLocking(const IntegerType seed)
double GetNormalVariate(const double mean=0.0, const double variance=1.0)
static constexpr double pi
Definition itkMath.h:63