!
      subroutine mksatcdf(data,locs,itime,nt,nz,nf,
     +                    flnm_cdf,f,h,ncid,imode,ier,satfile)
!
      use proc
      use hist,only: history,hdr
      use fields,only: field
      use input,only: mtimes
      implicit none
      include 'netcdf.inc'
!
! Args:
      integer,intent(in) :: nz,nt,nf,itime,imode
      integer,intent(inout) ::  ncid		! netcdf file id
      character(len=120),intent(in):: satfile
      type(field),intent(in) :: f(nf)
      type(history),intent(in) :: h
      real,intent(in) :: data(nz,nf),locs(9)
      character(len=120),intent(out) :: flnm_cdf	! netcdf file name
      integer,intent(out) :: ier
!
! Locals:
      character(len=6)   :: chnlev
      character(len=8)   :: charmtime,charmtime1,fname
      character(len=80)  :: char80
      character(len=120) :: char120
!
! Locals that must be saved for calls after 1st time:
      integer,save :: 
     +  id_time,id_mtime,		! dim ids for time and mtime
     +  idv_time,idv_mtime,idv_histvols,! time-dependent non-dim vars
     +  idv_date,idv_adflag,idv_angle,idv_yyddd,idv_hp,idv_cp,idv_by,
     +  idv_f107d,idv_f107a,idv_track,
     +  idv_lon,idv_lat,idv_dnflag,idv_lst
!
! Locals that need not be saved:
      integer :: 
     +  istat,i,ii,ipos,iier,		! status, loop vars, char pos
     +  id_lon,id_lat,id_lev,id_vert, 	! dimension id's
     +  id_hvol_len,
     +  ids1(1),ids2(2),ids4(4),  	! vectors of dim id's
!#     +  idv_lon,idv_lat,idv_lev,	! dim variable id's
     +  idv_lev,			! dim variable id's
     +  idv_flds(nf),			! field variable id's
     +  nflds,				! number of fields on netcdf file
     +  nnlev,				! temp for number of levels
!#     +  mins,				! total model time in minutes
     +  start_4d(4),count_4d(4),	! start loc and count for 4d vars
     +  start_2d(2),count_2d(2)		! start loc and count for 2d vars
!/      real,allocatable :: var(:,:,:)	! temp var
      real :: mins, lst
      character(len=8),allocatable :: field_names(:)	! field names
!
! Error init:
      ier = 0
!
! If this is not the 1st time, just grow the variables along the
! unlimited time dimension:

      if (imode > 0) goto 100
!
! If this is the 1st time iteration, create file, define variables, 
!   dimensions, and attributes:
!
! Construct file name: 
! If only 1 time, flnm_cdf = "1sthistvol_mtime.nc"
! If multiple times, flnm_cdr = "1sthistvol_1stmtime-lastmtime.nc"
!
      call tail(h%mssvol,flnm_cdf)
      write(charmtime,"(i3.3,i2.2,i2.2)") h%mtime
      if (ntimes == 1) then
        flnm_cdf = trim(flnm_cdf)//'_'//trim(charmtime)//'.nc'
      else
        write(charmtime1,"(i3.3,i2.2,i2.2)") mtimes(1:3,ntimes)
        flnm_cdf = trim(flnm_cdf)//'_'//trim(charmtime)//'-'
     +  //trim(charmtime1)//'.nc'
      endif
!
! Create a netcdf file (file id ncid is returned by nf_create):
!
      istat = nf_create(flnm_cdf,NF_CLOBBER,ncid)
      if (istat /= NF_NOERR) then
        write(char120,"('Error return from nf_create for netcdf ',
     +    'file ',a)") trim(flnm_cdf)
!       call handle_ncerr(istat,char120,ier)
        call handle_ncerr(istat,char120)
      else
        write(6,"('Created netcdf file ',a)") trim(flnm_cdf)
      endif
