calc_imp_acf Subroutine

   subroutine calc_imp_acf(long, larg, tau1, tau2, alpha, ang, vec_acf, vec_xy)
   !================================================================================================
   !<@note Function that returns the theoretical autocorrelation function in an array.<br/>
   !< The autocorrelation function is supposed to be obtained from a real surface which must be periodic
   !< or nearly periodic (because of the use of FFTs).
   !< In addition, the surface is supposed to be 0 mean and normalized (\(\sigma = 1 \)),
   !< therefore *acf* is zero-mean and normalized so that its max value is 1.<br/>
   !<
   !<@endnote
   !------------------------------------------------------------------------------------------------
   implicit none
   integer(kind=I4), intent(in) :: long   !! *surface acf width*
   integer(kind=I4), intent(in) :: larg   !! *surface acf height*
   real   (kind=R8), intent(in) :: tau1   !! *first correlation length*
   real   (kind=R8), intent(in) :: tau2   !! *surface second correlation length*
   real   (kind=R8), intent(in) :: alpha  !! *parameter that controls the expondential decrease*
   real   (kind=R8), intent(in) :: ang    !! *acf ellipsis angle*
   real   (kind=R8), dimension(1:long * larg),      intent(out) :: vec_acf  !! *resulting acf*
   real   (kind=R8), dimension(1:long * larg, 1:2), intent(out) :: vec_xy   !! *points coordinates*

      integer(kind=I4) :: i, j, k, long2, larg2
      real   (kind=R8) :: xi, xj, r

      long2 = long / 2 + 1
      larg2 = larg / 2 + 1

      k = 0
      do j = 1, larg
      do i = 1, long

         k = k + 1

         xi = real(i - long2, kind=R8) / long
         xj = real(j - larg2, kind=R8) / larg

         vec_xy(k, 1) = xi
         vec_xy(k, 2) = xj

         call random_number(r)

         r = 1. + 0.05 * (2 * (0.5 - r))

         vec_acf(k) = r * autocov_impo( xi    = xi,       &  ! IN
                                        xj    = xj,       &  ! IN
                                        tau1  = tau1,     &  ! IN
                                        tau2  = tau2,     &  ! IN
                                        alpha = alpha,    &  ! IN
                                        ang   = ang )        ! IN

      enddo
      enddo

   return
   endsubroutine calc_imp_acf