pyvib.filter module¶
-
pyvib.filter.differentiate(y, fs, order=3, cutoff=0.5, isnumeric=False)[source]¶ Differentiate y twice to get vel and acc
5-point stencil offers fairly good results in most cases [1]_
Differentiate signals should normally be avoided. Only use when distance sensors cannot be avoided, ie, rotating machinery.
In frequency domain, differentiation is given by Y(omega) = -omega**2*Y(omega) which shows that differentiation amplifies high-frequency noise
The normalized cutoff can be calculated as cutoff = 100 Hz nyq = 0.5 * fs normal_cutoff = cutoff / nyq
- cutoff : scalar (0-1)
- The freq where the gain drops to 1/sqrt(2) that of the passband (the “-3 dB point”). Normalized from 0 to 1, where 1 is the Nyquist frequency
Notes
. [1] K. Worden: Nonlinearity in structural dynamics. Appendix I
-
pyvib.filter.integrate(ddy, fs, lowcut=None, highcut=None, order=3, isnumeric=False)[source]¶ Integrate acceleration to get velocity and displacement.
Numerical integration is prone to low-frequency problems. Trapez-rule only have problems with low-freq, but is imprecise from 0.1*nyquist freq. Simpsons and Tick’s rule or give better precision, but have problems with high frequencies as well as low. They can(and should) be used until 1/4 of nyquist freq. See fig. I.9 p. 618 [1]_.
Instead of using band-pass filter, using low pass filter for measurement (acceleration) to remove noise and then high pass filter for integrated signals (velocity and position) for removing drift, yields better integrated signals.
Only filter ẏ (dy) after y have been calculated: In reality ẏ is zero-mean, but a finite sampling will cause it to be non-zero mean. If we now: ∫ ẏ(t) - ẏ_mean dt = y - ẏ_mean*t + k1 Ie. filtering ẏ before integration, introduce a linear trend in y (ẏ_mean*t) that then has to be removed somehow. k1 is the integration constant. We choose y(0) = 0, ie k1 = 0.
The cutoff frequencies are normalized with the Nyquist frequency (ie. 1 is the nyquist frequency). After filtering, they are the frequencies where the gain drops to 1/sqrt(2) of that of the passband (the “-3 dB point”).
- lowcut : float
- Cutoff frequency in Hz for the highpass filtering of integrated signals.
- Highcut : float
- Cutoff frequency in Hz for the lowpass filtering of acceleration signal before integrating.
- isnumeric : bool
- Only filter ‘real data’. Simulated data without noise should not be filtered
Notes
. [1] K. Worden: Nonlinearity in structural dynamics. Appendix I
-
pyvib.filter.resample(y, fs_in, fs_out, cutoff, order=3, isnumeric=False)[source]¶ Downsampling can be done like this:
- cutoff : int
- The cutoff frequency for lowpass filtering. Should be lower than the output nyquist frequency, ie. could be 98% of the output nyq
- y_out : ndarrray
- The resampled signal at fs_out Hz
In the simple case where your array’s size is divisible by the downsampling factor (R), you can reshape your array, and take the mean along the new axis: a = np.array([1.,2,6,2,1,7]) R = 3 a.reshape((-1, R)) => array([[ 1., 2., 6.],
[ 2., 1., 7.]])a.reshape((-1, R)).mean(axis=1) => array([ 3. , 3.33333333])
In the general case, you can pad your array with NaNs to a size divisible by R, and take the mean using scipy.nanmean: import math, scipy b = np.append(a, [ 4 ]) b.shape => (7,) pad_size = math.ceil(float(b.size)/R)*R - b.size b_padded = np.append(b, np.zeros(pad_size)*np.NaN) b_padded.shape => (9,) scipy.nanmean(b_padded.reshape(-1,R), axis=1) => array([ 3. , 3.33333333, 4.])