!
! Define dimensions (field vars will be (time,lev,lat,lon) in
! CDL notation): 
!
!/      istat = nf_def_dim(ncid,"lon",nlon,id_lon)
!/      if (istat /= NF_NOERR) call handle_ncerr(istat,
!/     +  'Error defining longitude dimension')
!/      istat = nf_def_dim(ncid,"lat",nlat,id_lat)
!/      if (istat /= NF_NOERR) call handle_ncerr(istat,
!/     +  'Error defining latitude dimension')
      istat = nf_def_dim(ncid,"lev",npress,id_lev)
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining levels dimension')
      istat = nf_def_dim(ncid,"time",NF_UNLIMITED,id_time)
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining time dimension')
!
! Define dimension variables and attributes:
!
! Longitudes:
!/      ids1(1) = id_lon
!/      istat = nf_def_var(ncid,"lon",NF_FLOAT,1,ids1,idv_lon)
!/      if (istat /= NF_NOERR) call handle_ncerr(istat,
!/     +  'Error defining longitude dimension variable')
!/      write(char80,"('geographic longitude (-west, +east)')")
!/      istat = nf_put_att_text(ncid,id_lon,"long_name",
!/     +  len_trim(char80),trim(char80))
!/      if (istat /= NF_NOERR) call handle_ncerr(istat,
!/     +  'Error defining long_name of longitude dimension variable')
!/      istat = nf_put_att_text(ncid,id_lon,"units",7,'degrees')
!/      if (istat /= NF_NOERR) call handle_ncerr(istat,
!/     +  'Error defining units of longitude dimension variable')
!
! Latitudes:
!/      ids1(1) = id_lat
!/      istat = nf_def_var(ncid,"lat",NF_FLOAT,1,ids1,idv_lat)
!/      if (istat /= NF_NOERR) call handle_ncerr(istat,
!/     +  'Error defining latitude dimension variable')
!/      write(char80,"('geographic latitude (-south, +north)')")
!/      istat = nf_put_att_text(ncid,id_lat,"long_name",
!/     +  len_trim(char80),trim(char80))
!/      if (istat /= NF_NOERR) call handle_ncerr(istat,
!/     +  'Error defining long_name of latitude dimension variable')
!/      istat = nf_put_att_text(ncid,id_lat,"units",7,'degrees')
!/      if (istat /= NF_NOERR) call handle_ncerr(istat,
!/     +  'Error defining units of latitude dimension variable')
!
! Pressures:
      ids1(1) = id_lev
      istat = nf_def_var(ncid,"LEV",NF_FLOAT,1,ids1,idv_lev)
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining pressure (lev) dimension variable')
      write(char80,"('LOG PRESSURE')")
      istat = nf_put_att_text(ncid,id_lev,"long_name",
     +  len_trim(char80),trim(char80))
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining long_name of pressure dimension variable')
      istat = nf_put_att_text(ncid,id_lev,"units",8,'LN(P0/P)')
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining units of pressure dimension variable')
!
! Time (total model time in hours):
      ids1(1) = id_time
      istat = nf_def_var(ncid,"TIME",NF_FLOAT,1,ids1,idv_time)
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining time dimension variable')
      write(char80,"('UNIVERSAL TIME')")
      istat = nf_put_att_text(ncid,id_time,"long_name",
     +  len_trim(char80),trim(char80))
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining long_name of time dimension variable')
      istat = nf_put_att_text(ncid,id_time,"units",5,'HOURS')
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining units of time dimension variable')
!
! Define local solar time (from model time)
      istat = nf_def_var(ncid,"LST",NF_FLOAT,1,ids1,idv_lst)
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining model lst variable')
      write(char80,"('LOCAL SOLAR TIME')")
      istat = nf_put_att_text(ncid,idv_lst,"long_name",
     +  len_trim(char80),trim(char80))
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining long_name attribute of lst variable')
      istat = nf_put_att_text(ncid,idv_lst,"units",5,'HOURS')
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining units of lst variable')
!
! Define time-dependent (non-dim) variables:
      ids1(1) = id_time	! for 1d time-unlimited vars
!
! Define model time var (3,unlimited):
      istat = nf_def_dim(ncid,"mtime",3,id_mtime)
      ids2(1) = id_mtime
      ids2(2) = id_time
      istat = nf_def_var(ncid,"MTIME",NF_SHORT,2,ids2,idv_mtime)
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining mtime variable')
      write(char80,"('MODEL TIME')")
      istat = nf_put_att_text(ncid,idv_mtime,"long_name",
     +  len_trim(char80),trim(char80))
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining long_name attribute of mtime variable')
      istat = nf_put_att_text(ncid,idv_mtime,"units",17,
     +  'DAY, HOUR, MINUTE')
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining units of mtime variable')
!
! Define history volumes var (unlimited):
! (must define a character position dimension for char var)
      istat = nf_def_dim(ncid,"hvol_len",len(h%mssvol),id_hvol_len)
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining id_hvol_len dimension')
      ids2(1) = id_hvol_len	! character position dimension
      ids2(2) = id_time
      istat = nf_def_var(ncid,"HISTVOLS",NF_CHAR,2,ids2,idv_histvols)
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining histvols variable')
      write(char80,"('mss paths to history volumes')")
      istat = nf_put_att_text(ncid,idv_histvols,"long_name",
     +  len_trim(char80),trim(char80))
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining long_name attribute of histvols variable')

