Преобразование Фурье в Matlab - Matlab - Страница 2

@vital792 · 29.10.2013, 16:07

Студворк — интернет-сервис помощи студентам

Сообщение от R2D2

это не разложение в ряд Фурье?

я не такой серьезный теоретик, чтобы правильно сформулировать разницу между "преобразованием" фурье и "разложением в ряд", тригонометрический разумеется, т.к. я писал ранее, таких разложений много Обобщенный ряд фурье. Я бы сказал, что это одно и то же. И все таки попытаюсь сформулировать

Математически разложение в тригонометрический ряд - это представление функции в виде суммы произведений коэффициентов разложения на сами гармоники. А преобразование фурье - операция, которая имеет результатом эти же самые коэффициенты разложения.
Ряд: функция = сумма( коэффициент * гармоника)
Преобразование: коэффициенты = f(функция)
Так что по поводу

Сообщение от R2D2

разложение сигнала на гармонические составляющие

оба определения тут уместны

@Tatyana Avs · 29.10.2013, 17:19 **[ТС]**

Спасибо всем, еще раз

Зосима,если я правильно понимаю, таким же образом выделить эти гармоники? Только как это можно сделать, аппроксимировать зависимость амплитуды от частоты?
Я себе представляю график с двумя синусоидами

Зосима · 29.10.2013, 18:11

Сообщение от Tatyana Avs

Только как это можно сделать, аппроксимировать зависимость амплитуды от частоты?

Танюша, но у тебя уже как бы есть эта зависимость!

или ты хочешь узнать на каких именно частотах имеются пики спектра? для этого есть функция findpeaks:

Matlab M
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
clear, clc
 
A = load('1Signal.txt'); % считываем данные
x = A(:,1); % первый столбец 
y = A(:,2); % второй столбец
figure(1)
subplot(2,1,1);% Выбор области окна для построения
plot(x,y) % рисуем
Tm=1200;% Длина сигнала (с)
Fd=1/4;% Частота дискретизации (Гц)
FftL=2048;% Количество линий Фурье спектра
 
% Спектральное представление сигнала
FftS=abs(fft(A(:,2),FftL));% Амплитуды преобразования Фурье сигнала
FftS=2*FftS./FftL;% Нормировка спектра по амплитуде
FftS(1)=FftS(1)/2;% Нормировка постоянной составляющей в спектре
 
F=0:Fd/FftL:Fd/2-Fd/FftL;% Массив частот вычисляемого спектра Фурье
 
subplot(2,1,2);% Выбор области окна для построения
Y = FftS(1:length(F));
[p i] = findpeaks(Y,'minpeakheight',4); % находим пики большие 4
plot(F,Y,'b', F(i),p,'or'),grid;% Построение спектра Фурье сигнала
title('Signal spectra');% Подпись графика
xlabel('Frequency (Hz)');% Подпись оси х графика
ylabel('Ampl');% Подпись оси у графика
disp('Частота  Амплитуда')
disp([F(i)' p'])

@Tatyana Avs · 30.10.2013, 11:05 **[ТС]**

Хммм, трудно объяснить, когда сама плохо понимаю)))
Значит так, мой исходный сигнал - это сумма двух синусоид (этих самых гармоник частоты которых нам пока неизвестны),
мы сделали преобразование Фурье и получили частоты этих гармоник, так?
Теперь надо построить сумму этих гармоник, на определенных нами частотах (тех самых пиках), я предполагаю, что рисунок что-то подобное должен представлять

Зосима · 30.10.2013, 16:41

Ну, получить сумму не очень сложно

сложно добиться, чтобы она напоминала входной сигнал

Matlab M
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
clear, clc
 
A = load('1Signal.txt'); % считываем данные
x = A(:,1); % первый столбец 
y = A(:,2); % второй столбец
figure(1)
subplot(2,1,1);% Выбор области окна для построения
plot(x,y) % рисуем
Tm=1200;% Длина сигнала (с)
Fd=1/4;% Частота дискретизации (Гц)
FftL=2048;% Количество линий Фурье спектра
 
