SinglePulse Class

The SinglePulse class handles an individual pulse data array.

class SinglePulse(data[, mpw=None, ipw=None, opw=None, prepare=False, align=None, period=None, windowsize=None, istemplate=False])
Parameters

  • data (list/numpy.ndarray) – Data array.

  • mpw (list/numpy.ndarray) – Indices of the main-pulse window

  • ipw (list/numpy.ndarray) – Indices of the inter-pulse window

  • opw (list/numpy.ndarray) – Indices of the off-pulse window

  • prepare (bool) – Run interpulse_align()

  • align (float) – Rotate the pulse by this amount in phase bin units.

  • period (float) – Provide pulse period.

  • windowsize (int) – Overrides mpw, ipw, opw. Defines an off-pulse window of a given size about the minimum-integral region.

  • istemplate (bool) – If true, run calcFT() to precalculate the Fourier Transform of the data. This is useful for speed-ups elsewhere throughout the code.

Usage:

sp = SinglePulse(data,windowsize=256) #will auto-calculate an offpulse region of length 256 bins
print(sp.getFWHM()) #prints the FWHM of the pulse
print(sp.getSN()) #prints a crude S/N of the pulse
print(sp.fitPulse(template_array)[5]) #prints a better S/N of the pulse using a template array

Methods

interpulse_align()

Rotate the pulse such that the main pulse is at phase=0.25 and the interpulse is at phase = 0.75.

Returns

None

center_align()

Rotate the pulse such that the peak is in the center.

Returns

None

normalize([area=False])

Normalize the pulse so that the peak has a value of 1.

Parameters

area (bool) – Normalize the pulse so that the area has a value of 1.

Returns

None

getFWHM([simple=False, timeunits=True])

Get the full width at half maximum of the main component of the pulse.

Returns

FWHM, float

getWeff([sumonly=False, timeunits=True])

Calculate the effective width of the pulse.

Returns

W_eff, float

getSN()

Calculate a crude signal-to-noise ratio as the maximum of the data array divided by the off-pulse RMS noise.

Returns

S/N, float

remove_baseline([save=True])

Subtract the baseline of the pulse so that it will be set to zero mean.

Parameters

save (bool) – Replace the data array with the rotated versions.

Returns

subtracted_array, numpy.ndarray

calcOffpulseWindow([windowsize=None])

Calculate an offpulse window automatically by finding where the area under the curve is minimized, i.e., where

\[\min \sum_{i=i_0}^{i_0+N} U(\phi_i)\]

is satisfied. The template (or profile) shape as a function of phase is \(U(\phi)\), \(i_0\) is the initial index and \(N\) is the length of the window. Returns the offpulse window indices.

Parameters

windowsize (int) – Window size (\(N\)) in bin units. The default is the length of the data array divided by 8.

Returns

opw, numpy.ndarray

calcFT()

Calculate the Fourier Transform and related frequencies Returns the Fourier frequencies, Fourier Transform Uses a real FFT

Returns

fs, np.ndarray, yfft, np.ndarray

getMainpulse()

Return the main-pulse intensities.

Returns

mpw, numpy.ndarray

getInterpulse()

Return the inter-pulse intensities.

Returns

mpw, numpy.ndarray

getOffpulse()

Return the off-pulse intensities.

Returns

mpw, numpy.ndarray

getAllpulse()

Return the main-pulse, inter-pulse, and off-pulse intensities.

Returns

mpw, numpy.ndarray, ipw, numpy.ndarray, opw, numpy.ndarray

getMainpulseACF()

Return the autocorrelation function of the main-pulse.

Returns

mpacf, nump.array

getInterpulseACF()

Return the autocorrelation function of the inter-pulse.

Returns

ipacf, nump.array

getOffpulseACF()

Return the autocorrelation function of the off-pulse.

Returns

opacf, nump.array

getAllpulseACF()

Return the main-pulse, inter-pulse, and off-pulse autocorrelation functions.

Returns

mpacf, numpy.ndarray, ipacf, numpy.ndarray, opacf, numpy.ndarray

getOffpulseNoise([mean=False, full=False])

Calculate the root-mean-square (RMS) of the off-pulse intensities.

Parameters

  • mean (bool) – Calculate the mean value of the off-pulse intensities instead of the root-mean-square.

  • full (bool) – Return both the off-pulse mean and RMS.

Returns

noise, float (if full is True then return two floats, the mean and the RMS)

getOffpulseZCT()

Perform a zero-crossing test of the offpulse noise. The expected number of zero crossings for \(N\) samples of white noise is \((N-1)/2 \pm \sqrt{(N-1)/2}\). If the number is within \(3\sigma\) of the expectation, return True, otherwise, return False.

Returns

success, test_statistic, count

fitPulse(template[, fixedphase=False, rms_baseline=None])

Perform the template-matching procedure of Taylor 1992 in the Fourier domain. Matched filtering assumes that the data profile \(I(t)\) is a scaled and shifted version of a template profile \(U(t)\), plus additive noise, i.e.,

\[I(t) = b U(t - \tau) + n(t)\]

where \(t\) is time (in pulse phase bin units), \(n(t)\) is the additive noise, and \(b\) and \(\tau\) are the scale factor and shift, respectively. PyPulse performs this fit in terms of phase bins but these numbers can be translated into absolute arrival times.

Parameters

  • template (list/numpy.ndarray/SinglePulse) – Template shape to fit to the data array.

  • fixedphase (bool) – Return only the signal-to-noise with the template shape fixed in its phase.

  • rms_baseline (float) – Provide a different value of the off-pulse root-mean-square noise.

Returns

