Jak ekstrapolować sygnał 1D?

Mam sygnał o pewnej długości, powiedzmy 1000 próbek. Chciałbym rozszerzyć ten sygnał na 5000 próbek, próbkowanych z tą samą prędkością co oryginał (tj. Chcę przewidzieć, jaki byłby sygnał, gdybym kontynuował próbkowanie przez dłuższy czas). Sygnał składa się z kilku połączonych ze sobą elementów sinusoidalnych.

Metodą, która po raz pierwszy do mnie przyszła, było wzięcie całego FFT i rozszerzenie go, ale pozostawia to bardzo silną nieciągłość na klatce 1001. Rozważyłem również użycie części widma w pobliżu pików i chociaż wydaje się nieco poprawić sygnał, nie wydaje mi się, aby faza gwarantowana była poprawna. Jaka jest najlepsza metoda przedłużenia tego sygnału?

Oto kod MATLAB pokazujący idealną metodę tego, czego chcę. Oczywiście nie będę wcześniej wiedział, że są dokładnie 3 składowe sinusoidalne, ani ich dokładna faza i częstotliwość. Chcę się upewnić, że funkcja jest ciągła, że ​​nie ma skoku, gdy przechodzimy do punktu 501,

vals = 1:50;
signal = 100+5*sin(vals/3.7+.3)+3*sin(vals/1.3+.1)+2*sin(vals/34.7+.7); % This is the measured signal
% Note, the real signal will have noise and not be known exactly.
output_vals = 1:200;
output_signal = 100+5*sin(output_vals/3.7+.3)+3*sin(output_vals/1.3+.1)+2*sin(output_vals/34.7+.7); % This is the output signal

figure;
plot(output_signal);
hold all;
plot(signal);

Zasadniczo, biorąc pod uwagę zieloną linię, chcę znaleźć niebieską linię.

W zależności od materiału źródłowego metoda interpolacji widmowej oparta na DCT opisana w następującym artykule wygląda obiecująco:

Łk, HG, Güler S. „Transformacja i interpolacja sygnału na podstawie zmodyfikowanej syntezy DCT” , Cyfrowe przetwarzanie sygnałów, Artykuł w prasie, 2011.

Oto jedna z liczb z artykułu pokazująca przykład interpolacji:

Zastosowanie tej techniki do odzyskiwania utraconych segmentów ( np. 4.2. Synteza utraconych okresów podziałki ) jest prawdopodobnie najbardziej odpowiednie dla ekstrapolacji. Innymi słowy, weź kawałek istniejącego materiału źródłowego, udawaj, że między krawędzią a wybranym przez ciebie segmentem jest zgubiony odcinek, a następnie „zrekonstruuj” „brakujący” fragment (i być może odrzuć fragment użyty na końcu). Wydaje się, że zastosowanie tej techniki byłoby jeszcze prostsze w celu płynnego zapętlenia materiału źródłowego ( np. 3.1. Interpolacja amplitudy widmowej ).

Myślę, że to, czego szukasz, to liniowe kodowanie predykcyjne (zwane inaczej automatyczną regresją średnią ruchomą ). LPC ekstrapoluje szereg czasowy, najpierw dopasowując model liniowy do szeregów czasowych, w których zakłada się, że każda próbka jest liniową kombinacją poprzednich próbek. Po dopasowaniu tego modelu do istniejących szeregów czasowych można go uruchomić do przodu, aby ekstrapolować dalsze wartości, zachowując stacjonarne (?) Spektrum mocy.

Oto mały przykład w Matlabie, wykorzystujący lpcfunkcję do oszacowania współczynników LPC.

N = 150;    % Order of LPC auto-regressive model
P = 500;    % Number of samples in the extrapolated time series
M = 150;    % Point at which to start predicting

t = 1:P;

x = 5*sin(t/3.7+.3)+3*sin(t/1.3+.1)+2*sin(t/34.7+.7); %This is the measured signal

a = lpc(x, N);

y = zeros(1, P);

% fill in the known part of the time series
y(1:M) = x(1:M);

% in reality, you would use `filter` instead of the for-loop
for ii=(M+1):P      
    y(ii) = -sum(a(2:end) .* y((ii-1):-1:(ii-N)));
end

plot(t, x, t, y);
l = line(M*[1 1], get(gca, 'ylim'));
set(l, 'color', [0,0,0]);
legend('actual signal', 'extrapolated signal', 'start of extrapolation');

Oczywiście w prawdziwym kodzie zastosowałbyś filterekstrapolację, wykorzystując współczynniki LPC ajako filtr IIR i wstępnie ładując znane wartości szeregów czasowych do stanu filtru; coś takiego:

% Run the initial timeseries through the filter to get the filter state 
[~, zf] = filter(-[0 a(2:end)], 1, x(1:M));     

% Now use the filter as an IIR to extrapolate
y((M+1):P) = filter([0 0], -a, zeros(1, P-M), zf); 

Oto wynik:

Wykonuje rozsądną robotę, choć z jakiegoś powodu przewidywania umierają z czasem.

Tak naprawdę niewiele wiem o modelach AR i chciałbym dowiedzieć się więcej.

-

EDYCJA: @china i @Emre mają rację, metoda Burg wydaje się działać znacznie lepiej niż LPC. Po prostu zmieniając lpcsię arburgw powyższych wydajnością kodu następujące wyniki:

Kod jest dostępny tutaj: https://gist.github.com/2843661

„Ekstrapolacja” 1-D jest dość prosta dzięki zastosowaniu metody BURG do oszacowania współczynników LP. Gdy współczynniki LP są dostępne, można łatwo obliczyć próbki czasu, stosując filtr. Próbki, które są przewidywane za pomocą Burga, są kolejnymi próbkami twojego wejściowego segmentu czasu.

Jeśli masz całkowitą pewność, że sygnał ma tylko kilka składowych częstotliwości, możesz wypróbować algorytm MUSIC , aby dowiedzieć się, jakie częstotliwości są zawarte w twoim sygnale i spróbować z niego pracować. Nie jestem do końca pewien, czy można to zrobić idealnie.

Ponadto, ponieważ Twoje dane są całkowicie deterministyczne, możesz spróbować zbudować jakiś nieliniowy predyktor, trenować go przy użyciu istniejącego zestawu danych i pozwolić ekstrapolować resztę.

Zasadniczo jest to problem ekstrapolacji, czy możesz chcieć Google coś w rodzaju ekstrapolacji ceny Fouriera .