wrfv2_fire/var/da/da_recursive_filter/da_recursive_filter_1d.inc

   1 subroutine da_recursive_filter_1d(pass, alpha, field, n)
   2
   3    !---------------------------------------------------------------------------
   4    ! Purpose: Perform one pass of recursive filter on 1D array.
   5    !
   6    ! Method:  Perform right-moving filter followed by left-moving filter.
   7    !---------------------------------------------------------------------------
   8
   9    implicit none
  10
  11    integer, intent(in)    :: pass           ! Current pass of filter.
  12    real,    intent(in)    :: alpha          ! Alpha coefficient for RF.
  13    real,    intent(inout) :: field(:)       ! Array to be filtered.
  14    integer, intent(in)    :: n              ! Size of field array.
  15
  16    integer :: j              ! Loop counter.
  17    real    :: one_alpha      ! 1 - alpha.
  18    real    :: a(1:n)         ! Input field.
  19    real    :: b(1:n)         ! Field after left-right pass.
  20    real    :: c(1:n)         ! Field after right-left pass.
  21
  22    if (trace_use_dull) call da_trace_entry("da_recursive_filter_1d")
  23
  24    !-------------------------------------------------------------------------
  25    ! [1.0] Initialise:
  26    !-------------------------------------------------------------------------
  27
  28    one_alpha = 1.0 - alpha
  29
  30    a(1:n) = field(1:n)
  31
  32    !-------------------------------------------------------------------------
  33    ! [2.0] Perform right-moving filter:
  34    !-------------------------------------------------------------------------
  35
  36    ! use turning conditions as in the appendix of Hayden & Purser (1995):
  37
  38    if (pass == 1) then
  39       b(1) = one_alpha * a(1)
  40    else if (pass == 2) then
  41       b(1) = a(1) / (1.0 + alpha)
  42    else
  43       b(1) = one_alpha * (a(1) - alpha**3 * a(2)) / (1.0 - alpha**2)**2
  44    end if
  45
  46    ! [2.2] Perform pass left to right:
  47
  48    do j = 2, n
  49       b(j) = alpha * b(j-1) + one_alpha * a(j)
  50    end do
  51
  52    !-------------------------------------------------------------------------
  53    ! [3.0] Perform left-moving filter:
  54    !-------------------------------------------------------------------------
  55
  56    ! use turning conditions as in the appendix of Hayden & Purser (1995):
  57
  58    if (pass == 1) then
  59       c(n) = b(n) / (1.0 + alpha)
  60    else
  61       c(n) = one_alpha * (b(n) - alpha**3 * b(n-1)) / (1.0 - alpha**2)**2
  62    end if
  63
  64    ! [3.2] Perform pass left to right:
  65
  66    do j = n-1, 1, -1
  67       c(j) = alpha * c(j+1) + one_alpha * b(j)
  68    end do
  69
  70    field(1:n) = c(1:n)
  71
  72    if (trace_use_dull) call da_trace_exit("da_recursive_filter_1d")
  73
  74 end subroutine da_recursive_filter_1d
  75
  76