hybridsynthesizerec.m

histogram based code for image retrival

    % a deterministic, reproduceable synthesis result

    % if (min_error == 0),
    %  disp('warning: hardcoded first patch!!!');
    %  bestpos = [128 128];        
    % end
    
    % -------------------------------------------------------------------------------------
    % ALGORITHM STEP: shift T by bestpos and compute the error surface S between I and T
    %                 mse -> mean squared error, ssd -> sum of square differences
    % -------------------------------------------------------------------------------------
    shiftT = circshift(data.T,-bestpos);
    [S,bestmatch_mse,bestmatch_ssd,numpixels] =...
        errorsurface(shiftT,grownpatchsize(1),grownpatchsize(2),J,I);

    % use mean squared error (or alternatively 'metric' error) to determine 
    % when to split patches
    patch_overlap_error = bestmatch_mse;
    % patch_overlap_error = min_error;
    
    % we currently stop splitting at patchsize defined in 'data' input structure
    finalsplit = 2*data.min_patchsize;
    if (patch_overlap_error > data.patch_delta_max & ...
            psize(1) >= finalsplit & psize(2) >= finalsplit),
        % ---------------------------------------------------------------------------------
        % ALGORITHM STEP: if patch_overlap_error too large, 
        %                 split (if possible), and recurse
        % ---------------------------------------------------------------------------------
        % make four patches out of this one patch by creating five new vertices and
        % four faces from the 9 vertices, see GenUniformQuadMeshWithOffset.m
        % NOTE: this implementation only works with regular quad patches!!!
        row_offset = patchlist.BB(pos,1) - 1;
        col_offset = patchlist.BB(pos,2) - 1;
        row_psize = psize(1)/2;
        col_psize = psize(2)/2;
        % also split pixel overlap for the next level of patches. NOTE:
        % we want a minimum of 3 lapping pixels (equivalent to a 7x7 nbhd for 1 pixel)
        overlap_new = max(3,ceil(overlap/2));

        [splitpatchlist.V,splitpatchlist.F,splitpatchlist.BB] = ...
            genquadmeshoffset(row_psize,col_psize,2,2,row_offset,col_offset);
        sub_error = hybridsynthesizerec(data,overlap_new,splitpatchlist);
        
        % accumulate the mse for this call
        error_accum = error_accum + sub_error;

        % some output
        if (data.verbatim),
            disp(sprintf('splitting. error: %f, patchsize: %d->%d, sub_error: %f',...
                bestmatch_mse, row_psize*2, row_psize, sub_error));
        end
    else
        
        % accumulate the error for this call
        error_accum = error_accum + bestmatch_mse;
        
        % ---------------------------------------------------------------------------------
        % ALGORITHM STEP: build alpha blending mask to feather valid pixels
        % ---------------------------------------------------------------------------------
        numalphasteps = overlap; % blend over entire overlap width
        if (data.feather == 0), numalphasteps = 0; end
        alphamask = buildalphamask(PS,J,grownpatchsize,numalphasteps);

        % ---------------------------------------------------------------------------------
        % ALGORITHM STEP: build validpixel mask and 'J OR PS' surface (=PSall)
        %                 PSall is original patch + area of grown patch which overlaps
        %                 already synthesized pixels
        % ---------------------------------------------------------------------------------
        % now we must use the entire patch plus the mask area as the lapping mask, set 
        % the validpixel-mask to 1 in errorvalid regions, and blackout invalid regions 
        % for later filling. note: a setting of 0 for data.pixel_delta_max results in 
        % simple feathering (if data.feather == 1), or abrupt borders
        PSall = PS;
        validpixels = PS;
        for jj=1:grownpatchsize(1),
            for ii=1:grownpatchsize(2),
                if (J(jj,ii) == 1 & PS(jj,ii) == 0),
                    PSall(jj,ii) = 1;
                    if (S(jj,ii) <= data.pixel_delta_max), validpixels(jj,ii) = 1; else
                        validpixels(jj,ii) = 0; end
                end
            end
        end
        
        % ---------------------------------------------------------------------------------
        % ALGORITHM STEP: build traversal map for the per-fragment resynthesis in J
        % ---------------------------------------------------------------------------------
        pixeltraversalmap = buildtraversalmap(GPS,PS,J,validpixels,grownpatchsize);
    
        % ---------------------------------------------------------------------------------
        % ALGORITHM STEP: compose the best match with the output and pre images
        %                 and prepare 'out' and 'pre' for per-fragment synthesis
        % ---------------------------------------------------------------------------------
        % use empty image ('bestpatch') of same size as output image, paste the best 
        % patch from texture example in upper left, circularly shift it by values stored in 
        % 'shift+deltashift' and add it to 'out' image (and 'pre' image)
    
        % build a copy of (flattened) indices into T for the source_map by
        % shifting a copy of the index_map so it's congruent with shiftT
        shift_index_map = circshift(data.index_map,-bestpos);
        best_src_indices = zeros(data.rows, data.cols);
        
        % first clear the 'bestpatch' image again
        bestpatch = zeros(data.rows, data.cols, 3); 
        % build best patch in 'bestpatch'. concurrently build the source indices.
        for jj=1:grownpatchsize(1),
            for ii=1:grownpatchsize(2),
                % if we have patch support here...
                if (PSall(jj,ii) == 1), 
                    % ...add pixel from shifted T to 'bestpatch'. 
                    % note: shiftT was already used to compute error surface S
                    bestpatch(jj,ii,1:3) = shiftT(jj,ii,1:3); 
                    % and add source index
                    best_src_indices(jj,ii) = shift_index_map(jj,ii);
                end
            end
        end
    
        % shift this best patch 'bestpatch' to final destination and compose with 
        % 'out' and 'pre' based on pixel validity and alpha values
        bestpatch = circshift(bestpatch,-(shift+deltashift));
        best_src_indices = circshift(best_src_indices,-(shift+deltashift));
        PSall = circshift(PSall,-(shift+deltashift));
        validpixels = circshift(validpixels,-(shift+deltashift));
        alphamask = circshift(alphamask,-(shift+deltashift));
        pixeltraversalmap = circshift(pixeltraversalmap,-(shift+deltashift));
  
        % for all pixels within support of PS or J (= PSall)
        for jj=1:data.rows,
            for ii=1:data.cols,
                % replace pixel values where validpixels is 1, set to value 
                % (PER_FRAGMENT_SYNTH 0 0) where validpixels is 0 (in pre)
                % reminder: PSall is ungrown patch PLUS mask support (overlap)
                if (PSall(jj,ii) == 1),
                    ispixelvalid = validpixels(jj,ii);
                    alpha = alphamask(jj,ii);
                    if (ispixelvalid == 1),
                        % VALID PIXEL
                        % compose final pixel into out and pre, using basic OVER compositing
                        out(jj,ii,1:3) = (1-alpha)*out(jj,ii,1:3) + alpha*bestpatch(jj,ii,1:3);
                        pre(jj,ii,1:3) = (1-alpha)*pre(jj,ii,1:3) + alpha*bestpatch(jj,ii,1:3); 
                        
                        % build source_map. if alpha > 0.5, use source_index from 
                        % new patch, otherwise retain old source value (= do nothing)
                        if (alpha > 0.5),
                            source_map(jj,ii) = best_src_indices(jj,ii);
                        end
                    
                    else
                        % INVALID PIXEL
                        % initialize out and pre for per-fragment synthesis
                        % pre now has (NOT_YET_SYNTHESIZED 0 0) where no synthesis has 
                        % taken place, and (PER_FRAGMENT_SYNTH 0 0) where pixels are to 
                        % be resynthesized per-fragment in this step 
                        out(jj,ii,1:3) = 0; 
                        % use (red) color gradient for displayed, intermediate result
                        factor = 1/pixeltraversalmap(jj,ii);
                        out(jj,ii,1) = min(factor*2.5,1);            
                        pre(jj,ii,2:3) = 0; pre(jj,ii,1) = PER_FRAGMENT_SYNTH;
                    end
                end
            end
        end
        
        % ---------------------------------------------------------------------------------
        % ALGORITHM STEP: fill invalid pixels using the traversal map
        % ---------------------------------------------------------------------------------
        % fill holes where invalid lapping pixels exist in traversal order 

        % if we pass in a knn-datastructure (globally defined, see hybridsynthesize.m)
        % for T, use it. otherwise default to exhaustive (fourier-based) search.
        if (knn ~= 0),
            % do ashikhmin/k-coherence search
            % using box-shaped neigborhood defined in 'data'
            overlap_resynthesis_knn(data,pixeltraversalmap);
        else
            % perform exhaustive, fft-based search 
            % using box-shaped neigborhood defined in 'data'
            overlap_resynthesis_exh(data,pixeltraversalmap);
        end

        % always display the current result (patch-progress)
        figure(1); imshow(out,'truesize')
    end
    
    % mse for this call: divide accumulated mse by the number of patches in this call
    error = error_accum/length(patchlist.F);
end