mathtools

Mathematical utility objects and functions.

_zero_function(x)[source]
_one_function(x)[source]
_fast_eval_interp1d_lin(x, xp0, dxp, fp)[source]

JIT fast linear interpolation.

_fast_eval_interp1d_cubic(x, xp0, dxp, fp, y2p)[source]

JIT fast cubic spline.

_spline_cubic_natural(fp, dxp)[source]

Constructs the second derivative array with natural boundary conditions.

As coded in Numerical Recipes in C Chaper 3.3, but specialized for equally spaced input data.

Parameters:
  • fp (ndarray) – 1D data points array to interpolate.

  • dxp (float) – Spacing of x points.

Returns:

y2p – Second derivative array. Same size as fp.

Return type:

ndarray

_spline_cubic_notaknot(fp, dxp)[source]

Constructs the second derivative array with not-a-knot boundary conditions:

\[y''_0 - 2 y''_1 + y''_2 = 0\]
Parameters:
  • fp (ndarray) – 1D data points array to interpolate.

  • dxp (float) – Spacing of x points.

Returns:

y2p – Second derivative array. Same size as fp.

Return type:

ndarray

mypoly1d(coef, x)[source]

My simple power series polynomial calculator.

Parameters:
  • coef (ndarray) – Coefficient array in increasing power starting with the constant.

  • x (ndarray) – Array to calculate polynomial.

Returns:

results – Polynomial calculated at x.

Return type:

ndarray

block_covariance_of_square(mean, var, cov)[source]

Return the block covariance of x^2, i.e. \(<x_i^2 x_j^2> - <x_i^2><x_j^2>\). Compatible with blockdim argument of SubsampleCov.get_mean_n_cov().

\[\begin{split}Cov[x_i^2, x_j^2] &= Var[x_i] Var[x_j] + Cov(x_i, x_j)^2 \\ &+ 2 Cov(x_i, x_j) E[x_i] E[x_j] \\ &+ Var[x_i] E[x_j]^2 + Var[x_j] E[x_i]^2\end{split}\]
Parameters:
  • mean (ndarray) – 1D array mean values.

  • var (ndarray) – 1D array variance values.

  • cov (ndarray) – 3D array covariance. Shape is (nblock, ndata, ndata).

Returns:

new_cov – 3D array propagated covariance.

Return type:

ndarray

fft_gaussian_smooth(x, sigma_pix=20, mode='edge')[source]

My Gaussian smoother using FFTs. Input array is padded with edge values at the boundary by default. Pad size is 3*sigma_pix.

Parameters:
  • x (ndarray) – 1D array to smooth.

  • sigma_pix (float, default: 20) – Smoothing Gaussian sigma in terms of number of pixels.

  • mode (str) – Padding method. See numpy.pad() for options.

Returns:

y – Smoothed x values. Same size as x.

Return type:

ndarray

get_median_outlier_mask(flux, ivar, esigma=3.5, nwindow=11)[source]

Calculates median[i] of flux and sigma = 1.4826 * MAD[i] based on moving median filter of size nwindow. Returns a boolean array to mask pixels if abs(f[i] - median[i]) > esigma * max(n[i], sigma[i]).

Parameters:
  • flux (ndarray) – Flux array.

  • ivar (ndarray) – Inverse variance array.

  • esigma (float, default: 3.5) – Sigma to identify outliers via MAD.

  • nwindow (int, default: 11) – Window size of the moving median filter.

Returns:

mask – Bool array. Mask if true.

Return type:

ndarray

get_smooth_ivar(ivar, sigma_pix=20, esigma=3.5)[source]

Smoothing ivar values to reduce signal-noise coupling.

Smoothing is done on error=1/sqrt(ivar), while replacing ivar=0 and outliers in error values with the median. These replaced values are put back in in the final result.

Parameters:
  • ivar (ndarray) – Inverse variance array.

  • sigma_pix (float, default: 20) – Smoothing Gaussian sigma in terms of number of pixels.

  • esigma (float, default: 3.5) – Sigma to identify outliers via MAD.

Returns:

ivar2 – Smoothed ivar values. Outliers and masked values are put back in.

Return type:

ndarray

class FastLinear1DInterp(xp0, dxp, fp, copy=False, ep=None)[source]

Bases: object

Fast interpolator class for equally spaced data. Out of domain points are linearly extrapolated without producing any warnings or errors.

Uses _fast_eval_interp1d_lin().

Example:

one_interp = FastLinear1DInterp(0., 1., np.ones(3))
one_interp(5) # = 1
Parameters:
  • xp0 (float) – Initial x point for interpolation data.

  • dxp (float) – Spacing of x points.

  • fp (ndarray) – Function calculated at interpolation points.

  • copy (bool, default: False) – Copy input data, specifically fp.

  • ep (ndarray, optional) – Error on fp points. Not used! Bookkeeping purposes only.

