<p><b>duda</b> 2012-04-10 16:13:32 -0600 (Tue, 10 Apr 2012)</p><p>BRANCH COMMIT<br>
<br>
Convert the mpas_io_input module to use the new IO layers.<br>
<br>
<br>
M    src/framework/mpas_io_input.F<br>
M    src/core_nhyd_atmos/mpas_atm_mpas_core.F<br>
M    src/registry/gen_inc.c<br>
</p><hr noshade><pre><font color="gray">Modified: branches/omp_blocks/io/src/core_nhyd_atmos/mpas_atm_mpas_core.F
===================================================================
--- branches/omp_blocks/io/src/core_nhyd_atmos/mpas_atm_mpas_core.F        2012-04-10 22:08:32 UTC (rev 1766)
+++ branches/omp_blocks/io/src/core_nhyd_atmos/mpas_atm_mpas_core.F        2012-04-10 22:13:32 UTC (rev 1767)
@@ -34,6 +34,7 @@
       type (MPAS_TimeInterval_type) :: timeDiff, minTimeDiff
       character(len=32) :: timeStamp
       integer :: i
+      integer :: ierr
 
       if (.not. config_do_restart) call atm_setup_test_case(domain)
 
@@ -58,48 +59,17 @@
          sfc_update_obj % filename = trim(config_sfc_update_name)
          sfc_update_obj % stream = STREAM_SFC
 
-         call mpas_io_input_init(sfc_update_obj, domain % dminfo)
+         call mpas_io_input_init(sfc_update_obj, domain % blocklist, domain % dminfo)
 
          !     
          ! We need to decide which time slice to read from the surface file - read the most recent time slice that falls before or on the start time
          !
-         sfc_update_obj % time = 1
-
-         if (sfc_update_obj % rdLocalTime &lt;= 0) then
-            write(0,*) 'Error: Couldn''t find any times in surface update file.'
+         sfc_update_obj % time = MPAS_seekStream(sfc_update_obj % io_stream, trim(config_start_time), MPAS_STREAM_LATEST_BEFORE, timeStamp, ierr)
+         if (ierr == MPAS_IO_ERR) then
+            write(0,*) 'Error: surface update file '//trim(sfc_update_obj % filename)//' did not contain any times at or before '//trim(config_start_time)
             call mpas_dmpar_abort(domain % dminfo)
          end if
-      
-         if (domain % dminfo % my_proc_id == IO_NODE) then
-            allocate(xtime % ioinfo)
-            xtime % ioinfo % start(1) = 1
-            xtime % ioinfo % count(1) = sfc_update_obj % rdLocalTime
-            allocate(xtime % array(sfc_update_obj % rdLocalTime))
 
-            xtime % ioinfo % fieldName = 'xtime'
-            call mpas_io_input_field(sfc_update_obj, xtime)
-
-            call mpas_set_timeInterval(interval=minTimeDiff, DD=10000)
-            call mpas_set_time(curr_time=startTime, dateTimeString=config_start_time)
-
-            do i=1,sfc_update_obj % rdLocalTime
-               call mpas_set_time(curr_time=sliceTime, dateTimeString=xtime % array(i))
-               timeDiff = abs(sliceTime - startTime)
-               if (sliceTime &lt;= startTime .and. timeDiff &lt; minTimeDiff) then
-                  minTimeDiff = timeDiff
-                  sfc_update_obj % time = i
-               end if
-            end do
-
-            timeStamp = xtime % array(sfc_update_obj % time)
-
-            deallocate(xtime % ioinfo)
-            deallocate(xtime % array)
-         end if
-
-         call mpas_dmpar_bcast_int(domain % dminfo, sfc_update_obj % time)
-         call mpas_dmpar_bcast_char(domain % dminfo, timeStamp)
-
          write(0,*) 'Starting model with surface time ', timeStamp
 
       end if
@@ -292,11 +262,7 @@
          if (mpas_is_alarm_ringing(clock, sfcAlarmID, ierr=ierr)) then
             call mpas_reset_clock_alarm(clock, sfcAlarmID, ierr=ierr)
 
-            call mpas_read_and_distribute_fields(domain % dminfo, sfc_update_obj, domain % blocklist, &amp;
-                                         readCellStart, nReadCells, readEdgeStart, nReadEdges, readVertexStart, nReadVertices, &amp;
-                                         readVertLevelStart, nReadVertLevels, &amp;
-                                         sendCellList, recvCellList, sendEdgeList, recvEdgeList, sendVertexList, recvVertexList, &amp;
-                                         sendVertLevelList, recvVertLevelList)
+            call mpas_read_and_distribute_fields(sfc_update_obj)
             sfc_update_obj % time = sfc_update_obj % time + 1
          end if
 

Modified: branches/omp_blocks/io/src/framework/mpas_io_input.F
===================================================================
--- branches/omp_blocks/io/src/framework/mpas_io_input.F        2012-04-10 22:08:32 UTC (rev 1766)
+++ branches/omp_blocks/io/src/framework/mpas_io_input.F        2012-04-10 22:13:32 UTC (rev 1767)
@@ -24,31 +24,9 @@
 
       type (MPAS_Stream_type) :: io_stream
 
-#include &quot;io_input_obj_decls.inc&quot;
    end type io_input_object
 
 
-   interface mpas_io_input_field
-      module procedure mpas_io_input_field0d_real
-      module procedure mpas_io_input_field1d_real
-      module procedure mpas_io_input_field2d_real
-      module procedure mpas_io_input_field3d_real
-      module procedure mpas_io_input_field1d_integer
-      module procedure mpas_io_input_field2d_integer
-      module procedure mpas_io_input_field0d_char
-      module procedure mpas_io_input_field1d_char
-   end interface mpas_io_input_field
-
-   interface mpas_io_input_field_time
-      module procedure mpas_io_input_field0d_real_time
-      module procedure mpas_io_input_field1d_real_time
-      module procedure mpas_io_input_field2d_real_time
-      module procedure mpas_io_input_field3d_real_time
-      module procedure mpas_io_input_field1d_integer_time
-      module procedure mpas_io_input_field0d_char_time
-      module procedure mpas_io_input_field1d_char_time
-   end interface mpas_io_input_field_time
-
    type (exchange_list), pointer :: sendCellList, recvCellList
    type (exchange_list), pointer :: sendEdgeList, recvEdgeList
    type (exchange_list), pointer :: sendVertexList, recvVertexList
@@ -75,6 +53,10 @@
 
       character (len=16) :: c_on_a_sphere
       real (kind=RKIND) :: r_sphere_radius
+
+      integer :: ierr
+      integer, dimension(:), pointer :: readIndices
+      type (MPAS_IO_Handle_type) :: inputHandle
    
       type (field1dInteger) :: indexToCellIDField
       type (field1dInteger) :: indexToEdgeIDField
@@ -167,7 +149,19 @@
          input_obj % filename = trim(config_input_name)
          input_obj % stream = STREAM_INPUT
       end if
