Viewing contents of file '../idllib/uit/pro/clustmod.pro'
pro clustmod, w, f, age, expon, log_Z, init_m, Z_evol=Z_evol, $
OldKurucz = oldkurucz, OldMaeder = oldmaeder, Up_mass = up_mass, $
Mass_Range= mass_range, NoInterp = NoInterp, Silent = Silent, $
Mdot = Mdot
;+
; NAME:
; CLUSTMOD
; PURPOSE:
; Return the flux of a model star cluster as a function of mass at a
; specified age and at specified wavelengths
;
; CALLING SEQUENCE:
; CLUSTMOD, w, f, age, expon, [ log_Z, init_M, Mass_Range = ,/Silent
; Mdot =, /NoInterp, /OldMaeder , /OldKurucz , Z_evol = ]
;
; INPUTS:
; W - wavelength(s) (A) at which to evaluate fluxes. Set W = -1 to
; evaluate Lyman continuum fluxes. Wavelengths may be
; modified upon output if the /NoInterp keyword is set.
; AGE - age (in Myr) of the the star cluster, scalar
;
; OPTIONAL INPUT:
; EXPON - power law exponent, uusally negative, scalar or vector
; The number of values in expon equals the number of different
; power-law components in the IMF
; A Saltpeter IMF has a scalar value of expon = -1.35
; LOG_Z - Log of the metallicity relative to solar. Values of
; LOG_Z = 0.3, 0, -0.4,-0.7, or -1.5 are preferred since these
; correspond exactly to available evolutionary tracks. User
; prompted if LOG_Z not supplied
;
; OPTIONAL INPUT-OUTPUT
; INIT_M - Masses at which to evaluate cluster flux. If not suppled
; or set to 0, then INIT_MASS is set equal to the masses used
; in the evolution grid of Schaller et al. plus a 20 element
; grid for the post main-sequence. Masses are truncated to
; those with lifetimes less than AGE.
;
; KEYWORD INPUT:
; MASS_RANGE - vector containing upper and lower limits of the IMF and
; masses where the IMF exponent changes. The number of values
; in mass_range should be one more than in expon. The values
; in mass_range should be monotonically increasing.
;
; OPTIONAL KEYWORD INPUTS:
; Z_Evol - scalar giving metal abundance to use with evolutionary tracks.
; Allowed values are 0.04, 0.02, 0.004, 0.008 or 0.001 for the
; new Schaller et al. tracks or 0.002, 0.005, 0.02, or 0.04 if
; the OldMaeder keyword is used. If not supplied, then the
; value of log_Z is used to select the track with the closest
; metallicity.
; Mdot - scalar either 1 (default) or 2. If Mdot = 2 then
; CLUSTMOD uses evolutionary tracks with the double the mass loss
; rate of the deJager parameterization whenever the winds are
; driven by radiation pressure (main-sequence + Wolf-Rayet phases)
; NoInterp - If this keyword set, then the input wavelengths are rounded
; off to the nearest Kurucz model wavelengths. This avoids
; interpolation in wavelength and saves considerable time.
; OldKurucz - If set and non-zero, then the old (1979) Kurucz models are
; used.
; OldMaeder - If set and non-zero, then the old (1989) Maeder tracks are
; are used.
; Silent - If set and non-zero, the informational messages are suppressed
;
; OUTPUT:
; F - Flux of star cluster at specified wavelengths and masses
; dimension NWAVE by NMASS. Flux units are erg s-1 A-1 cm-2
; per solar mass in the cluster. If W = -1 then F is the
; number of Lyman continuum photons (x 10(45)) emitted per
; solar mass in the cluster.
;
; EXAMPLE:
; Compute the flux at 1600 A for a solar metallicity 10Myr old model at
; the masses used in the Schaller et al. tracks. Use the IMF of
; Lequex (1981) which has an exponent of -0.6 between 0.007 and 1.8
; Msun, and an exponent of -1.7 between 1.8 and 120 Msun.
;
; IDL> clustmod, 1600, f, 10, [-.6,-1.7], 0, init_m, mass=[0.007,1.8,120]
;
; REVISION HISTORY:
; Written Wayne Landsman September, 1989
; Update for new Kurucz models July, 1992
; Added /NoInterp keyword September, 1992
;-
On_error,2
if N_params() LT 3 then begin
print,'Syntax - clustmod, w, f, age, [ expon, log_z, init_m, Mass_range ='
print,' Mdot = ,/Silent ,Z_evol = , /Oldkurucz, /OldMaeder ]'
return
endif
if ( N_elements(age) NE 1 ) then $
message,'A scalar star cluster age (third parameter) must be supplied'
rd_logz = N_elements(log_z) NE 1
Read_Z: if rd_logz then read, $
'Enter log of the abundance (0.3,0,-0.4,-0.7 or -1.5 preferred): ',log_z
if ( N_elements(expon) EQ 0 ) then expon = -1.35 ;Default IMF exponent
if not keyword_set( mass_range ) then $
if N_elements(expon) EQ 1 then mass_range = [0.1,110] $
else message,'ERROR - MASS_RANGE keyword not supplied
if not keyword_set( Mdot ) then mdot = 1
Ncomp = N_elements(expon)
if N_elements( mass_range) ne Ncomp+1 then message, $
'ERROR - MASS_RANGE vector must contain ' + strtrim(Ncomp,2) + ' elements'
up_mass = mass_range(Ncomp)
nwave = N_elements(w)
lyc = (w(0) EQ -1) ;Lyman continuum photons?
if N_elements(init_m) EQ 0 then set_mass = 1 else begin
if init_m(0) EQ 0 then set_mass = 1 else set_mass = 0
endelse
newkurucz = not keyword_set( OLDKURUCZ )
if not keyword_set(Z_evol) then begin
if keyword_set( OldMaeder) then zevol = [0.002,0.005,0.02,0.04] $
else zevol = [0.04, 0.02, 0.008, 0.004, 0.001]
tabinv, alog10( zevol/0.02), log_Z, i
Z_evol = zevol(nint(i))
endif
if set_mass then begin ;Get grid of masses for specified age
init_m = mass_grid( age(0),Z = Z_evol, Up_mass = up_mass, Mdot = mdot, $
OldMaeder = oldmaeder ) ;Fixed 10-2-92
if lyc then $ ;For Lyman continuum ignore low mass models
init_m= init_m(3:*)
endif
if newkurucz and (not lyc) and keyword_set( NoInterp) then begin
get_kur_wave, wmod
tabinv, wmod, w, W_index
W_index = nint(W_index)
w = wmod( W_index)
endif
; Compute evolutionary tracks for specified masses
starmod, age, init_m, m, l, log_teff, seq, $
Mdot = mdot, Z=Z_evol, OldMaeder = OldMaeder
age = age(0)
; Remove masses that have already supernovaed by the specified age
bad = where( init_m GT upper_mass( age,Z = Z_evol, MDot = mdot, $
Oldmaeder =oldmaeder), Nbad )
if ( Nbad GT 0 ) then begin
ngood = min(bad)
m = m(*,0:ngood)
log_teff = log_teff(*,0:ngood)
l = l(*,0:ngood)
init_m = init_m(0:ngood)
endif
if not keyword_set(SILENT) then begin
print,'CLUSTMOD - Current age '+ string(age,f='(f8.2)') + ' Million Years'
print,'CLUSTMOD - Current Mass Vector',init_m,f='(A,8F6.2,/,10(13f6.2,/))'
endif
teff = 10.^( log_teff )
r = star_rad( teff, l ) ;Get stellar radius
log_g, g, m, log_teff, l ;Get log surface gravity
frac = imf( init_m, expon, mass_range)/init_m ;Use IMF to obtain fractions
nmass = N_elements( init_m )
f = fltarr( nwave, nmass )
r2 = r^2
for j = 0,nmass-1 do begin
if teff(j) GT 100. then begin
if lyc then begin ;Lyman continuum flux
if newkurucz then kurucz_lyc, flux, teff(j), g(j), log_Z else $
kurucz_lyc_1979, flux, teff(j), g(j), log_Z
flux = flux/1e33
endif else $ ;Kurucz model flux
if newkurucz then $
kur_int, w, flux, teff(j), g(j), log_Z, W_index = w_index else $
kur_int_1979, w, flux, teff(j), g(j), log_Z
f(0,j) = flux*frac(j)*r2(j)
endif
endfor
f = f*4*!pi
if lyc then f = f/1e12
if N_elements( lun ) EQ 1 then free_lun,lun
return
end