!
! Define date var (unlimited):
      istat = nf_def_var(ncid,"DATE",NF_INT,1,ids1,idv_date)
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining date variable')
      write(char80,"('DATE')")
      istat = nf_put_att_text(ncid,idv_date,"long_name",
     +  len_trim(char80),trim(char80))
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining long_name attribute of date variable')
      istat = nf_put_att_text(ncid,idv_date,"units",7,'YYYYDDD')
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining units of date variable')
!
! Define adflag var (unlimited):
      istat = nf_def_var(ncid,"ADFLAG",NF_INT,1,ids1,idv_adflag)
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining adflag variable')
      write(char80,"('ORBIT NODE')")
      istat = nf_put_att_text(ncid,idv_adflag,"long_name",
     +  len_trim(char80),trim(char80))
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining long_name attribute of adflag variable')
      istat = nf_put_att_text(ncid,idv_adflag,"units",14,
     +  'ASC(0), DSC(1)')
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining units of adflag variable')
!
! Define dnflag var (unlimited):
      istat = nf_def_var(ncid,"DNFLAG",NF_INT,1,ids1,idv_dnflag)
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining dnflag variable')
      write(char80,"('DAY/NIGHT')")
      istat = nf_put_att_text(ncid,idv_dnflag,"long_name",
     +  len_trim(char80),trim(char80))
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining long_name attribute of dnflag variable')
      istat = nf_put_att_text(ncid,idv_dnflag,"units",16,
     +  'DAY(0), NIGHT(1)')
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining units of dnflag variable')
!
! Define angle var (unlimited):
      istat = nf_def_var(ncid,"WCFLAG",NF_INT,1,ids1,idv_angle)
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining wcflag variable')
      write(char80,"('TIDI WARM/COLD SIDE')")
      istat = nf_put_att_text(ncid,idv_angle,"long_name",
     +  len_trim(char80),trim(char80))
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining long_name attribute of wcflag variable')
      istat = nf_put_att_text(ncid,idv_angle,"units",16,
     +  'COLD(0), WARM(1)')
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining units of wcflag variable')
!
! Define track angle var (unlimited): TIDI only
      if(index(satfile, 'tidi') .gt. 0) then
      istat = nf_def_var(ncid,"TRACK",NF_FLOAT,1,ids1,idv_track)
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining track variable')
      write(char80,"('TIDI TRACK ANGLE')")
      istat = nf_put_att_text(ncid,idv_track,"long_name",
     +  len_trim(char80),trim(char80))
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining long_name attribute of track variable')
      istat = nf_put_att_text(ncid,idv_track,"units",3,
     +  'DEG')
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining units of track variable')  
      endif 
