【轉】 ISP-鏡頭陰影校正（LSC）

轉自：https://blog.csdn.net/xiaoyouck/article/details/77206505

介紹
鏡頭陰影校正（Lens Shading Correction）是為了解決由於lens的光學特性，由於鏡頭對於光學折射不均勻導致的鏡頭周圍出現陰影的情況。

shading可以細分為luma shading和color shading：

luma shading：
由於Lens的光學特性，Sensor影像區的邊緣區域接收的光強比中心小，所造成的中心和四角亮度不一致的現象。鏡頭本身就是一個凸透鏡，由於凸透鏡原理，中心的感光必然比周邊多。如圖所示：

chrom/color shading：
由於各種顏色的波長不同，經過了透鏡的折射，折射的角度也不一樣，因此會造成color shading的現象，這也是為什么太陽光經過三棱鏡可以呈現彩虹的效果。如圖所示：

 

此外，還有CRA的原因會導致shading現象的出現，這里不再贅述，這里推薦《What’s CRA》這篇文章，詳細講述了由於鏡頭的CRA帶來的shading。

影響
luma shading：會造成圖像邊角偏暗，就是所謂的暗角。

color shading：中心和四周顏色不一致，體現出來一般為中心或者四周偏色。如圖所示：

校正
lens shading的校正是分別對於bayer的四個通道進行校正，每個通道的校正過程是相對獨立的過程。

考慮到芯片設計的成本，因此一般情況下不會存儲整幅圖像的lut，目前主流的都是存儲128*128個點的增益，利用雙線性插值的方法計算每個pixel的增益。

算法
由於條件限制，圖像僅用於算法驗證，不做圖像質量評判標准
這里寫了一個shading的算法，將圖像分為16x16的方塊，求取每個交點的增益值，對平面進行四次方擬合，分別計算了luma shading 和 chrom shading，先計算出來一個lut用於存儲，校正的世行通過對這個lut進行雙線性插值得到每個pixel的值乘以原本像素點。

16x16的分塊並非固定，可以對塊的大小進行調整，比如中心塊偏大，靠近邊緣的方塊變小，這些都是可以自定義的，本算法由於做演示使用，故不做其他功能。如圖所示：

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

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

效果展示：
實驗條件有限，圖片有水波紋，僅用於理解算法

original image:

luma shading

 

chroma shading:

 

luma shading + chroma shading:

 

tuning
LSC的tuning一定要把校正圖采集好，一般情況下raw圖的G通道中心亮度在8bit的70%~80%之間，由於在不同色溫情況下是經過插值的，因此需要校正多個光源，一般情況下TL84、D65、A光源下進行校正。將得到的LUT寫入RAM中即可
注意：采集的raw圖不要有filcker。

LSC強度一般是可調的，由於圖像邊角的增益會很大，因此在高倍gain下，可以把強度給降低，防止圖像邊角噪聲壓不住的情況。

由於各個平台不同，這里不做詳細介紹，想到再補充。