Viewing contents of file '../idllib/astron/contrib/varosi/code/allpro/calc_effics.pro'
;+
; NAME:
; Calc_Effics
; PURPOSE:
; Compute efficiency of a grating using List of flux measurements.
; CALLING:
; Calc_Effics, eff_meas, effics
; INPUTS:
; eff_meas = array of structures, efficiency measurements
; collected by pro get_Eff_Data.
; KEYWORDS:
; REF_MIRROR = array of structures with mirror reflectivity,
; from function Read_Reflect (default = unity).
; REF_GRATING = array of structures with grating reflectivity,
; from function Read_Reflect (default = unity).
; OUTPUTS:
; effics = array of structures containing absolute & groove efficiency
; of grating and all measurements/errors used in computation.
; EXTERNAL CALLS:
; function unique
; function def_imhd_struct
; pro copy_struct_inx
; function spline
; function stdev
; COMMON BLOCKS:
; common Calc_Effics, efft ;structure template.
; PROCEDURE:
; Match wavelengths, average, interpolate reflectances, compute effic.
; HISTORY:
; Written: Frank Varosi NASA/GSFC 1994.
;-
pro Calc_Effics, eff_meas, effics, REF_MIRROR=ref_mirr, $
REF_GRATING=ref_grat, INITIALIZE=init
common Calc_Effics, efft
common Calc_Effic1, tension
if N_struct( efft ) NE 1 then begin
error = { EFF_ERROR_v2, WaveLength:0.0, $
flux_mirror:0.0, $
flux_grating:0.0, $
ref_mirror:0.0, $
ref_grating:0.0, $
absolute:0.0, $
groove:0.0 }
imh = def_imhd_struct()
efft ={ stime:"", $
grating:"", $
Lamp:"", $
object:"", $
environment:"", $
measurement:"", $
comment:"", $
filename:"", $
directory:"", $
flat_field:"", $
dev_axes: imh.dev_axes, $
dev_pos: imh.dev_pos, $
pbytes:0, $
psize:0, $
res:0, $
x0:0, y0:0, $
nx:0, ny:0, $
ydhms:intarr(5), $
counts:0L, $
dtime:0.0, $
WaveLength:0.0, $
flux_mirror:0.0, $
N_mirror:0, $
flux_grating:0.0, $
N_grating:0, $
ref_mirror:0.0, $
ref_m_file:"", $
ref_grating:0.0, $
ref_g_file:"", $
absolute:0.0, $
groove:0.0, $
errors:error, $
alpha_ang:0.0, $
gr_dens:0.0, $
gorder:1, $
saved:0 }
endif
if keyword_set( init ) then return
if N_struct( eff_meas ) LE 0 then begin
message,"no measurements in List",/INFO
return
endif
wg = where( eff_meas.mirror EQ 0, Ngrat )
wm = where( eff_meas.mirror EQ 1, Nmirr )
if (Ngrat LE 0) OR (Nmirr LE 0) then begin
if (Nmirr LE 0) then message,"no mirror measurements",/INFO
if (Ngrat LE 0) then message,"no grating measurements",/INFO
return
endif
if (Ngrat EQ 1) then begin
wgu = wg
Nwave = 1
endif else $
wgu = wg( unique( eff_meas(wg).wavelength, Nwave,/SORT,/ORIG ) )
effics = replicate( efft, Nwave )
copy_struct_inx, eff_meas, effics, INDEX_FROM=wgu, /RECUR_TANDEM
if (Nwave EQ Ngrat) then begin
effics.N_grating = 1
effics.flux_grating = eff_meas(wgu).flux
effics.errors.flux_grating = eff_meas(wgu).error
endif else begin
for iw = 0, Nwave-1 do begin
wavelen = effics(iw).wavelength
w = where( eff_meas(wg).wavelength EQ wavelen, nw )
effics(iw).N_grating = nw
if (nw GT 1) then begin
w = wg(w)
effics(iw).errors.flux_grating = $
sqrt( stdev( eff_meas(w).flux, fmean )^2 + $
total( eff_meas(w).error^2 )/nw )
effics(iw).flux_grating = fmean
endif else if (nw EQ 1) then begin
w = wg(w(0))
effics(iw).flux_grating = eff_meas(w).flux
effics(iw).errors.flux_grating=eff_meas(w).error
endif
endfor
endelse
for iw = 0, Nwave-1 do begin
wavelen = effics(iw).wavelength
w = where( eff_meas(wm).wavelength EQ wavelen, nw )
effics(iw).N_mirror = nw
if (nw GT 1) then begin
w = wm(w)
effics(iw).errors.flux_mirror = $
sqrt( stdev( eff_meas(w).flux, fmean )^2 + $
total( eff_meas(w).error^2 )/nw )
effics(iw).flux_mirror = fmean
endif else if (nw EQ 1) then begin
w = wm(w(0))
effics(iw).flux_mirror = eff_meas(w).flux
effics(iw).errors.flux_mirror = eff_meas(w).error
endif
endfor
if (Nwave GT 1) then effics = effics( sort( effics.wavelength ) )
wavl = effics.wavelength
if N_elements( tension ) NE 1 then tension = 20
BELL = string(7b)
errors = effics.errors
if N_struct( ref_mirr ) GT 0 then begin
message,"using mirror reflectance file: " + ref_mirr.file,/INFO
effics.ref_m_file = ref_mirr.file
wavr = ref_mirr.wavelength
if min( wavr ) GT min( wavl ) then $
message,"min wavelength not in range" +BELL,/INFO
if max( wavr ) LT max( wavl ) then $
message,"max wavelength not in range" +BELL,/INFO
effics.ref_mirror = $
spline( wavr, ref_mirr.reflectance, wavl, tension )
errors.ref_mirror = $
spline( wavr, ref_mirr.error, wavl, tension )
endif else begin
message,"default mirror reflectance = 1" + BELL,/INFO
effics.ref_mirror = 1
errors.ref_mirror = 0
effics.ref_m_file = ""
endelse
if N_struct( ref_grat ) GT 0 then begin
message,"using grating reflectance file: " + ref_grat.file,/INFO
effics.ref_g_file = ref_grat.file
wavr = ref_grat.wavelength
if min( wavr ) GT min( wavl ) then $
message,"min wavelength not in range" +BELL,/INFO
if max( wavr ) LT max( wavl ) then $
message,"max wavelength not in range" +BELL,/INFO
effics.ref_grating = $
spline( wavr, ref_grat.reflectance, wavl, tension )
errors.ref_grating = $
spline( wavr, ref_grat.error, wavl, tension )
endif else begin
message,"default grating reflectance = 1" + BELL,/INFO
effics.ref_grating = 1
errors.ref_grating = 0
effics.ref_g_file = ""
endelse
effics.absolute = $
( effics.flux_grating/effics.flux_mirror ) * effics.ref_mirror
errd2 = ( errors.flux_grating / effics.flux_grating )^2 + $
( errors.flux_mirror / effics.flux_mirror )^2 + $
( errors.ref_mirror / effics.ref_mirror )^2
errors.absolute = effics.absolute * sqrt( errd2 )
effics.groove = effics.absolute/effics.ref_grating
errors.groove = effics.groove * sqrt( errd2 + $
( errors.ref_grating / effics.ref_grating )^2 )
effics.errors = errors
end