1. 程式人生 > >傅立葉變換及其在opencv中影象去噪的實現

傅立葉變換及其在opencv中影象去噪的實現

前言

我保證這篇文章和你以前看過的所有文章都不同,這是12年還在果殼的時候寫的,但是當時沒有來得及寫 完就出國了……於是拖了兩年,嗯,我是拖延症患者……

這篇文章的核心思想就是:

要讓讀者在不看任何數學公式的情況下理解傅立葉分析。

傅立葉分析不僅僅是一個數學工具,更是一種可以徹底顛覆一個人以前世界觀的思維模式。但不幸的是,傅立葉分析的公式看起來太複雜了,所以很多大一新生上來就懵圈並從此對它深惡痛絕。老實說,這麼有意思的東西居然成了大學裡的殺手課程,不得不歸咎於編教材的人實在是太嚴肅了。(您把教材寫得好玩一點會死嗎?會死嗎?)所以我一直想寫一個有意思的文章來解釋傅立葉分析,有可能的話高中生都能看懂的那種。所以,不管讀到這裡的您從事何種工作,我保證您都能看懂,並且一定將體會到通過傅立葉分析看到世界另一個樣子時的快感。至於對於已經有一定基礎的朋友,也希望不要看到會的地方就急忙往後翻,仔細讀一定會有新的發現。

————以上是定場詩————

下面進入正題:

抱歉,還是要囉嗦一句:其實學習本來就不是易事,我寫這篇文章的初衷也是希望大家學習起來更加輕鬆,充滿樂趣。但是千萬!千萬不要把這篇文章收藏起來,或是存下地址,心裡想著:以後有時間再看。這樣的例子太多了,也許幾年後你都沒有再開啟這個頁面。無論如何,耐下心,讀下去。這篇文章要比讀課本要輕鬆、開心得多……

p.s.本文無論是cos還是sin,都統一用“正弦波”(Sine Wave)一詞來代表簡諧波。

一、什麼是頻域

從我們出生,我們看到的世界都以時間貫穿,股票的走勢、人的身高、汽車的軌跡都會隨著時間發生改變。這種以時間作為參照來觀察動態世界的方法我們稱其為時域分析。而我們也想當然的認為,世間萬物都在隨著時間不停的改變,並且永遠不會靜止下來。但如果我告訴你,用另一種方法來觀察世界的話,你會發現世界是永恆不變的,你會不會覺得我瘋了?我沒有瘋,這個靜止的世界就叫做頻域。

先舉一個公式上並非很恰當,但意義上再貼切不過的例子:

在你的理解中,一段音樂是什麼呢? 在這裡插入圖片描述 這是我們對音樂最普遍的理解,一個隨著時間變化的震動。但我相信對於樂器小能手們來說,音樂更直觀的理解是這樣的: 在這裡插入圖片描述 好的!下課,同學們再見。 是的,其實這一段寫到這裡已經可以結束了。上圖是音樂在時域的樣子,而下圖則是音樂在頻域的樣子。所以頻域這一概念對大家都從不陌生,只是從來沒意識到而已。

現在我們可以回過頭來重新看看一開始那句痴人說夢般的話:世界是永恆的。

將以上兩圖簡化:

時域: 在這裡插入圖片描述 頻域: 在這裡插入圖片描述 在時域,我們觀察到鋼琴的琴絃一會上一會下的擺動,就如同一支股票的走勢;而在頻域,只有那一個永恆的音符。

所以 你眼中看似落葉紛飛變化無常的世界,實際只是躺在上帝懷中一份早已譜好的樂章。

抱歉,這不是一句雞湯文,而是黑板上確鑿的公式:傅立葉同學告訴我們,任何周期函式,都可以看作是不同振幅,不同相位正弦波的疊加。在第一個例子裡我們可以理解為,利用對不同琴鍵不同力度,不同時間點的敲擊,可以組合出任何一首樂曲。

而貫穿時域與頻域的方法之一,就是傳中說的傅立葉分析。傅立葉分析可分為傅立葉級數(Fourier Serie)和傅立葉變換(Fourier Transformation),我們從簡單的開始談起。

二、傅立葉級數(Fourier Series)的頻譜

