anchor.diff:          call check( nf90_close(ncid) )
       end subroutine io_write2d
       
  -    ! simple error handling for common netcdf file errors
  +    !>------------------------------------------------------------
  +    !! Simple error handling for common netcdf file errors
  +    !!
  +    !! If status does not equal nf90_noerr, then print an error message and STOP 
  +    !! the entire program. 
  +    !!
  +    !! @param   status  integer return code from nc_* routines
  +    !! @param   extra   OPTIONAL string with extra context to print in case of an error
  +    !!
  +    !!------------------------------------------------------------
       subroutine check(status,extra)
           implicit none
           integer, intent ( in) :: status
           character(len=*), optional, intent(in) :: extra
  -    
  +        
  +        ! check for errors
           if(status /= nf90_noerr) then 
  +            ! print a useful message
               print *, trim(nf90_strerror(status))
               if(present(extra)) then
  +                ! print any optionally provided context
                   print*, trim(extra)
               endif
  +            ! STOP the program execute
               stop "Stopped"
           end if
       end subroutine check  
       
  -    ! Find an available file unit number.
  -    ! LUN_MIN and LUN_MAX define the range of possible LUNs to check.
  -    ! The UNIT value is returned by the function, and also by the optional
  -    ! argument. This allows the function to be used directly in an OPEN
  -    ! statement, and optionally save the result in a local variable.
  -    ! If no units are available, -1 is returned.
  -    ! Newer versions of fortran can do this automatically, but this keeps one thing
  -    ! a little more backwards compatible
  +    !>------------------------------------------------------------
  +    !! Find an available file unit number.
  +    !!
  +    !! LUN_MIN and LUN_MAX define the range of possible LUNs to check.
  +    !! The UNIT value is returned by the function, and also by the optional
  +    !! argument. This allows the function to be used directly in an OPEN
  +    !! statement, and optionally save the result in a local variable.
  +    !! If no units are available, -1 is returned.
  +    !! Newer versions of fortran can do this automatically, but this keeps one thing
  +    !! a little more backwards compatible
  +    !!
  +    !! @param[out]  unit    OPTIONAL integer to store the file logical unit number
  +    !! @retval      integer a file logical unit number
  +    !!
  +    !!------------------------------------------------------------
       integer function io_newunit(unit)
           implicit none
           integer, intent(out), optional :: unit
  
positive.code:         
      end subroutine io_read6d
  
      ! Reads in a variable from a netcdf file, allocating memory in data_in for it.  
      ! if extradim is provided specifies this index for any extra dimensions (dims>3)
      !   e.g. we may only want one time slice from a 3d variable
      subroutine io_read3d(filename,varname,data_in,extradim)
          implicit none
          ! This is the name of the data_in file and variable we will read. 
  
negative.code:         ! Read the dimension lengths
          call io_getdims(filename,varname,diminfo)
          allocate(data_in(diminfo(2),diminfo(3),diminfo(4),diminfo(5),diminfo(6),diminfo(7)))
          ! Open the file. NF90_NOWRITE tells netCDF we want read-only access to the file.
          call check(nf90_open(filename, NF90_NOWRITE, ncid),filename)
          ! Get the varid of the data_in variable, based on its name.
          call check(nf90_inq_varid(ncid, varname, varid),trim(filename)//":"//trim(varname))
  