-      call mpas_io_input_init(input_obj, domain % dminfo)
+      inputHandle = MPAS_io_open(trim(input_obj % filename), MPAS_IO_READ, MPAS_IO_PNETCDF, ierr)
+      if (ierr /= MPAS_IO_NOERR) then
+         write(0,*) ' '
+         if (input_obj % stream == STREAM_RESTART) then
+            write(0,*) 'Error opening restart file ''', trim(input_obj % filename), ''''
+         else if (input_obj % stream == STREAM_INPUT) then
+            write(0,*) 'Error opening input file ''', trim(input_obj % filename), ''''
+         else if (input_obj % stream == STREAM_SFC) then
+            write(0,*) 'Error opening sfc file ''', trim(input_obj % filename), ''''
+         end if
+         write(0,*) ' '
+         call mpas_dmpar_abort(domain % dminfo)
+      end if
    
 
       !
@@ -204,7 +198,13 @@
       indexToCellIDField % ioinfo % start(1) = readCellStart
       indexToCellIDField % ioinfo % count(1) = nReadCells
       allocate(indexToCellIDField % array(nReadCells))
-      call mpas_io_input_field(input_obj, indexToCellIDField)
+      allocate(readIndices(nReadCells))
+      do i=1,nReadCells
+         readIndices(i) = i + readCellStart - 1
+      end do
+      call MPAS_io_inq_var(inputHandle, 'indexToCellID', ierr=ierr)
+      call MPAS_io_set_var_indices(inputHandle, 'indexToCellID', readIndices, ierr=ierr)
+      call mpas_io_get_var(inputHandle, 'indexToCellID', indexToCellIDField % array, ierr)
    
 #ifdef HAVE_ZOLTAN
 #ifdef _MPI 
@@ -214,7 +214,9 @@
       xCellField % ioinfo % start(1) = readCellStart
       xCellField % ioinfo % count(1) = nReadCells
       allocate(xCellField % array(nReadCells))
-      call mpas_io_input_field(input_obj, xCellField)
+      call MPAS_io_inq_var(inputHandle, 'xCell', ierr=ierr)
+      call MPAS_io_set_var_indices(inputHandle, 'xCell', readIndices, ierr=ierr)
+      call mpas_io_get_var(inputHandle, 'xCell', xCellField % array, ierr)
 
       ! Cell y-coordinates (in 3d Cartesian space)
       allocate(yCellField % ioinfo)
@@ -222,7 +224,9 @@
       yCellField % ioinfo % start(1) = readCellStart
       yCellField % ioinfo % count(1) = nReadCells
       allocate(yCellField % array(nReadCells))
-      call mpas_io_input_field(input_obj, yCellField)
+      call MPAS_io_inq_var(inputHandle, 'yCell', ierr=ierr)
+      call MPAS_io_set_var_indices(inputHandle, 'yCell', readIndices, ierr=ierr)
+      call mpas_io_get_var(inputHandle, 'yCell', yCellField % array, ierr)
 
       ! Cell z-coordinates (in 3d Cartesian space)
       allocate(zCellField % ioinfo)
@@ -230,9 +234,12 @@
       zCellField % ioinfo % start(1) = readCellStart
       zCellField % ioinfo % count(1) = nReadCells
       allocate(zCellField % array(nReadCells))
-      call mpas_io_input_field(input_obj, zCellField)
+      call MPAS_io_inq_var(inputHandle, 'zCell', ierr=ierr)
+      call MPAS_io_set_var_indices(inputHandle, 'zCell', readIndices, ierr=ierr)
+      call mpas_io_get_var(inputHandle, 'zCell', zCellField % array, ierr)
 #endif
 #endif
+      deallocate(readIndices)
 
 
       ! Global edge indices
@@ -241,7 +248,14 @@
       indexToEdgeIDField % ioinfo % start(1) = readEdgeStart
       indexToEdgeIDField % ioinfo % count(1) = nReadEdges
       allocate(indexToEdgeIDField % array(nReadEdges))
-      call mpas_io_input_field(input_obj, indexToEdgeIDField)
+      allocate(indexToEdgeIDField % array(nReadEdges))
+      allocate(readIndices(nReadEdges))
+      do i=1,nReadEdges
+         readIndices(i) = i + readEdgeStart - 1
+      end do
+      call MPAS_io_inq_var(inputHandle, 'indexToEdgeID', ierr=ierr)
+      call MPAS_io_set_var_indices(inputHandle, 'indexToEdgeID', readIndices, ierr=ierr)
+      call mpas_io_get_var(inputHandle, 'indexToEdgeID', indexToEdgeIDField % array, ierr)
    
 #ifdef HAVE_ZOLTAN
 #ifdef _MPI 
@@ -251,7 +265,9 @@
       xEdgeField % ioinfo % start(1) = readEdgeStart
       xEdgeField % ioinfo % count(1) = nReadEdges
       allocate(xEdgeField % array(nReadEdges))
-      call mpas_io_input_field(input_obj, xEdgeField)
+      call MPAS_io_inq_var(inputHandle, 'xEdge', ierr=ierr)
+      call MPAS_io_set_var_indices(inputHandle, 'xEdge', readIndices, ierr=ierr)
+      call mpas_io_get_var(inputHandle, 'xEdge', xEdgeField % array, ierr)
 
       ! Edge y-coordinates (in 3d Cartesian space)
       allocate(yEdgeField % ioinfo)
@@ -259,7 +275,9 @@
       yEdgeField % ioinfo % start(1) = readEdgeStart
       yEdgeField % ioinfo % count(1) = nReadEdges
       allocate(yEdgeField % array(nReadEdges))
-      call mpas_io_input_field(input_obj, yEdgeField)
+      call MPAS_io_inq_var(inputHandle, 'yEdge', ierr=ierr)
+      call MPAS_io_set_var_indices(inputHandle, 'yEdge', readIndices, ierr=ierr)
+      call mpas_io_get_var(inputHandle, 'yEdge', yEdgeField % array, ierr)
 
       ! Edge z-coordinates (in 3d Cartesian space)
       allocate(zEdgeField % ioinfo)
@@ -267,17 +285,27 @@
       zEdgeField % ioinfo % start(1) = readEdgeStart
       zEdgeField % ioinfo % count(1) = nReadEdges
       allocate(zEdgeField % array(nReadEdges))
-      call mpas_io_input_field(input_obj, zEdgeField)
+      call MPAS_io_inq_var(inputHandle, 'zEdge', ierr=ierr)
+      call MPAS_io_set_var_indices(inputHandle, 'zEdge', readIndices, ierr=ierr)
+      call mpas_io_get_var(inputHandle, 'zEdge', zEdgeField % array, ierr)
 #endif
 #endif
+      deallocate(readIndices)
 
+
       ! Global vertex indices
       allocate(indexToVertexIDField % ioinfo)
       indexToVertexIDField % ioinfo % fieldName = 'indexToVertexID'
       indexToVertexIDField % ioinfo % start(1) = readVertexStart
       indexToVertexIDField % ioinfo % count(1) = nReadVertices
       allocate(indexToVertexIDField % array(nReadVertices))
-      call mpas_io_input_field(input_obj, indexToVertexIDField)
+      allocate(readIndices(nReadVertices))
+      do i=1,nReadVertices
+         readIndices(i) = i + readVertexStart - 1
+      end do
+      call MPAS_io_inq_var(inputHandle, 'indexToVertexID', ierr=ierr)
+      call MPAS_io_set_var_indices(inputHandle, 'indexToVertexID', readIndices, ierr=ierr)
+      call mpas_io_get_var(inputHandle, 'indexToVertexID', indexToVertexIDField % array, ierr)
    
 #ifdef HAVE_ZOLTAN
 #ifdef _MPI
@@ -287,7 +315,9 @@
       xVertexField % ioinfo % start(1) = readVertexStart
       xVertexField % ioinfo % count(1) = nReadVertices
       allocate(xVertexField % array(nReadVertices))
-      call mpas_io_input_field(input_obj, xVertexField)
+      call MPAS_io_inq_var(inputHandle, 'xVertex', ierr=ierr)
+      call MPAS_io_set_var_indices(inputHandle, 'xVertex', readIndices, ierr=ierr)
+      call mpas_io_get_var(inputHandle, 'xVertex', xVertexField % array, ierr)
 
       ! Vertex y-coordinates (in 3d Cartesian space)
       allocate(yVertexField % ioinfo)
@@ -295,7 +325,9 @@
       yVertexField % ioinfo % start(1) = readVertexStart
       yVertexField % ioinfo % count(1) = nReadVertices
       allocate(yVertexField % array(nReadVertices))
-      call mpas_io_input_field(input_obj, yVertexField)
+      call MPAS_io_inq_var(inputHandle, 'yVertex', ierr=ierr)
+      call MPAS_io_set_var_indices(inputHandle, 'yVertex', readIndices, ierr=ierr)
+      call mpas_io_get_var(inputHandle, 'yVertex', yVertexField % array, ierr)
 
       ! Vertex z-coordinates (in 3d Cartesian space)
       allocate(zVertexField % ioinfo)
@@ -303,9 +335,12 @@
       zVertexField % ioinfo % start(1) = readVertexStart
       zVertexField % ioinfo % count(1) = nReadVertices
       allocate(zVertexField % array(nReadVertices))
-      call mpas_io_input_field(input_obj, zVertexField)
+      call MPAS_io_inq_var(inputHandle, 'zVertex', ierr=ierr)
+      call MPAS_io_set_var_indices(inputHandle, 'zVertex', readIndices, ierr=ierr)
+      call mpas_io_get_var(inputHandle, 'zVertex', zVertexField % array, ierr)
 #endif
 #endif
+      deallocate(readIndices)
 
       ! Number of cell/edges/vertices adjacent to each cell
       allocate(nEdgesOnCellField % ioinfo)
@@ -313,7 +348,13 @@
       nEdgesOnCellField % ioinfo % start(1) = readCellStart
       nEdgesOnCellField % ioinfo % count(1) = nReadCells
       allocate(nEdgesOnCellField % array(nReadCells))
-      call mpas_io_input_field(input_obj, nEdgesOnCellField)
+      allocate(readIndices(nReadCells))
+      do i=1,nReadCells
+         readIndices(i) = i + readCellStart - 1
+      end do
+      call MPAS_io_inq_var(inputHandle, 'nEdgesOnCell', ierr=ierr)
+      call MPAS_io_set_var_indices(inputHandle, 'nEdgesOnCell', readIndices, ierr=ierr)
+      call mpas_io_get_var(inputHandle, 'nEdgesOnCell', nEdgesOnCellField % array, ierr)
    
       ! Global indices of cells adjacent to each cell
       allocate(cellsOnCellField % ioinfo)
@@ -323,7 +364,9 @@
       cellsOnCellField % ioinfo % count(1) = maxEdges
       cellsOnCellField % ioinfo % count(2) = nReadCells
       allocate(cellsOnCellField % array(maxEdges,nReadCells))
-      call mpas_io_input_field(input_obj, cellsOnCellField)
+      call MPAS_io_inq_var(inputHandle, 'cellsOnCell', ierr=ierr)
+      call MPAS_io_set_var_indices(inputHandle, 'cellsOnCell', readIndices, ierr=ierr)
+      call mpas_io_get_var(inputHandle, 'cellsOnCell', cellsOnCellField % array, ierr)
    
       ! Global indices of edges adjacent to each cell
       allocate(edgesOnCellField % ioinfo)
@@ -333,7 +376,9 @@
       edgesOnCellField % ioinfo % count(1) = maxEdges
       edgesOnCellField % ioinfo % count(2) = nReadCells
       allocate(edgesOnCellField % array(maxEdges,nReadCells))
-      call mpas_io_input_field(input_obj, edgesOnCellField)
+      call MPAS_io_inq_var(inputHandle, 'edgesOnCell', ierr=ierr)
+      call MPAS_io_set_var_indices(inputHandle, 'edgesOnCell', readIndices, ierr=ierr)
+      call mpas_io_get_var(inputHandle, 'edgesOnCell', edgesOnCellField % array, ierr)
    
       ! Global indices of vertices adjacent to each cell
       allocate(verticesOnCellField % ioinfo)
@@ -343,7 +388,10 @@
       verticesOnCellField % ioinfo % count(1) = maxEdges
       verticesOnCellField % ioinfo % count(2) = nReadCells
       allocate(verticesOnCellField % array(maxEdges,nReadCells))
-      call mpas_io_input_field(input_obj, verticesOnCellField)
+      call MPAS_io_inq_var(inputHandle, 'verticesOnCell', ierr=ierr)
+      call MPAS_io_set_var_indices(inputHandle, 'verticesOnCell', readIndices, ierr=ierr)
+      call mpas_io_get_var(inputHandle, 'verticesOnCell', verticesOnCellField % array, ierr)
+      deallocate(readIndices)
    
       ! Global indices of cells adjacent to each edge
       !    used for determining which edges are owned by a block, where 
@@ -355,7 +403,14 @@
       cellsOnEdgeField % ioinfo % count(1) = 2
       cellsOnEdgeField % ioinfo % count(2) = nReadEdges
       allocate(cellsOnEdgeField % array(2,nReadEdges))
-      call mpas_io_input_field(input_obj, cellsOnEdgeField)
+      allocate(readIndices(nReadEdges))
+      do i=1,nReadEdges
+         readIndices(i) = i + readEdgeStart - 1
+      end do
+      call MPAS_io_inq_var(inputHandle, 'cellsOnEdge', ierr=ierr)
+      call MPAS_io_set_var_indices(inputHandle, 'cellsOnEdge', readIndices, ierr=ierr)
+      call mpas_io_get_var(inputHandle, 'cellsOnEdge', cellsOnEdgeField % array, ierr)
+      deallocate(readIndices)
    
       ! Global indices of cells adjacent to each vertex
       !    used for determining which vertices are owned by a block, where 
@@ -367,7 +422,14 @@
       cellsOnVertexField % ioinfo % count(1) = vertexDegree
       cellsOnVertexField % ioinfo % count(2) = nReadVertices
       allocate(cellsOnVertexField % array(vertexDegree,nReadVertices))
-      call mpas_io_input_field(input_obj, cellsOnVertexField)
+      allocate(readIndices(nReadVertices))
+      do i=1,nReadVertices
+         readIndices(i) = i + readVertexStart - 1
+      end do
+      call MPAS_io_inq_var(inputHandle, 'cellsOnVertex', ierr=ierr)
+      call MPAS_io_set_var_indices(inputHandle, 'cellsOnVertex', readIndices, ierr=ierr)
+      call mpas_io_get_var(inputHandle, 'cellsOnVertex', cellsOnVertexField % array, ierr)
+      deallocate(readIndices)
    
    
       !
@@ -839,67 +901,69 @@
 #include &quot;dim_dummy_args.inc&quot;
                          )
 
+!!!!!!!!!!MGD HERE WE NEED TO READ IN indexTo*ID fields !!!!!!!!!!!!!!!!!
+      call MPAS_io_inq_var(inputHandle, 'indexToCellID', ierr=ierr)
+      call MPAS_io_set_var_indices(inputHandle, 'indexToCellID', local_cell_list(1:nOwnCells), ierr=ierr)
+      call mpas_io_get_var(inputHandle, 'indexToCellID', domain % blocklist % mesh % indexToCellID % array, ierr)
+
+      call MPAS_io_inq_var(inputHandle, 'indexToEdgeID', ierr=ierr)
+      call MPAS_io_set_var_indices(inputHandle, 'indexToEdgeID', local_edge_list(1:nOwnEdges), ierr=ierr)
+      call mpas_io_get_var(inputHandle, 'indexToEdgeID', domain % blocklist % mesh % indexToEdgeID % array, ierr)
+
+      call MPAS_io_inq_var(inputHandle, 'indexToVertexID', ierr=ierr)
+      call MPAS_io_set_var_indices(inputHandle, 'indexToVertexID', local_vertex_list(1:nOwnVertices), ierr=ierr)
+      call mpas_io_get_var(inputHandle, 'indexToVertexID', domain % blocklist % mesh % indexToVertexID % array, ierr)
+
+      domain % blocklist % mesh % nCellsSolve = nOwnCells
+      domain % blocklist % mesh % nEdgesSolve = nOwnEdges
+      domain % blocklist % mesh % nVerticesSolve = nOwnVertices
+      domain % blocklist % mesh % nVertLevelsSolve = domain % blocklist % mesh % nVertLevels   ! No vertical decomp yet...
+write(0,*) 'Setting nCellsSolve = ', domain % blocklist % mesh % nCellsSolve
+write(0,*) 'Setting nEdgesSolve = ', domain % blocklist % mesh % nEdgesSolve
+write(0,*) 'Setting nVerticesSolve = ', domain % blocklist % mesh % nVerticesSolve
+write(0,*) 'Setting nVertLevelsSolve = ', domain % blocklist % mesh % nVertLevelsSolve
+
+      call mpas_io_input_init(input_obj, domain % blocklist, domain % dminfo)
+
+
       !
       ! Read attributes
       !
-      call mpas_io_input_get_att_text(input_obj, 'on_a_sphere', c_on_a_sphere)
-      call mpas_io_input_get_att_real(input_obj, 'sphere_radius', r_sphere_radius)
-      if (index(c_on_a_sphere, 'YES') /= 0) then
+      call MPAS_readStreamAtt(input_obj % io_stream, 'sphere_radius', r_sphere_radius, ierr)
+      if (ierr /= MPAS_STREAM_NOERR) then
+         write(0,*) 'Warning: Attribute sphere_radius not found in '//trim(input_obj % filename)
+         write(0,*) '   Setting sphere_radius to 1.0'
+         domain % blocklist % mesh % sphere_radius = 1.0
+      else
+         domain % blocklist % mesh % sphere_radius = r_sphere_radius
+      end if
+
+      call MPAS_readStreamAtt(input_obj % io_stream, 'on_a_sphere', c_on_a_sphere, ierr)
+      if (ierr /= MPAS_STREAM_NOERR) then
+         write(0,*) 'Warning: Attribute on_a_sphere not found in '//trim(input_obj % filename)
+         write(0,*) '   Setting on_a_sphere to ''YES'''
          domain % blocklist % mesh % on_a_sphere = .true.
       else