還是舉個栗子並且有圖有真相才好理解。

如果我說我能用前面說的正弦曲線波疊加出一個帶90度角的矩形波來,你會相信嗎?你不會,就像當年的我一樣。但是看看下圖: 在這裡插入圖片描述 第一幅圖是一個鬱悶的正弦波cos(x)

第二幅圖是2個賣萌的正弦波的疊加cos(x)+a.cos(3x)

第三幅圖是4個發春的正弦波的疊加

第四幅圖是10個便祕的正弦波的疊加

隨著正弦波數量逐漸的增長,他們最終會疊加成一個標準的矩形,大家從中體會到了什麼道理?

(只要努力,彎的都能掰直!)

隨著疊加的遞增,所有正弦波中上升的部分逐漸讓原本緩慢增加的曲線不斷變陡,而所有正弦波中下降的部分又抵消了上升到最高處時繼續上升的部分使其變為水平線。一個矩形就這麼疊加而成了。但是要多少個正弦波疊加起來才能形成一個標準90度角的矩形波呢?不幸的告訴大家,答案是無窮多個。(上帝:我能讓你們猜著我?)

不僅僅是矩形,你能想到的任何波形都是可以如此方法用正弦波疊加起來的。這是沒 有接觸過傅立葉分析的人在直覺上的第一個難點,但是一旦接受了這樣的設定,遊戲就開始有意思起來了。

還是上圖的正弦波累加成矩形波,我們換一個角度來看看: 在這裡插入圖片描述 在這幾幅圖中,最前面黑色的線就是所有正弦波疊加而成的總和,也就是越來越接近矩形波的那個圖形。而後面依不同顏色排列而成的正弦波就是組合為矩形波的各個分量。這些正弦波按照頻率從低到高從前向後排列開來,而每一個波的振幅都是不同的。一定有細心的讀者發現了,每兩個正弦波之間都還有一條直線,那並不是分割線,而是振幅為0的正弦波!也就是說,為了組成特殊的曲線,有些正弦波成分是不需要的。

這裡,不同頻率的正弦波我們成為頻率分量。

好了,關鍵的地方來了!!

如果我們把第一個頻率最低的頻率分量看作“1”,我們就有了構建頻域的最基本單元。

對於我們最常見的有理數軸,數字“1”就是有理數軸的基本單元。

時域的基本單元就是“1秒”,如果我們將一個角頻率為w0w_0 的正弦波cos(w0tw_0t)看作基礎,那麼頻域的基本單元就是w0w_0

有了“1”,還要有“0”才能構成世界,那麼頻域的“0”是什麼呢?cos(0t)就是一個週期無限長的正弦波,也就是一條直線!所以在頻域,0頻率也被稱為直流分量,在傅立葉級數的疊加中,它僅僅影響全部波形相對於數軸整體向上或是向下而不改變波的形狀。

接下來,讓我們回到初中,回憶一下已經死去的八戒,啊不,已經死去的老師是怎麼定義正弦波的吧。 在這裡插入圖片描述 正弦波就是一個圓周運動在一條直線上的投影。所以頻域的基本單元也可以理解為一個始終在旋轉的圓: 在這裡插入圖片描述 介紹完了頻域的基本組成單元,我們就可以看一看一個矩形波,在頻域裡的另一個模樣了: 在這裡插入圖片描述 這是什麼奇怪的東西? 這就是矩形波在頻域的樣子,是不是完全認不出來了?教科書一般就給到這裡然後留給了讀者無窮的遐想,以及無窮的吐槽,其實教科書只要補一張圖就足夠了:頻域影象,也就是俗稱的頻譜,就是—— 在這裡插入圖片描述 再清楚一點: 在這裡插入圖片描述 可以發現,在頻譜中,偶數項的振幅都是0,也就對應了圖中的彩色直線。振幅為0的正弦波。 老實說,在我學傅立葉變換時,維基的這個圖還沒有出現,那時我就想到了這種表達方法,而且,後面還會加入維基沒有表示出來的另一個譜——相位譜。

