Viewing contents of file '../idllib/astron/contrib/varosi/vlibm/allpro/ep_homog.pro'
;+
; NAME:
;	EP_homog
; PURPOSE:
;	Approximate the radiative transfer of photons in a
;	homogenous medium (dust & gas) of spherical/disk geometry, by
;	calling EP_Rad_Trans, and then compute equilibrium dust temperatures
;	and infrared emission using function Dust_Emission.
; CALLING:
;	EP_homog, SED_wave, SED_emit, SED_esc, Lumins, dust_emiss
; INPUTS:
;	SED_wave = array of SED wavelengths, in microns.
;
;	SED_emit = emitted Spectral Energy Distribution (SED) vs. wavelength,
;		in solar units (so that integral gives solar luminosities).
; OUTPUTS:
;	SED_esc = escaping SED vs. wavelength (solar units).
;
;	Lumins = structure variable accounting of all Luminosities emitted,
;		escaped, and absorbed in dust/gas.
;
;	dust_emiss = structure variable containing infrared emission spectra of
;		and temperature distribution of each dust phase.
;
; KEYWORD INPUTS:
;
;	KAPPA_DUST = structure variable with dust absorption/scattering coeffs.
;
;	GEOMETRY = structure variable:
;		geometry.R = radius of disk/cylinder.
;		geometry.H = height of disk/cylinder.
;		geometry.Reff = effective radius used for escape probability.
;
;	MEDIUM = structure variable:
;
;	SOURCE = string (upper or lowercase) indicating source of photons:
;		"C" : do case of central source,
;		"X" : do case of uniformly illuminating external source,
;			default is "U": uniformly distributed internal source.
;
;	/PLOT : plot intermediate results.
;	/VERBOSE : print informational messages about results.
;	/FULL_STRUCT : include all results in output structure variable "Lumins"
;
; EXTERNAL CALLS:
;	pro EP_Rad_Trans
;	function Dust_Emission
; PROCEDURE:
; HISTORY:
;	Written: Frank Varosi, RSTX @ NASA/GSFC, 1997
;-

pro EP_homog, SED_wave, SED_emit, SED_esc, Lumins, dust_emiss, TMAX=Tmax, $
				SOURCE=srctyp, KAPPA_DUST=Dust_Kappa, $
				FULL_STRUCT=fulls, VERBOSE=verbose,PLOT=plotit,$
				POWER_LAW_TDIST=ptdist, TDIST_PARM=tdistp, $
				GEOMETRY=geometry, CLUMPY=clumpy, MEDIUM=medium

if( geometry.H LE 0 ) then begin	;case of just simple sphere:
	geometry.Volume = (4*!pi/3) * geometry.R^3
	geometry.Reff = geometry.R
 endif else begin
	geometry.Volume = geometry.H * !pi * geometry.R^2
	geometry.Reff = 1.5 * geometry.R * geometry.H /(geometry.R + geometry.H)
  endelse

dfac = 	!cv.pc * !cv.pc * (!cv.pc/!cv.msun)

medium.dust_mass = geometry.Volume * medium.dust_mass_dens * dfac

medium.dust_gas_ratio = medium.dust_mass_dens/(medium.gas_dens*1.4*!cv.mh)

medium.fsil = ( medium.fsil < 1 ) > 0
medium.fc = (1-medium.fsil)
dust_fracs = [ medium.fc, medium.fsil ]
dmass_phases = [ 1, 0 ]

if keyword_set( verbose ) then begin
	print," "
	print," SOURCE: ",srctyp
	print," "
	print_struct, geometry, FORM=['g','11','3']
	print," "
	print_struct, medium, TR=[0,6], FORM=['g','11','3']
	print," "
	print_struct, medium, TR=[7,20], FORM=['g','11','3']
   endif

EP_Rad_Trans, SED_wave, SED_emit, SED_esc, Lumins, SOURCE=srctyp, FULL=fulls, $
			RADIUS=geometry.Reff*!cv.pc, KAPPA_DUST=Dust_Kappa, $
			DUST_MASS_DEN = dust_fracs * medium.dust_mass_dens, $
			GAS_DEN=medium.gas_dens, GAS_TEMP=medium.gas_Temp, $
			PLOT=plotit, VERB=verbose

if N_struct( tdistp ) eq 1 then begin
	Lindex = [ tdistp.Lindex_C, tdistp.Lindex_Sil ]
	dindex = [ tdistp.dindex_ICM, tdistp.dindex_CL ]
   endif

dust_emiss = Dust_Emission( Lumins, KAPPA=Dust_Kappa, TMAX=Tmax, POWER=ptdist, $
				PLOT=plotit, FULL=fulls, VERB=verbose, $
				DUST_MASS= dust_fracs # dmass_phases, $
					LINDEX=Lindex, DINDEX=dindex )
end