% Спектральное представление сигнала
FftS=abs(fft(A(:,2),FftL));% Амплитуды преобразования Фурье сигнала
FftS=2*FftS./FftL;% Нормировка спектра по амплитуде
FftS(1)=FftS(1)/2;% Нормировка постоянной составляющей в спектре
 
F=0:Fd/FftL:Fd/2-Fd/FftL;% Массив частот вычисляемого спектра Фурье
 
subplot(2,1,2);% Выбор области окна для построения
Y = FftS(1:length(F));
[p i] = findpeaks(Y,'minpeakheight',2); % находим пики большие 4
plot(F,Y,'b', F(i),p,'or'),grid;% Построение спектра Фурье сигнала
title('Signal spectra');% Подпись графика
xlabel('Frequency (Hz)');% Подпись оси х графика
ylabel('Ampl');% Подпись оси у графика
disp('Часота  Амплитуда')
disp([F(i)' p'])
 
figure
Z = p*sin( 2*pi*F(i)'*x' ); % сумма гармоник (матричная магия)
plot(x,y, x,Z,'r')

@R2D2 · 30.10.2013, 17:42

Зосима, друг мой)) А как же фазы? А нормировка амплитуд?
ИМХО должно быть что то вроде этого:

Кликните здесь для просмотра всего текста

Matlab M
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
clear; clc; close all;
data = dlmread('1Signal.txt', '\t');
t = data(:,1);
y = data(:,2);
fs = 1/(t(2)-t(1));
L = length(data);
nfft = 2^nextpow2(L);
if L~=nfft
    data(end:nfft,1) = data(end,1):1/fs:(nfft-1)/fs;
    data(end:nfft,2) = 0;
    t = data(:,1); y = data(:,2);
end
%% Spectrum
S = fft(y, nfft);
freq = fs/2*linspace(-1, 1, nfft);
[A_harms, ind_harms] = findpeaks(abs(S/nfft), 'MINPEAKHEIGHT', 10,...
    'MINPEAKDISTANCE', 3);
mag_harms = zeros(size(S));
angle_harms = zeros(size(S));
mag_harms(ind_harms) = A_harms;
angle_harms(ind_harms) = angle(S(ind_harms));
S_harms = mag_harms.*exp(j*angle_harms);
%% plots
figure('units', 'normalized', 'position', [.3 .05 .4 .85])
% Сигналы во времени
subplot(311)
plot(t, y); grid on; hold all;
plot(t, 2*nfft*ifft(S_harms))
% Спектр амплитуд
subplot(312)
plot(freq, fftshift(abs(S/nfft))); grid on; hold all;
stem(freq, fftshift(mag_harms), 'or')
xlim([-fs/2, fs/2])
% Спектр фаз
subplot(313)
plot(freq, fftshift(angle(S))); grid on; hold all;
plot(freq, fftshift(angle_harms), 'or')
xlim([-fs/2, fs/2])

П.С. Юзайте тот файл, что во вложении. У ТС в исходном файле имеется небольшой баг. Надо будет - скажу какой

Зосима · 30.10.2013, 18:30

Не по теме:

R2D2, я работаю по методу простых итераций - даю неверное решение, его исправляют и получаем то что надо! :D

@R2D2 · 30.10.2013, 19:05

Не по теме:

Сообщение от Зосима

получаем то что надо

Надеюсь... Ждем отклик... :)

@Tatyana Avs · 04.11.2013, 17:36 **[ТС]**

Вновь ОГРОМНОЕ! спасибо!
Если я верно поняла, то мы получили зависимости 1) амплитуды (частоты)
2) Амплитуды (фазы)

Matlab M
1
2
3
4
5
[A_harms, ind_harms] = findpeaks(abs(S/nfft), 'MINPEAKHEIGHT', 10,...
'MINPEAKDISTANCE', 3);
mag_harms = zeros(size(S));
angle_harms = zeros(size(S));
mag_harms(ind_harms) = A_harms;