但是在講相位譜之前,我們先回顧一下剛剛的這個例子究竟意味著什麼。記得前面說過的那句“世界是靜止的”嗎?估計好多人對這句話都已經吐槽半天了。想象一下,世界上每一個看似混亂的表象,實際都是一條時間軸上不規則的曲線,但實際這些曲線都是由這些無窮無盡的正弦波組成。我們看似不規律的事情反而是規律的正弦波在時域上的投影,而正弦波又是一個旋轉的圓在直線上的投影。那麼你的腦海中會產生一個什麼畫面呢?

我們眼中的世界就像皮影戲的大幕布,幕布的後面有無數的齒輪,大齒輪帶動小齒輪,小齒輪再帶動更小的。在最外面的小齒輪上有一個小人——那就是我們自己。我們只看到這個小人毫無規律的在幕布前表演,卻無法預測他下一步會去哪。而幕布後面的齒輪卻永遠一直那樣不停的旋轉,永不停歇。這樣說來有些宿命論的感覺。說實話,這種對人生的描繪是我一個朋友在我們都是高中生的時候感嘆的,當時想想似懂非懂,直到有一天我學到了傅立葉級數……

三、傅立葉級數(Fourier Series)的相位譜

上一章的關鍵詞是:從側面看。這一章的關鍵詞是:從下面看。

在這一章最開始,我想先回答很多人的一個問題:傅立葉分析究竟是幹什麼用的?這段相對比較枯燥,已經知道了的同學可以直接跳到下一個分割線。

先說一個最直接的用途。無論聽廣播還是看電視,我們一定對一個詞不陌生——頻道。頻道頻道,就是頻率的通道,不同的頻道就是將不同的頻率作為一個通道來進行資訊傳輸。下面大家嘗試一件事:

先在紙上畫一個sin(x),不一定標準,意思差不多就行。不是很難吧。

好,接下去畫一個sin(3x)+sin(5x)的圖形。

別說標準不標準了,曲線什麼時候上升什麼時候下降你都不一定畫的對吧?

好,畫不出來不要緊,我把sin(3x)+sin(5x)的曲線給你,但是前提是你不知道這個曲線的方程式,現在需要你把sin(5x)給我從圖裡拿出去,看看剩下的是什麼。這基本是不可能做到的。

但是在頻域呢?則簡單的很,無非就是幾條豎線而已。

所以很多在時域看似不可能做到的數學操作,在頻域相反很容易。這就是需要傅立葉變換的地方。尤其是從某條曲線中去除一些特定的頻率成分,這在工程上稱為濾波,是訊號處理最重要的概念之一,只有在頻域才能輕鬆的做到。

再說一個更重要,但是稍微複雜一點的用途——求解微分方程。(這段有點難度,看不懂的可以直接跳過這段)微分方程的重要性不用我過多介紹了。各行各業都用的到。但是求解微分方程卻是一件相當麻煩的事情。因為除了要計算加減乘除,還要計算微分積分。而傅立葉變換則可以讓微分和積分在頻域中變為乘法和除法,大學數學瞬間變小學算術有沒有。

傅立葉分析當然還有其他更重要的用途,我們隨著講隨著提。

————————————————————————————————————

下面我們繼續說相位譜:

通過時域到頻域的變換,我們得到了一個從側面看的頻譜,但是這個頻譜並沒有包含時域中全部的資訊。因為頻譜只代表每一個對應的正弦波的振幅是多少,而沒有提到相位。基礎的正弦波A.sin(wt+θ)中,振幅,頻率,相位缺一不可,不同相位決定了波的位置,所以對於頻域分析,僅僅有頻譜(振幅譜)是不夠的,我們還需要一個相位譜。那麼這個相位譜在哪呢?我們看下圖,這次為了避免圖片太混論,我們用7個波疊加的圖。 在這裡插入圖片描述 鑑於正弦波是週期的,我們需要設定一個用來標記正弦波位置的東西。在圖中就是那些小紅點。小紅點是距離頻率軸最近的波峰,而這個波峰所處的位置離頻率軸有多遠呢?為了看的更清楚,我們將紅色的點投影到下平面,投影點我們用粉色點來表示。當然,這些粉色的點只標註了波峰距離頻率軸的距離,並不是相位。 在這裡插入圖片描述 這裡需要糾正一個概念:時間差並不是相位差。如果將全部週期看作2Pi或者360度的話,相位差則是時間差在一個週期中所佔的比例。我們將時間差除週期再乘2Pi,就得到了相位差。 在完整的立體圖中,我們將投影得到的時間差依次除以所在頻率的週期,就得到了最下面的相位譜。所以,頻譜是從側面看,相位譜是從下面看。下次偷看女生裙底被發現的話,可以告訴她:“對不起,我只是想看看你的相位譜。”

