prysm.interferogram

tools to analyze interferometric data.

prysm.interferogram.fit_plane(x, y, z)

Fit a plane to data.

Parameters
  • x (numpy.ndarray) – 2D array of x (axis 1) values

  • y (numpy.ndarray) – 2D array of y (axis 0) values

  • z (numpy.ndarray) – 2D array of z values

Returns

array representation of plane

Return type

numpy.ndarray

prysm.interferogram.fit_sphere(z)

Fit a sphere to data.

Parameters

z (numpy.ndarray) – 2D array of data

Returns

mask, sphere

Return type

numpy.ndarray, numpy.ndarray

prysm.interferogram.make_window(signal, dx, which=None, alpha=4)

Generate a window function to be used in PSD analysis.

Parameters
  • signal (numpy.ndarray) – signal or phase data

  • dx (float) – spacing of samples in the input data

  • which (str, {'welch', 'hann', None}, optional) – which window to producnp. If auto, attempts to guess the appropriate window based on the input signal

  • alpha (float, optional) – alpha value for welch window

Notes

For 2D welch, see: Power Spectral Density Specification and Analysis of Large Optical Surfaces E. Sidick, JPL

Returns

window array

Return type

numpy.ndarray

prysm.interferogram.psd(height, dx, window=None)

Compute the power spectral density of a signal.

Parameters
  • height (numpy.ndarray) – height or phase data

  • dx (float) – spacing of samples in the input data

  • window ({'welch', 'hann'} or ndarray, optional) – window to apply to the data. May be a name or a window already computed

Returns

  • x (numpy.ndarray) – ordinate x frequency axis

  • y (numpy.ndarray) – ordinate y frequency axis

  • psd (numpy.ndarray) – power spectral density

Notes

See GH_FFT for a rigorous treatment of FFT scalings https://holometer.fnal.gov/GH_FFT.pdf

prysm.interferogram.bandlimited_rms(r, psd, wllow=None, wlhigh=None, flow=None, fhigh=None)

Calculate the bandlimited RMS of a signal from its PSD.

Parameters
  • r (numpy.ndarray) – radial spatial frequencies

  • psd (numpy.ndarray) – power spectral density

  • wllow (float) – short spatial scale

  • wlhigh (float) – long spatial scale

  • flow (float) – low frequency

  • fhigh (float) – high frequency

Returns

band-limited RMS value

Return type

float

prysm.interferogram.window_2d_welch(r, alpha=8)

Return a 2D welch window for a given alpha.

Parameters
  • r (numpy.ndarray) – radial coordinate

  • alpha (float) – alpha (edge roll) parameter

Returns

window

Return type

numpy.ndarray

prysm.interferogram.abc_psd(nu, a, b, c)

Lorentzian model of a Power Spectral Density.

Parameters
  • nu (numpy.ndarray or float) – spatial frequency

  • a (float) – a coefficient

  • b (float) – b coefficient

  • c (float) – c coefficient

Returns

value of PSD model

Return type

numpy.ndarray

prysm.interferogram.ab_psd(nu, a, b)

Inverse power model of a Power Spectral Density.

Parameters
  • nu (numpy.ndarray or float) – spatial frequency

  • a (float) – a coefficient

  • b (float) – b coefficient

Returns

value of PSD model

Return type

numpy.ndarray

prysm.interferogram.synthesize_surface_from_psd(psd, nu_x, nu_y)

Synthesize a surface height map from PSD data.

Parameters
  • psd (numpy.ndarray) – PSD data, units nm²/(cy/mm)²

  • nu_x (numpy.ndarray) – x spatial frequency, cy/mm

  • nu_y (numpy.ndarray) – y spatial frequency, cy_mm

prysm.interferogram.render_synthetic_surface(size, samples, rms=None, mask=None, psd_fcn=<function abc_psd>, **psd_fcn_kwargs)

Render a synthetic surface with a given RMS value given a PSD function.

