pro vap_gts,dtr_prefix=dtr_prefix,tmp_prefix=tmp_prefix,year1,year2,$
    nor1=nor1,nor2=nor2,outprefix=outprefix
; calculates vapour pressure anomaly time series from dtr and tmp, using tmn 
;  as proxy for tdew
; this version also calculates the tdew normal from the vapour pressure normal
;  and then adjusts the tmin to have an average equal to tdew normal
; synthetic values are contrained to be between 0.1 and saturated vap
;  pressure at mean temperature
; Note that anomalies are in percent, whole numbers.
if n_params() lt 1 then begin
 print,' vap_gts,dtr_prefix=dtr_prefix,tmp_prefix=tmp_prefix,year1,year2,'
 print,'         nor1=nor1,nor2=nor2,outprefix=outprefix'
 return
endif
close,/all
if keyword_set(nor1) eq 0 then nor1=1961
if keyword_set(nor2) eq 0 then nor2=1990
if keyword_set(outprefix) eq 0 then outprefix='glo.vap.ano.'
outprefix2='glo.vap.'
if keyword_set(dtr_prefix) eq 0 then dtr_prefix='../../dtr/glo/glo.dtr.'
if keyword_set(tmp_prefix) eq 0 then tmp_prefix='../../tmp/glo/glo.tmp.'
rdbin,vapnor,'glo.vap.norm'
vapnor=vapnor/10.0
;
;
; calculate 1961-1990 tmn normal from monthly tmp and dtr values
print,'Calculating tmn normal'
for iy=nor1,nor2 do begin
 print,iy
 tmpfl=strip(string(tmp_prefix,iy))
 dtrfl=strip(string(dtr_prefix,iy))
 rdbin,tmpgrd,tmpfl
 rdbin,dtrgrd,dtrfl
 if iy eq 1961 then begin
  nsea=where(tmpgrd le -9999)
  nland=where(tmpgrd gt -9999)
  tmnnor=float(tmpgrd)*0.0
 endif
 tmnnor(nland)=tmnnor(nland)+(tmpgrd(nland)-(0.5*dtrgrd(nland)))/10.0
endfor
tmnnor(nland)=tmnnor(nland)/(nor2-nor1+1)
vapnor=vapnor(nland)>0.1
tdwnor=tvap(vapnor)
tadj=tdwnor-tmnnor(nland)

; calculate 1961-1990 synthetic normal from adjusted tmn
print,'Calculating synthetic vap normal'
for iy=nor1,nor2 do begin
 print,iy
 tmpfl=strip(string(tmp_prefix,iy))
 dtrfl=strip(string(dtr_prefix,iy))
 rdbin,tmpgrd,tmpfl
 rdbin,dtrgrd,dtrfl
 if iy eq 1961 then begin
  vapsyn=float(tmpgrd)*0.0
 endif
 tmn=(tmpgrd(nland)-(0.5*dtrgrd(nland)))/10.0 + tadj
 dtrgrd=1
 v=6.108*exp( (17.27*tmn) /(237.3+tmn) )>0.05
 v=v<esat(tmpgrd(nland)/10.0)
 vapsyn(nland)=vapsyn(nland)+v
 tmpgrd=1 
endfor
vapsyn(nland)=vapsyn(nland)/(nor2-nor1+1)>0.01
;
;
;  Calculate synthetic vap from tmin, convert to anomalies 
;  relative to synthetic mean vap, and apply to normal vap
third:
print,'Calculating synthetic anomalies'
for iy=year1,year2 do begin
 print,iy
 tmpfl=strip(string(tmp_prefix,iy))
 dtrfl=strip(string(dtr_prefix,iy))
 rdbin,tmpgrd,tmpfl
 rdbin,dtrgrd,dtrfl
 if iy eq year1 then begin
  vapgrd=float(tmpgrd)*0.0
  nsea=where(tmpgrd le -9999)
  nland=where(tmpgrd gt -9999)
 endif
 tmn=tmpgrd(nland)/10.0-0.5*dtrgrd(nland)/10.0 + tadj
 dtrgrd=1
 v=6.108*exp( (17.27*tmn) /(237.3+tmn) )>0.05
 v=v<esat(tmpgrd(nland)/10.0)
 tmpgrd=1
 vapgrd(nland)=v
 vapgrd(nland)=100*(vapgrd(nland)-vapsyn(nland))/vapsyn(nland)
 vapgrd(nland)=vapgrd(nland)>(-9998)
 vapgrd(nland)=vapgrd(nland)<9998
 vapgrd(nsea)=-9999
 wrbin,fix(round(vapgrd)),strip(string(outprefix,iy)),/compre  ; write anomalies
 ;
 vapgrd(nland)=(vapgrd(nland)/100.0)*(vapnor*10.0)+(vapnor*10.0)
 wrbin,fix(round(vapgrd)),strip(string(outprefix2,iy)),/compre  ; write absolutes
endfor



end

