what is the LSC?
lens shading 分为：Y-shading , color shading。
在讲LSC之前，我们先来理解一个重要的术语--CRA(Chief ray angle)。
CRA：分为lens cra , sensor cra两种。
1.1什么是lens CRA?
定义：最大像高处的主光线与光轴的倾角。备注：最大像高定义为 像素响应降低为零度角像素响应的80%的像素点。
因为Lens CRA 的存在，所以就出现了FOV （filed of view) 所谓的视角的概念。视角大小=2*CRA。
下图为lens CRA示意图。
1.2 什么是sensor的CRA
sensor的CRA 说白了是sensor与之配对的micro lens的CRA.
圈1 表示lens折射过来的光，圈2表示micro lens再次折射的光。可以看出，sensor的光电转换二极管接收到的光并不是 lens直接折射过来的，而是经过microlens再次折射过来的。所以你明白了什么呢？！ 这也就是为什么图像边缘会变暗 的原因，光电转换二极管并不会全部接收到光能。这是LensShading中Y-shading的主要原因。
一般情况下，当lens_cra < micro_cra的时候，在sensor传感器的边缘只会出现Y-shading。
倘若，lens_cra > micro_cra, sensor传感器可能会因为无法有效捕捉某频率的光，就会出现偏色的情况，即Color- shading。比如microlens 对R ,G,B的折射率不同，在R_pixel上接收了很小的能量，而在B_pixel上接受了相对较多的能 量。
如何处理或者叫减轻shading带来的图像不均问题呢？
1.lens cra要小于micro cra,这样既可以减缓y-shading,又可以大大避免color-shading.
2.因为镜头是圆的，画面的对角线是最接近镜头成像直径的，所以sensor有效区域最大尺寸是对角线尺寸
3.硬件矫正
4.软件矫正（下面会细讲）
how to LSC with software?
一般模组厂会根据矫正工具将图像划分为17*13(qualcomm),15*15(MTK),33*25(Samsung)等区间，计算每一块的四通道平均值，然后对每个通道与中心块的四通道平均值做除法，得到每个点的增益倍数，整个lut存储在otp中，手机端就可以根据lut进行插值计算每个pixel的增益。
一般用双线性差值或者多项式拟合的方法，确定每个pixel的增益。
原始图片的rgb通道(节选自网络)：
y-shading 增益之后，基本与中心数值在同一个平面(节选自网络)：
color-shading 增益之后，各区域的r/g，b/g均值相同(节选自网络)：
上面三张图片为了让大家理解，节选自网络。本人也模拟了shading矫正过程，但由于colorshading不严重所以矫正效果不明显。下面是我的矫正过程：
1.原始图：
2.三通道三维图
3.(17*13block)r,g,b平均
4.colorgain
5、lensgain=ygain+colorgain
6. 还原后三通道值
7.还原后图片
附上我的代码：
clc
clear all;
close all;
BLOCK_H = 17
BLOCK_V = 13
img = imread('1.bmp');
r= medfilt2(img(:,:,1));%中值滤波抑制边缘凸起噪声
g= medfilt2(img(:,:,2));
b= medfilt2(img(:,:,3));
figure(1)
[rows,cols,channels] = size(img);
width = 1:cols; height = 1:rows; [X1,Y1] = meshgrid(width,height);
mesh(X1,Y1,im2double(r)); hold on; mesh(X1,Y1,im2double(g)); hold on; mesh(X1,Y1,im2double(b));
subSampleRGB=ones(BLOCK_V+1,BLOCK_H+1,3);%6*6*3 gainlut
strideX=floor(cols/BLOCK_H)
strideY=floor(rows/BLOCK_V)
for i=1:BLOCK_H+1%cols
for j=1:BLOCK_V+1%rows
for k=1:3%channels
if i==1
if j==1%left_top
subSampleRGB(j,i,k)=mean2( img(1:floor(strideY*0.5), 1:floor(strideX*0.5),k) );
elseif j==BLOCK_V+1 %left_bottom
subSampleRGB(j,i,k)=mean2( img(floor((BLOCK_V-0.5)*strideY):BLOCK_V*strideY, 1:floor(strideX*0.5),k) );
else %left_middle
subSampleRGB(j,i,k)=mean2( img(floor((j-1-0.5)*strideY):floor((j-1+0.5)*strideY), 1:floor(strideX*0.5),k) );
end
elseif i==BLOCK_H+1
if j==1%right_top
subSampleRGB(j,i,k)=mean2( img(1:floor(strideY*0.5), floor(strideX*(BLOCK_H-0.5)):BLOCK_H*strideX,k) );
elseif j==BLOCK_V+1 %right_bottom
subSampleRGB(j,i,k)=mean2( img(floor((BLOCK_V-0.5)*strideY):BLOCK_V*strideY, floor(strideX*(BLOCK_H-0.5)):BLOCK_H*strideX,k) );
else %right_middle
subSampleRGB(j,i,k)=mean2( img(floor((j-1-0.5)*strideY):floor((j-1+0.5)*strideY), floor(strideX*(BLOCK_H-0.5)):BLOCK_H*strideX,k) );
end
elseif j==1%top without corner
if i~=1&&i~=BLOCK_H+1
subSampleRGB(j,i,k)=mean2( img(1:floor(0.5*strideY), floor(strideX*(i-1-0.5)):floor(strideX*(i-1+0.5)),k) );
end
elseif j==BLOCK_V+1%bottom without corner
if i~=1&&i~=BLOCK_H+1
subSampleRGB(j,i,k)=mean2( img(floor((BLOCK_V-0.5)*strideY):BLOCK_V*strideY, floor(strideX*(i-1-0.5)):floor(strideX*(i-1+0.5)),k));
end
else
if i~=1&&i~=BLOCK_H+1&&j~=1&&j~=BLOCK_V+1 %center area
subSampleRGB(j,i,k)=mean2( img(floor((j-1-0.5)*strideY):floor((j-1+0.5)*strideY), floor(strideX*(i-1-0.5)):floor(strideX*(i-1+0.5)),k));
end
end
end
end
end
figure(2)
[A,B] = meshgrid(1:BLOCK_H+1,1:BLOCK_V+1);
mesh(A,B,subSampleRGB(:,:,1));
hold on
mesh(A,B,subSampleRGB(:,:,2));
hold on
mesh(A,B,subSampleRGB(:,:,3));
%Gainlut luma gain (Y-shading gain)
gainlut=zeros(BLOCK_V+1,BLOCK_H+1,3);% gainlut
for i = 1:BLOCK_H+1
for j =1:BLOCK_V+1
for k = 1:3
gainlut(j,i,k)=mean2( subSampleRGB(floor(BLOCK_V/2+1) : floor(BLOCK_V/2+1)+1, floor(BLOCK_H/2+1) : floor(BLOCK_H/2+1)+1,k)) / subSampleRGB(j,i,k);
end
end
end
%Gainlut chroma gain (colorshading gain)
gainlut_c = zeros(BLOCK_V+1,BLOCK_H+1,2);
center_RoverG = mean2( subSampleRGB(floor(BLOCK_V/2+1) : floor(BLOCK_V/2+1)+1, floor(BLOCK_H/2+1) : floor(BLOCK_H/2+1)+1,1)) / mean2( subSampleRGB(floor(BLOCK_V/2+1) : floor(BLOCK_V/2+1)+1, floor(BLOCK_H/2+1) : floor(BLOCK_H/2+1)+1,2));
center_BoverG = mean2( subSampleRGB(floor(BLOCK_V/2+1) : floor(BLOCK_V/2+1)+1, floor(BLOCK_H/2+1) : floor(BLOCK_H/2+1)+1,3)) / mean2( subSampleRGB(floor(BLOCK_V/2+1) : floor(BLOCK_V/2+1)+1, floor(BLOCK_H/2+1) : floor(BLOCK_H/2+1)+1,2));
for i = 1:BLOCK_H+1
for j =1:BLOCK_V+1
RoverG = subSampleRGB(j,i,1)/subSampleRGB(j,i,2);
BoverG = subSampleRGB(j,i,3)/subSampleRGB(j,i,2);
gainlut_c(j,i,1) = center_RoverG/RoverG-1;
gainlut_c(j,i,2) = center_BoverG/BoverG-1;
end
end
%chroma+luma gain
gainlut_full = zeros(BLOCK_V+1,BLOCK_H+1,2);
for i = 1:BLOCK_H+1
for j =1:BLOCK_V+1
gainlut_full(:,:,1) = gainlut(:,:,1) + gainlut_c(j,i,1);
gainlut_full(:,:,2) = gainlut(:,:,2) + gainlut_c(j,i,2);
end
end
%(BLOCK_H+1)*(BLOCK_V+1) colorgain
figure(3)
[C,D] = meshgrid(1:BLOCK_H+1,1:BLOCK_V+1);
mesh(C,D,gainlut_c(:,:,1))
hold on
mesh(C,D,gainlut_c(:,:,2))
figure(4)
mesh(C,D,gainlut_full(:,:,1))
hold on
mesh(C,D,gainlut_full(:,:,2))
hold on
mesh(C,D,gainlut(:,:,2))
%双线性插值FullGainLut
x = 1:strideX:cols+1;
x(end) = cols;
y = 1:strideY:rows+1;
y(end) = rows;
xitp = 1:cols;
yitp = 1:rows;
[Xitp,Yitp]=meshgrid(xitp,yitp);
rgain=interp2(x,y,gainlut(:,:,1),Xitp,Yitp);
ggain=interp2(x,y,gainlut(:,:,2),Xitp,Yitp);
bgain=interp2(x,y,gainlut(:,:,3),Xitp,Yitp);
%shading r,b,g gain
figure(5)
r=im2uint8( (im2double(r).*rgain) );
mesh(X1,Y1,im2double(r));
b=im2uint8( (im2double(b).*bgain) );
hold on
mesh(X1,Y1,im2double(b));
g=im2uint8( (im2double(g).*ggain) );
hold on
mesh(X1,Y1,im2double(g));
%shading corrected img
figure(6)
lscgain = cat(3,rgain,ggain,bgain);%merge三通道gain
imgDoubleType = im2double(img);
imgDoubleType = imgDoubleType.*lscgain;
imgUint = im2uint8(imgDoubleType);
imshow(imgUint)