psmile_srch_nneigh_gauss2_real.F90

Go to the documentation of this file.

!-----------------------------------------------------------------------
! Copyright 2006-2010, NEC Europe Ltd., London, UK.
! All rights reserved. Use is subject to OASIS4 license terms.
!-----------------------------------------------------------------------
!BOP
!
!
! !ROUTINE: PSMILe_srch_nneigh_gauss2_real
!
subroutine psmile_srch_nneigh_gauss2_real (grid_id, nn_srch,    &
                arrays, search_mode, nref_3d,                   &
                sin_values_lon, cos_values_lon,                 &
                sin_values_lat, cos_values_lat,                 &
                grid_valid_shape,                               &
                z_coords, coords_shape,                         &
                neighbors_3d, nloc, num_neigh,                  &
                sin_search, cos_search, z_search,               &
                dist_real, dim1, extra_search, jbeg, jend,      &
                mask_array, mask_shape, mask_available,         &
                tol, ierror)
!
! !USES:
!
  use PRISM_constants
!
  use PSMILe ! , dummy_interface => PSMILe_srch_nneigh_gauss2_real

  implicit none
!
! !INPUT PARAMETERS:
!
  Integer, Intent (In)            :: grid_id
!
!    Info on the component in which the donor cells
!    should be searched.
!
      Type (Extra_search_nn), Intent (In) :: nn_srch
!
!    Info about the size of the coarse grid problem
!
      Type (Extra_search_real)        :: arrays
!
!    Memory for metric information
!
      Integer                         :: nref_3d
!
!    Something like the maximum number of neighbor to be searched for
!
      Integer, Intent (In)            :: search_mode
!
!     Specifies the search mode for nearest neigbours with
!        search_mode = 3 : Full 3d-search
!        search_mode = 2 : Search in 2d-hyperplane (1st 2nd direction lon, lat)
!        search_mode = 1 : Search in 1d-plane
!
      Integer, Intent (In)            :: grid_valid_shape (2, ndim_3d)
!
!     Specifies the valid block shape for the "ndim_3d"-dimensional block

      Real, Intent (In)               :: sin_values_lon (grid_valid_shape(1,1): 
                                                         grid_valid_shape(2,1))
      Real, Intent (In)               :: sin_values_lat (grid_valid_shape(1,1): 
                                                         grid_valid_shape(2,1))
!
!     Sin values of Longitudes and Latitudes (x_coords, y_coords)
!
      Real, Intent (In)               :: cos_values_lon (grid_valid_shape(1,1): 
                                                         grid_valid_shape(2,1))
      Real, Intent (In)               :: cos_values_lat (grid_valid_shape(1,1): 
                                                         grid_valid_shape(2,1))
!    
!     Cos values of Longitudes and Latitudes (x_coords, y_coords)
!
      Integer, Intent (In)            :: coords_shape (2, ndim_3d)

!     Dimension of coordinates (method) x_coords, ...

      Real, Intent (In)               :: z_coords(coords_shape(1,3): 
                                                  coords_shape(2,3))


      Integer, Intent (In)            :: nloc
!
!     Second dimension of neighbors array "neighbors_3d" and
!     Total number of locations to be transferred
!
      Integer, Intent (In)            :: num_neigh
!
!     Number of neighbors to be searched.
!     Last dimension of neighbors array "neighbors_3d" and
!     number of neighbors to be searched.
!
      Integer,          Intent (In)   :: jbeg, jend
!
!     Shape of input arrays
!
      Real, Intent (In)               :: sin_search (jbeg:jend, 2)
!
!     Sin values of the points for which the nearest neighbours
!     should be searched.
!
      Real, Intent (In)               :: cos_search (jbeg:jend, 2)
!
!     Cos values of the points for which the nearest neighbours
!     should be searched.
!
      Real, Intent (In)               :: z_search (jbeg:jend)
!
!     Z values of the points for which the nearest neighbours
!     should be searched.
!
      Integer, Intent (In)            :: dim1 (2)
!
!     Dimensions of "dist_real" in first direction 
!
      Type (Extra_search_info), Intent(InOut) :: extra_search
!
!     Indices of the extra points in entire list of all points to be
!     searched
!
      Integer, Intent (In)            :: mask_shape (2, ndim_3d)
!
!    Dimension of (method) mask array
!
      Logical, Intent (In)            :: mask_array (mask_shape (1,1): 
                                                     mask_shape (2,1), 
                                                     mask_shape (1,2): 
                                                     mask_shape (2,2), 
                                                     mask_shape (1,3): 
                                                     mask_shape (2,3))
!    Mask of the method
!
      Logical, Intent (In)            :: mask_available
!
!    Is mask specified in array "mask_array" ?
!    mask_available = .false. : Mask is not available
!                     .true.  : Mask is     available
!
      Real, Intent (In)                :: tol
!
!    Some tolerance may be used later to evaluate distances
!
! !INPUT/OUTPUT PARAMETERS:
!
      Real, Intent (InOut)             :: dist_real (dim1(1):dim1(2), num_neigh)
!
!    Initial distances and final distances
!
      Integer,           Intent (Out) :: neighbors_3d (ndim_3d, nloc, num_neigh)

