Документ взят из кэша поисковой машины. Адрес оригинального документа : http://hea-www.harvard.edu/PINTofALE/pro/flux_to_em.pro Дата изменения: Wed Jan 24 19:42:16 2001 Дата индексирования: Tue Oct 2 01:05:41 2012 Кодировка: Поисковые слова: apollo 11

function flux_to_em,emis,flux=flux,logT=logT,wvl=wvl,Z=Z,NH=NH,defEM=defEM,$
noph=noph,thresh=thresh, _extra=e
;+
;function flux_to_em
; return emission measures (-1 where undeterminable) at various
; temperatures consistent with given line fluxes and line emissivities.
;
; for each temperature, assume a delta-function emission measure,
; compute line fluxes seen through some instrument, account for
; interstellar absorption, and scale to match the observed flux.
;
;syntax
; em=flux_to_em(emis,flux=flux,logT=logT,wvl=wvl,Z=Z,NH=NH,$
; defEM=defEM,noph=noph,thresh=thresh, abund=abund,/temp,/ikev,$
; effar=effar,wvlar=wvlar,fh2=fh2,he1=he1,heII=heII,/fano)
;
;parameters
; emis [INPUT; required] 2D array of line emissivities, EMIS(LOGT,WVL)
; * if 1-D, assumed to be EMIS(LOGT)
;
;keywords
; flux [INPUT] observed fluxes
; * if not given, assumed to be 1 [(ph|erg)/s/cm^2]
; * size is forced to match 2nd dimension of EMIS
; if <, filled with 1s
; if >, extra elements are ignored
; * NOTE: the author realizes that there may be multiple IDs
; for the single observed line, and directs the inquirer
; to the function ID_TO_EMIS, which returns the appropriately
; concatenated emissivity table.
; logT [INPUT] log_10(Temperature [K]) at which EMIS is defined
; * if not given, assumed to go from 4 to 8 in even steps
; * interpolated if size does not match 1st dimension of EMIS
; wvl [INPUT] wavelengths at which EMIS is defined
; * if not given or if size does not match 2nd dimension of EMIS,
; (a) interstellar absorption correction is not made
; (b) all fluxes >gotta be< [ergs/s/cm^2]
; Z [INPUT] atomic numbers
; * if not given, assumed to be 1s (i.e., H. i.e., no
; abundance corrections are made)
; * size is forced to match 2nd dimension of EMIS
; if <, filled with 1s
; if >, extra elements are ignored
; NH [INPUT] H column density [cm^-2]
; * ignored if WVL is incorrect
; defEM [INPUT; default=1e14] default EM to use prior to scaling [cm^-5]
; noph [INPUT] if set, FLUX is assumed to be in [ergs/s/cm^2], and
; no attempts are made to convert things to [ph/...]
; * forcibly set if WVL is incorrect
; thresh [INPUT; default=1e-2] dynamic range in output curves,
; relative to maximum in each curve
; _extra [INPUT ONLY] use this to pass defined keywords to subroutines
; LINEFLX [ABUND,TEMP,IKEV,EFFAR,WVLAR]
; ISMTAU [FH2,HE1,HEII,FANO]
;
;restrictions
; * EFFAR and WVLAR must be correctly defined, or else the output
; will be garbage
; * requires subroutines
; LINEFLX
; WHEE
; ISMTAU
;
;history
; vinay kashyap (Oct98)
; bug correction re 0 EMIS (VK; Nov98)
; bug correction re 0 NH (VK; AugMM)
;-

; usage
sze=size(emis) & nsze=n_elements(sze) & n=sze(nsze-2)
if n eq 0 then begin
print,'Usage: em=flux_to_em(emis,flux=flux,logT=logT,wvl=wvl,Z=Z,NH=NH,$'
print,' defEM=defEM,noph=noph,thresh=thresh;'
print,' LINEFLX:abund,temp,ikev,effar,wvlar; ISMTAU:fh2,he1,heII,fano)'
return,-1L
endif

; consistency checks
edim=sze(0) ;dimension of EMIS
if edim eq 0 then begin & mt=1L & mw=1L & endif ;scalar
if edim eq 1 then begin & mt=sze(1) & mw=1L & endif ;EMIS=EMIS(LOGT)
if edim eq 2 then begin & mt=sze(1) & mw=sze(2) & endif ;EMIS=EMIS(LOGT,WVL)
if edim gt 2 then begin
message,'emissivity array is too complex to understand',/info
return,0*emis
endif
nf=n_elements(flux) & nt=n_elements(logT) & nw=n_elements(wvl) & nz=n_elements(Z)
;
fx=fltarr(mw)+1. ;fluxes
if nf gt 0 and nf lt mw then fx(0:nf-1)=flux
if nf ge mw then fx=flux(0:mw-1)
;
tlog=interpol([4.,8.],mt) ;temperatures
if nt gt 0 and nt ne mt then tlog=interpol(logt,mt)
if nt eq mt then tlog=logt
;
if keyword_set(NH) then setabs=1 else setabs=0 ;include ISM absorption?
;
ww=findgen(mw) ;wavelengths
if nw eq 0 or nw ne mw then begin
setnoph=1 & setabs=0
endif else ww=wvl
;
zz=intarr(mw)+1 ;atomic numbers
if nz gt 0 and nz lt mw then zz(0:nf-1)=Z
if nz ge mw then zz=Z(0:mw-1)
;
if keyword_set(noph) then setnoph=1 ;everything in [ergs/...]
;
if not keyword_set(defEM) then dem=1d14 else dem=double(defEM)
;default EM [cm^-5]
if not keyword_set(thresh) then thresh=1e-2 ;dynamic range

; get optical depths
tau=fltarr(mw)
if keyword_set(setabs) then begin
tau=ismtau(ww,NH=NH,fH2=0,/Fano,_extra=e) < 69.
endif

; get the predicted fluxes
ct=dblarr(nt,nw)
for iw=0,nw-1 do begin
if edim eq 0 then line=[emis]
if edim eq 1 then line=[emis(*)]
if edim eq 2 then line=reform(emis(*,iw))
for it=0,nt-1 do begin
tmp=lineflx(line(it),tlog(it),ww(iw),zz(iw),noph=setnoph,dem=dem, _extra=e)
ct(it,iw)=tmp*exp(-tau(iw))
endfor
endfor

; compare to observed and scale EM
em=dblarr(nt,nw)-1
for iw=0,nw-1 do begin
if edim eq 0 then begin
pct=[ct] & line=[emis]
endif
if edim eq 1 then begin
pct=[ct(*)] & line=[emis(*)]
endif
if edim eq 2 then begin
pct=reform(ct(*,iw)) & line=reform(emis(*,iw))
endif
oo=where(pct ge thresh*max(pct),moo)
if max(pct) le 0 then moo=0L
if moo gt 0 then begin
frac=dem*fx(iw)/pct(oo)
em(oo,iw)=frac(*)
endif
endfor

return,em
end