anchor.diff:          if (options%time_varying_z) then
               ! read in the updated vertical coordinate
               if (allocated(newbc%z)) deallocate(newbc%z)
  -            call io_read3d(file_list(curfile), options%zvar, newbc%z, curstep)
  +            call io_read(file_list(curfile), options%zvar, newbc%z, curstep)
               nx=size(newbc%z,1)
               ny=size(newbc%z,2)
               nz=size(newbc%z,3)
               if (trim(options%zbvar)/="") then
  -                call io_read3d(file_list(curfile),options%zbvar, zbase, curstep)
  +                call io_read(file_list(curfile),options%zbvar, zbase, curstep)
                   newbc%z=newbc%z+zbase
                   deallocate(zbase)
               endif
  
positive.code:         nz=size(domain%u,2)
  
  !       load low-res U data
          call io_read3d(filename,options%uvar,extra_data,curstep)
          domain%u=sum(extra_data(:,:,:nz))/size(extra_data(:,:,:nz))
          deallocate(extra_data)
  
  !       load low-res V data
          call io_read3d(filename,options%vvar,extra_data,curstep)
          domain%v=sum(extra_data(:,:,:nz))/size(extra_data(:,:,:nz))
          deallocate(extra_data)
  
  
negative.code:         real,dimension(:,:),allocatable :: inputdata
  
  !       Read the data in
          call io_read(filename,varname,inputdata,curstep)
  !       interpolate data onto the high resolution grid
          call geo_interp2d(highres,inputdata,geolut)
          deallocate(inputdata)
  