注意到,相位譜中的相位除了0,就是Pi。因為cost+Pi)=costcos(t+Pi)= -cos(t),所以實際上相位為Pi的波只是上下翻轉了而已。對於週期方波的傅立葉級數,這樣的相位譜已經是很簡單的了。另外值得注意的是,由於cost+2Pi=costcos(t+2Pi)=cos(t),所以相位差是週期的,pi和3pi,5pi,7pi都是相同的相位。人為定義相位譜的值域為(pipi](-pi,pi],所以圖中的相位差均為Pi。 最後來一張大集合: 在這裡插入圖片描述

四、傅立葉變換(Fourier Transformation)

相信通過前面三章,大家對頻域以及傅立葉級數都有了一個全新的認識。但是文章在一開始關於鋼琴琴譜的例子我曾說過,這個栗子是一個公式錯誤,但是概念典型的例子。所謂的公式錯誤在哪裡呢?

傅立葉級數的本質是將一個週期的訊號分解成無限多分開的(離散的)正弦波,但是宇宙似乎並不是週期的。曾經在學數字訊號處理的時候寫過一首打油詩: 往昔連續非週期,

回憶週期不連續,

任你ZT、DFT,

還原不回去。

(請無視我渣一樣的文學水平……)

在這個世界上,有的事情一期一會,永不再來,並且時間始終不曾停息地將那些刻骨銘心的往昔連續的標記在時間點上。但是這些事情往往又成為了我們格外寶貴的回憶,在我們大腦裡隔一段時間就會週期性的蹦出來一下,可惜這些回憶都是零散的片段,往往只有最幸福的回憶,而平淡的回憶則逐漸被我們忘卻。因為,往昔是一個連續的非週期訊號,而回憶是一個週期離散訊號。

是否有一種數學工具將連續非週期訊號變換為週期離散訊號呢?抱歉,真沒有。

比如傅立葉級數,在時域是一個週期且連續的函式,而在頻域是一個非週期離散的函式。這句話比較繞嘴,實在看著費事可以乾脆回憶第一章的圖片。

而在我們接下去要講的傅立葉變換,則是將一個時域非週期的連續訊號,轉換為一個在頻域非週期的連續訊號。

算了,還是上一張圖方便大家理解吧: 在這裡插入圖片描述 或者我們也可以換一個角度理解:傅立葉變換實際上是對一個週期無限大的函式進行傅立葉變換。 所以說,鋼琴譜其實並非一個連續的頻譜,而是很多在時間上離散的頻率,但是這樣的一個貼切的比喻真的是很難找出第二個來了。

因此在傅立葉變換在頻域上就從離散譜變成了連續譜。那麼連續譜是什麼樣子呢?

你見過大海麼?

為了方便大家對比,我們這次從另一個角度來看頻譜,還是傅立葉級數中用到最多的那幅圖,我們從頻率較高的方向看。 在這裡插入圖片描述 以上是離散譜,那麼連續譜是什麼樣子呢? 盡情的發揮你的想象,想象這些離散的正弦波離得越來越近,逐漸變得連續……

直到變得像波濤起伏的大海: 在這裡插入圖片描述

很抱歉,為了能讓這些波浪更清晰的看到,我沒有選用正確的計算引數,而是選擇了一些讓圖片更美觀的引數,不然這圖看起來就像屎一樣了。

不過通過這樣兩幅圖去比較,大家應該可以理解如何從離散譜變成了連續譜的了吧?原來離散譜的疊加,變成了連續譜的累積。所以在計算上也從求和符號變成了積分符號。

不過,這個故事還沒有講完,接下去,我保證讓你看到一幅比上圖更美麗壯觀的圖片,但是這裡需要介紹到一個數學工具才能然故事繼續,這個工具就是——

五、宇宙耍帥第一公式:尤拉公式

