texture_synth.m

this is a very very very nice code

%In is the input image, mask is the region to be synthesized, and prm is a
%parameter file. prm has three components, the similarity threshold,
%thresh, which defines how similar a sample must be to be considered a
%candidate for synthesization. W is a parameter that defines the size of
%the neighborhood to compare between a target region and a candidate.
%Larger values of W will allow to maintain large-scale textures, while
%small values of W will maintain small-scale textures. Sigma is the size of
%the gaussian kernel used when comparing neighborhoods, and skip is a value
%used to define how densely to sample the input image. To have a good guess
%for the value of skip computed automatically, simply define the value of
%skip to be 0.
%
%Note that mask is a binary image, where 0's represent the region to be
%synthesized. Within the algorithm though the reverse is true, 1's
%represent the regions to be synthesized, so the mask is inverted at the
%beginning of the algorithm. This is to make the algorithm follow the same
%convention used in the inpainting algorithm
function out = texture_synth(in, mask, prm)
if sum(class(prm) == 'struct') < 6
    prm = struct('thresh', 150, 'W', 17, 'sigma', 6, 'skip', 4);
end

%If you are unsure of what value of skip to use, plug in 0, and this will
%compute a reasonably good value for you
if prm.skip <= 0
    prm.skip = floor(sqrt(size(in,1)*size(in,2)/650))-2;
end

%Make sure we have just a single-channel image...
if size(in,3) > 1
    in = in(:,:,1);
end
in = double(in);

%I wrote this code originally using the inverted form of the mask, ie, ones
%specify where to inpaint, and 0's everywhere else. So I simply invert the
%mask to make to get this effect.
mask = ~mask;

%Create the full parameter file
param = struct('thresh', prm.thresh, 'W', prm.W, 'wh', 0, 'gauss', 0, 'sigma', prm.sigma, 'patches', 0, 'skip', prm.skip, 'fast', 1);
param.wh = (param.W-1)/2;

%Define the gaussian kernel
gauss = zeros(param.W, param.W);
sigs = param.sigma ^ 2;
wh = param.wh + 1;
for i = 1:param.W
    for j = 1:param.W
        gauss(i,j) = exp(-((i-wh)^2 + (j-wh)^2) / (2*sigs)) / sqrt(2*3.1415*sigs);
    end
end
gauss = gauss / sum(gauss(:));
param.gauss = gauss;

out = in;
M = size(out, 1);
N = size(out, 2);
wh = param.wh;
skip = param.skip;

%Sample the input image, these will be considered later for candidate regions
nm = 0;
for y = wh+1:skip:M-wh
    for x = wh+1:skip:N-wh
        if sum(sum(mask(y-wh:y+wh, x-wh:x+wh))) == 0
            nm = nm + 1;
            region = out(y-wh:y+wh, x-wh:x+wh);
            patches(nm, :) = region(:)';
        end
    end
end
param.patches = patches;
disp([num2str(size(patches,1)), ' sample patches collected'])

%For every pixel to by synthesized, starting from the border of the mask
%and eroding towards the center...
pix_set = get_pix_set(~mask);
for K = 1:size(pix_set,1)
    y = pix_set(K,1);
    x = pix_set(K,2);
    if ~(y < wh+1 | x < wh+1 | y > size(mask,1)-wh | x > size(mask,2)-wh) %Cant handle borders
        %Define the local mask of synthesized and unsynthesized pixels
        loc_mask = ~mask(y-wh:y+wh, x-wh:x+wh);
        %Gather the pixels into a target patch
        target = out(y-wh:y+wh, x-wh:x+wh);
        %Get a matching patch and fill the current pixel with the current pixel in that matching patch
        out(y, x) = get_match(target, loc_mask, param);
        mask(y, x) = 0;
    end
end