--- a/extra/nurbs/INDEX	Thu May 09 19:37:51 2013 +0000
+++ b/extra/nurbs/INDEX	Fri May 10 12:53:12 2013 +0000
@@ -44,6 +44,7 @@
  surfderiveval
 NURBS geometric entities and functions
  nrbbasisfun
+ nrbmeasure
  nrbbasisfunder
  nrbnumbasisfun
  nrbcrvderiveval

--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/extra/nurbs/inst/nrbmeasure.m	Fri May 10 12:53:12 2013 +0000
@@ -0,0 +1,90 @@
+%% Copyright (C) 2013 Carlo de Falco
+%% 
+%%   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 3 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 Octave; see the file COPYING.  If not, see
+%%   <http://www.gnu.org/licenses/>.  
+
+%% -*- texinfo -*- 
+%% @deftypefn {Function File}{[@var{dist}, @var{ddds}, @var{ddde}]=} 
+%% nrbmeasure (@var{line}, @var{start}, @var{end}, @var{tol}) 
+%% Compute the distance @{dist} between the points @var{start} and @{end} along the
+%% line @var{line}. 
+%% Use the tolerance @var{tol} when computing integrals. The output
+%% parameters  @var{ddds}, @var{ddde} are the derivatives of
+%% @var{dist} with respect to @var{start} and @var{end} respectively.
+%% @seealso{}
+%% @end deftypefn
+
+%% Author: Carlo de Falco 
+
+
+function [dist, ddistds, ddistde] = nrbmeasure (nrbline, s, e, tol)
+  
+  if (nargin < 4)
+    tol = 1e-6;
+    if (nargin < 3)
+      e = 1;
+      if (nargin < 2)
+        s = 0;
+      end
+    end
+  end
+
+  if (ismatrix (s) && isscalar (e))
+    e = repelems (e, size(s) .');
+  elseif (ismatrix (e) && isscalar (s))
+    s = repelems (s, size(e) .');
+  end
+
+  [ders, ders2] = nrbderiv (nrbline);
+  dist = quad (@(u) len (u, nrbline, ders), s, e, tol);
+  
+  if (nargout > 1)
+    ddistds = -len (s, nrbline, ders);
+    if (nargout > 2)
+      ddistde = +len (e, nrbline, ders);
+    end
+  end
+
+end
+
+function l = len (u, nrbline, ders)
+  [ignore, d] = nrbdeval (nrbline, ders, u);
+  f = d(1, :);
+  g = d(2, :);
+  h = d(3, :);
+  l = sqrt (f.^2 + g.^2 + h.^2);
+end
+
+%!test 
+%! c = nrbcirc (1, [0 0], 0,  pi/3);
+%! l = nrbmeasure(c, 0, 1, 1e-7);
+%! assert (l, pi/3, 1e-7)
+
+%!test 
+%! c = nrbcirc (1, [0 0], 0,  pi/2);
+%! s = zeros (1, 100); e = linspace (0, 1, 100);
+%! for ii = 1:100
+%!   l(ii) = nrbmeasure (c, s(ii), e(ii), 1e-7);
+%! endfor
+%! xx = nrbeval (c, e);
+%! theta = atan2 (xx(2,:), xx(1,:));
+%! assert (l, theta, 1e-7)
+
+%!test 
+%! c = nrbcirc (1, [0 0], 0,  pi/2);
+%! s = 0; e = linspace (0, 1, 100);
+%! for ii = 1:100
+%!   l(ii) = nrbmeasure (c, s, e(ii), 1e-7);
+%! endfor
+%! l2 = nrbmeasure (c, s, e, 1e-7);