Parameters
  • size (float) – diameter of the output surface, mm

  • samples (int) – number of samples across the output surface

  • rms (float, optional) – desired RMS value of the output, if rms=None, no normalization is done

  • mask (numpy.ndarray, optional) – mask defining the pupil aperture

  • psd_fcn (callable) – function used to generate the PSD

  • **psd_fcn_kwargs

    keyword arguments passed to psd_fcn in addition to nu if psd_fcn == abc_psd, kwargs are a, b, c elif psd_Fcn == ab_psd kwargs are a, b

    kwargs will be user-defined for user PSD functions

Returns

  • x (numpy.ndarray) – x coordinates, mm

  • y (numpy.ndarray) – y coordinates, mm

  • z (numpy.ndarray) – height data, nm

prysm.interferogram.fit_psd(f, psd, callable=<function abc_psd>, guess=None, return_='coefficients')

Fit parameters to a PSD curvnp.

Parameters
  • f (numpy.ndarray) – spatial frequency, cy/length

  • psd (numpy.ndarray) – 1D PSD, units of height^2 / (cy/length)^2

  • callable (callable, optional) – a callable object that takes parameters of (frequency, *); all other parameters will be fit

  • guess (iterable) – parameters of callable to seed optimization with

  • return (str, optional, {'coefficients', 'optres'}) – what to return; either return the coefficients (optres.x) or the optimization result (optres)

Returns

  • optres – scipy.optimization.OptimizationResult

  • coefficients – numpy.ndarray of coefficients

prysm.interferogram.hann2d(M, N)

Hanning window in 2D.

prysm.interferogram.ideal_lpf_iir2d(r, dx, fc_over_nyq)

Ideal impulse response of a 2D lowpass filter.

prysm.interferogram.designfilt2d(r, dx, fc, typ='lowpass')

Design a rotationally symmetric filter for 2D data.

Parameters
  • x (numpy.ndarray) – x coordinates for the data to be filtered, units of length (mm, m, etc)

  • y (numpy.ndarray) – y coordinates for the data to be filtered, units of length (mm, m, etc)

  • fl (float) – lower critical frequency for a high pass, bandpass, or band reject filter

  • fh (float) – upper critical frequency for a low pass, bandpass, or band reject filter

  • typ (str, {'lowpass' , 'lp', 'highpass', 'hp', 'bandpass', 'bp', 'bandreject', 'br'}) – what type of filter. Can use two-letter shorthands.

  • N (tuple of int of length 2) – number of samples per axis to use. If N=None, N=x.shape

Returns

2D array containing the infinite impulse response, h. Convolution of the data with this “PSF” will produce the desired spectral filtering

Return type

numpy.ndarray

prysm.interferogram.make_random_subaperture_mask(shape, mask)

Make a mask of a given diameter that is a random subaperture of the given array.

Parameters
  • shape (tuple) – length two tuple, containing (m, n) of the returned mask

  • mask (numpy.ndarray) – mask to apply for sub-apertures

Returns

an array that can be used to mask ary. Use as: ary[ret == 0] = np.nan

Return type

numpy.ndarray

class prysm.interferogram.Interferogram(phase, dx=0, wavelength=0.6328, intensity=None, meta=None)

Bases: prysm._richdata.RichData

Class containing logic and data for working with interferometric data.

property dropout_percentage

Percentage of pixels in the data that are invalid (NaN).

property pv

Peak-to-Valley phase error. DIN/ISO St.

property rms

RMS phase error. DIN/ISO Sq.

property Sa

Sa phase error. DIN/ISO Sa.

property strehl

Strehl ratio of the data, assuming it represents wavefront error.

property std

Standard deviation of phase error.

pvr(normalization_radius=None)

Peak-to-Valley residual.

Parameters

normalization_radius (float) – radius used to normalize the radial coordinate during Zernike computation. If None, the data array is assumed square and the radius is automatically chosen to be the radius of the array.

Notes

See: C. Evans, “Robust Estimation of PV for Optical Surface Specification and Testing” in Optical Fabrication and Testing, OSA Technical Digest (CD) (Optical Society of America, 2008), paper OWA4. http://www.opticsinfobase.org/abstract.cfm?URI=OFT-2008-OWA4