!     Indices of neighbor locations (interpolation bases)
!     Index array for points that have been found
!
! !OUTPUT PARAMETERS:
!
      Integer, Intent (Out)            :: ierror

!    Returns the error code of PSMILE_Neigh_nearx_sub_gauss2_real;
!            ierror = 0 : No error
!            ierror > 0 : Severe error!
!
! !DEFINED PARAMETERS:
!
!  real_earth = Earth radius in meter
!
  Real, Parameter                 :: real_earth = 6400000.0
!
  Real, Parameter                 :: acosp1 =  1.0
  Real, Parameter                 :: acosm1 = -1.0

  Integer, Parameter              :: lon = 1
  Integer, Parameter              :: lat = 2
!
! !LOCAL VARIABLES
!
  Real                            :: dist ( (grid_valid_shape(2,1)-grid_valid_shape(1,1)+1) )

  Real                            :: val, fac

  Integer, Pointer                :: indices (:)

  Integer                         :: leni
  Integer                         :: i, n
  Integer                         :: ii, kk
  Integer                         :: ind, jind
  Integer                         :: iloc(1)
#ifdef DEBUGX
  Integer                         :: m1, m2, m3, m4
#endif
!
! !DESCRIPTION:
!
! Subroutine "PSMILe_srch_nneigh_gauss2_real" searches the next "num_neigh"
! nearest neighbours on the method-grid with brut force but evaluating
! the distances to all source points.
!
! !REVISION HISTORY:
!
!   Date      Programmer   Description
! ----------  ----------   -----------
! 05/02/01    R. Redler    created
!
!EOP
!----------------------------------------------------------------------
!
!  Initialization

#ifdef VERBOSE
  print 9990, trim(ch_id)

  call psmile_flushstd
#endif /* VERBOSE */

  ierror = 0
  indices => extra_search%indices

!
! jind = loop index in compact vector "sin_search", "cos_search",
!        "z_search" and "dist_real"
!
! ind  = Corresponding index in "neighbors_3d"
!
  leni = grid_valid_shape(2,1) - grid_valid_shape(1,1) + 1

  do jind = jbeg, jend

     ! distance from source to target point

     fac = real_earth + z_search (jind)

        do i = grid_valid_shape(1,1), grid_valid_shape(2,1)

           val =    sin_values_lat(i) * sin_search(jind, lat)    &
                 +  cos_values_lat(i) * cos_search(jind, lat)    &
                 * (cos_values_lon(i) * cos_search(jind, lon)    &
                 +  sin_values_lon(i) * sin_search(jind, lon))

           val = min (acosp1, val)
           val = max (acosm1, val)

           if ( mask_available ) then
              if ( mask_array(i,1,1) ) then
                 dist ( i-grid_valid_shape(1,1) + 1) = fac * acos (val)
              else
                 dist ( i-grid_valid_shape(1,1) + 1) = huge(fac)
              endif
           else
              dist ( i-grid_valid_shape(1,1) + 1) = fac * acos (val)
           endif

     enddo

     ! To Do: Find closest level in k

     kk = 1

     ind = indices(jind)

     do n = 1, num_neigh

        iloc(:) = minloc(dist)

#ifdef MINLOCFIX
        if (iloc(1) == 0) iloc(1) = 1
#endif /* MINLOCFIX */
        ii =   iloc(1)
        
        neighbors_3d (1, ind, n) = ii
        neighbors_3d (2, ind, n) = 1
        neighbors_3d (3, ind, n) = kk

#ifdef DEBUGX
        print *, ' Latitude '
        print *, asin(sin_search(jind, lat)) / real_deg2rad
        print *, ' => ',  asin(sin_values_lat(ii)) / real_deg2rad 
        print *, ' Longitude '
        print *, acos(cos_search(jind, lon)) / real_deg2rad, &
                 asin(sin_search(jind, lon)) / real_deg2rad
        print *, ' => ', acos(cos_values_lon(ii)) / real_deg2rad , &
                         asin(sin_values_lon(ii)) / real_deg2rad
        print *, ' found neighbour ', ind,  neighbors_3d (1, ind, n)

        m1 = 1
        m3 = 0
        do m2 = grid_valid_shape(1,1), neighbors_3d (1, ind, n)
           m4 = Grids(grid_id)%partition(m1,1)
           m3 = m3 + 1
           if ( m3 == Grids(grid_id)%extent(m1,1) ) then
              m1 = m1 + 1
              m3 = 0
           endif
        enddo   
        print *, ' global index for found neighbour ', ind, m4+m3
#endif
        if (associated(extra_search%dist_real)) &
        extra_search%dist_real(jind,n) = dist(iloc(1))

     enddo

  end do ! jind
 
  !
  !===> All done
  !

#ifdef VERBOSE
  print 9980, trim(ch_id), ierror

  call psmile_flushstd
#endif /* VERBOSE */
  !
  !  Formats:
  !
9990 format (1x, a, ': psmile_srch_nneigh_gauss2_real')
9980 format (1x, a, ': psmile_srch_nneigh_gauss2_real: eof, ierror =', i3)

end subroutine PSMILe_srch_nneigh_gauss2_real