Анализ шума инерциального датчика с использованием дисперсии Аллана

Открыть сценарий

В этом примере показано, как использовать дисперсию Аллана для определения параметров шума гироскопа МЭМС. Эти параметры можно использовать для моделирования гироскопа при моделировании. Измерение гироскопа моделируют следующим образом:

$\Omega(t) = \Omega_{Ideal}(t) + Bias_N(t) + Bias_B(t) + Bias_K(t)$

Три параметра шума N (угловая случайная ходьба), K (скоростная случайная ходьба) и B (нестабильность смещения) оцениваются с использованием данных, записанных из стационарного гироскопа.

Фон

Дисперсию Аллана первоначально разработал Дэвид У. Аллан для измерения стабильности частоты прецизионных осцилляторов. Он также может использоваться для идентификации различных источников шума, присутствующих в стационарных измерениях гироскопа. Рассмотрим L образцов данных с гироскопа с временем выборки. $\tau_{0}$ Сформировать кластеры данных длительностью,,... и $\tau_{0}$ $2\tau_{0}$ $m\tau_{0}, (m < (L-1)/2)$ получить средние значения суммы точек данных, содержащихся в каждом кластере по длине кластера. Дисперсию Аллана определяют как дисперсию двух выборок средних значений кластера данных как функцию времени кластера. В этом примере используется перекрывающийся оценщик дисперсии Аллана. Это означает, что вычисляемые кластеры перекрываются. Оценщик работает лучше, чем неперекрывающиеся оценщики для больших значений L.

Расчет отклонения Аллана

Дисперсию Аллана рассчитывают следующим образом:

Log L стационарные образцы гироскопа с периодом образца. $\tau_{0}$ Пусть $\Omega$ это будут зарегистрированные образцы.

% Load logged data from one axis of a three-axis gyroscope. This recording
% was done over a six hour period with a 100 Hz sampling rate.
load('LoggedSingleAxisGyroscope', 'omega', 'Fs')
t0 = 1/Fs;

Для каждого образца рассчитайте выходной угол: $\theta$

