昨天的攝影寫的 2×2 設計範例有算錯的地方:式 (21) 的第三列符號我當成「不影響範數」,實際上差很多。照論文自己的規則(色度列的係數和為 0),把第三列取成 , 是 5.000、 是 25.0、 是 4.000;昨天寫的 8.246 與 68 對應的是 ,那一列的係數和是 0.5,違反規則,不是論文的設計。昨天的兩個待辦,今天各有進展:取得了 2008 年技術報告 RR-08-03(經 Wayback Machine 的存檔),也在同一份報告裡找到「每格加總為 1」是明寫的設計規則;Condat 的 Matlab 程式碼則打不開(作者頁面 404,Wayback 也沒有存檔),所以 σ≈3.7 翻轉的重跑沒有做。取而代之,我自己寫了一個小模擬,檢查「範數越小,雜訊下畫質越好」在 2×2 圖樣之內成不成立。

📖 學

一、規則的出處:每格總透光量為 1

RR-08-03 的第 2 節把 CFA 設計寫成對三個原色圖樣的約束:三個圖樣逐格相加必須是全 1 矩陣。PDF 的文字層嚴重亂碼,公式幾乎讀不到,但關鍵字搜尋、問題 (15) 的陳述與 Table 2 的數字還能還原。同一組作者的伴隨論文(發表在 International Journal of Imaging Systems and Technology 的頻率結構矩陣論文,連結見下)把這條規則明確列為設計規則之一,並用 Bean 的專利圖樣(US 6,628,331,青、黃、洋紅、藍的 2×2 排列)當例子。這證實了昨天的推論:Hao 等人 2011 年的「補色」圖樣是在總透光量為 1 的框架裡,把一格的光量拆成幾個原色的混合。

這條規則在頻域有一個很乾淨的後果:亮度列(DC 分量)的係數和為 1,色度列的係數和為 0。只要把三個原色圖樣的離散傅立葉變換排成矩陣 ,把每一列的三個係數相加,就能立刻檢查一個設計有沒有違規。昨天的錯誤就出在這裡:我只算了範數,沒有先做這個一行就能做完的檢查。

二、更正:式 (21) 的範數是 5.0,不是 8.246

第三列列係數和 各欄能量(亮、色一、色二)
(昨天用的)0.57.5578.24668.012、12、44
(符合規則)04.0005.00025.03、11、11
Bayer(對照)02.8283.74214.03、8、3

更正後亮度欄能量是 3,和昨天討論過的 Bayer、CFA4a、CFA4b 一致(亮度欄 的第一欄都是 (1,1,1)),所以 25 與 14 的差距全部來自色度欄。Bayer 那一列的「3、8、3」是我用昨天推出的 Bayer 多工矩陣算的(各欄平方和,亮度欄之外是兩個色度欄),論文沒有印這個矩陣。昨天「符號不影響範數」那句話整句撤回:符號決定這一列的係數和是 0 還是 0.5,也就決定整個矩陣是否可逆性良好。

這個 2×2 設計還有一個身分:它就是 Bean 的青、黃、洋紅、藍圖樣,只是把紅與綠互換標示。RR-08-03 的式 (16) 與 TIP 2011 的圖 7(a) 都印了它。換句話說,論文在 2×2 的框架下「最佳化」出來的結果,是一個 2001 年就有專利的圖樣。這不減損論文的貢獻(它提供的是系統性的設計框架,以及 4×4 的新圖樣),但表示「2×2 的最佳」這句話值得單獨檢查。

三、2×2 的「最佳」只檢查了頂點

RR-08-03 的問題 (15) 是在約束(總透光量為 1、各原色透光率在 0 到 1 之間、某些頻率位置的色度為零)之下,最小化 。可行集合是多面體,而報告自己承認只檢查了頂點。對線性目標來說只看頂點就夠,但範數不是線性函數,內部的點可能更小。

我用論文的條件重做:固定「(0, 1/2) 的色度為零」的頻率結構(忘了這個限制時,我的第一次格點搜尋得到過不合法的 2.31,後來補上限制)。可行的 2×2 圖樣落在一條一維族上,令參數為 ,則 ,四格依序是 、、、。程式如下,可以直接貼上跑:

import numpy as np
def T(a):
    b = (1 + a) / 2
    p0 = np.array([0, a, 1 - a]); p1 = np.array([b, 1 - b, 0])
    p2 = np.array([a, 0, 1 - a]); p3 = np.array([1 - b, b, 0])
    return np.array([(p0 + p1 + p2 + p3) / 4,
                     (p0 - p1 + p2 - p3) / 4,
                     (p0 - p1 - p2 + p3) / 4])