Если не сложно, не могли бы немного прокомментировать этот момент?
И что за график рисуется зеленым цветом?) Это и есть искомая сумма гармоник?

@R2D2 · 05.11.2013, 15:11

Tatyana Avs, спектр сигнала в общем случае есть комплексной функцией частоты. Есть такие хитрые сигналы у которых мнимая часть спекта = 0. Но в Вашем случае сигнал задан действительный и, соответственно спектр его будет комплексным.
И как у любой, уважающей себя, комплексной функции у него есть действительная и мнимая часть (или же модуль и аргумент).
В моих кругах принято второй вариант представления. Модуль (рис.2) - это зависимость амплитуд гармоник от частоты, а аргумент (рис.3) - зависимость фаз этих же гармоник от частоты.
И, да, на рис.1 зеленым - это сумма гармоник.
Про фрагмент кода - чуть позже отпишу.

Добавлено через 2 часа 13 минут

Matlab M
1
2
[A_harms, ind_harms] = findpeaks(abs(S/nfft), 'MINPEAKHEIGHT', 10,...
'MINPEAKDISTANCE', 3);

В этой строчке я нахожу пики в векторе abs(S/nfft), пареметры поиска: пороговое значение амплитуды = 10, минимальное расстояние между пиками (в единицах индексов вектора) 3 индекса (то есть пики, между которыми менее 3х значений игнорируются).

Matlab M
1
2
mag_harms = zeros(size(S));
angle_harms = zeros(size(S));

Здесь я задаю нулевые вектора для амплитуд и фаз гармоник. Эти вектора той же длины, что и вычисленный ранее спектр сигнала, это нужно для корректного отображения в плоте.

Matlab M
1
mag_harms(ind_harms) = A_harms;

Теперь я вписываю в вектор амплитуд гармоник полученные выше значения пиков.
Т.е. ненулевыми элементами вектора амплитуд гармоник стали только те, которые соответствуют пикам в спектре.

@Tatyana Avs · 13.11.2013, 14:01 **[ТС]**

Все работает, спасибо)
Как и было сказано Зосимой, надо добиваться соответствия экспериментальному сигналу. Я попробовала учитывать до 6 гармоник, выбирая частоты. Не знаете, при этом учитываются фазы? Или их надо тоже задавать и варьировать?

Зосима · 13.11.2013, 16:56

Tatyana Avs, а как же без фаз-то?

нужны! и они не задаются и не варьируются, а точно рассчитываются, как и амплитуды гармоник

@R2D2 · 13.11.2013, 19:15

Не по теме:

Сообщение от Tatyana Avs

Или их надо тоже задавать и варьировать?

"Не отвлекайся на босяцкий подгон! Только РАСЧЕТЫ! Только хардкор!"

Из спектра амплитуд ты берешь АМПЛИТУДУ гармоники (в пике допустим), а из спектра фаз - берешь соответствующую этой гармонике ФАЗУ!!!

@Tatyana Avs · 13.11.2013, 20:34 **[ТС]**

Да это я конечно поняла, спасибо)
Но при построении гармоник в написанной программе фазы уже тоже учтены? Или построение идет только с учетом амплитуд? Ну вроде как синусоида задается например y=a + bsin (cx+d),b-амплитуда, d - фаза. А тут не могу понять)
Насколько я поняла, вот в этом моменте считается сумма, а angle_harms это и есть фазы?

Matlab M
1
S_harms = mag_harms.*exp(j*angle_harms);

@R2D2 · 13.11.2013, 20:45

В этом моменте S_harms стала вектором комплексных чисел: модули этих чисел - это амплитуды синусоид, аргументы - фазы синусоид.

@Tatyana Avs · 22.11.2013, 18:50 **[ТС]**

Извините, что надоедаю снова.
Можно спросить почему используем функцию exp, а не sin, когда рассчитываем сумму гармоник? что такое j?
и как получить исходный массив частот (т.е. данные оси x)?

@R2D2 · 23.11.2013, 21:36

