Дискретное преобразование Фурье
dsp.FFT Система object™ вычисляет дискретное преобразование Фурье (DFT) входного сигнала с использованием быстрого преобразования Фурье (FFT). Объект использует один или несколько из следующих алгоритмов быстрого преобразования Фурье (FFT) в зависимости от сложности входа и того, находится ли выход в линейном или битовом обратном порядке:
Для вычисления DFT входного сигнала:
Создать dsp.FFT и задайте его свойства.
Вызовите объект с аргументами, как если бы это была функция.
Дополнительные сведения о работе системных объектов см. в разделе Что такое системные объекты?.
ft = dsp.FFTFFT объект, ft, которая вычисляет DFT N-D массива. Для векторов столбцов или многомерных массивов: FFT объект вычисляет DFT вдоль первого размера. Если вход является вектором строки, FFT объект вычисляет строку из одного образца DFT и выдает предупреждение.
ft = dsp.FFT(Name,Value)FFT для каждого указанного свойства задано заданное значение. Заключите каждое имя свойства в отдельные кавычки. Неопределенные свойства имеют значения по умолчанию.
Чтобы использовать функцию объекта, укажите объект System в качестве первого входного аргумента. Например, для освобождения системных ресурсов объекта System с именем obj, используйте следующий синтаксис:
release(obj)
Найти частотные составляющие сигнала в аддитивном шуме.
Примечание.Этот пример выполняется только в R2016b или более поздних версиях. При использовании более ранней версии замените каждый вызов функции эквивалентным step синтаксис. Например, myObject (x) становится шагом (myObject, x).
Fs = 800; L = 1000; t = (0:L-1)'/Fs; x = sin(2*pi*250*t) + 0.75*cos(2*pi*340*t); y = x + .5*randn(size(x)); % noisy signal ft = dsp.FFT('FFTLengthSource','Property', ... 'FFTLength',1024); Y = ft(y);
Постройте график одностороннего амплитудного спектра
plot(Fs/2*linspace(0,1,512), 2*abs(Y(1:512)/1024)) title('Single-sided amplitude spectrum of noisy signal y(t)') xlabel('Frequency (Hz)'); ylabel('|Y(f)|')
Вычислить БПФ шумного синусоидального входного сигнала. Энергия сигнала сохраняется как квадрат величины коэффициентов БПФ. Определите коэффициенты БПФ, которые занимают 99,99% энергии сигнала, и восстановите сигнал временной области, взяв IFFT этих коэффициентов. Сравните восстановленный сигнал с исходным сигналом.
Note: При использовании R2016a или более ранней версии замените каждый вызов объекта эквивалентным step синтаксис. Например, obj(x) становится step(obj(x)).
Рассмотрим сигнал временной области ], который определяется в течение конечного временного интервала 0≤n≤N-1. Энергия сигнала n] задается следующим уравнением:
Коэффициенты БПФ, ], считаются значениями сигнала в частотной области. Поэтому энергия сигнала n] в частотной области является суммой квадратов величины коэффициентов БПФ:
Согласно теореме Парсеваля, общая энергия сигнала во времени или частотной области одинакова.
Инициализация
Инициализация dsp.SineWave Система предназначена для генерации синусоидальной волны, дискретизированной на частоте 44,1 кГц и имеющей частоту 1000 Гц. Построить dsp.FFT и dsp.IFFT объекты для вычисления БПФ и ОБПФ входного сигнала.
'FFTLengthSource' свойство каждого из этих объектов преобразования имеет значение 'Auto'. Следовательно, длина БПФ рассматривается как размер входного кадра. Размер входного кадра в этом примере равен 1020, что не является мощностью 2, поэтому выберите 'FFTImplementation' как 'FFTW'.
L = 1020; Sineobject = dsp.SineWave('SamplesPerFrame',L,'PhaseOffset',10,... 'SampleRate',44100,'Frequency',1000); ft = dsp.FFT('FFTImplementation','FFTW'); ift = dsp.IFFT('FFTImplementation','FFTW','ConjugateSymmetricInput',true); rng(1);
Вытекание
Поток в шумном входном сигнале. Вычислите БПФ каждого кадра и определите коэффициенты, составляющие 99,99% энергии сигнала. Возьмите IFFT этих коэффициентов, чтобы восстановить сигнал временной области.
numIter = 1000; for Iter = 1:numIter Sinewave1 = Sineobject(); Input = Sinewave1 + 0.01*randn(size(Sinewave1)); FFTCoeff = ft(Input); FFTCoeffMagSq = abs(FFTCoeff).^2; EnergyFreqDomain = (1/L)*sum(FFTCoeffMagSq); [FFTCoeffSorted, ind] = sort(((1/L)*FFTCoeffMagSq),1,'descend'); CumFFTCoeffs = cumsum(FFTCoeffSorted); EnergyPercent = (CumFFTCoeffs/EnergyFreqDomain)*100; Vec = find(EnergyPercent > 99.99); FFTCoeffsModified = zeros(L,1); FFTCoeffsModified(ind(1:Vec(1))) = FFTCoeff(ind(1:Vec(1))); ReconstrSignal = ift(FFTCoeffsModified); end
99,99% энергии сигнала может быть представлено количеством коэффициентов БПФ, заданных Vec(1):
Vec(1)
ans = 296
Сигнал восстанавливается эффективно с использованием этих коэффициентов. Если сравнить последний кадр восстановленного сигнала с исходным сигналом во временной области, можно увидеть, что разница очень мала, и графики близко совпадают.
max(abs(Input-ReconstrSignal))
ans = 0.0431
plot(Input,'*'); hold on; plot(ReconstrSignal,'o'); hold off;
Этот объект реализует алгоритм, входы и выходы, описанные на странице ссылок блока БПФ. Свойства объекта соответствуют параметрам блока.
[1] FFTW (http://www.fftw.org)
[2] Фриго, М. и С. Г. Джонсон, «FFTW: адаптивная архитектура программного обеспечения для FFT», Труды Международной конференции по акустике, речи и обработке сигналов, том 3, 1998, стр. 1381-1384.