Преобразование Фурье в MATLAB
21.12.2020
Подскажите, пожалуйста, относительно Фурье преобразования в среде MATLAB. Сам несколько дней пытаюсь понять через гугл и периодические журналы... Есть набор ЭЭГ-сигналов: 500 значений за 2 сек (частот...
Подскажите, пожалуйста, относительно Фурье преобразования в среде MATLAB. Сам несколько дней пытаюсь понять через гугл и периодические журналы... Есть набор ЭЭГ-сигналов: 500 значений за 2 сек (частота 250 Гц). Мне надо просто разложить их в спектр (насколько понимаю от 0 до 125) и получить значения 2-х полос: от 0.5 до 14 Гц и от 14 до 40 Гц. Но даже разложение в спектр, полагаю с ошибкой:
% initfs = 250;
% loadda = da1_all;x = da(1:256,1);
% Построение графика спектральной плотностиfx = fft(x);sfx = fftshift(fx);f = linspace(-fs/2,fs/2,length(x));area(f,abs(sfx));title('Амплитудный спектр сигнала');xlabel('Частота(Гц)');ylabel('Амплитуда');
Комментарии
Может быть это поможет
https://www.mathworks.com/matlabcentral/fileexchange/50082-eeg-data-plotting-power-spectrum-spectrogram-frequency-spectrum-of-alpha-beta-delta-and-theta?s_tid=srchtitle
Благодарю. Но там самого кода преобразования в спектр не нашёл(фигурирует только fs/2).
1. sfx = fftshift(fx); - это зачем?
2. Расскажите хоть что-то о своем сигнале (форма, гармоники, и т.д.), чтоб можно было анализировать? Или там какие-то сейсмические данные о которых вы ничего не знаете? А то, то что у вас на картинке, от этого ни жарко, ни холодно, ибо информации входной нет. Ну или код выложите, при условии если он не огромных размеров :).
Для примера (тут 8 каналов). Спектральное представление нужно для каждого индивидуально. Также здесь 2 сек, а сейчас работаю с 4 сек. временными рядами. 250 отсчётов в 1 сек. По-хорошему, целевой сигнал находится в диапазоне частот 0,5 - 40 Гц (вне диапазона шумы электросети и др...). fftshit, как понял из периодики, нужен для сдвига нулевой частоты в начало координат (т. е. результат не был симметричным).
P.S. Да, сигнал я не синтезирую, а снимаю со своей головы (т.е. кода никакого нет).
fftshift думаю не нужен.
сохраните 1 сигнал в mat файле и выложите , проверить. или эксель, или в чем это у вас.
ps: а вобще в хелпе на fft всю жизнь был толковый пример.
psps: с головы? )))))))))
С поверхности головы, точнее не скажешь)
Приложил сигнал по одному каналу (1 из 8 как с графика выше).
А где в хелпе? Пробовал по примеру с офсайта https://se.mathworks.com/help/matlab/ref/fft.html
Так про него и говорю.
Вы наблюдаете ВЧ составляющие в вашем сигнале? Вот они и вылезают.
Не ясно следующее:
1. Вы прислали самый не показательный сигнал (нижнюю кривую) - она самая не изменяющаяся.
2. Увас на рисунке по оси ОХ написано - отсчеты. При этом максимальное значение 2е4. Т.е. 20000 - отсчетов. А в файле, что выложили, только 500.
Сигналы с графика и с mat-файла, действительно разные. Они сделаныв разное время, но при одинаковых условиях. Для графика нужна была показательность, чего нельзя сделать на коротком интервале. Задача перевода в спектр встала несколько позднее, тогда начал деалть отсчётыпо 2 сек. А прислал я не нижнюю кривую, а кривую 1 (наоборот, то ли верхнюю, то ли из центра).
Высокочастотные, да. Имеющий смысл диапазон находится до 40 Гц. Выше 50 Гц начинаются наводки электросети (поэтому собирался обрезать).
Вообще, спасибо за наглядные графики. Можете сориентировать, что скорректировать в моем коде(в идеале, конечно, ваш пример кода, по которому нижний график получился)?
Но у вас там есть "волны", а на рис.предыдущего поста - просто "прямая с шумом".
UPD: Быть не может, у вас там только одна кривая с отрицательными значениями.
C1 - из вашего файла.
Вот мне не жалко, честно, но кроме
- чтения вашего файла,
- функции Y = fft(C1) ,
- затем вывода графика,
. . . ничего больше нет. Никаких других хитрых функций и замысловатостей, которые у вас присутствуют. Попробуйте все лишнее выкинуть. Заодно Вам же понятнее станет. Если уж никак не получится, то конечно выложу. :)
PS: Но надо конечно "фильтрануть" всё лишнее, что вам не интересно, постояннное смещение особенно, иначе на фоне вы ничего полезного не увидите. Хотя я не знаю, если вам постоянка нужна, то графики, что нарисовал, не верные.
Добрый день!
Да, пример кода всё-таки не помешает. Попробовал упростить, вот что получилось:
Про постоянную составляющую: предполагаю, она всё-таки не нужна. С фильтрацией и удалением шумов, конечно, ещё предстоит работа...
Не знаю точно, почему на некоторых каналах разность потенциалов (амплитуда) отрицательная. Не стал углубляться, т.к. счёл ключевыми параметрыми свойства волны. Ещё раз оговорюсь, сигнал, котрый отображен на верхней картинке снимался ранее (был программынй сбой и проходили только 8 каналов). Насколько помню, его протяженность 80 сек. Сигнал из mat-файла по 1 из 16 каналов, с которыми работаю сейчас. Но на коротких промежутках он действительно выглядит как прямая с шумом:
Код у вас наверняка правильный, судя по тычку на 200 Гц, там ошибиться сложно. Удалите сперва постоянное смещение из сигналов (безо всяких фильтров). Самый простой деревенский способ:
Signal = Signal / mean(Signal);
Думаю сразу что-то да прояснится в спектре.