!
! Define yyddd var (unlimited):
!      istat = nf_def_var(ncid,"YYDDD",NF_INT,1,ids1,idv_yyddd)
!      if (istat /= NF_NOERR) call handle_ncerr(istat,
!     +  'Error defining yyddd variable')
!      write(char80,"('year and day of year (yyddd)')")
!      istat = nf_put_att_text(ncid,idv_yyddd,"long_name",
!     +  len_trim(char80),trim(char80))
!      if (istat /= NF_NOERR) call handle_ncerr(istat,
!     +  'Error defining long_name attribute of yyddd variable')
!
! Define longitude var (unlimited):
      istat = nf_def_var(ncid,"LON",NF_FLOAT,1,ids1,idv_lon)
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining lon variable')
      write(char80,"('LONGITUDE')")
      istat = nf_put_att_text(ncid,idv_lon,"long_name",
     +  len_trim(char80),trim(char80))
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining long_name attribute of lon variable')
      istat = nf_put_att_text(ncid,idv_lon,"units",3,'DEG')
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining units of lon variable')
!
! Define latitude var (unlimited):
      istat = nf_def_var(ncid,"LAT",NF_FLOAT,1,ids1,idv_lat)
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining lat variable')
      write(char80,"('LATITUDE')")
      istat = nf_put_att_text(ncid,idv_lat,"long_name",
     +  len_trim(char80),trim(char80))
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining long_name attribute of lat variable')
      istat = nf_put_att_text(ncid,idv_lat,"units",3,'DEG')
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error defining units of lat variable')
!
! Define f107d var (unlimited):
!      istat = nf_def_var(ncid,"F107D",NF_FLOAT,1,ids1,idv_f107d)
!      if (istat /= NF_NOERR) call handle_ncerr(istat,
!     +  'Error defining f107d variable')
!      write(char80,"('DAILY F10.7 CM SOLAR FLUX')")
!      istat = nf_put_att_text(ncid,idv_f107d,"long_name",
!     +  len_trim(char80),trim(char80))
!      if (istat /= NF_NOERR) call handle_ncerr(istat,
!     +  'Error defining long_name attribute of f107d variable')
!      istat = nf_put_att_text(ncid,idv_f107d,"units",1,' ')
!      if (istat /= NF_NOERR) call handle_ncerr(istat,
!     +  'Error defining units of f107d variable')
!
! Define f107a var (unlimited):
!      istat = nf_def_var(ncid,"F107A",NF_FLOAT,1,ids1,idv_f107a)
!      if (istat /= NF_NOERR) call handle_ncerr(istat,
!     +  'Error defining f107a variable')
!      write(char80,"('80-DAY AVERAGE F10.7 CM SOLAR FLUX')")
!      istat = nf_put_att_text(ncid,idv_f107a,"long_name",
!     +  len_trim(char80),trim(char80))
!      if (istat /= NF_NOERR) call handle_ncerr(istat,
!     +  'Error defining long_name attribute of f107a variable')
!      istat = nf_put_att_text(ncid,idv_f107a,"units",1,' ')
!      if (istat /= NF_NOERR) call handle_ncerr(istat,
!     +  'Error defining units of f107a variable')
!
! Define hp var (unlimited):
!      istat = nf_def_var(ncid,"HP",NF_FLOAT,1,ids1,idv_hp)
!      if (istat /= NF_NOERR) call handle_ncerr(istat,
!     +  'Error defining hp variable')
!      write(char80,"('HEMISPHERIC POWER')")
!      istat = nf_put_att_text(ncid,idv_hp,"long_name",
!     +  len_trim(char80),trim(char80))
!      if (istat /= NF_NOERR) call handle_ncerr(istat,
!     +  'Error defining long_name attribute of hp variable')
!      istat = nf_put_att_text(ncid,idv_hp,"units",2,'GW')
!      if (istat /= NF_NOERR) call handle_ncerr(istat,
!     +  'Error defining units of hp variable')
!
! Define cp var (unlimited):
!      istat = nf_def_var(ncid,"CP",NF_FLOAT,1,ids1,idv_cp)
!      if (istat /= NF_NOERR) call handle_ncerr(istat,
!     +  'Error defining cp variable')
!      write(char80,"('CROSS-TAIL POTENTIAL')")
!      istat = nf_put_att_text(ncid,idv_cp,"long_name",
!     +  len_trim(char80),trim(char80))
!      if (istat /= NF_NOERR) call handle_ncerr(istat,
!     +  'Error defining long_name attribute of cp variable')
!      istat = nf_put_att_text(ncid,idv_cp,"units",5,'VOLTS')
!      if (istat /= NF_NOERR) call handle_ncerr(istat,
!     +  'Error defining units of cp variable')
!
! Define byimf var (unlimited):
!      istat = nf_def_var(ncid,"BYIMF",NF_FLOAT,1,ids1,idv_by)
!      if (istat /= NF_NOERR) call handle_ncerr(istat,
!     +  'Error defining byimf variable')
!      write(char80,"('IMF BY COMPONENT')")
!      istat = nf_put_att_text(ncid,idv_by,"long_name",
!     +  len_trim(char80),trim(char80))
!      if (istat /= NF_NOERR) call handle_ncerr(istat,
!     +  'Error defining long_name attribute of byimf variable')
!      istat = nf_put_att_text(ncid,idv_by,"units",1,' ')
!      if (istat /= NF_NOERR) call handle_ncerr(istat,
!     +  'Error defining units of by variable')
!
! Define 4d field vars (always include Z):
!
!/      ids4(1) = id_lon
!/      ids4(2) = id_lat
!/      ids4(4) = id_time
      ids4(2) = id_time
      nnlev = 0
      nflds = 0
      floop1: do i=1,nf
        if ((.not.f(i)%requested.or..not.associated(f(i)%data)).and.
     +       f(i)%fname8 /= 'Z       ') cycle floop1
        nflds = nflds+1
        if (f(i)%nlev == h%nzp) then
          ids4(1) = id_lev
!
! field vert res different than history field (e.g., ht-independent,
! or derived oh):
!
        elseif (i==1.or.(i>1.and.f(i)%nlev /= f(i-1)%nlev)) then
          nnlev = nnlev+1
          write(chnlev,"('nlev',i2.2)") nnlev
          istat = nf_def_dim(ncid,chnlev,f(i)%nlev,id_vert)
          ids4(1) = id_vert
        endif
        fname = f(i)%fname8
!
! Netcdf does not allow '+' in a var field name, so change to 'P'
! if necessary (e.g., O2+ will be O2P on the netcdf file):
!
        ipos = index(fname,'+')
        if (ipos > 0) fname(ipos:ipos) = 'P'
        ipos = index(fname,'/')
        if (ipos > 0) fname(ipos:ipos) = '_'

! Define the field var:
!/        istat = nf_def_var(ncid,fname,NF_FLOAT,4,ids4,idv_flds(i))
        istat = nf_def_var(ncid,fname,NF_FLOAT,2,ids4,idv_flds(i))
        if (istat /= NF_NOERR) then
          write(char120,"('Error return from nf_def_var for ',
     +      'field var ',a,': i=',i2)") f(i)%fname8,i
          call handle_ncerr(istat,char120)
        endif
