一. 簡介
離散傅立葉、離散余弦和離散小波變換是圖像、音頻信號常用基礎操作,時域信號轉換到不同變換域以后,會導致不同程度的能量集中,信息隱藏利用這個原理在變換域選擇適當位置系數進行修改,嵌入信息,並確保圖像、音頻信號經處理后感官質量無明顯變化。
二. 數學公式
一維離散傅立葉變換對定義
一維離散傅里葉變換:
一維離散傅里葉逆變換:
一維離散余弦變換對定義
一維離散余弦正變換:
一維離散余弦反變換:
一維連續小波變換對定義
一維連續小波變換,其總h(t)是小波母函數:
一維連續小波逆變換:
二維離散傅立葉變換對定義
二維離散傅里葉變換:
二維離散傅里葉逆變換:
二維離散余弦變換對定義
二維離散余弦正變換:
二維離散余弦反變換:
三. 代碼實現
1. 用離散傅里葉變換分析合成音頻和圖像
(1)分析合成音頻文件
第一步:讀取音頻文件數據。
uigetfile() 是文件對話框函數,提供圖形界面供用戶選擇所需文件,返回目標的目錄名和文件名。
函數原型:y= wavread (FILE)
功能:讀取微軟音頻格式(wav)文件內容
輸入參數:file 表示音頻文件名,字符串
返回參數:y 表示音頻樣點,浮點型
第二步:一維離散傅立葉變換
fft 函數對輸入參數進行一維離散傅立葉變換並返回其系數,對應頻率從 0 到 fs(采樣頻率),使用 fftshift 將零頻對應系數移至中央。為了更好地觀察頻譜,可計算離散樣點對應的頻率值。
第三步:一維離散傅立葉逆變換
ifft 函數對輸入參數進行一維離散傅立葉逆變換並返回其系數。
第四步:觀察結果
figure(n)表示創建第 n 個圖形窗。
subplot 是子繪圖函數,第一、二個參數指明子圖像布局方式,例如,若參數為 2,3 則表示畫面共分為 2 行,每行有 3 個子圖像。第三個參數表明子圖像序號,排序順序為從左至右,從上至下。
plot 是繪圖函數,默認使用方式為 plot(y),參數 y 是要繪制的數據;如果需要指明圖像橫軸顯示序列,則命令行為 plot(x, y),默認方式等同於 plot ([0..len-1], y),len為序列y的長度。
例
len=40000; [fn,pn]=uigetfile('*.wav','請選擇音頻文件'); [x,fs]=wavread(strcat(pn,fn),len); %strcat()連接生成文件絕對路徑 %一維離散傅里葉變換 xf=fft(x); f1=[0:len-1]*fs/len; xff=fftshift(xf); % fftshift()將零頻對應系數移至中央 h1=floor(len/2); %取不大於len/2的最大整數 f2=[-h1:h1]*fs/len; %一維離散傅里葉逆變換 xsync=ifft(xf); figure; subplot(2,2,1);plot(x);title('原始圖像'); subplot(2,2,2);plot(xsync);title('synthesiaze audio'); subplot(2,2,3);plot(f1,abs(xf));title('fft coef. of audio'); %abs()求整數的絕對值 subplot(2,2,4);plot(f2(1:len),abs(xff));title('fftshift coef. of audio');
(2)分析合成圖像文件
第一步:讀取圖像文件數據
函數原型:A = imread(filename,fmt)
功能:讀取 fmt 指定格式的圖像文件內容
輸入參數:filename 表示圖像文件名,字符串。Fmt 表示圖像文件格式名,字符串,函數支持的圖像格式包括:JPEG,TIFF,GIF,BMP 等等,當參數中不包括文件格式名時,函數嘗試推斷出文件格式。
返回參數:A 表示圖像數據內容,整型。
rgb2gray 函數將 RGB 圖像轉換為灰度圖。
第二步:二維離散傅立葉變換
fft2 函數對輸入參數進行二維離散傅立葉變換並返回其系數,使用 fftshift 將零頻對應系數移至中央。
第三步:二維離散傅立葉逆變換
ifft2 函數對輸入參數進行二維離散傅立葉逆變換並返回其系數。
第四步:觀察結果
imshow 是二維數據繪圖函數,mesh 通過三維平面顯示數據。
例:
[fn,pn]=uigetfile('*.bmp','請選擇圖像文件'); [x,map]=imread(strcat(pn,fn),'bmp'); I=rgb2gray(x); %rgb2gray函數將 RGB 圖像轉換為灰度圖 %二維離散傅立葉變換 xf=fft2(I); xff=fftshift(xf); %二維離散傅立葉逆變換 xsync=ifft2(xf); %結果 figure; subplot(2,2,1);imshow(x);title('original image'); subplot(2,2,2);imshow(uint8(abs(xsync)));title('synthesize image'); subplot(2,2,3);mesh(abs(xf));title('fft coef. of image'); subplot(2,2,4);mesh(abs(xff));title('fftshift coef. of image');
2. 用離散余弦變換分析合成音頻和圖像
(1)分析合成音頻文件
第一步:讀取音頻文件數據。
第二步:一維離散余弦變換
dct 函數對輸入參數進行一維離散余弦變換並返回其系數,對應頻率從0到 fs(采樣頻率)。
第三步:一維離散余弦逆變換
idct 函數對輸入參數進行一維離散余弦逆變換並返回其系數。離散余弦變換常用於圖像壓縮,可以嘗試只使用部分系數重構語言,通過觀察可發現,原始音頻和合成后音頻兩者差別不大。
第四步:觀察結果
例:
len=40000; [fn,pn]=uigetfile('*.wav','請選擇音頻文件'); [x,fs]=wavread(strcat(pn,fn),len); %一維離散余弦變換 xf=dct(x); f1=[0:len-1]*fs/len; h1=floor(len/2); f2=[-h1:h1]*fs/len; %一維離散余弦逆變換 xsync=idct(xf); [row,col]=size(x); xff=zeros(row,col); xff(1:row,1:col)=xf(1:row,1:col); y=idct(xff); %結果 figure; subplot(2,2,1);plot(x);title('原始圖像'); subplot(2,2,2);plot(xsync);title('synthesiaze audio'); subplot(2,2,3);plot(f1,abs(xf));title('fft coef. of audio'); subplot(2,2,4);plot(f2(1:len),abs(xff));title('fftshift coef. of audio');
(2)分析合成圖像文件
第一步:讀取圖像文件數據
第二步:二維離散余弦變換
dct 函數對輸入參數進行二維離散余弦變換並返回其系數。
第三步:二維離散余弦逆變換
idct2 函數對輸入參數進行二維離散余弦逆變換並返回其系數。可以嘗試使用部分系數重構圖像,例如使用系數矩陣中4/5的數據,其它部分置零。
為了保證圖像能正確顯示,使用uint8 對重構圖像原始數據進行了數據類型轉換,確保其取值范圍在 0 到 255 之間。
第四步:觀察結果
請輸入命令顯示四個子圖,分別是原始圖像、使用全部系數恢復的圖像,使用部分系數恢復的圖像和用三維立體圖方式顯示系數。
例:
[fn,pn]=uigetfile('*.bmp','請選擇圖像文件'); [x,map]=imread(strcat(pn,fn),'bmp'); I=rgb2gray(x); %二維離散余弦變換 xf=dct2(I); %二維離散余弦逆變換 xsync=uint8(idct2(xf)); [row,col]=size(I); lenr=round(row*4/5); lenc=round(col*4/5); xff=zeros(row,col); xff(1:lenr,1:lenc)=xf(1:lenr,1:lenc); y=uint8(idct2(xff)); %結果 figure; subplot(2,2,1);imshow(x);title('original image'); subplot(2,2,2);imshow(uint8(abs(xsync)));title('synthesize image'); subplot(2,2,3);imshow(uint8(abs(y)));title('part synthesize image'); subplot(2,2,4);mesh(abs(xff));title('fftshift coef. of image');
3. 用離散小波變換分析合成音頻和圖像
(1)分析合成音頻文件
第一步:讀取音頻文件數據。
第二步:一維離散小波變換
wavedec 函數對輸入參數進行一維離散小波變換並返回其系數C 和各級系數長度L。第二個參數指明小波變換的級數,第三個參數指明小波變換使用的小波基名稱。
第三步:一維離散小波逆變換
waverec 函數對輸入參數進行一維離散小波逆變換並返回其系數。
appcoef 函數返回小波系數近似分量,第一個參數 C、第二個參數 L 是 wavedec 的返回參數,為各級小波系數和其長度,第三個參數指明小波基名稱,第四個參數指明級數。
detcoef 函數返回小波系數細節分量,第一個參數 C、第二個參數 L 是 wavedec 的返回參數,為各級小波系數和其長度,第三個參數指明級數。
第四步:觀察結果
例:
len=40000; [fn,pn]=uigetfile('*.wav','請選擇音頻文件'); [x,fs]=wavread(strcat(pn,fn),len); %一維離散小波變換 [C,L]=wavedec(x,2,'db4'); %一維離散小波逆變換 xsync=waverec(C,L,'db4'); cA2=appcoef(C,L,'db4',2); cD2=detcoef(C,L,2); cD1=detcoef(C,L,1); %結果 figure; subplot(2,3,1);plot(x);title('原始圖像'); subplot(2,3,2);plot(xsync);title('synthesiaze audio'); subplot(2,3,4);plot(cA2);title('app coef. of audio'); subplot(2,3,5);plot(cD2);title('det coef. of audio'); subplot(2,3,6);plot(cD1);title('det coef. of audio');
(2)分析合成圖像文件
第一步:讀取圖像文件數據
第二步:二維離散小波變換
dwt2函數對輸入參數進行二維一級離散小波變換並返回近似分量,水平細節分量,垂直細節分量和對角線細節分量。如果要對圖像進行多級小波分解,使用wavedec2函數。
第三步:二維離散小波逆變換
idwt2 函數對輸入參數進行二維離散小波逆變換並返回其系數。可以嘗試僅使用近似分量、水平細節分量、垂直細節分量或對角線細節分量重構圖像。
第四步:觀察結果
輸入命令顯示六個子圖,分別是原始圖像、使用全部系數恢復的圖像、小波系數近似分量、水平細節分量、垂直細節分量和對角線細節分量。
例:
[fn,pn]=uigetfile('*.bmp','請選擇bmp格式圖像文件'); [x,map]=imread(strcat(pn,fn),'bmp'); I=rgb2gray(x); %二維離散小波變換 sx=size(I); [cA1,cH1,cV1,cD1]=dwt2(I,'bior3.7'); %二維離散小波逆變換 xsync=uint8(idwt2(cA1,cH1,cV1,cD1,'bior3.7',sx)); A1=uint8(idwt2(cA1,[],[],[],'bior3.7',sx)); H1=uint8(idwt2([],cH1,[],[],'bior3.7',sx)); V1=uint8(idwt2([],[],cV1,[],'bior3.7',sx)); D1=uint8(idwt2([],[],[],cD1,'bior3.7',sx)); %結果 figure; subplot(2,3,1);imshow(x);title('original image'); subplot(2,3,2);imshow(uint8(abs(xsync)));title('synthesize image'); subplot(2,3,3);mesh(A1);title('app coef. of image'); subplot(2,3,4);mesh(H1);title('hor coef. of image'); subplot(2,3,5);mesh(V1);title('ver coef. of image'); subplot(2,3,6);mesh(D1);title('dia coef. of image');