DeltaVariance

class turbustat.statistics.DeltaVariance(img, header=None, weights=None, diam_ratio=1.5, lags=None, nlags=25, distance=None)[source]

Bases: BaseStatisticMixIn

The delta-variance technique as described in Ossenkopf et al. (2008).

Parameters:
imgnumpy.ndarray or astropy.io.fits.PrimaryHDU or astropy.io.fits.ImageHDU or spectral_cube.Projection or spectral_cube.Slice

The image calculate the delta-variance of.

headerFITS header, optional

Image header. Required when img is a ndarray.

weightsnumpy.ndarray or astropy.io.fits.PrimaryHDU or astropy.io.fits.ImageHDU or spectral_cube.Projection or spectral_cube.Slice

Weights to be used.

diam_ratiofloat, optional

The ratio between the kernel sizes.

lagsnumpy.ndarray or list, optional

The pixel scales to compute the delta-variance at.

nlagsint, optional

Number of lags to use.

distanceQuantity, optional

Physical distance to the region in the data.

Examples

>>> from turbustat.statistics import DeltaVariance
>>> from astropy.io import fits
>>> moment0 = fits.open("2D.fits") 
>>> delvar = DeltaVariance(moment0) 
>>> delvar.run(verbose=True) 

Attributes Summary

brk

Fitted break point.

brk_err

1-sigma on the break point in the segmented linear model.

convolve_arrays

convolve_weights

data

delta_var

Delta Variance values.

delta_var_error

1-sigma errors on the Delta variance values.

distance

fit_range

Range of lags used in the fit.

header

lags

Lag values.

need_header_flag

no_data_flag

slope

Fitted slope.

slope_err

Standard error on the fitted slope.

weights

Array of weights.

Methods Summary

compute_deltavar([allow_huge, boundary, ...])

Perform the convolution and calculate the delta variance at all lags.

fit_plaw([xlow, xhigh, brk, verbose, ...])

Fit a power-law to the Delta-variance spectrum.

fitted_model(xvals)

Computes the fitted power-law in log-log space using the given x values.

input_data_header(data, header[, need_copy])

Check if the header is given separately from the data type.

load_beam([beam])

Try loading the beam from the header or a given object.

load_results(pickle_file)

Load in a saved pickle file.

plot_fit([save_name, xunit, symbol, color, ...])

Plot the delta-variance curve and the fit.

run([show_progress, verbose, xunit, ...])

Compute the delta-variance.

save_results(output_name[, keep_data])

Save the results of the SCF to avoid re-computing.

Attributes Documentation

brk

Fitted break point.

brk_err

1-sigma on the break point in the segmented linear model.

convolve_arrays
convolve_weights
data
delta_var

Delta Variance values.

delta_var_error

1-sigma errors on the Delta variance values.

distance
fit_range

Range of lags used in the fit.

header
lags

Lag values.

need_header_flag = True
no_data_flag = False
slope

Fitted slope.

slope_err

Standard error on the fitted slope.

weights

Array of weights.

Methods Documentation

compute_deltavar(allow_huge=False, boundary='wrap', min_weight_frac=0.01, nan_treatment='fill', preserve_nan=False, use_pyfftw=False, threads=1, pyfftw_kwargs={}, show_progress=True, keep_convolve_arrays=False)[source]

Perform the convolution and calculate the delta variance at all lags.

Parameters:
allow_hugebool, optional

Passed to convolve_fft. Allows operations on images larger than 1 Gb.

boundary{“wrap”, “fill”}, optional

Use “wrap” for periodic boundaries, and “fill” for non-periodic.

min_weight_fracfloat, optional

Set the fraction of the peak of the weight array to mask below. Default is 0.01. This will remove most edge artifacts, but is not guaranteed to! Increase this value if artifacts are encountered (this typically results in large spikes in the delta-variance curve).

nan_treatmentbool, optional

Enable to interpolate over NaNs in the convolution. Default is True.

use_pyfftwbool, optional

Enable to use pyfftw, if it is installed.

threadsint, optional

Number of threads to use in FFT when using pyfftw.