Как уже было упомянуто, спектр - комплексая величина. В той строке, о которой идет речь, Вы видите вариант записи комплексного числа (массива чисел) в показательной форме. j - корень из -1 (мнимая единица). Что касается массива частот - не понял вопроса. По моему я в программе этот массив задал. Так что, не понятно в чем, собственно, вопрос.

@Tatyana Avs · 23.11.2013, 23:21 **[ТС]**

С этим разобралась, извините, уже поняла. Только почему использовали функцию exp, а не sin, где вычисляется сумма синусоид?
Или если считаем комплексные числа, используем exp, можно ли использовать sin в таком случае?
S_harms = mag_harms.*exp(j*angle_harms);

@R2D2 · 24.11.2013, 00:03

Посмотрите формы представления к.ч..

Сообщение от Tatyana Avs

можно ли использовать sin в таком случае?

Можно. Если использовать тригонометрическую форму записи.

UPD от Зосима: будет S_harms = mag_harms.*( cos(angle_harms) + j*sin(angle_harms) );

@Gidrobarbulyato · 11.05.2022, 16:49

Здравствуйте. Хотелось бы оживить тему. У меня есть сигнал с осциллографа. Там импульсное напряжение. Мне хотелось бы разложить этот сигнал по гармоникам. А далее представить это в тригонометрическом виде (как сумму синусов и косинусов). Зачем мне это нужно - программа, в которую я хочу загнать этот сигнал, понимает только тупую запись типа A1*sin(w1*t)+... Поэтому думаю сделать такое извращение, но не представляю как. Понимаю только, что на выходе у меня должна быть таблица амплитуд и таблица гармоник по номерам, а там я уж сам ручками запишу в тригонометрической форме. Думаю, записать будет достаточно первые 10-15 гармоник, мне 99,9999% точности не нужно. Сигнал прилагаю.

Новые блоги и статьи Все статьи Все блоги /
Почему SEO не начинается с ключевых слов: что проверить до написания текстов Neotwalker 01.08.2026 Когда владельцу сайта предлагают заняться SEO, первым шагом часто становится сбор запросов и написание текстов. Логика кажется понятной: 1. Находим ключевые слова. 2. Добавляем их на. . .	Знание — сила: Доктрина интенциональности знаний, углубление в формулу Hrethgir 01.08.2026 https:/ / www. cyberforum. ru/ blog_attachment. php?attachmentid=11957&stc=1&d=1785567302 Знаменитый афоризм Фрэнсиса Бэкона «Знание — сила» (Scientia potentia est) в массовой культуре принято понимать. . .	SUNO Ai - Река Без Дна zorxor 31.07.2026 Автор стихотворения - астрофизик Марина Катыс Ссылка на сгенерированную музыкальную композицию: https:/ / suno. com/ song/ 6f6e5464-b290-4650-be6c-44c85f8d8013 Я говорю, что Время- как вода течет. . .	Из невошедшего на форум (диалог с ИИ-гугла) zorxor 29.07.2026 А вот, что интересно, сказал мне ИИ-гугла: Этот текст — эмоциональный пост пользователя под ником zorxor на интернет-форуме (вероятно, посвященном мистике, непознанному или альтернативной науке). . . .
Был праздник вчера, а я и не знал. kumehtar 28.07.2026 27. 07. 2026г. Intel Core 2 Duo исполнилось 20 лет Салют, шампанское, овации! :drink:	Нейтральные знания, чистый код - бла-бла-бла-бла, на самом деле кликбейт и самореклама, плагиат, и вот почему Hrethgir 27.07.2026 То-есть отклонение такой публикации говорит само за себя, и пусть только возьмут на вооружение после отклонения публикации - это будет чистейшим актом плагиата. Отклонял Хабр. Дословно, отклонённая. . .	тв 16 бой ии anaschu 27.07.2026 Великий Перелом ИИ: Как уравнения ОДУ Radau дожали цензурные фильтры Алисы Фиксируем в мемофонде Теории Всего беспрецедентный факт в истории ИИ-зондирования. В затяжном многораундовом. . .	мв 15. непроверенное, возможно, глюк anaschu 27.07.2026 НАУЧНО-АНАЛИТИЧЕСКИЙ ОТЧЕТ. РАЗДЕЛ 1. 1: «НАУКА» (РАСШИРЕННАЯ СТЕХИОМЕТРИЧЕСКАЯ И ГЕНЕТИЧЕСКАЯ ВЕРСИЯ)Тема: Теоретическое обоснование инвариантности 19-мерного тензорного ядра непрерывных ОДУ и. . .

