Viewing contents of file '../idllib/deutsch/apo/dismkspcal.pro'
pro dismkspcal,rtnspec,calfile,h
;+
; No formal header yet. See the documentation in
; http//www.astro.washington.edu/deutsch/apoinfo.html
;-
if (n_params(0) ne 3) then begin
print,"Call> dismkspcal,rtnspec,calfile"
print,"e.g.> calfile='/astro/iraf/noao/lib/onedstds/spec50cal/hilt600.dat'"
print,"e.g.> disread,img,h,'n1.0015r,/spflat"
print,"e.g.> disspec,img,h,420,/tot,/wcal,/fitskyline,rtn=rtnspec"
print,"e.g.> dismkspcal,rtnspec,calfile,h"
return
endif
exptime=sxpar(h,'EXPOSURE')
airmass=sxpar(h,'AIRMASS')
openr,1,calfile
wave2=fltarr(1000) & mag2=fltarr(1000) & i=0 & lin=''
readf,1,lin & print,lin
while not EOF(1) do begin
readf,1,w2,f2,d2
wave2(i)=w2 & mag2(i)=f2 & i=i+1
endwhile
close,1
wave2=wave2(0:i-1) & mag2=mag2(0:i-1)
flux21=10^((mag2+48.59d)/(-2.5))
nAng=wave2*0+50
nHz=(3d18/(wave2-25)-3d18/(wave2+25))
flux2=flux21*nHz/nAng
; ----------------------------------------------------------------------------
openr,1,'/astro/iraf/noao/lib/onedstds/kpnoextinct.dat'
extinc=fltarr(2,81)
readf,1,extinc & close,1
; ----------------------------------------------------------------------------
s=size(rtnspec)
if (s(1) gt 790) then goto,REDPART
; ----------------------------------------------------------------------------
wave1=rtnspec(*,0)
cnts1=rtnspec(*,1)/exptime
plot,wave1,cnts1/max(cnts1),yr=[-0.05,1.05],ysty=1
oplot,wave2,flux2/max(flux2)
print,'Press any key....' & key=get_kbrd(1)
minwl=fix(min(wave1)/50)*50.0>3300-25 & maxwl=fix(max(wave1)/50)*50.0+25
wave5=findgen(maxwl-minwl+1)+minwl
cnts5=interpol(cnts1,wave1,wave5)
extincfac=interpol(10^(0.4*airmass*extinc(1,*)),extinc(0,*),wave5)
cnts5=cnts5*extincfac
wave6=fltarr(n_elements(wave5)/50-1) & cal6=wave6
for i=0,n_elements(wave6)-1 do begin
rng=where((wave5 ge i*50+minwl) and (wave5 le i*50+minwl+49))
totc=total(cnts5(rng))/50
wave6(i)=i*50+minwl+25
tmp1=where(wave2 eq i*50+minwl+25)
if (tmp1(0) ne -1) then cal6(i)=totc/flux2(tmp1)
endfor
plot,wave6,cal6/max(cal6),psym=-4
oplot,wave2,flux2/max(flux2)
print,'Press any key....' & key=get_kbrd(1)
print,'Writing bcal.dat....'
openw,1,'bcal.dat'
for i=0,n_elements(wave6)-1 do begin
if (cal6(i) ne 0) then printf,1,wave6(i),cal6(i)
endfor
close,1
return
; ----------------------------------------------------------------------------
REDPART:
rev=819-indgen(820)
wave1=reverse(rtnspec(*,0))
cnts1=reverse(rtnspec(*,1))/exptime
plot,wave1,cnts1/max(cnts1),yr=[-0.05,1.05],ysty=1
oplot,wave2,flux2/max(flux2)
print,'Press any key....' & key=get_kbrd(1)
minwl=fix(min(wave1)/50)*50.0>5100-25 & maxwl=fix(max(wave1)/50)*50.0<10200+25
wave5=findgen(maxwl-minwl+1)+minwl
cnts5=interpol(cnts1,wave1,wave5)
extincfac=interpol(10^(0.4*airmass*extinc(1,*)),extinc(0,*),wave5)
cnts5=cnts5*extincfac
wave6=fltarr(n_elements(wave5)/50-1) & cal6=wave6
for i=0,n_elements(wave6)-1 do begin
rng=where((wave5 ge i*50+minwl) and (wave5 le i*50+minwl+49))
totc=total(cnts5(rng))/50
wave6(i)=i*50+minwl+25
tmp1=where(wave2 eq i*50+minwl+25)
if (tmp1(0) ne -1) then cal6(i)=totc/flux2(tmp1)
endfor
plot,wave6,cal6/max(cal6),psym=-4
oplot,wave2,flux2/max(flux2)
print,'Press any key....' & key=get_kbrd(1)
print,'Writing rcal.dat....'
openw,1,'rcal.dat'
for i=0,n_elements(wave6)-1 do begin
if (cal6(i) ne 0) then printf,1,wave6(i),cal6(i)
endfor
close,1
return
end