pro frs_gts,dtr_prefix=dtr_prefix,tmp_prefix=tmp_prefix,year1,year2,$
    nor1=nor1,nor2=nor2,outprefix=outprefix,frssyn=frssyn
; calculates ground-frost anomaly time series from tmin, using dtr and tmp
;  to calculate tmin
; then uses anomalies with frost normal grid to calculate monthly frost freq
if n_params() lt 1 then begin
 print,' frs_gts,dtr_prefix=dtr_prefix,tmp_prefix=tmp_prefix,year1,year2,'
 print,'         nor1=nor1,nor2=nor2,outprefix=outprefix,frssyn=frssyn'
 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.frs.ano.'
outprefix2='glo.frs.'
if keyword_set(dtr_prefix) eq 0 then dtr_prefix='~/m1/gts/dtr/glo/glo.dtr.'
if keyword_set(tmp_prefix) eq 0 then tmp_prefix='~/m1/gts/tmp/glo/glo.tmp.'
rdbin,frsnor,'~/m1/gts/frs/glo/glo.frs.norm'
nland=where(frsnor gt -9999)
nsea=where(frsnor eq -9999)
frsnor=frsnor/10.0
frsnor=frsnor(nland)>0.0
;
; calculate 1961-1990 synthetic normal from adjusted tmn
print,'Calculating synthetic frs normal'
if n_elements(frssyn) gt 1 then goto,skipsyn
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
  frssyn=float(tmpgrd)*0.0
 endif
 tmn=(tmpgrd(nland)-(0.5*dtrgrd(nland)))/10.0
 dtrgrd=1
 tmpgrd=1 
 frssyn(nland)=frssyn(nland)+frscal(tmn)
endfor
frssyn(nland)=frssyn(nland)/(nor2-nor1+1)
for im=0,11 do begin
 temp=frssyn(*,*,im)
 nfin=where(temp gt 0)
 temp(nfin)=(temp(nfin)/100.0)*days(im)
 frssyn(*,*,im)=temp
endfor
tmn=0
frssyn(nsea)=-999.9
skipsyn:
;
;
;  Calculate synthetic frs from tmin, convert to anomalies 
;  relative to synthetic mean frs, and apply to normal frs
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
  frsgrd=float(tmpgrd)*0.0
 endif
 tmn=tmpgrd(nland)/10.0-0.5*dtrgrd(nland)/10.0 
 dtrgrd=1
 tmpgrd=1
 frsgrd(nland)=frscal(tmn)
 for im=0,11 do begin
  temp=frsgrd(*,*,im)
  nfin=where(temp gt 0)
  temp(nfin)=(temp(nfin)/100.0)*days(im)
  frsgrd(*,*,im)=temp
 endfor
 frsgrd(nland)=10*(frsgrd(nland)-frssyn(nland))
 frsgrd(nland)=frsgrd(nland)>(-9998)
 frsgrd(nland)=frsgrd(nland)<9998
 frsgrd(nsea)=-9999
 frsgrd=fix(round(frsgrd))
 tmn=0
 temp=0
 wrbin,frsgrd,strip(string(outprefix,iy)),/compre  ; write anomalies
 ;
 print,rnge(frsnor)
 frsgrd(nland)=frsgrd(nland)+(frsnor*10.0)
 frsgrd(nland)=frsgrd(nland)>0
 for im=0,11 do frsgrd(*,*,im)=frsgrd(*,*,im)<(days(im)*10)
 frsgrd=fix(round(frsgrd))
 wrbin,frsgrd,strip(string(outprefix2,iy)),/compre  ; write absolutes
endfor
end