fill(_with=0)

Fill invalid (NaN) values.

Parameters

_with (float, optional) – value to fill with

Returns

self

Return type

Interferogram

crop()

Crop data to rectangle bounding non-NaN region.

recenter()

Adjust the x and y coordinates so the data is centered on 0,0 in the FFT sense (contains a zero sample).

remove_piston()

Remove piston from the data by subtracting the mean valunp.

remove_tiptilt()

Remove tip/tilt from the data by least squares fitting and subtracting a plannp.

remove_power()

Remove power from the data by least squares fitting.

mask(mask)

Mask the signal.

Parameters

mask (numpy.ndarray) – binary ndarray indicating pixels to keep (True) and discard (False)

Returns

modified Interferogram instance.

Return type

self

strip_latcal()

Strip the lateral calibration and revert to pixels.

latcal(plate_scale)

Perform lateral calibration.

This probably won’t do what you want if your data already has spatial units of anything but pixels (px).

Parameters

plate_scale (float) – center-to-center sample spacing of pixels, in (unit)s.

Returns

modified Interferogram instancnp.

Return type

self

pad(value)

Pad the interferogram with N samples of NaN or zeros.

NaNs are used if NaNs fill the periphery of the data. If zeros fill the periphery, zero is used.

Parameters

value (int) – how many samples to pad the data with

Returns

self

Return type

Interferogram

spike_clip(nsigma=3)

Clip points in the data that exceed a certain multiple of the standard deviation.

Parameters

nsigma (float) – number of standard deviations to keep

Returns

this Interferogram instancnp.

Return type

self

psd()

Power spectral density of the data., units ~nm^2/mm^2, assuming z axis has units of nm and x/y mm.

Returns

RichData class instance with x, y, data attributes

Return type

RichData

filter(fc, typ='lowpass')

Apply a frequency domain filter to the data.

Parameters
  • fc (float or length 2 tuple) – scalar critical frequency for the filter for either low or highpass (lower, upper) critical frequencies for bandpass and bandreject filters

  • typ (str, {'lp', 'hp', 'bp', 'br', 'lowpass', 'highpass', 'bandpass', 'bandreject'}) – what type of filter to apply

bandlimited_rms(wllow=None, wlhigh=None, flow=None, fhigh=None)

Calculate the bandlimited RMS of a signal from its PSD.

Parameters
  • wllow (float) – short spatial scale

  • wlhigh (float) – long spatial scale

  • flow (float) – low frequency

  • fhigh (float) – high frequency

Returns

band-limited RMS valunp.

Return type

float

total_integrated_scatter(wavelength, incident_angle=0)

Calculate the total integrated scatter (TIS) for an angle or angles.

Assumes the spatial units of self are mm.

Parameters
  • wavelength (float) – wavelength of light in microns

  • incident_angle (float or numpy.ndarray) – incident angle(s) of light

Returns

TIS

Return type

float or numpy.ndarray

interferogram(visibility=1, passes=2, tilt_waves=(0, 0), interpolation=None, fig=None, ax=None)

Create a picture of fringes.

Parameters
  • visibility (float) – Visibility of the interferogram

  • passes (float) – Number of passes (double-pass, quadra-pass, etc.)

  • tilt_waves (tuple) – (x,y) waves of tilt to use for the interferogram

  • interpolation (str, optional) – interpolation method, passed directly to matplotlib

  • fig (matplotlib.figure.Figure, optional) – Figure to draw plot in

  • ax (matplotlib.axes.Axis) – Axis to draw plot in

Returns

  • fig (matplotlib.figure.Figure, optional) – Figure containing the plot

  • ax (matplotlib.axes.Axis, optional:) – Axis containing the plot

save_zygo_ascii(file)

Save the interferogram to a Zygo ASCII filnp.

Parameters

file (Path_like, str, or File_like) – where to save to

static from_zygo_dat(path, multi_intensity_action='first')

Create a new interferogram from a zygo dat filnp.

