Задать вопрос
sergiks
@sergiks
♬♬
математика

Как взять преобразование Фурье от дискретных точек?

Временные данные представляют из себя просто отсчёты времени, моменты, в которые прозошло событие. Напр. события отмечены в 1-ю секунду, 3, 5 и т.д.: t = [1, 3, 5, 7, 9].

Как правильно получить частотный спектр этих данных? Ожидаю, что должен вылезти пик на амплитуде "2".

Пробовал выбрать некую частоту дискретизации, напр. с шагом в 0.5 секунды, и сделать расчёску из "0", а там, где были события – "1". И от этого набора брать fft()

Но получается совсем не то, что ожидал:
f8371250cab64fe7a83fc75f3de84f5a.png

Мой код Matlab этой попытки:
A=[1,3,5,7,9];
subplot(2,1,1);
B = zeros(max(A),1);
B(A) = 1;
bar(B);

dx = 0.5;
X = zeros( round( max(A) / dx) - 1, 1);
X( round(A / dx)) = 1;
t = 1:length(X);
Y = fft(X);
P2 = abs(Y/length(X));
L = length(X);
P1 = P2(1:L/2+1);
P1(2:end-1) = 2*P1(2:end-1);
f = 0:L/2;
subplot(2,1,2);
plot(f,P1,'o-');
axis([0 10 0 0.5]);
  • Вопрос задан
  • 944 просмотра
Комментировать
Подписаться 2 Оценить Комментировать
Решения вопроса 1
longclaps
@longclaps
> Как правильно получить частотный спектр этих данных?
Спектр всегда частотный, а масло - масляное
> Ожидаю, что должен вылезти пик на амплитуде "2".
Вот это вообще непонятно что означает

Смотри: у тебя есть какие-то дискретные события, происходящие в произвольные моменты времени, например, ты сидишь под яблоней и, когда падает яблоко - засекаешь время. И ты хочешь что-то понять про периодичность, связанную с этими событиями.
Первое - временные рамки. Помимо времени падения первого и последнего яблока понадобятся моменты начала и конца наблюдений. Пусть это будут 00:00 и 24:00.
Второе. Дискретное преобразование фурье работает с периодическими отсчетами, например, сколько яблок упало за каждый час в течение суток (те 24 отсчета). БПФ - алгоритм, требующий, чтобы число отсчетов было степенью двойки (в моём примере возьму 32 отсчета по 45 мин).
Третье. Если яблок падает мало, а вытянуть из опыта хочется по-максимуму, можно попробовать интерполяцию, например, если в 12:30 упало 1 яблоко, это можно интерпретировать так: в 12:00 упала 1/3 яблока и в 12:45 - остальные 2/3. Для корректной интерполяции на краях временного диапазона неплохо бы прихватить по 45мин наблюдений до и после.
Ответ написан
Комментировать
Комментировать
Пригласить эксперта
Ответы на вопрос 2
@nirvimel
Не подскажу про Matlab, но на Python это делается так:
import numpy
import scipy.fftpack
x = numpy.array([0, 1, 0, 1, 0, 1, 0, 1, 0, 1])
numpy.abs(scipy.fftpack.fft(x)[:(len(x)/2)]) * (2./len(x))

На выходе: array([ 1., 0., 0., 0., 0.]) - тут виден единственный пик на самой высокой частоте n все более низкие частоты по нулям.
Если взять другой пример: 
x = numpy.array([0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 1])
, то получим array([ 0.5, 0. , 0. , 0.5, 0. , 0. ]) - тут два равных пика на n и на n/4, последний очевидно следует из данных, первый из того факта, что фигура исходного сигнала далека от синусоиды.
Ответ написан
Комментировать
Комментировать
begemot_sun
@begemot_sun
Программист в душе.
Я чего-то давно не практиковал Matlab. Лет 10 :) Все забыл.
Но скажу так -- преобразование фурье должно давать вам комплексное число, соотвественно пик вам нужно наблюдать на АЧХ (модуль комплексного числа).

Потом ваш сигнал - это меандр.
Может это поможет: https://habrahabr.ru/post/219337/
Ответ написан
Комментировать
Комментировать
Ваш ответ на вопрос

Войдите, чтобы написать ответ

Войти через центр авторизации
Похожие вопросы
Вакансии с Хабр Карьеры
Преподаватель по олимпиадному программированию
CODDY • Москва
от 20 000 до 40 000 ₽
Программист-разработчик систем САУ
Альбатрос • Москва
До 200 000 ₽
Программист С++
KeenTools • Ереван
от 150 000 ₽
Ещё вакансии
Заказы с Хабр Фриланса
Настройка self parser DLE
16 апр. 2024, в 12:48
1500 руб./за проект
Дописать софт на Python
16 апр. 2024, в 12:30
12000 руб./за проект
Доработать сайт на тильде
16 апр. 2024, в 12:13
10000 руб./за проект
Ещё заказы