語音學習筆記2------matlab實現傅立葉變換

Matlab是一個在很多科學和工程領域都非常有用的數學工具。傅立葉變換在訊號處理、物理、通訊、地質學、天文學、光學等很多領域都有應用。這個技術將一個函式或是一組資料從時域或是取樣域變換到頻域。這意味著，傅立葉變換可以展示一組時間序列資料的頻率分量。離散傅立葉變換是將取樣域的離散資料轉化到頻域。快速傅立葉變換是一種高效進行離散傅立葉變換的方法，並且存在很多種方法來完成快速傅立葉變換。Matlab 使用快速傅立葉變換來得到離散資料的頻域分量。下面是一個在 Matlab 中如何用快速傅立葉變換來分析音訊檔案的例子。這個例子中的檔案是記錄在 A5上的音叉錄音。這個展示了傅立葉變換如何進行和如何在 Matlab 中使用這項技術。

function test_20161214
[y,Fs]=audioread('C:\Users\wxq\Desktop\66666.wav')

Nsamps = length(y);
t = (1/Fs)*(1:Nsamps)          %Prepare time data for plot

%Do Fourier Transform
y_fft = abs(fft(y));            %Retain Magnitude
y_fft = y_fft(1:Nsamps/2);      %Discard Half of Points
f = Fs*(0:Nsamps/2-1)/Nsamps;   %Prepare freq data for plot

%Plot Sound File in Time Domain
figure
plot(t, y)
xlabel('Time (s)')
ylabel('Amplitude')
title('Tuning Fork A4 in Time Domain')

%Plot Sound File in Frequency Domain
figure
plot(f, y_fft)
xlim([0 1000])
xlabel('Frequency (Hz)')
ylabel('Amplitude')
title('Frequency Response of Tuning Fork A4')

快速傅立葉變換用 “fft” 函式執行。Matlab 沒有 “dft” 函式因為快速傅立葉變換實際上就是計算的離散傅立葉變換。儘管快速傅立葉變換的角度在很多應用中非常中央，但是隻有快速傅立葉變換的大小被儲存了。“fft” 函式允許指定一個快速傅立葉變換的輸出點數，但是在這個例子中，我們使用與輸入點數一樣數目的輸出點數。在下一行，快速傅立葉變換的一半資料被捨棄了。為了這個例子的目的所以這樣做了，但是在很多應用中，整個波譜都會用到（譯者注：我認為這裡捨棄一半的點是因為FFT是關於取樣頻率的一半對稱的，所以只要看一半就可以了）。接下來的一行，將會用於橫座標的資料通過使用取樣頻率和時遇取樣數量準備好了。這一步對於確定包含在音訊檔案的實際頻率是很重要的。

接下來，原始資料在時域上被畫了下來，快速傅立葉變換的資料也被畫了出來。為了展示在峰值的頻率上的更多詳情，在這個畫中，x軸被限定在了 [0,1000] 的範圍中。注意，在大約440Hz處，頻率響應有一個峰值，這個就是 “A5” 的頻率。在其他頻率也有一些很小東西，這些估計是音叉了。對於其他樂器比如吉他，在頻率響應的峰值的整數倍上都可以看見諧波。

傅立葉變換在很多不同的領域都是很有用的工具。二維的傅立葉變換也常常用在影象上。嘗試一下上面的程式碼，看看你能不能得到一樣的結果。

PS：有人問 FFT 結果的幅值問題。本質來說，FFT結果的幅值單位是什麼並不重要，只要你在分析過程中，需要分析的兩個幅值的單位是統一的就可以了。Matlab 的 FFT 最終結果的絕對值的確看上去並不好看（太大了），根據 Matlab 幫助檔案，FFT 的最終結果還需要進行一下小調整再來使用比較好（上面的圖未作調整），以下是我根據 Matlab 幫助文件編寫的一個計算 FFT 並繪製 FFT 結果幅值圖的函式，具體如何使用，請看函式內的註釋說明。

function [ frequency,fft_result ] = fft_plot( data,Fs,varargin )
% Calculate or plot directly fft results of data.
%
% [ frequency,fft_result ] = fft_plot( data,Fs,'plot' )
%
% inputs:
%   (1) data: data used to analysis. one row -> one data
%   (2) Fs: sample frequency
%   (3) 'plot': veriable input. if there is not this input, fft results will not be
%   ploted
% output:
%   (1) freqeuncy: frequency corresponding to the fft results
%   (2) fft_result: fft results

if nargin<2
    error('data and Fs must be given');
elseif nargin==2
    for k=1:size(data,1)
        size_data=size(data(k,:));
        if size_data(1)~=1 && size_data(2)~=1
            error('the length or the number of rows must be one.');
        end
        data(k,:)=data(k,:)-mean(data(k,:));
        L=length(data(k,:));
        NFFT=2^nextpow2(L);
        fft_result_temp=fft(data(k,:),NFFT)/length(data(k,:));
        fft_result(k,:)=fft_result_temp(k,1:NFFT/2+1);
        frequency(k,:)=Fs/2*linspace(0,1,NFFT/2+1);
    end
elseif nargin==3
    figure;
    title('FFT')
    for k=1:size(data,1)
        if strcmp(varargin,'plot')
             size_data=size(data(k,:));
            if size_data(1)~=1 && size_data(2)~=1
                error('the length or the number of rows must be one.');
            end
            data(k,:)=data(k,:)-mean(data(k,:));
            L=length(data(k,:));
            NFFT=2^nextpow2(L);
            fft_result_temp=fft(data(k,:),NFFT)/length(data(k,:));
            frequency(k,:)=Fs/2*linspace(0,1,NFFT/2+1);
            fft_result(k,:)=fft_result_temp(1:NFFT/2+1);
            subplot(size(data,1),1,k);
            plot(frequency(k,:),2*abs(fft_result(k,:)));
            xlabel('Frequency')
            ylabel('Amplitude')
        else
            error('variable input must be ''plot''');
        end
    end
elseif nargin>=3
    error('Too much inputs')
end
end

不懂的可以加我的QQ群：522869126（語音訊號處理） 歡迎你的

到來哦，看了博文給點腳印唄，謝謝啦~~