lans_gradient.m

来自「模式识别工具包」· M 代码 · 共 110 行

%	lans_gradient	- Compute the gradient of a principal curve
%
%	[dfdx,d2fdx2] = lans_gradient(psurf,para)
%
%	_____OUTPUT_____________________________________________________________
%	dfdx	approximated gradient				(col vectors)
%	d2fdx2	approximated 2nd derivative			(col vectors)
%
%	_____INPUT______________________________________________________________
%	psurf	current principal curve				(structure)
%		see lans_psurf.m
%	para	see lanspara.m paraget.m			(string)
%		-gradient				
%			'right'		right difference	(default)
%			'left'		left difference			
%			'avg'		average of left and right		
%			'wavg		weighted by segment length
%	
%	_____SEE ALSO___________________________________________________________
%	lans_plotgrad (for demo)
%
%	_____NOTES______________________________________________________________
%	- ONLY works for 1-D manifolds
%	- returns N gradients with the end point gradients approximated by that
%	  of the previous point
%	- vector difference less than a threshold will not be considered
%	- threshold currently set to
%		1e-4
%
%	(C) 1999.04.30 Kui-yu Chang
%	http://lans.ece.utexas.edu/~kuiyu

%	This program is free software; you can redistribute it and/or modify
%	it under the terms of the GNU General Public License as published by
%	the Free Software Foundation; either version 2 of the License, or
%	(at your option) any later version.
%
%	This program is distributed in the hope that it will be useful,
%	but WITHOUT ANY WARRANTY; without even the implied warranty of
%	MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
%	GNU General Public License for more details.
%
%	You should have received a copy of the GNU General Public License
%	along with this program; if not, write to the Free Software
%	Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307, USA
%	or check
%			http://www.gnu.org/

function [dfdx,d2fdx2] = lans_gradient(psurf,para)

[D M]	= size(psurf.f);

threshold	= 1e-4;

if nargin<2
	para	= [];
	if nargin<1
		clc;
		help lans_gradient;
		break;
	end
end

%-----	parameters
gtype		= paraget('-gradient',para);

f_diff		= diff(psurf.f')';

df_right	= [f_diff f_diff(:,end)];
df_left		= [f_diff(:,1) f_diff];

x_diff		= diff(psurf.x);

%_____	find ultra small gradient vector in EACH dimension of f and tag them
candidate	= sum(f_diff);
wzero		= find(candidate<threshold);
x_diff(wzero)	= ones(size(wzero));
f_diff(:,wzero)	= zeros(D,length(wzero));

dx_right	= [x_diff x_diff(end)];
dx_left		= [x_diff(1) x_diff];

dx_rightblock	= ones(D,1)*dx_right;
dx_leftblock	= ones(D,1)*dx_left;

switch(gtype)
	case 'left'
		dfdx		= df_left./dx_leftblock;

	case 'right'
		dfdx		= df_right./dx_rightblock;

	case 'avg'
		dfdx_left	= df_left./dx_leftblock;
		dfdx_right	= df_right./dx_rightblock;
		dfdx		= .5*(dfdx_left + dfdx_right);

	case 'wavg'		%weights prop. dist to end
		lambda		= dx_leftblock./(dx_leftblock+dx_rightblock);
		dfdx_left	= df_left./dx_leftblock;
		dfdx_right	= df_right./dx_rightblock;
		dfdx		= lambda.*dfdx_left + (1-lambda).*dfdx_right;
end

if nargout>1		% compute 2nd derivative
	gsurf.x	= psurf.x;
	gsurf.f	= dfdx;
	d2fdx2	= lans_gradient(gsurf,para);
end