ISP-镜头阴影校正（LSC）

偏执的太偏执、 2023-06-11 09:26 8阅读 0赞

本文链接： [https://blog.csdn.net/xiaoyouck/article/details/77206505][https_blog.csdn.net_xiaoyouck_article_details_77206505]

# 概述 #

## 介绍 ##

镜头阴影校正（Lens Shading Correction）是为了解决由于lens的光学特性，由于镜头对于光学折射不均匀导致的镜头周围出现阴影的情况。

shading可以细分为luma shading和color shading：

*  luma shading：  
    由于Lens的光学特性，Sensor影像区的边缘区域接收的光强比中心小，所造成的中心和四角亮度不一致的现象。镜头本身就是一个凸透镜，由于凸透镜原理，中心的感光必然比周边多。如图所示：  
    ![这里写图片描述][SouthEast]
 *  chrom/color shading：  
    由于各种颜色的波长不同，经过了透镜的折射，折射的角度也不一样，因此会造成color shading的现象，这也是为什么太阳光经过三棱镜可以呈现彩虹的效果。如图所示：  
    ![这里写图片描述][SouthEast 1]

此外，还有CRA的原因会导致shading现象的出现，这里不再赘述，这里推荐《[What’s CRA][What_s CRA]》这篇文章，详细讲述了由于镜头的CRA带来的shading。

## 影响 ##

luma shading：会造成图像边角偏暗，就是所谓的暗角。  
![这里写图片描述][SouthEast 2]

color shading：中心和四周颜色不一致，体现出来一般为中心或者四周偏色。如图所示：  
![这里写图片描述][SouthEast 3]

## 校正 ##

lens shading的校正是分别对于bayer的四个通道进行校正，每个通道的校正过程是相对独立的过程。

考虑到芯片设计的成本，因此一般情况下不会存储整幅图像的lut，目前主流的都是存储128\*128个点的增益，利用双线性插值的方法计算每个pixel的增益。

### 算法 ###

**由于条件限制，图像仅用于算法验证，不做图像质量评判标准**  
这里写了一个shading的算法，将图像分为16x16的方块，求取每个交点的增益值，对平面进行四次方拟合，分别计算了luma shading 和 chrom shading，先计算出来一个lut用于存储，校正的世行通过对这个lut进行双线性插值得到每个pixel的值乘以原本像素点。

16x16的分块并非固定，可以对块的大小进行调整，比如中心块偏大，靠近边缘的方块变小，这些都是可以自定义的，本算法由于做演示使用，故不做其他功能。如图所示：  
![这里写图片描述][SouthEast 4]

#### code ####

由于代码量较大，这里分别附上一部分算法

shading lut caculate:

function [image_r_gain, image_gr_gain, image_gb_gain, image_b_gain] = ...
    isp_lsc_lut(image_r, image_gr, image_gb, image_b, side_num)
    [height, width] = size(image_r);
    side_y = floor(height/side_num);
    side_x = floor(width/side_num);
    
    % figure,imshow(image_r);
    % hold on;
    % for k=0:side_num
    %     line_x = side_x * k;
    %     line_y = side_y * k;
    %     if(k==side_num && line_y ~= width) line_y = height;end
    %     if(k==side_num && line_x ~= width) line_x = width;end
    %     line([line_x,line_x],[0,height],'Color','red');
    %     line([0,width], [line_y, line_y],'Color','red');
    % %     line(Xd,Yd,'Color','red');
    % end
    % hold off
    
    %% compress resolution
    image_point = zeros(side_num,side_num);
    for i = 0:side_num
        for j = 0:side_num
            x_clip = floor([j*side_x - side_x/2, j*side_x + side_x/2]);
            y_clip = floor([i*side_y - side_y/2, i*side_y + side_y/2]);
            if(i==side_num && y_clip(2) ~= height) y_clip(2) = height;end
            if(j==side_num && x_clip(2) ~= width) x_clip(2) = width;end
            x_clip(x_clip<1) = 1;x_clip(x_clip>width) = width;
            y_clip(y_clip<1) = 1;y_clip(y_clip>height) = height;
            data_r_in = image_r(y_clip(1):y_clip(2), x_clip(1):x_clip(2));
            image_r_point(i+1,j+1) = mean(mean(data_r_in));
            data_gr_in = image_gr(y_clip(1):y_clip(2), x_clip(1):x_clip(2));
            image_gr_point(i+1,j+1) = mean(mean(data_gr_in));
            data_gb_in = image_gb(y_clip(1):y_clip(2), x_clip(1):x_clip(2));
            image_gb_point(i+1,j+1) = mean(mean(data_gb_in));
            data_b_in = image_b(y_clip(1):y_clip(2), x_clip(1):x_clip(2));
            image_b_point(i+1,j+1) = mean(mean(data_b_in));
        end
    end
    
    % figure,imshow(uint8(image_r_point));
    %% caculate lsc luma gain
    for i = 1:side_num+1
        for j = 1:side_num+1
            image_r_luma_gain_point(i,j) = mean2(image_r_point(uint8(side_num/2)-1:uint8(side_num/2)+1, uint8(side_num/2)-1:uint8(side_num/2)+1)) / image_r_point(i,j);
            image_gr_luma_gain_point(i,j) = mean2(image_gr_point(uint8(side_num/2)-1:uint8(side_num/2)+1, uint8(side_num/2)-1:uint8(side_num/2)+1)) / image_gr_point(i,j);
            image_gb_luma_gain_point(i,j) = mean2(image_gb_point(uint8(side_num/2)-1:uint8(side_num/2)+1, uint8(side_num/2)-1:uint8(side_num/2)+1)) / image_gb_point(i,j);
            image_b_luma_gain_point(i,j) = mean2(image_b_point(uint8(side_num/2)-1:uint8(side_num/2)+1, uint8(side_num/2)-1:uint8(side_num/2)+1)) / image_b_point(i,j);
        end
    end
    
     
      
     
      
       1 
       2 
       3 
       4 
       5 
       6 
       7 
       8 
       9 
       10 
       11 
       12 
       13 
       14 
       15 
       16 
       17 
       18 
       19 
       20 
       21 
       22 
       23 
       24 
       25 
       26 
       27 
       28 
       29 
       30 
       31 
       32 
       33 
       34 
       35 
       36 
       37 
       38 
       39 
       40 
       41 
       42 
       43 
       44 
       45 
       46 
       47 
       48 
       49 
       50

bilinear interpolation:

image_r_luma_gain_reshape = reshape(image_r_luma_gain_point, [], 1);
    image_gr_luma_gain_reshape = reshape(image_gr_luma_gain_point, [], 1);
    image_gb_luma_gain_reshape = reshape(image_gb_luma_gain_point, [], 1);
    image_b_luma_gain_reshape = reshape(image_b_luma_gain_point, [], 1);
    for i = 1:17
        for j = 1:17
            x((i-1)*17+j) = i;
            y((i-1)*17+j) = j;
        end
    end
    x=x';
    y=y';
    % scatter3(x,y,image_r_luma_gain_reshape)
    % hold on
    Z=[ones(length(x),1),x,y,x.^2,x.*y,y.^2,x.^3,x.^2.*y,x.*y.^2,y.^3];
    [x y]=meshgrid(1:17,1:17);
    A=Z\image_r_luma_gain_reshape;
    image_r_luma_gain=A(1)+A(2)*x+A(3)*y+A(4)*x.^2+A(5)*x.*y+A(6)*y.^2+A(7)*x.^3+A(8)*x.^2.*y+A(9)*x.*y.^2+A(10)*y.^3;
    A=Z\image_gr_luma_gain_reshape;
    image_gr_luma_gain=A(1)+A(2)*x+A(3)*y+A(4)*x.^2+A(5)*x.*y+A(6)*y.^2+A(7)*x.^3+A(8)*x.^2.*y+A(9)*x.*y.^2+A(10)*y.^3;
    A=Z\image_gb_luma_gain_reshape;
    image_gb_luma_gain=A(1)+A(2)*x+A(3)*y+A(4)*x.^2+A(5)*x.*y+A(6)*y.^2+A(7)*x.^3+A(8)*x.^2.*y+A(9)*x.*y.^2+A(10)*y.^3;
    A=Z\image_b_luma_gain_reshape;
    image_b_luma_gain=A(1)+A(2)*x+A(3)*y+A(4)*x.^2+A(5)*x.*y+A(6)*y.^2+A(7)*x.^3+A(8)*x.^2.*y+A(9)*x.*y.^2+A(10)*y.^3;
    % surf(x,y,image_r_luma_gain)
    % hold on 
    % surf(x,y,image_r_luma_gain_point)
    
    
    %% calulate lsc chroma gain
    for i = 1:side_num+1
        for j = 1:side_num+1
            image_r_chroma_gain(i,j) = image_r_luma_gain(i,j) - image_r_luma_gain_point(i,j);
            image_gr_chroma_gain(i,j) = image_gr_luma_gain(i,j) - image_gr_luma_gain_point(i,j);
            image_gb_chroma_gain(i,j) = image_gb_luma_gain(i,j) - image_gb_luma_gain_point(i,j);
            image_b_chroma_gain(i,j) = image_b_luma_gain(i,j) - image_b_luma_gain_point(i,j);
        end
    end
    %% caculate lsc result gain
    image_r_gain = image_r_luma_gain - image_r_chroma_gain;
    image_gr_gain = image_gr_luma_gain - image_gr_chroma_gain;
    image_gb_gain = image_gb_luma_gain - image_gb_chroma_gain;
    image_b_gain = image_b_luma_gain - image_b_chroma_gain;
    
    
    
    function image_gain_lut = lsc_data_gain_interpolation(image_gain, height, width, side_num)
    side_y_ori = floor(height/side_num);
    side_x_ori = floor(width/side_num);
    k = 0;
    l = 0;
    [gain_height, gain_width] = size(image_gain);
    for i = 1:gain_height-1
        for j = 1:gain_width-1
            data_gain_11 = image_gain(i, j);
            data_gain_12 = image_gain(i, j+1);
            data_gain_21 = image_gain(i+1, j);
            data_gain_22 = image_gain(i+1, j+1);
            if(j == gain_width-1 && ((j-1)*side_x + l) ~= width) 
                side_x = width - (j-1)*side_x_ori;
            else
                side_x = side_x_ori;
            end
    
            if(i == gain_width-1 && ((i-1)*side_y + k) ~= width)
                side_y = height - (i-1)*side_y_ori;
            else
                side_y = side_y_ori;
            end
    
            for k = 1:side_y
                for l = 1:side_x
                    label_y1 = 1;
                    label_x1 = 1;
                    label_y2 = side_y;
                    label_x2 = side_x;
                    image_gain_lut((i-1)*side_y_ori + k, (j-1)*side_x_ori + l) = ...
                        data_gain_22/(label_x2-label_x1)/(label_y2-label_y1)* ...
                        (l - label_x1) * (k - label_y1) + ...
                        data_gain_21/(label_x2-label_x1)/(label_y2-label_y1)* ...
                        (label_x2 - l) * (k - label_y1) + ...
                        data_gain_12/(label_x2-label_x1)/(label_y2-label_y1)* ...
                        (l - label_x1) * (label_y2 - k) + ...
                        data_gain_11/(label_x2-label_x1)/(label_y2-label_y1)* ...
                        (label_x2 - l) * (label_y2 - k);
                end
            end
    
        end
    end
    end
    
     
      
     
      
       1 
       2 
       3 
       4 
       5 
       6 
       7 
       8 
       9 
       10 
       11 
       12 
       13 
       14 
       15 
       16 
       17 
       18 
       19 
       20 
       21 
       22 
       23 
       24 
       25 
       26 
       27 
       28 
       29 
       30 
       31 
       32 
       33 
       34 
       35 
       36 
       37 
       38 
       39 
       40 
       41 
       42 
       43 
       44 
       45 
       46 
       47 
       48 
       49 
       50 
       51 
       52 
       53 
       54 
       55 
       56 
       57 
       58 
       59 
       60 
       61 
       62 
       63 
       64 
       65 
       66 
       67 
       68 
       69 
       70 
       71 
       72 
       73 
       74 
       75 
       76 
       77 
       78 
       79 
       80 
       81 
       82 
       83 
       84 
       85 
       86 
       87 
       88 
       89 
       90 
       91

#### 效果展示： ####

**实验条件有限，图片有水波纹，仅用于理解算法**

original image:  
![这里写图片描述][SouthEast 5]

luma shading

![这里写图片描述][SouthEast 6]  
![这里写图片描述][SouthEast 7]  
![这里写图片描述][SouthEast 8]

chroma shading  
![这里写图片描述][SouthEast 9]  
![这里写图片描述][SouthEast 10]  
![这里写图片描述][SouthEast 11]

luma shading + chroma shading:  
![这里写图片描述][SouthEast 12]  
![这里写图片描述][SouthEast 13]  
![这里写图片描述][SouthEast 14]

### tuning ###

LSC的tuning一定要把校正图采集好，一般情况下raw图的G通道中心亮度在8bit的70%~80%之间，由于在不同色温情况下是经过插值的，因此需要校正多个光源，一般情况下TL84、D65、A光源下进行校正。将得到的LUT写入RAM中即可  
**注意：**采集的raw图不要有filcker。

LSC强度一般是可调的，由于图像边角的增益会很大，因此在高倍gain下，可以把强度给降低，防止图像边角噪声压不住的情况。

由于各个平台不同，这里不做详细介绍，想到再补充。

[Link 1]: http://creativecommons.org/licenses/by-sa/4.0/
[https_blog.csdn.net_xiaoyouck_article_details_77206505]: https://blog.csdn.net/xiaoyouck/article/details/77206505
[SouthEast]: /images/20230531/7d5cf45aca4640d9814f71f3a955984d.png
[SouthEast 1]: /images/20230531/a182aa93f263457087d09f2f7635341b.png
[What_s CRA]: http://blog.csdn.net/huddheaven/article/details/52300262
[SouthEast 2]: /images/20230531/1188fc625c344e6697b8c4af879e461b.png
[SouthEast 3]: /images/20230531/f6897afea46c40e29dd6df3b36e8f3a4.png
[SouthEast 4]: /images/20230531/3c72877ecac044b89ccd9796bce73fd4.png
[SouthEast 5]: /images/20230531/55f80ba4bbf1441cb6ae28b97fe991eb.png
[SouthEast 6]: /images/20230531/42d4ad7b506f4051a15b4c0e1125964c.png
[SouthEast 7]: /images/20230531/c96208bdac5943b3ba619bb8989f92fb.png
[SouthEast 8]: /images/20230531/04c3d4d3211042fcae5da3c7dc51bb9f.png
[SouthEast 9]: /images/20230531/cfcbbd35946a4652a5127d5c5611fc8f.png
[SouthEast 10]: /images/20230531/b5b8bde9efa0431e9b17678acf75dce3.png
[SouthEast 11]: /images/20230531/4bd7b02b10864ec19ffc1ddad441a1b4.png
[SouthEast 12]: /images/20230531/54bd1284f1eb47ccbf696a6b38541c66.png
[SouthEast 13]: /images/20230531/671a69c8d5e24f53aaba65eb7287b40f.png
[SouthEast 14]: /images/20230531/2fd7daaed8f64bb8a9c07779e63a1808.png