close all;
clear all;
cd /local/class/psych221/tutorials
load macbethChart; 
greenChip = 11; redChip = 15; whiteChip = 19; grayChip = 20;
macbethImage = flipud(rot90(reshape(1:24,6,4)));

load D65; 
lgt = D65;

load cMatch/monitor; 			% This contains the monitor phosphors
load cMatch/monitorGam;  % This contains monitorInvGamma
load XYZ;
cd ~/proj;

% Scale factor used to color cast image.
factor=1.5;

% The gamma we found for the monitor
myGamma=2.7;

macbethXYZ = XYZ'* diag(lgt)* macbethChart;
xyz2rgb = inv(XYZ'*phosphors);
macbethLinearRGB = xyz2rgb*macbethXYZ;

% Darken surface 20
%macbethLinearRGB(:,20)=0.5*macbethLinearRGB(:,20);

macbethLinearRGB(1,:) = factor*macbethLinearRGB(1,:);
macbethLinearRGB(2,:) = factor*macbethLinearRGB(2,:);
%macbethLinearRGB(3,:) = factor*macbethLinearRGB(3,:);

wht = macbethLinearRGB(:,whiteChip);

wht(1) = wht(1)/factor;
wht(2) = wht(2)/factor;
%wht(3) = wht(3)/factor;

macbethLinearRGB = diag(1./wht)*macbethLinearRGB;

%macbethLinearRGB(:,16)=0;
%macbethLinearRGB(:,19)=0;

macbethLinearRGB = scale(macbethLinearRGB,0,1);
macbethRGB = (macbethLinearRGB).^(1/myGamma);
figure; colormap(macbethRGB'); 
image(macbethImage);
title('Original Macbeth Image');

% Perform Gray World assumption color balancing on filtered image.
newLinearRGB = GWbalance(macbethLinearRGB');
newLinearRGB = scale(newLinearRGB,0,1);
newRGB = (newLinearRGB).^(1/myGamma);
newColorMap = newRGB;

figure; colormap(newColorMap)
image(macbethImage);
title('Forced Gray World Assumption Image');
