var/da/da_define_structures/da_gauss_noise.inc

   1 subroutine da_gauss_noise( z)
   2
   3    !-----------------------------------------------------------------------
   4    ! Purpose: TBD
   5    !-----------------------------------------------------------------------
   6
   7    implicit none
   8
   9    real,intent(out)     :: z
  10    real                 :: x, y, r, coeff
  11
  12    if (trace_use) call da_trace_entry("da_gauss_noise")
  13
  14    ! [2.1] Get two uniform variate random numbers in range 0 to 1:
  15
  16    do
  17       call random_number( x)
  18       call random_number( y)
  19
  20       ! [2.2] Transform to range -1 to 1 and calculate sum of squares:
  21
  22       x = 2.0 * x - 1.0
  23       y = 2.0 * y - 1.0
  24       r = x * x + y * y
  25
  26       if (r > 0.0 .and. r < 1.0) exit
  27    end do
  28
  29    ! [2.3] use Box-Muller transformation to get normal deviates:
  30
  31    coeff = sqrt( -2.0 * log(r) / r)
  32    z = coeff * x
  33
  34    if (trace_use) call da_trace_exit("da_gauss_noise")
  35
  36 end subroutine da_gauss_noise
  37
  38