digi_fil Subroutine

private subroutine digi_fil()

Note

Function that applies the digital filter to the random heights

Arguments

None

Calls

proc~~digi_fil~~CallsGraph proc~digi_fil digi_fil calc_fftw3_real_bwd calc_fftw3_real_bwd proc~digi_fil->calc_fftw3_real_bwd calc_fftw3_real_fwd calc_fftw3_real_fwd proc~digi_fil->calc_fftw3_real_fwd calc_moments calc_moments proc~digi_fil->calc_moments

Called by

proc~~digi_fil~~CalledByGraph proc~digi_fil digi_fil proc~read_job read_job proc~read_job->proc~digi_fil proc~prg_surf prg_surf proc~prg_surf->proc~read_job program~main main program~main->proc~prg_surf

Source Code

   subroutine digi_fil()
   !================================================================================================
   !<@note Function that applies the digital filter to the random heights
   !<
   !<@endnote
   !------------------------------------------------------------------------------------------------
   implicit none

      integer(kind=I4) :: w, h

      complex(kind=R8), dimension(:,:), allocatable :: cmple
      complex(kind=R8), dimension(:,:), allocatable :: ftab         ! FFT(tab_prf)

      type(MOMENT_STAT) :: m_res

      write(SPY,*) 'digi_fil -> 1 - extends then windows PARAM%surf to calculate its FFT ftab'
      write(SPY,*) 'digi_fil -> 2 - multiplies ftab and the digital filter PARAM%fhi, then FFT-1'
      write(SPY,*) 'digi_fil -> 3 - retrieves PARAM%surf by removing the padded extension'

      w = PARAM%width
      h = PARAM%height

      allocate( cmple(1:w, 1:h) )     !

      allocate(  ftab(1:w, 1:h) )     !

      cmple(1:w, 1:h) = cmplx( PARAM%surf(1:w, 1:h), 0, kind = R8 )

      call calc_fftw3_real_fwd( tab_in = PARAM%surf(1:w, 1:h),     &  ! IN
                                tab_ou = cmple(1:w, 1:h),          &  ! OUT
                                  long = w,                        &  ! IN
                                  larg = h )                          ! IN

      ftab(1:w, 1:h) = cmple(1:w, 1:h)

      where( abs(cmple(1:w, 1:h)) > 100 * EPS_R8 )

         ftab(1:w, 1:h) = cmple(1:w, 1:h) / abs( cmple(1:w, 1:h) )

      elsewhere

         ftab(1:w, 1:h) = cmplx( UN, 0._R8, kind=R8 )

      endwhere

      write(SPY,*) 'ftab normalized'

      !--------------------------------------------------------------------
      ! FFT of the filter * FFT of the random heights
      !--------------------------------------------------------------------

      cmple(1:w, 1:h) = PARAM%fhi(1:w, 1:h) * ftab(1:w, 1:h)

      call calc_fftw3_real_bwd( tab_in = cmple(1:w, 1:h),        &  ! IN
                                tab_ou = PARAM%surf(1:w, 1:h),   &  ! OUT
                                  long = w,                      &  ! IN
                                  larg = h  )                       ! IN

      ! signal centré et normé
      call calc_moments(    tab = PARAM%surf(1:w, 1:h),   &  ! IN
                             mx = m_res,                  &  ! OUT
                         nb_mom = 4 )                        ! IN

      PARAM%surf(1:w, 1:h) = ( PARAM%surf(1:w, 1:h) - m_res%mu ) / m_res%si

      write(SPY,*) 'sk ku fin ', m_res%sk, m_res%ku

      deallocate( cmple )
      deallocate( ftab )

   return
   endsubroutine digi_fil