This page documents the Korg API for users of the package. Low-level functions that only people digging into the internals of Korg will be interested in can be found in the Developer documentation.
Top-level functions
If you are synthesize a spectrum Korg, these are the functions you will call. These functions are exported, so if you do using Korg
, you can call them unqualified (i.e. synthesize
instead of Korg.synthesize
).
Korg.synth
— Functionsynth(kwargs...)
This function creates a synthetic spectrum. It's easier to use than synthesize
, but it gives you less control. Unlike synthesize
, it returns a tuple of (wavelengths, rectified_flux, cntm)
(Wavelength in Å, rectified flux as a unitless number between 0 and 1, and continuum in erg/s/cm^5). Korg.synth
also provides shortcuts for some ways you might want to post-process the spectrum (applying a LSF, rotation, etc).
Keyword arguments
Teff
: effective temperature in K (default: 5000)logg
: surface gravity in cgs units (default: 4.5)m_H
: metallicity, [metals/H], (default: 0.0) (Seeformat_A_X
for precisely how this is interpreted.)alpha_H
: alpha enhancement, [α/H], (default:m_H
) (Seeformat_A_X
for precisely how this is interpreted.)- Any atomic symbol (e.g.
Fe
orC
) can be used to to specify a (solar relative, [X/H]) abundance. These overridem_H
andalpha_H
. Specifying an individual abundance means that the true metallicity and alpha will not correspond precisely to the values ofm_H
andalpha_H
. Seeformat_A_X
for details. linelist
: a linelist, (default:get_VALD_solar_linelist()
). See alsoread_linelist
.wavelengths
: a tuple of the start and end wavelengths (default: (5000, 6000)), or a vector of(λstart, λstop)
pairs. SeeWavelengths
for all the ways the wavelengths can be specified.rectify
: whether to rectify (continuum normalize) the spectrum (default: true)R
: resolution (default:Inf
, no LSF applied).R
can be a scalar, or a function from wavelength (in Å) to resolving power. Seeapply_LSF
for details on how to do this manually.vsini
: projected rotational velocity in km/s (default: 0). Seeapply_rotation
for details on how to do this manually.vmic
: microturbulent velocity in km/s (default: 1.0).synthesize_kwargs
: additional keyword arguments pass tosynthesize
.format_A_X_kwargs
: additional keyword arguments pass toformat_A_X
.
Korg.synthesize
— Functionsynthesize(atm, linelist, A_X, (λ_start, λ_stop); kwargs... )
synthesize(atm, linelist, A_X, wavelength_ranges; kwargs... )
Compute a synthetic spectrum. Returns a SynthesisResult
.
Arguments
atm
: the model atmosphere (seeread_model_atmosphere
)linelist
: A vector ofLine
s (seeread_linelist
,get_APOGEE_DR17_linelist
, andget_VALD_solar_linelist
).A_X
: a vector containing the A(X) abundances (log(X/H) + 12) for elements from hydrogen to uranium. (seeformat_A_X
)- The wavelengths at which to synthesize the spectrum. They can be specified either as a pair
(λstart, λstop)
, or as a list of pairs[(λstart1, λstop1), (λstart2, λstop2), ...]
. λ_start
: the lower bound (in Å) of the region you wish to synthesize.λ_stop
: the upper bound (in Å) of the region you wish to synthesize.λ_step
(default: 0.01): the (approximate) step size to take (in Å).
Example
to synthesize a spectrum between 5000 Å and 5100 Å, with all metal abundances set to 0.5 dex less than the solar value except carbon, except carbon, which we set to [C/H]=-0.25:
atm = read_model_atmosphere("path/to/atmosphere.mod")
linelist = read_linelist("path/to/linelist.vald")
A_X = format_A_X(-0.5, Dict("C" => -0.25))
result = synthesize(atm, linelist, A_X, 5000, 5100)
Optional arguments:
vmic
(default: 0) is the microturbulent velocity, $\xi$, in km/s.line_buffer
(default: 10): the farthest (in Å) any line can be from the provided wavelength range before it is discarded. If the edge of your window is near a strong line, you may have to turn this up.cntm_step
(default 1): the distance (in Å) between point at which the continuum opacity is calculated.hydrogen_lines
(default:true
): whether or not to include H lines in the synthesis.use_MHD_for_hydrogen_lines
(default:true
): whether or not to use the MHD occupation probability formalism for hydrogen lines. (MHD is always used for hydrogen bound-free absorption.)hydrogen_line_window_size
(default: 150): the maximum distance (in Å) from each hydrogen line center at which to calculate its contribution to the total absorption coefficient.mu_values
(default: 20): the number of μ values at which to calculate the surface flux, or a vector of the specific values to use when doing transfer in spherical geometry. Ifmu_points
is an integer, the values are chosen per Gauss-Legendre integration. If they are specified directly, the trapezoid rule is used for the astrophysical flux. The default values is sufficient for accuracy at the 10^-3 level. Note that if you are using the default radiative transfer scheme, with a plane-parallel model atmosphere, the integral over μ is exact, so this parameter has no effect. The points and weights are returned in themu_grid
field of the output.line_cutoff_threshold
(default:3e-4
): the fraction of the continuum absorption coefficient at which line profiles are truncated. This has major performance impacts, since line absorption calculations dominate more syntheses. Turn it down for more precision at the expense of runtime. The default value should effect final spectra below the 10^-3 level.electron_number_density_warn_threshold
(default:Inf
): if the relative difference between the calculated electron number density and the input electron number density is greater than this value, a warning is printed. By default, this warning is suppress (threshold isInf
) because it is very easily raised in cases where it is of no observable consequence. See alsoelectron_number_density_warn_min_value
, below.electron_number_density_warn_min_value
(default:1e-4
): The minimum value of the ratio of the electron number density to the total number density at which a warning is printed.return_cntm
(default:true
): whether or not to return the continuum at each wavelength. If this is false,solution.cntm
will benothing
.ionization_energies
, aDict
mappingSpecies
to their first three ionization energies, defaults toKorg.ionization_energies
.partition_funcs
, aDict
mappingSpecies
to partition functions (in terms of ln(T)). Defaults to data from Barklem & Collet 2016,Korg.default_partition_funcs
.equilibrium_constants
, aDict
mappingSpecies
representing diatomic molecules to the base-10 log of their molecular equilibrium constants in partial pressure form. Defaults to data from Barklem and Collet 2016,Korg.default_log_equilibrium_constants
.use_chemical_equilibrium_from
(default:nothing
): Takes another solution returned bysynthesize
. When provided, the chemical equilibrium solution will be taken from this object, rather than being recomputed. This is physically self-consistent only when the abundances,A_X
, and model atmosphere,atm
, are unchanged.molecular_cross_sections
(default:[]
): A vector of precomputed molecular cross-sections. SeeMolecularCrossSection
for how to generate these. If you are using the default radiative transfer scheme, your molecular cross-sections should cover 5000 Å only if your linelist does.tau_scheme
(default: "linear"): how to compute the optical depth. Options are "linear" and "bezier" (testing only–not recommended).I_scheme
(default:"linear_flux_only"
): how to compute the intensity. Options are"linear"
,"linear_flux_only"
, and"bezier"
."linear_flux_only"
is the fastest, but does not return the intensity values anywhere except at the top of the atmosphere. "linear" performs an equivalent calculation, but stores the intensity at every layer."bezier"
is for testing and not recommended.verbose
(default:false
): Whether or not to print information about progress, etc.
Korg.read_linelist
— Functionread_linelist(filename; format="vald", isotopic_abundances=Korg.isotopic_abundances)
Parse a linelist file, returning a vector of Line
s.
The format
keyword argument can be used to specify one of these linelist formats (default: "vald"
):
"vald"
for a VALD linelist. These can be either "short" or "long" format, "extract all" or "extract stellar". Air wavelengths will automatically be converted into vacuum wavelengths, and energy levels will be automatically converted from cm$^{-1}$ to eV."kurucz"
for an atomic or molecular Kurucz linelist (format=kurucz_vac if it uses vacuum wavelengths; be warned that Korg will not assume that wavelengths are vacuum below 2000 Å),"moog"
for a MOOG linelist (doesn't support broadening parameters or dissociation energies, assumed to be in vacuum wavelengths)."moog_air"
for a MOOG linelist in air wavelengths."turbospectrum"
for a Turbospectrum linelist in air wavelengths. Note that Korg doesn't make use of the (optional) orbital angular momentum quantum number, l, for the upper or lower levels, so it won't fall back on generic ABO recipes when the ABO parameters are not available. Korg's interpretation of thefdamp
parameter is also slightly different from Turbospectrum's. See the documentation of thevdW
parameter ofLine
for details. Korg will error if encounters an Unsoeld fudge factor, which it does not support.- "turbospectrum_vac" for a Turbospectrum linelist in vacuum wavelengths.
- "korg" for a Korg linelist (saved with hdf5). If the filename ends in
.h5
, this will be used by default.
For VALD and Turbospectrum linelists with isotope information available, Korg will scale log gf values by isotopic abundance (unless VALD has already pre-scaled them), using isotopic abundances from NIST ([Korg.isotopicabundances]). To use custom isotopic abundances, just pass `isotopicabundances` with the same structure: a dict mapping atomic number to a dict mapping from atomic weight to abundance.
Be warned that for linelists which are pre-scaled for isotopic abundance, the estimation of radiative broadening from log(gf) is not accurate.
See also: load_ExoMol_linelist
, save_linelist
.
Korg.load_ExoMol_linelist
— Functionload_ExoMol_linelist(specs, states_file, transitions_file, upper_wavelength, lower_wavelength)
Load a linelist from ExoMol. Returns a vector of Line
s, the same as read_linelist
.
Arguments
spec
: the species, i.e. the molecule that the linelist is forstates_file
: the path to the ExoMol states filetransitions_file
: the path to the ExoMol, transitions fileupper_wavelength
: the upper limit of the wavelength range to load (Å)lower_wavelength
: the lower limit of the wavelength range to load (Å)
Keyword Arguments
line_strength_cutoff
: the cutoff for the line strength (default: -15) used to filter the linelist. Seeapproximate_line_strength
for more information.T_line_strength
: the temperature (K) at which to evaluate the line strength (default: 3500.0)
This functionality is in beta.
Korg.read_model_atmosphere
— Functionread_model_atmosphere(filename)
Parse the provided model atmosphere file in MARCS ".mod" format. Returns either a PlanarAtmosphere
or a ShellAtmosphere
.
Korg.interpolate_marcs
— Functioninterpolate_marcs(Teff, logg, A_X; kwargs...)
interpolate_marcs(Teff, logg, m_H=0, alpha_m=0, C_m=0; kwargs...)
Returns a model atmosphere computed by interpolating models from MARCS ((Gustafsson+ 2008)[https://ui.adsabs.harvard.edu/abs/2008A&A...486..951G/abstract]). Along with Teff
and logg
, the atmosphere is specified by m_H
, alpha_m
, and C_m
, which can be automatically determined from an A_X
abundance vector (the recommended method, see format_A_X
). Note that the MARCS atmosphere models were constructed with the Grevesse+ 2007 solar abundances (Korg.grevesse_2007_solar_abundances
). This is handled automatically when A_X
is provided.
interpolate_marcs
uses three different interpolation schemes for different stellar parameter regimes. In the standard case the model atmosphere grid is the one generated for SDSS, transformed and linearly interpolated. For cool dwarfs (Teff
≤ 4000 K, logg
≥ 3.5), the grid is resampled onto unchanging tau_5000
values and interpolated with a cubic spline. For low-metallicity models (-5 ≤ m_H
< -2.5), a grid of standard composition (i.e. fixed alpha and C) atmospheres is used. (The microturbulence is 1km/s for dwarfs and 2km/s for giants and the mass for spherical models is 1 solar mass.) The interpolation method is the same as in the standard case. See Wheeler+ 2024 for more details and a discussion of errors introduced by model atmosphere interpolation. (Note that the cubic scheme for cool dwarfs is referred to as not-yet-implemented in the paper but is now available.)
keyword arguments
spherical
: whether or not to return a ShellAtmosphere (as opposed to a PlanarAtmosphere). By default true whenlogg
< 3.5.solar_abundances
: (default:grevesse_2007_solar_abundances
) The solar abundances to use whenA_X
is provided instead ofM_H
,alpha_M
, andC_M
. The default is chosen to match that of the atmosphere grid, and if you change it you are likely trying to do something else.clamp_abundances
: (default:false
) allowed only when specifyingA_X
. Whether or not to clamp the abundance parameters to be within range to avoid throwing an out of bounds error.perturb_at_grid_values
(default:true
): whether or not to add or a subtract a very small number to each parameter which is exactly at a grid value. This prevents null derivatives, which can cause problems for minimizers.resampled_cubic_for_cool_dwarfs
(default:true
): whether or not to used specialized method for cool dwarfs.archives
: A tuple containing the atmosphere grids to use. For testing purposes.
Korg.format_A_X
— Functionformat_A_X(default_metals_H, default_alpha_H, abundances; kwargs... )
Returns a 92 element vector containing abundances in $A(X)$ ($\log_{10}(X/H) + 12$) format for elements from hydrogen to uranium.
Arguments
You can specify abundance with these positional arguments. All are optional, but if default_alpha_H
is provided, default_metals_H
must be as well.
default_metals_H
(default: 0), i.e. [metals/H] is the $\log_{10}$ solar-relative abundance of elements heavier than He. It is overridden bydefault_alpha
andabundances
on a per-element basis.default_alpha_H
(default: same asdefault_metals_H
), i.e. [alpha/H] is the $\log_{10}$ solar-relative abundance of the alpha elements (Seealpha_elements
, below). It is overridden byabundances
on a per-element basis.abundances
is aDict
mapping atomic numbers or symbols to [$X$/H] abundances. (Setsolar_relative=false
to use $A(X)$ abundances instead.) These overridedefault_metals_H
. This is the only way to specify an abundance of He that is non-solar.
Keyword arguments
solar_relative
(default: true): When true, interpret abundances as being in [$X$/H] ($\log_{10}$ solar-relative) format. When false, interpret them as $A(X)$ abundances, i.e. $A(x) = \log_{10}(n_X/n_\mathrm{H}) + 12$, where $n_X$ is the number density of $X$. Note that abundances not specified default to the solar value still depend on the solar value, as they are set according todefault_metals_H
anddefault_alpha_H
.solar_abundances
(default:Korg.asplund_2020_solar_abundances
) is the set of solar abundances to use, as a vector indexed by atomic number.Korg.asplund_2009_solar_abundances
andKorg.grevesse_2007_solar_abundances
are also provided for convenience.alpha_elements
(default:Korg.default_alpha_elements
): vector of atomic numbers of the alpha elements. (Useful since conventions vary.)
Built-in linelists
Korg.get_APOGEE_DR17_linelist
— Functionget_APOGEE_DR17_linelist(; include_water=true)
The APOGEE DR 17 linelist. It ranges from roughly 15,000 Å to 17,000 Å. It is nearly the same at the DR 16 linelist described in Smith+ 2021.
Korg.get_GALAH_DR3_linelist
— Functionget_GALAH_DR3_linelist()
The GALAH DR 3 linelist (also used for DR 4). It ranges from roughly 4,675 Å to 7,930 Å. This linelist is based on, but distinct from Heiter 2021. See Buder et al. 2021 for details.
Korg.get_GES_linelist
— Functionget_GES_linelist()
The Gaia-ESO survey linelist from Heiter et al. 2021. This linelist contains > 15 million lines, which means that it can take a while to synthesize spectra. If you don't need molecular lines, you can set include_molecules=false
to speed things up.
Keyword Arguments
include_molecules
(default:true
): whether to include molecular lines.
Fitting
Korg.Fit.fit_spectrum
— Functionfit_spectrum(obs_wls, obs_flux, obs_err, linelist, initial_guesses, fixed_params; kwargs...)
Find the parameters and abundances that best match a rectified observed spectrum.
Arguments:
obs_wls
: the wavelengths of the observed spectrum in Å. These must be vacuum wavelengths.obs_flux
: the rectified flux of the observed spectrumobs_err
: uncertainty influx
linelist
: a linelist to use for the synthesisinitial_guesses
: a NamedTuple specifying initial guesses for the parameters to be fit. See "Specifying parameters" below.fixed_params
: a NamedTuple specifying parameters to be held fixed. See "Specifying parameters" below.
initial_guesses
and fixed_params
can also be specified as Dicts instead of NamedTuples, which is more convenient when calling Korg from python.
Specifying parameters
Parameters are specified as named tuples or dictionaries. Named tuples look like this: (Teff=5000, logg=4.5, m_H=0.0)
. Single-element named tuples require a semicolon: (; Teff=5000)
.
Required parameters
Teff
and logg
must be specified in either initial_guesses
or fixed_params
.
Optional Parameters
These can be specified in either initial_guesses
or fixed_params
, but if they are not default values are used.
m_H
: the metallicity of the star, in dex. Default:0.0
alpha_H
: the alpha enhancement of the star, in dex. Default:m_H
. Note that, because of the parameter range supported byKorg.interpolate_marcs
, only values within ±1 ofm_H
are supported.vmic
: the microturbulence velocity, in km/s. Default:1.0
vsini
: the projected rotational velocity of the star, in km/s. Default:0.0
. SeeKorg.apply_rotation
for details.epsilon
: the linear limb-darkening coefficient. Default:0.6
. Used for applying rotational broadening only. SeeKorg.apply_rotation
for details.- Individual elements, e.g.
Na
, specify the solar-relative ([X/H]) abundance of that element.
Keyword arguments
R
, the resolution of the observed spectrum. This is required. It can be specified as a function of wavelength, in which case it will be evaluated at the observed wavelengths.windows
is a vector of wavelength pairs, each of which specifies a wavelength "window" to synthesize and contribute to the total χ². If not specified, the entire spectrum is used. Overlapping windows are automatically merged.adjust_continuum
(default:false
) if true, adjust the continuum with the best-fit linear correction within each window, minimizing the chi-squared between data and model at every step of the optimization.wl_buffer
is the number of Å to add to each side of the synthesis range for each window.time_limit
is the maximum number of seconds to spend in the optimizer. (default:10_000
). The optimizer will only checks against the time limit after each step, so the actual wall time may exceed this limit.precision
specifies the tolerance for the solver to accept a solution. The solver operates on transformed parameters, soprecision
doesn't translate straightforwardly to Teff, logg, etc, but the default value,1e-4
, provides a theoretical worst-case tolerance of about 0.15 K inTeff
, 0.0002 inlogg
, 0.0001 inm_H
, and 0.0004 in detailed abundances. In practice the precision achieved by the optimizer is about 10x bigger than this.postprocess
can be used to arbitrarilly transform the synthesized (and LSF-convolved) spectrum before calculating the chi2. It should take the formpostprocess(flux, data, err)
and write its changes in-place to the flux array.LSF_matrix
: this can be provedided along withsynthesis_wls
in place of specifyingR
if you have a precomputed custom LSF matrix.synthesis_wls
: seeLSF_matrix
above. This can be a Korg.Wavelengths object or any arguments that can be passed to its constructor, e.g. a range or vector of ranges. Wavelengths are in Å.- Any additional keyword arguments will be passed to
Korg.synthesize
when synthesizing the spectra for the fit.
Returns
A NamedTuple with the following fields:
best_fit_params
: the best-fit parametersbest_fit_flux
: the best-fit flux, with LSF applied, resampled, and rectified.obs_wl_mask
: a bitmask forobs_wls
which selects the wavelengths used in the fit (i.e. those in thewindows
)solver_result
: the result object fromOptim.jl
trace
: a vector of NamedTuples, each of which contains the parameters at each step of the optimization. This is empty for single parameter fits, because the underlying solver doesn't supply it.covariance
: a pair(params, Σ)
whereparams
is vector of parameter name (providing an order), andΣ
is an estimate of the covariance matrix of the parameters. It is the approximate inverse hessian of the log likelihood at the best-fit parameter calculated by the BGFS algorithm, and should be interpreted with caution.
This function takes a long time to compile the first time it is called. Compilation performance is significantly better on Julia 1.10+ than previous versions, so if you are using an older version of Julia, you may want to upgrade.
Korg.Fit.ews_to_abundances
— Functionews_to_abundances(atm, linelist, A_X, measured_EWs; kwargs... )
Functions using equivalent widths have been dissabled while we fix some bugs.
Compute per-line abundances on the linear part of the curve of growth given a model atmosphere and a list of lines with equivalent widths.
Arguments:
atm
: the model atmosphere (seeKorg.read_model_atmosphere
andKorg.interpolate_marcs
).linelist
: A vector ofKorg.Line
s (seeKorg.read_linelist
). The lines must be sorted by wavelength.A_X
: a vector containing the A(X) abundances (log(nX/nH) + 12) for elements from hydrogen to uranium (seeKorg.format_A_X
). All syntheses are done with these abundances, so if the resulting abundances deviate significantly from these, you may wish to iterate.measured_EWs
: a vector of equivalent widths (in mÅ)
Returns
A vector of abundances (A(X) = log10(n_X/n_H) + 12
format) for each line in linelist
.
Optional arguments:
wl_step
(default: 0.01) is the resolution in Å at which to synthesize the spectrum around each line.ew_window_size
(default: 2): the farthest (in Å) to consider equivalent width contributions for each line. It's very important that this is large enough to include each line entirely.blend_warn_threshold
(default: 0.01) is the minimum absorption between two lines allowed before triggering a warning that they may be blended. All other keyword arguments are passed toKorg.synthesize
when synthesizing each line.
Korg.Fit.ews_to_stellar_parameters
— Functionews_to_stellar_parameters(linelist, measured_EWs, [measured_EW_err]; kwargs...)
Functions using equivalent widths have been dissabled while we fix some bugs.
Find stellar parameters from equivalent widths the "old fashioned" way. This function finds the values of $T_\mathrm{eff}$, $\log g$, $v_{mic}$, and [m/H] which satisfy the following conditions (using a Newton-Raphson solver):
- The slope of the abundances of neutral lines with respect to lower excitation potential is zero.
- The difference between the mean abundances of neutral and ionized lines is zero.
- The slope of the abundances of neutral lines with respect to reduced equivalent width is zero.
- The difference between the mean abundances of all lines and the model-atmosphere input [m/H] is zero. Here the "slope" refers to the slope of a linear fit to the abundances of the lines in question.
Arguments:
linelist
: A vector ofKorg.Line
objects (seeKorg.read_linelist
). The lines must be sorted by wavelength.measured_EWs
: a vector of equivalent widths (in mÅ).measured_EW_err
(optional): the uncertainty inmeasured_EWs
. If not specified, all lines are assumed to have the same uncertainty. These uncertainties are used when evaluating the equations above, and are propagated to provide uncertainties in the resulting parameters.
Returns:
A tuple containing:
- the best-fit parameters:
[Teff, logg, vmic, [m/H]]
as a vector - the statistical uncertainties in the parameters, propagated from the uncertainties in the equivalent widths. This is zero when EW uncertainties are not specified.
- the systematic uncertainties in the parameters, estimated from the scatter the abundances computed from each line not accounted for by the EW uncertainties.
Keyword arguments:
Teff0
(default: 5000.0) is the starting guess for Tefflogg0
(default: 3.5) is the starting guess for loggvmic0
(default: 1.0) is the starting guess for vmic. Note that this must be nonzero in order to avoid null derivatives. Very small values are fine.m_H0
(default: 0.0) is the starting guess for [m/H]tolerances
(default:[1e-3, 1e-3, 1e-3, 1e-3]
) is the tolerance for the residuals each equation listed above. The solver stops when all residuals are less than the corresponding tolerance.max_step_sizes
(default:[1000.0, 1.0, 0.3, 0.5]
) is the maximum step size to take in each parameter direction. This is used to prevent the solver from taking too large of a step and missing the solution. Be particularly cautious with the vmic (third) parameter, as the unadjusted step size is often too large.parameter_ranges
(default:[(2800.0, 8000.0), (-0.5, 5.5), (1e-3, 10.0), (-2.5, 1.0)]
) is the allowed range for each parameter. This is used to prevent the solver from wandering into unphysical parameter space, or outside the range of the MARCS grid supported byKorg.interpolate_marcs
. The default ranges $T_\mathrm{eff}$, $\log g$, and [m/H] are the widest supported by the MARCS grid. Note that vmic must be nonzero in order to avoid null derivatives.callback
: is a function which is called at each step of the optimization. It is passed three arguments:- the current values of the parameters
- the residuals of each equation being solved
- the abundances of each line computed with the current parameters. You can pass a callback function, to e.g. make a plot of the residuals at each step.
max_iterations
(default: 30) is the maximum number of iterations to allow before stopping the optimization.
Secondary functions
You don't need use these to synthesize spectra, but they might be relevant depending on what you are doing.
Postprocessing
These are used to transform observational or synthetic spectra.
Korg.apply_LSF
— Functionapply_LSF(flux, wls, R; window_size=4)
Applies a gaussian line spread function the the spectrum with flux vector flux
and wavelengths wls
in any format accepted by synthesize, e.g. as a pair (λstart, λstop)
) with constant spectral resolution (, $R = \lambda/\Delta\lambda$, where $\Delta\lambda$ is the LSF FWHM. The window_size
argument specifies how far out to extend the convolution kernel in standard deviations.
For the best match to data, your wavelength range should extend a couple $\Delta\lambda$ outside the region you are going to compare.
If you are convolving many spectra defined on the same wavelenths to observational resolution, you will get much better performance using compute_LSF_matrix
.
Keyword Arguments
window_size
(default: 4): how far out to extend the convolution kernel in units of the LSF width (σ, not HWHM)renormalize_edge
(default:true
): whether or not to renormalize the LSF at the edge of the wl range. This doen't matter as long assynth_wls
extends to large and small enough wavelengths.
This is a naive, slow implementation. Do not use it when performance matters.
apply_LSF
will have weird behavior if your wavelength grid is not locally linearly-spaced. It is intended to be run on a fine wavelength grid, then downsampled to the observational (or otherwise desired) grid.
Korg.compute_LSF_matrix
— Functioncompute_LSF_matrix(synth_wls, obs_wls, R; kwargs...)
Construct a sparse matrix, which when multiplied with a flux vector defined over wavelenths synth_wls
, applies a gaussian line spead function (LSF) and resamples to the wavelenths obswls
.
Arguments
synth_wls
: the synthesis wavelengths in any form recognized bysynthesize
, e.g. a pair containing a lower and upper bound in Å, a vector of pairs, or a vector of Julia range objects.obs_wls
: the wavelengths of the observed spectrumR
: the resolving power, $R = \lambda/\Delta\lambda$
Keyword Arguments
window_size
(default: 4): how far out to extend the convolution kernel in units of the LSF width (σ, not HWHM)verbose
(default:true
): whether or not to emit warnings and information to stdout/stderr.step_tolerance
: the maximum difference between adjacent wavelengths insynth_wls
for them to be considered linearly spaced. This is only used ifsynth_wls
is a vector of wavelengths rather than a range or vector or ranges.
For the best match to data, your wavelength range should extend a couple $\Delta\lambda$ outside the region you are going to compare.
Korg.apply_LSF
can apply an LSF to a single flux vector efficiently. This function is relatively slow, but one the LSF matrix is constructed, convolving spectra to observational resolution via matrix multiplication is fast.
Korg.air_to_vacuum
— Functionair_to_vacuum(λ; cgs=λ<1)
Convert λ from an air to vacuum. λ is assumed to be in Å if it is ⩾ 1, in cm otherwise. Formula from Birch and Downs (1994) via the VALD website.
See also: vacuum_to_air
.
Korg.vacuum_to_air
— Functionvacuum_to_air(λ; cgs=λ<1)
convert λ from a vacuum to air. λ is assumed to be in Å if it is ⩾ 1, in cm otherwise. Formula from Birch and Downs (1994) via the VALD website.
See also: air_to_vacuum
.
Misc
Korg.blackbody
— Functionblackbody(T, λ)
The value of the Planck blackbody function for temperature T
at wavelength λ
[cm].
Korg.prune_linelist
— Functionprune_linelist(atm, linelist, A_X, wls...; threshold=1.0, sort=true, synthesis_kwargs...)
Return the vector containing the strongest lines in linelist
, (optionally) sorted by approximate equivalent width.
Arguments
atm
: the atmosphere modellinelist
: the linelist (a vector ofLine
objects)A_X
: the abundance of each element (seeformat_A_X
)wls...
: the wavelength ranges to synthesize over. These are specified the same way as thewls
forsynthesize
.
Keyword Arguments
threshold=0.1
: The threshold ratio in the line center absorption to the continuum absorption computed at the photosphere for a line to be included in the returned list.0.1
is a reasonable default for getting a sense of what might be measurable in a high-res, high-quality spectrum, but it should not be used to create a linelist for synthesis.sort_by_EW=true
: Iftrue
, the returned linelist will be sorted by approximate reduced equivalent width. Iffalse
, the linelist will be in wavelength order. Leaving the list in wavelength order is much faster, but sorting by strength is useful for visualizing the strongest lines.verbose=true
: Iftrue
, a progress bar will be displayed while measuring the EWs. All other kwargs are passed to internal calls tosynthesize
.max_distance=0.0
, how far fromwls
lines can be (in Å) before they are excluded from the returned list.
While this function can be used to prune a linelist for synthesis, the default behavior too aggressive for this purpose. Set a much lower threshold (e.g. threshold=1e-4
) and use sort_by_EW=false
if you are pruning the linelist to speedup synthesis. Note that Korg will dynamically choose which lines to include even if you use a large linelist (see the line_cutoff_threshold
keyword argument to synthesize
).
See also merge_close_lines
if you are using this for plotting.
Korg.merge_close_lines
— Functionmerge_close_lines(linelist; merge_distance=0.2)
Produce a list of the species and wavelengths of the lines in linelist
, merging lines of the same species that are within merge_distance
(default: 0.2 Å). This is useful for labeling lines in a plot after running prune_linelist
.
Arguments
linelist
: the linelist (a vector ofLine
objects)
Keyword Arguments
merge_distance=0.2
: The maximum distance in Å between lines of the same species to be merged into a single entry
Returns
A vector of tuples (wl, wl_low, wl_high, species)
where wl
is gf-weighted wavelength of each set of merged lines (Å), wl_low
and wl_high
are their highest and lowest wavelength, and species
is a string (not a Korg.Species
) identifying the species of the line. These will be in wavelength order.