-         domain % blocklist % mesh % on_a_sphere = .false.
+         if (index(c_on_a_sphere, 'YES') /= 0) then
+            domain % blocklist % mesh % on_a_sphere = .true.
+         else
+            domain % blocklist % mesh % on_a_sphere = .false.
+         end if
       end if
-      domain % blocklist % mesh % sphere_radius = r_sphere_radius
 
       if (.not. config_do_restart) then
          input_obj % time = 1
       else
-         input_obj % time = 1
-
          !
          ! If doing a restart, we need to decide which time slice to read from the 
          !   restart file
          !
-         if (input_obj % rdLocalTime &lt;= 0) then
-            write(0,*) 'Error: Couldn''t find any times in restart file.'
+         input_obj % time = MPAS_seekStream(input_obj % io_stream, config_start_time, MPAS_STREAM_EXACT_TIME, timeStamp, ierr)
+         if (ierr == MPAS_IO_ERR) then
+            write(0,*) 'Error: restart file '//trim(filename)//' did not contain time '//trim(config_start_time)
             call mpas_dmpar_abort(domain % dminfo)
          end if
-         if (domain % dminfo % my_proc_id == IO_NODE) then
-            allocate(xtime % ioinfo)
-            xtime % ioinfo % start(1) = 1
-            xtime % ioinfo % count(1) = input_obj % rdLocalTime
-            allocate(xtime % array(input_obj % rdLocalTime))
-
-            xtime % ioinfo % fieldName = 'xtime'
-            call mpas_io_input_field(input_obj, xtime)
-
-            call mpas_set_timeInterval(interval=minTimeDiff, DD=10000)
-            call mpas_set_time(curr_time=startTime, dateTimeString=config_start_time)
-
-            do i=1,input_obj % rdLocalTime
-               call mpas_set_time(curr_time=sliceTime, dateTimeString=xtime % array(i))
-               timeDiff = abs(sliceTime - startTime)
-               if (timeDiff &lt; minTimeDiff) then
-                  minTimeDiff = timeDiff
-                  input_obj % time = i
-               end if
-            end do
-
-            ! require restart time to exactly match start time (this error should never be reached as we have by this point opened the restart file with a name containing the startTime)
-            if(sliceTime /= startTime) then
-               write(0,*) &quot;Error: restart file &quot;, filename, &quot; did not contain time &quot;, config_start_time
-               call mpas_dmpar_abort(domain % dminfo)
-            end if
-
-            timeStamp = xtime % array(input_obj % time)
-
-            deallocate(xtime % ioinfo)
-            deallocate(xtime % array)
-         end if
-
-         call mpas_dmpar_bcast_int(domain % dminfo, input_obj % time)
-         call mpas_dmpar_bcast_char(domain % dminfo, timeStamp)
-
+write(0,*) 'MGD DEBUGGING time = ', input_obj % time
          write(0,*) 'Restarting model from time ', timeStamp
 
       end if
@@ -915,17 +979,134 @@
       !      processes that own those indices based on 
       !      {send,recv}{Cell,Edge,Vertex,VertLevel}List
       !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 
-      call mpas_read_and_distribute_fields(domain % dminfo, input_obj, domain % blocklist, &amp;
-                                      readCellStart, nReadCells, readEdgeStart, nReadEdges, readVertexStart, nReadVertices, &amp;
-                                      readVertLevelStart, nReadVertLevels, &amp;
-                                      sendCellList, recvCellList, sendEdgeList, recvEdgeList, sendVertexList, recvVertexList, &amp;
-                                      sendVertLevelList, recvVertLevelList) 
+      call mpas_read_and_distribute_fields(input_obj)
 
-
       call mpas_io_input_finalize(input_obj, domain % dminfo)
 
+      call MPAS_io_close(inputHandle, ierr)
+
    
       !
