本文是isp算法学习的第五篇章,主要是讲述DPC原理和算法。

DPC-Defect Pixel Correction,即为坏点检测

比较便宜的Sensor,暗处拍摄,亮点会比较多。

1产生原因

通常都会是工艺方面为主,使用时间较长也会使得像素点异常

  • 感光元件芯片自身工艺技术瑕疵造成
  • 光线采集存在缺陷
  • 制造商产品差异

2校正方法

2.1离线方式

静态校正

通常由Sensor产商在生产后进行标定,把所有坏点的坐标位置记录下来,然后校正的时候直接通过查表方式找到坏点进行校正。

2.2在线校正

动态校正

在ISP算法中通过特殊算法判断一个点是否为坏点,如果是坏点就进行校正处理,其他不做处理。

2.3关于静态校正局限性

  • 这种标定方式需要消耗人力物力,而且存储这些数据对硬件也有一定要求
  • 随着产品使用时间地增加,坏点地数目会增加,而厂商地数据没法标定这些新增地坏点
  • 有些坏点表现为正常情况下的光电特性正常,只是在长时间使用或者高ISO的情况下才会变为坏点,这部分坏点正常使用时不用当作坏点处理,毕竟插值的数据和真是数据还是有一定的差异

基于以上几点,静态坏点检测方法不进行讨论,而且实际应用的时候很多时候静态坏点检测也是disable了。以下讨论几种常见的检测算法。

3方法实现

3.1PINTO算法

主要思想是坏点往往是在一个邻域中的极亮或极暗点,所以以一个5x5的邻域为检测区域

如图所示5x5的邻域内同意颜色通道相对于中心像素都有8个临近像素。那么矫正按以下步骤操作:

  1. 计算中心像素与周围八个像素值的差;
  2. 判断八个差值是否都为正值或者都为负值;
  3. 如果有的为正有的为负,那么就为正常值,否则进行下一步;
  4. 设置一个阈值,如果八个差值的绝对值都查过阈值,那么就判断为坏点;
  5. 判断为坏点后就用八个临近的像素值的中位值来替换当前的像素值;

3.2梯度法

针对三个通道都使用这一种窗口进行检测,具体步骤如下:

  1. 计算四个方向的梯度:

    水平方向三个二阶梯度:

    Dh1 = |P1+P3-2P2|

    Dh2 = |P4+P5-2Pc|

    Dh3 = |P6+P8-2P7|

    竖直方向三个二阶梯度:

    Dv1 = |P1+P6-2P4|

    Dv2 = |P2+P7-2Pc|

    Dv3 = |P3+P8-2P5|

    45°三个二阶梯度:

    D45_1 = 2|P4-P2|

    D45_2 =|P3+P6-2Pc|

    D45_3 = 2|P7-P5|;

    135°三个二阶梯度

    D135_1 = 2|P2-P5|

    D135_2 =|P1+P8-2Pc|

    D135_3 = 2|P7-P4|

  2. 取出各个方向梯度绝对值的中值 $$ median_D_h = median(D_{h1},D_{h2},D_{h3}) $$ 同理求出其他三个方向的中值 $$ median_D_v = median(D_{v1},D_{v2},D_{v3}) \
    median_D_{45} = median(D_{45_1},D_{45_2},D_{45_3}) \
    median_D_{135} = median(D_{135_1},D_{135_2},D_{135_3})
    $$

  3. 求出四个中值的最小值作为边缘方向: $$ min_{grad} = min(median_D_h,median_D_v,median_D_{45},median_D_{135}); $$

  4. 坏点判定方法1(水平和竖直方向)

    如果最小梯度方向为水平或者竖直,若过Pc那个梯度的绝对值大于同方向的另外两个梯度绝对值和的4倍,则Pc为坏点;

    如果是水平方向 $$ |D_{h2}| > 4(|D_{h1}|+|D_{h3}|) $$ 如果是竖直方向 $$ |D_{v2}| > 4(|D_{v1}|+|D_{v3}|) $$

  5. 坏点判定方法2(45°和135°方向)

    如果是45°方向,则计算135°三个梯度绝对值两两之差的绝对值的和 $$ D_{135_sum} = |D_{135_1} - D_{135_2}| \+ |D_{135_1} - D_{135_3}| + |D_{135_2} - D_{135_3}| $$ 1)如果D135_sum小于100,Pc为坏点 $$ D_{45_2}>3(D_{45_1}+D_{45_3}) && \ D_{135_2}>3(D_{135_1}+D_{135_3}) $$ 2)如果D135_sum大于100,Pc为坏点 $$ D_{45_2}>3(D_{45_1}+D_{45_3}) $$

  6. 135°方向和45°相反的方向计算和判断即可;

  7. 坏点判定方法3

    为减少漏判,当Pc小于15且周围点都大于Pc40以上,则Pc为坏点。如果Pc大于230,且周围的点都小于Pc30以下,则该点为坏点

  8. 水平竖直校正。边缘为水平方向,且判断为坏点,如过|P4-Pc|<|Pc-P5|则Pc更靠近P4,根据同一颜色通道亮度的渐变性可以推导出 $$ ouput=P_4+(P_2+P_7-P_1-P_6)/2 $$

    否则 $$ ouput=P_6+(P_2+P_7-P_3-P_9)/2 $$

  9. 如果为竖直方向可以参考水平方向求出

  10. 45°和135°校正。边缘为45°,如果|P3-Pc|<|P6-Pc|则根据同一原则 $$ output=P_3+(P_4+P_7-P_2-P_5)/2 $$

    否则为 $$ output=P_7+(P_2+P_5-P_7-P_3)/2 $$

  11. 边缘为135°则按照45°的方式反过来计算即可。