TOA from cross-correlation function (tauccf), TOA from template matching procedure (tauhat), scale factor (bhat), error on TOA (sigma_Tau), error on scale factor (sigma_b), signal-to-noise ratio (snr), cross-correlation coefficient (rho)

Todo

include a goodness-of-fit flag as a measure of the residuals

shiftit(shift[, save=False])

Shift the pulse by some number of phase bins

Parameters

  • shift (float) – Rotate the pulse in bin units.

  • save (bool) – Replace the data array with the rotated versions.

spline_smoothing([lam=None, **kwargs])

Iterative cubic spline interpolation for smooth template generation. The spline-generation algorithm is defined by D.S.G. Pollock 1999 (Smoothing with Cubic Splines), in which the cubic spline preserves slope and curvature through each of the control points, i.e., it minimizes the function:

\[L = \lambda \chi^2 + (1-\lambda) \int \left[S^{''}(\phi)\right]^2 d\phi\]

where \(\chi^2\) is the sum of the residuals (“chi-squared”), the integral is the curvature of the spline over all of pulse phase \(\phi\), and \(\lambda\) represents the paramter (“lam”) that weights the importance of the two components.

Control points (knots) are generated iteratively. Polynomials up to quartic order are fit over the data array. If a quadratic well fits the data, then return without adding a knot, else subdivide the data in half and return the midpoint as a knot. The data array is continuously divided into two halves up to some minimum separation (default=16). In locations \(\phi_i\) where a knot is decided to be placed, the phase bins closest to that knot location are fit by a cubic \(f(\phi\) and the knot is placed at \((\phi_i,f(\phi))\). After the initial knots locations are chosen, run the cubic spline interpolation as described before. If there is a location where a new knot will help the fit (the residuals must be greater than 4 times the RMS and the location in phase must be 3 times the off-pulse RMS, i.e., the current spline model fails there and it must be a location where actual pulse intensity is present), include that knot and re-run the cubic spline interpolation until the residuals RMS is low.

Parameters

  • lam (float) – The lambda parameter that determines the relative weight of fitting through the knots and minimizing the curvature of the spline. Can take values from (0,1], where 1 = fit exactly through the knots and 0 = minimized curvature of the spline. If not provided, it is calculated as 1-(off-pulse rms)^2.

  • **kwargs – Additional arguments to pass to utils.subdivide().

Returns

template, numpy.ndarray

component_fitting([mode='gaussian', nmax=10, full=False, minamp=None, alpha=0.05, allownegative=False, verbose=False, save=False])

Iteratively fits a total of nmax components to a pulse profile. The iteration will end when the \(\chi_r^2\) is minimized.

If mode == ‘gaussian’, will fit components of the form

\[f(\phi|A, \mu, \sigma) = A \exp\left[-\frac{1}{2}\left(\frac{\phi-\mu}{\sigma}\right)^2\right]\]

where \(\phi\) is the pulse phase, \(A\) is the component amplitude, \(\mu\) is the mean pulse phase, and \(\sigma\) is the component width.

If mode == ‘vonmises’, will fit components of the form

\[f(\phi|A, \mu, \kappa) = A \frac{\exp\left[\kappa \cos(\phi-\mu)\right]}{2\pi I_0(\kappa)}\]

where again \(\phi\) is the pulse phase and \(A\) is the component amplitude, \(\mu\) is the mean pulse phase, and \(1/\kappa\) is the analogous to the component width (squared) of a Gaussian, \(\sigma^2\). The function \(I_0\) is the modified Bessel function of order 0.

Parameters

  • mode (str) – The function to fit, can choose: ‘gaussian’, ‘vonmises’.

  • nmax (int) – Maximum number of components to test.

  • full (bool) – Return the template shape as well as the array of the best-fit parameters and the number of parameters.

  • minamp (float) – If all of the profile residuals are less than the fraction of the maximum set by minamp, then return. For example, if minamp is 0.05, once all profile residuals are less than 5% of the maximum, return.

  • alpha (float) – F-test significance value. For example, if alpha is 0.05, then new components are kept if 1 - 0.05 = 95% significant.

  • allownegative (bool) – If true, allow negative-amplitude components in the fitting

  • verbose (bool) – Print when each component is fit

  • plot (bool) – If true, show a plot of the fit, with the components

  • filename (str) – If provided, save the plot above to the filename.

  • save (bool) – If true, replace the data array with the smoothed profile.

Returns

template, numpy.ndarray (if full is True then the template, the fit of the parameters, and the number of parameters)

gaussian_smoothing([nmax=10])

Runs component_fitting() with mode = ‘gaussian’

Parameters

nmax (int) – Maximum number of components to test.

Returns

template, numpy.ndarray

vonmises_smoothing([nmax=10])

Runs component_fitting() with mode = ‘vonmises’

Parameters

nmax (int) – Maximum number of components to test.

Returns

template, numpy.ndarray

savgol_smoothing([window_length=11, polyorder=5, save=False])

Runs scipy’s signal.savgol_filter()

Parameters

  • window_length (int) – Averaging length to fit the polynomial of order polyorder

  • polyorder (int) – Polynomial order to fit local data segments.

  • save (bool) – If true, replace the data array with the smoothed profile.

Returns

template, numpy.ndarray

getPeriod()

Return the pulse period.

Returns

period, float

getNbin()

Return the number of phase bins.

Returns

nbins, int

getTbin()

Return the time per phase bin

Returns

t, float

plot([show=True])

Simple plot of the pulse.

Parameters

show (bool) – Show the plot.

plot_windows([show=True])

Diagnostic plot of the main-, inter-, and off-pulse regions.

Parameters

show (bool) – Show the plot.