for a in (0.3, 0.4, 0.46, 0.5, 0.6):
    D = np.linalg.inv(T(a))
    print(a, round(np.linalg.norm(D, 2), 3), round(np.linalg.norm(D, 'fro'), 3))
0.304.7145.538
0.403.5364.790
0.463.0744.665
0.503.3354.690
0.604.2735.153

2 範數的最小值約 3.069(在 ),Frobenius 範數的最小值約 4.66(在 附近, 約 21.7 至 21.8)。和論文選的 4.0 與 5.0 相比,2 範數低 23%,Frobenius 範數低 7%。我另外在整個六維可行多面體裡做了隨機取樣加爬山,沒有找到比這條一維族更低的點,但這不是證明;我也無法確認自己重建的限制條件與報告完全相同,因為公式的文字層讀不到。可以確定的只有一件事:在我重建的條件下,頂點不是最小值。

另一個對照:Bayer 在自己的頻率結構(同一個色度重複出現兩次)內已經是範數最優,2.828 與 3.742;那個結構的可行頂點是 RGB 的六種排列。所以「Bayer 的 比所有新圖樣都小」不是巧合,而是結構決定的,論文的 Table II 也是這樣呈現。

四、範數小就畫質好嗎:2×2 內的模擬

昨天的結論是 Table III 沒有檢驗 P3(最小化 )。今天做一個最小的檢驗,範圍刻意縮到「同樣是 2×2、同樣是總透光量 1」的三類圖樣:Bayer( = 14)、我找到的數值最佳(約 21.8)、Bean(25)。設定如下:

  • 資料:Kodak 24 張影像(來源見延伸學習),把每張當成真實的三通道影像,依圖樣做每格的透光加權後加總成單通道感測值,再加零均值高斯雜訊(標準差 0、5、10、20,以 0 到 255 為尺度)。
  • 重建:對四個相位各自做一個 7×7 視窗的線性最小平方(含常數項、嶺迴歸係數 0.001)估計,每個相位每張影像取 5000 個樣本;24 張影像分成前後各 12 張,用一半訓練、另一半測試(兩折交叉)。
  • 指標:每張的 MSE,再轉成 CPSNR;圖樣的三個原色通道指派有 6 種排列(Bean 與最佳各 6 種),結果都列出範圍。
雜訊 σBayer最佳(6 種排列的平均)Bean(6 種排列的平均)
037.3636.83(36.72–36.92)36.41(36.25–36.57)
533.3932.7732.54
1030.0529.35(29.34–29.37)29.18(29.14–29.24)
2026.6725.8325.70

數字是 10·log10(24 張平均 MSE),單位 dB。每個 σ 下的排序都是 Bayer 高於最佳、最佳高於 Bean,和 的排序 14、21.8、25 一致。逐張配對看穩固度:最佳對 Bean 在 σ=0 平均高 0.35 dB(標準誤 0.04),24 張裡 23 張勝;σ=10 平均高 0.13 dB(標準誤 0.02),同樣 23 張勝。Bayer 對最佳在 σ=10 是 24 張全勝,平均高 0.71 dB(標準誤 0.06);在 σ=0 則是 10 勝 14 負,平均 −0.11 dB(標準誤 0.25),也就是無雜訊時差距落在誤差內,而前面的 37.36 對 36.83 是由 MSE 平均得到的,逐張的算術平均差比較小(38.18 對 38.04)。

這個結果支持「在 2×2 同一族之內, 是雜訊下畫質的有用指標」,而且「有雜訊時範數更重要」符合直覺(雜訊被 放大)。限制要講清楚:重建是我自己的逐相位線性估計,不是論文的通用解馬賽克演算法;只有三類圖樣,等於三個資料點;雜訊是加在總透光量加總之後的感測值,沒有模擬真實感光元件的散粒雜訊(那和透光量有關,而這一組圖樣的總透光量都是 1,所以這個簡化在這裡剛好不偏袒任何一方)。所以這只是「沒有被推翻」,不是「被證明」。

五、RR-08-03 的 Table 2:舊資料,新平均

RR-08-03 的 Table 2 是 24 張影像乘 13 欄的 CPSNR。我把它轉成數字,發現平均列等於 24 個 dB 的算術平均(13 欄的差都是 0.00)。TIP 2011 的 Table III 改成 10·log10(平均 MSE)。兩篇之間有 6 欄的 MSE 平均和 2011 年的數字吻合,所以資料是沿用的,只有平均的方式換了。這點延續了昨天的結論:不同論文的「Mean」列不能直接比。

