The spectrum of Prony
The Proni spectrum is a method for estimating the power spectrum of a random process based on an autoregressive model (AR). The method was proposed in 1969 by French electrical engineer Rene Proni. The basic idea is to approximate the time series as a sum of sinusoidal components and noise using autoregression parameters. Spectral analysis by the Proni method is particularly useful for isolating frequency components in noisy signals and is actively used in signal processing and radar systems.
Engee does not have built-in functions for constructing spectra using the Prony transform, however, in Julia you can implement the algorithm yourself or use available libraries to solve similar problems.
The purpose of this demonstration is to implement such an algorithm. The following is an example of the implementation of an algorithm for estimating the power spectrum using the Proni method.
# Установите необходимые пакеты
Pkg.add("LinearAlgebra")
using LinearAlgebra
Next, we implement the prony_spectrum function. This function performs a Prony conversion for a given signal. It takes three arguments: the signal itself (signal), the order of the model (order), and the sampling frequency (fs). The function returns two values: frequencies and their corresponding amplitudes, which represent the components of the signal spectrum. The process includes forming an autocorrelation matrix, solving a system of linear equations to obtain the coefficients of an autoregressive model, finding the roots of a characteristic equation, and calculating frequencies and amplitudes.
# Функция для преобразования Прони
function prony_spectrum(signal::Vector{Float64}, order::Int, fs::Float64)
N = length(signal)
if order >= N
error("Порядок модели должен быть меньше длины сигнала.")
end
# Формирование вектора автокорреляции
X = zeros(Float64, N - order, order)
for i in 1:(N - order)
X[i, :] = signal[i:(i + order - 1)]
end
y = signal[(order + 1):N]
a = X \ y # Решение задачи линейной регрессии для коэффициентов AR-модели
a = vcat(1.0, -a) # Коэффициенты AR модели
roots_a = roots(a) # Нахождение корней характеристического уравнения
roots_a = roots_a[abs.(roots_a) .< 1.0] # Оставляем только корни внутри единичной окружности
# Рассчитываем частоты и амплитуды
frequencies = angle.(roots_a) * fs / (2 * π)
amplitudes = abs.(roots_a)
return frequencies, amplitudes
end
Next, we will set the noisy input signal.
fs = 1000.0 # Частота дискретизации
t = 0:1/fs:1-1/fs # Временная шкала
# Сигнал из нескольких синусоид
signal = sin.(2π * 50 * t) + 0.5 * sin.(2π * 120 * t)
# Добавление шума
noisy_signal = signal + 0.05 * randn(length(signal))
plot(t, noisy_signal)
Now let's perform the calculation and construction of the Proni spectrum.
order = 10 # Порядок модели
frequencies, amplitudes = prony_spectrum(noisy_signal, order, fs)
# Визуализация
plot(frequencies, amplitudes, seriestype=:stem, yaxis=:log, xlabel="Частота (Гц)", ylabel="Амплитуда",
title="Спектр Прони", marker=:circle)
Conclusion
In this example, we have demonstrated how to build a spectrum in Engee using a custom spectrum calculation function.