$\theta(t) = \int^{t}\Omega(t')dt'$

Для дискретных выборок $\tau_{0}$ можно использовать суммарную сумму, умноженную на.

theta = cumsum(omega, 1)*t0;

Затем вычислите дисперсию Аллана:

$\sigma^2(\tau) =
\frac{1}{2\tau^2}<(\theta_{k+2m}-2\theta_{k+m}+\theta_{k})^2>$

где $\tau = m\tau_{0}$ и $<>$ - среднее значение ансамбля.

Среднее значение ансамбля можно расширить до:

$\sigma^2(\tau) =
\frac{1}{2\tau^2(L-2m)}\sum_{k=1}^{L-2m}(\theta_{k+2m} - 2\theta_{k+m}
+ \theta_{k})^2$

maxNumM = 100;
L = size(theta, 1);
maxM = 2.^floor(log2(L/2));
m = logspace(log10(1), log10(maxM), maxNumM).';
m = ceil(m); % m must be an integer.
m = unique(m); % Remove duplicates.

tau = m*t0;

avar = zeros(numel(m), 1);
for i = 1:numel(m)
    mi = m(i);
    avar(i,:) = sum( ...
        (theta(1+2*mi:L) - 2*theta(1+mi:L-mi) + theta(1:L-2*mi)).^2, 1);
end
avar = avar ./ (2*tau.^2 .* (L - 2*m));

Наконец, отклонение Аллана $\sigma(t) = \sqrt{\sigma^2(t)}$ используется для определения параметров шума гироскопа.

adev = sqrt(avar);

figure
loglog(tau, adev)
title('Allan Deviation')
xlabel('\tau');
ylabel('\sigma(\tau)')
grid on
axis equal

Отклонение Аллана также может быть вычислено с помощью allanvar функция.

[avarFromFunc, tauFromFunc] = allanvar(omega, m, Fs);
adevFromFunc = sqrt(avarFromFunc);

figure
loglog(tau, adev, tauFromFunc, adevFromFunc);
title('Allan Deviations')
xlabel('\tau')
ylabel('\sigma(\tau)')
legend('Manual Calculation', 'allanvar Function')
grid on
axis equal

Идентификация параметров шума

Чтобы получить параметры шума для гироскопа, используйте следующую зависимость между дисперсией Аллана и двусторонней спектральной плотностью мощности (PSD) параметров шума в исходном наборе данных. $\Omega$ Отношения таковы:

$\sigma^2(\tau) = 4\int_{0}^{\infty}S_\Omega(f)
\frac{sin^4(\pi f\tau)}{(\pi f\tau)^2}df$

Из приведенного выше уравнения дисперсия Аллана пропорциональна суммарной мощности шума гироскопа при прохождении через фильтр с передаточной функцией. $sin^4(x)/(x)^2$ Эта передаточная функция возникает в результате операций, выполняемых для создания кластеров и управления ими.

Использование этой интерпретации передаточной функции зависит от полосы пропускания фильтра. $\tau$ Это означает, что различные параметры шума могут быть идентифицированы путем изменения полосы пропускания фильтра или изменения. $\tau$

Случайная прогулка под углом

Угловая случайная ходьба характеризуется спектром белого шума выхода гироскопа. PSD представлен следующим образом:

$S_\Omega(f) = N^2$

где

N = коэффициент случайного обхода под углом

Подстановка в исходное уравнение PSD и выполнение интеграционных выходов:

$\sigma^2(\tau) = \frac{N^2}{\tau}$

Вышеприведенное уравнение представляет собой линию с наклоном -1/2 при построении на логарифмическом графике $\sigma(\tau)$ зависимости. $\tau$ Значение N можно считать непосредственно за пределами этой строки. $\tau = 1$

% Find the index where the slope of the log-scaled Allan deviation is equal
% to the slope specified.
slope = -0.5;
logtau = log10(tau);
logadev = log10(adev);
dlogadev = diff(logadev) ./ diff(logtau);
[~, i] = min(abs(dlogadev - slope));

% Find the y-intercept of the line.
b = logadev(i) - slope*logtau(i);

% Determine the angle random walk coefficient from the line.
logN = slope*log(1) + b;
N = 10^logN

% Plot the results.
tauN = 1;
lineN = N ./ sqrt(tau);
figure
loglog(tau, adev, tau, lineN, '--', tauN, N, 'o')
title('Allan Deviation with Angle Random Walk')
xlabel('\tau')
ylabel('\sigma(\tau)')
legend('\sigma', '\sigma_N')
text(tauN, N, 'N')
grid on
axis equal

N =

    0.0126

Оценить случайную прогулку

Скорость случайной ходьбы характеризуется спектром красного шума (броуновского шума) выхода гироскопа. PSD представлен следующим образом:

$S_\Omega(f) = (\frac{K}{2\pi})^2\frac{1}{f^2}$

где

K = коэффициент случайного обхода

Подстановка в исходное уравнение PSD и выполнение интеграционных выходов:

$\sigma^2(\tau) = \frac{K^2\tau}{3}$

Вышеприведенное уравнение представляет собой линию с наклоном 1/2 при построении на логарифмическом графике $\sigma(\tau)$ зависимости. $\tau$ Значение K можно считать непосредственно за пределами этой строки в. $\tau = 3$

% Find the index where the slope of the log-scaled Allan deviation is equal
% to the slope specified.
slope = 0.5;
logtau = log10(tau);
logadev = log10(adev);
dlogadev = diff(logadev) ./ diff(logtau);
[~, i] = min(abs(dlogadev - slope));

% Find the y-intercept of the line.
b = logadev(i) - slope*logtau(i);

% Determine the rate random walk coefficient from the line.
logK = slope*log10(3) + b;
K = 10^logK

% Plot the results.
tauK = 3;
lineK = K .* sqrt(tau/3);
figure
loglog(tau, adev, tau, lineK, '--', tauK, K, 'o')
title('Allan Deviation with Rate Random Walk')
xlabel('\tau')
ylabel('\sigma(\tau)')
legend('\sigma', '\sigma_K')
text(tauK, K, 'K')
grid on
axis equal

K =

   9.0679e-05

Неустойчивость смещения

Нестабильность смещения характеризуется спектром розового шума (мерцающего шума) выхода гироскопа. PSD представлен следующим образом:

$S_{\Omega}(f) = \left\{\begin{array}{lr}(\frac{B^2}{2\pi})\frac{1}{f}
& : f \leq f_0\\0 & : f > f_0\end{array}\right.$

где

B = коэффициент нестабильности смещения

$f_0$ = частота отсечки

Подстановка в исходное уравнение PSD и выполнение интеграционных выходов:

$\sigma^2(\tau) = \frac{2B^2}{\pi}[\ln{2} + \\
-\frac{sin^3x}{2x^2}(sinx + 4xcosx) + Ci(2x) - Ci(4x)]$

где

$x = \pi f_0\tau$

Ci = косинусно-интегральная функция

Когда $\tau$ это намного длиннее, чем обратная частота отсечки, уравнение PSD равно:

$\sigma^2(\tau) = \frac{2B^2}{\pi}\ln{2}$

Вышеприведенное уравнение представляет собой линию с наклоном 0 при построении на логарифмическом графике $\sigma(\tau)$ зависимости. $\tau$ Значение B можно считывать непосредственно с этой строки с масштабированием. $\sqrt{\frac{2\ln{2}}{\pi}} \approx 0.664$

% Find the index where the slope of the log-scaled Allan deviation is equal
% to the slope specified.
slope = 0;
logtau = log10(tau);
logadev = log10(adev);
dlogadev = diff(logadev) ./ diff(logtau);
[~, i] = min(abs(dlogadev - slope));

% Find the y-intercept of the line.
b = logadev(i) - slope*logtau(i);

% Determine the bias instability coefficient from the line.
scfB = sqrt(2*log(2)/pi);
logB = b - log10(scfB);
B = 10^logB

% Plot the results.
tauB = tau(i);
lineB = B * scfB * ones(size(tau));
figure
loglog(tau, adev, tau, lineB, '--', tauB, scfB*B, 'o')
title('Allan Deviation with Bias Instability')
xlabel('\tau')
ylabel('\sigma(\tau)')
legend('\sigma', '\sigma_B')
text(tauB, scfB*B, '0.664B')
grid on
axis equal

B =

    0.0020

Теперь, когда все параметры шума были рассчитаны, постройте график отклонения Аллана со всеми линиями, используемыми для количественной оценки параметров.

tauParams = [tauN, tauK, tauB];
params = [N, K, scfB*B];
figure
loglog(tau, adev, tau, [lineN, lineK, lineB], '--', ...
    tauParams, params, 'o')
title('Allan Deviation with Noise Parameters')
xlabel('\tau')
ylabel('\sigma(\tau)')
legend('\sigma', '\sigma_N', '\sigma_K', '\sigma_B')
text(tauParams, params, {'N', 'K', '0.664B'})
grid on
axis equal

Моделирование гироскопа

Используйте imuSensor изобретение позволяет моделировать измерения гироскопа шумовыми параметрами, определенными выше.

% Simulating the gyroscope measurements takes some time. To avoid this, the
% measurements were generated and saved to a MAT-file. By default, this
% example uses the MAT-file. To generate the measurements instead, change
% this logical variable to true.
generateSimulatedData = false;

if generateSimulatedData
    % Set the gyroscope parameters to the noise parameters determined
    % above.
    gyro = gyroparams('NoiseDensity', N, 'RandomWalk', K, ...
        'BiasInstability', B);
    omegaSim = helperAllanVarianceExample(L, Fs, gyro);
else
    load('SimulatedSingleAxisGyroscope', 'omegaSim')
end

Рассчитайте смоделированное отклонение Аллана и сравните его с записанными данными.

[avarSim, tauSim] = allanvar(omegaSim, 'octave', Fs);
adevSim = sqrt(avarSim);
adevSim = mean(adevSim, 2); % Use the mean of the simulations.

figure
loglog(tau, adev, tauSim, adevSim, '--')
title('Allan Deviation of HW and Simulation')
xlabel('\tau');
ylabel('\sigma(\tau)')
legend('HW', 'SIM')
grid on
axis equal

На графике показано, что модель гироскопа, созданная из imuSensor генерирует измерения со сходным с записанными данными отклонением Аллана. Измерения модели содержат немного меньше шума, поскольку параметры квантования и температуры не задаются с использованием gyroparams. Модель гироскопа может быть использована для генерации измерений с использованием движений, которые не могут быть легко зафиксированы аппаратными средствами.

Ссылки

Стандарт IEEE. 647-2006 Руководство по формату стандартных спецификаций IEEE и процедура испытаний одноосных лазерных гироскопов

Документация

Анализ шума инерциального датчика с использованием дисперсии Аллана

Фон

Расчет отклонения Аллана

Идентификация параметров шума

Случайная прогулка под углом

Оценить случайную прогулку

Неустойчивость смещения

Моделирование гироскопа

Ссылки

Документация по комплекту инструментов для слияния и отслеживания датчиков

Поддержка