--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/scripts/geometry/delaunayn.m	Fri Aug 24 08:27:29 2007 +0000
@@ -0,0 +1,76 @@
+## Copyright (C) 2007  David Bateman
+##
+## This file is part of Octave.
+##
+## Octave 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, or (at your option)
+## any later version.
+##
+## Octave 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 Octave; see the file COPYING.  If not, write to the Free
+## Software Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
+## 02110-1301, USA.
+
+## -*- texinfo -*-
+## @deftypefn {Loadable Function} {@var{T} =} delaunayn (@var{P})
+## @deftypefnx {Loadable Function} {@var{T} =} delaunayn (@var{P}, @var{opt})
+## Form the Delaunay triangulation for a set of points.
+## The Delaunay triangulation is a tessellation of the convex hull of the
+## points such that no n-sphere defined by the n-triangles contains
+## any other points from the set.\n
+## The input matrix @var{P} of size [n, dim] contains n points in a space
+## of dimension dim. The return matrix @var{T} has the size [m, dim+1]. It
+## contains for each row a set of indices to the points, which describes a
+## simplex of dimension dim.  For example, a 2d simplex is a triangle and 3d
+## simplex is a tetrahedron.
+## 
+## Extra options for the underlying Qhull command can be specified by the
+## second argument. This argument is a cell array of strings. The default
+## options depend on the dimension of the input: 
+## 
+## @itemize 
+## @item  2D and 3D: @var{opt} = @{'Qt','Qbb','Qc'@}  
+## @item  4D and higher: @var{opt} = @{'Qt','Qbb','Qc','Qz'@} 
+## @end itemize 
+## 
+## If @var{opt} is [], then the default arguments are used. If @var{opt}
+## is @{'@w{}'@}, then none of the default arguments are used by Qhull. 
+## See the Qhull documentation for the available options. 
+## 
+## All options can also be specified as single string, for example
+## 'Qt Qbb Qc Qz'. 
+## 
+## @end deftypefn
+
+function t = delaunayn (x, varargin)
+  if (nargin < 1)
+    print_usage()
+  endif
+
+  t = __delaunayn__ (x, varargin {:});
+
+  ## Try to remove the zero volume simplices. The volume of the i-th simplex is
+  ## given by abs(det(x(t(i,1:end-1),:)-x(t(i,2:end),:)))/prod(1:n) 
+  ## (reference http://en.wikipedia.org/wiki/Simplex). Any simplex with a 
+  ## relative volume less than some arbitrary criteria is rejected. The 
+  ## criteria we use is the volume of the simplex corresponding to an 
+  ## orthogonal simplex is equal edge length all equal to the edge length of 
+  ## the original simplex. If the relative volume is 1e3*eps then the simplex
+  ## is rejected. Note division of the two volumes means that the factor 
+  ## prod(1:n) is dropped.
+  idx = [];
+  [nt, n] = size(t);
+  for i = 1 : nt
+    X = x(t(i,1:end-1),:) - x(t(i,2:end),:);
+    if (abs (det (X)) /  sqrt (sum (X .^ 2, 2)) < 1e3 * eps)
+     idx = [idx, i];
+    endif
+  endfor
+  t(idx,:) = [];
+endfunction