;
pro deflatdata,lat,data
;
; Define lat.data(nlon,nlev) from processed var in data.
;
field = lat.field
nlon = field.nlon & nlat = field.nlat & nlev = field.nlev
;print,'deflatdata: Field ',field.name,' nlon=',nlon,' nlev=',nlev
;
; Check dimensions (return if field is only 2d lat x lon)
;
dims = size(data,/dimensions)
ndims = n_elements(dims)
error = 0
if (dims[0] ne field.nlon) then begin
  print,'>>> WARNING deflatdata: 1st dimension is not nlons: dims[0]=',$
    dims[0],' but nlon=',field.nlon
  error = 1 
endif
if (dims[1] ne field.nlat) then begin
  print,'>>> WARNING deflatdata: 2nd dimension is not nlats: dims[1]=',$
    dims[1],' but nlat=',field.nlat
  error = 1 
endif
if (dims[2] ne field.nlev) then begin
  print,'>>> WARNING deflatdata: 3rd dimension is not nlats: dims[2]=',$
    dims[2],' but nlev=',field.nlev
  error = 1 
endif
if nlev le 0 then begin
  print,'>>> deflatdata: nlev=',nlev,' -- cannot make vertical slice.'
  error = 1 
endif
if error gt 0 then begin
  size = size(lat.data)
  latdata = fltarr(size[0],size[2])       ; latdata(nlon,nlev)
  latdata[*,*] = float(lat.missing_value)
  lat.data = ptr_new(latdata)
  return
endif
lons = *field.lons
lats = *field.lats
fieldlevs = *field.levs
nfieldlevs = n_elements(fieldlevs)
;print,'deflatdata: nfieldlevs=',nfieldlevs,' fieldlevs=' & print,fieldlevs
;
; If axis limits levmin,max are set, reset levs accordingly:
;
if lat.levmin lt lat.levmax then begin ; set axes limits levs:
; print,'deflatdata: lat.levmin,max=',lat.levmin,lat.levmax
  nlev = 0
  for k=0,nfieldlevs-1 do $
    if fieldlevs[k] ge lat.levmin and fieldlevs[k] le lat.levmax then nlev = nlev+1
  if nlev lt 2 then print,'>>> WARNING deflatdata: nlev point for y-axis must be',$
    ' at least 2: nlev=',nlev
  levs = fltarr(nlev)
  kk = 0
  for k=0,nfieldlevs-1 do $
    if fieldlevs[k] ge lat.levmin and fieldlevs[k] le lat.levmax then begin
      levs[kk] = fieldlevs[k]
      kk = kk+1
    endif
endif else levs = fieldlevs
nlev = n_elements(levs)
;
;print,'deflatdata: nfieldlevs=',nfieldlevs,' nlev=',nlev,' levs=' & print,levs
;
; Define (nlon,nlev) at lat.slon:
;
latdata = fltarr(nlon,nlev)
latdata[*,*] = 0.
ilat = ixfind_nearest(lats,lat.slat)
if ilat lt 0 then begin
  print,'>>> deflatdata: cannot find latitude ',lat.slat,': nlat=',nlat,$
    ' ilat=',ilat,' lats=' & print,lats
  lat.data = ptr_new(latdata)
  return
endif
;
if lat.levmin lt lat.levmax then begin 
  kk = 0
  for k=0,nfieldlevs-1 do $
    if fieldlevs[k] ge lat.levmin and fieldlevs[k] le lat.levmax then begin
      latdata[*,kk] = float(data[*,ilat,k,lat.imtime])
      kk = kk+1
    endif
endif else begin
  latdata[*,*] = float(data[*,ilat,*,lon.imtime])
endelse
;
; Take log10 of densities:
;
if (lat.log10 eq 'density fields' and isdensity(field.name) gt 0) or $
   (lat.log10 eq 'this frame only') or $
   (lat.log10 eq 'all fields') then begin
; print,'deflatdata: taking log10 of field ',field.name,'...'
  missing_value = float(lat.missing_value)
  log10f,latdata,missing_value
endif
;
; Update lon structure:
;
lat.slat = lats[ilat]
lat.data = ptr_new(latdata)
lat.lons = ptr_new(lons)
lat.levs = ptr_new(levs)
end