!
! Assign attributes to field variables:
!
! Short and long field names, and units:
        istat = nf_put_att_text(ncid,idv_flds(i),"long_name",
     +    len_trim(f(i)%fname56),trim(f(i)%fname56))
        istat = nf_put_att_text(ncid,idv_flds(i),"units",
     +    len_trim(f(i)%units),trim(f(i)%units))
        if (istat /= NF_NOERR) then
          write(char120,"('Error return from nf_put_att_text',
     +      ' of long_name and units for field ',a)") f(i)%fname8
          call handle_ncerr(istat,char120)
        endif
!
! Field "type" (DENSITY, VECTOR, SCALAR, etc):
        istat = nf_put_att_text(ncid,idv_flds(i),"type",
     +    len(f(i)%type),f(i)%type)
        if (istat /= NF_NOERR) then
          write(char120,"('Error return from nf_put_att_text',
     +      ' of type ',a,' for field ',a)") f(i)%type,f(i)%fname8
          call handle_ncerr(istat,char120)
        endif
!
! Vertical coords for each field (usually zp as in model, but
! certain fields may have their own height coords, e.g. OH emissions, 
! etc) (vtype indicated pressure or height coords)
        istat = nf_put_att_int(ncid,idv_flds(i),"nlev",NF_INT,
     +    1,f(i)%nlev)
        istat = nf_put_att_real(ncid,idv_flds(i),"blev",NF_FLOAT,
     +    1,f(i)%lev(1))
        istat = nf_put_att_real(ncid,idv_flds(i),"dlev",NF_FLOAT,
     +    1,f(i)%dlev)
        istat = nf_put_att_text(ncid,idv_flds(i),"vtype",
     +    len(f(i)%vtype),f(i)%vtype)
        if (istat /= NF_NOERR) then
          write(char120,"('Error preparing vertical coord attributes',
     +      'for field ',a,': nlev=',i3,' blev=',f8.1,' dlev=',f8.1)")
     +      f(i)%fname8,f(i)%nlev,f(i)%lev(1),f(i)%dlev
          call handle_ncerr(istat,char120)
        endif
      enddo floop1
!
! Define global attributes:
!
! Model version string:
      istat = nf_put_att_text(ncid,NF_GLOBAL,"model_version",
     +  len_trim(h%version),trim(h%version)) 
      if (istat /= NF_NOERR) then
        write(char120,"('Error return from nf_put_att_text ',
     +    'for version global attribute')")
        call handle_ncerr(istat,char120)
      endif

! Satellite sampling
      if(index(satfile, 'tidi') .gt. 0) then
      write(char80,"(a,' ',i7)") 
     +  trim('tidi'),ifix(locs(4))
      istat = nf_put_att_text(ncid,NF_GLOBAL,"satellite_sampling",
     +  len_trim(char80),trim(char80))
      else if(index(satfile, 'saber') .gt. 0) then
       write(char80,"(a,' ',i7)") 
     +  trim('saber'),ifix(locs(4))
      istat = nf_put_att_text(ncid,NF_GLOBAL,"satellite_sampling",
     +  len_trim(char80),trim(char80))
      else if(index(satfile, 'champ') .gt. 0) then 
       write(char80,"(a,' ',i7)")
     +  trim('champ'),ifix(locs(4))
      istat = nf_put_att_text(ncid,NF_GLOBAL,"satellite_sampling",
     +  len_trim(char80),trim(char80))
      else if(index(satfile, 'guvi') .gt. 0) then 
       write(char80,"(a,' ',i7)")
     +  trim('guvi'),ifix(locs(4))
      istat = nf_put_att_text(ncid,NF_GLOBAL,"satellite_sampling",
     +  len_trim(char80),trim(char80))
      else if(index(satfile, 'cedar') .gt. 0) then
       write(char80,"(a,' ',i7)") 
     +  trim('cedar'),ifix(locs(4))
      istat = nf_put_att_text(ncid,NF_GLOBAL,"satellite_sampling",
     +  len_trim(char80),trim(char80))      
      else
      write(char80,"(a,' ',i7)") 
     +  trim('unknown'),ifix(locs(4))
      istat = nf_put_att_text(ncid,NF_GLOBAL,"satellite_sampling",
     +  len_trim(char80),trim(char80))
      endif
      if (istat /= NF_NOERR) then
        write(char120,"('Error return from nf_put_att_text ',
     +    'for satfile global attribute')")
        call handle_ncerr(istat,char120)
      endif

!      istat = nf_put_att_text(ncid,NF_GLOBAL,"satellite_sampling",
!     +        4,"test")
! Missing value:
      istat = nf_put_att_real(ncid,NF_GLOBAL,"missing_value",NF_FLOAT,
     +  1,spval)
      if (istat /= NF_NOERR) then
        write(char120,"('Error return from nf_put_att_real ',
     +    'for missing_value global attribute')")
        call handle_ncerr(istat,char120)
      endif
