function Aest=solver(imageSize,queryLimit)
Aest=128+zeros(imageSize);
if queryLimit==0,return;end
V=zeros(queryLimit,1);N=ones(queryLimit,1);C=ones(queryLimit,1);D=ones(queryLimit,1);K=zeros(queryLimit,4);
K(1,:)=[1,1,imageSize,imageSize];V(1)=queryImage(true(imageSize));N(1)=imageSize*imageSize;
for n1=2:queryLimit,
[maxC,idxC]=max(C(1:n1-1));
if ~maxC,queryLimit=n1-1;break;end
V0=V(idxC);N0=N(idxC);D0=D(idxC);K0=K(idxC,:);
if K0(3)>K0(4),K([idxC,n1],:)=[K0(1)+floor(K0(3)/2),K0(2),ceil(K0(3)/2),K0(4);K0(1:2),floor(K0(3)/2),K0(4)];
else K([idxC,n1],:)=[K0(1),K0(2)+floor(K0(4)/2),K0(3),ceil(K0(4)/2);K0(1:3),floor(K0(4)/2)]; end
M=false(imageSize);M(K(n1,1)+(0:K(n1,3)-1),K(n1,2)+(0:K(n1,4)-1))=true;
V(n1)=queryImage(M);N(n1)=K(n1,3)*K(n1,4);
V(idxC)=V0-V(n1);N(idxC)=N0-N(n1);
[C([idxC,n1]),D([idxC,n1])]=cost(V(idxC),V(n1),N(idxC),N(n1),D0);
end
X=V./N;I=zeros(imageSize);
for n1=1:queryLimit,
I(K(n1,1)+(0:K(n1,3)-1),K(n1,2)+(0:K(n1,4)-1))=n1;
end
Aest=X(I);
Aest2=Aest;
[x,y]=ndgrid(min(0:2*size(Aest2,1)-1,2*size(Aest2,1):-1:1));z=[1:size(Aest2,1),size(Aest2,1):-1:1];h=exp(-(x.^2+y.^2)/imageSize/imageSize);
for lazy=1:20,
fAest2=fft2([Aest2(z,z)]);
fAest2=real(ifft2(fAest2.*h));
d1=fAest2(1:end/2,1:end/2);
d1=max(1,min(255, d1*(d1(:)'*Aest2(:))/(d1(:)'*d1(:)) ));
X2=accumarray(I(:),Aest2(:),[queryLimit,1])./N;
d2=X2(I);
Aest2=max(1,min(255, Aest2+.95*(d1-Aest2+Aest-d2) ));
end
Aest=Aest2;
end
function [x,y]=cost(v1,v2,m1,m2,d0)
p1=v1/m1/255;p2=v2/m2/255;
y=d0/2+abs(p1-p2);
x=y.*[(m1-1).*p1.*(1-p1),(m2-1).*p2.*(1-p2)];
end
|
