1f419f759735f33721a9506230d9444fb3dce5024Martin Josefsson/*
2f419f759735f33721a9506230d9444fb3dce5024Martin Josefsson *  Copyright (c) 2015 The WebRTC project authors. All Rights Reserved.
3f419f759735f33721a9506230d9444fb3dce5024Martin Josefsson *
45d9678ad3eabc34ac40dfe055d7f6a8e44445a5aJan Engelhardt *  Use of this source code is governed by a BSD-style license
5e62f426c7ead7c0025d15860df97426db6509942Patrick McHardy *  that can be found in the LICENSE file in the root of the source
6a2a7f2b531cc582ab6cc3c2b73715ed1d58b9eabJan Engelhardt *  tree. An additional intellectual property rights grant can be found
7f419f759735f33721a9506230d9444fb3dce5024Martin Josefsson *  in the file PATENTS.  All contributing project authors may
88d14aeb8c4c3dc8ce9264b04b97f2e8634c1f381Jan Engelhardt *  be found in the AUTHORS file in the root of the source tree.
98d14aeb8c4c3dc8ce9264b04b97f2e8634c1f381Jan Engelhardt */
108d14aeb8c4c3dc8ce9264b04b97f2e8634c1f381Jan Engelhardt#include "webrtc/base/random.h"
118d14aeb8c4c3dc8ce9264b04b97f2e8634c1f381Jan Engelhardt
123964023f8640b60456373825b326b91badd7a058Jan Engelhardt#include <math.h>
133964023f8640b60456373825b326b91badd7a058Jan Engelhardt
148d14aeb8c4c3dc8ce9264b04b97f2e8634c1f381Jan Engelhardt#include "webrtc/base/checks.h"
158d14aeb8c4c3dc8ce9264b04b97f2e8634c1f381Jan Engelhardt
161d5b63d12984d12c8d87242179855e17657be16dJan Engelhardtnamespace webrtc {
17f419f759735f33721a9506230d9444fb3dce5024Martin Josefsson
18f419f759735f33721a9506230d9444fb3dce5024Martin JosefssonRandom::Random(uint64_t seed) {
198b7c64d6ba156a99008fcd810cba874c73294333Jan Engelhardt  RTC_DCHECK(seed != 0x0ull);
2018f1aff721e19486d87342abb594831b08b1083eHarald Welte  state_ = seed;
21cf655eb194951a93e4e1371747273c12466c1952Harald Welte}
2205e0b01bd1cd4035893c33c7084164bd8fab37c8Harald Welte
2305e0b01bd1cd4035893c33c7084164bd8fab37c8Harald Welteuint32_t Random::Rand(uint32_t t) {
24cf655eb194951a93e4e1371747273c12466c1952Harald Welte  // Casting the output to 32 bits will give an almost uniform number.
25cf655eb194951a93e4e1371747273c12466c1952Harald Welte  // Pr[x=0] = (2^32-1) / (2^64-1)
26ae4b0b3aa70c67f2eff303a3e75834e45c3794a7Eric Leblond  // Pr[x=k] = 2^32 / (2^64-1) for k!=0
27ae4b0b3aa70c67f2eff303a3e75834e45c3794a7Eric Leblond  // Uniform would be Pr[x=k] = 2^32 / 2^64 for all 32-bit integers k.
288b7c64d6ba156a99008fcd810cba874c73294333Jan Engelhardt  uint32_t x = NextOutput();
29f419f759735f33721a9506230d9444fb3dce5024Martin Josefsson  // If x / 2^32 is uniform on [0,1), then x / 2^32 * (t+1) is uniform on
30f419f759735f33721a9506230d9444fb3dce5024Martin Josefsson  // the interval [0,t+1), so the integer part is uniform on [0,t].
318d14aeb8c4c3dc8ce9264b04b97f2e8634c1f381Jan Engelhardt  uint64_t result = x * (static_cast<uint64_t>(t) + 1);
328d14aeb8c4c3dc8ce9264b04b97f2e8634c1f381Jan Engelhardt  result >>= 32;
338d14aeb8c4c3dc8ce9264b04b97f2e8634c1f381Jan Engelhardt  return result;
348d14aeb8c4c3dc8ce9264b04b97f2e8634c1f381Jan Engelhardt}
358d14aeb8c4c3dc8ce9264b04b97f2e8634c1f381Jan Engelhardt
368d14aeb8c4c3dc8ce9264b04b97f2e8634c1f381Jan Engelhardtuint32_t Random::Rand(uint32_t low, uint32_t high) {
37f419f759735f33721a9506230d9444fb3dce5024Martin Josefsson  RTC_DCHECK(low <= high);
38f419f759735f33721a9506230d9444fb3dce5024Martin Josefsson  return Rand(high - low) + low;
39f419f759735f33721a9506230d9444fb3dce5024Martin Josefsson}
40e62f426c7ead7c0025d15860df97426db6509942Patrick McHardy
41f419f759735f33721a9506230d9444fb3dce5024Martin Josefssonint32_t Random::Rand(int32_t low, int32_t high) {
428d14aeb8c4c3dc8ce9264b04b97f2e8634c1f381Jan Engelhardt  RTC_DCHECK(low <= high);
43bd9438420d92c41a5cf20a53b7a18d3ddea4216dJan Engelhardt  // We rely on subtraction (and addition) to be the same for signed and
44f419f759735f33721a9506230d9444fb3dce5024Martin Josefsson  // unsigned numbers in two-complement representation. Thus, although
458d14aeb8c4c3dc8ce9264b04b97f2e8634c1f381Jan Engelhardt  // high - low might be negative as an int, it is the correct difference
468d14aeb8c4c3dc8ce9264b04b97f2e8634c1f381Jan Engelhardt  // when interpreted as an unsigned.
478d14aeb8c4c3dc8ce9264b04b97f2e8634c1f381Jan Engelhardt  return Rand(high - low) + low;
48e62f426c7ead7c0025d15860df97426db6509942Patrick McHardy}
49f419f759735f33721a9506230d9444fb3dce5024Martin Josefsson
5005e0b01bd1cd4035893c33c7084164bd8fab37c8Harald Weltetemplate <>
51f419f759735f33721a9506230d9444fb3dce5024Martin Josefssonfloat Random::Rand<float>() {
52f419f759735f33721a9506230d9444fb3dce5024Martin Josefsson  double result = NextOutput() - 1;
53f419f759735f33721a9506230d9444fb3dce5024Martin Josefsson  result = result / 0xFFFFFFFFFFFFFFFEull;
541e01b0b82f70b0b11dcfbced485dbe7aeac4fb8cJan Engelhardt  return static_cast<float>(result);
55f419f759735f33721a9506230d9444fb3dce5024Martin Josefsson}
561829ed482efbc8b390cc760d012b3a4450494e1aJan Engelhardt
57f419f759735f33721a9506230d9444fb3dce5024Martin Josefssontemplate <>
58f419f759735f33721a9506230d9444fb3dce5024Martin Josefssondouble Random::Rand<double>() {
5905e0b01bd1cd4035893c33c7084164bd8fab37c8Harald Welte  double result = NextOutput() - 1;
6005e0b01bd1cd4035893c33c7084164bd8fab37c8Harald Welte  result = result / 0xFFFFFFFFFFFFFFFEull;
611e01b0b82f70b0b11dcfbced485dbe7aeac4fb8cJan Engelhardt  return result;
6205e0b01bd1cd4035893c33c7084164bd8fab37c8Harald Welte}
631829ed482efbc8b390cc760d012b3a4450494e1aJan Engelhardt
6405e0b01bd1cd4035893c33c7084164bd8fab37c8Harald Weltetemplate <>
6505e0b01bd1cd4035893c33c7084164bd8fab37c8Harald Weltebool Random::Rand<bool>() {
66f419f759735f33721a9506230d9444fb3dce5024Martin Josefsson  return Rand(0, 1) == 1;
6705e0b01bd1cd4035893c33c7084164bd8fab37c8Harald Welte}
6805e0b01bd1cd4035893c33c7084164bd8fab37c8Harald Welte
691829ed482efbc8b390cc760d012b3a4450494e1aJan Engelhardtdouble Random::Gaussian(double mean, double standard_deviation) {
7005e0b01bd1cd4035893c33c7084164bd8fab37c8Harald Welte  // Creating a Normal distribution variable from two independent uniform
718d14aeb8c4c3dc8ce9264b04b97f2e8634c1f381Jan Engelhardt  // variables based on the Box-Muller transform, which is defined on the
72f419f759735f33721a9506230d9444fb3dce5024Martin Josefsson  // interval (0, 1]. Note that we rely on NextOutput to generate integers
73f419f759735f33721a9506230d9444fb3dce5024Martin Josefsson  // in the range [1, 2^64-1]. Normally this behavior is a bit frustrating,
748d14aeb8c4c3dc8ce9264b04b97f2e8634c1f381Jan Engelhardt  // but here it is exactly what we need.
75f419f759735f33721a9506230d9444fb3dce5024Martin Josefsson  const double kPi = 3.14159265358979323846;
768d14aeb8c4c3dc8ce9264b04b97f2e8634c1f381Jan Engelhardt  double u1 = static_cast<double>(NextOutput()) / 0xFFFFFFFFFFFFFFFFull;
77e62f426c7ead7c0025d15860df97426db6509942Patrick McHardy  double u2 = static_cast<double>(NextOutput()) / 0xFFFFFFFFFFFFFFFFull;
78f419f759735f33721a9506230d9444fb3dce5024Martin Josefsson  return mean + standard_deviation * sqrt(-2 * log(u1)) * cos(2 * kPi * u2);
798d14aeb8c4c3dc8ce9264b04b97f2e8634c1f381Jan Engelhardt}
808d14aeb8c4c3dc8ce9264b04b97f2e8634c1f381Jan Engelhardt
818d14aeb8c4c3dc8ce9264b04b97f2e8634c1f381Jan Engelhardtdouble Random::Exponential(double lambda) {
8205e0b01bd1cd4035893c33c7084164bd8fab37c8Harald Welte  double uniform = Rand<double>();
831829ed482efbc8b390cc760d012b3a4450494e1aJan Engelhardt  return -log(uniform) / lambda;
8405e0b01bd1cd4035893c33c7084164bd8fab37c8Harald Welte}
8505e0b01bd1cd4035893c33c7084164bd8fab37c8Harald Welte
8605e0b01bd1cd4035893c33c7084164bd8fab37c8Harald Welte}  // namespace webrtc
878d14aeb8c4c3dc8ce9264b04b97f2e8634c1f381Jan Engelhardt