Geogrid doesn't find geog tiles that contains certain grid points

[mgduda]

When a WRF grid point maps to a half-integer value in the projection space of a geographical dataset, and when that point is a distance of 0.5 outside of the nominal coverage of the tile, geogrid may incorrectly determine that the point lies within the tile, while interpolation schemes will think that the point lies outside of the valid coverage of the tile.

This issue can lead to points that are assigned missing values along the boundaries of source tiles.

The following namelist will reproduce the issue in the WPS v4.0 pre-release 1 code:

&share
 wrf_core = 'ARW',
 max_dom = 5,
 active_grid = f, f, f, f, t
/

&geogrid
 parent_id         = 1,1,2,3,4,
 parent_grid_ratio = 1,3,3,3,11,
 i_parent_start    = 1,20,32,25,28,
 j_parent_start    = 1,28,45,35,41,
 e_we          = 70,100,100,151,1002,
 e_sn          = 70,100,100,151,1013,
 geog_data_res = '2m','2m','30s','30s','30s',
 dx = 27000, 
 dy = 27000,
 map_proj =  'polar',
 ref_lat   = 66.5,
 ref_lon   = 20,
 truelat1  = 70,
 stand_lon = 20,
 geog_data_path = '/glade/u/home/wrfhelp/WPS_GEOG/',
/

Specifically, in the geo_em.d05.nc file, look for bad values at, e.g., (ix,jy)=(349,615).

[mgduda]

In geogrid/src/proc_point_module.F, in the function is_point_in_tile, we use the following logic to determine whether an (rx,ry) point is within a tile for the purposes of interpolation:

      if (src_min_x+src_npts_bdr <= floor(rx+0.5) .and. ceiling(rx-0.5) <= src_max_x-src_npts_bdr .and. &
          src_min_y+src_npts_bdr <= floor(ry+0.5) .and. ceiling(ry-0.5) <= src_max_y-src_npts_bdr) then
         is_point_in_tile = .true.
      else
         is_point_in_tile = .false.
      end if

This is inconsistent with interpolation schemes, e.g., nearest_neighbor in geogrid/src/interp_module.F:

      ix = nint(xx)
      iy = nint(yy)

As a temporary fix that needs more careful consideration, we can use the following logic in is_point_in_tile:

      if (src_min_x+src_npts_bdr <= nint(rx) .and. nint(rx) <= src_max_x-src_npts_bdr .and. &
          src_min_y+src_npts_bdr <= nint(ry) .and. nint(ry) <= src_max_y-src_npts_bdr) then
         is_point_in_tile = .true.
      else
         is_point_in_tile = .false.
      end if