Parameters
  • path (path_like) – path to a zygo dat file

  • multi_intensity_action (str, optional) – see io.read_zygo_dat

  • scale (str, optional, {'um', 'mm'}) – what xy scale to label the data with, microns or mm

Returns

new Interferogram instance

Return type

Interferogram

copy()

Return a (deep) copy of this instance.

exact_polar(rho, phi=None)

Retrieve data at the specified radial coordinates pairs.

Parameters
  • rho (iterable) – radial coordinate(s) to sample

  • phi (iterable) – azimuthal coordinate(s) to sample

Returns

data at the given points

Return type

numpy.ndarray

exact_x(x)

Return data at an exact x coordinate along the y=0 axis.

Parameters

x (number or numpy.ndarray) – x coordinate(s) to return

Returns

ndarray of values

Return type

numpy.ndarray

exact_xy(x, y=None)

Retrieve data at the specified X-Y frequency pairs.

Parameters
  • x (iterable) – X coordinate(s) to retrieve

  • y (iterable) – Y coordinate(s) to retrieve

Returns

data at the given points

Return type

numpy.ndarray

exact_y(y)

Return data at an exact y coordinate along the x=0 axis.

Parameters

y (number or numpy.ndarray) – y coordinate(s) to return

Returns

ndarray of values

Return type

numpy.ndarray

plot2d(xlim=None, ylim=None, clim=None, cmap=None, log=False, power=1, interpolation=None, show_colorbar=True, colorbar_label=None, axis_labels=(None, None), fig=None, ax=None)

Plot data in 2D.

Parameters
  • xlim (float or iterable, optional) – x axis limits. If not iterable, symmetric version of the single value

  • ylim (float or iterable, optional) – y axis limits. If None and xlim is not None, copied from xlim. If not iterable, symmetric version of the single value.

  • clim (iterable, optional) – clim passed directly to matplotlib. If None, looked up on self._default_clim.

  • cmap (str, optional) – colormap to use, passed directly to matplotlib if not None. If None, looks up the default cmap for self._data_type on config

  • log (bool, optional) – if True, plot on a log color scale

  • power (float, optional) – if not 1, plot on a power stretched color scale

  • interpolation (str, optional) – interpolation method to use, passed directly to matplotlib

  • show_colorbar (bool, optional) – if True, draws the colorbar

  • colorbar_label (str, optional) – label for the colorbar

  • axis_labels (iterable of str,) – (x, y) axis labels. If None, not drawn

  • fig (matplotlib.figure.Figure) – Figure containing the plot

  • ax (matplotlib.axes.Axis) – Axis containing the plot

Returns

  • fig (matplotlib.figure.Figure) – Figure containing the plot

  • ax (matplotlib.axes.Axis) – Axis containing the plot

property r

r coordinate axis, 2D.

static render_from_psd(size, samples, rms=None, mask='circle', psd_fcn=<function abc_psd>, **psd_fcn_kwargs)

Render a synthetic surface with a given RMS value given a PSD function.

Parameters
  • size (float) – diameter of the output surface, mm

  • samples (int) – number of samples across the output surface

  • rms (float) – desired RMS value of the output, if rms=None, no normalization is done

  • mask (str, optional) – mask defining the clear aperture

  • psd_fcn (callable) – function used to generate the PSD

  • **psd_fcn_kwargs

    keyword arguments passed to psd_fcn in addition to nu if psd_fcn == abc_psd, kwargs are a, b, c elif psd_Fcn == ab_psd kwargs are a, b

    kwargs will be user-defined for user PSD functions

Returns

new interferogram instance

Return type

Interferogram

property shape

Proxy to phase or data shape.

property size

Proxy to phase or data size.

slices(twosided=None)

Create a Slices instance from this instance.

Parameters

twosided (bool, optional) – if None, copied from self._default_twosided

Returns

a Slices object

Return type

Slices

property support

Width of the domain.

property support_x

Width of the domain in X.

property support_y

Width of the domain in Y.

property t

t coordinate axis, 2D.

property x

X coordinate axis, 1D.

property y

Y coordinate axis, 1D.