netCDF: Better get_file_fields; hack get_file_info

Two changes, one an iterative improvement, the other a hack to get
things running.

1. The get_file_fields has been rewritten to fetch both axes and fields,
and there is some minimal effort to distinguish between them when
sweeping through variables.  The results are then populated into the
linked lists of the handle.

Two caveats:
    - it's untested, probably not even working
    - we don't even use this function anymore, see change #2

2. get_file_info now just returns -1 for nvars.

Why?  because previously FMS would initialize its nvar to -1, and if you
query a newly opened file then this is what you see.  That was an
indicator to "recreate" the file (actually its metadata, which pulls
from existing templates and now matches the opened file).

To spoof this, I just always return -1 and trigger this "recreate"
effect.  I can do this because the file metadata is static across runs.

Not robust, but it makes sense in a perverse way...

I will try to make this work honestly, so that get_file_fields will
return the correct content.  But for now this feels sufficient.

src/framework/MOM_io_file.F90

@@ -1498,7 +1498,9 @@ subroutine get_file_info_nc(handle, ndim, nvar, ntime)
 
   ! MOM I/O follows legacy FMS behavior and excludes axes from field count
   if (present(ndim)) ndim = ndim_nc
-  if (present(nvar)) nvar = nvar_nc - ndim_nc
+  !if (present(nvar)) nvar = nvar_nc - ndim_nc
+  ! just trying this...
+  if (present(nvar)) nvar = -1
 end subroutine get_file_info_nc
 
 
@@ -1507,21 +1509,22 @@ subroutine get_file_fields_nc(handle, fields)
   class(MOM_netcdf_file), intent(inout) :: handle
   type(MOM_field), intent(inout) :: fields(:)
 
+  type(netcdf_axis), allocatable :: axes_nc(:)
   type(netcdf_field), allocatable :: fields_nc(:)
   integer :: i
-  character(len=64) :: label
 
-  allocate(fields_nc(size(fields)))
-  call get_netcdf_fields(handle%handle_nc, fields_nc)
+  call get_netcdf_fields(handle%handle_nc, axes_nc, fields_nc)
+  if (size(fields) /= size(fields_nc)) error stop
+  ! TODO: call MOM_error
 
-  do i = 1, size(fields)
-    ! TODO
-    !call get_field_atts(fields_nc(i), name=label)
-    !call handle%fields%append(fields_nc(i), trim(label))
-    !fields(i)%label = trim(label)
+  do i = 1, size(axes_nc)
+    call handle%axes%append(axes_nc(i), axes_nc(i)%label)
   enddo
 
-  !call get_file_fields(handle%handle_nc, fields)
+  do i = 1, size(fields)
+    fields(i)%label = trim(fields_nc(i)%label)
+    call handle%fields%append(fields_nc(i), fields_nc(i)%label)
+  enddo
 end subroutine get_file_fields_nc
 
 

src/framework/MOM_netcdf.F90

@@ -11,11 +11,12 @@ module MOM_netcdf
 use netcdf, only : nf90_enddef
 use netcdf, only : nf90_def_dim, nf90_def_var
 use netcdf, only : NF90_UNLIMITED
+use netcdf, only : nf90_get_var
 use netcdf, only : nf90_put_var, nf90_put_att
 use netcdf, only : NF90_FLOAT, NF90_DOUBLE
 use netcdf, only : nf90_strerror, NF90_NOERR
 use netcdf, only : NF90_GLOBAL
-use netcdf, only : nf90_inquire, nf90_inquire_dimension
+use netcdf, only : nf90_inquire, nf90_inquire_dimension, nf90_inquire_variable
 use netcdf, only : nf90_inq_dimids, nf90_inq_varids
 
 use MOM_error_handler, only : MOM_error, FATAL
@@ -534,48 +535,102 @@ subroutine get_netcdf_size(handle, ndims, nvars, nsteps)
 end subroutine get_netcdf_size
 
 
-! TODO: Can this be a function?
-subroutine get_netcdf_fields(handle, fields)
+subroutine get_netcdf_fields(handle, axes, fields)
   type(netcdf_file_type), intent(in) :: handle