虛數i這個概念大家在高中就接觸過,但那時我們只知道它是-1的平方根,可是它真正的意義是什麼呢? 在這裡插入圖片描述 這裡有一條數軸,在數軸上有一個紅色的線段,它的長度是1。當它乘以3的時候,它的長度發生了變化,變成了藍色的線段,而當它乘以-1的時候,就變成了綠色的線段,或者說線段在數軸上圍繞原點旋轉了180度。 我們知道乘-1其實就是乘了兩次 i使線段旋轉了180度,那麼乘一次 i 呢——答案很簡單——旋轉了90度。 在這裡插入圖片描述 同時,我們獲得了一個垂直的虛數軸。實數軸與虛數軸共同構成了一個複數的平面,也稱複平面。這樣我們就瞭解到,乘虛數i的一個功能——旋轉。 現在,就有請宇宙第一耍帥公式尤拉公式隆重登場——

在這裡插入圖片描述 這個公式在數學領域的意義要遠大於傅立葉分析,但是乘它為宇宙第一耍帥公式是因為它的特殊形式——當x等於Pi的時候。 在這裡插入圖片描述

經常有理工科的學生為了跟妹子表現自己的學術功底,用這個公式來給妹子解釋數學之美:”石榴姐你看,這個公式裡既有自然底數e,自然數1和0,虛數i還有圓周率pi,它是這麼簡潔,這麼美麗啊!“但是姑娘們心裡往往只有一句話:”臭屌絲……“ 這個公式關鍵的作用,是將正弦波統一成了簡單的指數形式。我們來看看影象上的涵義: 在這裡插入圖片描述 尤拉公式所描繪的,是一個隨著時間變化,在複平面上做圓周運動的點,隨著時間的改變,在時間軸上就成了一條螺旋線。如果只看它的實數部分,也就是螺旋線在左側的投影,就是一個最基礎的餘弦函式。而右側的投影則是一個正弦函式。 關於複數更深的理解,大家可以參考:

複數的物理意義是什麼?

這裡不需要講的太複雜,足夠讓大家理解後面的內容就可以了。

六、指數形式的傅立葉變換

有了尤拉公式的幫助,我們便知道:正弦波的疊加,也可以理解為螺旋線的疊加在實數空間的投影。而螺旋線的疊加如果用一個形象的栗子來理解是什麼呢?

光波

高中時我們就學過,自然光是由不同顏色的光疊加而成的,而最著名的實驗就是牛頓師傅的三稜鏡實驗: 在這裡插入圖片描述 所以其實我們在很早就接觸到了光的頻譜,只是並沒有瞭解頻譜更重要的意義。 但不同的是,傅立葉變換出來的頻譜不僅僅是可見光這樣頻率範圍有限的疊加,而是頻率從0到無窮所有頻率的組合。

這裡,我們可以用兩種方法來理解正弦波:

第一種前面已經講過了,就是螺旋線在實軸的投影。

另一種需要藉助尤拉公式的另一種形式去理解: 在這裡插入圖片描述 我們剛才講過,eite^{it}可以理解為一條逆時針旋轉的螺旋線,那麼eite^{-it}則可以理解為一條順時針旋轉的螺旋線。而cos(t)則是這兩條旋轉方向不同的螺旋線疊加的一半,因為這兩條螺旋線的虛數部分相互抵消掉了!

舉個例子的話,就是極化方向不同的兩束光波,磁場抵消,電場加倍。

這裡,逆時針旋轉的我們稱為正頻率,而順時針旋轉的我們稱為負頻率(注意不是複頻率)。

好了,剛才我們已經看到了大海——連續的傅立葉變換頻譜,現在想一想,連續的螺旋線會是什麼樣子:

想象一下再往下翻: 在這裡插入圖片描述 是不是很漂亮? 你猜猜,這個圖形在時域是什麼樣子? 在這裡插入圖片描述 哈哈,是不是覺得被狠狠扇了一個耳光。數學就是這麼一個把簡單的問題搞得很複雜的東西。 順便說一句,那個像大海螺一樣的圖,為了方便觀看,我僅僅展示了其中正頻率的部分,負頻率的部分沒有顯示出來。

