;
function yaxpos,x,vec0,vec1,nvec,vb,vt
;             zp   ht
  bracket, x,vec1,nvec,1,k0,k1,ier
  f1 = (x-vec1(k0)) / (vec1(k1)-vec1(k0))
  xvec0 = f1*vec0(k1) + (1.-f1)*vec0(k0)
  ypos = vb + (vt-vb)*(xvec0-vec0(0))/(vec0(nvec-1)-vec0(0))
; print,'yaxpos returning: x=',x,' ypos=',ypos
  return,ypos
end
;-----------------------------------------------------------------------
;
pro bracket, x,xx,nx,inc,n1,n2,ier
;
; Bracket x in xx(nx), returning lower index in n1 and upper index in n2
; If inc > 0 -> array increases from bottom to top,
; If inc <= 0 -> array increases from top to bottom
;
ier = 0
;
; Array increases from bottom to top:
;
if (inc gt 0) then begin
  if (x lt xx(0)) then begin
    n1 = 0
    n2 = 1
    ier = 1
    return
  endif
  if (x gt xx(nx-1)) then begin
    n1 = nx-2
    n2 = nx-1
    ier = 2
    return
  endif
  for i=0,nx-2 do begin
    if (x ge xx(i) and x le xx(i+1)) then begin
      n1 = i
      n2 = i+1
      return
    endif
  endfor
;   
; Array increases from top to bottom:
;
endif else begin
  if (x gt xx(0)) then begin
    n1 = 0
    n2 = 1
    ier = 3
    return
  endif
  if (x lt xx(nx-1)) then begin
    n1 = nx-2
    n2 = nx-1
    ier = 4
    return
  endif
  for i=1,nx-2 do begin
    if (x le xx(i) and x ge xx(i+1)) then begin
      n1 = i
      n2 = i+1
      return
    endif
  endfor
endelse
return
end
;-----------------------------------------------------------------------
function rnd, val,rval,dir=dir
;
; Round val to the nearest rval (there must be a better way...)
; if dir == 0 (or dir not defined) then find nearest rval
; if dir > 0 then round up, if dir < 0 round down
;
ival = round(val) & irval = round(rval)
ivalp = ival
while ivalp mod irval ne 0 do ivalp = ivalp+1
ivalm = ival
while ivalm mod irval ne 0 do ivalm = ivalm-1
ans = ivalp
if n_elements(dir) > 0 then begin
  if dir lt 0 then ans = ivalm else return,ans
endif else begin
  if ival-ivalm lt ivalp-ival then ans = ivalm
endelse
return,ans
end
;-----------------------------------------------------------------------
pro altyax,info,slice,slicetype,vl,vr,vb,vt,charsize
zp = *slice.levs 
nlev = n_elements(zp)
;print,'altyax: vl,r,b,t=',vl,vr,vb,vt
;print,'altyax: zp min,max=',min(zp),max(zp)
fields = *info.fields
;
; Get heights from Z field:
;
ixz = -1
for i=0,n_elements(fields)-1 do begin
  if strcompress(fields[i].name,/remove_all) eq 'Z' then ixz = i
endfor
if ixz eq -1 then begin
  print,'>>> WARNING altyax: Could not find Z field -- extra right-hand',$
    ' axis NOT added.'
  return
endif
;
; Read Z field if necessary:
;
if not ptr_valid(fields[ixz].data) then begin ; field has not been read
; print,'WARNING altyax: Z has not been read..'
  print,'Reading Z field to get heights for right-hand axes..'
  widget_control,/hourglass
  varget,info,fields[ixz],zdata
  field = fields[ixz]
  procfield,info,zdata,field,info.z_hist
  fields[ixz] = field
  fmin = min(*fields[ixz].data) & fmax = max(*fields[ixz].data)
  rpt_minmax,info,fields[ixz],fmin,fmax
  *info.fields = fields
endif
fielddata = *fields[ixz].data
slicedata = *slice.data       ; save
log10save = slice.log10       ; save
slice.log10 = 'none'          ; want linear heights
if slicetype eq 'lon' then begin
; print,'altyax calling deflondata for Z..'
  deflondata,slice,fielddata   ; get for lon slices (lat,lev)
endif else $
  deflatdata,slice,fielddata   ; get for lat slices (lon,lev)
