
load /local/class/psych221/tutorials/cMatch/monitor;
load XYZ;
rgb2xyz = XYZ'*phosphors;
xyz2rgb = inv(rgb2xyz);
clear XYZ phosphors comment wavelength

%xyz = 0.5*ones(3,256);
%xyz(1,:) = (0:1/255:1);
%rgb = xyz2rgb * xyz;
%r = ones(64,1)*rgb(1,:);
%g = ones(64,1)*rgb(2,:);
%b = ones(64,1)*rgb(3,:);

if (exist('bg1')~=1) bg1 = [0.1 0.1 0.1]; end
if (exist('bg2')~=1) bg2 = [0.8 0.8 0.8]; end
if (exist('fg1')~=1) fg1 = [0.5 0.6 0.3]; end

xyz1 = rgb2xyz*fg1';
xyzw = rgb2xyz*[1 1 1]';
tmp = rgb2xyz*bg1'; Yb1 = tmp(2)+1; 
tmp = rgb2xyz*bg2'; Yb2 = tmp(2)+1; 

LA = 20;  	% white is approx 100 cd/m^2

% avg surround, samples subtending 4o
% c = 0.69;	Nc = 1.0;	FLL = 0.0; 	F = 1.0;
% average surround
c = 0.69;	Nc = 1.0;	FLL = 1.0;	F = 1.0;
% dim surround
% c = 0.59;	Nc = 1.1;	FLL = 1.0;	F = 0.9;
% dark surround
% c = 0.525;	Nc = 0.8;	FLL = 1.0;	F = 0.9;
% cut sheet
% c = 0.41;	Nc = 0.8;	FLL = 1.0;	F = 0.9;

D = 1.0;	% discounting the illuminant completely (complete adaptation)

jch = ciecam_jc(xyz1,xyzw,LA,c,Nc,FLL,F,Yb1,D);
xyz2 = inv_ciecam_jc(jch,xyzw,LA,c,Nc,FLL,F,Yb2,D);

fg2 = (xyz2rgb * xyz2)';

indximg = [ 1 1 1 1 1 3 3 3 3 3; ...
            1 1 1 1 1 3 3 3 3 3; ...
            1 1 2 1 1 3 3 4 3 3; ...
            1 1 1 1 1 3 3 3 3 3; ...
            1 1 1 1 1 3 3 3 3 3];

cmap = [bg1;fg1;bg2;fg2];
displayImage(indximg);
imshow(indximg,cmap);
axis('image');
axis('off');
title(sprintf('[%6.4f %6.4f %6.4f] foregrounds [%6.4f %6.4f %6.4f]', ...
              fg1,fg2));

