; --- Make the new ancil -----------------------------------------------------

dir='/data/beowulf2/gcm/um/vn*/ancil/atmos/cl*/hadcm3/'

;
; First: convert a normal ancil file into pp format
;
spawn,'wconv.pl '+dir+'waterfix.ancil.4 waterfix.ancil.32.pp'

; 
; Now read in the ancil-as-pp, and modify it
;
a=get_field('waterfix.ancil.32.pp',/nat)
y=yc(a,/arr)
x=xc(a,/arr)

i=where(y lt -53 and a.data ne a.bmdi,c)
print,'Selecting ',c,' points'
a1=a
; a1.data(i)=3.888e-6       
a1.data(i)=3.888e-5       

;
; Write out the pp field (in pp format)
;
pp_write,a1,'glen-waterfix.ancil.new5.32.pp',/nat

;
; Now insert the new pp-field into an ancil file
;
spawn,'wunconv2.pl '+dir+'waterfix.ancil.4 glen-waterfix.ancil.new5.32.pp glen-waterfix.new5.ancil.4'

;
; --- Now just plot etc -----------------------------------------------------
;

a=get_field('waterfix.ancil.32.pp',/nat)
b=get_field('glen-waterfix.ancil.new5.32.pp',/nat)

; lev=(1+indgen(10))*4e-7
lev=(1+indgen(10))*4e-6
pp_plot,a,/blo,la=-45,ymarg=[2,0],xmarg=[0,0],lev=lev    
pp_plot,b,/blo,la=-45,ymarg=[2,0],xmarg=[0,0],lev=lev
pp_zonal_mean,a
pp_zonal_mean,b

print,pp_area_total(a,la=-10)
print,pp_area_total(b,la=-10)

gettwogifs,out='glen5'