+      ! Work out halo exchange lists for cells, edges, and vertices
+      ! NB: The next pointer in each element of, e.g., cellsToSend, acts as the head pointer of
+      !     the list, since Fortran does not allow arrays of pointers
+      !
+
+      !--------- Create Cell Exchange Lists ---------!
+
+      ! pass in neededList of ownedCells and halo layer 1 cells
+      call mpas_dmpar_get_owner_list(domain % dminfo, &amp;
+                                nOwnCells, nCellsCumulative(2), &amp;
+                                block_graph_2Halo % vertexID(1:nOwnCells), block_graph_2Halo % vertexID(1 : nCellsCumulative(2)), &amp;
+                                domain % blocklist % parinfo % cellsToSend(1) % next, domain % blocklist % parinfo % cellsToRecv(1) % next)
+
+      ! pass in neededList of ownedCells and halo layer 2 cells; offset of number of halo 1 cells is required
+      offset = nCellsHalo(1)
+      nTempIDs = nOwnCells + nCellsHalo(2)
+      allocate(tempIDs(nTempIDs))
+      tempIDs(1:nOwnCells) = block_graph_2Halo % vertexID(1:nOwnCells)
+      tempIDs(nOwnCells+1:nTempIDs) = block_graph_2Halo % vertexID(nCellsCumulative(2)+1 : nCellsCumulative(3))
+      call mpas_dmpar_get_owner_list(domain % dminfo, &amp;
+                                nOwnCells, nTempIDs, &amp;
+                                block_graph_2Halo % vertexID(1:nOwnCells), tempIDs, &amp;
+                                domain % blocklist % parinfo % cellsToSend(2) % next, domain % blocklist % parinfo % cellsToRecv(2) % next, &amp;
+                                offset)
+      deallocate(tempIDs)
+
+
+      !--------- Create Edge Exchange Lists ---------!
+
+      ! pass in neededList of ownedEdges and ownedCell perimeter edges
+      call mpas_dmpar_get_owner_list(domain % dminfo, &amp;
+                                nOwnEdges, nEdgesCumulative(2), &amp;
+                                local_edge_list(1:nOwnEdges), local_edge_list(1 : nEdgesCumulative(2)), &amp;
+                                domain % blocklist % parinfo % edgesToSend(1) % next, domain % blocklist % parinfo % edgesToRecv(1) % next)
+
+      ! pass in neededList of owned edges and yet-to-be-included edges from halo 1 cells; offset of number of ownedCell perimeter edges is required
+      offset = nEdgesHalo(1)
+      nTempIDs = nOwnEdges + nEdgesHalo(2)
+      allocate(tempIDs(nTempIDs))
+      tempIDs(1:nOwnEdges) = local_edge_list(1:nOwnEdges)
+      tempIDs(nOwnEdges+1:nTempIDs) = local_edge_list(nEdgesCumulative(2)+1 : nEdgesCumulative(3))
+      call mpas_dmpar_get_owner_list(domain % dminfo, &amp;
+                                nOwnEdges, nTempIDs, &amp;
+                                local_edge_list(1:nOwnEdges), tempIDs, &amp;  
+                                domain % blocklist % parinfo % edgesToSend(2) % next, domain % blocklist % parinfo % edgesToRecv(2) % next, &amp;
+                                offset)
+      deallocate(tempIDs)
+
+      ! pass in neededList of owned edges and yet-to-be-included edges from halo 2 cells; offset of number of ownedCell perimeter edges and halo 1 edges is required
+      offset = nEdgesHalo(1) + nEdgesHalo(2)
+      nTempIDs = nOwnEdges + nEdgesHalo(3)
+      allocate(tempIDs(nTempIDs))
+      tempIDs(1:nOwnEdges) = local_edge_list(1:nOwnEdges)
+      tempIDs(nOwnEdges+1:nTempIDs) = local_edge_list(nEdgesCumulative(3)+1 : nEdgesCumulative(4))
+      call mpas_dmpar_get_owner_list(domain % dminfo, &amp;
+                                nOwnEdges, nTempIDs, &amp;
+                                local_edge_list(1:nOwnEdges), tempIDs, &amp;  
+                                domain % blocklist % parinfo % edgesToSend(3) % next, domain % blocklist % parinfo % edgesToRecv(3) % next, &amp;
+                                offset)
+      deallocate(tempIDs)
+
+
+      !--------- Create Vertex Exchange Lists ---------!
+
+
+      ! pass in neededList of ownedVertices and ownedCell perimeter vertices
+      call mpas_dmpar_get_owner_list(domain % dminfo, &amp;
+                                nOwnVertices, nVerticesCumulative(2), &amp;
+                                local_vertex_list(1:nOwnVertices), local_vertex_list(1 : nVerticesCumulative(2)), &amp;
+                                domain % blocklist % parinfo % verticesToSend(1) % next, domain % blocklist % parinfo % verticesToRecv(1) % next)
+
+      ! pass in neededList of owned vertices and yet-to-be-included vertices from halo 1 cells; offset of number of ownedCell perimeter vertices is required
+      offset = nVerticesHalo(1)
+      nTempIDs = nOwnVertices + nVerticesHalo(2)
+      allocate(tempIDs(nTempIDs))
+      tempIDs(1:nOwnVertices) = local_vertex_list(1:nOwnVertices)
+      tempIDs(nOwnVertices+1:nTempIDs) = local_vertex_list(nVerticesCumulative(2)+1 : nVerticesCumulative(3))
+      call mpas_dmpar_get_owner_list(domain % dminfo, &amp;
+                                nOwnVertices, nTempIDs, &amp;
+                                local_vertex_list(1:nOwnVertices), tempIDs, &amp;  
+                                domain % blocklist % parinfo % verticesToSend(2) % next, domain % blocklist % parinfo % verticesToRecv(2) % next, &amp;
+                                offset)
+      deallocate(tempIDs)
+
+      ! pass in neededList of owned vertices and yet-to-be-included vertices from halo 2 cells; offset of number of ownedCell perimeter vertices and halo 1 vertices is required
+      offset = nVerticesHalo(1) + nVerticesHalo(2)
+      nTempIDs = nOwnVertices + nVerticesHalo(3)
+      allocate(tempIDs(nTempIDs))
+      tempIDs(1:nOwnVertices) = local_vertex_list(1:nOwnVertices)
+      tempIDs(nOwnVertices+1:nTempIDs) = local_vertex_list(nVerticesCumulative(3)+1 : nVerticesCumulative(4))
+      call mpas_dmpar_get_owner_list(domain % dminfo, &amp;
+                                nOwnVertices, nTempIDs, &amp;
+                                local_vertex_list(1:nOwnVertices), tempIDs, &amp;  
+                                domain % blocklist % parinfo % verticesToSend(3) % next, domain % blocklist % parinfo % verticesToRecv(3) % next, &amp;
+                                offset)
+      deallocate(tempIDs)
+
+
+      domain % blocklist % mesh % nCellsSolve = nOwnCells
+      domain % blocklist % mesh % nEdgesSolve = nOwnEdges
+      domain % blocklist % mesh % nVerticesSolve = ghostVertexStart-1
+      domain % blocklist % mesh % nVertLevelsSolve = domain % blocklist % mesh % nVertLevels   ! No vertical decomp yet...
+write(0,*) 'Setting nCellsSolve = ', domain % blocklist % mesh % nCellsSolve
+write(0,*) 'Setting nEdgesSolve = ', domain % blocklist % mesh % nEdgesSolve
+write(0,*) 'Setting nVerticesSolve = ', domain % blocklist % mesh % nVerticesSolve
+write(0,*) 'Setting nVertLevelsSolve = ', domain % blocklist % mesh % nVertLevelsSolve
+
+      ! Link the sendList and recvList pointers in each field type to the appropriate lists 
+      !   in parinfo, e.g., cellsToSend and cellsToRecv; in future, it can also be extended to 
+      !   link blocks of fields to eachother
+      call mpas_create_field_links(domain % blocklist)
+
+
+      !
+      ! Exchange halos for all of the fields that were read from the input file
+      !
+      call mpas_exch_input_field_halos(domain, input_obj)
+
+   
+      !
       ! Rename vertices in cellsOnCell, edgesOnCell, etc. to local indices
       !
       allocate(cellIDSorted(2,domain % blocklist % mesh % nCells))
@@ -960,7 +1141,6 @@
                domain % blocklist % mesh % cellsOnCell % array(j,i) = cellIDSorted(2,k)
             else
                domain % blocklist % mesh % cellsOnCell % array(j,i) = domain % blocklist % mesh % nCells + 1
