Contenido principal

Esta página se ha traducido mediante traducción automática. Haga clic aquí para ver la última versión en inglés.

Análisis de ruido del sensor inercial utilizando la varianza de Allan

Este ejemplo muestra cómo utilizar la varianza de Allan para determinar los parámetros de ruido de un giroscopio MEMS. Estos parámetros se pueden utilizar para modelar el giroscopio en simulación. La medida del giroscopio se modela como:

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

Los tres parámetros de ruido N (recorrido aleatorio del ángulo), K (recorrido aleatorio de velocidad) y B (inestabilidad de sesgo) se estiman utilizando datos registrados desde un giróscopo estacionario.

Fondo

La varianza de Allan fue desarrollada originalmente por David W. Allan para medir la estabilidad de frecuencia de osciladores de precisión. También se puede utilizar para identificar diversas fuentes de ruido presentes en mediciones de giroscopio estacionario. Considere L muestras de datos de un giroscopio con un tiempo de muestreo de $\tau_{0}$. Forme clústeres de datos con duraciones $\tau_{0}$, $2\tau_{0}$, ..., $m\tau_{0}, (m < (L-1)/2)$ y obtenga los promedios de la suma de los puntos de datos contenidos en cada clúster a lo largo de la longitud del clúster. La varianza de Allan se define como la varianza de dos muestras de los promedios del grupo de datos en función del tiempo del grupo. Este ejemplo utiliza el estimador de varianza de Allan superpuesto. Esto significa que los grupos calculados se superponen. El estimador funciona mejor que los estimadores no superpuestos para valores mayores de L.

Cálculo de la varianza de Allan

La varianza de Allan se calcula de la siguiente manera:

Registra L muestras de giroscopio estacionario con un período de muestra $\tau_{0}$. Sean $\Omega$ las muestras registradas.

% 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. The unit of
% the gyroscope data is in radians per second.
load('LoggedSingleAxisGyroscope', 'omega', 'Fs')
t0 = 1/Fs;

Para cada muestra, calcule el ángulo de salida $\theta$:

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

Para muestras discretas, se puede utilizar la suma acumulada multiplicada por $\tau_{0}$.

theta = cumsum(omega, 1)*t0;

A continuación, calcule la varianza de Allan:

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

donde $\tau = m\tau_{0}$ y $<>$ es el promedio del conjunto.

El promedio del ensemble se puede ampliar a:

$$\sigma^2(\tau) =&#10;\frac{1}{2\tau^2(L-2m)}\sum_{k=1}^{L-2m}(\theta_{k+2m} - 2\theta_{k+m}&#10;+ \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));

Finalmente, la desviación de Allan $\sigma(t) = \sqrt{\sigma^2(t)}$ se utiliza para determinar los parámetros de ruido del giroscopio.

adev = sqrt(avar);

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

La varianza de Allan también se puede calcular utilizando la función 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

Identificación de parámetros de ruido

Para obtener los parámetros de ruido del giroscopio, utilice la siguiente relación entre la varianza de Allan y la densidad espectral de potencia (PSD) bilateral de los parámetros de ruido en el conjunto de datos original $\Omega$. La relación es:

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

De la ecuación anterior, la varianza de Allan es proporcional a la potencia de ruido total del giroscopio cuando pasa a través de un filtro con una función de transferencia de $sin^4(x)/(x)^2$. Esta función de transferencia surge de las operaciones realizadas para crear y operar sobre los clusters.

Utilizando esta interpretación de la función de transferencia, la paso banda del filtro depende de $\tau$. Esto significa que se pueden identificar diferentes parámetros de ruido cambiando la paso banda del filtro o variando $\tau$.

Paseo aleatorio en ángulo

La caminata aleatoria en ángulo se caracteriza por el espectro de ruido blanco de la salida del giroscopio. El PSD está representado por:

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

dónde

N = coeficiente del recorrido aleatorio del ángulo

Sustituyendo en la ecuación PSD original y realizando la integración se obtiene:

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

La ecuación anterior es una línea con una pendiente de -1/2 cuando se traza en un gráfico logarítmico de $\sigma(\tau)$ versus $\tau$. El valor de N se puede leer directamente desde esta línea en $\tau = 1$. Las unidades de N son $(rad/s)/\sqrt{Hz}$.

% 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

Caminata aleatoria de la velocidad

La caminata aleatoria de velocidad se caracteriza por el espectro de ruido rojo (ruido browniano) de la salida del giroscopio. El PSD está representado por:

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

dónde

K = coeficiente de caminata aleatoria de velocidad

Sustituyendo en la ecuación PSD original y realizando la integración se obtiene:

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

La ecuación anterior es una línea con una pendiente de 1/2 cuando se traza en un gráfico logarítmico de $\sigma(\tau)$ versus $\tau$. El valor de K se puede leer directamente desde esta línea en $\tau = 3$. Las unidades de K son $(rad/s)\sqrt{Hz}$.

% 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

Inestabilidad del sesgo

La inestabilidad del sesgo se caracteriza por el espectro de ruido rosa (ruido de parpadeo) de la salida del giroscopio. El PSD está representado por:

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

dónde

B = coeficiente de inestabilidad de sesgo

$f_0$ = frecuencia de corte

Sustituyendo en la ecuación PSD original y realizando la integración se obtiene:

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

dónde

$x = \pi f_0\tau$

Ci = función integral del coseno

Cuando $\tau$ es mucho más largo que la inversa de la frecuencia de corte, la ecuación PSD es:

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

La ecuación anterior es una línea con una pendiente de 0 cuando se traza en un gráfico logarítmico de $\sigma(\tau)$ versus $\tau$. El valor de B se puede leer directamente desde esta línea con una escala de $\sqrt{\frac{2\ln{2}}{\pi}} \approx 0.664$. Las unidades de B son $rad/s$.

% 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

Ahora que se han calculado todos los parámetros de ruido, trace la desviación de Allan con todas las líneas utilizadas para cuantificar los parámetros.

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 (rad/s)$', '$\sigma_N ((rad/s)/\sqrt{Hz})$', ...
    '$\sigma_K ((rad/s)\sqrt{Hz})$', '$\sigma_B (rad/s)$', 'Interpreter', 'latex')
text(tauParams, params, {'N', 'K', '0.664B'})
grid on
axis equal

Simulación de giroscopio

Utilice el objeto imuSensor para simular mediciones de giroscopio con los parámetros de ruido identificados anteriormente.

% 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. Use single-sided noise to match the parameter scaling and 500
    % poles for the bias instability model to match the hardware output.
    numpoles = 500;
    gyro = gyroparams('NoiseDensity', N, 'RandomWalk', K, ...
        'BiasInstability', B, 'NoiseType', 'single-sided', ...
        'BiasInstabilityCoefficients', fractalcoef(numpoles));
    omegaSim = helperAllanVarianceExample(L, Fs, gyro);
else
    load('SimulatedSingleAxisGyroscope', 'omegaSim')
end

Calcule la desviación de Allan simulada y compárela con los datos registrados.

[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

La gráfica muestra que el modelo de giroscopio creado a partir de imuSensor genera mediciones con una desviación de Allan similar a los datos registrados. Las mediciones del modelo contienen un poco menos de ruido ya que los parámetros de cuantificación y relacionados con la temperatura no se establecen utilizando gyroparams. El modelo de giroscopio se puede utilizar para generar mediciones utilizando movimientos que no se capturan fácilmente con hardware.

Referencias

  • Estándar IEEE. Guía de formato de especificación estándar IEEE 647-2006 y procedimiento de prueba para giroscopios láser de un solo eje