# Diffing "Revenge of Antichaff" and "Same people different day"

 Title: Revenge of Antichaff Same people different day Author: Nicholas Howe Mike Russell Submitted: 2011-11-04 15:34:08 UTC 2011-11-04 17:02:58 UTC Status: Passed Passed Score: 671223.0 1035290.0 Result: 67086812 (cyc: 33, node: 2645) 65841499 (cyc: 33, node: 8284) CPU Time: 81.089 142.901 Code: ```function [moves, vine] = nochaff4(board, limit) function [moves, vine] = solver(board, limit) [moves,vine] = nick(board,limit); [moves1,vine1] = alfi(board,limit); score1 = grade(moves,vine,board,limit); score2 = grade(moves1,vine1,board,limit); if score1 > score2 moves = moves1; vine = vine1; elseif score1 == score2 && length(moves1) < length(moves) moves = moves1; vine = vine1; end end function [moves, vine] = nick(board, limit) % Added filter to crosseyed % % N. Howe [vine1,scr1] = trickle2(board); moves1 = []; [moves2,vine2,scr2] = pusher(board,limit); [nrow,ncol] = size(board); ntot = nrow*ncol; [moves3,vine3,scr3] = pusher(board',limit); [moves4,vine4,scr4] = pusher(fliplr(board),limit); [moves5,vine5,scr5] = pusher(fliplr(board'),limit); [moves6,vine6,scr6] = dofilter(board,limit); if (scr1==max([scr1 scr2 scr3 scr4 scr5 scr6])) moves = moves1; vine = vine1; elseif (scr2==max([scr1 scr2 scr3 scr4 scr5 scr6])) moves = moves2; vine = vine2; elseif (scr3==max([scr1 scr2 scr3 scr4 scr5 scr6])) id = reshape(1:ntot,nrow,ncol); id3 = id'; moves = id3(moves3); vine = id3(vine3); elseif (scr4==max([scr1 scr2 scr3 scr4 scr5 scr6])) id = reshape(1:ntot,nrow,ncol); id4 = fliplr(id); moves = id4(moves4); vine = id4(vine4); elseif (scr5==max([scr1 scr2 scr3 scr4 scr5 scr6])) id = reshape(1:ntot,nrow,ncol); id5 = fliplr(id'); moves = id5(moves5); vine = id5(vine5); else moves = moves6; vine = vine6; end; end function [moves, vine, scr] = pusher(board, limit) moves = zeros(limit,2); [nrow,ncol] = size(board); if (limit < nrow) scr = -inf; moves = []; vine = []; return; end; ntot = nrow*ncol; xg = repmat(1:ncol,nrow,1); yg = repmat((1:nrow)',1,ncol); pos = reshape(1:ntot,nrow,ncol); fitloc = reshape(1:ntot,nrow,ncol); fitloc(:,2:2:end) = flipud(fitloc(:,2:2:end)); [s,r] = sort(board(:),'descend'); rank = zeros(nrow,ncol); rank(r) = 1:ntot; remain = true(1,ntot); nextfit = 1; steps = 0; cboard = board; for i = 1:ntot if remain(i) [cy,cx] = find(rank==i); tx = xg(fitloc(nextfit)); ty = yg(fitloc(nextfit)); while (cx~=tx)|(cy~=ty) dx = sign(tx-cx); dy = sign(ty-cy); if (dx~=0)&(dy~=0) if cboard(cy,cx+dx)<=cboard(cy+dy,cx) dy = 0; else dx = 0; end; end; if (rank(cy+dy,cx+dx)>0) vx = cx+dx; vy = cy+dy; if (dx==0) if (vx>1)&(cboard(vy,vx)>cboard(vy,vx-1))&((vx==ncol)|(cboard(vy,vx+1)>cboard(vy,vx-1))) % move left steps = steps+1; if (steps > limit) break; end; moves(steps,1) = pos(vy,vx); moves(steps,2) = pos(vy,vx-1); cboard(vy,vx-1) = cboard(vy,vx); cboard(vy,vx) = 0; if (rank(vy,vx-1)>0) remain(rank(vy,vx-1)) = 0; end; rank(vy,vx-1) = rank(vy,vx); rank(vy,vx) = 0; %imagesc(cboard); %axis image; %drawnow; elseif (vx < ncol)&(cboard(vy,vx)>cboard(vy,vx+1)) % move right steps = steps+1; if (steps > limit) break; end; moves(steps,1) = pos(vy,vx); moves(steps,2) = pos(vy,vx+1); cboard(vy,vx+1) = cboard(vy,vx); cboard(vy,vx) = 0; if (rank(vy,vx+1)>0) remain(rank(vy,vx+1)) = 0; end; rank(vy,vx+1) = rank(vy,vx); rank(vy,vx) = 0; %imagesc(cboard); %axis image; %drawnow; else remain(rank(cy+dy,cx+dx)) = 0; end; else if (vy>1)&(cboard(vy,vx)>cboard(vy-1,vx))&((vy==nrow)|(cboard(vy+1,vx)>cboard(vy-1,vx))) % move north steps = steps+1; if (steps > limit) break; end; moves(steps,1) = pos(vy,vx); moves(steps,2) = pos(vy-1,vx); cboard(vy-1,vx) = cboard(vy,vx); cboard(vy,vx) = 0; if (rank(vy-1,vx)>0) remain(rank(vy-1,vx)) = 0; end; rank(vy-1,vx) = rank(vy,vx); rank(vy,vx) = 0; %imagesc(cboard); %axis image; %drawnow; elseif (vy < nrow)&(cboard(vy,vx)>cboard(vy+1,vx)) % move right steps = steps+1; if (steps > limit) break; end; moves(steps,1) = pos(vy,vx); moves(steps,2) = pos(vy+1,vx); cboard(vy+1,vx) = cboard(vy,vx); cboard(vy,vx) = 0; if (rank(vy+1,vx)>0) remain(rank(vy+1,vx)) = 0; end; rank(vy+1,vx) = rank(vy,vx); rank(vy,vx) = 0; %imagesc(cboard); %axis image; %drawnow; else remain(rank(cy+dy,cx+dx)) = 0; end; end; end; steps = steps+1; if (steps > limit) break; end; moves(steps,1) = pos(cy,cx); moves(steps,2) = pos(cy+dy,cx+dx); cboard(cy+dy,cx+dx) = s(i); %cboard(cy,cx); cboard(cy,cx) = 0; %rank(cy+dy,cx+dx) = i; %cboard(cy,cx); rank(cy,cx) = 0; cx = cx+dx; cy = cy+dy; %imagesc(cboard); %axis image; %drawnow; end; if (steps > limit) break; end; rank(ty,tx) = i; nextfit = nextfit+1; end; end; moves = moves(1:min(steps,limit),:); [vine,scr] = trickle2(cboard); end function [vine,scr] = trickle2(board) % Second effort % Somewhat better handling of ties. % %N. Howe moves = []; [nrow,ncol] = size(board); xg = repmat(1:ncol,nrow,1); yg = repmat((1:nrow)',1,ncol); limbo = nrow*ncol+1; vals = [board(:)' -inf]; %up = [zeros(1,ncol); board(1:end-1,:)]; %down = [board(2:end,:);zeros(1,ncol)]; %left = [zeros(nrow,1); board(:,1:end-1)]; %right = [board(:,2:end);zeros(nrow,1)]; ind = reshape(1:nrow*ncol,nrow,ncol); up = [repmat(limbo,1,ncol); ind(1:end-1,:)]; down = [ind(2:end,:);repmat(limbo,1,ncol)]; left = [repmat(limbo,nrow,1) ind(:,1:end-1)]; right = [ind(:,2:end) repmat(limbo,nrow,1)]; nbr = [[up(:) down(:) left(:) right(:)]; limbo limbo limbo limbo]; %vnbr = vals(nbr); %big = repmat(board(:),1,4)<=vnbr; cum = vals; next = zeros(1,nrow*ncol+1); next(limbo) = -1; levels = unique(vals); done(nrow*ncol+1) = true; for lvl = levels(end:-1:2) lvli = find(board==lvl); lvln = numel(lvli); for ii = 1:lvln moves = (vals(nbr(lvli,:))>=lvl).*done(nbr(lvli,:)); [heur,jheur] = max(cum(nbr(lvli,:)).*moves,[],2); heur(done(lvli)) = 0; [mheur,mheuri] = max(heur); if (mheur==0) done(lvli) = true; break; end; i = lvli(mheuri); cum(i) = cum(i)+mheur; next(i) = nbr(i,jheur(mheuri)); done(i) = true; end; end; vine = zeros(1,nrow*ncol); [scr,vine(1)] = max(cum); vlen = 1; while (next(vine(vlen))>0) vine(vlen+1) = next(vine(vlen)); vlen = vlen+1; end; % add on at end i = vine(vlen); done(i) = false; lvl = board(i); mv = (vals(nbr(i,:))==lvl).*done(nbr(i,:)); while any(mv) vine(vlen+1) = nbr(i,find(mv,1)); vlen = vlen+1; i = vine(vlen); done(i) = false; lvl = board(i); mv = (vals(nbr(i,:))==lvl).*done(nbr(i,:)); end; vine = vine(1:vlen); % now do augmentation j = 1; used = false(nrow,ncol); used(vine) = true; while (j < vlen) i1 = vine(j); i2 = vine(j+1); x1 = xg(i1); x2 = xg(i2); y1 = yg(i1); %y2 = yg(i2); if (x1==x2) % vertical move if (x1>1)&&(~used(i1-nrow))&&(~used(i2-nrow))&&(board(i1)<=board(i1-nrow))&&(board(i1-nrow)<=board(i2-nrow))&&(board(i2-nrow)<=board(i2)) % augment to the left vlen = vlen+2; vine = [vine(1:j) i1-nrow i2-nrow vine(j+1:end)]; used(i1-nrow) = true; used(i2-nrow) = true; elseif (x11)&&(~used(i1-1))&&(~used(i2-1))&&(board(i1)<=board(i1-1))&&(board(i1-1)<=board(i2-1))&&(board(i2-1)<=board(i2)) % augment to the north vlen = vlen+2; vine = [vine(1:j) i1-1 i2-1 vine(j+1:end)]; used(i1-1) = true; used(i2-1) = true; elseif (y10.67*ncol)||(sum(sign(diff(m2)))>0.67*nrow)||(sum(sign(diff(m2)))<-0.67*nrow)||(sum(sign(diff(m2)))<-0.67*nrow) d = diff(board,1,2); mb = median(board,2); md = median(d,2); rmb = repmat(mb,1,ncol); rmd = repmat(md,1,ncol); guess = rmb+kron(md,-(ncol-1)/2:(ncol-1)/2); bad = (abs(board-guess)>20*abs(rmd)); for i = nrow:-1:1 offset = 0; for j = ceil(ncol/2):ncol if bad(i,j) offset = offset+1; else moves(nmove+1:nmove+offset,:) = [pos(i,j:-1:j-offset+1)' pos(i,j-1:-1:j-offset)']; nmove = nmove+offset; if (nmove>=limit) break; end; end; end; end; offset = 0; for j = ceil(ncol/2):-1:1 if bad(i,j) offset = offset+1; else moves(nmove+1:nmove+offset,:) = [pos(i,j:1:j+offset-1)' pos(i,j+1:1:j+offset)']; nmove = nmove+offset; if (nmove>=limit) break; end; end; end; end; if (nmove>=limit) break; end; end; % b2 = board; % b2(bad) = guess(bad); % % sd = sign(d); % rs = median(sd,2); % rrs = repmat(rs,1,ncol); % for i = 1:nrow % end; % % spd = [rs sd rs]; % bad = ((spd(:,1:end-1)~=gd)&(gd>0))|((spd(:,2:end)~=gd)&(gd<0)); % b2 = board; % b2(bad) = guess(bad); % for i = 1:nrow % for j = 1:ncol % % end; % end; moves = moves(1:min(nmove,limit),:); for i = 1:size(moves,1) board(moves(i,2)) = board(moves(i,1)); board(moves(i,1)) = 0; %imagesc(board); %axis image; %drawnow; end; [vine,scr] = trickle2(board); else scr = -inf; moves = []; vine = []; end; end function [moves, vine] = alfi(board, limit) [moves, vine, gain,newboard] = solver1(board, limit); if 1,%limit>.1*numel(board) [board_preproc,moves_preproc]=preprocess(board,ceil(limit*.125)); %%% limit_preproc=limit-size(moves_preproc,1); [moves2, vine2, gain2,newboard2] = solver1(board_preproc, limit_preproc); if gain2>gain moves=cat(1,moves_preproc,moves2); vine=vine2; gain=gain2; newboard=newboard2; end end % disp(gain) % imagesc(newboard);hold on;[i1,i2]=ind2sub(size(board),vine);plot(i2,i1,'k.-'); hold off; end function [moves, vine, gain, board] = solver1(board, limit) moves=[]; ab=accumarray(1+board(:),1); cab=flipud(cumsum(flipud(ab)))-ab(end); [m,n]=size(board); board=[zeros(1,n+2);[zeros(m,1),board,zeros(m,1)];zeros(1,n+2)]; [m,n]=size(board); boardab=cab(1+board); if max(ab(2:end))>1, % compute same-valued clusters [BoardCluster,ClusterValue,IdxList,IdxSegments,Nclusters]=connected(board); ClusterSize=diff(IdxSegments); ClusterNeighb=neighb(BoardCluster,board,Nclusters); % search between-clusters ClustersOrder = bellman(ClusterNeighb,ClusterValue'.*sqrt(ClusterSize')); % search within-clusters iC=bellman_postprocess(ClustersOrder,IdxList,IdxSegments,BoardCluster); else % search between-pixels iC = bellman_pixel(board,limit,boardab); end % post-processing moves if limit>0 [board,iC,moves,limit] = postprocess(board,iC,moves,limit,boardab); [iC1,iC2]=ind2sub([m,n],moves); moves=sub2ind([m-2,n-2],iC1-1,iC2-1); end [iC1,iC2]=ind2sub([m,n],iC); vine=sub2ind([m-2,n-2],iC1-1,iC2-1); vine=vine(end:-1:1); board=board(2:end-1,2:end-1); gain=sum(board(vine)); end function [board,moves]=preprocess(board,limit) [m,n]=size(board); moves=[]; minmaxm=[1,m]; minmaxn=[1,n]; tn=n; tm=m; tboard=board; ok=1; while ok if max(tm,tn)<=1, break; end if tm>1 [minboard,minidx]=min(tboard,[],1); if mean(minidx)=(minidx(n1)-1) tboard(2:minidx(n1),n1)=tboard(1:minidx(n1)-1,n1); tboard(1,n1)=0; moves=[moves; [(minmaxn(1)-1+n1-1)*m+(minmaxm(1)-1+(minidx(n1)-1:-1:1)'),(minmaxn(1)-1+n1-1)*m+(minmaxm(1)-1+(minidx(n1):-1:2)')]]; limit=limit-(minidx(n1)-1); else ok=0; end end board(minmaxm(1)-1+(1:tm),minmaxn(1)-1+(1:tn))=tboard; minmaxm(1)=minmaxm(1)+1; tboard=tboard(2:tm,:); tm=tm-1; else for n1=1:tn, if limit>=(tm-minidx(n1)) tboard(minidx(n1):tm-1,n1)=tboard(minidx(n1)+1:tm,n1); tboard(tm,n1)=0; moves=[moves; [(minmaxn(1)-1+n1-1)*m+(minmaxm(1)-1+(minidx(n1)+1:tm)'),(minmaxn(1)-1+n1-1)*m+(minmaxm(1)-1+(minidx(n1):tm-1)')]]; limit=limit-(tm-minidx(n1)); else ok=0; end end board(minmaxm(1)-1+(1:tm),minmaxn(1)-1+(1:tn))=tboard; minmaxm(2)=minmaxm(2)-1; tboard=tboard(1:tm-1,:); tm=tm-1; end end if tn>1 [minboard,minidx]=min(tboard,[],2); if mean(minidx)=(minidx(n1)-1) tboard(n1,2:minidx(n1))=tboard(n1,1:minidx(n1)-1); tboard(n1,1)=0; moves=[moves; [(minmaxn(1)-1+(minidx(n1)-1:-1:1)'-1)*m+(minmaxm(1)-1+n1),(minmaxn(1)-1+(minidx(n1):-1:2)'-1)*m+(minmaxm(1)-1+n1)]]; limit=limit-(minidx(n1)-1); else ok=0; end end board(minmaxm(1)-1+(1:tm),minmaxn(1)-1+(1:tn))=tboard; minmaxn(1)=minmaxn(1)+1; tboard=tboard(:,2:tn); tn=tn-1; else for n1=1:tm, if limit>=(tn-minidx(n1)) tboard(n1,minidx(n1):tn-1)=tboard(n1,minidx(n1)+1:tn); tboard(n1,tn)=0; moves=[moves; [(minmaxn(1)-1+(minidx(n1)+1:tn)'-1)*m+(minmaxm(1)-1+n1),(minmaxn(1)-1+(minidx(n1):tn-1)'-1)*m+(minmaxm(1)-1+n1)]]; limit=limit-(tn-minidx(n1)); else ok=0; end end board(minmaxm(1)-1+(1:tm),minmaxn(1)-1+(1:tn))=tboard; minmaxn(2)=minmaxn(2)-1; tboard=tboard(:,1:tn-1); tn=tn-1; end end end end function [board,tiC,moves,limit] = postprocess(board,tiC,moves,limit,boardab) [m,n]=size(board); nC=numel(tiC); if nC>1&&limit>0, %%% % grow laterally d=abs(diff(tiC))==1; E=[m*~d+d; m*d+~d]; BorderUp={repmat(1+(0:n-1)*m,[m,1]),repmat((1:m)',[1,n])}; BorderDown={repmat(m+(0:n-1)*m,[m,1]),repmat(m*(n-1)+(1:m)',[1,n])}; tP=board>0; tP(tiC)=false; for n1=1:1e3, [board,moves,tiC,tP,E,limit,ok]=postprocess_lateral(board,moves,tiC,tP,E,limit,BorderDown,BorderUp,boardab); if ~ok, break; end end end if limit>0 % grow end-points tP=board>0; tP(tiC)=false; for n1=1:1e3, [board,moves,tiC,tP,limit,ok]=postprocess_endpoint(board,moves,tiC,tP,limit,1); if ~ok, break; end end for n1=1:1e3, [board,moves,tiC,tP,limit,ok]=postprocess_endpoint(board,moves,tiC,tP,limit,2); if ~ok, break; end end end end function [board,moves,tiC,tP,E,limit,ok]=postprocess_lateral(board,moves,tiC,tP,E,limit,BorderDown,BorderUp,boardab) [m,n]=size(board); ok=0; % D=Dijkstra(tP,tiC(1)); % maskout=tP&~isinf(D)&board>board(tiC(1)); % nmaskout=sum(D(maskout)); for ndir=1:4 K=size(tiC,2); switch(ndir) case 1, idx1=find(tP(tiC(1:2:K-1)+E(2,1:2:K-1))&tP(tiC(2:2:K)+E(2,1:2:K-1))); case 2, idx1=find(tP(tiC(1:2:K-1)-E(2,1:2:K-1))&tP(tiC(2:2:K)-E(2,1:2:K-1))); case 3, idx1=find(tP(tiC(2:2:K-1)+E(2,2:2:K-1))&tP(tiC(3:2:K)+E(2,2:2:K-1))); case 4, idx1=find(tP(tiC(2:2:K-1)-E(2,2:2:K-1))&tP(tiC(3:2:K)-E(2,2:2:K-1))); end %for n2=1:numel(idx1), for n2=numel(idx1):-1:1, %%% switch(ndir) case 1, tempA=tiC(2*idx1(n2)-1)+E(2,2*idx1(n2)-1):E(2,2*idx1(n2)-1):BorderDown{1+(E(2,2*idx1(n2)-1)==m)}(tiC(2*idx1(n2)-1)); tempB=tiC(2*idx1(n2))+E(2,2*idx1(n2)-1):E(2,2*idx1(n2)-1):BorderDown{1+(E(2,2*idx1(n2)-1)==m)}(tiC(2*idx1(n2))); case 2, tempA=tiC(2*idx1(n2)-1)-E(2,2*idx1(n2)-1):-E(2,2*idx1(n2)-1):BorderUp{1+(E(2,2*idx1(n2)-1)==m)}(tiC(2*idx1(n2)-1)); tempB=tiC(2*idx1(n2))-E(2,2*idx1(n2)-1):-E(2,2*idx1(n2)-1):BorderUp{1+(E(2,2*idx1(n2)-1)==m)}(tiC(2*idx1(n2))); case 3, tempA=tiC(2*idx1(n2))+E(2,2*idx1(n2)):E(2,2*idx1(n2)):BorderDown{1+(E(2,2*idx1(n2))==m)}(tiC(2*idx1(n2))); tempB=tiC(2*idx1(n2)+1)+E(2,2*idx1(n2)):E(2,2*idx1(n2)):BorderDown{1+(E(2,2*idx1(n2))==m)}(tiC(2*idx1(n2)+1)); case 4, tempA=tiC(2*idx1(n2))-E(2,2*idx1(n2)):-E(2,2*idx1(n2)):BorderUp{1+(E(2,2*idx1(n2))==m)}(tiC(2*idx1(n2))); tempB=tiC(2*idx1(n2)+1)-E(2,2*idx1(n2)):-E(2,2*idx1(n2)):BorderUp{1+(E(2,2*idx1(n2))==m)}(tiC(2*idx1(n2)+1)); end idxZ=find(~tP(tempA),1,'first'); %idxZ=find(~tP(tempA)|maskout(tempA),1,'first'); switch(ndir) case {1,2}, %idxA=find(board(tiC(2*idx1(n2)-1))>=board(tempA)&board(tempA)>=board(tiC(2*idx1(n2))),1,'first'); [nill,idxA]=max((1./(abs(board(tiC(2*idx1(n2)-1))-board(tempA))+1)).*(board(tiC(2*idx1(n2)-1))>=board(tempA)&board(tempA)>=board(tiC(2*idx1(n2))))); if ~nill,idxA=[]; end case {3,4}, %idxA=find(board(tiC(2*idx1(n2)))>=board(tempA)&board(tempA)>=board(tiC(2*idx1(n2)+1)),1,'first'); [nill,idxA]=max((1./(abs(board(tiC(2*idx1(n2)))-board(tempA))+1)).*(board(tiC(2*idx1(n2)))>=board(tempA)&board(tempA)>=board(tiC(2*idx1(n2)+1)))); if ~nill,idxA=[]; end end if ~isempty(idxA)&&idxA=board(tempB)&board(tempB)>=board(tiC(2*idx1(n2))),1,'first'); case {3,4}, idxB=find(board(tempA(idxA))>=board(tempB)&board(tempB)>=board(tiC(2*idx1(n2)+1)),1,'first'); end if ~isempty(idxB)&&idxB=board(tiC(1))); if isempty(idx1), mask=tP; D=Dijkstra(mask,tiC(1)); idx1=find(mask&D-1<=limit&board>=board(tiC(1))); if isempty(idx1),ok=0;end end if ok [nill,idx2]=min(board(idx1)+1e-3*D(idx1)); %%%t3 idx0=idx1(idx2); limit=limit-(D(idx0)-1); for n2=D(idx0)-1:-1:1 if D(idx0+1)==n2, tidx0=idx0+1; elseif D(idx0-1)==n2, tidx0=idx0-1; elseif D(idx0+m)==n2, tidx0=idx0+m; elseif D(idx0-m)==n2, tidx0=idx0-m; end moves=[moves; idx0, tidx0]; board(tidx0)=board(idx0); board(idx0)=0; idx0=tidx0; end tiC=[idx0,tiC]; tP(idx0)=false; end case 2 mask=tP; D=Dijkstra(mask,tiC(end)); idx1=find(mask&D-1<=limit&board<=board(tiC(end))&board>0); if isempty(idx1),ok=0;end if ok [nill,idx2]=min(-board(idx1)+1e-3*D(idx1)); %%%t3 idx0=idx1(idx2); limit=limit-(D(idx0)-1); for n2=D(idx0)-1:-1:1 if D(idx0+1)==n2, tidx0=idx0+1; elseif D(idx0-1)==n2, tidx0=idx0-1; elseif D(idx0+m)==n2, tidx0=idx0+m; elseif D(idx0-m)==n2, tidx0=idx0-m; end moves=[moves; idx0, tidx0]; board(tidx0)=board(idx0); board(idx0)=0; idx0=tidx0; end tiC=[tiC,idx0]; tP(idx0)=false; end end end function iC = bellman_pixel(board,limit,boardab) [m,n]=size(board); cD=board; IDX=zeros(m,n); touched=board>0; valid=touched; while any(touched(:)) touch=find(touched); touched(touch)=false; for dx=[1,-1,m,-m], cDnew=board(touch+dx)+cD(touch)+0*min(limit,boardab(touch+dx)*(n+m)/2).*board(touch+dx)/((n+m)/2); %%% idx=find(board(touch+dx)>board(touch)&cDnew>cD(touch+dx)); touched(touch(idx)+dx)=valid(touch(idx)+dx); cD(touch(idx)+dx)=cDnew(idx); IDX(touch(idx)+dx)=touch(idx); end end [nill,idx]=max(cD(:)); iC=zeros(1,m*n); for n1=1:m*n, if idx>0 iC(n1)=idx; idx=IDX(idx); else break; end end iC=iC(1:n1-1); end function iC = bellman(C,D) N=size(C,1); c=cell(1,N); for n1=1:N,c{n1}=find(C(:,n1)); end %C=full(C); IDX=zeros(N,1); cD=D; touched=true(1,N); while any(touched) for touch=find(touched), touched(touch)=false; idx=c{touch};%find(C(:,touch)); cDnew=D(idx)+cD(touch); idx2=find(cD(idx)0 iC(n1)=idx; idx=IDX(idx); else break; end end iC=iC(1:n1-1); end function iC=bellman_postprocess(ClustersOrder,IdxList,IdxSegments,BoardCluster) [m,n]=size(BoardCluster); iC=[]; for n1=1:numel(ClustersOrder) idx=IdxList(IdxSegments(ClustersOrder(n1)):IdxSegments(ClustersOrder(n1)+1)-1); if numel(idx)==1, iC=[iC,idx(1)]; else % longest shortest-path ThisCluster=false(m,n); ThisCluster(idx)=true; D=Dijkstra(ThisCluster,iC); if n1==numel(ClustersOrder) temp=D(ThisCluster); else temp=D(ThisCluster).*(BoardCluster([false(1,n);ThisCluster(1:m-1,:)])==ClustersOrder(n1+1)|BoardCluster([ThisCluster(2:m,:);false(1,n)])==ClustersOrder(n1+1)|BoardCluster([false(m,1),ThisCluster(:,1:n-1)])==ClustersOrder(n1+1)|BoardCluster([ThisCluster(:,2:n),false(m,1)])==ClustersOrder(n1+1)); end [nill,idx0]=max(temp(:)); idx0=idx(idx0); E=zeros(2,D(idx0)-1); tiC=zeros(1,D(idx0)); tiC(end)=idx0; for n2=D(idx0)-1:-1:1 if D(idx0+1)==n2, idx0=idx0+1; E(1,n2)=1;E(2,n2)=m; elseif D(idx0-1)==n2, idx0=idx0-1; E(1,n2)=1;E(2,n2)=m; elseif D(idx0+m)==n2, idx0=idx0+m; E(1,n2)=m;E(2,n2)=1; elseif D(idx0-m)==n2, idx0=idx0-m; E(1,n2)=m;E(2,n2)=1; end tiC(n2)=idx0; end % now grow it tP=ThisCluster; tP(tiC)=false; ok=1; while ok ok=0; K=size(tiC,2); idx1=find(tP(tiC(1:2:K-1)+E(2,1:2:K-1))&tP(tiC(2:2:K)+E(2,1:2:K-1))); for n2=numel(idx1):-1:1 ok=1; tiC=[tiC(1:2*idx1(n2)-1),tiC(2*idx1(n2)-1)+E(2,2*idx1(n2)-1), tiC(2*idx1(n2))+E(2,2*idx1(n2)-1), tiC(2*idx1(n2):end)]; E=[E(:,1:2*idx1(n2)-2),E([2,1],2*idx1(n2)-1), E([1,2],2*idx1(n2)-1), E([2,1],2*idx1(n2)-1), E(:,2*idx1(n2):end)]; tP(tiC)=false; end K=size(tiC,2); idx1=find(tP(tiC(1:2:K-1)-E(2,1:2:K-1))&tP(tiC(2:2:K)-E(2,1:2:K-1))); for n2=numel(idx1):-1:1 ok=1; tiC=[tiC(1:2*idx1(n2)-1),tiC(2*idx1(n2)-1)-E(2,2*idx1(n2)-1), tiC(2*idx1(n2))-E(2,2*idx1(n2)-1), tiC(2*idx1(n2):end)]; E=[E(:,1:2*idx1(n2)-2),E([2,1],2*idx1(n2)-1), E([1,2],2*idx1(n2)-1), E([2,1],2*idx1(n2)-1), E(:,2*idx1(n2):end)]; tP(tiC)=false; end K=size(tiC,2); idx1=find(tP(tiC(2:2:K-1)+E(2,2:2:K-1))&tP(tiC(3:2:K)+E(2,2:2:K-1))); for n2=numel(idx1):-1:1 ok=1; tiC=[tiC(1:2*idx1(n2)),tiC(2*idx1(n2))+E(2,2*idx1(n2)), tiC(2*idx1(n2)+1)+E(2,2*idx1(n2)), tiC(2*idx1(n2)+1:end)]; E=[E(:,1:2*idx1(n2)-1),E([2,1],2*idx1(n2)), E([1,2],2*idx1(n2)), E([2,1],2*idx1(n2)), E(:,2*idx1(n2)+1:end)]; tP(tiC)=false; end K=size(tiC,2); idx1=find(tP(tiC(2:2:K-1)-E(2,2:2:K-1))&tP(tiC(3:2:K)-E(2,2:2:K-1))); for n2=numel(idx1):-1:1 ok=1; tiC=[tiC(1:2*idx1(n2)),tiC(2*idx1(n2))-E(2,2*idx1(n2)), tiC(2*idx1(n2)+1)-E(2,2*idx1(n2)), tiC(2*idx1(n2)+1:end)]; E=[E(:,1:2*idx1(n2)-1),E([2,1],2*idx1(n2)), E([1,2],2*idx1(n2)), E([2,1],2*idx1(n2)), E(:,2*idx1(n2)+1:end)]; tP(tiC)=false; end end iC=[iC,tiC]; end end end function B = neighb(A,V,nA) [m,n] = size(A); B=sparse(A(:,1:n-1),A(:,2:n),double(V(:,1:n-1)>V(:,2:n)),nA,nA)+sparse(A(:,2:n),A(:,1:n-1),double(V(:,2:n)>V(:,1:n-1)),nA,nA)+sparse(A(1:m-1,:),A(2:m,:),double(V(1:m-1,:)>V(2:m,:)),nA,nA)+sparse(A(2:m,:),A(1:m-1,:),double(V(2:m,:)>V(1:m-1,:)),nA,nA); B=B>0; end function [B,C,p,r,nc] = connected(A) [m,n] = size(A); N = m*n ; K = reshape (1:N, m, n) ; V = A(:,1:n-1) == A(:,2:n); H = A(1:m-1,:) == A(2:m,:); G = sparse(K([V,false(m,1)]),K([false(m,1),V]),1,N,N)+sparse(K([H; false(1,n)]),K([false(1,n); H]),1,N,N); G = G + G' + speye(N); [p q r] = dmperm(G); nc = numel(r) - 1; C = A(p(r(1:nc))); B = ones(m, n); for a = 2:nc B(p(r(a):r(a+1)-1)) = a; end end function D = Dijkstra(A,i1) [m,n]=size(A); mn=m*n; D=inf(m,n); D(~A)=0; if isempty(i1), i1=find(A,1,'first'); mode=0; else i1=i1(end); mode=1; % if A(i1+1); i1=i1+1; % elseif A(i1-1); i1=i1-1; % elseif A(i1+m); i1=i1+m; % elseif A(i1-m); i1=i1-m; % end end D(i1)=1; for n1=2:mn, idx1=D(i1+1)>n1; idx2=D(i1-1)>n1; idx3=D(i1+m)>n1; idx4=D(i1-m)>n1; %i1=unique([i1(idx1)+1,i1(idx2)-1,i1(idx3)+m,i1(idx4)-m]); X=false(m,n);X(i1(idx1)+1)=true;X(i1(idx2)-1)=true;X(i1(idx3)+m)=true;X(i1(idx4)-m)=true; i1=find(X); if isempty(i1), break; end D(i1)=n1; end if mode, D(D>0)=D(D>0)-1; end end```