-!               domain % blocklist % mesh % cellsOnCell % array(j,i) = 0
             end if
 
             k = mpas_binary_search(edgeIDSorted, 2, 1, domain % blocklist % mesh % nEdges, &amp;
@@ -969,7 +1149,6 @@
                domain % blocklist % mesh % edgesOnCell % array(j,i) = edgeIDSorted(2,k)
             else
                domain % blocklist % mesh % edgesOnCell % array(j,i) = domain % blocklist % mesh % nEdges + 1
-!               domain % blocklist % mesh % edgesOnCell % array(j,i) = 0
             end if
 
             k = mpas_binary_search(vertexIDSorted, 2, 1, domain % blocklist % mesh % nVertices, &amp;
@@ -978,7 +1157,6 @@
                domain % blocklist % mesh % verticesOnCell % array(j,i) = vertexIDSorted(2,k)
             else
                domain % blocklist % mesh % verticesOnCell % array(j,i) = domain % blocklist % mesh % nVertices + 1
-!               domain % blocklist % mesh % verticesOnCell % array(j,i) = 0
             end if
 
          end do
@@ -993,7 +1171,6 @@
                domain % blocklist % mesh % cellsOnEdge % array(j,i) = cellIDSorted(2,k)
             else
                domain % blocklist % mesh % cellsOnEdge % array(j,i) = domain % blocklist % mesh % nCells + 1
-!               domain % blocklist % mesh % cellsOnEdge % array(j,i) = 0
             end if
 
             k = mpas_binary_search(vertexIDSorted, 2, 1, domain % blocklist % mesh % nVertices, &amp;
@@ -1002,7 +1179,6 @@
                domain % blocklist % mesh % verticesOnEdge % array(j,i) = vertexIDSorted(2,k)
             else
                domain % blocklist % mesh % verticesOnEdge % array(j,i) = domain % blocklist % mesh % nVertices + 1
-!               domain % blocklist % mesh % verticesOnEdge % array(j,i) = 0
             end if
 
          end do
@@ -1015,7 +1191,6 @@
                domain % blocklist % mesh % edgesOnEdge % array(j,i) = edgeIDSorted(2,k)
             else
                domain % blocklist % mesh % edgesOnEdge % array(j,i) = domain % blocklist % mesh % nEdges + 1
-!               domain % blocklist % mesh % edgesOnEdge % array(j,i) = 0
             end if
 
          end do
@@ -1030,7 +1205,6 @@
                domain % blocklist % mesh % cellsOnVertex % array(j,i) = cellIDSorted(2,k)
             else
                domain % blocklist % mesh % cellsOnVertex % array(j,i) = domain % blocklist % mesh % nCells + 1
-!               domain % blocklist % mesh % cellsOnVertex % array(j,i) = 0
             end if
 
             k = mpas_binary_search(edgeIDSorted, 2, 1, domain % blocklist % mesh % nEdges, &amp;
@@ -1039,7 +1213,6 @@
                domain % blocklist % mesh % edgesOnVertex % array(j,i) = edgeIDSorted(2,k)
             else
                domain % blocklist % mesh % edgesOnVertex % array(j,i) = domain % blocklist % mesh % nEdges + 1
-!               domain % blocklist % mesh % edgesOnVertex % array(j,i) = 0
             end if
 
          end do
@@ -1049,116 +1222,6 @@
       deallocate(edgeIDSorted)
       deallocate(vertexIDSorted)
 
-
-      !
-      ! Work out halo exchange lists for cells, edges, and vertices
-      ! NB: The next pointer in each element of, e.g., cellsToSend, acts as the head pointer of
-      !     the list, since Fortran does not allow arrays of pointers
-      !
-
-      !--------- Create Cell Exchange Lists ---------!
-
-      ! pass in neededList of ownedCells and halo layer 1 cells
-      call mpas_dmpar_get_owner_list(domain % dminfo, &amp;
-                                nOwnCells, nCellsCumulative(2), &amp;
-                                block_graph_2Halo % vertexID(1:nOwnCells), block_graph_2Halo % vertexID(1 : nCellsCumulative(2)), &amp;
-                                domain % blocklist % parinfo % cellsToSend(1) % next, domain % blocklist % parinfo % cellsToRecv(1) % next)
-
-      ! pass in neededList of ownedCells and halo layer 2 cells; offset of number of halo 1 cells is required
-      offset = nCellsHalo(1)
-      nTempIDs = nOwnCells + nCellsHalo(2)
-      allocate(tempIDs(nTempIDs))
-      tempIDs(1:nOwnCells) = block_graph_2Halo % vertexID(1:nOwnCells)
-      tempIDs(nOwnCells+1:nTempIDs) = block_graph_2Halo % vertexID(nCellsCumulative(2)+1 : nCellsCumulative(3))
-      call mpas_dmpar_get_owner_list(domain % dminfo, &amp;
-                                nOwnCells, nTempIDs, &amp;
-                                block_graph_2Halo % vertexID(1:nOwnCells), tempIDs, &amp;
-                                domain % blocklist % parinfo % cellsToSend(2) % next, domain % blocklist % parinfo % cellsToRecv(2) % next, &amp;
-                                offset)
-      deallocate(tempIDs)
-
-
-      !--------- Create Edge Exchange Lists ---------!
-
-      ! pass in neededList of ownedEdges and ownedCell perimeter edges
-      call mpas_dmpar_get_owner_list(domain % dminfo, &amp;
-                                nOwnEdges, nEdgesCumulative(2), &amp;
-                                local_edge_list(1:nOwnEdges), local_edge_list(1 : nEdgesCumulative(2)), &amp;
-                                domain % blocklist % parinfo % edgesToSend(1) % next, domain % blocklist % parinfo % edgesToRecv(1) % next)
-
-      ! pass in neededList of owned edges and yet-to-be-included edges from halo 1 cells; offset of number of ownedCell perimeter edges is required
-      offset = nEdgesHalo(1)
-      nTempIDs = nOwnEdges + nEdgesHalo(2)
-      allocate(tempIDs(nTempIDs))
-      tempIDs(1:nOwnEdges) = local_edge_list(1:nOwnEdges)
-      tempIDs(nOwnEdges+1:nTempIDs) = local_edge_list(nEdgesCumulative(2)+1 : nEdgesCumulative(3))
-      call mpas_dmpar_get_owner_list(domain % dminfo, &amp;
-                                nOwnEdges, nTempIDs, &amp;
-                                local_edge_list(1:nOwnEdges), tempIDs, &amp;  
-                                domain % blocklist % parinfo % edgesToSend(2) % next, domain % blocklist % parinfo % edgesToRecv(2) % next, &amp;
-                                offset)
-      deallocate(tempIDs)
-
-      ! pass in neededList of owned edges and yet-to-be-included edges from halo 2 cells; offset of number of ownedCell perimeter edges and halo 1 edges is required
-      offset = nEdgesHalo(1) + nEdgesHalo(2)
-      nTempIDs = nOwnEdges + nEdgesHalo(3)
-      allocate(tempIDs(nTempIDs))
-      tempIDs(1:nOwnEdges) = local_edge_list(1:nOwnEdges)
-      tempIDs(nOwnEdges+1:nTempIDs) = local_edge_list(nEdgesCumulative(3)+1 : nEdgesCumulative(4))
-      call mpas_dmpar_get_owner_list(domain % dminfo, &amp;
-                                nOwnEdges, nTempIDs, &amp;
-                                local_edge_list(1:nOwnEdges), tempIDs, &amp;  
-                                domain % blocklist % parinfo % edgesToSend(3) % next, domain % blocklist % parinfo % edgesToRecv(3) % next, &amp;
-                                offset)
-      deallocate(tempIDs)
-
-
-      !--------- Create Vertex Exchange Lists ---------!
-
-
-      ! pass in neededList of ownedVertices and ownedCell perimeter vertices
-      call mpas_dmpar_get_owner_list(domain % dminfo, &amp;
-                                nOwnVertices, nVerticesCumulative(2), &amp;
-                                local_vertex_list(1:nOwnVertices), local_vertex_list(1 : nVerticesCumulative(2)), &amp;
-                                domain % blocklist % parinfo % verticesToSend(1) % next, domain % blocklist % parinfo % verticesToRecv(1) % next)
-
-      ! pass in neededList of owned vertices and yet-to-be-included vertices from halo 1 cells; offset of number of ownedCell perimeter vertices is required
-      offset = nVerticesHalo(1)
-      nTempIDs = nOwnVertices + nVerticesHalo(2)
-      allocate(tempIDs(nTempIDs))
-      tempIDs(1:nOwnVertices) = local_vertex_list(1:nOwnVertices)
-      tempIDs(nOwnVertices+1:nTempIDs) = local_vertex_list(nVerticesCumulative(2)+1 : nVerticesCumulative(3))
-      call mpas_dmpar_get_owner_list(domain % dminfo, &amp;
-                                nOwnVertices, nTempIDs, &amp;
-                                local_vertex_list(1:nOwnVertices), tempIDs, &amp;  
-                                domain % blocklist % parinfo % verticesToSend(2) % next, domain % blocklist % parinfo % verticesToRecv(2) % next, &amp;
-                                offset)
-      deallocate(tempIDs)
-
-      ! pass in neededList of owned vertices and yet-to-be-included vertices from halo 2 cells; offset of number of ownedCell perimeter vertices and halo 1 vertices is required
-      offset = nVerticesHalo(1) + nVerticesHalo(2)
-      nTempIDs = nOwnVertices + nVerticesHalo(3)
-      allocate(tempIDs(nTempIDs))
-      tempIDs(1:nOwnVertices) = local_vertex_list(1:nOwnVertices)
-      tempIDs(nOwnVertices+1:nTempIDs) = local_vertex_list(nVerticesCumulative(3)+1 : nVerticesCumulative(4))
-      call mpas_dmpar_get_owner_list(domain % dminfo, &amp;
-                                nOwnVertices, nTempIDs, &amp;
-                                local_vertex_list(1:nOwnVertices), tempIDs, &amp;  
-                                domain % blocklist % parinfo % verticesToSend(3) % next, domain % blocklist % parinfo % verticesToRecv(3) % next, &amp;
-                                offset)
-      deallocate(tempIDs)
-
-
-      domain % blocklist % mesh % nCellsSolve = nOwnCells
-      domain % blocklist % mesh % nEdgesSolve = nOwnEdges
-      domain % blocklist % mesh % nVerticesSolve = ghostVertexStart-1
-      domain % blocklist % mesh % nVertLevelsSolve = domain % blocklist % mesh % nVertLevels   ! No vertical decomp yet...
-
-      ! Link the sendList and recvList pointers in each field type to the appropriate lists 
-      !   in parinfo, e.g., cellsToSend and cellsToRecv; in future, it can also be extended to 
-      !   link blocks of fields to eachother
-      call mpas_create_field_links(domain % blocklist)
-
    
       !
       ! Deallocate fields, graphs, and other memory
@@ -1237,89 +1300,34 @@
    end subroutine mpas_insert_string_suffix
 
 
-   subroutine mpas_read_and_distribute_fields(dminfo, input_obj, block, &amp;
-                                     readCellsStart, readCellsCount, &amp;
-                                     readEdgesStart, readEdgesCount, &amp;
-                                     readVerticesStart, readVerticesCount, &amp;
-                                     readVertLevelsStart, readVertLevelsCount, &amp;
-                                     sendCellsList, recvCellsList, &amp;
-                                     sendEdgesList, recvEdgesList, &amp;
-                                     sendVerticesList, recvVerticesList, &amp;
-                                     sendVertLevelsList, recvVertLevelsList)
+   subroutine mpas_read_and_distribute_fields(input_obj)
       
       implicit none
 
-      type (dm_info), intent(in) :: dminfo
       type (io_input_object), intent(in) :: input_obj
-      type (block_type), intent(inout) :: block
-      integer, intent(in) :: readCellsStart, readCellsCount, readEdgesStart, readEdgesCount, readVerticesStart, readVerticesCount
-      integer, intent(in) :: readVertLevelsStart, readVertLevelsCount
-      type (exchange_list), pointer :: sendCellsList, recvCellsList
-      type (exchange_list), pointer :: sendEdgesList, recvEdgesList
-      type (exchange_list), pointer :: sendVerticesList, recvVerticesList
-      type (exchange_list), pointer :: sendVertLevelsList, recvVertLevelsList
 
-      type (field1dInteger) :: int1d
-      type (field2dInteger) :: int2d
-      type (field0dReal) :: real0d
-      type (field1dReal) :: real1d
-      type (field2dReal) :: real2d
-      type (field3dReal) :: real3d
-      type (field0dChar) :: char0d
-      type (field1dChar) :: char1d
+      integer :: ierr
 
-      integer :: i1, i2, i3, i4
 
-      integer, dimension(:), pointer :: super_int1d
-      integer, dimension(:,:), pointer :: super_int2d
-      real (kind=RKIND) :: super_real0d
-      real (kind=RKIND), dimension(:), pointer :: super_real1d
-      real (kind=RKIND), dimension(:,:), pointer :: super_real2d
-      real (kind=RKIND), dimension(:,:,:), pointer :: super_real3d
-      character (len=64) :: super_char0d
-      character (len=64), dimension(:), pointer :: super_char1d
+      call MPAS_readStream(input_obj % io_stream, 1, ierr)
 
-      integer :: i, k
 
-#include &quot;nondecomp_dims.inc&quot;
-
-      allocate(int1d % ioinfo)
-      allocate(int2d % ioinfo)
-      allocate(real0d % ioinfo)
-      allocate(real1d % ioinfo)
-      allocate(real2d % ioinfo)
-      allocate(real3d % ioinfo)
-      allocate(char0d % ioinfo)
-      allocate(char1d % ioinfo)
-
-
-#include &quot;io_input_fields.inc&quot;
-
-#include &quot;nondecomp_dims_dealloc.inc&quot;
-
    end subroutine mpas_read_and_distribute_fields
 
 
 
-   subroutine mpas_io_input_init(input_obj, dminfo)
+   subroutine mpas_io_input_init(input_obj, blocklist, dminfo)
  
       implicit none
 
       type (io_input_object), intent(inout) :: input_obj
+      type (block_type), intent(in) :: blocklist
       type (dm_info), intent(in) :: dminfo
  
-      include 'netcdf.inc'

       integer :: nferr
  

-#ifdef OFFSET64BIT
-      nferr = nf_open(trim(input_obj % filename), ior(NF_SHARE,NF_64BIT_OFFSET), input_obj % rd_ncid)
-#else
-      nferr = nf_open(trim(input_obj % filename), NF_SHARE, input_obj % rd_ncid)
-#endif
-
-      if (nferr /= NF_NOERR) then
+      call MPAS_createStream(input_obj % io_stream, trim(input_obj % filename), MPAS_IO_PNETCDF, MPAS_IO_READ, 1, nferr)
+      if (nferr /= MPAS_STREAM_NOERR) then
          write(0,*) ' '
          if (input_obj % stream == STREAM_RESTART) then
             write(0,*) 'Error opening restart file ''', trim(input_obj % filename), ''''
@@ -1331,9 +1339,9 @@
          write(0,*) ' '
          call mpas_dmpar_abort(dminfo)
       end if

-#include &quot;netcdf_read_ids.inc&quot;
 
+#include &quot;add_input_fields.inc&quot;
+
    end subroutine mpas_io_input_init
 
   
@@ -1345,7 +1353,7 @@
       character (len=*), intent(in) :: dimname
       integer, intent(out) :: dimsize
 
-#include &quot;get_dimension_by_name.inc&quot;
+!include &quot;get_dimension_by_name.inc&quot;
 
    end subroutine mpas_io_input_get_dimension
 
@@ -1358,24 +1366,8 @@
       character (len=*), intent(in) :: attname
       real (kind=RKIND), intent(out) :: attvalue
 
-      include 'netcdf.inc'
-
       integer :: nferr
 
-      if (RKIND == 8) then
-         nferr = nf_get_att_double(input_obj % rd_ncid, NF_GLOBAL, attname, attvalue)
-      else
-         nferr = nf_get_att_real(input_obj % rd_ncid, NF_GLOBAL, attname, attvalue)
-      end if
-      if (nferr /= NF_NOERR) then
-         write(0,*) 'Warning: Attribute '//trim(attname)//&amp;
-           ' not found in '//trim(input_obj % filename)
-         if (index(attname, 'sphere_radius') /= 0) then
-            write(0,*) '   Setting '//trim(attname)//' to 1.0'
-            attvalue = 1.0
-         end if
-      end if
-
    end subroutine mpas_io_input_get_att_real
 
    
@@ -1387,448 +1379,23 @@
       character (len=*), intent(in) :: attname
       character (len=*), intent(out) :: attvalue
 
-      include 'netcdf.inc'
-
       integer :: nferr
 
-      nferr = nf_get_att_text(input_obj % rd_ncid, NF_GLOBAL, attname, attvalue)
-      if (nferr /= NF_NOERR) then
-         write(0,*) 'Warning: Attribute '//trim(attname)//&amp;
-            ' not found in '//trim(input_obj % filename)
-         if (index(attname, 'on_a_sphere') /= 0) then
-            write(0,*) '   Setting '//trim(attname)//' to ''YES'''
-            attvalue = 'YES'
-         end if
-      end if
-
    end subroutine mpas_io_input_get_att_text
 
 
-   subroutine mpas_io_input_field0d_real(input_obj, field)

-      implicit none
+   subroutine mpas_exch_input_field_halos(domain, input_obj)
 
-      type (io_input_object), intent(in) :: input_obj      
-      type (field0dReal), intent(inout) :: field

-      include 'netcdf.inc'

-      integer :: nferr
-      integer :: varID
-      integer, dimension(1) :: start1, count1

-      start1(1) = 1
-      count1(1) = 1
-
-#include &quot;input_field0dreal.inc&quot;
-
-#if SINGLE_PRECISION
-      nferr = nf_get_vara_real(input_obj % rd_ncid, varID, start1, count1, field % scalar)
-#else
-      nferr = nf_get_vara_double(input_obj % rd_ncid, varID, start1, count1, field % scalar)
-#endif

-   end subroutine mpas_io_input_field0d_real
-
-
-   subroutine mpas_io_input_field1d_real(input_obj, field)

       implicit none
 
-      type (io_input_object), intent(in) :: input_obj      
-      type (field1dReal), intent(inout) :: field

-      include 'netcdf.inc'

-      integer :: nferr
-      integer :: varID
-      integer, dimension(1) :: start1, count1

-      start1(1) = field % ioinfo % start(1)
-      count1(1) = field % ioinfo % count(1)
+      type (domain_type), intent(inout) :: domain
+      type (io_input_object), intent(inout) :: input_obj
 
-      !
-      ! Special case: we may want to read the xtime variable across the
-      !   time dimension as a 1d array.
-      !
-      if (trim(field % ioinfo % fieldName) == 'xtime') then
-         varID = input_obj % rdVarIDxtime
-      end if

-#include &quot;input_field1dreal.inc&quot;
+#include &quot;exchange_input_field_halos.inc&quot;
 
-#if SINGLE_PRECISION
-      nferr = nf_get_vara_real(input_obj % rd_ncid, varID, start1, count1, field % array)
-#else
-      nferr = nf_get_vara_double(input_obj % rd_ncid, varID, start1, count1, field % array)
-#endif

-   end subroutine mpas_io_input_field1d_real
+   end subroutine mpas_exch_input_field_halos
 
 
-   subroutine mpas_io_input_field2d_real(input_obj, field)

-      implicit none
-
-      type (io_input_object), intent(in) :: input_obj      
-      type (field2dReal), intent(inout) :: field

-      include 'netcdf.inc'

-      integer :: nferr
-      integer :: varID
-      integer, dimension(2) :: start2, count2

-      start2(1) = field % ioinfo % start(1)
-      start2(2) = field % ioinfo % start(2)
-      count2(1) = field % ioinfo % count(1)
-      count2(2) = field % ioinfo % count(2)

-#include &quot;input_field2dreal.inc&quot;
-
-#if SINGLE_PRECISION
-      nferr = nf_get_vara_real(input_obj % rd_ncid, varID, start2, count2, field % array)
-#else
-      nferr = nf_get_vara_double(input_obj % rd_ncid, varID, start2, count2, field % array)
-#endif
-
-   end subroutine mpas_io_input_field2d_real
-
-
-   subroutine mpas_io_input_field3d_real(input_obj, field)

-      implicit none
-
-      type (io_input_object), intent(in) :: input_obj      
-      type (field3dReal), intent(inout) :: field

-      include 'netcdf.inc'

-      integer :: nferr
-      integer :: varID
-      integer, dimension(3) :: start3, count3

-      start3(1) = field % ioinfo % start(1)
-      start3(2) = field % ioinfo % start(2)
-      start3(3) = field % ioinfo % start(3)
-      count3(1) = field % ioinfo % count(1)
-      count3(2) = field % ioinfo % count(2)
-      count3(3) = field % ioinfo % count(3)

-#include &quot;input_field3dreal.inc&quot;
-
-#if SINGLE_PRECISION
-      nferr = nf_get_vara_real(input_obj % rd_ncid, varID, start3, count3, field % array)
-#else
-      nferr = nf_get_vara_double(input_obj % rd_ncid, varID, start3, count3, field % array)
-#endif
-
-   end subroutine mpas_io_input_field3d_real
-
-
-   subroutine mpas_io_input_field0d_real_time(input_obj, field)

-      implicit none
-
-      type (io_input_object), intent(in) :: input_obj      
-      type (field0dReal), intent(inout) :: field

-      include 'netcdf.inc'

-      integer :: nferr
-      integer :: varID
-      integer, dimension(1) :: start1, count1

-      start1(1) = input_obj % time
-      count1(1) = 1

-#include &quot;input_field0dreal_time.inc&quot;
-
-#if SINGLE_PRECISION
-      nferr = nf_get_vara_real(input_obj % rd_ncid, varID, start1, count1, field % scalar)
-#else
-      nferr = nf_get_vara_double(input_obj % rd_ncid, varID, start1, count1, field % scalar)
-#endif
-
-   end subroutine mpas_io_input_field0d_real_time
-
-
-   subroutine mpas_io_input_field1d_real_time(input_obj, field)

-      implicit none
-
-      type (io_input_object), intent(in) :: input_obj      
-      type (field1dReal), intent(inout) :: field

-      include 'netcdf.inc'

-      integer :: nferr
-      integer :: varID
-      integer, dimension(2) :: start2, count2

-      start2(1) = field % ioinfo % start(1)
-      start2(2) = input_obj % time
-      count2(1) = field % ioinfo % count(1)
-      count2(2) = 1

-#include &quot;input_field1dreal_time.inc&quot;
-
-#if SINGLE_PRECISION
-      nferr = nf_get_vara_real(input_obj % rd_ncid, varID, start2, count2, field % array)
-#else
-      nferr = nf_get_vara_double(input_obj % rd_ncid, varID, start2, count2, field % array)
-#endif
-
-   end subroutine mpas_io_input_field1d_real_time
-
-
-   subroutine mpas_io_input_field2d_real_time(input_obj, field)

-      implicit none
-
-      type (io_input_object), intent(in) :: input_obj      
-      type (field2dReal), intent(inout) :: field

-      include 'netcdf.inc'

-      integer :: nferr
-      integer :: varID
-      integer, dimension(3) :: start3, count3

-      start3(1) = field % ioinfo % start(1)
-      start3(2) = field % ioinfo % start(2)
-      start3(3) = input_obj % time
-      count3(1) = field % ioinfo % count(1)
-      count3(2) = field % ioinfo % count(2)
-      count3(3) = 1

-#include &quot;input_field2dreal_time.inc&quot;
-
-#if SINGLE_PRECISION
-      nferr = nf_get_vara_real(input_obj % rd_ncid, varID, start3, count3, field % array)
-#else
-      nferr = nf_get_vara_double(input_obj % rd_ncid, varID, start3, count3, field % array)
-#endif
-
-   end subroutine mpas_io_input_field2d_real_time
-
-
-   subroutine mpas_io_input_field3d_real_time(input_obj, field)

-      implicit none
-
-      type (io_input_object), intent(in) :: input_obj      
-      type (field3dReal), intent(inout) :: field

-      include 'netcdf.inc'

-      integer :: nferr
-      integer :: varID
-      integer, dimension(4) :: start4, count4

-      start4(1) = field % ioinfo % start(1)
-      start4(2) = field % ioinfo % start(2)
-      start4(3) = field % ioinfo % start(3)
-      start4(4) = input_obj % time
-      count4(1) = field % ioinfo % count(1)
-      count4(2) = field % ioinfo % count(2)
-      count4(3) = field % ioinfo % count(3)
-      count4(4) = 1

-#include &quot;input_field3dreal_time.inc&quot;
-
-#if SINGLE_PRECISION
-      nferr = nf_get_vara_real(input_obj % rd_ncid, varID, start4, count4, field % array)
-#else
-      nferr = nf_get_vara_double(input_obj % rd_ncid, varID, start4, count4, field % array)
-#endif
-
-   end subroutine mpas_io_input_field3d_real_time
-
-
-   subroutine mpas_io_input_field1d_integer(input_obj, field)

-      implicit none
-
-      type (io_input_object), intent(in) :: input_obj      
-      type (field1dInteger), intent(inout) :: field

-      include 'netcdf.inc'

-      integer :: nferr
-      integer :: varID
-      integer, dimension(1) :: start1, count1

-      start1(1) = field % ioinfo % start(1)
-      count1(1) = field % ioinfo % count(1)
-      
-#include &quot;input_field1dinteger.inc&quot;
-
-      nferr = nf_get_vara_int(input_obj % rd_ncid, varID, start1, count1, field % array)

-   end subroutine mpas_io_input_field1d_integer
-
-
-   subroutine mpas_io_input_field2d_integer(input_obj, field)

-      implicit none
-
-      type (io_input_object), intent(in) :: input_obj      
-      type (field2dInteger), intent(inout) :: field

-      include 'netcdf.inc'

-      integer :: nferr
-      integer :: varID
-      integer, dimension(2) :: start2, count2

-      start2(1) = field % ioinfo % start(1)
-      start2(2) = field % ioinfo % start(2)
-      count2(1) = field % ioinfo % count(1)
-      count2(2) = field % ioinfo % count(2)
-
-#include &quot;input_field2dinteger.inc&quot;
-
-      nferr = nf_get_vara_int(input_obj % rd_ncid, varID, start2, count2, field % array)
-
-   end subroutine mpas_io_input_field2d_integer
-
-
-   subroutine mpas_io_input_field1d_integer_time(input_obj, field)

-      implicit none
-
-      type (io_input_object), intent(in) :: input_obj      
-      type (field1dInteger), intent(inout) :: field

-      include 'netcdf.inc'

-      integer :: nferr
-      integer :: varID
-      integer, dimension(2) :: start2, count2

-      start2(1) = field % ioinfo % start(1)
-      start2(2) = input_obj % time
-      count2(1) = field % ioinfo % count(1)
-      count2(2) = 1

-#include &quot;input_field1dinteger_time.inc&quot;
-
-      nferr = nf_get_vara_int(input_obj % rd_ncid, varID, start2, count2, field % array)
-
-   end subroutine mpas_io_input_field1d_integer_time
-
-
-   subroutine mpas_io_input_field0d_char_time(input_obj, field)

-      implicit none
-
-      type (io_input_object), intent(in) :: input_obj      
-      type (field0dChar), intent(inout) :: field

-      include 'netcdf.inc'

-      integer :: nferr
-      integer :: varID
-      integer, dimension(2) :: start1, count1

-      start1(1) = 1
-      count1(1) = 64
-      start1(2) = input_obj % time
-      count1(2) = 1

-#include &quot;input_field0dchar_time.inc&quot;
-
-      nferr = nf_get_vara_text(input_obj % rd_ncid, varID, start1, count1, field % scalar)
-
-   end subroutine mpas_io_input_field0d_char_time
-
-
-   subroutine mpas_io_input_field1d_char_time(input_obj, field)

-      implicit none
-
-      type (io_input_object), intent(in) :: input_obj      
-      type (field1dChar), intent(inout) :: field

-      include 'netcdf.inc'

-      integer :: nferr
-      integer :: varID
-      integer, dimension(3) :: start2, count2

-      start2(1) = 1
-      start2(2) = field % ioinfo % start(1)
-      start2(3) = input_obj % time
-      count2(1) = 64
-      count2(2) = field % ioinfo % count(1)
-      count2(3) = 1

-#include &quot;input_field1dchar_time.inc&quot;
-
-      nferr = nf_get_vara_text(input_obj % rd_ncid, varID, start2, count2, field % array)
-
-   end subroutine mpas_io_input_field1d_char_time
-
-
-   subroutine mpas_io_input_field0d_char(input_obj, field)

-      implicit none
-
-      type (io_input_object), intent(in) :: input_obj      
-      type (field0dChar), intent(inout) :: field

-      include 'netcdf.inc'

-      integer :: nferr
-      integer :: varID
-      integer, dimension(2) :: start1, count1

-      start1(1) = 1
-      count1(1) = 64
-      start1(2) = 1
-      count1(2) = 1
-
-#include &quot;input_field0dchar.inc&quot;
-
-      nferr = nf_get_vara_text(input_obj % rd_ncid, varID, start1, count1, field % scalar)

-   end subroutine mpas_io_input_field0d_char
-
-
-   subroutine mpas_io_input_field1d_char(input_obj, field)

-      implicit none
-
-      type (io_input_object), intent(in) :: input_obj      
-      type (field1dChar), intent(inout) :: field

-      include 'netcdf.inc'

-      integer :: nferr
-      integer :: varID
-      integer, dimension(2) :: start1, count1

-      start1(1) = 1
-      count1(1) = 64
-      start1(2) = field % ioinfo % start(1)
-      count1(2) = field % ioinfo % count(1)
-
-      !
-      ! Special case: we may want to read the xtime variable across the
-      !   time dimension as a 1d array.
-      !
-      if (trim(field % ioinfo % fieldName) == 'xtime') then
-         varID = input_obj % rdVarIDxtime
-      end if

-#include &quot;input_field1dchar.inc&quot;
-
-      nferr = nf_get_vara_text(input_obj % rd_ncid, varID, start1, count1, field % array)

-   end subroutine mpas_io_input_field1d_char
-
-
    subroutine mpas_io_input_finalize(input_obj, dminfo)
  
       implicit none
@@ -1836,11 +1403,9 @@
       type (io_input_object), intent(inout) :: input_obj
       type (dm_info), intent(in) :: dminfo
 
-      include 'netcdf.inc'

       integer :: nferr
  
-      nferr = nf_close(input_obj % rd_ncid)
+      call MPAS_closeStream(input_obj % io_stream, nferr)
  
    end subroutine mpas_io_input_finalize
  

Modified: branches/omp_blocks/io/src/registry/gen_inc.c
===================================================================
--- branches/omp_blocks/io/src/registry/gen_inc.c        2012-04-10 22:08:32 UTC (rev 1766)
+++ branches/omp_blocks/io/src/registry/gen_inc.c        2012-04-10 22:13:32 UTC (rev 1767)
@@ -275,8 +275,8 @@
    fd = fopen(&quot;read_dims.inc&quot;, &quot;w&quot;);
    dim_ptr = dims;
    while (dim_ptr) {
-      if (dim_ptr-&gt;constant_value &lt; 0 &amp;&amp; !dim_ptr-&gt;namelist_defined &amp;&amp; !is_derived_dim(dim_ptr-&gt;name_in_code)) fortprintf(fd, &quot;      call mpas_io_input_get_dimension(input_obj, \'%s\', %s)</font>
<font color="red">&quot;, dim_ptr-&gt;name_in_file, dim_ptr-&gt;name_in_code);
-      else if (dim_ptr-&gt;constant_value &lt; 0 &amp;&amp; dim_ptr-&gt;namelist_defined &amp;&amp; !is_derived_dim(dim_ptr-&gt;name_in_code)) fortprintf(fd, &quot;      call mpas_io_input_get_dimension(input_obj, \'%s\', %s)</font>
<font color="blue">&quot;, dim_ptr-&gt;name_in_file, dim_ptr-&gt;name_in_file);
+      if (dim_ptr-&gt;constant_value &lt; 0 &amp;&amp; !dim_ptr-&gt;namelist_defined &amp;&amp; !is_derived_dim(dim_ptr-&gt;name_in_code)) fortprintf(fd, &quot;      call MPAS_io_inq_dim(inputHandle, \'%s\', %s, ierr)</font>
<font color="blue">&quot;, dim_ptr-&gt;name_in_file, dim_ptr-&gt;name_in_code);
+      else if (dim_ptr-&gt;constant_value &lt; 0 &amp;&amp; dim_ptr-&gt;namelist_defined &amp;&amp; !is_derived_dim(dim_ptr-&gt;name_in_code)) fortprintf(fd, &quot;      %s = %s</font>
<font color="gray">&quot;, dim_ptr-&gt;name_in_file, dim_ptr-&gt;name_in_code);
       dim_ptr = dim_ptr-&gt;next;
    }
 
@@ -1052,12 +1052,14 @@
    FILE * fd;
    char vtype[5];
    char fname[32];
+   char super_array[1024];
    char struct_deref[1024];
    char * cp1, * cp2;
    int i, j;
    int ivtype;
 
 
+#ifdef LEGACY_CODE
    /*
     *  Generate declarations of IDs belonging in io_input_object
     */
@@ -1087,6 +1089,7 @@
 
    fclose(fd);
 
+
    /*
     *  Definitions of read bounds and exchange lists for non-decomposed fields
     */
@@ -1273,7 +1276,6 @@
    fclose(fd);
 
 
-#ifdef LEGACY_CODE
    /*
     *  Definitions of read bounds and exchange lists for non-decomposed fields
     */
@@ -1477,7 +1479,6 @@
    }
 
    fclose(fd);
-#endif
    
 
    /*
@@ -1746,6 +1747,7 @@
    }
 
    fclose(fd);
+#endif
 
 
    /*
@@ -1759,14 +1761,93 @@
       while (var_list_ptr) {
          var_ptr = var_list_ptr-&gt;var;
 
-         if (group_ptr-&gt;vlist-&gt;var-&gt;ntime_levs &gt; 1)
-            snprintf(struct_deref, 1024, &quot;block %% %s %% time_levs(1) %% %s&quot;, group_ptr-&gt;name, group_ptr-&gt;name);
+         if (var_ptr-&gt;ntime_levs &gt; 1)
+            snprintf(struct_deref, 1024, &quot;blocklist %% %s %% time_levs(1) %% %s&quot;, group_ptr-&gt;name, group_ptr-&gt;name);
          else
-            snprintf(struct_deref, 1024, &quot;block %% %s&quot;, group_ptr-&gt;name);
+            snprintf(struct_deref, 1024, &quot;blocklist %% %s&quot;, group_ptr-&gt;name);
+         
+         if (strncmp(var_ptr-&gt;super_array, &quot;-&quot;, 1024) != 0) {
+            fortprintf(fd, &quot;      if ((%s %% %s %% ioinfo %% input .and. input_obj %% stream == STREAM_INPUT) .or. &amp;</font>
<font color="blue">&quot;, struct_deref, var_ptr-&gt;super_array);
+            fortprintf(fd, &quot;          (%s %% %s %% ioinfo %% restart .and. input_obj %% stream == STREAM_RESTART) .or. &amp;</font>
<font color="blue">&quot;, struct_deref, var_ptr-&gt;super_array);
+            fortprintf(fd, &quot;          (%s %% %s %% ioinfo %% sfc .and. input_obj %% stream == STREAM_SFC)) then</font>
<font color="blue">&quot;, struct_deref, var_ptr-&gt;super_array);
+            memcpy(super_array, var_ptr-&gt;super_array, 1024);
+            fortprintf(fd, &quot;         write(0,*) \'adding input field %s\'</font>
<font color="blue">&quot;, var_ptr-&gt;super_array);
+            fortprintf(fd, &quot;         call MPAS_streamAddField(input_obj %% io_stream, %s %% %s, nferr)</font>
<font color="blue">&quot;, struct_deref, var_ptr-&gt;super_array);
+            while (var_list_ptr &amp;&amp; strncmp(super_array, var_list_ptr-&gt;var-&gt;super_array, 1024) == 0) {
+               var_list_ptr = var_list_ptr-&gt;next;
+            }
+         }
+         else {
+            fortprintf(fd, &quot;      if ((%s %% %s %% ioinfo %% input .and. input_obj %% stream == STREAM_INPUT) .or. &amp;</font>
<font color="blue">&quot;, struct_deref, var_ptr-&gt;name_in_code);
+            fortprintf(fd, &quot;          (%s %% %s %% ioinfo %% restart .and. input_obj %% stream == STREAM_RESTART) .or. &amp;</font>
<font color="blue">&quot;, struct_deref, var_ptr-&gt;name_in_code);
+            fortprintf(fd, &quot;          (%s %% %s %% ioinfo %% sfc .and. input_obj %% stream == STREAM_SFC)) then</font>
<font color="blue">&quot;, struct_deref, var_ptr-&gt;name_in_code);
+            fortprintf(fd, &quot;         write(0,*) \'adding input field %s\'</font>
<font color="blue">&quot;, var_ptr-&gt;name_in_code);
+            fortprintf(fd, &quot;         call MPAS_streamAddField(input_obj %% io_stream, %s %% %s, nferr)</font>
<font color="red">&quot;, struct_deref, var_ptr-&gt;name_in_code);
+         }
    
-         fortprintf(fd, &quot;   call MPAS_streamAddField(input_obj %% io_stream, %s %% %s, ierr)</font>
<font color="blue">&quot;, struct_deref, var_ptr-&gt;name_in_code);
+         fortprintf(fd, &quot;      end if</font>
<font color="black"></font>
<font color="red">&quot;);
+
+         if (var_list_ptr) var_list_ptr = var_list_ptr-&gt;next;
+      }
+      group_ptr = group_ptr-&gt;next;
+   }
+
+   fclose(fd);
+
+
+   /*
+    * MGD NEW CODE
+    */
+   fd = fopen(&quot;exchange_input_field_halos.inc&quot;, &quot;w&quot;);
+
+   group_ptr = groups;
+   while (group_ptr) {
+      var_list_ptr = group_ptr-&gt;vlist;
+      while (var_list_ptr) {
+         var_ptr = var_list_ptr-&gt;var;
+
+         dimlist_ptr = var_ptr-&gt;dimlist;
+         i = 1;
+         while (dimlist_ptr) {
+            if (i == var_ptr-&gt;ndims) { 
+               if (!strncmp(dimlist_ptr-&gt;dim-&gt;name_in_file, &quot;nCells&quot;, 1024) ||
+                   !strncmp(dimlist_ptr-&gt;dim-&gt;name_in_file, &quot;nEdges&quot;, 1024) ||
+                   !strncmp(dimlist_ptr-&gt;dim-&gt;name_in_file, &quot;nVertices&quot;, 1024)) {
    
-         var_list_ptr = var_list_ptr-&gt;next;
+                  if (var_ptr-&gt;ntime_levs &gt; 1)
+                     snprintf(struct_deref, 1024, &quot;domain %% blocklist %% %s %% time_levs(1) %% %s&quot;, group_ptr-&gt;name, group_ptr-&gt;name);
+                  else
+                     snprintf(struct_deref, 1024, &quot;domain %% blocklist %% %s&quot;, group_ptr-&gt;name);
+                  
+                  if (strncmp(var_ptr-&gt;super_array, &quot;-&quot;, 1024) != 0) {
+                     fortprintf(fd, &quot;      if ((%s %% %s %% ioinfo %% input .and. input_obj %% stream == STREAM_INPUT) .or. &amp;</font>
<font color="blue">&quot;, struct_deref, var_ptr-&gt;super_array);
+                     fortprintf(fd, &quot;          (%s %% %s %% ioinfo %% restart .and. input_obj %% stream == STREAM_RESTART) .or. &amp;</font>
<font color="blue">&quot;, struct_deref, var_ptr-&gt;super_array);
+                     fortprintf(fd, &quot;          (%s %% %s %% ioinfo %% sfc .and. input_obj %% stream == STREAM_SFC)) then</font>
<font color="blue">&quot;, struct_deref, var_ptr-&gt;super_array);
+                     memcpy(super_array, var_ptr-&gt;super_array, 1024);
+                     fortprintf(fd, &quot;         write(0,*) \'exchange halo for %s\'</font>
<font color="blue">&quot;, var_ptr-&gt;super_array);
+                     fortprintf(fd, &quot;         call mpas_dmpar_exch_halo_field(%s %% %s)</font>
<font color="blue">&quot;, struct_deref, var_ptr-&gt;super_array);
+                     while (var_list_ptr &amp;&amp; strncmp(super_array, var_list_ptr-&gt;var-&gt;super_array, 1024) == 0) {
+                        var_list_ptr = var_list_ptr-&gt;next;
+                     }
+                  }
+                  else {
+                     fortprintf(fd, &quot;      if ((%s %% %s %% ioinfo %% input .and. input_obj %% stream == STREAM_INPUT) .or. &amp;</font>
<font color="blue">&quot;, struct_deref, var_ptr-&gt;name_in_code);
+                     fortprintf(fd, &quot;          (%s %% %s %% ioinfo %% restart .and. input_obj %% stream == STREAM_RESTART) .or. &amp;</font>
<font color="blue">&quot;, struct_deref, var_ptr-&gt;name_in_code);
+                     fortprintf(fd, &quot;          (%s %% %s %% ioinfo %% sfc .and. input_obj %% stream == STREAM_SFC)) then</font>
<font color="blue">&quot;, struct_deref, var_ptr-&gt;name_in_code);
+                     fortprintf(fd, &quot;         write(0,*) \'exchange halo for %s\'</font>
<font color="blue">&quot;, var_ptr-&gt;name_in_code);
+                     fortprintf(fd, &quot;         call mpas_dmpar_exch_halo_field(%s %% %s)</font>
<font color="blue">&quot;, struct_deref, var_ptr-&gt;name_in_code);
+                  }
+            
+                  fortprintf(fd, &quot;      end if</font>
<font color="black"></font>
<font color="gray">&quot;);
+   
+               }
+            }
+   
+            i++;
+            dimlist_ptr = dimlist_ptr -&gt; next;
+         }
+
+         if (var_list_ptr) var_list_ptr = var_list_ptr-&gt;next;
       }
       group_ptr = group_ptr-&gt;next;
    }
@@ -1774,6 +1855,7 @@
    fclose(fd);
 
 
+#ifdef LEGACY_CODE
    /*
     *  Generate NetCDF reads of dimension and variable IDs
     */
@@ -1936,6 +2018,7 @@
    
       fclose(fd);
    } 
+#endif
    
 }
 

</font>
</pre>