!
! Number of fields:
      istat = nf_put_att_int(ncid,NF_GLOBAL,"number_of_fields",
     +  NF_INT,1,nflds)
      if (istat /= NF_NOERR) then
        write(char120,"('Error return from nf_put_att_int ',
     +    'for nf global attribute')")
        call handle_ncerr(istat,char120)
      endif
!
! Field names:
      if (allocated(field_names)) deallocate(field_names)
      allocate(field_names(nflds),stat=iier)
      if (iier /= 0) call allocerr(ier,"mkcdf allocate field_names")
      ii = 0
      do i=1,nf
        if ((f(i)%requested.and.associated(f(i)%data)).or.
     +       f(i)%fname8 == 'Z       ') then
          ii = ii+1
          field_names(ii) = f(i)%fname8
        endif
      enddo
      istat = nf_put_att_text(ncid,NF_GLOBAL,"field_names",
     +  nflds*8,field_names)
      if (istat /= NF_NOERR) then
        write(char120,"('Error return from nf_put_att_text ',
     +    'for field_names')")
        call handle_ncerr(istat,char120)
      endif
      if (allocated(field_names)) deallocate(field_names)

!
! Number of times:
      istat = nf_put_att_int(ncid,NF_GLOBAL,"number_of_times",
     +  NF_INT,1,ntimes)
      if (istat /= NF_NOERR) then
        write(char120,"('Error return from nf_put_att_int ',
     +    'for ntimes global attribute: ntimes=',i4)") ntimes
        call handle_ncerr(istat,char120)
      endif
!
! First history (mssvol and mtime):
! (note last history is added when imode>1 below)
      write(char80,"(a,', ',i3.3,':',i2.2,':',i2.2)") 
     +  trim(h%mssvol),h%mtime
      istat = nf_put_att_text(ncid,NF_GLOBAL,"first_history",
     +  len_trim(char80),trim(char80))
      if (istat /= NF_NOERR) then
        write(char120,"('Error return from nf_put_att_text ',
     +    'for first_history global attribute: ',a)") trim(char80)
        call handle_ncerr(istat,char120)
      endif

! Exit define mode:
      istat = nf_enddef(ncid)
      if (istat /= NF_NOERR) call handle_ncerr(istat,
     +  'Error return from nf_enddef')
!
! Give values to dimension variables:
!#      print *, ' << MKSATCDF >> : gcmlon = ', gcmlon
      istat = nf_put_var_real(ncid,id_lon,gcmlon) 
      istat = nf_put_var_real(ncid,id_lat,gcmlat) 
      istat = nf_put_var_real(ncid,id_lev,gcmlev) 
      if (istat /= NF_NOERR) then
        write(char120,"('Error return from nf_put_var_real ',
     +    ' to assign values to dimension vars')")
        call handle_ncerr(istat,char120)
      endif
!
! End imode==0
!
      if (ier /= 0) return
      if (imode==0) return
 100  continue
!
! Give value to time var (total model HOURS) for current model time:
!      mins = h%mtime(2)*1.+h%mtime(3)/60.
      mins = locs(1)/60.
      istat = nf_put_var1_real(ncid,idv_time,itime,mins)
      if (istat /= NF_NOERR) then
        write(char120,"('Error from nf_put_var1_int defining ',
     +    'time at itime=',i3,': hours=',i10)") itime,mins
        call handle_ncerr(istat,char120)
      endif
! Give value to lst var for current model time:
      istat = nf_put_var1_real(ncid,idv_lst,itime,locs(8))
      if (istat /= NF_NOERR) then
        write(char120,"('Error from nf_put_var1_real defining ',
     +    'lst at itime=',i3,': lst=',i8)")
     +    itime,locs(8)
        call handle_ncerr(istat,char120)
      endif
!
!
! Give values to model time var:
      start_2d(1) = 1
      start_2d(2) = itime
      count_2d(1) = 3
      count_2d(2) = 1
      istat = nf_put_vara_int(ncid,idv_mtime,start_2d,
     +  count_2d,h%mtime)
      if (istat /= NF_NOERR) then
        write(char120,"('Error defining mtime at itime=',i2,
     +    ' h%mtime=',3i4,' istat=',i5)") itime,h%mtime,istat
        call handle_ncerr(istat,char120)
      endif
!
! Give value to histvols var for current model time:
! (note special treatment of character variable)
      start_2d(1) = 1
      start_2d(2) = itime
      count_2d(1) = len_trim(h%mssvol)+1 ! extra for null terminate
      count_2d(2) = 1
      char80 = h%mssvol
      ipos = len_trim(h%mssvol)+1
      char80(ipos:ipos) = char(0)	! null terminate
      istat = nf_put_vara_text(ncid,idv_histvols,start_2d,count_2d,
     +  char80,80)
      if (istat /= NF_NOERR) then
        write(char120,"('Error from nf_put_vara_text defining ',
     +    'histvol at itime=',i3,': h%mssvol=',a)") itime,trim(h%mssvol)
        call handle_ncerr(istat,char120)
      endif
