Line data    Source code

       1             : /* Copyright (c) 2003, Roger Dingledine
       2             :  * Copyright (c) 2004-2006, Roger Dingledine, Nick Mathewson.
       3             :  * Copyright (c) 2007-2021, The Tor Project, Inc. */
       4             : /* See LICENSE for licensing information */
       5             : 
       6             : /**
       7             :  * \file laplace.c
       8             :  *
       9             :  * \brief Implements a Laplace distribution, used for adding noise to things.
      10             :  **/
      11             : 
      12             : #include "orconfig.h"
      13             : #include "lib/math/laplace.h"
      14             : #include "lib/math/fp.h"
      15             : 
      16             : #include "lib/log/util_bug.h"
      17             : 
      18             : #include <math.h>
      19             : #include <stdlib.h>
      20             : 
      21             : /** Transform a random value <b>p</b> from the uniform distribution in
      22             :  * [0.0, 1.0[ into a Laplace distributed value with location parameter
      23             :  * <b>mu</b> and scale parameter <b>b</b>. Truncate the final result
      24             :  * to be an integer in [INT64_MIN, INT64_MAX]. */
      25             : int64_t
      26          36 : sample_laplace_distribution(double mu, double b, double p)
      27             : {
      28          36 :   double result;
      29          36 :   tor_assert(p >= 0.0 && p < 1.0);
      30             : 
      31             :   /* This is the "inverse cumulative distribution function" from:
      32             :    * https://en.wikipedia.org/wiki/Laplace_distribution */
      33          36 :   if (p <= 0.0) {
      34             :     /* Avoid taking log(0.0) == -INFINITY, as some processors or compiler
      35             :      * options can cause the program to trap. */
      36             :     return INT64_MIN;
      37             :   }
      38             : 
      39          27 :   result = mu - b * (p > 0.5 ? 1.0 : -1.0)
      40          27 :                   * tor_mathlog(1.0 - 2.0 * fabs(p - 0.5));
      41             : 
      42          27 :   return clamp_double_to_int64(result);
      43             : }
      44             : 
      45             : /** Add random noise between INT64_MIN and INT64_MAX coming from a Laplace
      46             :  * distribution with mu = 0 and b = <b>delta_f</b>/<b>epsilon</b> to
      47             :  * <b>signal</b> based on the provided <b>random</b> value in [0.0, 1.0[.
      48             :  * The epsilon value must be between ]0.0, 1.0]. delta_f must be greater
      49             :  * than 0. */
      50             : int64_t
      51          31 : add_laplace_noise(int64_t signal_, double random_, double delta_f,
      52             :                   double epsilon)
      53             : {
      54          31 :   int64_t noise;
      55             : 
      56             :   /* epsilon MUST be between ]0.0, 1.0] */
      57          31 :   tor_assert(epsilon > 0.0 && epsilon <= 1.0);
      58             :   /* delta_f MUST be greater than 0. */
      59          31 :   tor_assert(delta_f > 0.0);
      60             : 
      61             :   /* Just add noise, no further signal */
      62          31 :   noise = sample_laplace_distribution(0.0,
      63             :                                       delta_f / epsilon,
      64             :                                       random_);
      65             : 
      66             :   /* Clip (signal + noise) to [INT64_MIN, INT64_MAX] */
      67          31 :   if (noise > 0 && INT64_MAX - noise < signal_)
      68             :     return INT64_MAX;
      69          27 :   else if (noise < 0 && INT64_MIN - noise > signal_)
      70             :     return INT64_MIN;
      71             :   else
      72          22 :     return signal_ + noise;
      73             : }