www.gusucode.com > matlab非局部均值工具箱 > matlab非局部均值工具箱/matlab非局部均值工具箱/toolbox_nlmeans/toolbox/compute_patch_library.m
function [H,X,Y] = compute_patch_library(M,w,options) % [H,X,Y] = compute_patch_library(M,w,options); % % M is the texture % w is the half-width of the patch, so that each patch % has size (2*w+1,2*w+1,s) where s is the number of colors. % % H(:,:,:,i) is the ith patch (can be a color patch). % X(i) is the x location of H(:,:,:,i) in the image M. % Y(i) is the y location of H(:,:,:,i) in the image M. % % options.sampling can be set to 'random' or 'uniform' % If options.sampling=='random' then you can define % options.nbr as the number of random patches and/or % option.locations_x and option.locations_y. % % If w is a vector of integer size, then the algorithm % compute a set of library, one for each size, % and use the same location X,Y for each patch. % H{k} is the k-th library. % % You can define options.wmax to avoid too big patch. % In this case, if w>wmax, the texture M will be filtered, and the patches % will be subsampled to match the size. % % Copyright (c) 2006 Gabriel Peyr options.null = 0; if length(w)>1 % process the library of patches X = []; Y = []; for i=1:length(w) if ~isempty(X) option.locations_x = X; option.locations_y = Y; [H{i},X,Y] = compute_patch_library(M,w(i),options); end end return; end if isfield(options, 'sampling') sampling = options.sampling; else sampling = 'random'; end if isfield(options, 'wmax') wmax = options.wmax; else wmax = 9; end wmax = min(w,wmax); m = size(M,1); n = size(M,2); s = size(M,3); if w>wmax % perform some smoothing to avoid aliasing when subsampling sigma = w/wmax; % in pixel size n1 = max( 16, round(n/8)*2+1 ); h = compute_gaussian_filter( [1 1]*n1,sigma/(2*n),[1 1]*n); M = perform_convolution(M,h); end % do some padding to avoid boundary problems M = symmetric_extension(M,w); ww = 2*w+1; ww1 = 2*wmax+1; if isfield(options, 'nbr') && strcmp(sampling, 'random') p = options.nbr; else if strcmp(sampling, 'random') p = 2000; else p = n*m; end end if strcmp(sampling, 'random') if ~isfield(options, 'locations_x') X = w + 1 + floor(rand(p,1)*m); Y = w + 1 + floor(rand(p,1)*n); else X = options.locations_x(:)+w; Y = options.locations_y(:)+w; p = length(X); end else [Y,X] = meshgrid(w+1:w+n, w+1:w+m); X = X(:); Y = Y(:); end H = zeros(ww1,ww1,s,p); B = zeros(ww1,ww1,s); if mod(w/wmax,1)==0 %%% in this case, a fast sampling can be used %%% [dY,dX] = meshgrid(-w:w/wmax:w,-w:w/wmax:w); Xp = repmat( reshape(X,[1,1,1,p]) ,[ww1 ww1 s 1]) + repmat(dX,[1 1 s p]); Yp = repmat( reshape(Y,[1,1,1,p]) ,[ww1 ww1 s 1]) + repmat(dY,[1 1 s p]); Cp = repmat( reshape(1:s,[1 1 s]), [ww1 ww1 1 p]); I = sub2ind(size(M), Xp,Yp,Cp); H = M(I); % remove padding in sampling location X = X-w; Y = Y-w; return; end %%% in this case, interpolation is needed %%% xi = linspace(1,ww,ww1); [Yi,Xi] = meshgrid(xi,xi); h = waitbar(0,['Extracting patches of size ' num2str(w) ' ...']); for i=1:p waitbar(i/p,h); x = X(i); y = Y(i); A = M( x-w:x+w, y-w:y+w,: ); % need to perform interpolation for t=1:s B(:,:,t) = interp2(A(:,:,t),Xi,Yi)'; end % perform sub-sampling H(:,:,:,i) = B; end close(h); % remove padding in sampling location X = X-w; Y = Y-w;