anchor.diff:  if (GRD%var%missing) then
     if (GRD%var%missing_isnan) then
       where(isnan(F))
  -      F = zero
  +      F = spval
       elsewhere
         F = GRD%var%add_offset + GRD%var%scale_factor*F
       endwhere
     else
       where((F.eq.GRD%var%missing_value))
  -      F = zero
  +      F = spval
       elsewhere
         F = GRD%var%add_offset + GRD%var%scale_factor*F
       endwhere
  
positive.code: ! ...
  ! =====================================================================
  ! ...
  function grid_read_point(GRD,coords) result(F)
  
  class(type_grid), intent(in)           :: GRD
  integer, dimension(:), intent(in)      :: coords
  real(dp)                               :: F
  
  ! ... Local variables
  
negative.code: ! ...
  ! =====================================================================
  ! ...
  function grid_read2D(GRD,layer,step,missing_value) result(F)
  
  class(type_grid), intent(in)           :: GRD
  integer, intent(in), optional          :: layer 
  integer, intent(in), optional          :: step 
  real(dp), intent(in), optional         :: missing_value
  real(dp), dimension(GRD%ni,GRD%nj)     :: F
  
  ! ... Local variables
  ! ...
  logical ZZ,TT
  integer k,l,err
  real(dp) spval
  
  if (present(missing_value)) then
    spval = missing_value
  else
    spval = 0.0D0
  endif
  
  
  STATUS_ERROR = 0; STATUS_TEXT  = ""
  