ConjugateSpectrum¶
- class screens.conjspec.ConjugateSpectrum(conjspec, tau, fd, noise=None, d_eff=None, mu_eff=None, theta=None, magnification=None)[source] [edit on github]¶
Bases:
object
Conjugate spectrum and methods to fit it.
The code is meant to be agnostic to which axes are which, but some may assume a shape of
(..., doppler_axis, delay_axis)
.- Parameters:
- conjspec
ndarray
Fourier transform of a dynamic spectrum.
- fd
Quantity
Doppler factors of the conjugate spectrum. Normally time conjugate but can be arbitrary (e.g., conjugate of
f*t
). Should have the the proper shape to broadcast withconjspec
.- tau
Quantity
Delays of the conjugate spectrum. Should have the proper shape to broadcast with
dynspec
.- noisefloat
The uncertainty in the real and imaginary components of the conjugate spectrum.
- d_eff
Quantity
Assumed effective distance. This is used throughout and not fit, but can be treated as a scaling parameters.
- mu_eff
Quantity
, optional Initial guess for the effective proper motion,
v_eff/d_eff
.- theta
Quantity
, optional Grid of theta angles to use for modelling the dynamic spectrum. Probably more usefully calculated later using
theta_grid
.- magnification
Quantity
, optional Magnifications at each
theta
. More typically inferred by fitting the secondary spectrum.
- conjspec
Attributes Summary
Secondary spectrum, i.e., the power of the conjugate spectrum.
Methods Summary
from_dynamic_spectrum
(dynspec[, normalization])Create a conjugate spectrum from a dynamic one.
locate_mu_eff
([mu_eff_trials, power, verbose])Try reproducing the secondary spectrum for a range of proper motion.
model
([magnification, mu_eff, conserve])Calculate the expected conjugate spectrum for given parameters.
theta_grid
([oversample_tau, oversample_fd])Calculate a grid of theta for modelling the secondary spectrum.
theta_theta
([mu_eff, conserve, theta_grid])Project the secondary spectrum into theta-theta space.
Attributes Documentation
- secspec¶
Secondary spectrum, i.e., the power of the conjugate spectrum.
Methods Documentation
- classmethod from_dynamic_spectrum(dynspec, normalization='mean', **kwargs)[source] [edit on github]¶
Create a conjugate spectrum from a dynamic one.
Easiest if the input is a
DynamicSpectrum
instance.By passing in an explicit time axis using
t
, one can get a different delay factor conjugate. Particularly useful witht=f*t
, which takes into account the frequency dependence of the time variation of scintles.Note that the dynamic spectrum is assumed to have shape
(..., time_axis, frequency_axis)
.- Parameters:
- dynspecarray_like or
DynamicSpectrum
Input dynamic spectrum for which the fourier transform will be calculated. If it has attributes
f
,t
,d_eff
,theta
,magnification
, andnoise
, those will be used as default inputs. TODO:noise
is likely wrong!- normalization‘mean’ or None
Normalize such that the 0, 0 element equals the mean of the dynamic spectrum.
- **kwargs
Other arguments to initialize the conjugate spectrum.
- dynspecarray_like or
- locate_mu_eff(mu_eff_trials=None, power=True, verbose=False)[source] [edit on github]¶
Try reproducing the secondary spectrum for a range of proper motion.
For each proper motion, construct a theta-theta array, calculte the largest eigenvalue and use the corresponding eigenvector as the model one-dimensional screen.
- Parameters:
- mu_eff_trials
Quantity
Proper motions to try.
- powerbool
Whether to just decompose the powers or the full complex secondary spectrum. The former is faster and varies less quickly with proper motion, so can be used to find the global minimum before re-running on the full complex values.
- verbosebool
Whether or not to give summary statistics for each trial.
- mu_eff_trials
- Returns:
- curvature
QTable
Table with the following columns: -
mu_eff
: Input proper motions. -theta
: grid in theta used. -w
: Largest eigenvalue. -recovered
: corresponding eigenvector, i.e., magnifications. -th_ms
: Mean-square residual in theta-theta space. -ndof
: degrees of freedomn_dynspec - n_theta - 2
. -redchi2
: reduced chi2((dynspec - model)/noise)**2/ndof
.
- curvature
Notes
The resulting table is also stored on the instance, as
curvature
.Note that the noise in the secondary spectrum is currently ignored.
- model(magnification=None, mu_eff=None, conserve=False)[source] [edit on github]¶
Calculate the expected conjugate spectrum for given parameters.
- theta_grid(oversample_tau=2, oversample_fd=4, **kwargs)[source] [edit on github]¶
Calculate a grid of theta for modelling the secondary spectrum.
Wraps
screens.fields.theta_grid
with defaults from the class. See that function for details. Oversampling is set relatively high here, since covariances between the theta are not as important as for fitting a dynamic spectrum.Note that this does not set the angles on the class, so typical usage is
ds.theta = ds.theta_grid()
.
- theta_theta(mu_eff=None, conserve=False, theta_grid=True, **kwargs)[source] [edit on github]¶
Project the secondary spectrum into theta-theta space.
For a given grid in
theta
(possibly calculated) and a set of pairs found usingscreens.fields.theta_theta_indices
, interpolate in the secondary spectrum to find the theta-theta arrays.- Parameters:
- mu_eff
Quantity
, optional Effective proper motion to use. Defaults to that stored on the instance. Will update the instance if given.
- conservebool
Whether to conserve flux per surface area. Doing so reduces sensitivity to points near the axes, but means one cannot directly use any eigenvectors directly in constructing dynamic spectra.
- theta_gridbool, optional
Whether to calculate a new theta grid, or use the one stored on the instance. By default, calculate it only if
mu_eff
is passed in. IfTrue
, this will update the grid stored on the instance.- **kwargs
Any further arguments are passed on to
theta_grid
- mu_eff