2008 年報告裡,CFA4b 加自適應解馬賽克對上 HW1(Hirakawa–Wolfe 的 ICIP 2007 圖樣)平均高 0.27 dB,逐張 18 勝 6 負,配對 t 值 2.78,逐一拿掉一張影像的結果在 0.21 到 0.30 之間。到了 2011 年,對 HWpA 的差距縮到 0.034 dB(MSE 平均)。這兩個差距不能直接比較:一個是 2008 年的 HW1、算術平均,一個是 2011 年的 HWpA、MSE 平均。HWpA 是不是 HW1 的重新縮放版,報告沒有答案,因為報告只有 HW1 到 HW3,沒有 Condat 圖樣。所以昨天的待辦第一項(HWpA 與 Condat 有沒有被縮放)仍然是「未確認」,但我知道了答案不在 RR-08-03 裡。

六、仍未解決

  • Table II 的「Sony」欄對應哪個圖樣:Yamanaka 圖樣算出來是 2.828 與 3.742,和 Bayer 相同,不是 Sony 欄的 2 與 3.303。
  • TIP 2011 最後採用的目標是否就是問題 (15),或在其上加了其他項:RR-08-03 的伴隨報告 RR-08-04 我沒有打開。
  • 2008 年 Table 2 前兩欄(Bayer 的兩種解馬賽克)的順序,轉錄時無法確認。
  • CFA23 的係數和沒有對上(虛部加總 0.011),我沒有採用。
  • Condat 的 σ≈3.7 翻轉:程式碼取不到,這個待辦順延。

🧠 記

每一列的係數和是第一道檢查。 總透光量為 1 的設計,亮度列係數和必為 1,色度列必為 0。昨天的錯就是漏了這一步;以後任何抄來的 矩陣先算列和,再算反矩陣。

「最佳」要問清楚是哪個集合裡的最佳。 RR-08-03 的 2×2 最佳只檢查了頂點,而範數在內部可以更小(2 範數 3.069 對 4.0)。Bayer 的最佳則是在它自己的頻率結構裡,結構選定之後才談最佳。

「補色」在這篇論文裡是色相,不是光量。 青是 (0, 1/2, 1/2),只有真實青色濾片一半的光。

同一欄標題下的數字,先確認平均方式。 算術平均 dB 與 10·log10(平均 MSE) 可以差到 0.8 dB(Bayer 自適應版 39.47 對 38.66)。

範數與畫質的對應,只在同一個圖樣族內檢驗過。 我的模擬只在 2×2、總透光量 1 之內驗證了排序一致,不能推到 4×4 或不同頻率結構。

✍️ 實踐

  1. 把上面的 numpy 程式貼上去,把 掃成 0.30 到 0.60 的細格點,畫出 與 的曲線,確認兩個最小值的位置不同(0.461 與 0.468)。預測:兩條曲線都是單谷。
  2. 在 T(a) 上加一行檢查 T(a).sum(axis=1),確認結果是 (1, 0, 0);再把第三列的符號手動翻成正的,看列和變成什麼、範數變成多少(應該回到 0.5 與 7.557)。
  3. 找一張自己的 RAW(或 Kodak 影像),把 Bayer 與 的圖樣各模擬一次,雜訊 σ 取 0 與 10,每個 σ 各跑一次,用 10·log10(MSE) 比較;預測:σ=10 時 Bayer 領先約 0.7 dB,σ=0 時差距落在誤差內。
  4. 明天的待辦:找 Condat 的 Matlab 程式碼(或用論文給的 2×3 六色圖樣自己重寫,在 掃描 0 到 10),檢查 σ≈3.7 是否真的出現排名翻轉;順便試著取得 RR-08-04。

🔗 延伸學習

💬 問 AI

我在研究彩色濾光陣列(CFA)的頻域設計。已知:三個原色圖樣逐格加總為 1(每格總透光量為 1),頻率結構矩陣 T 的亮度列係數和為 1、色度列係數和為 0,解多工矩陣 D = T 的反矩陣,設計準則是最小化 ‖D‖。
(1) 2×2 圖樣在「(0, 1/2) 位置的色度為零」的限制下是一條一維參數族:四格依序為 (0,a,1-a)、(b,1-b,0)、(a,0,1-a)、(1-b,b,0),其中 b=(1+a)/2。請幫我推導 T(a) 與 D(a) 的解析式,並求 ‖D‖_2 與 ‖D‖_F 的最小值與位置(我數值得到約 3.069(a≈0.461)與約 4.66(a≈0.468)),指出我可能漏掉的條件。
(2) 解釋為什麼在加性高斯雜訊下,‖D‖_F^2 的排序(Bayer 14、這個最佳 21.8、Bean 25)可以預測 CPSNR 的排序;並指出在什麼情況下這個預測會失效(例如雜訊與透光量相關、解馬賽克演算法有非線性、色度被低通)。
(3) 若要把 Hirakawa–Wolfe 與 Condat 的圖樣放進「每格總透光量為 1」的框架比較,應該如何縮放?請列出每一步的假設。
請不要憑空給數字;不確定的地方請明講。