Implemented statistics functions
Low-level access to the statistics functions
The deviation functions are generally of the form:
(tau_out, adev, adeverr, n) = allantools.adev(data, rate=1.0, data_type="phase", taus=None)
Inputs:
data = list of phase measurements in seconds, or list of fractional frequency measurements (nondimensional)
rate = sample rate of data in Hz , i.e. interval between phase measurements is 1/rate seconds.
data_type= = either “phase” or “freq”
taus = list of tau-values for ADEV computation. The keywords “all”, “octave”, or “decade” can also be used.
Outputs
tau_out = list of tau-values for which deviations were computed
adev = list of ADEV (or another statistic) deviations
adeverr = list of estimated errors of allan deviations. some functions instead return a confidence interval (err_l, err_h)
n = list of number of pairs in allan computation. standard error is adeverr = adev/sqrt(n)
Statistics
- allantools.adev(data, rate=1.0, data_type='phase', taus=None)
- Allan deviation.
Classic - use only if required - relatively poor confidence.
\[\sigma^2_{ADEV}(\tau) = { 1 \over 2 \tau^2 } \langle ( {x}_{n+2} - 2x_{n+1} + x_{n} )^2 \rangle = { 1 \over 2 (N-2) \tau^2 } \sum_{n=1}^{N-2} ( {x}_{n+2} - 2x_{n+1} + x_{n} )^2\]where \(x_n\) is the time-series of phase observations, spaced by the measurement interval \(\tau\), and with length \(N\).
Or alternatively calculated from a time-series of fractional frequency:
\[\sigma^{2}_{ADEV}(\tau) = { 1 \over 2 } \langle ( \bar{y}_{n+1} - \bar{y}_n )^2 \rangle\]where \(\bar{y}_n\) is the time-series of fractional frequency at averaging time \(\tau\)
Parameters
- data: np.array
Input data. Provide either phase or frequency (fractional, adimensional).
- rate: float
The sampling rate for data, in Hz. Defaults to 1.0
- data_type: {‘phase’, ‘freq’}
Data type, i.e. phase or frequency. Defaults to “phase”.
- taus: np.array
Array of tau values, in seconds, for which to compute statistic. Optionally set taus=[“all”|”octave”|”decade”] for automatic tau-list generation.
Returns
- (taus2, ad, ade, ns): tuple
Tuple of values
- taus2: np.array
Tau values for which td computed
- ad: np.array
Computed adev for each tau value
- ade: np.array
adev errors
- ns: np.array
Values of N used in each adev calculation
References
NIST [SP1065] eqn (6) and (7), pages 14 and 15.
- allantools.oadev(data, rate=1.0, data_type='phase', taus=None)
Overlapping Allan deviation. General purpose - most widely used - first choice.
\[\sigma^2_{OADEV}(m\tau_0) = { 1 \over 2 (m \tau_0 )^2 (N-2m) } \sum_{n=1}^{N-2m} ( {x}_{n+2m} - 2x_{n+1m} + x_{n} )^2\]where \(\sigma_{OADEV}(m\tau_0)\) is the overlapping Allan deviation at an averaging time of \(\tau=m\tau_0\), and \(x_n\) is the time-series of phase observations, spaced by the measurement interval \(\tau_0\), with length \(N\).
Parameters
- data: np.array
Input data. Provide either phase or frequency (fractional, adimensional).
- rate: float
The sampling rate for data, in Hz. Defaults to 1.0
- data_type: {‘phase’, ‘freq’}
Data type, i.e. phase or frequency. Defaults to “phase”.
- taus: np.array
Array of tau values, in seconds, for which to compute statistic. Optionally set taus=[“all”|”octave”|”decade”] for automatic tau-list generation.
Returns
- (taus2, ad, ade, ns): tuple
Tuple of values
- taus2: np.array
Tau values for which td computed
- ad: np.array
Computed oadev for each tau value
- ade: np.array
oadev errors
- ns: np.array
Values of N used in each oadev calculation
References
NIST [SP1065] eqn (11), page 16.
- allantools.mdev(data, rate=1.0, data_type='phase', taus=None)
Modified Allan deviation. Used to distinguish between White and Flicker Phase Modulation.
\[\sigma^2_{MDEV}(m\tau_0) = { 1 \over 2 (m \tau_0 )^2 (N-3m+1) } \sum_{j=1}^{N-3m+1} \left[ \sum_{i=j}^{j+m-1} {x}_{i+2m} - 2x_{i+m} + x_{i} \right]^2\]Parameters
- data: np.array
Input data. Provide either phase or frequency (fractional, adimensional).
- rate: float
The sampling rate for data, in Hz. Defaults to 1.0
- data_type: {‘phase’, ‘freq’}
Data type, i.e. phase or frequency. Defaults to “phase”.
- taus: np.array
Array of tau values, in seconds, for which to compute statistic. Optionally set taus=[“all”|”octave”|”decade”] for automatic tau-list generation.
Returns
- (taus2, md, mde, ns): tuple
Tuple of values
- taus2: np.array
Tau values for which td computed
- md: np.array
Computed mdev for each tau value
- mde: np.array
mdev errors
- ns: np.array
Values of N used in each mdev calculation
References
NIST [SP1065] eqn (14), page 17.
- allantools.pdev(data, rate=1.0, data_type='phase', taus=None)
Parabolic deviation.
Use for evaluating uncertainty of omega-average of frequency.
\[\sigma^2_{PDEV}(m\tau_0) = { 72 \over (N-2m) (m \tau_0 )^2 } \sum_{i=0}^{N-2m-1} \left[ \sum_{k=0}^{m-1} \left( { m-1 \over 2} - k \right) {x}_{i+k} - x_{i+k+m} \right]^2\]for \(m>1\) and for an averaging-factor of \(m=1\) PDEV equals ADEV/MDEV: \(\sigma_{PDEV}(\tau_0)=\sigma_{ADEV}(\tau_0)\).
Parameters
- data: np.array
Input data. Provide either phase or frequency (fractional, adimensional).
- rate: float
The sampling rate for data, in Hz. Defaults to 1.0
- data_type: {‘phase’, ‘freq’}
Data type, i.e. phase or frequency. Defaults to “phase”.
- taus: np.array
Array of tau values, in seconds, for which to compute statistic. Optionally set taus=[“all”|”octave”|”decade”] for automatic tau-list generation.
Returns
- (taus2, ad, ade, ns): tuple
Tuple of values
- taus2: np.array
Tau values for which td computed
- ad: np.array
Computed adev for each tau value
- ade: np.array
adev errors
- ns: np.array
Values of N used in each adev calculation
References
- allantools.hdev(data, rate=1.0, data_type='phase', taus=None)
- Hadamard deviation.
Rejects frequency drift, and handles divergent noise.
\[\sigma^2_{HDEV}( \tau ) = { 1 \over 6 \tau^2 (N-3) } \sum_{i=1}^{N-3} ( {x}_{i+3} - 3x_{i+2} + 3x_{i+1} - x_{i} )^2\]where \(x_i\) is the time-series of phase observations, spaced by the measurement interval \(\tau\), and with length \(N\).
Parameters
- datanp.array
Input data. Provide either phase or frequency (fractional, adimensional).
- ratefloat
The sampling rate for data, in Hz. Defaults to 1.0
- data_typestring, {‘phase’, ‘freq’}
Data type, i.e. phase or frequency. Defaults to “phase”.
- tausnp.array
Array of tau values, in seconds, for which to compute statistic. Optionally set taus=[“all”|”octave”|”decade”] for automatic tau-list generation.
References
NIST [SP1065] eqn (17) and (18), page 20
- allantools.ohdev(data, rate=1.0, data_type='phase', taus=None)
- Overlapping Hadamard deviation.
Better confidence than normal Hadamard.
\[\sigma^2_{OHDEV}(m\tau_0) = { 1 \over 6 (m \tau_0 )^2 (N-3m) } \sum_{i=1}^{N-3m} ( {x}_{i+3m} - 3x_{i+2m} + 3x_{i+m} - x_{i} )^2\]where \(x_i\) is the time-series of phase observations, spaced by the measurement interval \(\tau_0\), and with length \(N\).
Parameters
- data: np.array
Input data. Provide either phase or frequency (fractional, adimensional).
- rate: float
The sampling rate for data, in Hz. Defaults to 1.0
- data_type: {‘phase’, ‘freq’}
Data type, i.e. phase or frequency. Defaults to “phase”.
- taus: np.array
Array of tau values, in seconds, for which to compute statistic. Optionally set taus=[“all”|”octave”|”decade”] for automatic tau-list generation.
Returns
- (taus2, hd, hde, ns): tuple
Tuple of values
- taus2: np.array
Tau values for which td computed
- hd: np.array
Computed hdev for each tau value
- hde: np.array
hdev errors
- ns: np.array
Values of N used in each hdev calculation
References
NIST [SP1065] eqn (20), page 21
- allantools.tdev(data, rate=1.0, data_type='phase', taus=None)
Time deviation.
Based on modified Allan variance.
\[\sigma^2_{TDEV}( \tau ) = { \tau^2 \over 3 } \sigma^2_{MDEV}( \tau )\]Note that TDEV has a unit of seconds.
Parameters
- data: np.array
Input data. Provide either phase or frequency (fractional, adimensional).
- rate: float
The sampling rate for data, in Hz. Defaults to 1.0
- data_type: {‘phase’, ‘freq’}
Data type, i.e. phase or frequency. Defaults to “phase”.
- taus: np.array
Array of tau values, in seconds, for which to compute statistic. Optionally set taus=[“all”|”octave”|”decade”] for automatic tau-list generation.
Returns
- (taus, tdev, tdev_error, ns): tuple
Tuple of values
- taus: np.array
Tau values for which td computed
- tdev: np.array
Computed time deviations (in seconds) for each tau value
- tdev_errors: np.array
Time deviation errors
- ns: np.array
Values of N used in mdev_phase()
References
NIST [SP1065] eqn (15), page 18.
- allantools.totdev(data, rate=1.0, data_type='phase', taus=None)
- Total deviation.
Better confidence at long averages for Allan deviation.
\[\sigma^2_{TOTDEV}( m\tau_0 ) = { 1 \over 2 (m\tau_0)^2 (N-2) } \sum_{i=2}^{N-1} ( {x}^*_{i-m} - 2x^*_{i} + x^*_{i+m} )^2\]Where \(x^*_i\) is a new time-series of length \(3N-4\) derived from the original phase time-series \(x_n\) of length \(N\) by reflection at both ends. The original data \(x_n\) is in the center of \(x^*\):
\[ \begin{align}\begin{aligned}x^*_{1-j} = 2x_1 - x_{1+j} \quad \text{for} j=1..N-2\\x^*_i = x_i \quad \text{for} i=1..N\\x^*_{N+j} = 2x_N - x_{N-j} \quad \text{for} j=1..N-2\end{aligned}\end{align} \]FIXME: bias correction http://www.wriley.com/CI2.pdf page 5
Parameters
- phase: np.array
Phase data in seconds. Provide either phase or frequency.
- frequency: np.array
Fractional frequency data (nondimensional). Provide either frequency or phase.
- rate: float
The sampling rate for phase or frequency, in Hz
- taus: np.array
Array of tau values for which to compute measurement
References
NIST [SP1065] eqn (25) page 23
- allantools.mtotdev(data, rate=1.0, data_type='phase', taus=None)
Modified Total deviation.
Better confidence at long averages for modified Allan
FIXME: bias-correction http://www.wriley.com/CI2.pdf page 6
The variance is scaled up (divided by this number) based on the noise-type identified.
noise type
bias correction
WPM
0.94
FPM
0.83
WFM
0.73
FFM
0.70
RWFM
0.69
Parameters
- data: np.array
Input data. Provide either phase or frequency (fractional, adimensional).
- rate: float
The sampling rate for data, in Hz. Defaults to 1.0
- data_type: {‘phase’, ‘freq’}
Data type, i.e. phase or frequency. Defaults to “phase”.
- taus: np.array
Array of tau values, in seconds, for which to compute statistic. Optionally set taus=[“all”|”octave”|”decade”] for automatic tau-list generation.
References
NIST [SP1065] eqn (27) page 25
- allantools.ttotdev(data, rate=1.0, data_type='phase', taus=None)
Time Total Deviation
Modified total variance scaled by \(\tau^2 / 3\)
\[\sigma^2_{TTOTDEV}( \tau ) = { \tau^2 \over 3 } \sigma^2_{MTOTDEV}( \tau )\]Note that [SP1065] erroneously has tau-cubed here (!).
References
NIST [SP1065] eqn (28) page 26.
- allantools.htotdev(data, rate=1.0, data_type='phase', taus=None)
Hadamard Total deviation.
Better confidence at long averages for Hadamard deviation
PRELIMINARY - REQUIRES FURTHER TESTING.
Computed for N fractional frequency points y_i with sampling period tau0, analyzed at tau = m*tau0 1. remove linear trend by averaging first and last half, and dividing by interval 2. extend sequence by uninverted even reflection 3. compute Hadamard for extended, length 9m, sequence.
FIXME: bias corrections from http://www.wriley.com/CI2.pdf W FM 0.995 alpha= 0 F FM 0.851 alpha=-1 RW FM 0.771 alpha=-2 FW FM 0.717 alpha=-3 RR FM 0.679 alpha=-4
Parameters
- data: np.array
Input data. Provide either phase or frequency (fractional, adimensional).
- rate: float
The sampling rate for data, in Hz. Defaults to 1.0
- data_type: {‘phase’, ‘freq’}
Data type, i.e. phase or frequency. Defaults to “phase”.
- taus: np.array
Array of tau values, in seconds, for which to compute statistic. Optionally set taus=[“all”|”octave”|”decade”] for automatic tau-list generation.
- allantools.theo1(data, rate=1.0, data_type='phase', taus=None)
Theo1 is a two-sample variance with improved confidence and extended averaging factor range.
PRELIMINARY - REQUIRES FURTHER TESTING.
\[\sigma^2_{THEO1}(m\tau_0) = { 1 \over (m \tau_0 )^2 (N-m) } \sum_{i=1}^{N-m} \sum_{\delta=0}^{m/2-1} {1\over m/2-\delta}\lbrace ({x}_{i} - x_{i-\delta +m/2}) + (x_{i+m}- x_{i+\delta +m/2}) \rbrace^2\]Where \(10<=m<=N-1\) is even.
FIXME: bias correction
Parameters
- data: np.array
Input data. Provide either phase or frequency (fractional, adimensional).
- rate: float
The sampling rate for data, in Hz. Defaults to 1.0
- data_type: {‘phase’, ‘freq’}
Data type, i.e. phase or frequency. Defaults to “phase”.
- taus: np.array
Array of tau values, in seconds, for which to compute statistic. Optionally set taus=[“all”|”octave”|”decade”] for automatic tau-list generation.
References
NIST [SP1065] eq (30) page 29
- allantools.mtie(data, rate=1.0, data_type='phase', taus=None)
Maximum Time Interval Error.
Parameters
- data: np.array
Input data. Provide either phase or frequency (fractional, adimensional).
- rate: float
The sampling rate for data, in Hz. Defaults to 1.0
- data_type: {‘phase’, ‘freq’}
Data type, i.e. phase or frequency. Defaults to “phase”.
- taus: np.array
Array of tau values, in seconds, for which to compute statistic. Optionally set taus=[“all”|”octave”|”decade”] for automatic tau-list generation.
Notes
this seems to correspond to Stable32 setting “Fast(u)” Stable32 also has “Decade” and “Octave” modes where the dataset is extended somehow?
- allantools.tierms(data, rate=1.0, data_type='phase', taus=None)
Time Interval Error RMS.
Parameters
- data: np.array
Input data. Provide either phase or frequency (fractional, adimensional).
- rate: float
The sampling rate for data, in Hz. Defaults to 1.0
- data_type: {‘phase’, ‘freq’}
Data type, i.e. phase or frequency. Defaults to “phase”.
- taus: np.array
Array of tau values, in seconds, for which to compute statistic. Optionally set taus=[“all”|”octave”|”decade”] for automatic tau-list generation.
- allantools.gcodev(data_1, data_2, rate=1.0, data_type='phase', taus=None)
Groslambert codeviation a.k.a. Allan Covariance
Similarly to the three-cornered hat method, we consider three uncorrelated oscillators A, B, C. The Groslambert codeviation estimates the noise of one oscillator (e.g. B), given two synchronous measurements AB and BC. Unlike three-cordenred hat, Gcodev is not affected by the (uncorrelated) noise of the measurement devices (time-interval or frequency counter) used for the measurements AB and BC.
Parameters
- data_1: np.array
Oscillator 1 input data. Provide either phase or frequency
- data_2: np.array
Oscillator 2 input data. Provide either phase or frequency
- rate: float
The sampling rate for data, in Hz. Defaults to 1.0
- data_type: {‘phase’, ‘freq’}
Data type, i.e. phase or frequency. Defaults to “phase”.
- taus: np.array
Array of tau values, in seconds, for which to compute statistic. Optionally set taus=[“all”|”octave”|”decade”] for automatic tau-list generation.
Returns
- (taus, gd): tuple
Tuple of values
- taus: np.array
Tau values for which gcodev computed
- gd: np.array
Computed gcodev for each tau value
References
Real-Time Statistics
- class allantools.oadev_realtime(afs=[1], tau0=1.0, auto_afs=False, pts_per_decade=4)
Overlapping Allan deviation in real-time from a stream of phase/frequency samples.
Reference [Dobrogowski2007]
- add_phase(xnew)
add new phase point, in units of seconds
- update_S(idx)
update S, sum-of-squares
- class allantools.tdev_realtime(afs=[1], tau0=1.0, auto_afs=False, pts_per_decade=4)
Time deviation and Modified Allan deviation in real-time from a stream of phase/frequency samples.
Reference [Dobrogowski2007]
- add_phase(xnew)
add new phase point
- mdev()
scale tdev to output mdev
- update_S(idx)
update S, sum-of-squares
- update_S3n(idx)
eqn (13) of paper
- class allantools.ohdev_realtime(afs=[1], tau0=1.0, auto_afs=False, pts_per_decade=4)
Overlapping Hadamard deviation in real-time from a stream of phase/frequency samples.
Reference [Dobrogowski2007]
- add_phase(xnew)
add new phase point
- update_S(idx)
update S, sum-of-squares
Noise Generation
- class allantools.noise_kasdin.Noise(nr=2, qd=1, b=0)
Generate discrete colored noise
Python / Numpy implementation of [Kasdin1992] Kasdin, N.J., Walter, T., “Discrete simulation of power law noise [for oscillator stability evaluation],” Frequency Control Symposium, 1992. 46th., Proceedings of the 1992 IEEE, pp.274,283, 27-29 May 1992 http://dx.doi.org/10.1109/FREQ.1992.270003
Parameters
- nr: integer
length of generated time-series must be power of two
- qd: float
discrete variance
- b: float
noise type
b
noise type
0
White Phase Modulation (WPM)
-1
Flicker Phase Modulation (FPM)
-2
White Frequency Modulation (WFM)
-3
Flicker Frequency Modulation (FFM)
-4
Random Walk Frequency Modulation (RWFM)
Returns
- Noise()
A Noise() instance
- Example:
import numpy as np noise = allantools.Noise(nr=2*8, qd=1.0e-20, b=-1) noise.generateNoise() print noise.time_series
- adev(tau0, tau)
return predicted ADEV of noise-type at given tau
- adev_from_qd(tau0=1.0, tau=1.0)
prefactor for Allan deviation for noise type defined by (qd, b, tau0)
Colored noise generated with (qd, b, tau0) parameters will show an Allan variance of:
\[AVAR = prefactor \cdot h_a \cdot \tau^c\]where \(a = b + 2\) is the slope of the frequency PSD. and \(h_a\) is the frequency PSD prefactor \(S_y(f) = h_a f^a\)
The relation between a, b, c is:
a
b
c(AVAR)
c(MVAR)
-2
-4
1
1
-1
-3
0
0
0
-2
-1
-1
+1
-1
-2
-2
+2
0
-2
-3
Coefficients from [Dawkins2007].
Vernotte2015 Table I
- c_avar()
return tau exponent “c” for noise type. AVAR = prefactor * h_a * tau^c
- c_mvar()
return tau exponent “c” for noise type. MVAR = prefactor * h_a * tau^c
- frequency_psd_from_qd(tau0=1.0)
return frequency power spectral density coefficient \(h_a\) for the noise type defined by (qd, b, tau0)
Colored noise generated with (qd, b, tau0) parameters will show a frequency power spectral density of
\[S_y(f) = Frequency_{PSD}(f) = h_a f^a\]where the slope \(a\) comes from the phase PSD slope \(b\): \(a = b + 2\)
[Kasdin1992] eqn (39)
- generateNoise()
Generate noise time series based on input parameters
Returns
- time_series: np.array
Time series with colored noise. len(time_series) == nr
- mdev(tau0, tau)
return predicted MDEV of noise-type at given tau
- mdev_from_qd(tau0=1.0, tau=1.0)
prefactor for Modified Allan deviation for noise type defined by (qd, b, tau0)
Colored noise generated with (qd, b, tau0) parameters will show an Modified Allan variance of:
\[MVAR = prefactor \cdot h_a \cdot \tau^c\]where \(a = b + 2\) is the slope of the frequency PSD. and \(h_a\) is the frequency PSD prefactor \(S_y(f) = h_a f^a\)
- pdev_from_qd(tau0=1.0, tau=1.0)
prefactor for Parabolic Allan deviation for noise type defined by (qd, b, tau0)
Colored noise generated with (qd, b, tau0) parameters will show an Parabolic Allan variance of:
\[PVAR = prefactor \cdot h_a \cdot \tau^c\]where \(a = b + 2\) is the slope of the frequency PSD. and \(h_a\) is the frequency PSD prefactor \(S_y(f) = h_a f^a\)
- phase_psd_from_qd(tau0=1.0)
return phase power spectral density coefficient \(g_b\) for noise-type defined by (qd, b, tau0) where tau0 is the interval between data points
Colored noise generated with (qd, b, tau0) parameters will show a phase power spectral density of \(S_x(f) = Phase_{PSD}(f) = g_b f^b\)
[Kasdin1992] eqn (39)
- set_input(nr=2, qd=1, b=0)
Set inputs after initialization
Parameters
- nr: integer
length of generated time-series number must be power of two
- qd: float
discrete variance
- b: float
noise type
b
noise type
0
White Phase Modulation (WPM)
-1
Flicker Phase Modulation (FPM)
-2
White Frequency Modulation (WFM)
-3
Flicker Frequency Modulation (FFM)
-4
Random Walk Frequency Modulation (RWFM)
- allantools.noise.white(num_points=1024, b0=1.0, fs=1.0)
White noise generator
Generate time series with white noise that has constant PSD = b0, up to the nyquist frequency fs/2.
The PSD is at ‘height’ b0 and extends from 0 Hz up to the nyquist frequency fs/2 (prefactor math.sqrt(b0*fs/2.0))
Parameters
- num_points: int, optional
number of samples
- b0: float, optional
desired power-spectral density in [X^2/Hz] where X is the unit of x
- fs: float, optional
sampling frequency, i.e. 1/fs is the time-interval between datapoints
Returns
White noise sample: numpy.array
- allantools.noise.brown(num_points=1024, b_minus2=1.0, fs=1.0)
Brownian or random walk (diffusion) noise with 1/f^2 PSD
Not really a color… rather Brownian or random-walk. Obtained by integrating white-noise.
Parameters
- num_points: int, optional
number of samples
- b_minus2: float, optional
desired power-spectral density is b2*f^-2
- fs: float, optional
sampling frequency, i.e. 1/fs is the time-interval between datapoints
Returns
Random walk sample: numpy.array
- allantools.noise.violet(num_points=1024, b2=1, fs=1)
Violet noise with f^2 PSD
Obtained by differentiating white noise
Parameters
- num_points: int, optional
number of samples
- b2: float, optional
desired power-spectral density is b2*f^2
- fs: float, optional
sampling frequency, i.e. 1/fs is the time-interval between datapoints
Returns
Violet noise sample: numpy.array
- allantools.noise.pink(num_points=1024, depth=80)
Pink noise (approximation) with 1/f PSD
Fills a sample with results from a pink noise generator from http://pydoc.net/Python/lmj.sound/0.1.1/lmj.sound.noise/, based on the Voss-McCartney algorithm, discussion and code examples at http://www.firstpr.com.au/dsp/pink-noise/
Parameters
- num_points: int, optional
number of samples
- depth: int, optional
number of iteration for each point. High numbers are slower but generates a more correct spectrum on low-frequencies end.
Returns
Pink noise sample: numpy.array
Utilities
- allantools.frequency2phase(freqdata, rate)
integrate fractional frequency data and output phase data
Parameters
- freqdata: np.array
Data array of fractional frequency measurements (nondimensional)
- rate: float
The sampling rate for phase or frequency, in Hz
Returns
- phasedata: np.array
Time integral of fractional frequency data, i.e. phase (time) data in units of seconds. For phase in units of radians, see phase2radians()
- allantools.phase2frequency(phase, rate)
Convert phase in seconds to fractional frequency
Parameters
- phase: np.array
Data array of phase in seconds, length N
- rate: float
The sampling rate for phase, in Hz
Returns
- y: np.array
Data array of fractional frequency, length N-1
- allantools.phase2radians(phasedata, v0)
Convert phase in seconds to phase in radians
Parameters
- phasedata: np.array
Data array of phase in seconds
- v0: float
Nominal oscillator frequency in Hz
Returns
- fi: np.array
phase data in radians
- allantools.psd2allan(S_y, f=1.0, kind='adev', base=2)
- Convert a given (one-sided) power spectral density \(S_y(f)\) to Allan
deviation or modified Allan deviation
For ergodic noise, the Allan variance or modified Allan variance is related to the power spectral density \(S_y\) of the fractional frequency deviation:
\[\sigma^2_y(\tau) = 2 \int_0^\infty S_y(f) \left| \sin(\pi f \tau)^{(k+1)} \over (\pi f \tau)^k \right|^2 df,\]where \(f\) is the Fourier frequency and \(\tau\) the averaging time. The exponent \(k\) is 1 for the Allan variance and 2 for the modified Allan variance.
psd2allan() implements the integral by discrete numerical integration via a sum.
Parameters
- S_y: np.array
Single-sided power spectral density (PSD) of fractional frequency deviation S_y in 1/Hz^2. First element is S_y(f=0).
- f: np.array or scalar numeric (float or int)
if np.array: Spectral frequency vector in Hz if numeric scalar: Spectral frequency step in Hz default: Spectral frequency step 1 Hz
- kind: {‘adev’, ‘mdev’}
Which kind of Allan deviation to compute. Defaults to ‘adev’
- base: float
Base for logarithmic spacing of tau values. E.g. base= 10: decade, base= 2: octave, base <= 1: all
Returns
- (taus_used, ad): tuple
Tuple of 2 values
- taus_used: np.array
tau values for which ad computed
- ad: np.array
Computed Allan deviation of requested kind for each tau value
References
NIST [SP1065] eqs (65-66), page 73.
[Benkler2015] eqn (23).
- allantools.tau_generator(data, rate, taus=None, v=False, even=False, maximum_m=-1)
pre-processing of the tau-list given by the user (Helper function)
Does sanity checks, sorts data, removes duplicates and invalid values. Generates a tau-list based on keywords ‘all’, ‘decade’, ‘octave’. Uses ‘octave’ by default if no taus= argument is given.
Parameters
- data: np.array
data array
- rate: float
Sample rate of data in Hz. Time interval between measurements is 1/rate seconds.
- taus: np.array
Array of tau values for which to compute measurement. Alternatively one of the keywords: “all”, “octave”, “decade”. Defaults to “octave” if omitted.
keyword
averaging-factors
“all”
1, 2, 3, 4, …, len(data)
“octave”
1, 2, 4, 8, 16, 32, …
“decade”
1, 2, 4, 10, 20, 40, 100, …
“log10”
approx. 10 points per decade
- v: bool
verbose output if True
- even: bool
require even m, where tau=m*tau0, for Theo1 statistic
- maximum_m: int
limit m, where tau=m*tau0, to this value. used by mtotdev() and htotdev() to limit maximum tau.
Returns
- (data, m, taus): tuple
List of computed values
- data: np.array
Data
- m: np.array
Tau in units of data points
- taus: np.array
Cleaned up list of tau values
- allantools.edf_simple(N, m, alpha)
Equivalent degrees of freedom. Simple approximate formulae.
Parameters
- Nint
the number of phase samples
- mint
averaging factor, tau = m * tau0
- alpha: int
exponent of f for the frequency PSD: ‘wp’ returns white phase noise. alpha=+2 ‘wf’ returns white frequency noise. alpha= 0 ‘fp’ returns flicker phase noise. alpha=+1 ‘ff’ returns flicker frequency noise. alpha=-1 ‘rf’ returns random walk frequency noise. alpha=-2 If the input is not recognized, it defaults to idealized, uncorrelated noise with (N-1) degrees of freedom.
Notes
See [Stein1985].
Returns
- edffloat
Equivalent degrees of freedom
- allantools.edf_greenhall(alpha, d, m, N, overlapping=False, modified=False, verbose=False)
returns Equivalent degrees of freedom
Parameters
- alpha: int
noise type, +2…-4
- d: int
1 first-difference variance 2 Allan variance 3 Hadamard variance require alpha+2*d>1
- m: int
averaging factor tau = m*tau0 = m*(1/rate)
- N: int
number of phase observations (length of time-series)
- overlapping: bool
True for oadev, ohdev
- modified: bool
True for mdev, tdev
Returns
- edf: float
Equivalent degrees of freedom
Notes
Reference [Greenhall2004].
Used for the following deviations (see [Riley_CI] page 8) adev() oadev() mdev() tdev() hdev() ohdev()
- allantools.edf_totdev(N, m, alpha)
Equivalent degrees of freedom for Total Deviation FIXME: what is the right behavior for alpha outside 0,-1,-2?
NIST [SP1065] page 41, Table 7
- allantools.edf_mtotdev(N, m, alpha)
Equivalent degrees of freedom for Modified Total Deviation
NIST [SP1065] page 41, Table 8
- allantools.three_cornered_hat_phase(phasedata_ab, phasedata_bc, phasedata_ca, rate, taus, function)
Three Cornered Hat Method
Given three clocks A, B, C, we seek to find their variances \(\sigma^2_A\), \(\sigma^2_B\), \(\sigma^2_C\). We measure three phase differences, assuming no correlation between the clocks, the measurements have variances:
\[ \begin{align}\begin{aligned}\sigma^2_{AB} = \sigma^2_{A} + \sigma^2_{B}\\\sigma^2_{BC} = \sigma^2_{B} + \sigma^2_{C}\\\sigma^2_{CA} = \sigma^2_{C} + \sigma^2_{A}\end{aligned}\end{align} \]Which allows solving for the variance of one clock as:
\[\sigma^2_{A} = {1 \over 2} ( \sigma^2_{AB} + \sigma^2_{CA} - \sigma^2_{BC} )\]and similarly cyclic permutations for \(\sigma^2_B\) and \(\sigma^2_C\)
Parameters
- phasedata_ab: np.array
phase measurements between clock A and B, in seconds
- phasedata_bc: np.array
phase measurements between clock B and C, in seconds
- phasedata_ca: np.array
phase measurements between clock C and A, in seconds
- rate: float
The sampling rate for phase, in Hz
- taus: np.array
The tau values for deviations, in seconds
- function: allantools deviation function
The type of statistic to compute, e.g. allantools.oadev
Returns
- tau_ab: np.array
Tau values corresponding to output deviations
- dev_a: np.array
List of computed values for clock A
References