與半二次分裂快速反卷積)
簡(jiǎn)介這份資源是圖像去模糊方向經(jīng)典論文《Fast Image Deconvolution using Hyper-Laplacian Priors》的配套Matlab實(shí)現(xiàn)面向圖像處理研究者、計(jì)算機(jī)視覺方向?qū)W生及需要復(fù)現(xiàn)去卷積算法的開發(fā)者。它針對(duì)模糊、噪聲等成像缺陷導(dǎo)致的圖像質(zhì)量下降問題利用超拉普拉斯先驗(yàn)對(duì)圖像邊緣與細(xì)節(jié)的稀疏建模能力在盲去卷積場(chǎng)景下同時(shí)估計(jì)模糊核與清晰圖像并借助迭代優(yōu)化提升恢復(fù)速度與質(zhì)量。壓縮包共7個(gè)文件約2.05MB包含4個(gè)m腳本、1個(gè)mat數(shù)據(jù)文件、1張jpg測(cè)試圖及1份readme說明腳本覆蓋主流程、圖像求解與信噪比評(píng)估等環(huán)節(jié)可直接運(yùn)行驗(yàn)證算法效果。目前已有1126人學(xué)習(xí)下載適合作為理解超拉普拉斯先驗(yàn)、復(fù)現(xiàn)論文結(jié)果或搭建自身去模糊項(xiàng)目的起點(diǎn)也可在此基礎(chǔ)上對(duì)比PSNR、SSIM等指標(biāo)并探索后續(xù)改進(jìn)。1. 圖像去模糊的反直覺真相為什么銳化濾鏡救不了糊圖很多人第一次處理模糊照片第一反應(yīng)是拉高銳化。結(jié)果呢噪點(diǎn)被放大邊緣出現(xiàn)白邊糊的地方還是糊。這不是濾鏡不夠強(qiáng)而是方向從一開始就錯(cuò)了——模糊在數(shù)學(xué)上是一個(gè)卷積過程想還原它必須做反卷積deconvolution而不是在像素上做局部增強(qiáng)。Fast Image Deconvolution using Hyper-Laplacian Priors 這份資源講的正是反卷積里最核心的一環(huán)先驗(yàn)怎么選。傳統(tǒng)做法假設(shè)圖像梯度服從高斯分布解出來又快又平滑但邊緣全被抹掉換成超拉普拉斯先驗(yàn)邊緣能保住代價(jià)是求解變成非凸優(yōu)化慢得讓人想砸鍵盤。這份工作要解決的就是既要邊緣銳利、又要跑得快這個(gè)看似矛盾的命題。搞懂它你就明白為什么有些去模糊算法一秒出圖還干凈有些跑十分鐘還是一團(tuán)漿糊。2. 從卷積到反卷積模糊到底把什么信息弄丟了2.1 模糊的數(shù)學(xué)模型一次卷積信息被攤平一張清晰圖像 $x$ 經(jīng)過模糊核 $k$ 卷積再加上噪聲 $n$得到觀測(cè)圖像 $y$$$y k \otimes x n$$這里的 $\otimes$ 是卷積。直觀理解模糊核就像一個(gè)涂抹刷子把每個(gè)像素的能量按權(quán)重?cái)偟街車?。運(yùn)動(dòng)模糊的核是一條線失焦模糊的核是一個(gè)圓盤。卷積的致命之處在于——它是一個(gè)低通濾波高頻信息也就是邊緣、紋理被壓制甚至歸零。反卷積就是已知 $y$ 和 $k$求 $x$。但問題來了高頻被壓到接近零反解時(shí)就要除以一個(gè)接近零的數(shù)噪聲會(huì)被無限放大。這就是為什么直接做逆濾波必然翻車——數(shù)學(xué)上叫病態(tài)問題ill-posed。2.2 為什么必須加先驗(yàn)沒有約束解有無窮多個(gè)病態(tài)意味著解不唯一。同一個(gè)模糊圖 $y$可以對(duì)應(yīng)無數(shù)張清晰圖 $x$它們卷積后都能得到 $y$。要挑出最像真實(shí)照片的那一個(gè)就得引入先驗(yàn)——也就是我們對(duì)自然圖像的統(tǒng)計(jì)認(rèn)知。最常見的先驗(yàn)是高斯先驗(yàn)等價(jià)于在優(yōu)化里加一個(gè) $L_2$ 正則項(xiàng)$$\min_x |k \otimes x - y|^2 \lambda | \nabla x |^2$$這個(gè)式子有閉式解用 FFT 幾行代碼就能算完快得飛起。但它假設(shè)圖像梯度服從高斯分布而真實(shí)自然圖像的梯度分布是重尾的——大部分地方梯度接近零平坦區(qū)域少數(shù)地方梯度極大邊緣。高斯先驗(yàn)壓不住這個(gè)重尾結(jié)果就是邊緣被當(dāng)成異常值懲罰掉圖像變糊。2.3 超拉普拉斯先驗(yàn)用 $L_p$$0p1$逼近真實(shí)梯度分布自然圖像梯度的真實(shí)分布用 $L_p$ 范數(shù)擬合時(shí)$p$ 通常在 0.5 到 0.8 之間這就是超拉普拉斯hyper-Laplacian。對(duì)應(yīng)的優(yōu)化目標(biāo)變成$$\min_x |k \otimes x - y|^2 \lambda \sum_i |\nabla x_i|^p, \quad 0 p 1$$$p1$ 時(shí)$| \cdot |^p$ 在零點(diǎn)不可導(dǎo)、且是非凸的。好處是它對(duì)小梯度懲罰極重平坦區(qū)域被強(qiáng)力平滑對(duì)大梯度懲罰相對(duì)輕邊緣被保留。壞處是——非凸、非光滑沒有閉式解傳統(tǒng)梯度下降收斂慢且容易卡住。這就是整份資源要啃的硬骨頭怎么在 $p1$ 的非凸目標(biāo)上做到接近 $L_2$ 的速度。3. 半二次分裂把非凸問題拆成兩個(gè)能解的凸子問題3.1 變量分裂的核心思路引入輔助變量解耦直接優(yōu)化 $| \nabla x |^p$ 很難因?yàn)榉稊?shù)和卷積耦合在一起。半二次分裂Half-Quadratic Splitting, HQS的做法是引入輔助變量 $w$讓它去逼近 $\nabla x$把目標(biāo)改寫成$$\min_{x,w} |k \otimes x - y|^2 \lambda \sum_i |w_i|^p \beta \sum_i (w_i - \nabla x_i)^2$$$\beta$ 是懲罰系數(shù)。當(dāng) $\beta \to \infty$ 時(shí)$w$ 被迫等于 $\nabla x$原問題等價(jià)。但實(shí)際計(jì)算時(shí) $\beta$ 取有限值交替優(yōu)化 $x$ 和 $w$固定 $w$對(duì) $x$ 求最小這是一個(gè)最小二乘問題因?yàn)橹皇?$L_2$ 項(xiàng)可以在頻域用 FFT 一步解出。固定 $x$對(duì) $w$ 求最小每個(gè)像素獨(dú)立變成一個(gè)一維標(biāo)量?jī)?yōu)化問題。這個(gè)拆分的精妙之處在于原本耦合的非凸問題被拆成了頻域可閉式求解和逐像素可查表求解兩塊每塊都極快。3.2 $x$ 子問題頻域閉式解為什么能一步到位固定 $w$ 后$x$ 子問題是$$\min_x |k \otimes x - y|^2 \beta |w - \nabla x|^2$$卷積在空域是復(fù)雜的積分但轉(zhuǎn)到頻域就是逐點(diǎn)相乘。對(duì)兩邊做 FFT卷積變成乘法梯度算子 $\nabla$ 也有對(duì)應(yīng)的頻域形式。設(shè) $K$、$Y$、$W$ 分別是 $k$、$y$、$w$ 的傅里葉變換$D_x$、$D_y$ 是梯度算子的頻域表示則最優(yōu)解為$$X \frac{\overline{K} Y \beta (\overline{D_x} W_x \overline{D_y} W_y)}{|K|^2 \beta (|D_x|^2 |D_y|^2)}$$橫線表示復(fù)共軛。分母是逐點(diǎn)實(shí)數(shù)分子是逐點(diǎn)復(fù)數(shù)整個(gè)求解就是幾次 FFT 加逐點(diǎn)除法。這就是為什么它能一步到位——卷積在頻域沒有耦合每個(gè)頻率點(diǎn)獨(dú)立。import numpy as np def solve_x_subproblem(y, k, wx, wy, beta): 頻域閉式求解 x 子問題 y: 觀測(cè)模糊圖 (H, W) k: 模糊核需補(bǔ)零到圖像尺寸 wx, wy: 輔助變量 w 的兩個(gè)分量 beta: 懲罰系數(shù) H, W y.shape K np.fft.fft2(k) # 模糊核頻域 Y np.fft.fft2(y) Wx np.fft.fft2(wx) Wy np.fft.fft2(wy) # 梯度算子的頻域表示前向差分 dx np.zeros((H, W)); dx[0, 0] -1; dx[0, -1] 1 dy np.zeros((H, W)); dy[0, 0] -1; dy[-1, 0] 1 Dx np.fft.fft2(dx) Dy np.fft.fft2(dy) # 逐點(diǎn)求解分母加 eps 防止除零 numerator np.conj(K) * Y beta * (np.conj(Dx) * Wx np.conj(Dy) * Wy) denominator np.abs(K)**2 beta * (np.abs(Dx)**2 np.abs(Dy)**2) 1e-8 X numerator / denominator return np.real(np.fft.ifft2(X))這段代碼里np.conj(K) * Y對(duì)應(yīng)數(shù)據(jù)保真項(xiàng)beta * (...)對(duì)應(yīng)先驗(yàn)耦合項(xiàng)分母的1e-8是數(shù)值穩(wěn)定用的。注意dx、dy的構(gòu)造方式?jīng)Q定了邊界條件——這里用的是循環(huán)邊界實(shí)際工程中若邊界偽影明顯需要改用其他邊界處理。3.3 $w$ 子問題逐像素查表為什么比迭代快幾十倍固定 $x$ 后$w$ 子問題對(duì)每個(gè)像素獨(dú)立$$\min_{w_i} |w_i|^p \beta (w_i - v_i)^2$$其中 $v_i \nabla x_i$ 是已知量。這個(gè)一維問題的解 $w_i^$ 只依賴于 $v_i$ 和 $\beta$。由于 $p1$ 時(shí)沒有解析解常見做法是預(yù)先建表把 $v_i$ 的取值范圍離散化對(duì)每個(gè)格點(diǎn)用牛頓法離線解出 $w_i^$運(yùn)行時(shí)直接查表加插值。這就是速度的關(guān)鍵——把在線迭代變成離線查表。一張 512×512 的圖有 26 萬個(gè)像素每個(gè)像素做一次查表比每個(gè)像素跑十幾次牛頓迭代快一個(gè)數(shù)量級(jí)以上。查表的精度取決于離散粒度粒度太粗會(huì)有階梯效應(yīng)太細(xì)則內(nèi)存吃緊通常取幾千個(gè)格點(diǎn)就能在精度和速度間取得平衡。4. 避坑與排查超拉普拉斯去模糊最容易翻車的五個(gè)地方4.1 現(xiàn)象結(jié)果出現(xiàn)振鈴邊緣一圈圈波紋原因頻域求解默認(rèn)循環(huán)邊界圖像左右、上下邊界被當(dāng)成相連卷積核在邊界處繞回來造成不連續(xù)。解決對(duì)圖像做邊界延拓如鏡像填充或者改用基于空間卷積的求解方式。簡(jiǎn)單驗(yàn)證方法是把圖像四周各裁掉 20 像素再看如果振鈴消失基本就是邊界問題。4.2 現(xiàn)象迭代不收斂$x$ 和 $w$ 來回震蕩原因$\beta$ 初始值太大或太小。$\beta$ 太大$w$ 被強(qiáng)行拉向 $\nabla x$$x$ 子問題條件數(shù)惡化$\beta$ 太小$w$ 和 $\nabla x$ 脫節(jié)等價(jià)性不成立。解決采用連續(xù)化策略$\beta$ 從小到大逐步增大每次增大后熱啟動(dòng)上一輪的解。典型序列是 $\beta_01$每輪乘 2跑 5 到 10 輪。4.3 現(xiàn)象平坦區(qū)域出現(xiàn)塊狀偽影原因$w$ 子問題查表粒度太粗或者 $p$ 取值過小如 0.3導(dǎo)致小梯度被過度懲罰平坦區(qū)被切成一塊塊。解決$p$ 一般取 0.5 到 0.8查表格點(diǎn)數(shù)不低于 2000。如果偽影仍在檢查 $\lambda$ 是否過大——$\lambda$ 控制先驗(yàn)強(qiáng)度過大會(huì)讓圖像過度平滑。4.4 現(xiàn)象噪聲被放大暗部出現(xiàn)彩色斑點(diǎn)原因反卷積本質(zhì)是放大高頻噪聲也在高頻。如果數(shù)據(jù)保真項(xiàng)權(quán)重相對(duì)先驗(yàn)太強(qiáng)噪聲就會(huì)壓過先驗(yàn)約束。解決適當(dāng)增大 $\lambda$或在反卷積前先做一次輕度降噪。注意降噪不能太狠否則邊緣信息先丟了反卷積也救不回來。4.5 現(xiàn)象模糊核估計(jì)不準(zhǔn)去模糊后出現(xiàn)鬼影原因這份資源假設(shè)模糊核已知。如果核估計(jì)有偏差反卷積會(huì)把偏差放大成重影。解決先用盲去模糊方法估核或者用一張標(biāo)定板拍下核形狀。核的尺寸也要匹配——核估小了邊緣去不干凈核估大了會(huì)引入額外振鈴。5. 設(shè)計(jì)取舍為什么快和好能同時(shí)成立超拉普拉斯先驗(yàn)這套方法最值得琢磨的地方是它沒有在快和好之間做二選一而是通過問題結(jié)構(gòu)找到了兩全的路徑。核心洞察有三層。第一層是先驗(yàn)選擇的取舍。高斯先驗(yàn)快但糊超拉普拉斯好但慢這是統(tǒng)計(jì)建模層面的矛盾。第二層是求解策略的取舍。變量分裂把耦合問題拆開代價(jià)是引入輔助變量和懲罰系數(shù)需要調(diào)參和迭代。第三層是計(jì)算方式的取舍。頻域閉式解快但受邊界條件限制查表快但精度受離散粒度限制。這三層取舍疊在一起最終達(dá)成的效果是在 $p1$ 的非凸目標(biāo)上每輪迭代只花幾次 FFT 加一次查表整體速度比通用非凸優(yōu)化器快一到兩個(gè)數(shù)量級(jí)同時(shí)邊緣保留能力遠(yuǎn)超 $L_2$ 方法。這不是靠某個(gè)單點(diǎn)技巧而是靠把每個(gè)子問題都推到它能被最快求解的形式。方法先驗(yàn)求解方式邊緣保留單輪速度逆濾波無頻域除法極差極快維納濾波高斯頻域閉式差極快TV 正則$L_1$迭代中中超拉普拉斯$L_p, p1$分裂查表好快驗(yàn)證這套方法是否真的 work最直接的辦法是造一張已知清晰圖和已知核的合成模糊圖跑完反卷積后算 PSNR 和 SSIM再和維納濾波對(duì)比。如果 PSNR 高出 2dB 以上且邊緣目視更銳說明實(shí)現(xiàn)方向?qū)α恕_M(jìn)階方向可以往盲去模糊走——把核也當(dāng)成未知量交替優(yōu)化或者把先驗(yàn)從梯度域擴(kuò)展到更復(fù)雜的圖像塊先驗(yàn)。我自己的習(xí)慣是每次實(shí)現(xiàn)這類反卷積算法先用 64×64 的小圖跑通全流程確認(rèn)頻域求解和查表都沒寫錯(cuò)再上大圖。因?yàn)轭l域里一個(gè)共軛符號(hào)寫反小圖上可能只是輕微模糊大圖上就是滿屏鬼影排查起來極其痛苦。從那以后我每次寫頻域代碼都強(qiáng)制先跑一個(gè) delta 函數(shù)的 FFT 驗(yàn)證對(duì)稱性這個(gè)習(xí)慣幫我省了無數(shù)時(shí)間。希望幫到你。本文還有配套的精品資源點(diǎn)擊獲取