pro test_v2p

vs='Icarusr6r6'
dir='/misc/mcc20/common/ongocart/output/'+vs+'/2010/monthly/tavg3d_chm_v/'
fn=vs+'.tavg3d_chm_v.monthly.201008.nc4'
readncdf,dir,fn,'lon',lon
readncdf,dir,fn,'lat',lat
readncdf,dir,fn,'lev',lev
nlon=n_elements(lon)
nlat=n_elements(lat)
nlev=n_elements(lev)

READNCDF,dir,fn,'ps',ps
READNCDF,dir,fn,'delp',delp
READNCDF,dir,fn,'AIRDENS',airdens
READNCDF,dir,fn,'CO',co
;delz=delp/(airdens*9.81) ; delz in m
ps=ps*0.01
delp=delp*0.01

psel=150.
; Call V2P to get the level index ll (nlon,nlat) for selected psel presure
V2P,nlon,nlat,nlev,ps,delp,psel,ll
print,'min,max ll ',min(ll(where(ll ge 0))),max(ll(where(ll ge 0)))

cop=fltarr(nlon,nlat) & cop(*,*)=-999.9 ; co at the selected pressure level
for i=0,nlon-1 do begin
for j=0,nlat-1 do begin
 lp=ll(i,j)
 if lp ge 0 then cop(i,j)=co(i,j,lp)*1.e9
endfor
endfor
print,'min,max CO at ',psel,'hPa',min(cop),max(cop)

;fign='test_v2p'
;lvl=findgen(11)*max(cop)/10
;pmap,cop,lon,lat,fign='test_v2p',levels=lvl
;print,fign+'.ps'

stop
end