Fitting a user-defined model to a spectrumΒΆ
damped_lya_profile is the function that generates the model and damped_lya_fitter creates the fitter class.
def damped_lya_profile(wave_to_fit,vs_n,am_n,fw_n,vs_b,am_b,fw_b,h1_col,h1_b,
h1_vel,d2h,resolution,single_component_flux=False,
return_components=False,
return_hyperfine_components=False):
"""
Computes a damped Lyman-alpha profile (by calling the functions
lya_intrinsic_profile_func and total_tau_profile_func) and convolves it
to the proper resolution.
"""
lya_intrinsic_profile = lya_intrinsic_profile_func(wave_to_fit,vs_n,am_n,fw_n,
vs_b,am_b,fw_b,single_component_flux=single_component_flux)
total_tau_profile = total_tau_profile_func(wave_to_fit,h1_col,h1_b,h1_vel,d2h)
lya_obs_high = lya_intrinsic_profile * total_tau_profile
## Convolving the data ##
fw = lya_rest/resolution
aa = make_kernel(grid=wave_to_fit,fwhm=fw)
lyman_fit = np.convolve(lya_obs_high,aa,mode='same')
return lyman_fit*1e14
def damped_lya_fitter(multisingle='multi'):
"""
Generator for Damped LyA fitter class
"""
myclass = pyspeckit.models.model.SpectralModel(damped_lya_profile,
11,
parnames=['vs_n','am_n','fw_n','vs_b','am_b','fw_b','h1_col',
'h1_b','h1_vel','d2h','resolution'],
parlimited=[(False,False),(True,False),(True,False),(False,False),
(True,False),(True,False),(True,True),(True,True),
(False,False),(True,True),(True,False)],
parlimits=[(14.4,14.6), (1e-16,0), (50.,0),(0,0), (1e-16,0),(50.,0),
(17.0,19.5), (5.,20.), (0,0),(1.0e-5,2.5e-5),(10000.,0)])
myclass.__name__ = "damped_lya"
return myclass
## Load the spectrum, register the fitter, and fit
spec = pyspeckit.Spectrum(xarr=wave_to_fit, data=flux_to_fit*1e14,
error=error_to_fit*1e14, doplot=False,
header=spec_header)
spec.Registry.add_fitter('damped_lya',damped_lya_fitter(),11)
spec.specfit(fittype='damped_lya',
guesses=[vs_n, am_n, fw_n, vs_b, am_b, fw_b, h1_col, h1_b,
h1_vel, d2h, resolution],
quiet=False, fixed=[False, False, False, False, False,
False, False, False, False, True, True])
## Check out the results
spec.specfit.parinfo.values
spec.specfit.parinfo.errors
reduced_chi2 = spec.specfit.chi2/spec.specfit.dof