diff --git a/sw/lib/ocaml/ocaml_tools.ml b/sw/lib/ocaml/ocaml_tools.ml index b32435ac91..66cdbfdf0f 100644 --- a/sw/lib/ocaml/ocaml_tools.ml +++ b/sw/lib/ocaml/ocaml_tools.ml @@ -24,7 +24,7 @@ * *) - +let pi = 3.14159265358979323846;; let open_compress file = if Filename.check_suffix file "gz" or Filename.check_suffix file "Z" or @@ -54,3 +54,16 @@ let affine_transform = fun format -> [a;b] -> float_of_string a, float_of_string b | [a] -> float_of_string a, 0. | _ -> 1., 0. + +(* Box-Muller transform to generate a normal distribution from a uniform one + http://en.wikipedia.org/wiki/Normal_distribution *) +let normal = fun mu sigma -> + let u1 = Random.float 1. + and u2 = Random.float 1. in + mu +. sigma *. sqrt (-2. *. log u1) *. cos (2. *. pi *. u2) + +let make_1st_order_noise_generator = fun ?(init = 0.) k sigma -> + let x = ref init in + fun () -> + x := k *. !x +. normal 0. sigma; + !x diff --git a/sw/lib/ocaml/ocaml_tools.mli b/sw/lib/ocaml/ocaml_tools.mli index d8ab9e9677..3b945bf46d 100644 --- a/sw/lib/ocaml/ocaml_tools.mli +++ b/sw/lib/ocaml/ocaml_tools.mli @@ -37,3 +37,12 @@ val find_file : string list -> string -> string val affine_transform : string -> float * float (* [affine_transform format] Parses [format] as a+b and returns (a,b). Returns a if +b is not present. Returns 1. if the parsing fails *) + +val normal : float -> float -> float +(* [normal mean stdev] Normal distribution *) + +val make_1st_order_noise_generator : + ?init:float -> float -> float -> (unit -> float) +(* [make_1st_order_noise_generator ?init k sigma] Returns a generator + initialized to [init], damped by [k] with a variation of stdev [sigma]. + 0 < [k] < 1 . x_n+1 <- k x_n + normal 0 sigma . *)