function [B,N,u,s,v]=svdConv2(A,K,N,shape,u,s,v)
% svdConv2 <fast svd 2D convolution>
% Usage:: [B,[N]|[u,s,v]]=svdConv2(A,[K],N[1],shape['same'],[u,s,v])
% Usage:: B=svdConv2(A,K,[N[1]|tol],shape{'full'|['same']|'valid'})
% Uses [u,s,v]=svd(K) to calculate approximate 2D convolution
% if N is int from [0:rank[K]] then N is # of u'v terms
% if tol is not int then calc N=find(s<tol);
% K is optional if u,s,v specified
% 
%
% conv2(A,K)~=sum[n=1:N] conv2(conv2(A,s(n)*u(:,n)),v(:,n)')


% revision history
% 08/05/1992 Mark D. Shattuck <mds> fastConv2.m from fastConv2D.c
% 08/22/1992 mds add proper svd to remove minimization step 
% 06/21/1996 mds frk fastConv2=>svdConv2.m and fftConv2.m
%                    remove fft version
% 04/23/1997 mds add default args
% 06/17/2001 mds add N vs tol N=find(s<tol)
% 10/12/2003 mds add proper shape support
% 04/04/2004 mds rmv rank calc 
% 12/26/2015 mds add some comments to top
% 01/22/2016 mds add N|{u,s.v} output (Not finished!)
%

%% Parse Input
if(~exist('N','var') || isempty(N))
	N=1;
end

if(~exist('shape','var') || isempty(shape))
	shape='same';
end

if(~exist('u','var') || isempty(u))
	[u,s,v]=svd(K);
end

if N~=fix(N)
  tol=N;
  N=find(diag(s)<tol*s(1),1);
end

%% Main

B=0;
for n=1:N
  B=B+conv2(conv2(A,s(n,n)*u(:,n),shape),v(:,n)',shape);
end 