pyfftw_kwargsPassed to

See here for a list of accepted kwargs.

show_progressbool, optional

Show a progress bar while convolving the image at each lag.

keep_convolve_arraysbool, optional

Keeps the convolved arrays at each lag. Disabled by default to minimize memory usage.

fit_plaw(xlow=None, xhigh=None, brk=None, verbose=False, bootstrap=False, bootstrap_kwargs={}, **fit_kwargs)[source]

Fit a power-law to the Delta-variance spectrum.

Parameters:
xlowQuantity, optional

Lower lag value to consider in the fit.

xhighQuantity, optional

Upper lag value to consider in the fit.

brkQuantity, optional

Give an initial guess for a break point. This enables fitting with a turbustat.statistics.Lm_Seg.

bootstrapbool, optional

Bootstrap using the model residuals to estimate the standard errors.

bootstrap_kwargsdict, optional

Pass keyword arguments to residual_bootstrap.

verbosebool, optional

Show fit summary when enabled.

fitted_model(xvals)[source]

Computes the fitted power-law in log-log space using the given x values.

Parameters:
xvalsndarray

Values of log(lags) to compute the model at (base 10 log).

Returns:
model_valuesndarray

Values of the model at the given values.

input_data_header(data, header, need_copy=False)

Check if the header is given separately from the data type.

load_beam(beam=None)

Try loading the beam from the header or a given object.

Parameters:
beamBeam, optional

The beam.

static load_results(pickle_file)

Load in a saved pickle file.

Parameters:
pickle_filestr

Name of filename to load in.

Returns:
selfSave statistic class

Statistic instance with saved results.

Examples

Load saved results. >>> stat = Statistic.load_results(“stat_saved.pkl”) # doctest: +SKIP

plot_fit(save_name=None, xunit=Unit('pix'), symbol='o', color='r', fit_color='k', label=None, show_residual=True)[source]

Plot the delta-variance curve and the fit.

Parameters:
save_namestr,optional

Save the figure when a file name is given.

xunitu.Unit, optional

The unit to show the x-axis in.

symbolstr, optional

Shape to plot the data points with.

color{str, RGB tuple}, optional

Color to show the delta-variance curve in.

fit_color{str, RGB tuple}, optional

Color of the fitted line. Defaults to color when no input is given.

labelstr, optional

Label to later be used in a legend.

show_residualbool, optional

Plot the fit residuals.

run(show_progress=True, verbose=False, xunit=Unit('pix'), nan_treatment='fill', preserve_nan=False, allow_huge=False, boundary='wrap', use_pyfftw=False, threads=1, pyfftw_kwargs={}, xlow=None, xhigh=None, brk=None, fit_kwargs={}, save_name=None)[source]

Compute the delta-variance.

Parameters:
show_progressbool, optional

Show a progress bar during the creation of the covariance matrix.

verbosebool, optional

Plot delta-variance transform.

xunitu.Unit, optional

The unit to show the x-axis in.

allow_hugebool, optional

See do_convolutions.

nan_treatmentbool, optional

Enable to interpolate over NaNs in the convolution. Default is True.

boundary{“wrap”, “fill”}, optional

Use “wrap” for periodic boundaries, and “cut” for non-periodic.

use_pyfftwbool, optional

Enable to use pyfftw, if it is installed.

threadsint, optional

Number of threads to use in FFT when using pyfftw.

pyfftw_kwargsPassed to

See here for a list of accepted kwargs.

xlowQuantity, optional

Lower lag value to consider in the fit.

xhighQuantity, optional

Upper lag value to consider in the fit.

brkQuantity, optional

Give an initial break point guess. Enables fitting a segmented linear model.

fit_kwargsdict, optional

Passed to fit_model when using a broken linear fit.

save_namestr,optional

Save the figure when a file name is given.

save_results(output_name, keep_data=False)

Save the results of the SCF to avoid re-computing. The pickled file will not include the data cube by default.

Parameters:
output_namestr

Name of the outputted pickle file.

keep_databool, optional

Save the data cube in the pickle file when enabled.