请教高手谁能把这段程序读懂（这是一段联合图像降噪和插值研究的程序）
function [Y]= TLSdemosaic(X,k0,k1)
% G1 R
% B G2
% With respect to the center pixel, the colors are named as follows:
% 2 3 2
% 4 1 4
% 2 3 2
% for example, if the center pixel is B, then 1=B, 2=R, 3,4=G
% control variables
m = 25; % size of the search area (must be odd)
n = 5 ; % size of the image patch (must be odd)
% initialization
TLSindex; % script file that generates indeces
X = TLSpad(X,numM);
Y = zeros(size(X)); % output preallocation
X = TLSlimit(X,11,4);
X = sum(X,3); % input in Bayer pattern
tic
for j=numM+1:size(X,2)-numM % for each column...
if mod(j,2)==1 % GBGBGBGB
for i=numM+1:size(X,1)-numM
x = X(i+indexM,j+indexM);
if mod(i,2)==1 % 1=G1, 2=G2, 3=B, 4=R
Y(i,j,:) = TLSestimate(x,k0,k1,n,M4,M12,M3,IC0,IC,I4,I12,I3);
else % 1=B, 2=R, 3=G1, 4=G2
Y(i,j,:) = TLSestimate(x,k0,k1,n,M2,M34,M1,IC0,IC,I2,I34,I1);
end
end
else % RGRGRGRG
for i=numM+1:size(X,1)-numM
x = X(i+indexM,j+indexM);
if mod(i,2)==0 % 1=G2, 2=G1, 3=R, 4=B
Y(i,j,:) = TLSestimate(x,k0,k1,n,M3,M12,M4,IC0,IC,I3,I12,I4);
else % 1=R, 2=B, 3=G2, 4=G1
Y(i,j,:) = TLSestimate(x,k0,k1,n,M1,M34,M2,IC0,IC,I1,I34,I2);
end
end
end
j
imwrite(Y/255,'tempY.bmp');
end
Y = TLSpad(Y,-numM);
toc
% index_gen.m
% m size of the search area (must be odd)
% n size of the image patch (must be odd)
k = 10;
% indeces for navigating around Bayer pattern
numM = (m-1)/2; % maximum displacement for search
numN = (n-1)/2; % size of the image patch, max distance from center pixel
indexM = -numM:numM; % search area index
indexN = -numN:numN; % image patch index
% indeces used repeatedly in TLSestimate functions
%
% nxn image patch indexing
J = reshape(1:n^2,n,n);
if mod(numN,2)==0
shuffle1 = J(1:2:end,1:2:end); % color 1 downsample
shuffle2 = J(2:2:end,2:2:end); % color 2 downsample
shuffle3 = J(2:2:end,1:2:end); % color 3 downsample
shuffle4 = J(1:2:end,2:2:end); % color 4 downsample
else
shuffle1 = J(2:2:end,2:2:end); % color 1 downsample
shuffle2 = J(1:2:end,1:2:end); % color 2 downsample
shuffle3 = J(1:2:end,2:2:end); % color 3 downsample
shuffle4 = J(2:2:end,1:2:end); % color 4 downsample
end
shuffle = [shuffle1(:);shuffle2(:);shuffle3(:);shuffle4(:)];
% We use x(shuffle,:) to group the same colors together.
% search area indexing (for X(i+indexM,j+indexM))
I = reshape(1:m^2,m,m);
% cross-color image patch
t = I(1+numM+indexN,1+numM+indexN);
IC0 = t(:);
IC = [];
temp = min(numM-numN,8);
if mod(numM,2)==0
for j = -temp:2:temp
for i = -temp:2:temp
t = I(i+1+numM+indexN,j+1+numM+indexN);
IC = [IC t(:)];
end
end
else
for j = -(numM-numN):2:(numM-numN)
for i = -(numM-numN):2:(numM-numN)
t = I(2+numM+i+indexN,2+numM+j+indexN);
IC = [IC t(:)];
end
end
end
% single-color image patches
if mod(numM,2)==0
C1 = I(1:2:end,1:2:end); % color 1 downsample
C2 = I(2:2:end,2:2:end); % color 2 downsample
C3 = I(2:2:end,1:2:end); % color 3 downsample
C4 = I(1:2:end,2:2:end); % color 4 downsample
else
C1 = I(2:2:end,2:2:end); % color 1 downsample
C2 = I(1:2:end,1:2:end); % color 2 downsample
C3 = I(1:2:end,2:2:end); % color 3 downsample
C4 = I(2:2:end,1:2:end); % color 4 downsample
end
I1 = im2col(C1,[n,n]);
I2 = im2col(C2,[n,n]);
I3 = im2col(C3,[n,n]);
I4 = im2col(C4,[n,n]);
I12 = [I1 I2];
I34 = [I3 I4];
clear C1 C2 C3 C4
% distance weights
H = exp(-(-numM:numM).^2/k);
H = H'*H;
H1 = H(I1((end+1)/2,:));
H2 = H(I2((end+1)/2,:));
H3 = H(I3((end+1)/2,:));
H4 = H(I4((end+1)/2,:));
H12 = H(I12((end+1)/2,:));
H34 = H(I34((end+1)/2,:));
clear H
% TLS constraint matrix
len1 = length(shuffle1(:));
len2 = length(shuffle2(:));
len3 = length(shuffle3(:));
len4 = length(shuffle4(:));
M1 = [ eye(len1,len1) zeros(len1,len2-1) zeros(len1,len3) zeros(len1,len4-1)
zeros(len2,len1) [eye(len2-1);-ones(1,len2-1)] zeros(len2,len3) zeros(len2,len4-1)
zeros(len3,len1) zeros(len3,len2-1) eye(len3,len3) zeros(len3,len4-1)
zeros(len4,len1) zeros(len4,len2-1) [zeros(len4-1,len3);-ones(1,len3)] [eye(len4-1);-ones(1,len4-1)]];
M2 = [[eye(len1-1);-ones(1,len1-1)] zeros(len1,len2) zeros(len1,len3) zeros(len1,len4-1)
zeros(len2,len1-1) eye(len2,len2) zeros(len2,len3) zeros(len2,len4-1)
zeros(len3,len1-1) zeros(len3,len2) eye(len3,len3) zeros(len3,len4-1)
zeros(len4,len1-1) zeros(len4,len2) [zeros(len4-1,len3);-ones(1,len3)] [eye(len4-1);-ones(1,len4-1)]];
M3 = [ eye(len1,len1) zeros(len1,len2-1) zeros(len1,len3) zeros(len1,len4-1)
[zeros(len2-1,len1);-ones(1,len1)] [eye(len2-1);-ones(1,len2-1)] zeros(len2,len3) zeros(len2,len4-1)
zeros(len3,len1) zeros(len3,len2-1) eye(len3,len3) zeros(len3,len4-1)
zeros(len4,len1) zeros(len4,len2-1) zeros(len4,len3) [eye(len4-1);-ones(1,len4-1)]];
M4 = [ eye(len1,len1) zeros(len1,len2-1) zeros(len1,len3-1) zeros(len1,len4)
[zeros(len2-1,len1);-ones(1,len1)] [eye(len2-1);-ones(1,len2-1)] zeros(len2,len3-1) zeros(len2,len4)
zeros(len3,len1) zeros(len3,len2-1) [eye(len3-1);-ones(1,len3-1)] zeros(len3,len4)
zeros(len4,len1) zeros(len4,len2-1) zeros(len4,len3-1) eye(len4,len4) ];
M12 = [ eye(len1,len1) zeros(len1,len2) zeros(len1,len3-1) zeros(len1,len4-1)
zeros(len2,len1) eye(len2,len2) zeros(len2,len3-1) zeros(len2,len4-1)
zeros(len3,len1) zeros(len3,len2) [eye(len3-1);-ones(1,len3-1)] zeros(len3,len4-1)
zeros(len4,len1) zeros(len4,len2) zeros(len4,len3-1) [eye(len4-1);-ones(1,len4-1)]];
M34 = [[eye(len1-1);-ones(1,len1-1)] zeros(len1,len2-1) zeros(len1,len3) zeros(len1,len4)
zeros(len2,len1-1) [eye(len2-1);-ones(1,len2-1)] zeros(len2,len3) zeros(len2,len4)
zeros(len3,len1-1) zeros(len3,len2-1) eye(len3,len3) zeros(len3,len4)
zeros(len4,len1-1) zeros(len4,len2-1) zeros(len4,len3) eye(len4,len4) ];
M1(shuffle,:)=M1;
M2(shuffle,:)=M2;
M3(shuffle,:)=M3;
M4(shuffle,:)=M4;
M12(shuffle,:)=M12;
M34(shuffle,:)=M34;
function [X]=TLSlimit(X,m,n)
H = zeros(m,m);
K = zeros(m,m);
if mod((m-1)/2,2)==0
H(1:2:end,1:2:end) = 1;
else
H(2:2:end,2:2:end) = 1;
end
K(1:2:end,1:2:end) = 1;
K(2:2:end,2:2:end) = 1;
Xmin(:,:,1) = ordfilt2(X(:,:,1),n,H);
Xmin(:,:,2) = ordfilt2(X(:,:,2),n*2,K);
Xmin(:,:,3) = ordfilt2(X(:,:,3),n,H);
Xmax(:,:,1) = ordfilt2(X(:,:,1),sum(H(:))-n+1,H);
Xmax(:,:,2) = ordfilt2(X(:,:,2),sum(K(:))-n*2+1,K);
Xmax(:,:,3) = ordfilt2(X(:,:,3),sum(H(:))-n+1,H);
X = max(min(X,Xmax),Xmin);
function [Y] = TLSpad(X,n)
if n>0
if mod(n,2)==1 % if odd
n = n+1;
end
Y = [X(n+1:-1:2,n+1:-1:2,:) X(n+1:-1:2,:,:) X(n+1:-1:2,end-1:-1:end-n,:)
X(:,n+1:-1:2,:) X X(:,end-1:-1:end-n,:)
X(end-1:-1:end-n,n+1:-1:2,:) X(end-1:-1:end-n,:,:) X(end-1:-1:end-n,end-1:-1:end-n,:)];
end
if n<0
n=-n;
if mod(n,2)==1
n=n+1;
end
Y = X(n+1:end-n,n+1:end-n,:);
end
