pro open_files,stateptr
  state = *stateptr
  files = *state.files
  nfiles = n_elements(files) ; number of currently selected files
  mxtimes = 10000            ; max number of total times (all files)
  mtimes = intarr(3,mxtimes) ; mtimes to plot
  mxtime_file = 1000         ; max number of times per file
  mtime_files = intarr(3,mxtime_file,nfiles)  ; mtimes on each file
  ntime_files = intarr(nfiles)                ; number of times on each file
  ntimes = 0
  mxfields = 100
; fields = strarr(mxfields)
  fields = replicate({field_struct},mxfields)
;
; Open files and get number of times (histories):
;
  for ifile=0,nfiles-1 do begin
    ncid = ncdf_open(files[ifile],/nowrite)
    id = ncdf_dimid(ncid,'time')
    if (id eq -1) then begin
      print,'>>> Open_files: Could not find id of time dimension on file ',files[ifile]
      return
    endif
    name = ''
    ncdf_diminq,ncid,id,name,ntime
    if (ntimes+ntime gt mxtimes) then begin
      print,'>>> Open_files: ntimes gt mxtimes: ntimes=',ntimes,' mxtimes=',mxtimes
      return
    endif 
    if (ntime gt mxtime_file) then begin
      print,'>>> Open_files: ntime gt mxtime_file: ntime=',ntime,' mxtime_file=',$
        mxtime_file
      return
    endif
;
; Get model times and add to total mtimes:
;
    id = ncdf_varid(ncid,'mtime')
    if (id eq -1) then begin
      print,'>>> Open_files: Could not find id of mtime var on file ',files[ifile]
      return
    endif
    ncdf_varget,ncid,id,mtime
    mtimes[*,ntimes:ntimes+ntime-1] = mtime[*,*]
    ntimes = ntimes+ntime
;
; Get list of 4-d fields from first file:
;
    if (ifile eq 0) then begin
      fileinfo = ncdf_inquire(ncid)
      ndims = fileinfo.ndims
      nvars = fileinfo.nvars
      nfields = 0
      for ivar = 0,nvars-1 do begin
        vinfo = ncdf_varinq(ncid,ivar)
;
; For now, assume that all 4-d fields are [lat,lon,lev,time] 
; (where lev could be lev or ilev)
; To Do: should exclude fields on mag grid
;
        if (vinfo.ndims eq 4) then begin
          fields[nfields].name = vinfo.name
          fields[nfields].model = 'TGCM'
;
; Set density species if same as any msis density sp:
          if (vinfo.name eq 'O2' or vinfo.name eq 'O1' or vinfo.name eq 'H' or $
              vinfo.name eq 'HE' or vinfo.name eq 'N2' or vinfo.name eq 'N' or $ 
              vinfo.name eq 'AR') then fields[nfields].density = 1
          nfields = nfields+1
        endif
;
; Get grid locations:
;
        case vinfo.name of
          'lon': begin
            ncdf_varget,ncid,ivar,lons
;           print,format="('First file: lons = ',/,(8f8.2))",lons
            *state.lons = lons
            state.horiz_resolution = lons[1]-lons[0] ; assume regular grid, and that dlat==dlon
          end
          'lat': begin
            ncdf_varget,ncid,ivar,lats
;           print,format="('First file: lats = ',/,(8f8.2))",lats
            *state.lats = lats
          end
          'lev': begin
            ncdf_varget,ncid,ivar,levs
;           print,format="('First file: levs = ',/,(8f8.2))",levs
            *state.levs = levs
            state.vert_resolution = levs[1]-levs[0] ; assume regular grid
          end
          'ilev': begin
            ncdf_varget,ncid,ivar,ilevs
;           print,format="('First file: ilevs = ',/,(8f8.2))",ilevs
            *state.ilevs = ilevs
          end
          else: 
        endcase ; vinfo.name
      endfor ; ivar=0,nvars-1
;     print,format="('First file: nfields=',i4,' fields = ',/,(8a8))",$
;       nfields,fields[0:nfields-1].name
    endif
;
    if (ntime gt 1) then begin
      for it=1,ntime-1 do begin
        deltamin = mtime_to_mins(mtime[*,it])-mtime_to_mins(mtime[*,it-1])
        if deltamin ne 0 then break
      endfor
    endif else deltamin = 0.
    mtime_files[*,0:ntime-1,ifile] = mtime[*,*]
    ntime_files[ifile] = ntime
    print,format="('File ',a,' ntime=',i4,' (',3(i3,','),' to ',3(i3,','),' by ',i4,')')",$
      file_basename(files[ifile]),ntime,mtime[*,0],mtime[*,ntime-1],deltamin
    ncdf_close,ncid
  endfor ; ifile=0,nfiles-1
  print,format=$
    "('Total ntimes = ',i5,' (',i3,',',i3,',',i3,' to ',i3,',',i3,',',i3,')')",$
    ntimes,mtimes[*,0],mtimes[*,ntimes-1]
  if (ntimes le 1) then begin
    print,'>>> WARNING open_files: ntimes=',ntimes,$
      ': must have more than 1 time to make plots.'
  endif
;
; "old" histories do not have vert coord at interface (ilev):
;
  if (n_elements(*state.ilevs) eq 0) then begin
    print,'Note open_files: Did not read ilevs. Will assign ilevs <- levs'
    *state.ilevs = *state.levs
  endif
; print,'open_files: n_elements(*state.ilevs)=',n_elements(*state.ilevs)
;
; Update state:
;
  state.ntimes = ntimes
  *state.mtimes = mtimes[*,0:ntimes-1]
  state.mtime_beg = mtimes[*,0]
  state.mtime_end = mtimes[*,ntimes-1]
  *state.fields_tgcm = fields[0:nfields-1]
  *state.mtime_files = mtime_files
  *state.ntime_files = ntime_files
  *stateptr = state
end