polarimetricImage( x, y, S11, S12, S21, S22, coding )

%% function polim = polarimetricImage( x, y, S11, S12, S21, S22, coding )
%
% polarimetricImage.m produces an RGB image from polarimetric S-parameters.  
% The x and y parameters the 2-dim coordinates and may be empty.
% The four channels are labelled S11, S12, S21 and S22.  They may be empty.
% The input value for 'coding' must be a string and can be 'lexicographic' 
% or 'pauli'.  The default value is 'lexicographic'.
% The function returns rgb-matrix, which can be viewed with image.m
%
% author:  [e-mail]

function polim = polarimetricImage( x, y, S11, S12, S21, S22, coding )

if(nargin<7)
	coding = 'lexicographic';
end

polmatSize = [ size(S11); size(S12); size(S21); size(S22) ];
nr = max(polmatSize(:,1));
nc = max(polmatSize(:,2));
nel = prod(polmatSize,2);

if( find(and((nel<nr*nc),(nel>0))) )
	error( '\npolarimetricCoding() requires channels of equal size or empty.\n');
end

% Handle one or more co-pol channels missing
if(nel(1)==0)
	S11 = zeros(nr,nc);
end
if(nel(4)==0) 
	S12 = zeros(nr,nc);
end
if(nel(1)==0 && nel(4)==nr*nc) 
	S11 = S22;
	display( 'Missing S11 in polarimetricCoding():  assert S11 = S22' );
	display( 'To over-ride the assertion, use a matrix of zeros for S11.');
end
if(nel(1)==nr*nc && nel(4)==0)
	S22=S11;
	display( 'Missing S22 in polarimetricCoding():  assert S22 = S11' );
	display( 'To over-ride the assertion, use a matrix of zeros for S22.');
end

% Handle one or more cross-pol channels missing
if(nel(2)==0)
	S12 = zeros(nr,nc);
end
if(nel(3)==0) 
	S21 = zeros(nr,nc);
end
if(nel(2)==0 && nel(3)==nr*nc) 
	S12=S21;
end
if(nel(2)==nr*nc && nel(3)==0)
	S21=S12;
end

switch coding
	case 'lexicographic'
		polim(:,:,3) = abs(S11);
		polim(:,:,2) = abs((S12+S21));
		polim(:,:,1) = abs(S22);
	case 'pauli'
		polim(:,:,3) = abs(S11+S22);
		polim(:,:,2) = abs(S12+S21);
		polim(:,:,1) = abs(S11-S22);
end

polim = polim/max(polim(:));