reset(fp, copy=False, ep=None)[source]
class FastCubic1DInterp(xp0, dxp, fp, copy=False, ep=None, bc_type='not-a-knot')[source]

Bases: object

Fast cubic spline for equally spaced data. Out of domain points are extrapolated without producing any warnings or errors.

This class supports natural, not-a-knot boundary conditions. Uses _spline_cubic_natural(), _spline_cubic_notaknot(), and _fast_eval_interp1d_cubic().

Parameters:
  • xp0 (float) – Initial x point for interpolation data.

  • dxp (float) – Spacing of x points.

  • fp (ndarray) – Function calculated at interpolation points.

  • copy (bool, default: False) – Copy input data, specifically fp.

  • ep (ndarray, optional) – Error on fp points. Not used! Bookkeeping purposes only.

  • bc_type (str, default: 'not-a-knot') – Boundary condition type. Other option is ‘natural’. See scipy.interpolate.CubicSpline.

reset(fp, copy=False, ep=None)[source]
class SubsampleCov(ndata, nsamples, istart=0)[source]

Bases: object

Utility class to store all subsamples with weights and calculate the delete-one Jackknife covariance matrix.

Usage:

subsampler = SubsampleCov(ndata, nsamples)

for measurement, weights in SomeDataAfterFunction:
    # You can do this more than nsamples times
    subsampler.add_measurement(measurement, weights)

mean, covariance = subsampler.get_mean_n_cov()

Warning

You cannot call add_measurement() after calling get_mean(), get_mean_n_cov(), or get_mean_n_var().

Parameters:
  • ndata (int or tuple(int)) – Size or shape of the data vector. If tuple, it should be (nset, size1d). For example, 3 quantities share the same weights, data vector shape should pass ndata=(3, size1d).

  • nsamples (int) – Number of samples. You can add more measurements then this.

  • istart (int, default: 0) – Start index for the subsampling array

_istart

Sampler initial index.

Type:

int

_isample

Sample counter. Wraps around nsamples

Type:

int

_is_normalized

If the weights are normalized. Keeps track if _normalize() is called.

Type:

bool

all_measurements

3D array of shape (nsamples, nset, ndata).

Type:

ndarray

all_weights

3D array of shape (nsamples, nset, ndata).

Type:

ndarray

mean

Mean. 2D array of shape (nset, ndata)

Type:

ndarray or None

variance

Variance. 2D array of shape (nset, ndata)

Type:

ndarray or None

covariance

Covariance. 2D arrays of shape (ndata, ndata) or 3D arrays of shape (nblock, blockdim, blockdim).

Type:

list(ndarray) or None

add_measurement(xvec, wvec)[source]

Adds a measurement to the sample.

You can call this function more then nsamples times. The provided measurement should be weighted, but unnormalized, i.e. xvec=wi*xi. After a mean or covariance is obtained, you cannot add more measurements.

Parameters:
  • xvec (ndarray) – Data (measurement) vector.

  • wvec (ndarray) – 1D weight vector.

Raises:

RuntimeError – If the object is normalized by calling _normalize, get_mean, get_mean_n_cov and get_mean_n_var.

allreduce(comm, inplace)[source]

Sums statistics from all MPI process.

Note

Call this with inplace=MPI.IN_PLACE.

Parameters:
  • comm (MPI.COMM_WORLD) – MPI comm object for Allreduce

  • inplace (BufSpec) – MPI.IN_PLACE

_normalize()[source]
get_mean()[source]

Get the mean of all subsamples.

Warning

You cannot call add_measurement() after calling this unless you reset().

Returns:

mean – Mean.

Return type:

ndarray

_get_xdiff(mean_xvec, bias_correct=False)[source]
_get_block_covariance(x, blockdim)[source]
get_mean_n_cov(indices=None, blockdim=None, bias_correct=False)[source]

Get the mean and covariance of the mean using delete-one Jackknife.

Also sets mean and covariance.

Warning

You cannot call add_measurement() after calling this unless you reset().

Parameters:
  • indices (list(int), default: None) – Data set indices to estimate the covariance.

  • blockdim (int, default: None) – Calculate covariance by this block size instead of the full space.

  • bias_correct (bool, default: False) – Jackknife bias correction term for the mean.

Returns:

  • mean (ndarray) – Mean.

  • cov (list(ndarray)) – Covariances of the mean.

get_mean_n_var(bias_correct=False)[source]

Get the mean and variance of themean (i.e. diagonal of the covariance) using delete-one Jackknife.

Warning

You cannot call add_measurement() after calling this unless you reset().

Parameters:

bias_correct (bool, default: False) – Jackknife bias correction term for the mean.

Returns:

  • mean_xvec (ndarray) – Mean.

  • var_xvec (ndarray) – Variance of the mean. 1D array

reset()[source]