3.3DPC和demosaic结合法

这种方法就是将DPC算法和去马赛则算法融合到一块而

思路

  • 先对bayer图像插值成全彩图像

  • 对每个点进行坏点检测,检测的时候只用同一同到的像素检测,而且每个像素点用的颜色通道即为该点raw图对应的颜色通道

    换句话说就是RGGB图像,对于这四个点,第一个点我只通过周围插值出的R来检测R是否为坏点,其他两个同样这么处理,这样就保证和raw处理的时候一样,只用处理MxN个像素点就行,

  • 如果该点是坏点就对其进行矫正并重新对齐进行颜色插值,那么执行完这整个流程,坏点被矫正了而且直接输出了全彩图,

3.4行检测法

上图第一行保存的是上一行的信息点,是否为候选点;第二行为当前的信息点。

这类算法检测和矫正的时候并不使用整帧图像,而是通过几行数据来处理,这样对硬件的buffer要求就没那么高,而且CMOS也是行扫面的,一次正好手机一行数据。具体的操作如下:

  • 缓存一行数据,在同一行中通过比较同一通道相邻的数据的插值,如果待检测的点同时大于相邻点一个阈值,或者同时小于相邻点一个阈值,那么这个点就是候选坏点;
  • 有的算法还会利用buffer缓存上一行的数据,然后判断这个坏点在当前行周围有没有候选坏点,并判断上一行对应位置的点是否为候选坏点,如果都不是,那么当前点就是真实坏点,就通过周围点矫正,如果周围有候选坏点,那么就判断为非真是坏点,不用矫正;
  • 在上一条提到的缓存的上一行数据,有的算法中不是缓存上一行每个像素的信息,而是上一行经过判断后的属性值,比如上一行存储的是每个点是否为候选坏点或者坏点,那么每个点就可以用0或者1来表示,那么上一行每个像素只需要1bit的数据量来存储,这样就能进一步降低对硬件的要求了。

4代码

4.1文件信息配置并读取

 1% --------raw parameters-------
 2filePath = 'images/HisiRAW_4208x3120_8bits_RGGB.raw';
 3bayerFormat = 'RGGB';
 4bayerBits = 8;
 5row = 4208;
 6col = 3120;
 7% -----------------------------
 8
 9rawData = readRaw(filePath, bayerBits, row, col);