zheights = *slice.data
*slice.data = slicedata ; restore
slice.log10 = log10save ; restore
;
; Average over first dimension (lat or lon):
;
dims = size(zheights,/dimensions)
heights = fltarr(dims[1])
for k=0,dims[1]-1 do begin
  heights[k] = calcmean(zheights[*,k],dims[0],0,1.e-20,$
    float(slice.missing_value),0)
endfor
;help,heights
;print,'Averaged heights=' & print,heights
;print,'Range of averaged heights=',heights[nlev-1]-heights[0]
 dht_axis = heights[nlev-1]-heights[0]
;
offtoax  = 0.    ; x offset from existing from existing to new axes
offtonum = .02   ; x offset from existing axis to numberic label on new axis
offtolab = .07   ; x offset from existing axis to new axis string label
offylab  = .009  ; y offset to center numerics next to tic marks
ticlen   = .012  ; major tic len (minor tic len = .5*ticlen)
if charsize eq 0. then chsize = 1.2 else chsize = charsize
;
; Draw new vertical axis offset to right of right axis:
;
xpos = vr+offtoax
plots,[xpos,xpos],[vb,vt],/norm
xyouts,vr+offtolab,.5*(vt+vb),'AVE HEIGHT (KM)',orientation=-90.,/norm,$
  align=.5,charsize=chsize
;
; Draw major tic marks and labels:
;
rndval = 10.
if dht_axis lt 25. then rndval = 5.
if dht_axis lt 20. then rndval = 3.
if dht_axis lt 10. then rndval = 2.
;print,'dht_axis=',dht_axis,' rndval=',rndval
botrnd = rnd(heights[0],rndval,dir=1)       ; round up
toprnd = rnd(heights[nlev-1],rndval,dir=-1) ; round down
;print,'botrnd=',botrnd,' toprnd=',toprnd
nhtlabmin = 6                   ; minimum number of major tics
del = (toprnd-botrnd)/nhtlabmin ; delta
delrnd = rnd(del,rndval)
if delrnd le 0 then delrnd = 1
ht = botrnd
while ht le toprnd do begin
  ypos = yaxpos(ht,zp,heights,nlev,vb,vt)
; print,'ht=',ht,' delrnd=',delrnd,' toprnd=',toprnd,' ypos=',ypos,$
;   ' vb=',vb,' vt=',vt
  if ypos le vt and ypos ge vb then begin
    plots,[xpos,xpos+ticlen],[ypos,ypos],/norm
    altlab = strcompress(string(format="(i4)",fix(ht+.5)),/remove_all)
    xyouts,vr+offtonum,ypos-offylab,altlab,/norm,charsize=chsize
;   print,'after label altlab=',altlab,' dht_axis=',dht_axis
  endif
;
; Draw minor tic marks:
;
  if dht_axis ge 10 then begin
    htnext = ht+delrnd                              ; height of next major
    yposnext = yaxpos(htnext,zp,heights,nlev,vb,vt) ; axis position of next major
    nminor = 3 & if yposnext-ypos gt 0.5*(vt-vb) then nminor = 6
    while (htnext-ht) mod nminor ne 0 do nminor=nminor+1
    delminor = (htnext-ht)/nminor
    if delminor le 0 then delminor = 1
;   print,'delminor=',delminor
    htminor = ht-delrnd
;
; Draw minor tics below bottom major tic:
    while htminor le ht do begin
      yposminor = yaxpos(htminor,zp,heights,nlev,vb,vt)
      if yposminor gt vb then $
        plots,[xpos,xpos+0.5*ticlen],[yposminor,yposminor],/norm
      htminor = htminor+delminor
    end 
    htminor = ht+delminor
;
; Draw minor tics between current and next majors:
    while htminor le htnext do begin
      yposminor = yaxpos(htminor,zp,heights,nlev,vb,vt)
      if yposminor le vt and yposminor ge vb then begin
        plots,[xpos,xpos+0.5*ticlen],[yposminor,yposminor],/norm
      endif 
      if yposnext-ypos gt 0.5*(vt-vb) and htminor eq ht+3.*delminor then begin
        altlab='('+strcompress(string(format="(i4)",fix(ht+.5)),/remove_all)+')'
;       xyouts,vr+offtonum,yposminor-offylab,altlab,/norm,charsize=chsize
      endif
      htminor = htminor+delminor
    end
  endif ; dht_axis ge 10
;
; Increment major tic:
  ht = ht+delrnd
end ; ht le toprnd
end