Your assignment it to write code in matlab that takes an ima

Your assignment it to write code in matlab that takes an image and removes shadow from the image. Write the whole code yourself, without using any built-in function. I would like the code to be as much long as you can please

1.call ratio

function ratio = cal_ratio(seg,pair,hard_map,soft_map,im)
hard_map_2=conv2(hard_map,fspecial(\'average\',20),\'same\');
imcanny=edge(hard_map_2,\'canny\');
%imcanny=edge(hard_map,\'canny\');
[x,y,t]=bdry_extract_3(imcanny);

ps = 4;
[height width] = size(double(imcanny));
npt = length(x);

pair_ratio=[];

nonpair_ratio=[];

for n=1:npt
skip=false;
x0 = x(n); y0 = y(n); t0 = t(n);
x1 = round(x0 + 1.5*ps*cos(t0+pi/2));
y1 = round(y0 + 1.5*ps*sin(t0+pi/2));
x2 = round(x0 + 1.5*ps*cos(t0-pi/2));
y2 = round(y0 + 1.5*ps*sin(t0-pi/2));
%hist{1} = ones(3, 1); hist{2} = ones(3, 1);
val{1} = zeros(3, 1); val{2} = zeros(3, 1);
k=zeros(2,1);
%cnt1 = 0; cnt2 = 0;
ix = max(1, x1 - ps):min(width, x1 + ps);
iy = max(1, y1 - ps):min(height, y1 + ps);
for ch=1:3
temp = im(iy, ix, ch);
val{1}(ch) = mean(temp(:));

%% find the label for this tmp area

if isnan(val{1}(ch))
skip=true;
end
end

label=seg(iy, ix);
label1=mode(double(label(:)));

temp=soft_map(iy,ix);
k(1)=mean(temp(:));
if isnan(k(1))
skip=true;
end

temp=hard_map(iy,ix);
if mean(temp(:))<0.8 && mean(temp(:))>0.2
skip=true;
end

%region1=sum(sum(hard_map(iy,ix)));

ix = max(1, x2 - ps):min(width, x2 + ps);
iy = max(1, y2 - ps):min(height, y2 + ps);
for ch=1:3
temp = im(iy, ix, ch);
val{2}(ch) = mean(temp(:));
if isnan(val{2}(ch))
skip=true;
end
end

label=seg(iy, ix);
label2=mode(double(label(:)));

temp=soft_map(iy,ix);
k(2)=mean(temp(:));
if isnan(k(2))
skip=true;
end

temp=hard_map(iy,ix);
if mean(temp(:))<0.8 && mean(temp(:))>0.2
skip=true;
end

if abs(k(1)-k(2))<0.5
skip = true;
end

tmp_ratio = [];

if ~skip
for ch=1:3
tmp_ratio=[tmp_ratio,(val{2}(ch)-val{1}(ch))/(val{1}(ch)*k(2)-val{2}(ch)*k(1))];
if isnan(tmp_ratio(ch))
skip = true;
end
end
end

if ~skip
% test to see if label1 and label2 make a pair
if ~isempty(pair) && sum((pair(:,1) == label1 & pair(:,2) == label2) | (pair(:,2) == label2 & pair(:,2) == label1))~=0
pair_ratio = [pair_ratio; tmp_ratio];
else
nonpair_ratio = [nonpair_ratio; tmp_ratio];
end
end

end

if ~isempty(pair_ratio)
pair_ratio(pair_ratio(:,1)<=0|pair_ratio(:,2)<=0|pair_ratio(:,3)<=0,:)=[];
end
if ~isempty(nonpair_ratio)
nonpair_ratio(nonpair_ratio(:,1)<=0|nonpair_ratio(:,2)<=0|nonpair_ratio(:,3)<=0,:)=[];
end
if size(pair_ratio,1) == 0 && size(nonpair_ratio,1)>0
ratio=nonpair_ratio;
else
ratio = pair_ratio;
end

step_size = 0.1;

pt_num = size(ratio,1);

hash_mode = 99991;
hash_map = zeros(4, hash_mode);


%% voting

for i = 1:pt_num
if mod(i, 100)==0
fprintf(1, \'%d\ \', i);
end
ir = floor(ratio(i,1)/step_size); ir = [ir ir+1];
ig = floor(ratio(i,2)/step_size); ig = [ig ig+1];
ib = floor(ratio(i,3)/step_size); ib = [ib ib+1];
if isinf(ir) | isinf(ig) | isinf(ib)
continue
end
for i1=ir
for i2=ig
for i3=ib
hashins([i1, i2, i3]);
end
end
end

end

[dummy,ind] = max(hash_map(1,:));

ratio = hash_map(2:4,ind);
ratio = ratio*step_size;

function ind = hashins(entry)
h = round(entry(1)*1e2 + entry(2)*10 + entry(3));
hm = h;
cnt = 0;
while 1
hm = mod(hm, hash_mode) + 1;
if hash_map(1, hm) == 0
hash_map(1, hm) = 1;
hash_map(2:4, hm) = entry;
break;
elseif sum ( hash_map(2:4, hm) == entry\' ) == 3
hash_map(1, hm) = hash_map(1, hm) + 1;

break;
else
hm = hm+2*cnt+1;
cnt = cnt + 1;
end
end
ind = hm;
end

end

2.deshadow

addpath(\'meanshift\');
addpath(\'used_desc\');
addpath(\'svmlight_mex\');
addpath(\'utils\');
addpath(genpath(\'UGM\'));

unix(\'rm data/binary/* data/mask/* data/cache/* data/matting/* data/removal/* data/unary/*\');
opt = {};
opt.dir = \'data/\';
fnlist1 = {\'DSC01536.jpg\', \'DSC01611.jpg\', \'DSC01508.jpg\', \'DSC01603.jpg\'};
fnlist2 = {\'p21_1.jpg\', \'p14_1.jpg\', \'siebel.bmp\'};
for i=1:length(fnlist1)
opt.fn = fnlist1{i};
opt.save = 1;
opt.binaryClassifier = \'models/model_pair.mat\';
opt.unaryClassifier = \'models/unary_model.mat\';
opt.resize = 1;
opt.adjecent = 0;
opt.pairwiseonly = 0;
opt.linearize = 0;
h = findshadowcut_cross(opt);
end

for i=1:length(fnlist2)
opt.fn = fnlist2{i};
opt.save = 1;
opt.binaryClassifier = \'models/model_pair_our.mat\';
opt.unaryClassifier = \'models/unary_model_our.mat\';
opt.resize = 1;
opt.adjecent = 0;
opt.pairwiseonly = 0;
opt.linearize = 0;
h = findshadowcut_cross(opt);
end

testfnlist = {}; cnt = 0;
for i=1:length(fnlist1)
cnt = cnt+1;
testfnlist{cnt} = fnlist1{i}(1:end-4);
end
for i=1:length(fnlist2)
cnt = cnt+1;
testfnlist{cnt} = fnlist2{i}(1:end-4);
end

runmatting;

3.find shadow cut

function [hardmap] = findshadow(opt)
if ~isfield(opt, \'save\')
opt.save = 0;
end
if ~isfield(opt, \'adjecent\')
opt.adjecent = 0;
end
if ~isfield(opt, \'pairwiseonly\')
opt.pairwiseonly = 0;
end
basename = opt.fn(1:end-4);
im=imread([opt.dir \'original/\' opt.fn]);

%im = (double(im)/255).^(.45).*255;

oim = im;
org_siz = [size(im, 1) size(im, 2)];

if opt.resize == 0
resiz_ratio = 1;
else
resiz_ratio = 640/org_siz(2);
end
im = imresize(oim, resiz_ratio);
if opt.linearize
gamma_im = (double(im)/255).^(1/.45);
im = gamma_im;
end
%im = (double(im)/255);
%im = (double(im)/255).^(1/.45);

try
load([opt.dir \'cache/\' basename \'_seg.mat\']);
catch exp1
disp \'Segmenting\'
[dummy seg] = edison_wrapper(im, @RGB2Luv, ...
\'SpatialBandWidth\', 9, \'RangeBandWidth\', 15, ...
\'MinimumRegionArea\', 200);
seg = seg + 1;
save([opt.dir \'cache/\' basename \'_seg.mat\'], \'seg\');
end

numlabel = length(unique(seg(:)));

try
load([opt.dir \'cache/\' basename \'_single.mat\']);
catch exp1
disp \'Single region feature\'

%load(\'cache/model_our.mat\', \'model\');
labhist = calcLabHist(im, seg, numlabel);
texthist = calcTextonHistNoInv(im, seg, numlabel);
testdata=[labhist, texthist];
testlabel = ones(numlabel, 1);
save([opt.dir \'cache/\' basename \'_single.mat\'], \'testdata\', \'testlabel\');
end
load(opt.unaryClassifier, \'model\');
disp \'Single region classification\'
[err, s]=svmclassify(testdata, testlabel, model);
ss = double(sign(s));
ssmap=(1-ss(seg))/2;

try
load([opt.dir \'cache/\' basename \'_pair.mat\']);
catch exp1

disp \'Pair region feature\'

shapemean = calcShapeMean(seg, numlabel);
area = shapemean.area;
centroids = shapemean.center;
rgbmean = calcRGBMean(im, seg, numlabel);
texthist = calcTextonHist(im, seg, numlabel);
labhist = calcLabHist(im, seg, numlabel);
d_text=dist_chi2(texthist\', texthist\');
d_lab =dist_chi2(labhist\', labhist\');

ind = [];
finalvector = zeros(numlabel,numlabel, 8);

for i=1:numlabel
for j=1:numlabel
if i==j , continue; end
r = rgbmean(i,:)./rgbmean(j,:);
dist = norm(centroids(i,:) - centroids(j,:))/sqrt(sqrt(area(i)*area(j)));
finalfeature = [d_text(i,j), d_lab(i,j), ...
r, dist, abs(r(1)./r(2)-1)*10, abs(r(2)./r(3)-1)*10];
finalvector(j,i,:)=finalfeature;
end
end
finalvector = reshape(finalvector, [numlabel*numlabel 8]);
g_diff = zeros(numlabel, numlabel);
g_same = zeros(numlabel, numlabel);
save([opt.dir \'cache/\' basename \'_pair.mat\'], \'finalvector\', \'shapemean\');
end
load(opt.binaryClassifier, \'diffmodel\', \'samemodel\');
disp \'Pair region classification\'
[err, s1]=svmclassify(finalvector, zeros(numlabel*numlabel, 1) , diffmodel);
[err, s2]=svmclassify(finalvector, zeros(numlabel*numlabel, 1), samemodel);
s1 = reshape(s1, [numlabel numlabel]);
s2 = reshape(s2, [numlabel numlabel]);
k1=100; k2=200;
t1 = sort(s1(:));t1(isnan(t1))=[];
t2 = sort(s2(:));t2(isnan(t2))=[];
thresh1 = t1(max(1, length(t1(:))-k1)); thresh1 = max(.6, thresh1);
thresh2 = t2(max(1, length(t2(:))-k2)); thresh2 = max(.6, thresh2);
% thresh1 = t1(max(1, length(t1(:))-k1)); thresh1 = max(0, thresh1);
% thresh2 = t2(max(1, length(t2(:))-k2)); thresh2 = max(0, thresh2);
% thresh1 = 0;
% thresh2 = 0;
%FIXME!!!!!!!!!!!!
if strcmp(\'models/model_pair_our.mat\', opt.binaryClassifier)
for i=1:numlabel
for j=1:numlabel
if i==j , continue; end
w = sqrt(shapemean.area(i)*shapemean.area(j));
g_diff(i,j) = w*s1(j,i);
g_same(i,j) = w*s2(j,i);
end
end
else
for i=1:numlabel
for j=1:numlabel
if i==j , continue; end
w = sqrt(shapemean.area(i)*shapemean.area(j));
g_diff(i,j) = w*s1(i,j);
g_same(i,j) = w*s2(i,j);
end
end
end

nNodes = numlabel;
nStates = 2;
adj1 = logical(sparse(nNodes,nNodes));
adj2 = logical(sparse(nNodes,nNodes));

m=buildAdjacentMatrix(seg, numlabel);

for i=1:numlabel
for j=1:numlabel
if opt.adjecent,
if ~m(i,j) continue; end;
end;
if s1(i,j)>thresh1
adj1(i,j)=1;
end
if s2(i,j)>thresh2
adj2(i,j)=1;
end
end
end

nodePot = zeros(nNodes,nStates);
w1 = 1;
if ~opt.pairwiseonly
for i=1:numlabel
wi = w1 * shapemean.area(i);
nodePot(i,1) = -s(i)*wi;
nodePot(i,2) = s(i)*wi;
end
end

if 1
sc = shapemean.center;
nim = im;
[gx gy] = gradient(double(seg));
eim = (gx.^2+gy.^2)>1e-10;

t = nim(:,:,1); t(eim)=0; nim(:,:,1)=t;
t = nim(:,:,2); t(eim)=0; nim(:,:,2)=t;
t = nim(:,:,3); t(eim)=0; nim(:,:,3)=t;
f3 = figure(3);
clf;
imshow(nim);
hold on
for i=1:numlabel
for j=1:numlabel
if ~adj1(i,j) && ~adj2(i,j), continue; end
if s1(i,j)>thresh1
plot([sc(i,1) sc(j,1)], [sc(i,2) sc(j,2)], \'b\');
plot(sc(i,1), sc(i,2), \'bo\')
elseif s2(i,j)>thresh2
plot([sc(i,1) sc(j,1)], [sc(i,2) sc(j,2)], \'r\');
end
end
end
%print(f3, \'-dpsc\', [opt.dir \'cache/\' basename \'_graph.eps\']);
if ~opt.save
figure(3)
imshow(nim);
hold on
for i=1:numlabel
for j=1:numlabel
if ~adj1(i,j) && ~adj2(i,j), continue; end
if s1(i,j)>thresh1
plot([sc(i,1) sc(j,1)], [sc(i,2) sc(j,2)], \'b\');
plot(sc(i,1), sc(i,2), \'bo\')
elseif s2(i,j)>thresh2
plot([sc(i,1) sc(j,1)], [sc(i,2) sc(j,2)], \'r\');
end
end
end
figure(6), imagesc(s(seg))
figure(4), imagesc(seg)
end
end

edgeStruct1 = UGM_makeEdgeStruct_directed(adj1,nStates);
edgeStruct2 = UGM_makeEdgeStruct_directed(adj2,nStates);
edgePot = [];
w3=1;
w2=2;
for e = 1:edgeStruct1.nEdges
n1 = edgeStruct1.edgeEnds(e,1);
n2 = edgeStruct1.edgeEnds(e,2);
%nodePot(n1,1) = nodePot(n1,1)+ w2*(g_diff(n1, n2)-g_diff(n2, n1));
%nodePot(n2,2) = nodePot(n2,2)+ w2*(g_diff(n1, n2)-g_diff(n2, n1));
nodePot(n1,1) = nodePot(n1,1)+ w2*g_diff(n1, n2);
nodePot(n1,2) = nodePot(n1,2)- w2*g_diff(n1, n2);
nodePot(n2,1) = nodePot(n2,1)- w2*g_diff(n1, n2);
nodePot(n2,2) = nodePot(n2,2)+ w2*g_diff(n1, n2);
end

for e = 1:edgeStruct2.nEdges
n1 = edgeStruct2.edgeEnds(e,1);
n2 = edgeStruct2.edgeEnds(e,2);

edgePot(:,:,e) = [g_same(n1, n2) 0;...
0, g_same(n1, n2) ].*[w3 1; 1 w3];
end

if ~isempty(edgePot)
Decoding = UGM_Decode_ModifiedCut(nodePot,edgePot,edgeStruct2);
else
Decoding = double(sign(s))+1;
end

hardmap = Decoding(seg)-1;

if ~opt.save
figure(4)
ss = double(sign(s));
ssmap=(1-ss(seg))/2;
imshow(double(ssmap))

figure(5)
imshow(1-double(hardmap))

figure(1);
imshow(im);
else
imwrite(hardmap,[opt.dir \'binary/\' basename \'_binary.png\']);
ss = double(sign(s));
ssmap=(1-ss(seg))/2;
imwrite(1-ssmap,[opt.dir \'unary/\' basename \'_unary.png\']);
save([opt.dir \'cache/\' basename \'_detect.mat\'],\'hardmap\',\'ssmap\');

shadowim = hardmap;

tmp = ones(size(im));
tmp2 = cat(3, zeros([size(shadowim) 2]), 0.5*shadowim);

mask = logical(repmat(shadowim, [1 1 3]));
tmp(mask) = tmp2(mask);
%tmp = shadowim.*ones(size(shadowim));
%tmp = cat(3, zeros([size(shadowim) 2]), shadowim);

grayim = im2double(repmat(rgb2gray(im), [1 1 3]));
im2 = (grayim+tmp)/2;
imwrite(im2, [opt.dir \'mask/\' basename \'_mask.png\']);
end

RegionScore = s;
DiffScore = s1;
SameScore = s2;
DiffAdj = adj1;
SameAdj = adj2;
save([opt.dir \'cache/\' basename \'_everything.mat\'], \'RegionScore\', \'DiffScore\', \'SameScore\', \'DiffAdj\', \'SameAdj\', \'seg\', \'im\', \'hardmap\', \'ssmap\');

% get the pairing information consistent with the final detection.
[shadow, non_shadow] = find(adj1 == 1);

pair = [];

num_pair = numel(non_shadow);
for i = 1:num_pair
if Decoding(non_shadow(i))==1 && Decoding(shadow(i)) == 2
pair = [pair; non_shadow(i) shadow(i)];
end
end

save([opt.dir \'cache/\' basename \'_finalpair.mat\'], \'pair\');
end

4run matting

addpath(\'matting\');
addpath(\'utils\');

if (~exist(\'thr_alpha\',\'var\'))
thr_alpha=[];
end
if (~exist(\'epsilon\',\'var\'))
epsilon=[];
end
if (~exist(\'win_size\',\'var\'))
win_size=[];
end

if (~exist(\'levels_num\',\'var\'))
levels_num=1;
end
if (~exist(\'active_levels_num\',\'var\'))
active_levels_num=1;
end

for i =1:numel(testfnlist)
fprintf(\'processing %d out of %d\ \', i, numel(testfnlist));

load([\'data/cache/\' testfnlist{i} \'_finalpair.mat\']);

load([\'data/cache/\' testfnlist{i} \'_everything.mat\']);
img = im;
I=double(img)/255;

load([\'data/cache/\' testfnlist{i} \'_detect.mat\']);

hard_map=hardmap;

hard_map_1=conv2(hard_map,fspecial(\'average\',5),\'same\');
hard_map_2=conv2(hard_map,fspecial(\'average\',5),\'same\');

consts_map=hard_map_2>1-1e-8 | hard_map_1==0;

consts_vals=hard_map_1==0;

alpha=solveAlpha(I,consts_map,consts_vals,epsilon,win_size);

soft_map=alpha;

ratio = cal_ratio(seg,pair,hard_map,soft_map,I);

red=I(:,:,1);
green=I(:,:,2);
blue=I(:,:,3);

red_val=ratio(1);
green_val=ratio(2);
blue_val=ratio(3);

red_recover=(ones(size(soft_map))*(red_val+1))./(soft_map*red_val+1).*red;
green_recover=(ones(size(soft_map))*(green_val+1))./(soft_map*green_val+1).*green;
blue_recover=(ones(size(soft_map))*(blue_val+1))./(soft_map*blue_val+1).*blue;

shadow_free=zeros(size(I));

shadow_free(:,:,1)=red_recover;
shadow_free(:,:,2)=green_recover;
shadow_free(:,:,3)=blue_recover;

shadow_free(shadow_free(:)>1) = 1;

imwrite(soft_map, [\'data/matting/\' testfnlist{i} \'_soft.jpg\'], \'jpg\');
imwrite(shadow_free, [\'data/removal/\' testfnlist{i} \'_recovery.jpg\'], \'jpg\');

end

5UGM decode