+  type(netcdf_axis), intent(inout), allocatable :: axes(:)
   type(netcdf_field), intent(inout), allocatable :: fields(:)
 
   integer :: ndims
     ! Number of netCDF dimensions
   integer :: nvars
-    ! Number of netCDF variables (including axes)
+    ! Number of netCDF dimensions
   integer :: nfields
-    ! Number of fields (i.e. non-axis netCDF variables)
+    ! Number of fields in the file (i.e. non-axis variables)
+  integer, allocatable :: dimids(:)
+    ! netCDF dimension IDs of file
+  integer, allocatable :: varids(:)
+    ! netCDF variable IDs of file
+  character(len=:), allocatable :: label
+    ! Current dimension or variable label
+  integer :: len
+    ! Current dimension length
   integer :: rc
     ! netCDF return code
-  integer, allocatable :: dim_ids(:)
-    ! netCDF dimension IDs
-  integer, allocatable :: field_ids(:)
-    ! netCDF field variable IDs
   integer :: grp_ndims, grp_nvars
     ! Group-based counts for nf90_inq_* (unused)
-  integer :: i
+  logical :: is_axis
+    ! True if the current variable is an axis
+  integer :: i, j, n
 
   integer, save :: no_parent_groups = 0
     ! Flag indicating exclusion of parent groups in netCDF file
     ! NOTE: This must be passed as a variable, and cannot be declared as a
-    !   parameter, so must be declared as a variable.
+    !   parameter.
 
   rc = nf90_inquire(handle%ncid, nDimensions=ndims, nVariables=nvars)
   call check_netcdf_call(rc, 'get_netcdf_fields', &
       'File "' // handle%filename // '"')
 
-  ! Gather the axes, to be removed from the field list
-  rc = nf90_inq_dimids(handle%ncid, grp_ndims, dim_ids, no_parent_groups)
+  allocate(dimids(ndims))
+  rc = nf90_inq_dimids(handle%ncid, grp_ndims, dimids, no_parent_groups)
   call check_netcdf_call(rc, 'get_netcdf_fields', &
       'File "' // handle%filename // '"')
 
-  allocate(field_ids(nvars))
-  rc = nf90_inq_varids(handle%ncid, grp_nvars, field_ids)
+  allocate(varids(nvars))
+  rc = nf90_inq_varids(handle%ncid, grp_nvars, varids)
   call check_netcdf_call(rc, 'get_netcdf_fields', &
       'File "' // handle%filename // '"')
 
+  do i = 1, ndims
+    rc = nf90_inquire_dimension(handle%ncid, dimids(i), name=label, len=len)
+    call check_netcdf_call(rc, 'get_netcdf_fields', &
+        'File "' // handle%filename // '"')
+
+    axes(i)%dimid = dimids(i)
+    axes(i)%label = label
+    allocate(axes(i)%points(len))
+    ! NOTE: points are populated in the variable loop
+
+    ! TODO: Attributes?
+  enddo
+
   nfields = nvars - ndims
   allocate(fields(nfields))
+
+  n = 0
+  do i = 1, nvars
+    rc = nf90_inquire_variable(handle%ncid, varids(i), name=label)
+    call check_netcdf_call(rc, 'get_netcdf_fields', &
+        'File "' // handle%filename // '"')
+
+    ! Check if variable is an axis
+    is_axis = .false.
+    do j = 1, ndims
+      if (label == axes(j)%label) then
+        rc = nf90_get_var(handle%ncid, varids(i), axes(j)%points)
+        call check_netcdf_call(rc, 'get_netcdf_fields', &
+            'File "' // handle%filename // '"')
+        axes(j)%varid = varids(i)
+
+        is_axis = .true.
+        exit
+      endif
+    enddo
+    if (is_axis) continue
+
+    n = n + 1
+    fields(n)%label = label
+    fields(n)%varid = varids(i)
+    fields(n)%xtype = NF90_DOUBLE   ! TODO: Generalize this
+    ! TODO: This is important
+    ! fields(n)%axes
+
+    ! TODO: Probably not important but maybe?
+    ! fields(n)%longname
+    ! fields(n)%units
+    ! fields(n)%standard_name
+    ! fields(n)%checksum
+  enddo
 end subroutine get_netcdf_fields