!
! Give value to date var for current model time:
      istat = nf_put_var1_int(ncid,idv_date,itime,ifix(locs(4)))
      if (istat /= NF_NOERR) then
        write(char120,"('Error from nf_put_var1_int defining ',
     +    'date at itime=',i3,': locs(4)=',i8)") itime,ifix(locs(4))
        call handle_ncerr(istat,char120)
      endif
!
! Give value to adflag var for current model time:
      istat = nf_put_var1_int(ncid,idv_adflag,itime,ifix(locs(5)))
      if (istat /= NF_NOERR) then
        write(char120,"('Error from nf_put_var1_int defining ',
     +    'adflag at itime=',i3,': locs(5)=',i8)") itime,ifix(locs(5))
        call handle_ncerr(istat,char120)
      endif
!
! Give value to dnflag var for current model time:
      istat = nf_put_var1_int(ncid,idv_dnflag,itime,ifix(locs(7)))
      if (istat /= NF_NOERR) then
        write(char120,"('Error from nf_put_var1_int defining ',
     +    'dnflag at itime=',i3,': locs(7)=',i8)") itime,ifix(locs(7))
        call handle_ncerr(istat,char120)
      endif
!
! Give value to angle var for current model time:
      istat = nf_put_var1_int(ncid,idv_angle,itime,ifix(locs(6)))
      if (istat /= NF_NOERR) then
        write(char120,"('Error from nf_put_var1_int defining ',
     +    'wcflag at itime=',i3,': locs(6)=',i8)") itime,ifix(locs(6))
        call handle_ncerr(istat,char120)
      endif
!
! Give value to track var for current model time:
      if(index(satfile, 'tidi') .gt. 0) then
      istat = nf_put_var1_real(ncid,idv_track,itime,locs(9))
      if (istat /= NF_NOERR) then
        write(char120,"('Error from nf_put_var1_real defining ',
     +   'track at itime=',i3,': locs(9)=',f8.2)") itime,locs(9)
        call handle_ncerr(istat,char120)
      endif
      endif 
!
! Give value to yyddd var for current model time:
!      istat = nf_put_var1_int(ncid,idv_yyddd,itime,h%iyd)
!      if (istat /= NF_NOERR) then
!        write(char120,"('Error from nf_put_var1_int defining ',
!     +    'yyddd at itime=',i3,': h%iyd=',i8)") itime,h%iyd
!        call handle_ncerr(istat,char120)
!      endif
!
! Give value to longitude var for current model time:
      istat = nf_put_var1_real(ncid,idv_lon,itime,locs(2))
      if (istat /= NF_NOERR) then
        write(char120,"('Error from nf_put_var1_real defining ',
     +    'lon at itime=',i3,': locs(2)=',i8)")
     +    itime,locs(2)
        call handle_ncerr(istat,char120)
      endif
!
! Give value to lat var for current model time:
      istat = nf_put_var1_real(ncid,idv_lat,itime,locs(3))
      if (istat /= NF_NOERR) then
        write(char120,"('Error from nf_put_var1_int defining ',
     +    'lat at itime=',i3,': locs(3)=',i8)")
     +    itime,locs(3)
        call handle_ncerr(istat,char120)
      endif
!
! Give value to f107d var for current model time:
!      istat = nf_put_var1_real(ncid,idv_f107d,itime,hdr%f107d)
!      if (istat /= NF_NOERR) then
!        write(char120,"('Error from nf_put_var1_real defining ',
!     +    'f107d at itime=',i3,': hdr%f107d=',f8.2)") itime,hdr%f107d
!        call handle_ncerr(istat,char120)
!      endif
!
! Give value to f107a var for current model time:
!      istat = nf_put_var1_real(ncid,idv_f107a,itime,hdr%f107a)
!      if (istat /= NF_NOERR) then
!        write(char120,"('Error from nf_put_var1_real defining ',
!     +    'f107a at itime=',i3,': hdr%f107a=',f8.2)") itime,hdr%f107a
!        call handle_ncerr(istat,char120)
!      endif
!
! Give value to hp var for current model time:
!      istat = nf_put_var1_real(ncid,idv_hp,itime,hdr%hp)
!      if (istat /= NF_NOERR) then
!        write(char120,"('Error from nf_put_var1_real defining ',
!     +    'hp at itime=',i3,': hdr%hp=',f8.2)") itime,hdr%hp
!        call handle_ncerr(istat,char120)
!      endif
!
! Give value to cp var for current model time:
!      istat = nf_put_var1_real(ncid,idv_cp,itime,hdr%cp)
!      if (istat /= NF_NOERR) then
!        write(char120,"('Error from nf_put_var1_real defining ',
!     +    'cp at itime=',i3,': hdr%cp=',f8.2)") itime,hdr%cp
!        call handle_ncerr(istat,char120)
!      endif
!
! Give value to byimf var for current model time:
!      istat = nf_put_var1_real(ncid,idv_by,itime,hdr%byimf)
!      if (istat /= NF_NOERR) then
!        write(char120,"('Error from nf_put_var1_real defining ',
!     +    'byimf at itime=',i3,': hdr%byimf=',f8.2)") itime,hdr%byimf
!        call handle_ncerr(istat,char120)
!      endif
!
! Give values to field variables for current model time:
!
!/      start_4d(1:3) = 1
!/      start_4d(4) = itime
!/      count_4d(1) = nlon
!/      count_4d(2) = nlat
      start_4d(1) = 1
      start_4d(2) = itime
      count_4d(1) = npress
      count_4d(2) = 1
      floop2: do i=1,nf
        if ((.not.f(i)%requested.or..not.associated(f(i)%data)).and.
     +       f(i)%fname8 /= 'Z       ') cycle floop2
