c/*@@
c  @file      AHFinder_mask.F
c  @date      November 1998
c  @author    Miguel Alcubierre
c  @desc 
c             Set up horizon mask for black hole excision.
c  @enddesc 
c@@*/

c Note that including cactus.h will also include AHFinder.h
!#include "cactus.h"
#include "cctk.h"
#include "cctk_parameters.h"
#include "cctk_arguments.h"

      subroutine AHFinder_mask(CCTK_FARGUMENTS)

!     Set up the horizon mask.  All this routine does is
!     set up the value of the grid function "ahmask" to
!     one for points outside the horizon, and to zero
!     for points inside the horizon.  This is only necessary
!     if we want to use black hole excision.

      use AHFinder_dat

      implicit none

      DECLARE_CCTK_FARGUMENTS
      DECLARE_CCTK_PARAMETERS

      integer i,j,k

      CCTK_REAL zero,one
      CCTK_REAL shrink


!     *****************************
!     ***   DEFINE {zero,one}   ***
!     *****************************

      zero = 0.0D0
      one  = 1.0D0


!     *************************
!     ***   START ROUTINE   ***
!     *************************

      call AHFinder_fun(CCTK_FARGUMENTS)

!     Find out how big is the buffer zone.  Notice that
!     a negative value for shrink puts the masked region
!     safely inside the horizon, while a positive value
!     moves the masked region outside the horizon.

!      shrink = getr8("ahf_maskshrink")*dx
      shrink = ahf_maskshrink*dx

!     Mask points.

      do k=1,nz
         do j=1,ny
            do i=1,nx

               if (ahfgrid(i,j,k).ge.shrink) then
                  ahmask(i,j,k) = one
               else
                  ahmask(i,j,k) = zero
               end if

            end do
         end do
      end do


!     ***************
!     ***   END   ***
!     ***************

      end subroutine AHFinder_mask