10[height, width, channel] = size(rawData);

4.2对边缘进行处理

因为处理边缘的时候没有邻域了,需要将图片扩充一下,以便对边缘进行处理

1img_expand = expandRaw(rawData, expandNum);

4.3判断坏点

判断周边8个点差值的绝对值是否大于阈值TH,则为坏点,用八个点中位值校正。

 1function correctP = judgeDefectPixel(aroundP, currentP, Th)
 2% judgeDefectPixel.m    correct the curren pixel
 3%   Input:
 4%       aroundP    	the pixel around the current pixel 
 5%       currentP    the value of current pixel
 6%       Th          the threshold of the defect pixel
 7%   Output:
 8%       correctP    the corrected value of the pixel
 9%   Instructions:
10%       author:     wtzhu
11%       e-mail:     wtzhu_13@163.com
12% Last Modified by wtzhu v1.0 2021-07-16
13% Note: 
14    % get the median value of the around list
15    medianV = median(aroundP);
16    % get the difference between the around pixel and the current pixel
17    diff = aroundP - ones(1, numel(aroundP)) * currentP;
18    % if all difference bigger than 0 or all smaller than 0 and all abs of the diff are bigger than Th, that pixel is
19    % a defect pixel and replace it with the median;
20    if (nnz(diff > 0) ==  numel(aroundP)) || (nnz(diff < 0) ==  numel(aroundP))
21        if length(find((abs(diff)>Th)==1)) == numel(aroundP)
22            correctP = medianV;
23        else
24            correctP = currentP;
25        end
26    else
27        correctP = currentP;
28    end
29end

4.4归一化四通道遍历

 1disImg = zeros(height, width);
 2for i = expandNum+1 : 2 : height+expandNum
 3    for j = expandNum+1 : 2 : width+expandNum
 4        % R
 5        % get the pixel around the current R pixel
 6        around_R_pixel = [img_expand(i-2, j-2) img_expand(i-2, j) img_expand(i-2, j+2) img_expand(i, j-2) img_expand(i, j+2) img_expand(i+2, j-2) img_expand(i+2, j) img_expand(i+2, j+2)];
 7        disImg(i-expandNum, j-expandNum) = judgeDefectPixel(around_R_pixel, img_expand(i, j), Th);
 8        % Gr
 9        % get the pixel around the current Gr pixel
10        around_Gr_pixel = [img_expand(i-1, j) img_expand(i-2, j+1) img_expand(i-1, j+2)  img_expand(i, j-1) img_expand(i, j+3) img_expand(i+1, j) img_expand(i+2, j+1) img_expand(i+1, j+2)];
11        disImg(i-expandNum, j-expandNum+1) = judgeDefectPixel(around_Gr_pixel, img_expand(i, j+1), Th);
12        % B
13        % get the pixel around the current B pixel
14        around_B_pixel = [img_expand(i-1, j-1) img_expand(i-1, j+1) img_expand(i-1, j+3) img_expand(i+1, j-1) img_expand(i+1, j+3) img_expand(i+3, j-1) img_expand(i+3, j+1) img_expand(i+3, j+3)];
15        disImg(i-expandNum+1, j-expandNum+1) = judgeDefectPixel(around_B_pixel, img_expand(i+1, j+1), Th);
16        % Gb
17        % get the pixel around the current Gb pixel
18        around_Gb_pixel = [img_expand(i, j-1) img_expand(i-1, j) img_expand(i, j+1) img_expand(i+1, j-2) img_expand(i+1, j+2) img_expand(i+2, j-1) img_expand(i+3, j) img_expand(i+2, j+1)];
19        disImg(i-expandNum+1, j-expandNum) = judgeDefectPixel(around_Gb_pixel, img_expand(i+1, j), Th);
20    end
21end
22figure();
23imshow(rawData);title('org');
24figure();
25imshow(uint8(disImg));title('corrected');
26disp(['cost time£º',num2str(toc)])