!
! Get var id from name (actually only necessary when itime > 1)
        fname = f(i)%fname8
        ipos = index(fname,'+')
        if (ipos > 0) fname(ipos:ipos) = 'P'
        ipos = index(fname,'/')
        if (ipos > 0) fname(ipos:ipos) = '_'
        istat = nf_inq_varid(ncid,fname,idv_flds(i))
        if (istat /= NF_NOERR) then
          write(char120,"('Error getting id of field var ',a)") 
     +      f(i)%fname8
          call handle_ncerr(istat,char120)
        endif
!
! Netcdf doesn't like the pointer f(i)%data to be passed directly
!  to nf_put_vara_real, so transfer to tmp array var first. 
! Note alloc/dealloc of var must be inside loop because f(i)%nlev 
!   may vary.
!
!#        if (allocated(var)) deallocate(var)
!#        allocate(var(nlon,nlat,f(i)%nlev),stat=iier)
!#        if (iier /= 0) call allocerr(ier,"mkcdf allocate tmp var")
!#        var(:,:,:) = data(:,:,:)
!#        istat = nf_put_vara_real(ncid,idv_flds(i),start_4d,count_4d,var)
        istat = nf_put_vara_real(ncid,idv_flds(i),
     +                           start_4d,count_4d,data(:,i))
        if (istat /= NF_NOERR) then
          write(char120,"('Error return from nf_put_vara_real',
     +      ' for field var ',a,' itime=',i2)") f(i)%fname8,itime
          call handle_ncerr(istat,char120)
        endif
      enddo floop2
!#      if (allocated(var)) deallocate(var)
      if (ier /= 0) return
!#	print *, ' << MKSATCDF >> ', h%mtime, locs
!#      write(6,"('Wrote to netcdf file ',a,' at mtime=',
!#     +  i3,':',i2,':',i2,' itime=',i3)") trim(flnm_cdf),h%mtime,itime
!
! If last time iteration, add last_history attribute and close the file:
      if (imode > 1) then
!
! Add last history (mssvol and mtime) as global attribute:
        istat = nf_redef(ncid)
        if (istat /= NF_NOERR) call handle_ncerr(istat,
     +      'Error reentering define mode at last history')
        write(char80,"(a,', ',i3.3,':',i2.2,':',i2.2)") 
     +    trim(h%mssvol),h%mtime
        istat = nf_put_att_text(ncid,NF_GLOBAL,"last_history",
     +    len_trim(char80),trim(char80))
        if (istat /= NF_NOERR) then
          write(char120,"('Error return from nf_put_att_text ',
     +      'for last_history global attribute: ',a)") trim(char80)
          call handle_ncerr(istat,char120)
        endif
        istat = nf_enddef(ncid)
        if (istat /= NF_NOERR) call handle_ncerr(istat,
     +    'Error return from nf_enddef')
!
! Close netcdf file:
!/        istat = nf_close(ncid)
!/        if (istat /= NF_NOERR) then
!/          call handle_ncerr(istat,
!/     +      'Error return from nf_close to close netcdf file')
!/        else
!/          write(6,"('Closed netcdf file ',a)") trim(flnm_cdf)
!/        endif
      endif
      return
      end
!-------------------------------------------------------------------
!     subroutine handle_ncerr(istat,msg,ier)
!     include 'netcdf.inc'
!
! Handle a netcdf lib error:
!
!     integer,intent(in) :: istat
!     integer,intent(inout) :: ier
!     character(len=*),intent(in) :: msg
!
!     ier = ier+1
!     write(6,"(/72('-'))")
!     write(6,"('>>> Error from netcdf library:')")
!     write(6,"(a)") trim(msg)
!     write(6,"('istat=',i5)") istat
!     write(6,"(a)") nf_strerror(istat)
!     write(6,"(72('-')/)")
!     return
!     end