@Tatyana Avs 0 / 0 / 0 Регистрация: 24.10.2013 Сообщений: 22
	29.10.2013, 17:19 [ТС]
	Спасибо всем, еще раз Зосима,если я правильно понимаю, таким же образом выделить эти гармоники? Только как это можно сделать, аппроксимировать зависимость амплитуды от частоты? Я себе представляю график с двумя синусоидами 0

@Tatyana Avs 0 / 0 / 0 Регистрация: 24.10.2013 Сообщений: 22
	30.10.2013, 11:05 [ТС]
	Хммм, трудно объяснить, когда сама плохо понимаю))) Значит так, мой исходный сигнал - это сумма двух синусоид (этих самых гармоник частоты которых нам пока неизвестны), мы сделали преобразование Фурье и получили частоты этих гармоник, так? Теперь надо построить сумму этих гармоник, на определенных нами частотах (тех самых пиках), я предполагаю, что рисунок что-то подобное должен представлять 0

Зосима
	30.10.2013, 18:30
	Не по теме: R2D2, я работаю по методу простых итераций - даю неверное решение, его исправляют и получаем то что надо! :D 0

@Tatyana Avs 0 / 0 / 0 Регистрация: 24.10.2013 Сообщений: 22
	13.11.2013, 14:01 [ТС]
	Все работает, спасибо) Как и было сказано Зосимой, надо добиваться соответствия экспериментальному сигналу. Я попробовала учитывать до 6 гармоник, выбирая частоты. Не знаете, при этом учитываются фазы? Или их надо тоже задавать и варьировать? 0

Зосима 5246 / 3574 / 379 Регистрация: 02.04.2012 Сообщений: 6,477 Записей в блоге: 18
	13.11.2013, 16:56
	Tatyana Avs, а как же без фаз-то? нужны! и они не задаются и не варьируются, а точно рассчитываются, как и амплитуды гармоник 1

@R2D2 939 / 839 / 121 Регистрация: 23.11.2012 Сообщений: 2,488
	13.11.2013, 20:45
	В этом моменте S_harms стала вектором комплексных чисел: модули этих чисел - это амплитуды синусоид, аргументы - фазы синусоид. 1

@Tatyana Avs 0 / 0 / 0 Регистрация: 24.10.2013 Сообщений: 22
	22.11.2013, 18:50 [ТС]
	Извините, что надоедаю снова. Можно спросить почему используем функцию exp, а не sin, когда рассчитываем сумму гармоник? что такое j? и как получить исходный массив частот (т.е. данные оси x)? 0

@R2D2 939 / 839 / 121 Регистрация: 23.11.2012 Сообщений: 2,488
	23.11.2013, 21:36
	Как уже было упомянуто, спектр - комплексая величина. В той строке, о которой идет речь, Вы видите вариант записи комплексного числа (массива чисел) в показательной форме. j - корень из -1 (мнимая единица). Что касается массива частот - не понял вопроса. По моему я в программе этот массив задал. Так что, не понятно в чем, собственно, вопрос. 0

@Tatyana Avs 0 / 0 / 0 Регистрация: 24.10.2013 Сообщений: 22
	23.11.2013, 23:21 [ТС]
	С этим разобралась, извините, уже поняла. Только почему использовали функцию exp, а не sin, где вычисляется сумма синусоид? Или если считаем комплексные числа, используем exp, можно ли использовать sin в таком случае? S_harms = mag_harms.exp(jangle_harms); 0