如果你認真去看,海螺圖上的每一條螺旋線都是可以清楚的看到的,每一條螺旋線都有著不同的振幅(旋轉半徑),頻率(旋轉週期)以及相位。而將所有螺旋線連成平面,就是這幅海螺圖了。

好了,講到這裡,相信大家對傅立葉變換以及傅立葉級數都有了一個形象的理解了,我們最後用一張圖來總結一下: 在這裡插入圖片描述

七、基於opencv的影象傅立葉變換

#include <opencv2/core/core.hpp>
#include <opencv2/highgui/highgui.hpp>
#include <opencv2/imgproc/imgproc.hpp>
#include <iostream>
#include <Windows.h>

using namespace cv;
using namespace std;

static void help()
{
	printf("\nThis program demonstrated the use of the discrete Fourier transform (dft)\n"
		"The dft of an image is taken and it's power spectrum is displayed.\n"
		"Usage:\n"
		"./dft [image_name -- default lena.jpg]\n");
}

int main()
{
	SetConsoleTextAttribute(GetStdHandle(STD_OUTPUT_HANDLE), FOREGROUND_INTENSITY | FOREGROUND_GREEN);
	help();
	Mat src = imread("C:\\Users\\18301\\Desktop\\1.jpg");
	if (src.empty())
	{
		return -1;
	}

	Mat srcGray;
	cvtColor(src, srcGray, CV_RGB2GRAY); //灰度影象做傅立葉變換

	int m = getOptimalDFTSize(srcGray.rows); //2,3,5的倍數有更高效率的傅立葉變換
	int n = getOptimalDFTSize(srcGray.cols);
	cout << m << endl;
	cout << n << endl;
	Mat padded;
	//把灰度影象放在左上角,在右邊和下邊擴充套件影象,擴充套件部分填充為0;
	copyMakeBorder(srcGray, padded, 0, m - srcGray.rows, 0, n - srcGray.cols, BORDER_CONSTANT, Scalar::all(0));
	cout << padded.size() << endl;
	//這裡是獲取了兩個Mat,一個用於存放dft變換的實部,一個用於存放虛部,初始的時候,實部就是影象本身,虛部全為零
	Mat planes[] = { Mat_<float>(padded), Mat::zeros(padded.size(), CV_32F) };//在這裡是宣告兩個Mat矩陣,分別存放實部和虛部
	Mat planes_true = Mat_<float>(padded);

	Mat complexImg;
	//將幾個單通道的mat融合成一個多通道的mat,這裡融合的complexImg既有實部又有虛部
	merge(planes, 2, complexImg);

	//對上邊合成的mat進行傅立葉變換,***支援原地操作***,傅立葉變換結果為複數.通道1存的是實部,通道二存的是虛部
	dft(complexImg, complexImg);
	//把變換後的結果分割到兩個mat,一個實部,一個虛部,方便後續操作
	split(complexImg, planes);
	//這一部分是為了計算dft變換後的幅值,傅立葉變換的幅度值範圍大到不適合在螢幕上顯示。高值在螢幕上顯示為白點,而低值為黑點,高低值的變化無法有效分辨。為了在螢幕上凸顯出高低變化的連續性,我們可以用對數尺度來替換線性尺度,以便於顯示幅值,計算公式如下:
	//=> log(1 + sqrt(Re(DFT(I))^2 +Im(DFT(I))^2))
	//計算實數和虛數的幅值,並把幅值儲存到palnes[0],允許線上覆蓋
	magnitude(planes[0], planes[1], planes_true);
	Mat A = planes[0];
	Mat B = planes[1];
	Mat mag = planes_true;
	//對幅值加1
	mag += Scalar::all(1);
	log(mag, mag);//計算出的幅值一般很大,達到10^4,通常沒有辦法在影象中顯示出來,需要對其進行log求解。
	//在這裡進行一次低通濾波
	for (int i = 0; i < mag.rows; i++)
	{
		for (int j = 0; j < mag.cols; j++)
		{
			float num = mag.at<float>(i, j);
			if (num>13.5)//13.5為幅度閾值,這裡是低通濾波,高通濾波只需要改成<就好
			{
				planes[0].at<float>(i, j) = 0;
				planes[1].at<float>(i, j) = 0;
			}
		}
	}
	merge(planes, 2, complexImg);
	//crop the spectrum, if it has an odd number of rows or columns
	//修剪頻譜,如果影象的行或者列是奇數的話,那其頻譜是不對稱的,因此要修剪
	//這裡為什麼要用  &-2這個操作,我會在程式碼後面的 注2 說明
	//我們知道x&-2代表x與 - 2按位相與,而 - 2的二進位制形式是2的二進位制取反加一的結果(這是補碼的問題)。2 的二進位制結果是(假設
	//8位表示,實際整型是32位,但是描述方式是一樣的,為便於描述,用8位表示)0000 0010,則 - 2的二進位制形式為:1111 1110,
	//x與 - 2按位相與後,不管x是奇數還是偶數,最後x都會變成一個偶數。 
	//就是說dft這個函式雖然對於輸入mat的尺寸不做要求,但是如果其行數和列數可以分解為2、3、5的乘積,那麼對於dft運算的速度會加快很多。
	mag = mag(Rect(0, 0, mag.cols & -2, mag.rows & -2));
	Mat _magI = mag.clone();
	//這一步的目的仍然是為了顯示,但是幅度值仍然超過可顯示範圍[0,1],我們使用 normalize() 函式將幅度歸一化到可顯示範圍。
	normalize(_magI, _magI, 0, 1, CV_MINMAX);
	namedWindow("before rearrange", 0);
	imshow("before rearrange", _magI);

	//rearrange the quadrants of Fourier image
	//so that the origin is at the image center
	//重新分配象限,使(0,0)移動到影象中心,  
	//在《數字影象處理》中,傅立葉變換之前要對源影象乘以(-1)^(x+y)進行中心化。  
	//這是是對傅立葉變換結果進行中心化
	int cx = mag.cols / 2;
	int cy = mag.rows / 2;

	//這裡是以中心為標準,把mag影象分成四部分
	Mat tmp;
	Mat q0(mag, Rect(0, 0, cx, cy));   //Top-Left - Create a ROI per quadrant
	Mat q1(mag, Rect(cx, 0, cx, cy));  //Top-Right
	Mat q2(mag, Rect(0, cy, cx, cy));  //Bottom-Left
	Mat q3(mag, Rect(cx, cy, cx, cy)); //Bottom-Right

	//swap quadrants(Top-Left with Bottom-Right),交換象限
	q0.copyTo(tmp);
	q3.copyTo(q0);
	tmp.copyTo(q3);

	// swap quadrant (Top-Rightwith Bottom-Left),交換象限
	q1.copyTo(tmp);
	q2.copyTo(q1);
	tmp.copyTo(q2);


	normalize(mag, mag, 0, 1, CV_MINMAX);
	namedWindow("Input Image", 0);
	namedWindow("spectrum magnitude", 0);
	imshow("Input Image", srcGray);
	imshow("spectrum magnitude", mag);

	//傅立葉的逆變換
	Mat ifft;
	//傅立葉逆變換
	idft(complexImg, ifft, DFT_REAL_OUTPUT);
	normalize(ifft, ifft, 0, 1, CV_MINMAX);
	namedWindow("inverse fft", 0);
	imshow("inverse fft", ifft);
	waitKey();
	return 0;
}

在這裡插入圖片描述 在這裡插入圖片描述 在這裡插入圖片描述 傅立葉轉換的結果是複數,這也顯示出了傅立葉變換是一副實數影象(real image)和虛數影象(complex image)疊加或者是幅度影象(magitude image)和相點陣圖像(phase image)疊加的結果。在實際的影象處理演算法中僅有幅度影象(magnitude image)影象能夠用到,因為幅度影象包含了我們所需要的所有影象幾何結構的資訊。但是,如果想通過修改幅度影象或者相點陣圖像來間接修改原空間影象,需要保留幅度影象和相點陣圖像來進行傅立葉逆變換,從而得到修改後影象。

注意:通過修改幅度影象或者相點陣圖像來間接修改原空間影象。
幅度影象(改過) + 相點陣圖像 = 修改後影象(通過傅立葉逆變換)
影象 = 幅度影象 + 相點陣圖像(傅立葉變換)