Mercurial > forge
changeset 6910:97806581afd6 octave-forge
New version compatible with Matlab
line wrap: on
line diff
--- a/extra/nurbs/DESCRIPTION Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/DESCRIPTION Mon Mar 22 18:13:09 2010 +0000 @@ -1,7 +1,7 @@ Name: Nurbs -Version: 1.0.3 -Date: 2010-02-09 -Author: Mark Spink, Daniel Claxton, Carlo de Falco +Version: 1.1.0 +Date: 2010-03-19 +Author: Mark Spink, Daniel Claxton, Carlo de Falco, Rafael Vazquez Maintainer: Carlo de Falco Title: Nurbs. Description: Collection of routines for the creation, and manipulation of Non-Uniform Rational B-Splines (NURBS). @@ -12,4 +12,3 @@ Url: http://octave.sf.net SVNRelease: 6861 -
--- a/extra/nurbs/INDEX Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/INDEX Mon Mar 22 18:13:09 2010 +0000 @@ -1,41 +1,46 @@ nurbs >> Nurbs -NURBS curves and surfaces +Basic operations for NURBS curves, surfaces and volumes nrbmak - nrbtform nrbkntins nrbdegelev nrbderiv nrbdeval - nrbkntmult + nrbeval +Operations for constructing NURBS curves and surfaces + nrbtform nrbreverse nrbtransp nrbline nrbcirc nrbrect nrb4surf - nrbeval + nrbcylind nrbextrude nrbrevolve nrbruled nrbcoons nrbplot - nrbsrfgradient -NURBS basis functions - nrbbasisfun - nrbbasisfunder - nrbnumbasisfun - nrbbasisfungradient - nrbsurfderiveval + nrbtestcrv + nrbtestsrf B-Spline functions bspeval bspderiv bspkntins bspdegelev + basisfun basisfunder findspan + numbasisfun + tbasisfun +B-splines geometric entities curvederivcpts surfderivcpts - tbasisfun + surfderiveval +NURBS geometric entities and functions + nrbbasisfun + nrbbasisfunder + nrbnumbasisfun + nrbsurfderiveval Vector and Transformation Utilities vecnorm vecmag
--- a/extra/nurbs/inst/basisfun.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/basisfun.m Mon Mar 22 18:13:09 2010 +0000 @@ -1,18 +1,3 @@ -%% Copyright (C) 2003 Mark Spink, 2007 Daniel Claxton, 2009 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 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, see <http://www.gnu.org/licenses/>. - function B = basisfun (iv, uv, p, U) % BASISFUN: Basis function for B-Spline @@ -23,20 +8,35 @@ % % Calling Sequence: % -% N = basisfun(i,u,p,U) +% N = basisfun(iv,uv,p,U) % % INPUT: % -% i - knot span ( from FindSpan() ) -% u - parametric point -% p - spline degree -% U - knot sequence +% iv - knot span ( from FindSpan() ) +% uv - parametric points +% p - spline degree +% U - knot sequence % % OUTPUT: % % N - Basis functions vector(numel(u)x(p+1)) % -% Algorithm A2.2 from 'The NURBS BOOK' pg70. +% Adapted from Algorithm A2.2 from 'The NURBS BOOK' pg70. +% +% Copyright (C) 2000 Mark Spink, 2007 Daniel Claxton, 2009 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 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, see <http://www.gnu.org/licenses/>. B = zeros(numel(uv), p+1); @@ -67,6 +67,7 @@ end +end %!test %! n = 3;
--- a/extra/nurbs/inst/basisfunder.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/basisfunder.m Mon Mar 22 18:13:09 2010 +0000 @@ -1,31 +1,16 @@ -%% Copyright (C) 2009 Rafael Vazquez -%% -%% 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, see <http://www.gnu.org/licenses/>. - function dersv = basisfunder (ii, pl, uu, u_knotl, nders) % BASISFUNDER: B-Spline Basis function derivatives. % % Calling Sequence: % -% ders = basisfunder (i, pl, u, k, nd) +% ders = basisfunder (ii, pl, uu, k, nd) % % INPUT: % -% i - knot span +% ii - knot span % pl - degree of curve -% u - parametric points +% uu - parametric points % k - knot vector % nd - number of derivatives to compute % @@ -33,6 +18,22 @@ % % ders - ders(n, i, :) (i-1)-th derivative at n-th point % +% Adapted from Algorithm A2.3 from 'The NURBS BOOK' pg72. +% +% Copyright (C) 2009 Rafael Vazquez +% +% 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, see <http://www.gnu.org/licenses/>. for jj = 1:numel(uu) @@ -108,6 +109,8 @@ end +end + %!test %! k = [0 0 0 0 1 1 1 1]; %! p = 3; @@ -138,7 +141,7 @@ %! sumders = sum (squeeze(ders(ii,:,:)), 2); %! assert (sumders(1), 1, 1e-15); %! assert (sumders(2:end), zeros(rows(squeeze(ders(ii,:,:)))-1, 1), 1e-13); -%! endfor +%! end %! assert (ders(:, (p+2):end, :), zeros(numel(u), 8-p-1, p+1), 1e-13) %! assert (all(all(ders(:, 1, :) <= 1)), true)
--- a/extra/nurbs/inst/bspdegelev.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/bspdegelev.m Mon Mar 22 18:13:09 2010 +0000 @@ -1,18 +1,3 @@ -%% Copyright (C) 2003 Mark Spink, 2007 Daniel Claxton -%% -%% 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, see <http://www.gnu.org/licenses/>. - function [ic,ik] = bspdegelev(d,c,k,t) % BSPDEGELEV: Degree elevate a univariate B-Spline. @@ -21,16 +6,32 @@ % % [ic,ik] = bspdegelev(d,c,k,t) % -% Parameters: +% INPUT: % % d - Degree of the B-Spline. % c - Control points, matrix of size (dim,nc). % k - Knot sequence, row vector of size nk. % t - Raise the B-Spline degree t times. % +% OUTPUT: +% % ic - Control points of the new B-Spline. % ik - Knot vector of the new B-Spline. % +% Copyright (C) 2000 Mark Spink, 2007 Daniel Claxton +% +% 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, see <http://www.gnu.org/licenses/>. [mc,nc] = size(c); % @@ -269,7 +270,7 @@ % mxFree(alfs); % % return(ierr); - % } +end % } function b = bincoeff(n,k) @@ -286,9 +287,10 @@ % double bincoeff(int n, int k) % { b = floor(0.5+exp(factln(n)-factln(k)-factln(n-k))); % return floor(0.5+exp(factln(n)-factln(k)-factln(n-k))); - % } +end % } function f = factln(n) % computes ln(n!) if n <= 1, f = 0; return, end -f = gammaln(n+1); %log(factorial(n)); \ No newline at end of file +f = gammaln(n+1); %log(factorial(n)); +end \ No newline at end of file
--- a/extra/nurbs/inst/bspderiv.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/bspderiv.m Mon Mar 22 18:13:09 2010 +0000 @@ -1,18 +1,3 @@ -%% Copyright (C) 2003 Mark Spink, 2007 Daniel Claxton -%% -%% 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, see <http://www.gnu.org/licenses/>. - function [dc,dk] = bspderiv(d,c,k) % BSPDERIV: B-Spline derivative. @@ -33,6 +18,21 @@ % dk - knot sequence of the derivative double vector(nk) % % Modified version of Algorithm A3.3 from 'The NURBS BOOK' pg98. +% +% Copyright (C) 2000 Mark Spink, 2007 Daniel Claxton +% +% 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, see <http://www.gnu.org/licenses/>. [mc,nc] = size(c); nk = numel(k); @@ -64,4 +64,4 @@ % freevec2mat(ctrl); % % return ierr; - % } +end % }
--- a/extra/nurbs/inst/bspeval.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/bspeval.m Mon Mar 22 18:13:09 2010 +0000 @@ -1,18 +1,3 @@ -%% Copyright (C) 2003 Mark Spink, 2007 Daniel Claxton -%% -%% 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, see <http://www.gnu.org/licenses/>. - function p = bspeval(d,c,k,u) % BSPEVAL: Evaluate B-Spline at parametric points. @@ -32,6 +17,20 @@ % % p - Evaluated points, matrix of size (dim,nu) % +% Copyright (C) 2000 Mark Spink, 2007 Daniel Claxton +% +% 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, see <http://www.gnu.org/licenses/>. nu = numel(u); [mc,nc] = size(c); @@ -70,4 +69,4 @@ % freevec2mat(ctrl); % % return ierr; - % } +end % }
--- a/extra/nurbs/inst/bspkntins.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/bspkntins.m Mon Mar 22 18:13:09 2010 +0000 @@ -1,18 +1,3 @@ -%% Copyright (C) 2003 Mark Spink, 2007 Daniel Claxton -%% -%% 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, see <http://www.gnu.org/licenses/>. - function [ic,ik] = bspkntins(d,c,k,u) % BSPKNTINS: Insert knots into a B-Spline @@ -35,8 +20,23 @@ % % Modified version of Algorithm A5.4 from 'The NURBS BOOK' pg164. % +% Copyright (C) 2000 Mark Spink, 2007 Daniel Claxton, 2010 Rafael Vazquez +% +% 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, see <http://www.gnu.org/licenses/>. [mc,nc] = size(c); +u = sort(u); nu = numel(u); nk = numel(k); % @@ -108,5 +108,5 @@ % freevec2mat(ictrl); % % return ierr; - % } +end % }
--- a/extra/nurbs/inst/curvederivcpts.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/curvederivcpts.m Mon Mar 22 18:13:09 2010 +0000 @@ -1,55 +1,54 @@ -%% Copyright (C) 2009 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 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, see <http://www.gnu.org/licenses/>. +function pk = curvederivcpts (n, p, U, P, d, r1, r2) +% +% usage: pk = curvederivcpts (n, p, U, P, d) +% +% INPUT: +% n+1 = number of control points +% p = degree of the spline +% d = maximum derivative order (d<=p) +% U = knots +% P = control points +% OUTPUT: +% pk(k,i) = i-th control point (k-1)-th derivative +% +% Adaptation of algorithm A3.3 from the NURBS book, pg98. +% +% Copyright (C) 2009 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 2 of the License, or +% (at your option) any later version. -%% usage: pk = curvederivcpts (n, p, U, P, d) -%% -%% INPUT: -%% n+1 = number of control points -%% p = degree of the spline -%% d = maximum derivative order (d<=p) -%% U = knots -%% P = control points -%% OUTPUT: -%% pk(k,i) = i-th control point (k-1)-th derivative -%% -%% Adaptation of algorithm A3.3 from the NURBS book - -function pk = curvederivcpts (n, p, U, P, d, r1, r2) +% 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, see <http://www.gnu.org/licenses/>. if (nargin <= 5) r1 = 0; r2 = n; - endif + end - r = r2 - r1, d + r = r2 - r1; for i=0:r pk(1, i+1) = P(r1+i+1); - endfor + end for k=1:d tmp = p - k + 1; for i=0:r-k pk (k+1, i+1) = tmp * (pk(k,i+2)-pk(k,i+1)) / ... (U(r1+i+p+2)-U(r1+i+k+1)); - endfor - endfor - size(pk) -endfunction + end + end +end %!test %! line = nrbmak([0.0 1.5; 0.0 3.0],[0.0 0.0 1.0 1.0]); -%! pk = curvederivcpts (line.number-1, line.order-1, line.knots, +%! pk = curvederivcpts (line.number-1, line.order-1, line.knots,... %! line.coefs(1,:), 2); %! assert (pk, [0 3/2; 3/2 0], 100*eps);
--- a/extra/nurbs/inst/deg2rad.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/deg2rad.m Mon Mar 22 18:13:09 2010 +0000 @@ -6,9 +6,11 @@ % % rad = deg2rad(deg); % -% Parameters: +% INPUT: % % deg : Angle in degrees. +% +% OUTPUT: % % rad : Angle in radians. % @@ -21,9 +23,22 @@ % % // Convert 35 degrees to radians % rad = deg2rad(35); +% +% Copyright (C) 2000 Mark Spink +% +% 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. -% D.M. Spink -% Copyright (c) 2000. +% 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, see <http://www.gnu.org/licenses/>. rad = pi*deg/180.0; +end \ No newline at end of file
--- a/extra/nurbs/inst/demo4surf.m Sun Mar 21 18:40:41 2010 +0000 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,11 +0,0 @@ -function demo4surf -% Demonstration of a bilinear surface. -% - -% D.M. Spink -% Copyright (c) 2000 - -srf = nrb4surf([0.0 0.0 0.5],[1.0 0.0 -0.5],[0.0 1.0 -0.5],[1.0 1.0 0.5]); -nrbplot(srf,[10,10]); -title('Construction of a bilinear surface.'); -
--- a/extra/nurbs/inst/democirc.m Sun Mar 21 18:40:41 2010 +0000 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,16 +0,0 @@ -function democirc -% Demonstration of a circle and arcs in the x-y plane. -% - -% D.M. Spink -% Copyright (c) 2000 - -for r = 1:9 - crv = nrbcirc(r,[],deg2rad(45),deg2rad(315)); - nrbplot(crv,50); - hold on; -end -hold off; -axis equal; -title('NURBS construction of a 2D circle and arc.'); -
--- a/extra/nurbs/inst/democoons.m Sun Mar 21 18:40:41 2010 +0000 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,35 +0,0 @@ -% Construction of a bilinearly blended Coons surface -% - -% D.M. Spink -% Copyright (c) 2000. - -% Boundary curve 1 -pnts = [ 0.0 3.0 4.5 6.5 8.0 10.0; - 0.0 0.0 0.0 0.0 0.0 0.0; - 2.0 2.0 7.0 4.0 7.0 9.0]; -crv1 = nrbmak(pnts, [0 0 0 1/3 0.5 2/3 1 1 1]); - -% Boundary curve 2 -pnts= [ 0.0 3.0 5.0 8.0 10.0; - 10.0 10.0 10.0 10.0 10.0; - 3.0 5.0 8.0 6.0 10.0]; -crv2 = nrbmak(pnts, [0 0 0 1/3 2/3 1 1 1]); - -% Boundary curve 3 -pnts= [ 0.0 0.0 0.0 0.0; - 0.0 3.0 8.0 10.0; - 2.0 0.0 5.0 3.0]; -crv3 = nrbmak(pnts, [0 0 0 0.5 1 1 1]); - -% Boundary curve 4 -pnts= [ 10.0 10.0 10.0 10.0 10.0; - 0.0 3.0 5.0 8.0 10.0; - 9.0 7.0 7.0 10.0 10.0]; -crv4 = nrbmak(pnts, [0 0 0 0.25 0.75 1 1 1]); - -srf = nrbcoons(crv1, crv2, crv3, crv4); - -% Draw the surface -nrbplot(srf,[20 20]); -title('Construction of a bilinearly blended Coons surface.');
--- a/extra/nurbs/inst/democurve.m Sun Mar 21 18:40:41 2010 +0000 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,44 +0,0 @@ -%% Copyright (C) 2003 Mark Spink -%% -%% 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, see <http://www.gnu.org/licenses/>. - -function democurve - -% Shows a simple test curve. -% - - -crv = nrbtestcrv; - -% plot the control points -plot(crv.coefs(1,:),crv.coefs(2,:),'ro'); -title('Arbitrary Test 2D Curve.'); -hold on; -plot(crv.coefs(1,:),crv.coefs(2,:),'r--'); - -% plot the nurbs curve -nrbplot(crv,48); -hold off; - -crv.knots(4)=0.1; -figure -% plot the control points -plot(crv.coefs(1,:),crv.coefs(2,:),'ro'); -title('Arbitrary Test 2D Curve.'); -hold on; -plot(crv.coefs(1,:),crv.coefs(2,:),'r--'); - -% plot the nurbs curve -nrbplot(crv,48); -hold off;
--- a/extra/nurbs/inst/democylind.m Sun Mar 21 18:40:41 2010 +0000 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,12 +0,0 @@ -function democylind -% Demonstration of the construction of a cylinder. -% - -% D.M. Spink -% Copyright (c) 2000 - -srf = nrbcylind(3,1,[],deg2rad(270),deg2rad(180)); -nrbplot(srf,[20,20]); -axis equal; -title('Cylinderical section by extrusion of a circular arc.'); -
--- a/extra/nurbs/inst/demodegelev.m Sun Mar 21 18:40:41 2010 +0000 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,22 +0,0 @@ -% Demonstration of the degree elevation algorithm. -% - -crv = nrbtestcrv; - -% plot the control points -plot(crv.coefs(1,:),crv.coefs(2,:),'bo') -title('Degree elevation of test curve by 1'); -hold on; -plot(crv.coefs(1,:),crv.coefs(2,:),'b--'); - -% draw nurbs curve -nrbplot(crv,48); - -% degree elevate the curve by 1 -icrv = nrbdegelev(crv, 1); - -% insert new knots and plot new control points -plot(icrv.coefs(1,:),icrv.coefs(2,:),'ro') -plot(icrv.coefs(1,:),icrv.coefs(2,:),'r--'); - -hold off;
--- a/extra/nurbs/inst/demodercrv.m Sun Mar 21 18:40:41 2010 +0000 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,48 +0,0 @@ -%% Copyright (C) 2003 Mark Spink -%% -%% 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, see <http://www.gnu.org/licenses/>. - -function demodercrv - -% Demonstrates the construction of a general -% curve and determine of the derivative. -% - - - - -% make and draw nurbs test curve -crv = nrbtestcrv; -nrbplot(crv,48); -title('First derivatives along a test curve.'); - -npts = 9; -tt = linspace(0.0,1.0,npts); - -dcrv = nrbderiv(crv); - -% first derivative -[p1, dp] = nrbdeval(crv,dcrv,tt); - -p2 = vecnorm(dp); - -hold on; -plot(p1(1,:),p1(2,:),'ro'); -h = quiver(p1(1,:),p1(2,:),p2(1,:),p2(2,:),0); -set(h,'Color','black'); -hold off; - - - -
--- a/extra/nurbs/inst/demodersrf.m Sun Mar 21 18:40:41 2010 +0000 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,34 +0,0 @@ -function demodersrf -% Demonstrates the construction of a general -% surface derivatives. -% - -% D.M. Spink -% Copyright (c) 2000 - -% make and draw a test surface -srf = nrbtestsrf; -p = nrbeval(srf,{linspace(0.0,1.0,20) linspace(0.0,1.0,20)}); -h = surf(squeeze(p(1,:,:)),squeeze(p(2,:,:)),squeeze(p(3,:,:))); -set(h,'FaceColor','blue','EdgeColor','blue'); -title('First derivatives over a test surface.'); - -npts = 5; -tt = linspace(0.0,1.0,npts); - -dsrf = nrbderiv(srf); - -[p1, dp] = nrbdeval(srf, dsrf, {tt, tt}); - -up2 = vecnorm(dp{1}); -vp2 = vecnorm(dp{2}); - -hold on; -plot3(p1(1,:),p1(2,:),p1(3,:),'ro'); -h1 = quiver3(p1(1,:),p1(2,:),p1(3,:),up2(1,:),up2(2,:),up2(3,:)); -h2 = quiver3(p1(1,:),p1(2,:),p1(3,:),vp2(1,:),vp2(2,:),vp2(3,:)); -set(h1,'Color','black'); -set(h2,'Color','black'); - -hold off; -
--- a/extra/nurbs/inst/demoellip.m Sun Mar 21 18:40:41 2010 +0000 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,14 +0,0 @@ -function demoellip -% Demonstration of a unit circle transformed to a inclined ellipse -% by first scaling, then rotating and finally translating. -% - -% D.M. Spink -% Copyright (c) 2000 - -xx = vectrans([2.0 1.0])*vecroty(pi/8)*vecrotx(pi/4)*vecscale([1.0 2.0]); -c0 = nrbtform(nrbcirc, xx); -nrbplot(c0,50); -title('Construction of an ellipse by transforming a unit circle.'); -grid on; -
--- a/extra/nurbs/inst/demoextrude.m Sun Mar 21 18:40:41 2010 +0000 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,6 +0,0 @@ -% Demonstration of surface construction by extrusion. -% - -srf = nrbextrude(nrbtestcrv,[0 0 5]); -nrbplot(srf,[40 10]); -title('Extrusion of a test curve along the z-axis');
--- a/extra/nurbs/inst/demogeom.m Sun Mar 21 18:40:41 2010 +0000 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,32 +0,0 @@ -function demogeom -% Demonstration of how to construct a 2D geometric -% shape from a piece-wise set of NURBSs -% - -% D.M. Spink -% Copyright (c) 2000. - -coefs = [0.0 7.5 15.0 25.0 35.0 30.0 27.5 30.0; - 0.0 2.5 0.0 -5.0 5.0 15.0 22.5 30.0]; -knots = [0.0 0.0 0.0 1/6 1/3 1/2 2/3 5/6 1.0 1.0 1.0]; - -% Geometric boundary curve -geom = [ -nrbmak(coefs,knots) -nrbline([30.0 30.0],[20.0 30.0]) -nrbline([20.0 30.0],[20.0 20.0]) -nrbcirc(10.0,[10.0 20.0],1.5*pi,0.0) -nrbline([10.0 10.0],[0.0 10.0]) -nrbline([0.0 10.0],[0.0 0.0]) -nrbcirc(5.0,[22.5 7.5]) -]; - -ng = length(geom); -for i = 1:ng - nrbplot(geom(i),50); - hold on; -end -hold off; -axis equal; -title('2D Geometry formed by a series of NURBS curves'); -
--- a/extra/nurbs/inst/demohelix.m Sun Mar 21 18:40:41 2010 +0000 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,26 +0,0 @@ -function demohelix -% Demonstration of a 3D helical curve -% - -% D.M. Spink -% Copyright (c) 2000 - -coefs =[ 6.0 0.0 6.0 1; - -5.5 0.5 5.5 1; - -5.0 1.0 -5.0 1; - 4.5 1.5 -4.5 1; - 4.0 2.0 4.0 1; - -3.5 2.5 3.5 1; - -3.0 3.0 -3.0 1; - 2.5 3.5 -2.5 1; - 2.0 4.0 2.0 1; - -1.5 4.5 1.5 1; - -1.0 5.0 -1.0 1; - 0.5 5.5 -0.5 1; - 0.0 6.0 0.0 1]'; -knots = [0 0 0 0 .1 .2 .3 .4 .5 .6 .7 .8 .9 1 1 1 1]; - -crv = nrbmak(coefs,knots); -nrbplot(crv,100); -title('3D helical curve.'); -grid on;
--- a/extra/nurbs/inst/demokntins.m Sun Mar 21 18:40:41 2010 +0000 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,27 +0,0 @@ -function demokntins -% Demonstration of the knot insertion algorithm. -% - -% D.M. Spink -% Copyright (c) 2000. - -crv = nrbtestcrv; - -% plot the control points -plot(crv.coefs(1,:),crv.coefs(2,:),'bo') -title('Knot insertion along test curve.'); -hold on; -plot(crv.coefs(1,:),crv.coefs(2,:),'b--'); - -% draw nurbs curve -nrbplot(crv,48); - -% insert new knots and plot new control points -icrv = nrbkntins(crv,[0.125 0.375 0.625 0.875] ); -plot(icrv.coefs(1,:),icrv.coefs(2,:),'ro') -plot(icrv.coefs(1,:),icrv.coefs(2,:),'r--'); - -hold off; - - -
--- a/extra/nurbs/inst/demoline.m Sun Mar 21 18:40:41 2010 +0000 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,11 +0,0 @@ -function demoline -% Demonstration of a 3D straight line -% - -% D.M. Spink -% Copyright (c) 2000. - -crv = nrbline([0.0 0.0 0.0]',[5.0 4.0 2.0]'); -nrbplot(crv,1); -title('3D straight line.'); -grid on;
--- a/extra/nurbs/inst/demorect.m Sun Mar 21 18:40:41 2010 +0000 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,11 +0,0 @@ -function demorect -% Demonstrate of rectangluar curve -% - -% D.M. Spink -% Copyright (c) 2000 - -crv = nrbtform(nrbrect(2,1), vecrotz(deg2rad(35))); -nrbplot(crv,4); -title('Construction and rotation of a rectangle.'); -axis equal;
--- a/extra/nurbs/inst/demorevolve.m Sun Mar 21 18:40:41 2010 +0000 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,31 +0,0 @@ -function demorevolve -% Demonstration of surface construction by revolving a -% profile curve. - -% D.M. Spink -% Copyright (c) 2000 - -% Construct a test profile in the x-z plane -pnts = [3.0 5.5 5.5 1.5 1.5 4.0 4.5; - 0.0 0.0 0.0 0.0 0.0 0.0 0.0; - 0.5 1.5 4.5 3.0 7.5 6.0 8.5]; -crv = nrbmak(pnts,[0 0 0 1/4 1/2 3/4 3/4 1 1 1]); - -% rotate and vectrans by some arbitrary amounts. -xx = vecrotz(deg2rad(25))*vecroty(deg2rad(15))*vecrotx(deg2rad(20)); -nrb = nrbtform(crv,vectrans([5 5])*xx); - -% define axes of rotation -pnt = [5 5 0]'; -vec = xx*[0 0 1 1]'; -srf = nrbrevolve(nrb,pnt,vec(1:3)); - -% make and draw nurbs curve -p = nrbeval(srf,{linspace(0.0,1.0,20) linspace(0.0,1.0,20)}); -surfl(squeeze(p(1,:,:)),squeeze(p(2,:,:)),squeeze(p(3,:,:))); -title('Construct of a 3D surface by revolution of a curve.'); -shading interp; -colormap(copper); -axis equal; - -
--- a/extra/nurbs/inst/demoruled.m Sun Mar 21 18:40:41 2010 +0000 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,14 +0,0 @@ -% Demonstration of ruled surface construction. -% - -% D.M. Spink -% Copyright (c) 2000. - -% clear all recorded graphics and the current window -title('Ruled surface construction from two NURBS curves.'); - -crv1 = nrbtestcrv; -crv2 = nrbtform(nrbcirc(4,[4.5;0],pi,0.0),vectrans([0.0 4.0 -4.0])); -srf = nrbruled(crv1,crv2); -nrbplot(srf,[40 20]); -title('Ruled surface construction from two NURBS curves.');
--- a/extra/nurbs/inst/demosurf.m Sun Mar 21 18:40:41 2010 +0000 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,10 +0,0 @@ -function demosurf -% DEMOSURF: Shows a simple test surface -% - -% D.M. Spink -% Copyright (c) 2000 - -nrbplot(nrbtestsrf,[20 20]); -title('Arbitrary test surface.'); -
--- a/extra/nurbs/inst/demotorus.m Sun Mar 21 18:40:41 2010 +0000 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,33 +0,0 @@ -%% Copyright (C) 2003 Mark Spink -%% -%% 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, see <http://www.gnu.org/licenses/>. - -function demotorus - -% DEMOTORUS: A second demonstration of surface construction -% by revolution. -% - -% D.M. Spink - - -sphere = nrbrevolve(nrbcirc(1,[],0.0,pi),[0.0 0.0 0.0],[1.0 0.0 0.0]); -nrbplot(sphere,[20 20],'light','on'); -title('Ball and torus - surface construction by revolution'); -hold on; -torus = nrbrevolve(nrbcirc(0.2,[0.9 1.0]),[0.0 0.0 0.0],[1.0 0.0 0.0]); -nrbplot(torus,[20 10],'light','on'); -nrbplot(nrbtform(torus,vectrans([-1.8])),[20 10],'light','on'); -hold off; -
--- a/extra/nurbs/inst/findspan.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/findspan.m Mon Mar 22 18:13:09 2010 +0000 @@ -1,56 +1,64 @@ -%% Copyright (C) 2009 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 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, see <http://www.gnu.org/licenses/>. - -function sv = findspan (n, p, uv, U) - -% FINDSPAN: Find the span of a B-Spline knot vector at a parametric point -% -% Calling Sequence: -% -% s = findspan(n,p,u,U) -% -% INPUT: -% -% n - number of control points - 1 -% p - spline degree -% u - parametric point -% U - knot sequence -% -% U(1) <= u <= U(end) -% RETURN: -% -% s - knot span -% - - -% This should be equivalent to -% Algorithm A2.1 from 'The NURBS BOOK' pg68 - -sv = min (lookup (U, uv, "lr") - 1, n); - -%!test -%! n = 3; -%! U = [0 0 0 1/2 1 1 1]; -%! p = 2; -%! u = linspace(0, 1, 10); -%! s = findspan (n, p, u, U); -%! assert (s, [2*ones(1, 5) 3*ones(1, 5)]); - -%!test -%! p = 2; m = 7; n = m - p - 1; -%! U = [zeros(1,p) linspace(0,1,m+1-2*p) ones(1,p)]; -%! u = [ 0 0.11880 0.55118 0.93141 0.40068 0.35492 0.44392 0.88360 0.35414 0.92186 0.83085 1]; -%! s = [2 2 3 4 3 3 3 4 3 4 4 4]; -%! assert (findspan (n, p, u, U), s, 1e-10); \ No newline at end of file +function s = findspan(n,p,u,U) +% FINDSPAN Find the span of a B-Spline knot vector at a parametric point +% ------------------------------------------------------------------------- +% ADAPTATION of FINDSPAN from C +% ------------------------------------------------------------------------- +% +% Calling Sequence: +% +% s = findspan(n,p,u,U) +% +% INPUT: +% +% n - number of control points - 1 +% p - spline degree +% u - parametric point +% U - knot sequence +% +% OUTPUT: +% +% s - knot span +% +% Modification of Algorithm A2.1 from 'The NURBS BOOK' pg68 +% +% Copyright (C) 2010 Rafael Vazquez +% +% 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, see <http://www.gnu.org/licenses/>. + +if max(u) > U(end) || min(u) < U(1) + error('Some value is outside the knot span') +end + +s = zeros(size(u)); +for j = 1:numel(u) + if (u(j)==U(n+2)), s(j)=n; continue, end + s(j) = find(u(j) >= U,1,'last')-1; +end + +end + +%!test +%! n = 3; +%! U = [0 0 0 1/2 1 1 1]; +%! p = 2; +%! u = linspace(0, 1, 10); +%! s = findspan (n, p, u, U); +%! assert (s, [2*ones(1, 5) 3*ones(1, 5)]); + +%!test +%! p = 2; m = 7; n = m - p - 1; +%! U = [zeros(1,p) linspace(0,1,m+1-2*p) ones(1,p)]; +%! u = [ 0 0.11880 0.55118 0.93141 0.40068 0.35492 0.44392 0.88360 0.35414 0.92186 0.83085 1]; +%! s = [2 2 3 4 3 3 3 4 3 4 4 4]; +%! assert (findspan (n, p, u, U), s, 1e-10);
--- a/extra/nurbs/inst/findspanc.m Sun Mar 21 18:40:41 2010 +0000 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,53 +0,0 @@ -%% Copyright (C) 2009 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 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, see <http://www.gnu.org/licenses/>. - -function sv = findspanc (uv, U) - -% FINDSPANC: Find the span of a B-Spline knot vector at a parametric point -% -% Calling Sequence: -% -% s = findspanc(u, U) -% -% INPUT: -% -% u - parametric point -% U - knot sequence -% -% U(1) <= u <= U(end) -% RETURN: -% -% s - knot span -% -% NOTE: the difference between findspan and findspanc is that -% the latter works for non-open knot vectors - - -sv = lookup (U, uv, "lr") - 1; - -%!test -%! n = 3; -%! U = [0 0 0 1/2 1 1 1]; -%! p = 2; -%! u = linspace(0, 1, 10)(2:end-1); -%! s = findspanc (u, U); -%! assert (s, [2*ones(1, 4) 3*ones(1, 4)]); - -%!test -%! p = 2; m = 7; n = m - p - 1; -%! U = [zeros(1,p) linspace(0,1,m+1-2*p) ones(1,p)]; -%! u = [ 0.11880 0.55118 0.93141 0.40068 0.35492 0.44392 0.88360 0.35414 0.92186 0.83085 ]; -%! s = [2 3 4 3 3 3 4 3 4 4]; -%! assert (findspanc (u, U), s, 1e-10); \ No newline at end of file
--- a/extra/nurbs/inst/nrb4surf.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/nrb4surf.m Mon Mar 22 18:13:09 2010 +0000 @@ -6,7 +6,7 @@ % % srf = nrb4surf(p11,p12,p21,p22) % -% Parameters: +% INPUT: % % p11 : Cartesian coordinate of the lhs bottom corner point. % @@ -15,6 +15,8 @@ % p21 : Cartesian coordinate of the lhs top corner point. % % p22 : Cartesian coordinate of the rhs top corner point. +% +% OUTPUT: % % srf : NURBS bilinear surface, see nrbmak. % @@ -34,11 +36,23 @@ % |p11 p12| % -------------------> U direction % +% +% Copyright (C) 2000 Mark Spink +% +% 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. -% D.M. Spink -% Copyright (c) 2000. +% 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, see <http://www.gnu.org/licenses/>. -if nargin < 4 +if nargin ~= 4 error('Four corner points must be defined'); end @@ -50,4 +64,11 @@ knots = {[0 0 1 1] [0 0 1 1]}; srf = nrbmak(coefs, knots); - + +end + +%!demo +%! srf = nrb4surf([0.0 0.0 0.5],[1.0 0.0 -0.5],[0.0 1.0 -0.5],[1.0 1.0 0.5]); +%! nrbplot(srf,[10,10]); +%! title('Construction of a bilinear surface.'); +%! hold off \ No newline at end of file
--- a/extra/nurbs/inst/nrb_srf_basisfun__.m Sun Mar 21 18:40:41 2010 +0000 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,51 +0,0 @@ -%% Copyright (C) 2009 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 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, see <http://www.gnu.org/licenses/>. - -function [B, N] = nrb_srf_basisfun__ (points, nrb); - - %% __NRB_SRF_BASISFUN__: Undocumented internal function - - m = size (nrb.coefs, 2) -1; - n = size (nrb.coefs, 3) -1; - - p = nrb.order(1) -1; - q = nrb.order(2) -1; - - u = points(1,:); - v = points(2,:); - npt = length(u); - - U = nrb.knots{1}; - V = nrb.knots{2}; - - w = squeeze(nrb.coefs(4,:,:)); - - spu = findspan (m, p, u, U); - spv = findspan (n, q, v, V); - NuIkuk = basisfun (spu, u, p, U); - NvJkvk = basisfun (spv, v, q, V); - - indIkJk = nrbnumbasisfun (points, nrb); - - for k=1:npt - wIkaJkb(1:p+1, 1:q+1) = reshape (w(indIkJk(k, :)), p+1, q+1); - NuIkukaNvJkvk(1:p+1, 1:q+1) = (NuIkuk(k, :).' * NvJkvk(k, :)); - RIkJk(k, :) = reshape((NuIkukaNvJkvk .* wIkaJkb ./ sum(sum(NuIkukaNvJkvk .* wIkaJkb))),1,[]); - end - - B = RIkJk; - N = indIkJk; - - end \ No newline at end of file
--- a/extra/nurbs/inst/nrb_srf_basisfun_der__.m Sun Mar 21 18:40:41 2010 +0000 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,62 +0,0 @@ -%% Copyright (C) 2009 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 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, see <http://www.gnu.org/licenses/>. - -function [Bu, Bv, N] = nrb_srf_basisfun_der__ (points, nrb); - - %% __NRB_SRF_BASISFUN_DER__: Undocumented internal function - - m = size (nrb.coefs, 2) -1; - n = size (nrb.coefs, 3) -1; - - p = nrb.order(1) -1; - q = nrb.order(2) -1; - - u = points(1,:); - v = points(2,:); - npt = length(u); - - U = nrb.knots{1}; - V = nrb.knots{2}; - - w = squeeze(nrb.coefs(4,:,:)); - - spu = findspan (m, p, u, U); - spv = findspan (n, q, v, V); - N = nrbnumbasisfun (points, nrb); - - NuIkuk = basisfun (spu, u, p, U); - NvJkvk = basisfun (spv, v, q, V); - - NuIkukprime = basisfunder (spu, p, u, U, 1); - NuIkukprime = reshape (NuIkukprime(:,2,:), npt, []); - - NvJkvkprime = basisfunder (spv, q, v, V, 1); - NvJkvkprime = reshape (NvJkvkprime(:,2,:), npt, []); - - for k=1:npt - wIkaJkb(1:p+1, 1:q+1) = reshape (w(N(k, :)), p+1, q+1); - - Num = (NuIkuk(k, :).' * NvJkvk(k, :)) .* wIkaJkb; - Num_du = (NuIkukprime(k, :).' * NvJkvk(k, :)) .* wIkaJkb; - Num_dv = (NuIkuk(k, :).' * NvJkvkprime(k, :)) .* wIkaJkb; - Denom = sum(sum(Num)); - Denom_du = sum(sum(Num_du)); - Denom_dv = sum(sum(Num_dv)); - - Bu(k, :) = reshape((Num_du/Denom - Denom_du.*Num/Denom.^2),1,[]); - Bv(k, :) = reshape((Num_dv/Denom - Denom_dv.*Num/Denom.^2),1,[]); - end - -end \ No newline at end of file
--- a/extra/nurbs/inst/nrbbasisfun.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/nrbbasisfun.m Mon Mar 22 18:13:09 2010 +0000 @@ -1,26 +1,11 @@ -%% Copyright (C) 2009 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 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, see <http://www.gnu.org/licenses/>. - function [B, id] = nrbbasisfun (points, nrb) % NRBBASISFUN: Basis functions for NURBS % % Calling Sequence: % -% B = nrbbasisfun (u, crv) -% B = nrbbasisfun ({u, v}, srf) +% B = nrbbasisfun (u, crv) +% B = nrbbasisfun ({u, v}, srf) % [B, N] = nrbbasisfun ({u, v}, srf) % [B, N] = nrbbasisfun (p, srf) % @@ -33,13 +18,28 @@ % % OUTPUT: % -% B - Basis functions +% B - Value of the basis functions at the points % size(B)=[numel(u),(p+1)] for curves % or [numel(u)*numel(v), (p+1)*(q+1)] for surfaces % % N - Indices of the basis functions that are nonvanishing at each % point. size(N) == size(B) % +% +% Copyright (C) 2009 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 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, see <http://www.gnu.org/licenses/>. if ( (nargin<2) ... || (nargout>2) ... @@ -48,15 +48,14 @@ || (~iscell(points) && iscell(nrb.knots) && (size(points,1)~=2)) ... || (~iscell(nrb.knots) && (nargout>1)) ... ) - print_usage(); + error('Incorrect input arguments in nrbbasisfun'); end - if (~iscell(nrb.knots)) %% NURBS curve + if (~iscell(nrb.knots)) %% NURBS curve [B, id] = nrb_crv_basisfun__ (points, nrb); - else %% NURBS surface - + elseif size(nrb.knots,2) == 2 %% NURBS surface if (iscell(points)) [v, u] = meshgrid(points{2}, points{1}); p = [u(:), v(:)]'; @@ -66,24 +65,23 @@ [B, id] = nrb_srf_basisfun__ (p, nrb); + else %% NURBS volume + error('The function nrbbasisfun is not yet ready for volumes') end -end +end - function [B, nbfu] = nrb_crv_basisfun__ (points, nrb); - n = size (nrb.coefs, 2) -1; - p = nrb.order -1; - u = points; - U = nrb.knots; - w = nrb.coefs(4,:); - - spu = findspan (n, p, u, U); - nbfu = numbasisfun (spu, u, p, U); - - N = w(nbfu+1) .* basisfun (spu, u, p, U); - B = bsxfun (@(x,y) x./y, N, sum (N,2)); - - end - +%!demo +%! U = [0 0 0 0 1 1 1 1]; +%! x = [0 1/3 2/3 1] ; +%! y = [0 0 0 0]; +%! w = [1 1 1 1]; +%! nrb = nrbmak ([x;y;y;w], U); +%! u = linspace(0, 1, 30); +%! B = nrbbasisfun (u, nrb); +%! xplot = sum(bsxfun(@(x,y) x.*y, B, x),2); +%! plot(xplot, B) +%! title('Cubic Bernstein polynomials') +%! hold off %!test %! U = [0 0 0 0 1 1 1 1]; @@ -97,30 +95,35 @@ %! %! yy = y; yy(1) = 1; %! nrb2 = nrbmak ([x.*w;yy;y;w], U); -%! #figure, plot(xplot, B(:,1), nrbeval(nrb2, u)(1,:).', w(1)*nrbeval(nrb2, u)(2,:).') -%! assert(B(:,1), w(1)*nrbeval(nrb2, u)(2,:).', 1e-6) +%! aux = nrbeval(nrb2,u); +%! %figure, plot(xplot, B(:,1), aux(1,:).', w(1)*aux(2,:).') +%! assert(B(:,1), w(1)*aux(2,:).', 1e-6) %! %! yy = y; yy(2) = 1; %! nrb2 = nrbmak ([x.*w;yy;y;w], U); -%! #figure, plot(xplot, B(:,2), nrbeval(nrb2, u)(1,:).', w(2)*nrbeval(nrb2, u)(2,:).') -%! assert(B(:,2), w(2)*nrbeval(nrb2, u)(2,:).', 1e-6) +%! aux = nrbeval(nrb2, u); +%! %figure, plot(xplot, B(:,2), aux(1,:).', w(2)*aux(2,:).') +%! assert(B(:,2), w(2)*aux(2,:).', 1e-6) %! %! yy = y; yy(3) = 1; %! nrb2 = nrbmak ([x.*w;yy;y;w], U); -%! #figure, plot(xplot, B(:,3), nrbeval(nrb2, u)(1,:).', w(3)*nrbeval(nrb2, u)(2,:).') -%! assert(B(:,3), w(3)*nrbeval(nrb2, u)(2,:).', 1e-6) +%! aux = nrbeval(nrb2,u); +%! %figure, plot(xplot, B(:,3), aux(1,:).', w(3)*aux(2,:).') +%! assert(B(:,3), w(3)*aux(2,:).', 1e-6) %! %! yy = y; yy(4) = 1; %! nrb2 = nrbmak ([x.*w;yy;y;w], U); -%! #figure, plot(xplot, B(:,4), nrbeval(nrb2, u)(1,:).', w(4)*nrbeval(nrb2, u)(2,:).') -%! assert(B(:,4), w(4)*nrbeval(nrb2, u)(2,:).', 1e-6) +%! aux = nrbeval(nrb2,u); +%! %figure, plot(xplot, B(:,4), aux(1,:).', w(4)*aux(2,:).') +%! assert(B(:,4), w(4)*aux(2,:).', 1e-6) %!test %! p = 2; q = 3; m = 4; n = 5; %! Lx = 1; Ly = 1; %! nrb = nrb4surf ([0 0], [1 0], [0 1], [1 1]); %! nrb = nrbdegelev (nrb, [p-1, q-1]); -%! nrb = nrbkntins (nrb, {linspace(0,1,m)(2:end-1), linspace(0,1,n)(2:end-1)}); +%! aux1 = linspace(0,1,m); aux2 = linspace(0,1,n); +%! nrb = nrbkntins (nrb, {aux1(2:end-1), aux2(2:end-1)}); %! u = rand (1, 30); v = rand (1, 10); %! [B, N] = nrbbasisfun ({u, v}, nrb); %! assert (sum(B, 2), ones(300, 1), 1e-6)
--- a/extra/nurbs/inst/nrbbasisfunder.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/nrbbasisfunder.m Mon Mar 22 18:13:09 2010 +0000 @@ -1,18 +1,3 @@ -%% Copyright (C) 2009 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 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, see <http://www.gnu.org/licenses/>. - function varargout = nrbbasisfunder (points, nrb) % NRBBASISFUNDER: NURBS basis functions derivatives @@ -45,6 +30,20 @@ % N - Indices of the basis functions that are nonvanishing at each % point. size(N) == size(B) % +% Copyright (C) 2009 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 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, see <http://www.gnu.org/licenses/>. if ( (nargin<2) ... @@ -54,14 +53,14 @@ || (~iscell(points) && iscell(nrb.knots) && (size(points,1)~=2)) ... || (~iscell(nrb.knots) && (nargout>2)) ... ) - print_usage(); + error('Incorrect input arguments in nrbbasisfun'); end - if (~iscell(nrb.knots)) %% NURBS curve + if (~iscell(nrb.knots)) %% NURBS curve [varargout{1}, varargout{2}] = nrb_crv_basisfun_der__ (points, nrb); - else %% NURBS surface + elseif size(nrb.knots,2) == 2 %% NURBS surface if (iscell(points)) [v, u] = meshgrid(points{2}, points{1}); @@ -72,39 +71,22 @@ [varargout{1}, varargout{2}, varargout{3}] = nrb_srf_basisfun_der__ (p, nrb); + else %% NURBS volume + error('The function nrbbasisfunder is not yet ready for volumes') end end - - function [Bu, nbfu] = nrb_crv_basisfun_der__ (points, nrb); - n = size (nrb.coefs, 2) -1; - p = nrb.order -1; - u = points; - U = nrb.knots; - w = nrb.coefs(4,:); - - spu = findspan (n, p, u, U); - nbfu = numbasisfun (spu, u, p, U); - - N = basisfun (spu, u, p, U); - - Nprime = basisfunder (spu, p, u, U, 1); - Nprime = squeeze(Nprime(:,2,:)); - - - [Dpc, Dpk] = bspderiv (p, w, U); - D = bspeval (p, w, U, u); - Dprime = bspeval (p-1, Dpc, Dpk, u); - - - Bu1 = bsxfun (@(np, d) np/d , Nprime.', D); - Bu2 = bsxfun (@(n, dp) n*dp, N.', Dprime./D.^2); - Bu = w(nbfu+1) .* (Bu1 - Bu2).'; - - end - - - +%!demo +%! U = [0 0 0 0 1 1 1 1]; +%! x = [0 1/3 2/3 1] ; +%! y = [0 0 0 0]; +%! w = [1 1 1 1]; +%! nrb = nrbmak ([x;y;y;w], U); +%! u = linspace(0, 1, 30); +%! [Bu, id] = nrbbasisfunder (u, nrb); +%! plot(u, Bu) +%! title('Derivatives of the cubic Bernstein polynomials') +%! hold off %!test %! U = [0 0 0 0 1 1 1 1]; @@ -132,9 +114,10 @@ %! Lx = 1; Ly = 1; %! nrb = nrb4surf ([0 0], [1 0], [0 1], [1 1]); %! nrb = nrbdegelev (nrb, [p-1, q-1]); -%! nrb = nrbkntins (nrb, {linspace(0,1,m)(2:end-1), linspace(0,1,n)(2:end-1)}); -%! nrb.coefs (4,:,:) += rand (size (nrb.coefs (4,:,:))); -%! tic(); [Bu, Bv, N] = nrbbasisfunder ({rand(1, 20), rand(1, 20)}, nrb); toc() +%! aux1 = linspace(0,1,m); aux2 = linspace(0,1,n); +%! nrb = nrbkntins (nrb, {aux1(2:end-1), aux2(2:end-1)}); +%! nrb.coefs (4,:,:) = nrb.coefs(4,:,:) + rand (size (nrb.coefs (4,:,:))); +%! [Bu, Bv, N] = nrbbasisfunder ({rand(1, 20), rand(1, 20)}, nrb); %! #plot3(squeeze(u(1,:,:)), squeeze(u(2,:,:)), reshape(Bu(:,10), 20, 20),'o') %! assert (sum (Bu, 2), zeros(20^2, 1), 1e-10)
--- a/extra/nurbs/inst/nrbbasisfungradient.m Sun Mar 21 18:40:41 2010 +0000 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,103 +0,0 @@ -%% Copyright (C) 2009 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 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., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA - -%% -*- texinfo -*- -%% @deftypefn {Function File} {[@var{dRdx}, @var{dRdy}]=} nrbbasisfungradient (@{@var{dzdu}, @var{dzdv}@}, @{@var{dxdu}, @var{dydu}, @var{dxdv}, @var{dydv}@}) -%% @deftypefnx {Function File} {[@var{dRdx}]=} nrbbasisfungradient (@var{dzdu}, @var{dxdu}, @var{nrb}) -%% NRBBASISFUNGRADIENT Compute the gradient of the basis functions of a NURBS surface at the -%% specified parametric points. -%% -%% INPUT: -%% @itemize @minus -%% @item @var{dzdu}, @var{dzdv}: basis function derivatives with respect -%% to parameters @code{u} and @code{v} -%% @item @var{dxdu}, @var{dydu}, @var{dxdv}, -%% @var{dydv}: NURBS geometry map derivatives -%% @var{nrb}: NURBS structure describing the surface -%% @end itemize -%% -%% OUTPUT: -%% @itemize @minus -%% @item @var{dRdx} derivative of the basis functions with respect to -%% the @code{x} direction -%% @item @var{dRdy} derivative of the basis functions with respect to -%% the @code{y} direction -%% @end itemize -%% -%% @seealso{nrbbasisfunder,nrbbasisfun,nrbderiv} -%% @end deftypefn - -%% Author: Carlo de Falco <carlo@guglielmo.local> -%% Created: 2009-04-29 - -function [varargout] = nrbbasisfungradient (dz, N, nrb) - - if (iscell(dz)) - - dzdu=dz{1}; - dzdv=dz{2}; - - X = squeeze(nrb.coefs(1,:,:)./nrb.coefs(4,:,:)); - Y = squeeze(nrb.coefs(2,:,:)./nrb.coefs(4,:,:)); - - nfun = size(N,2); - tmp = sum(X(N).*dzdu, 2); - dxdu= tmp(:, ones(nfun,1)); - tmp = sum(Y(N).*dzdu, 2); - dydu= tmp(:, ones(nfun,1)); - tmp = sum(X(N).*dzdv, 2); - dxdv= tmp(:, ones(nfun,1)); - tmp = sum(Y(N).*dzdv, 2); - dydv= tmp(:, ones(nfun,1)); - - detjac = dxdu.*dydv - dydu.*dxdv; - - varargout{1} = ( dydv .* dzdu - dydu .*dzdv)./detjac; - varargout{2} = (-dxdv .* dzdu + dxdu .*dzdv)./detjac; - %%keyboard - elseif (~iscell(dz) && (nargout==1)) - - varargout{1} = dzdu ./ dxdu; - - else - - print_usage(); - - end - -end - -%!test -%! p = 2; q = 3; m = 4; n = 5; -%! Lx = 1; Ly = 1; -%! nrb = nrb4surf ([0 0], [Lx 0], [0 Ly], [Lx Ly]); -%! nrb = nrbdegelev (nrb, [p-1, q-1]); -%! nrb = nrbkntins (nrb, {linspace(0,1,m)(2:end-1), linspace(0,1,n)(2:end-1)}); -%! nrb.coefs (4,:,:) += rand (size (nrb.coefs (4,:,:))); -%! -%! [u(2,:,:), u(1,:,:)] = meshgrid(rand (1, 20), rand (1, 20)); -%! -%! -%! uv (1,:) = u(1,:,:)(:)'; -%! uv (2,:) = u(2,:,:)(:)'; -%! -%! [dzdu, dzdv, connect] = nrbbasisfunder (uv, nrb); -%! nd = nrbderiv(nrb); -%! [ndp, dp] = nrbdeval(nrb, nd, uv); -%! -%! [dzdx, dzdy]= nrbbasisfungradient ({dzdu, dzdv}, connect, nrb); -%! assert(norm(sum(dzdx, 2), inf), 0, 1e-10) -%! assert(norm(sum(dzdy, 2), inf), 0, 1e-10) \ No newline at end of file
--- a/extra/nurbs/inst/nrbcirc.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/nrbcirc.m Mon Mar 22 18:13:09 2010 +0000 @@ -4,12 +4,12 @@ % % Calling Sequence: % -% crv = nrbarc() -% crv = nrbarc(radius) -% crv = nrbarc(radius,center) -% crv = nrbarc(radius,center,sang,eang) +% crv = nrbcirc() +% crv = nrbcirc(radius) +% crv = nrbcirc(radius,center) +% crv = nrbcirc(radius,center,sang,eang) % -% Parameters: +% INPUT: % % radius : Radius of the circle, default 1.0 % @@ -19,6 +19,8 @@ % % eang : End angle 360 degrees % +% OUTPUT: +% % crv : NURBS curve for a circular arc. % % Description: @@ -28,9 +30,21 @@ % constructed. % % Angles are defined as positive in the anti-clockwise direction. +% +% Copyright (C) 2000 Mark Spink +% +% 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. -% D.M. Spink -% Copyright (c) 2000. +% 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, see <http://www.gnu.org/licenses/>. if nargin < 1 radius = 1; @@ -95,3 +109,15 @@ end curve = nrbmak(coefs,knots); + +end + +%!demo +%! for r = 1:9 +%! crv = nrbcirc(r,[],deg2rad(45),deg2rad(315)); +%! nrbplot(crv,50); +%! hold on; +%! end +%! hold off; +%! axis equal; +%! title('NURBS construction of several 2D arcs.');
--- a/extra/nurbs/inst/nrbcoons.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/nrbcoons.m Mon Mar 22 18:13:09 2010 +0000 @@ -6,7 +6,7 @@ % % srf = nrbcoons(ucrv1, ucrv2, vcrv1, vcrv2) % -% Parameters: +% INPUT: % % ucrv1 : NURBS curve defining the bottom U direction boundary of % the constructed NURBS surface. @@ -19,6 +19,8 @@ % % vcrv2 : NURBS curve defining the top V direction boundary of % the constructed NURBS surface. +% +% OUTPUT: % % srf : Coons NURBS surface patch. % @@ -67,9 +69,21 @@ % % srf = nrbcoons(crv1, crv2, crv3, crv4); % nrbplot(srf,[20 20],220,45); +% +% Copyright (C) 2000 Mark Spink, 2010 Rafael Vazquez +% +% 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. -% D.M. Spink -% Copyright (c) 2000. +% 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, see <http://www.gnu.org/licenses/>. if nargin ~= 4 error('Incorrect number of input arguments'); @@ -77,7 +91,7 @@ r1 = nrbruled(u1, u2); r2 = nrbtransp(nrbruled(v1, v2)); -t = nrb4surf(u1.coefs(:,1), u1.coefs(:,end), u2.coefs(:,1), u2.coefs(:,end)); +t = nrb4surf(u1.coefs(:,1), u2.coefs(:,1), u1.coefs(:,end), u2.coefs(:,end)); % raise all surfaces to a common degree du = max([r1.order(1), r2.order(1), t.order(1)]); @@ -137,3 +151,31 @@ coefs(4,:,:) = r1.coefs(4,:,:) + r2.coefs(4,:,:) - t.coefs(4,:,:); srf = nrbmak(coefs, r1.knots); +end + +%!demo +%! pnts = [ 0.0 3.0 4.5 6.5 8.0 10.0; +%! 0.0 0.0 0.0 0.0 0.0 0.0; +%! 2.0 2.0 7.0 4.0 7.0 9.0]; +%! crv1 = nrbmak(pnts, [0 0 0 1/3 0.5 2/3 1 1 1]); +%! +%! pnts= [ 0.0 3.0 5.0 8.0 10.0; +%! 10.0 10.0 10.0 10.0 10.0; +%! 3.0 5.0 8.0 6.0 10.0]; +%! crv2 = nrbmak(pnts, [0 0 0 1/3 2/3 1 1 1]); +%! +%! pnts= [ 0.0 0.0 0.0 0.0; +%! 0.0 3.0 8.0 10.0; +%! 2.0 0.0 5.0 3.0]; +%! crv3 = nrbmak(pnts, [0 0 0 0.5 1 1 1]); +%! +%! pnts= [ 10.0 10.0 10.0 10.0 10.0; +%! 0.0 3.0 5.0 8.0 10.0; +%! 9.0 7.0 7.0 10.0 10.0]; +%! crv4 = nrbmak(pnts, [0 0 0 0.25 0.75 1 1 1]); +%! +%! srf = nrbcoons(crv1, crv2, crv3, crv4); +%! +%! nrbplot(srf,[20 20]); +%! title('Construction of a bilinearly blended Coons surface.'); +%! hold off \ No newline at end of file
--- a/extra/nurbs/inst/nrbcylind.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/nrbcylind.m Mon Mar 22 18:13:09 2010 +0000 @@ -10,7 +10,7 @@ % srf = nrbcylind(height,radius,center) % srf = nrbcylind(height,radius,center,sang,eang) % -% Parameters: +% INPUT: % % height : Height of the cylinder along the axis, default 1.0 % @@ -18,17 +18,32 @@ % % center : Center of the cylinder, default (0,0,0) % -% sang : Start angle relative to the origin, default 0. +% sang : Start angle relative to the origin, default 0. % % eang : End angle relative to the origin, default 2*pi. -% +% +% OUTPUT: +% +% srf : cylindrical surface patch % % Description: % % Construct a cylinder or cylindrical patch by extruding a circular arc. +% +% Copyright (C) 2000 Mark Spink +% +% 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. -% D.M. Spink -% Copyright (c) 2000. +% 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, see <http://www.gnu.org/licenses/>. if nargin < 1 height = 1; @@ -49,3 +64,11 @@ surf = nrbextrude(nrbcirc(radius,center,sang,eang),[0.0 0.0 height]); +end + +%!demo +%! srf = nrbcylind(3,1,[],deg2rad(270),deg2rad(180)); +%! nrbplot(srf,[20,20]); +%! axis equal; +%! title('Cylinderical section by extrusion of a circular arc.'); +%! hold off \ No newline at end of file
--- a/extra/nurbs/inst/nrbdegelev.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/nrbdegelev.m Mon Mar 22 18:13:09 2010 +0000 @@ -1,96 +1,173 @@ -function inurbs = nrbdegelev(nurbs, ntimes) -% -% NRBDEGELEV: Elevate the degree of the NURBS curve or surface. -% -% Calling Sequence: -% -% ecrv = nrbdegelev(crv,utimes); -% esrf = nrbdegelev(srf,{utimes,vtimes}); -% -% Parameters: -% -% crv : NURBS curve, see nrbmak. -% -% srf : NURBS surface, see nrbmak. -% -% utimes : Increase the degree along U direction utimes. -% -% vtimes : Increase the degree along V direction vtimes. -% -% ecrv : new NURBS structure for a curve with degree elevated. -% -% esrf : new NURBS structure for a surface with degree elevated. -% -% -% Description: -% -% Degree elevates the NURBS curve or surface. This function uses the -% B-Spline function bspdegelev, which interface to an internal 'C' -% routine. -% -% Examples: -% -% Increase the NURBS surface twice along the V direction. -% esrf = nrbdegelev(srf, 0, 1); -% -% See: -% -% bspdegelev - -% D.M. Spink -% Copyright (c) 2000. - -if nargin < 2 - error('Input argument must include the NURBS and degree increment.'); -end - -if ~isstruct(nurbs) - error('NURBS representation is not structure!'); -end - -if ~strcmp(nurbs.form,'B-NURBS') - error('Not a recognised NURBS representation'); -end - -degree = nurbs.order-1; - -if iscell(nurbs.knots) - % NURBS represents a surface - [dim,num1,num2] = size(nurbs.coefs); - - % Degree elevate along the v direction - if ntimes(2) == 0 - coefs = nurbs.coefs; - knots{2} = nurbs.knots{2}; - else - coefs = reshape(nurbs.coefs,4*num1,num2); - [coefs,knots{2}] = bspdegelev(degree(2),coefs,nurbs.knots{2},ntimes(2)); - num2 = size(coefs,2); - coefs = reshape(coefs,[4 num1 num2]); - end - - % Degree elevate along the u direction - if ntimes(1) == 0 - knots{1} = nurbs.knots{1}; - else - coefs = permute(coefs,[1 3 2]); - coefs = reshape(coefs,4*num2,num1); - [coefs,knots{1}] = bspdegelev(degree(1),coefs,nurbs.knots{1},ntimes(1)); - coefs = reshape(coefs,[4 num2 size(coefs,2)]); - coefs = permute(coefs,[1 3 2]); - end - -else - - % NURBS represents a curve - if isempty(ntimes) - coefs = nurbs.coefs; - knots = nurbs.knots; - else - [coefs,knots] = bspdegelev(degree,nurbs.coefs,nurbs.knots,ntimes); - end - -end - -% construct new NURBS -inurbs = nrbmak(coefs,knots); +function inurbs = nrbdegelev(nurbs, ntimes) +% +% NRBDEGELEV: Elevate the degree of the NURBS curve or surface. +% +% Calling Sequence: +% +% ecrv = nrbdegelev(crv,utimes); +% esrf = nrbdegelev(srf,{utimes,vtimes}); +% evol = nrbdegelev(vol,{utimes,vtimes,wtimes}); +% +% INPUT: +% +% crv : NURBS curve, see nrbmak. +% +% srf : NURBS surface, see nrbmak. +% +% vol : NURBS volume, see nrbmak. +% +% utimes : Increase the degree along U direction utimes. +% +% vtimes : Increase the degree along V direction vtimes. +% +% wtimes : Increase the degree along W direction vtimes. +% +% OUTPUT: +% +% ecrv : new NURBS structure for a curve with degree elevated. +% +% esrf : new NURBS structure for a surface with degree elevated. +% +% evol : new NURBS structure for a volume with degree elevated. +% +% +% Description: +% +% Degree elevates the NURBS curve or surface. This function uses the +% B-Spline function bspdegelev, which interface to an internal 'C' +% routine. +% +% Examples: +% +% Increase the NURBS surface twice along the V direction. +% esrf = nrbdegelev(srf, [0, 2]); +% +% See also: +% +% bspdegelev +% +% Copyright (C) 2000 Mark Spink, 2010 Rafel Vazquez +% +% 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, see <http://www.gnu.org/licenses/>. + +if nargin < 2 + error('Input argument must include the NURBS and degree increment.'); +end + +if ~isstruct(nurbs) + error('NURBS representation is not structure!'); +end + +if ~strcmp(nurbs.form,'B-NURBS') + error('Not a recognised NURBS representation'); +end + +degree = nurbs.order-1; + +if iscell(nurbs.knots) + if size(nurbs.knots,2) == 3 + % NURBS represents a volume + [dim,num1,num2,num3] = size(nurbs.coefs); + + % Degree elevate along the w direction + if ntimes(3) == 0 + coefs = nurbs.coefs; + knots{3} = nurbs.knots{3}; + else + coefs = reshape(nurbs.coefs,4*num1*num2,num3); + [coefs,knots{3}] = bspdegelev(degree(3),coefs,nurbs.knots{3},ntimes(3)); + num3 = size(coefs,2); + coefs = reshape(coefs,[4 num1 num2 num3]); + end + + % Degree elevate along the v direction + if ntimes(2) == 0 + knots{2} = nurbs.knots{2}; + else + coefs = permute(coefs,[1 2 4 3]); + coefs = reshape(coefs,4*num1*num3,num2); + [coefs,knots{2}] = bspdegelev(degree(2),coefs,nurbs.knots{2},ntimes(2)); + num2 = size(coefs,2); + coefs = reshape(coefs,[4 num1 num3 num2]); + coefs = permute(coefs,[1 2 4 3]); + end + + % Degree elevate along the u direction + if ntimes(1) == 0 + knots{1} = nurbs.knots{1}; + else + coefs = permute(coefs,[1 3 4 2]); + coefs = reshape(coefs,4*num2*num3,num1); + [coefs,knots{1}] = bspdegelev(degree(1),coefs,nurbs.knots{1},ntimes(1)); + coefs = reshape(coefs,[4 num2 num3 size(coefs,2)]); + coefs = permute(coefs,[1 4 2 3]); + end + + elseif size(nurbs.knots,2) == 2 + % NURBS represents a surface + [dim,num1,num2] = size(nurbs.coefs); + + % Degree elevate along the v direction + if ntimes(2) == 0 + coefs = nurbs.coefs; + knots{2} = nurbs.knots{2}; + else + coefs = reshape(nurbs.coefs,4*num1,num2); + [coefs,knots{2}] = bspdegelev(degree(2),coefs,nurbs.knots{2},ntimes(2)); + num2 = size(coefs,2); + coefs = reshape(coefs,[4 num1 num2]); + end + + % Degree elevate along the u direction + if ntimes(1) == 0 + knots{1} = nurbs.knots{1}; + else + coefs = permute(coefs,[1 3 2]); + coefs = reshape(coefs,4*num2,num1); + [coefs,knots{1}] = bspdegelev(degree(1),coefs,nurbs.knots{1},ntimes(1)); + coefs = reshape(coefs,[4 num2 size(coefs,2)]); + coefs = permute(coefs,[1 3 2]); + end + end +else + + % NURBS represents a curve + if isempty(ntimes) + coefs = nurbs.coefs; + knots = nurbs.knots; + else + [coefs,knots] = bspdegelev(degree,nurbs.coefs,nurbs.knots,ntimes); + end + +end + +% construct new NURBS +inurbs = nrbmak(coefs,knots); + +end + +%!demo +%! crv = nrbtestcrv; +%! plot(crv.coefs(1,:),crv.coefs(2,:),'bo') +%! title('Degree elevation along test curve: curve and control polygons.'); +%! hold on; +%! plot(crv.coefs(1,:),crv.coefs(2,:),'b--'); +%! nrbplot(crv,48); +%! +%! icrv = nrbdegelev(crv, 1); +%! +%! plot(icrv.coefs(1,:),icrv.coefs(2,:),'ro') +%! plot(icrv.coefs(1,:),icrv.coefs(2,:),'r--'); +%! +%! hold off;
--- a/extra/nurbs/inst/nrbdemo.m Sun Mar 21 18:40:41 2010 +0000 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,43 +0,0 @@ -%% Copyright (C) 2003 Mark Spink -%% -%% 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, see <http://www.gnu.org/licenses/>. - -%% nrbdemo () -%% demos for the nurbs package - -function nrbdemo() - -lst = {'3D Line', 'demoline', 'Rectangle','demorect', ... - 'Circle and Arc','democirc','Helix','demohelix', ... - 'Ellipse', 'demoellip', ... - 'Test Curve','democurve', ... - 'Test Surface','demosurf', ... - 'Bilinear Surface','demo4surf', ... - 'Knot Insertion','demokntins', ... - 'Degree Elevation','demodegelev', ... - 'Curve derivatives','demodercrv', ... - 'Surface derivatives','demodersrf' ... - 'Cylinderical arc','democylind', ... - 'Ruled Surface','demoruled', ... - 'Coons Surface','democoons', ... - 'Test Extrusion','demoextrude', ... - 'Test Revolution - Cup','demorevolve', ... - 'Test Revolution - Ball & Torus','demotorus', ... - '2D Geomtery','demogeom'}; - -for ii=1:2:length(lst) - printf("Demo number %d: %s\n", (ii+1)/2, lst{ii}) - eval(lst{ii+1}) - pause -endfor
--- a/extra/nurbs/inst/nrbderiv.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/nrbderiv.m Mon Mar 22 18:13:09 2010 +0000 @@ -1,18 +1,20 @@ function dnurbs = nrbderiv(nurbs) % % NRBDERIV: Construct the first derivative representation of a -% NURBS curve or surface. +% NURBS curve, surface or volume. % % Calling Sequence: % % ders = nrbderiv(nrb); % -% Parameters: +% INPUT: % % nrb : NURBS data structure, see nrbmak. +% +% OUTPUT: % % ders : A data structure that represents the first -% derivatives of a NURBS curve or surface. +% derivatives of a NURBS curve, surface or volume. % % Description: % @@ -32,12 +34,24 @@ % % See the function nrbdeval for an example. % -% See: +% See also: % % nrbdeval +% +% Copyright (C) 2000 Mark Spink, 2010 Rafael Vazquez +% +% 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. -% D.M. Spink -% Copyright (c) 2000. +% 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, see <http://www.gnu.org/licenses/>. if ~isstruct(nurbs) error('NURBS representation is not structure!'); @@ -50,26 +64,56 @@ degree = nurbs.order - 1; if iscell(nurbs.knots) - % NURBS structure represents a surface - - num1 = nurbs.number(1); - num2 = nurbs.number(2); + if size(nurbs.knots,2) == 3 + % NURBS structure represents a volume + num1 = nurbs.number(1); + num2 = nurbs.number(2); + num3 = nurbs.number(3); % taking derivatives along the u direction - dknots = nurbs.knots; - dcoefs = permute(nurbs.coefs,[1 3 2]); - dcoefs = reshape(dcoefs,4*num2,num1); - [dcoefs,dknots{1}] = bspderiv(degree(1),dcoefs,nurbs.knots{1}); - dcoefs = permute(reshape(dcoefs,[4 num2 size(dcoefs,2)]),[1 3 2]); - dnurbs{1} = nrbmak(dcoefs, dknots); + dknots = nurbs.knots; + dcoefs = permute(nurbs.coefs,[1 3 4 2]); + dcoefs = reshape(dcoefs,4*num2*num3,num1); + [dcoefs,dknots{1}] = bspderiv(degree(1),dcoefs,nurbs.knots{1}); + dcoefs = permute(reshape(dcoefs,[4 num2 num3 size(dcoefs,2)]),[1 4 2 3]); + dnurbs{1} = nrbmak(dcoefs, dknots); % taking derivatives along the v direction - dknots = nurbs.knots; - dcoefs = reshape(nurbs.coefs,4*num1,num2); - [dcoefs,dknots{2}] = bspderiv(degree(2),dcoefs,nurbs.knots{2}); - dcoefs = reshape(dcoefs,[4 num1 size(dcoefs,2)]); - dnurbs{2} = nrbmak(dcoefs, dknots); + dknots = nurbs.knots; + dcoefs = permute(nurbs.coefs,[1 2 4 3]); + dcoefs = reshape(dcoefs,4*num1*num3,num2); + [dcoefs,dknots{2}] = bspderiv(degree(2),dcoefs,nurbs.knots{2}); + dcoefs = permute(reshape(dcoefs,[4 num1 num3 size(dcoefs,2)]),[1 2 4 3]); + dnurbs{2} = nrbmak(dcoefs, dknots); + + % taking derivatives along the w direction + dknots = nurbs.knots; + dcoefs = reshape(nurbs.coefs,4*num1*num2,num3); + [dcoefs,dknots{3}] = bspderiv(degree(3),dcoefs,nurbs.knots{3}); + dcoefs = reshape(dcoefs,[4 num1 num2 size(dcoefs,2)]); + dnurbs{3} = nrbmak(dcoefs, dknots); + + elseif size(nurbs.knots,2) == 2 + % NURBS structure represents a surface + num1 = nurbs.number(1); + num2 = nurbs.number(2); + + % taking derivatives along the u direction + dknots = nurbs.knots; + dcoefs = permute(nurbs.coefs,[1 3 2]); + dcoefs = reshape(dcoefs,4*num2,num1); + [dcoefs,dknots{1}] = bspderiv(degree(1),dcoefs,nurbs.knots{1}); + dcoefs = permute(reshape(dcoefs,[4 num2 size(dcoefs,2)]),[1 3 2]); + dnurbs{1} = nrbmak(dcoefs, dknots); + + % taking derivatives along the v direction + dknots = nurbs.knots; + dcoefs = reshape(nurbs.coefs,4*num1,num2); + [dcoefs,dknots{2}] = bspderiv(degree(2),dcoefs,nurbs.knots{2}); + dcoefs = reshape(dcoefs,[4 num1 size(dcoefs,2)]); + dnurbs{2} = nrbmak(dcoefs, dknots); + end else % NURBS structure represents a curve @@ -78,7 +122,48 @@ end +end +%!demo +%! crv = nrbtestcrv; +%! nrbplot(crv,48); +%! title('First derivatives along a test curve.'); +%! +%! tt = linspace(0.0,1.0,9); +%! +%! dcrv = nrbderiv(crv); +%! +%! [p1, dp] = nrbdeval(crv,dcrv,tt); +%! +%! p2 = vecnorm(dp); +%! +%! hold on; +%! plot(p1(1,:),p1(2,:),'ro'); +%! h = quiver(p1(1,:),p1(2,:),p2(1,:),p2(2,:),0); +%! set(h,'Color','black'); +%! hold off; - - +%!demo +%! srf = nrbtestsrf; +%! p = nrbeval(srf,{linspace(0.0,1.0,20) linspace(0.0,1.0,20)}); +%! h = surf(squeeze(p(1,:,:)),squeeze(p(2,:,:)),squeeze(p(3,:,:))); +%! set(h,'FaceColor','blue','EdgeColor','blue'); +%! title('First derivatives over a test surface.'); +%! +%! npts = 5; +%! tt = linspace(0.0,1.0,npts); +%! dsrf = nrbderiv(srf); +%! +%! [p1, dp] = nrbdeval(srf, dsrf, {tt, tt}); +%! +%! up2 = vecnorm(dp{1}); +%! vp2 = vecnorm(dp{2}); +%! +%! hold on; +%! plot3(p1(1,:),p1(2,:),p1(3,:),'ro'); +%! h1 = quiver3(p1(1,:),p1(2,:),p1(3,:),up2(1,:),up2(2,:),up2(3,:)); +%! h2 = quiver3(p1(1,:),p1(2,:),p1(3,:),vp2(1,:),vp2(2,:),vp2(3,:)); +%! set(h1,'Color','black'); +%! set(h2,'Color','black'); +%! +%! hold off;
--- a/extra/nurbs/inst/nrbdeval.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/nrbdeval.m Mon Mar 22 18:13:09 2010 +0000 @@ -1,9 +1,10 @@ function [pnt,jac] = nrbdeval(nurbs, dnurbs, tt) -% NRBDEVAL: Evaluation of the derivative NURBS curve or surface. +% NRBDEVAL: Evaluation of the derivative NURBS curve, surface or volume. % % [pnt, jac] = nrbdeval(crv, dcrv, tt) % [pnt, jac] = nrbdeval(srf, dsrf, {tu tv}) +% [pnt, vol] = nrbdeval(vol, dvol, {tu tv tw}) % % INPUTS: % @@ -11,13 +12,19 @@ % % srf - original NUBRS surface % -% dcrv - NURBS derivative represention of crv +% vol - original NUBRS volume +% +% dcrv - NURBS derivative representation of crv % -% dsrf - NURBS derivative represention of surface +% dsrf - NURBS derivative representation of surface +% +% dvol - NURBS derivative representation of volume % % tt - parametric evaluation points -% If the nurbs is a surface then tt is a cell -% {tu, tv} are the parametric coordinates +% If the nurbs is a surface or a volume then tt is a cell +% {tu, tv} or {tu, tv, tw} are the parametric coordinates +% +% OUTPUT: % % pnt - evaluated points. % jac - evaluated first derivatives (Jacobian). @@ -29,9 +36,21 @@ % tt = linspace(0.0, 1.0, 9); % dcrv = nrbderiv(crv); % [pnts,jac] = nrbdeval(crv, dcrv, tt); +% +% Copyright (C) 2000 Mark Spink, 2010 Rafael Vazquez +% +% 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. -% D.M. Spink -% Copyright (c) 2000. +% 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, see <http://www.gnu.org/licenses/>. if ~isstruct(nurbs) error('NURBS representation is not structure!'); @@ -44,19 +63,37 @@ [cp,cw] = nrbeval(nurbs, tt); if iscell(nurbs.knots) - - % NURBS structure represents a surface - temp = cw(ones(3,1),:,:); - pnt = cp./temp; + if size(nurbs.knots,2) == 3 + % NURBS structure represents a volume + temp = cw(ones(3,1),:,:,:); + pnt = cp./temp; + + [cup,cuw] = nrbeval(dnurbs{1}, tt); + tempu = cuw(ones(3,1),:,:,:); + jac{1} = (cup-tempu.*pnt)./temp; - [cup,cuw] = nrbeval(dnurbs{1}, tt); - tempu = cuw(ones(3,1),:,:); - jac{1} = (cup-tempu.*pnt)./temp; + [cvp,cvw] = nrbeval(dnurbs{2}, tt); + tempv = cvw(ones(3,1),:,:,:); + jac{2} = (cvp-tempv.*pnt)./temp; + + [cwp,cww] = nrbeval(dnurbs{3}, tt); + tempw = cww(ones(3,1),:,:,:); + jac{3} = (cwp-tempw.*pnt)./temp; + + elseif size(nurbs.knots,2) == 2 + % NURBS structure represents a surface + temp = cw(ones(3,1),:,:); + pnt = cp./temp; - [cvp,cvw] = nrbeval(dnurbs{2}, tt); - tempv = cvw(ones(3,1),:,:); - jac{2} = (cvp-tempv.*pnt)./temp; + [cup,cuw] = nrbeval(dnurbs{1}, tt); + tempu = cuw(ones(3,1),:,:); + jac{1} = (cup-tempu.*pnt)./temp; + + [cvp,cvw] = nrbeval(dnurbs{2}, tt); + tempv = cvw(ones(3,1),:,:); + jac{2} = (cvp-tempv.*pnt)./temp; + end else % NURBS is a curve @@ -70,4 +107,91 @@ end +end +%!demo +%! crv = nrbtestcrv; +%! nrbplot(crv,48); +%! title('First derivatives along a test curve.'); +%! +%! tt = linspace(0.0,1.0,9); +%! +%! dcrv = nrbderiv(crv); +%! +%! [p1, dp] = nrbdeval(crv,dcrv,tt); +%! +%! p2 = vecnorm(dp); +%! +%! hold on; +%! plot(p1(1,:),p1(2,:),'ro'); +%! h = quiver(p1(1,:),p1(2,:),p2(1,:),p2(2,:),0); +%! set(h,'Color','black'); +%! hold off; + +%!demo +%! srf = nrbtestsrf; +%! p = nrbeval(srf,{linspace(0.0,1.0,20) linspace(0.0,1.0,20)}); +%! h = surf(squeeze(p(1,:,:)),squeeze(p(2,:,:)),squeeze(p(3,:,:))); +%! set(h,'FaceColor','blue','EdgeColor','blue'); +%! title('First derivatives over a test surface.'); +%! +%! npts = 5; +%! tt = linspace(0.0,1.0,npts); +%! dsrf = nrbderiv(srf); +%! +%! [p1, dp] = nrbdeval(srf, dsrf, {tt, tt}); +%! +%! up2 = vecnorm(dp{1}); +%! vp2 = vecnorm(dp{2}); +%! +%! hold on; +%! plot3(p1(1,:),p1(2,:),p1(3,:),'ro'); +%! h1 = quiver3(p1(1,:),p1(2,:),p1(3,:),up2(1,:),up2(2,:),up2(3,:)); +%! h2 = quiver3(p1(1,:),p1(2,:),p1(3,:),vp2(1,:),vp2(2,:),vp2(3,:)); +%! set(h1,'Color','black'); +%! set(h2,'Color','black'); +%! +%! hold off; + +%!test +%! knots{1} = [0 0 0 1 1 1]; +%! knots{2} = [0 0 0 .5 1 1 1]; +%! knots{3} = [0 0 0 0 1 1 1 1]; +%! cx = [0 0.5 1]; nx = length(cx); +%! cy = [0 0.25 0.75 1]; ny = length(cy); +%! cz = [0 1/3 2/3 1]; nz = length(cz); +%! coefs(1,:,:,:) = repmat(reshape(cx,nx,1,1),[1 ny nz]); +%! coefs(2,:,:,:) = repmat(reshape(cy,1,ny,1),[nx 1 nz]); +%! coefs(3,:,:,:) = repmat(reshape(cz,1,1,nz),[nx ny 1]); +%! coefs(4,:,:,:) = 1; +%! nurbs = nrbmak(coefs, knots); +%! x = rand(5,1); y = rand(5,1); z = rand(5,1); +%! tt = [x y z]'; +%! ders = nrbderiv(nurbs); +%! [points,jac] = nrbdeval(nurbs,ders,tt); +%! assert(points,tt,1e-10) +%! assert(jac{1}(1,:,:),ones(size(jac{1}(1,:,:))),1e-12) +%! assert(jac{2}(2,:,:),ones(size(jac{2}(2,:,:))),1e-12) +%! assert(jac{3}(3,:,:),ones(size(jac{3}(3,:,:))),1e-12) +%! +%!test +%! knots{1} = [0 0 0 1 1 1]; +%! knots{2} = [0 0 0 0 1 1 1 1]; +%! knots{3} = [0 0 0 1 1 1]; +%! cx = [0 0 1]; nx = length(cx); +%! cy = [0 0 0 1]; ny = length(cy); +%! cz = [0 0.5 1]; nz = length(cz); +%! coefs(1,:,:,:) = repmat(reshape(cx,nx,1,1),[1 ny nz]); +%! coefs(2,:,:,:) = repmat(reshape(cy,1,ny,1),[nx 1 nz]); +%! coefs(3,:,:,:) = repmat(reshape(cz,1,1,nz),[nx ny 1]); +%! coefs(4,:,:,:) = 1; +%! coefs = coefs([2 1 3 4],:,:,:); +%! nurbs = nrbmak(coefs, knots); +%! x = rand(5,1); y = rand(5,1); z = rand(5,1); +%! tt = [x y z]'; +%! dnurbs = nrbderiv(nurbs); +%! [points, jac] = nrbdeval(nurbs,dnurbs,tt); +%! assert(points,[y.^3 x.^2 z]',1e-10); +%! assert(jac{2}(1,:,:),3*y'.^2,1e-12) +%! assert(jac{1}(2,:,:),2*x',1e-12) +%! assert(jac{3}(3,:,:),ones(size(z')),1e-12) \ No newline at end of file
--- a/extra/nurbs/inst/nrbeval.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/nrbeval.m Mon Mar 22 18:13:09 2010 +0000 @@ -6,20 +6,27 @@ % % [p,w] = nrbeval(crv,ut) % [p,w] = nrbeval(srf,{ut,vt}) +% [p,w] = nrbeval(srf,{ut,vt,wt}) % -% Parameters: +% INPUT: % % crv : NURBS curve, see nrbmak. % % srf : NURBS surface, see nrbmak. +% +% vol : NURBS volume, see nrbmak. % % ut : Parametric evaluation points along U direction. % % vt : Parametric evaluation points along V direction. % -% p : Evaluated points on the NURBS curve or surface as cartesian -% coordinates (x,y,z). If w is included on the lhs argument list -% the points are returned as homogeneous coordinates (wx,wy,wz). +% wt : Parametric evaluation points along W direction. +% +% OUTPUT: +% +% p : Evaluated points on the NURBS curve, surface or volume as +% Cartesian coordinates (x,y,z). If w is included on the lhs argument +% list the points are returned as homogeneous coordinates (wx,wy,wz). % % w : Weights of the homogeneous coordinates of the evaluated % points. Note inclusion of this argument changes the type @@ -27,9 +34,9 @@ % % Description: % -% Evaluation of NURBS curves or surfaces at parametric points along the -% U and V directions. Either homogeneous coordinates are returned if the -% weights are requested in the lhs arguments, or as cartesian coordinates. +% Evaluation of NURBS curves, surfaces or volume at parametric points along +% the U, V and W directions. Either homogeneous coordinates are returned +% if the weights are requested in the lhs arguments, or as Cartesian coordinates. % This function utilises the 'C' interface bspeval. % % Examples: @@ -40,14 +47,25 @@ % ut = linspace(0.0,1.0,20); % p = nrbeval(nrb,ut); % -% See: +% See also: % % bspeval % +% Copyright (C) 2000 Mark Spink, 2010 Rafael Vazquez +% +% 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. -% D.M. Spink -% Copyright (c) 2000. -ctrpt=1; +% 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, see <http://www.gnu.org/licenses/>. + if nargin < 2 error('Not enough input arguments'); end @@ -66,70 +84,117 @@ end if iscell(nurbs.knots) - % NURBS structure represents a surface - - num1 = nurbs.number(1); - num2 = nurbs.number(2); - degree = nurbs.order-1; - - if iscell(tt) - % Evaluate over a [u,v] grid - % tt{1} represents the u direction - % tt{2} represents the v direction + if size(nurbs.knots,2) ==3 + %% NURBS structure represents a volume - nt1 = length(tt{1}); - nt2 = length(tt{2}); - - % Evaluate along the v direction - val = reshape(nurbs.coefs,4*num1,num2); - val = bspeval(degree(2),val,nurbs.knots{2},tt{2}); - val = reshape(val,[4 num1 nt2]); - - % Evaluate along the u direction - val = permute(val,[1 3 2]); - val = reshape(val,4*nt2,num1); - val = bspeval(degree(1),val,nurbs.knots{1},tt{1}); - val = reshape(val,[4 nt2 nt1]); - val = permute(val,[1 3 2]); + num1 = nurbs.number(1); + num2 = nurbs.number(2); + num3 = nurbs.number(3); + degree = nurbs.order-1; - w = val(4,:,:); - p = val(1:3,:,:); - if foption - p = p./repmat(w,[3 1 1]); + if iscell(tt) + tt = meshgrid(tt{3},tt{2},tt{1}); + tt = reshape(tt,3,[]); end - else + %% Evaluate at scattered points + %% tt(1,:) represents the u direction + %% tt(2,:) represents the v direction + %% tt(3,:) represents the w direction - % Evaluate at scattered points - % tt(1,:) represents the u direction - % tt(2,:) represents the v direction - + %% evaluate along the w direction nt = size(tt,2); + val = reshape(nurbs.coefs,4*num1*num2,num3); + val = bspeval(degree(3),val,nurbs.knots{3},tt(3,:)); + val = reshape(val,[4 num1 num2 nt]); - val = reshape(nurbs.coefs,4*num1,num2); - val = bspeval(degree(2),val,nurbs.knots{2},tt(2,:)); - val = reshape(val,[4 num1 nt]); + %% evaluate along the v direction + val2 = zeros(4*num1,nt); + for v = 1:nt + coefs = reshape(val(:,:,:,v),4*num1,num2); + val2(:,v) = bspeval(degree(2),coefs,nurbs.knots{2},tt(2,v)); + end + val2 = reshape(val2,[4 num1 nt]); - - % evaluate along the u direction + %% evaluate along the u direction pnts = zeros(4,nt); for v = 1:nt - coefs = squeeze(val(:,:,v)); + coefs = squeeze(val2(:,:,v)); pnts(:,v) = bspeval(degree(1),coefs,nurbs.knots{1},tt(1,v)); end w = pnts(4,:); p = pnts(1:3,:); if foption - p = p./repmat(w,[3, 1]); + p = p./repmat(w,[3, 1]); end + + elseif size(nurbs.knots,2) ==2 + %% NURBS structure represents a surface + + num1 = nurbs.number(1); + num2 = nurbs.number(2); + degree = nurbs.order-1; + + if iscell(tt) + %% Evaluate over a [u,v] grid + %% tt{1} represents the u direction + %% tt{2} represents the v direction + + nt1 = length(tt{1}); + nt2 = length(tt{2}); + + %% Evaluate along the v direction + val = reshape(nurbs.coefs,4*num1,num2); + val = bspeval(degree(2),val,nurbs.knots{2},tt{2}); + val = reshape(val,[4 num1 nt2]); + + %% Evaluate along the u direction + val = permute(val,[1 3 2]); + val = reshape(val,4*nt2,num1); + val = bspeval(degree(1),val,nurbs.knots{1},tt{1}); + val = reshape(val,[4 nt2 nt1]); + val = permute(val,[1 3 2]); + + w = val(4,:,:); + p = val(1:3,:,:); + if foption + p = p./repmat(w,[3 1 1]); + end + + else + + %% Evaluate at scattered points + %% tt(1,:) represents the u direction + %% tt(2,:) represents the v direction + + nt = size(tt,2); + + val = reshape(nurbs.coefs,4*num1,num2); + val = bspeval(degree(2),val,nurbs.knots{2},tt(2,:)); + val = reshape(val,[4 num1 nt]); + + + %% evaluate along the u direction + pnts = zeros(4,nt); + for v = 1:nt + coefs = squeeze(val(:,:,v)); + pnts(:,v) = bspeval(degree(1),coefs,nurbs.knots{1},tt(1,v)); + end + + w = pnts(4,:); + p = pnts(1:3,:); + if foption + p = p./repmat(w,[3, 1]); + end + end + end - else - % NURBS structure represents a curve - % tt represent a vector of parametric points in the u direction + %% NURBS structure represents a curve + %% tt represent a vector of parametric points in the u direction val = bspeval(nurbs.order-1,nurbs.coefs,nurbs.knots,tt); @@ -141,5 +206,65 @@ end +end +%!demo +%! srf = nrbtestsrf; +%! p = nrbeval(srf,{linspace(0.0,1.0,20) linspace(0.0,1.0,20)}); +%! h = surf(squeeze(p(1,:,:)),squeeze(p(2,:,:)),squeeze(p(3,:,:))); +%! title('Test surface.'); +%! hold off +%!test +%! knots{1} = [0 0 0 1 1 1]; +%! knots{2} = [0 0 0 .5 1 1 1]; +%! knots{3} = [0 0 0 0 1 1 1 1]; +%! cx = [0 0.5 1]; nx = length(cx); +%! cy = [0 0.25 0.75 1]; ny = length(cy); +%! cz = [0 1/3 2/3 1]; nz = length(cz); +%! coefs(1,:,:,:) = repmat(reshape(cx,nx,1,1),[1 ny nz]); +%! coefs(2,:,:,:) = repmat(reshape(cy,1,ny,1),[nx 1 nz]); +%! coefs(3,:,:,:) = repmat(reshape(cz,1,1,nz),[nx ny 1]); +%! coefs(4,:,:,:) = 1; +%! nurbs = nrbmak(coefs, knots); +%! x = rand(5,1); y = rand(5,1); z = rand(5,1); +%! tt = [x y z]'; +%! points = nrbeval(nurbs,tt); +%! +%! assert(points,tt,1e-10) +%! +%!test +%! knots{1} = [0 0 0 1 1 1]; +%! knots{2} = [0 0 0 0 1 1 1 1]; +%! knots{3} = [0 0 1 1]; +%! cx = [0 0 1]; nx = length(cx); +%! cy = [0 0 0 1]; ny = length(cy); +%! cz = [0 1]; nz = length(cz); +%! coefs(1,:,:,:) = repmat(reshape(cx,nx,1,1),[1 ny nz]); +%! coefs(2,:,:,:) = repmat(reshape(cy,1,ny,1),[nx 1 nz]); +%! coefs(3,:,:,:) = repmat(reshape(cz,1,1,nz),[nx ny 1]); +%! coefs(4,:,:,:) = 1; +%! nurbs = nrbmak(coefs, knots); +%! x = rand(5,1); y = rand(5,1); z = rand(5,1); +%! tt = [x y z]'; +%! points = nrbeval(nurbs,tt); +%! assert(points,[x.^2 y.^3 z]',1e-10); +%! +%!test +%! knots{1} = [0 0 0 1 1 1]; +%! knots{2} = [0 0 0 0 1 1 1 1]; +%! knots{3} = [0 0 1 1]; +%! cx = [0 0 1]; nx = length(cx); +%! cy = [0 0 0 1]; ny = length(cy); +%! cz = [0 1]; nz = length(cz); +%! coefs(1,:,:,:) = repmat(reshape(cx,nx,1,1),[1 ny nz]); +%! coefs(2,:,:,:) = repmat(reshape(cy,1,ny,1),[nx 1 nz]); +%! coefs(3,:,:,:) = repmat(reshape(cz,1,1,nz),[nx ny 1]); +%! coefs(4,:,:,:) = 1; +%! coefs = coefs([2 1 3 4],:,:,:); +%! nurbs = nrbmak(coefs, knots); +%! x = rand(5,1); y = rand(5,1); z = rand(5,1); +%! tt = [x y z]'; +%! points = nrbeval(nurbs,tt); +%! [y.^3 x.^2 z]'; +%! assert(points,[y.^3 x.^2 z]',1e-10);
--- a/extra/nurbs/inst/nrbextrude.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/nrbextrude.m Mon Mar 22 18:13:09 2010 +0000 @@ -1,18 +1,3 @@ -%% Copyright (C) 2003 Mark Spink -%% -%% 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, see <http://www.gnu.org/licenses/>. - function srf = nrbextrude(curve,vector) % @@ -22,11 +7,13 @@ % % srf = nrbextrude(crv,vec); % -% Parameters: +% INPUT: % % crv : NURBS curve to extrude, see nrbmak. % % vec : Vector along which the curve is extruded. +% +% OUTPUT: % % srf : NURBS surface constructed. % @@ -42,9 +29,21 @@ % Form a hollow cylinder by extruding a circle along the z-axis. % % srf = nrbextrude(nrbcirc, [0,0,1]); +% +% Copyright (C) 2000 Mark Spink +% +% 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. -% D.M. Spink -% Copyright (c) 2000. +% 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, see <http://www.gnu.org/licenses/>. if iscell(curve.knots) error('Nurb surfaces cannot be extruded!'); @@ -61,3 +60,11 @@ coefs = cat(3,curve.coefs,vectrans(vector)*curve.coefs); srf = nrbmak(coefs,{curve.knots, [0 0 1 1]}); +end + +%!demo +%! crv = nrbtestcrv; +%! srf = nrbextrude(crv,[0 0 5]); +%! nrbplot(srf,[40 10]); +%! title('Extrusion of a test curve along the z-axis'); +%! hold off \ No newline at end of file
--- a/extra/nurbs/inst/nrbkntins.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/nrbkntins.m Mon Mar 22 18:13:09 2010 +0000 @@ -1,118 +1,190 @@ -function inurbs = nrbkntins(nurbs,iknots) -% -% NRBKNTINS: Insert a single or multiple knots into a NURBS curve or -% surface. -% -% Calling Sequence: -% -% icrv = nrbkntins(crv,iuknots); -% isrf = nrbkntins(srf,{iuknots ivknots}); -% -% Parameters: -% -% crv : NURBS curve, see nrbmak. -% -% srf : NURBS surface, see nrbmak. -% -% iuknots : Knots to be inserted along U direction. -% -% ivknots : Knots to be inserted along V direction. -% -% icrv : new NURBS structure for a curve with knots inserted. -% -% isrf : new NURBS structure for a surface with knots inserted. -% -% Description: -% -% Inserts knots into the NURBS data structure, these can be knots at -% new positions or at the location of existing knots to increase the -% multiplicity. Note that the knot multiplicity cannot be increased -% beyond the order of the spline. Knots along the V direction can only -% inserted into NURBS surfaces, not curve that are always defined along -% the U direction. This function use the B-Spline function bspkntins, -% which interfaces to an internal 'C' routine. -% -% Examples: -% -% Insert two knots into a curve, one at 0.3 and another -% twice at 0.4 -% -% icrv = nrbkntins(crv, [0.3 0.4 0.4]) -% -% Insert into a surface two knots as (1) into the U knot -% sequence and one knot into the V knot sequence at 0.5. -% -% isrf = nrbkntins(srf, {[0.3 0.4 0.4] [0.5]}) -% -% Also see: -% -% bspkntins -% -% Note: -% -% No knot multiplicity will be increased beyond the order of the spline. - -% D.M. Spink -% Copyright (c) 2000. - -if nargin < 2 - error('Input argument must include the NURBS and knots to be inserted'); -end - -if ~isstruct(nurbs) - error('NURBS representation is not structure!'); -end - -if ~strcmp(nurbs.form,'B-NURBS') - error('Not a recognised NURBS representation'); -end - -degree = nurbs.order-1; - -if iscell(nurbs.knots) - % NURBS represents a surface - num1 = nurbs.number(1); - num2 = nurbs.number(2); - - % Insert knots along the v direction - if isempty(iknots{2}) - coefs = nurbs.coefs; - knots{2} = nurbs.knots{2}; - else - coefs = reshape(nurbs.coefs,4*num1,num2); - [coefs,knots{2}] = bspkntins(degree(2),coefs,nurbs.knots{2},iknots{2}); - num2 = size(coefs,2); - coefs = reshape(coefs,[4 num1 num2]); - end - - % Insert knots along the u direction - if isempty(iknots{1}) - knots{1} = nurbs.knots{1}; - else - coefs = permute(coefs,[1 3 2]); - coefs = reshape(coefs,4*num2,num1); - [coefs,knots{1}] = bspkntins(degree(1),coefs,nurbs.knots{1},iknots{1}); - coefs = reshape(coefs,[4 num2 size(coefs,2)]); - coefs = permute(coefs,[1 3 2]); - end -else - - % NURBS represents a curve - if isempty(iknots) - coefs = nurbs.coefs; - knots = nurbs.knots; - else - [coefs,knots] = bspkntins(degree,nurbs.coefs,nurbs.knots,iknots); - end - -end - -% construct new NURBS -inurbs = nrbmak(coefs,knots); - - - - - - - +function inurbs = nrbkntins(nurbs,iknots) +% +% NRBKNTINS: Insert a single or multiple knots into a NURBS curve, +% surface or volume. +% +% Calling Sequence: +% +% icrv = nrbkntins(crv,iuknots); +% isrf = nrbkntins(srf,{iuknots ivknots}); +% ivol = nrbkntins(vol,{iuknots ivknots iwknots}); +% +% INPUT: +% +% crv : NURBS curve, see nrbmak. +% +% srf : NURBS surface, see nrbmak. +% +% srf : NURBS volume, see nrbmak. +% +% iuknots : Knots to be inserted along U direction. +% +% ivknots : Knots to be inserted along V direction. +% +% iwknots : Knots to be inserted along W direction. +% +% OUTPUT: +% +% icrv : new NURBS structure for a curve with knots inserted. +% +% isrf : new NURBS structure for a surface with knots inserted. +% +% ivol : new NURBS structure for a volume with knots inserted. +% +% Description: +% +% Inserts knots into the NURBS data structure, these can be knots at +% new positions or at the location of existing knots to increase the +% multiplicity. Note that the knot multiplicity cannot be increased +% beyond the order of the spline. Knots along the V direction can only +% inserted into NURBS surfaces, not curve that are always defined along +% the U direction. This function use the B-Spline function bspkntins, +% which interfaces to an internal 'C' routine. +% +% Examples: +% +% Insert two knots into a curve, one at 0.3 and another +% twice at 0.4 +% +% icrv = nrbkntins(crv, [0.3 0.4 0.4]) +% +% Insert into a surface two knots as (1) into the U knot +% sequence and one knot into the V knot sequence at 0.5. +% +% isrf = nrbkntins(srf, {[0.3 0.4 0.4] [0.5]}) +% +% See also: +% +% bspkntins +% +% Note: +% +% No knot multiplicity will be increased beyond the order of the spline. +% +% Copyright (C) 2000 Mark Spink, 2010 Rafael Vazquez +% +% 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, see <http://www.gnu.org/licenses/>. + +if nargin < 2 + error('Input argument must include the NURBS and knots to be inserted'); +end + +if ~isstruct(nurbs) + error('NURBS representation is not structure!'); +end + +if ~strcmp(nurbs.form,'B-NURBS') + error('Not a recognised NURBS representation'); +end + +degree = nurbs.order-1; + +if iscell(nurbs.knots) + if size(nurbs.knots,2)==3 + % NURBS represents a volume + num1 = nurbs.number(1); + num2 = nurbs.number(2); + num3 = nurbs.number(3); + + % Insert knots along the w direction + if isempty(iknots{3}) + coefs = nurbs.coefs; + knots{3} = nurbs.knots{3}; + else + coefs = reshape(nurbs.coefs,4*num1*num2,num3); + [coefs,knots{3}] = bspkntins(degree(3),coefs,nurbs.knots{3},iknots{3}); + num3 = size(coefs,2); + coefs = reshape(coefs,[4 num1 num2 num3]); + end + + % Insert knots along the v direction + if isempty(iknots{2}) + knots{2} = nurbs.knots{2}; + else + coefs = permute(coefs,[1 2 4 3]); + coefs = reshape(coefs,4*num1*num3,num2); + [coefs,knots{2}] = bspkntins(degree(2),coefs,nurbs.knots{2},iknots{2}); + num2 = size(coefs,2); + coefs = reshape(coefs,[4 num1 num3 num2]); + coefs = permute(coefs,[1 2 4 3]); + end + + % Insert knots along the u direction + if isempty(iknots{1}) + knots{1} = nurbs.knots{1}; + else + coefs = permute(coefs,[1 3 4 2]); + coefs = reshape(coefs,4*num2*num3,num1); + [coefs,knots{1}] = bspkntins(degree(1),coefs,nurbs.knots{1},iknots{1}); + coefs = reshape(coefs,[4 num2 num3 size(coefs,2)]); + coefs = permute(coefs,[1 4 2 3]); + end + + elseif size(nurbs.knots,2)==2 + % NURBS represents a surface + num1 = nurbs.number(1); + num2 = nurbs.number(2); + + % Insert knots along the v direction + if isempty(iknots{2}) + coefs = nurbs.coefs; + knots{2} = nurbs.knots{2}; + else + coefs = reshape(nurbs.coefs,4*num1,num2); + [coefs,knots{2}] = bspkntins(degree(2),coefs,nurbs.knots{2},iknots{2}); + num2 = size(coefs,2); + coefs = reshape(coefs,[4 num1 num2]); + end + + % Insert knots along the u direction + if isempty(iknots{1}) + knots{1} = nurbs.knots{1}; + else + coefs = permute(coefs,[1 3 2]); + coefs = reshape(coefs,4*num2,num1); + [coefs,knots{1}] = bspkntins(degree(1),coefs,nurbs.knots{1},iknots{1}); + coefs = reshape(coefs,[4 num2 size(coefs,2)]); + coefs = permute(coefs,[1 3 2]); + end + end +else + + % NURBS represents a curve + if isempty(iknots) + coefs = nurbs.coefs; + knots = nurbs.knots; + else + [coefs,knots] = bspkntins(degree,nurbs.coefs,nurbs.knots,iknots); + end + +end + +% construct new NURBS +inurbs = nrbmak(coefs,knots); + +end + +%!demo +%! crv = nrbtestcrv; +%! plot(crv.coefs(1,:),crv.coefs(2,:),'bo') +%! title('Knot insertion along test curve: curve and control polygons.'); +%! hold on; +%! plot(crv.coefs(1,:),crv.coefs(2,:),'b--'); +%! +%! nrbplot(crv,48); +%! +%! icrv = nrbkntins(crv,[0.125 0.375 0.625 0.875] ); +%! plot(icrv.coefs(1,:),icrv.coefs(2,:),'ro') +%! plot(icrv.coefs(1,:),icrv.coefs(2,:),'r--'); +%! hold off \ No newline at end of file
--- a/extra/nurbs/inst/nrbline.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/nrbline.m Mon Mar 22 18:13:09 2010 +0000 @@ -7,11 +7,13 @@ % crv = nrbline() % crv = nrbline(p1,p2) % -% Parameters: +% INPUT: % % p1 : 2D or 3D cartesian coordinate of the start point. % % p2 : 2D or 3D cartesian coordinate of the end point. +% +% OUTPUT: % % crv : NURBS curve for a straight line. % @@ -20,9 +22,21 @@ % Constructs NURBS data structure for a straight line. If no rhs % coordinates are included the function returns a unit straight % line along the x-axis. +% +% Copyright (C) 2000 Mark Spink +% +% 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. -% D.M. Spink -% Copyright (c) 2000. +% 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, see <http://www.gnu.org/licenses/>. coefs = [zeros(3,2); ones(1,2)]; @@ -35,3 +49,11 @@ curve = nrbmak(coefs, [0 0 1 1]); +end + +%!demo +%! crv = nrbline([0.0 0.0 0.0]',[5.0 4.0 2.0]'); +%! nrbplot(crv,1); +%! grid on; +%! title('3D straight line.'); +%! hold off
--- a/extra/nurbs/inst/nrbmak.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/nrbmak.m Mon Mar 22 18:13:09 2010 +0000 @@ -1,132 +1,245 @@ -function nurbs = nrbmak(coefs,knots) -% -% NRBMAK: Construct the NURBS structure given the control points -% and the knots. -% -% Calling Sequence: -% -% nurbs = nrbmak(cntrl,knots); -% -% Parameters: -% -% cntrl : Control points, these can be either Cartesian or -% homogeneous coordinates. -% -% For a curve the control points are represented by a -% matrix of size (dim,nu) and for a surface a multidimensional -% array of size (dim,nu,nv). Where nu is number of points along -% the parametric U direction, and nv the number of points -% along the V direction. Dim is the dimension valid options -% are -% 2 .... (x,y) 2D Cartesian coordinates -% 3 .... (x,y,z) 3D Cartesian coordinates -% 4 .... (wx,wy,wz,w) 4D homogeneous coordinates -% -% knots : Non-decreasing knot sequence spanning the interval -% [0.0,1.0]. It's assumed that the curves and surfaces -% are clamped to the start and end control points by knot -% multiplicities equal to the spline order. -% For curve knots form a vector and for a surface the knot -% are stored by two vectors for U and V in a cell structure. -% {uknots vknots} -% -% nurbs : Data structure for representing a NURBS curve. -% -% NURBS Structure: -% -% Both curves and surfaces are represented by a structure that is -% compatible with the Spline Toolbox from Mathworks -% -% nurbs.form .... Type name 'B-NURBS' -% nurbs.dim .... Dimension of the control points -% nurbs.number .... Number of Control points -% nurbs.coefs .... Control Points -% nurbs.order .... Order of the spline -% nurbs.knots .... Knot sequence -% -% Note: the control points are always converted and stored within the -% NURBS structure as 4D homogeneous coordinates. A curve is always stored -% along the U direction, and the vknots element is an empty matrix. For -% a surface the spline degree is a vector [du,dv] containing the degree -% along the U and V directions respectively. -% -% Description: -% -% This function is used as a convenient means of constructing the NURBS -% data structure. Many of the other functions in the toolbox rely on the -% NURBS structure been correctly defined as shown above. The nrbmak not -% only constructs the proper structure, but also checks for consistency. -% The user is still free to build his own structure, in fact a few -% functions in the toolbox do this for convenience. -% -% Examples: -% -% Construct a 2D line from (0.0,0.0) to (1.5,3.0). -% For a straight line a spline of order 2 is required. -% Note that the knot sequence has a multiplicity of 2 at the -% start (0.0,0.0) and end (1.0 1.0) in order to clamp the ends. -% -% line = nrbmak([0.0 1.5; 0.0 3.0],[0.0 0.0 1.0 1.0]); -% nrbplot(line, 2); -% -% Construct a surface in the x-y plane i.e -% -% ^ (0.0,1.0) ------------ (1.0,1.0) -% | | | -% | V | | -% | | Surface | -% | | | -% | | | -% | (0.0,0.0) ------------ (1.0,0.0) -% | -% |------------------------------------> -% U -% -% coefs = cat(3,[0 0; 0 1],[1 1; 0 1]); -% knots = {[0 0 1 1] [0 0 1 1]} -% plane = nrbmak(coefs,knots); -% nrbplot(plane, [2 2]); - -% D.M. Spink -% Copyright (c) 2000. - -nurbs.form = 'B-NURBS'; -nurbs.dim = 4; -np = size(coefs); -dim = np(1); -if iscell(knots) - % constructing a surface - nurbs.number = np(2:3); - if (dim < 4) - nurbs.coefs = repmat([0.0 0.0 0.0 1.0]',[1 np(2:3)]); - nurbs.coefs(1:dim,:,:) = coefs; - else - nurbs.coefs = coefs; - end - uorder = size(knots{1},2)-np(2); - vorder = size(knots{2},2)-np(3); - uknots = sort(knots{1}); - vknots = sort(knots{2}); - uknots = (uknots-uknots(1))/(uknots(end)-uknots(1)); - vknots = (vknots-vknots(1))/(vknots(end)-vknots(1)); - nurbs.knots = {uknots vknots}; - nurbs.order = [uorder vorder]; - -else - - % constructing a curve - nurbs.number = np(2); - if (dim < 4) - nurbs.coefs = repmat([0.0 0.0 0.0 1.0]',[1 np(2)]); - nurbs.coefs(1:dim,:) = coefs; - else - nurbs.coefs = coefs; - end - nurbs.order = size(knots,2)-np(2); - knots = sort(knots); - nurbs.knots = (knots-knots(1))/(knots(end)-knots(1)); - -end - - - +function nurbs = nrbmak(coefs,knots) +% +% NRBMAK: Construct the NURBS structure given the control points +% and the knots. +% +% Calling Sequence: +% +% nurbs = nrbmak(cntrl,knots); +% +% INPUT: +% +% cntrl : Control points, these can be either Cartesian or +% homogeneous coordinates. +% +% For a curve the control points are represented by a +% matrix of size (dim,nu), for a surface a multidimensional +% array of size (dim,nu,nv), for a volume a multidimensional array +% of size (dim,nu,nv,nw). Where nu is number of points along +% the parametric U direction, nv the number of points along +% the V direction and nw the number of points along the W direction. +% dim is the dimension. Valid options +% are +% 2 .... (x,y) 2D Cartesian coordinates +% 3 .... (x,y,z) 3D Cartesian coordinates +% 4 .... (wx,wy,wz,w) 4D homogeneous coordinates +% +% knots : Non-decreasing knot sequence spanning the interval +% [0.0,1.0]. It's assumed that the geometric entities +% are clamped to the start and end control points by knot +% multiplicities equal to the spline order (open knot vector). +% For curve knots form a vector and for surfaces (volumes) +% the knots are stored by two (three) vectors for U and V (and W) +% in a cell structure {uknots vknots} ({uknots vknots wknots}). +% +% OUTPUT: +% +% nurbs : Data structure for representing a NURBS entity +% +% NURBS Structure: +% +% Both curves and surfaces are represented by a structure that is +% compatible with the Spline Toolbox from Mathworks +% +% nurbs.form .... Type name 'B-NURBS' +% nurbs.dim .... Dimension of the control points +% nurbs.number .... Number of Control points +% nurbs.coefs .... Control Points +% nurbs.order .... Order of the spline +% nurbs.knots .... Knot sequence +% +% Note: the control points are always converted and stored within the +% NURBS structure as 4D homogeneous coordinates. A curve is always stored +% along the U direction, and the vknots element is an empty matrix. For +% a surface the spline order is a vector [du,dv] containing the order +% along the U and V directions respectively. For a volume the order is +% a vector [du dv dw]. Recall that order = degree + 1. +% +% Description: +% +% This function is used as a convenient means of constructing the NURBS +% data structure. Many of the other functions in the toolbox rely on the +% NURBS structure been correctly defined as shown above. The nrbmak not +% only constructs the proper structure, but also checks for consistency. +% The user is still free to build his own structure, in fact a few +% functions in the toolbox do this for convenience. +% +% Examples: +% +% Construct a 2D line from (0.0,0.0) to (1.5,3.0). +% For a straight line a spline of order 2 is required. +% Note that the knot sequence has a multiplicity of 2 at the +% start (0.0,0.0) and end (1.0 1.0) in order to clamp the ends. +% +% line = nrbmak([0.0 1.5; 0.0 3.0],[0.0 0.0 1.0 1.0]); +% nrbplot(line, 2); +% +% Construct a surface in the x-y plane i.e +% +% ^ (0.0,1.0) ------------ (1.0,1.0) +% | | | +% | V | | +% | | Surface | +% | | | +% | | | +% | (0.0,0.0) ------------ (1.0,0.0) +% | +% |------------------------------------> +% U +% +% coefs = cat(3,[0 0; 0 1],[1 1; 0 1]); +% knots = {[0 0 1 1] [0 0 1 1]} +% plane = nrbmak(coefs,knots); +% nrbplot(plane, [2 2]); +% +% Copyright (C) 2000 Mark Spink, 2010 Rafael Vazquez +% +% 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, see <http://www.gnu.org/licenses/>. + +nurbs.form = 'B-NURBS'; +nurbs.dim = 4; +np = size(coefs); +dim = np(1); +if iscell(knots) + if size(knots,2) == 3 + if (numel(np) == 3) + np(4) = 1; + elseif (numel(np)==2) + np(3:4) = 1; + end + % constructing a volume + nurbs.number = np(2:4); + if (dim < 4) + nurbs.coefs = repmat([0.0 0.0 0.0 1.0]',[1 np(2:4)]); + nurbs.coefs(1:dim,:,:) = coefs; + else + nurbs.coefs = coefs; + end + uorder = size(knots{1},2)-np(2); + vorder = size(knots{2},2)-np(3); + worder = size(knots{3},2)-np(4); + uknots = sort(knots{1}); + vknots = sort(knots{2}); + wknots = sort(knots{3}); + uknots = (uknots-uknots(1))/(uknots(end)-uknots(1)); + vknots = (vknots-vknots(1))/(vknots(end)-vknots(1)); + wknots = (wknots-wknots(1))/(wknots(end)-wknots(1)); + nurbs.knots = {uknots vknots wknots}; + nurbs.order = [uorder vorder worder]; + + elseif size(knots,2) == 2 + if (numel(np)==2); np(3) = 1; end + % constructing a surface + nurbs.number = np(2:3); + if (dim < 4) + nurbs.coefs = repmat([0.0 0.0 0.0 1.0]',[1 np(2:3)]); + nurbs.coefs(1:dim,:,:) = coefs; + else + nurbs.coefs = coefs; + end + uorder = size(knots{1},2)-np(2); + vorder = size(knots{2},2)-np(3); + uknots = sort(knots{1}); + vknots = sort(knots{2}); + uknots = (uknots-uknots(1))/(uknots(end)-uknots(1)); + vknots = (vknots-vknots(1))/(vknots(end)-vknots(1)); + nurbs.knots = {uknots vknots}; + nurbs.order = [uorder vorder]; + + end + +else + + % constructing a curve + nurbs.number = np(2); + if (dim < 4) + nurbs.coefs = repmat([0.0 0.0 0.0 1.0]',[1 np(2)]); + nurbs.coefs(1:dim,:) = coefs; + else + nurbs.coefs = coefs; + end + nurbs.order = size(knots,2)-np(2); + knots = sort(knots); + nurbs.knots = (knots-knots(1))/(knots(end)-knots(1)); + +end + +end + +%!demo +%! pnts = [0.5 1.5 4.5 3.0 7.5 6.0 8.5; +%! 3.0 5.5 5.5 1.5 1.5 4.0 4.5; +%! 0.0 0.0 0.0 0.0 0.0 0.0 0.0]; +%! crv = nrbmak(pnts,[0 0 0 1/4 1/2 3/4 3/4 1 1 1]); +%! nrbplot(crv,100) +%! title('Test curve') +%! hold off + +%!demo +%! pnts = [0.5 1.5 4.5 3.0 7.5 6.0 8.5; +%! 3.0 5.5 5.5 1.5 1.5 4.0 4.5; +%! 0.0 0.0 0.0 0.0 0.0 0.0 0.0]; +%! crv = nrbmak(pnts,[0 0 0 0.1 1/2 3/4 3/4 1 1 1]); +%! nrbplot(crv,100) +%! title('Test curve with a slight variation of the knot vector') +%! hold off + +%!demo +%! pnts = zeros(3,5,5); +%! pnts(:,:,1) = [ 0.0 3.0 5.0 8.0 10.0; +%! 0.0 0.0 0.0 0.0 0.0; +%! 2.0 2.0 7.0 7.0 8.0]; +%! pnts(:,:,2) = [ 0.0 3.0 5.0 8.0 10.0; +%! 3.0 3.0 3.0 3.0 3.0; +%! 0.0 0.0 5.0 5.0 7.0]; +%! pnts(:,:,3) = [ 0.0 3.0 5.0 8.0 10.0; +%! 5.0 5.0 5.0 5.0 5.0; +%! 0.0 0.0 5.0 5.0 7.0]; +%! pnts(:,:,4) = [ 0.0 3.0 5.0 8.0 10.0; +%! 8.0 8.0 8.0 8.0 8.0; +%! 5.0 5.0 8.0 8.0 10.0]; +%! pnts(:,:,5) = [ 0.0 3.0 5.0 8.0 10.0; +%! 10.0 10.0 10.0 10.0 10.0; +%! 5.0 5.0 8.0 8.0 10.0]; +%! +%! knots{1} = [0 0 0 1/3 2/3 1 1 1]; +%! knots{2} = [0 0 0 1/3 2/3 1 1 1]; +%! +%! srf = nrbmak(pnts,knots); +%! nrbplot(srf,[20 20]) +%! title('Test surface') +%! hold off + +%!demo +%! coefs =[ 6.0 0.0 6.0 1; +%! -5.5 0.5 5.5 1; +%! -5.0 1.0 -5.0 1; +%! 4.5 1.5 -4.5 1; +%! 4.0 2.0 4.0 1; +%! -3.5 2.5 3.5 1; +%! -3.0 3.0 -3.0 1; +%! 2.5 3.5 -2.5 1; +%! 2.0 4.0 2.0 1; +%! -1.5 4.5 1.5 1; +%! -1.0 5.0 -1.0 1; +%! 0.5 5.5 -0.5 1; +%! 0.0 6.0 0.0 1]'; +%! knots = [0 0 0 0 .1 .2 .3 .4 .5 .6 .7 .8 .9 1 1 1 1]; +%! +%! crv = nrbmak(coefs,knots); +%! nrbplot(crv,100); +%! grid on; +%! title('3D helical curve.'); +%! hold off +
--- a/extra/nurbs/inst/nrbnumbasisfun.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/nrbnumbasisfun.m Mon Mar 22 18:13:09 2010 +0000 @@ -1,27 +1,12 @@ -%% Copyright (C) 2009 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 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, see <http://www.gnu.org/licenses/>. - function idx = nrbnumbasisfun (points, nrb) - +% % NRBNUMBASISFUN: Numbering of basis functions for NURBS % % Calling Sequence: % -% N = nrbbasisfun (u, crv) -% N = nrbbasisfun ({u, v}, srf) -% N = nrbbasisfun (p, srf) +% N = nrbnumbasisfun (u, crv) +% N = nrbnumbasisfun ({u, v}, srf) +% N = nrbnumbasisfun (p, srf) % % INPUT: % @@ -34,7 +19,21 @@ % % N - Indices of the basis functions that are nonvanishing at each % point. size(N) == size(B) -% +% +% Copyright (C) 2009 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 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, see <http://www.gnu.org/licenses/>. if ( (nargin<2) ... || (nargout>1) ... @@ -42,16 +41,16 @@ || (iscell(points) && ~iscell(nrb.knots)) ... || (~iscell(points) && iscell(nrb.knots) && (size(points,1)~=2)) ... ) - print_usage(); + error('Incorrect input arguments in nrbnumbasisfun'); end - if (~iscell(nrb.knots)) %% NURBS curve + if (~iscell(nrb.knots)) %% NURBS curve iv = findspan (nrb.number-1, nrb.order-1, points, nrb.knots); idx = numbasisfun (iv, points, nrb.order-1, nrb.knots); - else %% NURBS surface + elseif size(nrb.knots,2) == 2 %% NURBS surface if (iscell(points)) [v, u] = meshgrid(points{2}, points{1}); @@ -63,46 +62,21 @@ end idx = nrb_srf_numbasisfun__ (p, nrb); - + else + error('The function nrbnumbasisfun is not yet ready for volumes') end end -function idx = nrb_srf_numbasisfun__ (points, nrb) - - m = nrb.number(1)-1; - n = nrb.number(2)-1; - - npt = size(points,2); - u = points(1,:); - v = points(2,:); - - U = nrb.knots{1}; - V = nrb.knots{2}; - - p = nrb.order(1)-1; - q = nrb.order(2)-1; - - spu = findspan (m, p, u, U); - Ik = numbasisfun (spu, u, p, U); - - spv = findspan (n, q, v, V); - Jk = numbasisfun (spv, v, q, V); - - for k=1:npt - [Jkb, Ika] = meshgrid(Jk(k, :), Ik(k, :)); - idx(k, :) = sub2ind([m+1, n+1], Ika(:)+1, Jkb(:)+1); - end - -end %!test %! p = 2; q = 3; m = 4; n = 5; %! Lx = 1; Ly = 1; %! nrb = nrb4surf ([0 0], [1 0], [0 1], [1 1]); %! nrb = nrbdegelev (nrb, [p-1, q-1]); -%! nrb = nrbkntins (nrb, {linspace(0,1,m)(2:end-1), linspace(0,1,n)(2:end-1)}); -%! nrb.coefs (4,:,:) += rand (size (nrb.coefs (4,:,:))); +%! ikx = linspace(0,1,m); iky = linspace(0,1,n); +%! nrb = nrbkntins (nrb, {ikx(2:end-1), iky(2:end-1)}); +%! nrb.coefs (4,:,:) = nrb.coefs (4,:,:) + rand (size (nrb.coefs (4,:,:))); %! u = rand (1, 30); v = rand (1, 10); %! u = (u-min (u))/max (u-min (u)); %! v = (v-min (v))/max (v-min (v));
--- a/extra/nurbs/inst/nrbplot.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/nrbplot.m Mon Mar 22 18:13:09 2010 +0000 @@ -1,20 +1,4 @@ -%% Copyright (C) 2003 Mark Spink -%% -%% 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, see <http://www.gnu.org/licenses/>. - function nrbplot(nurbs,subd,p1,v1) - % % NRBPLOT: Plot a NURBS curve or surface. % @@ -23,8 +7,7 @@ % nrbplot(nrb,subd) % nrbplot(nrb,subd,p,v) % -% Parameters: -% +% INPUT: % % nrb : NURBS curve or surface, see nrbmak. % @@ -47,9 +30,21 @@ % and 30 along the V direction % % plot(nrbtestsrf, [20 30]) +% +% Copyright (C) 2000 Mark Spink +% +% 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. -% D.M. Spink -% Copyright (c) 2000. +% 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, see <http://www.gnu.org/licenses/>. nargs = nargin; if nargs < 2 @@ -99,7 +94,7 @@ % plot the curve or surface if iscell(nurbs.knots) - % plot a NURBS surface + if size(nurbs.knots,2) == 2 % plot a NURBS surface p = nrbeval(nurbs,{linspace(0.0,1.0,subd(1)) linspace(0.0,1.0,subd(2))}); if strcmp(light,'on') % light surface @@ -110,6 +105,9 @@ surf(squeeze(p(1,:,:)),squeeze(p(2,:,:)),squeeze(p(3,:,:))); shading faceted; end + else + error('The function nrbplot is not yet ready for volumes') + end else % plot a NURBS curve p = nrbeval(nurbs,linspace(0.0,1.0,subd)); @@ -125,8 +123,48 @@ end axis equal; +end + % plot the control surface % hold on; % mesh(squeeze(pnts(1,:,:)),squeeze(pnts(2,:,:)),squeeze(pnts(3,:,:))); % hold off; +%!demo +%! crv = nrbtestcrv; +%! nrbplot(crv,100) +%! title('Test curve') +%! hold off + +%!demo +%! coefs = [0.0 7.5 15.0 25.0 35.0 30.0 27.5 30.0; +%! 0.0 2.5 0.0 -5.0 5.0 15.0 22.5 30.0]; +%! knots = [0.0 0.0 0.0 1/6 1/3 1/2 2/3 5/6 1.0 1.0 1.0]; +%! +%! geom = [ +%! nrbmak(coefs,knots) +%! nrbline([30.0 30.0],[20.0 30.0]) +%! nrbline([20.0 30.0],[20.0 20.0]) +%! nrbcirc(10.0,[10.0 20.0],1.5*pi,0.0) +%! nrbline([10.0 10.0],[0.0 10.0]) +%! nrbline([0.0 10.0],[0.0 0.0]) +%! nrbcirc(5.0,[22.5 7.5]) +%! ]; +%! +%! ng = length(geom); +%! for i = 1:ng +%! nrbplot(geom(i),500); +%! hold on; +%! end +%! hold off; +%! axis equal; +%! title('2D Geometry formed by a series of NURBS curves'); + +%!demo +%! sphere = nrbrevolve(nrbcirc(1,[],0.0,pi),[0.0 0.0 0.0],[1.0 0.0 0.0]); +%! nrbplot(sphere,[40 40],'light','on'); +%! title('Ball and torus - surface construction by revolution'); +%! hold on; +%! torus = nrbrevolve(nrbcirc(0.2,[0.9 1.0]),[0.0 0.0 0.0],[1.0 0.0 0.0]); +%! nrbplot(torus,[40 40],'light','on'); +%! hold off
--- a/extra/nurbs/inst/nrbrect.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/nrbrect.m Mon Mar 22 18:13:09 2010 +0000 @@ -1,6 +1,6 @@ function curve = nrbrect(w,h) % -% NRBRECT: Construct NURBS representation of a rectangle. +% NRBRECT: Construct NURBS representation of a rectangular curve. % % Calling Sequence: % @@ -8,14 +8,16 @@ % crv = nrbrect(size) % crv = nrbrect(width, height) % -% Parameters: +% INPUT: % % size : Size of the square (width = height). % % width : Width of the rectangle (along x-axis). % % height : Height of the rectangle (along y-axis). -% +% +% OUTPUT: +% % crv : NURBS curve, see nrbmak. % % @@ -24,9 +26,21 @@ % Construct a rectangle or square in the x-y plane with the bottom % lhs corner at (0,0,0). If no rhs arguments provided the function % constructs a unit square. +% +% Copyright (C) 2000 Mark Spink +% +% 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. -% D.M. Spink -% Copyright (c) 2000. +% 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, see <http://www.gnu.org/licenses/>. if nargin < 1 w = 1; @@ -45,3 +59,12 @@ knots = [0 0 0.25 0.25 0.5 0.5 0.75 0.75 1 1]; curve = nrbmak(coefs, knots); + +end + +%!demo +%! crv = nrbtform(nrbrect(2,1), vecrotz(deg2rad(35))); +%! nrbplot(crv,4); +%! axis equal +%! title('Construction and rotation of a rectangular curve.'); +%! hold off \ No newline at end of file
--- a/extra/nurbs/inst/nrbreverse.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/nrbreverse.m Mon Mar 22 18:13:09 2010 +0000 @@ -6,9 +6,11 @@ % % rnrb = nrbreverse(nrb); % -% Parameters: +% INPUT: % % nrb : NURBS data structure, see nrbmak. +% +% OUTPUT: % % rnrb : Reversed NURBS. % @@ -16,20 +18,35 @@ % % Utility function to reverse the evaluation direction of a NURBS % curve or surface. +% +% Copyright (C) 2000 Mark Spink +% +% 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. -% D.M. Spink -% Copyright (c) 2000. +% 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, see <http://www.gnu.org/licenses/>. if nargin ~= 1 error('Incorrect number of input arguments'); end if iscell(nrb.knots) - + if size(nrb.knots,2) == 3 + error('The function nrbreverse is not yet ready for volumes') + else % reverse a NURBS surface coefs = nrb.coefs(:,:,end:-1:1); rnrb = nrbmak(coefs(:,end:-1:1,:), {1.0-fliplr(nrb.knots{1}),... - 1.0-fliplr(nrb.knots{2})}); + 1.0-fliplr(nrb.knots{2})}); + end else @@ -37,3 +54,25 @@ rnrb = nrbmak(fliplr(nrb.coefs), 1.0-fliplr(nrb.knots)); end + +end + +%!demo +%! pnts = [0.5 1.5 3.0 7.5 8.5; +%! 3.0 5.5 1.5 4.0 4.5; +%! 0.0 0.0 0.0 0.0 0.0]; +%! crv1 = nrbmak(pnts,[0 0 0 1/2 3/4 1 1 1]); +%! crv2 = nrbreverse(crv1); +%! fprintf('Knots of the original curve\n') +%! disp(crv1.knots) +%! fprintf('Knots of the reversed curve\n') +%! disp(crv2.knots) +%! fprintf('Control points of the original curve\n') +%! disp(crv1.coefs(1:2,:)) +%! fprintf('Control points of the reversed curve\n') +%! disp(crv2.coefs(1:2,:)) +%! nrbplot(crv1,100) +%! hold on +%! nrbplot(crv2,100) +%! title('The curve and its reverse are the same') +%! hold off \ No newline at end of file
--- a/extra/nurbs/inst/nrbrevolve.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/nrbrevolve.m Mon Mar 22 18:13:09 2010 +0000 @@ -1,18 +1,3 @@ -%% Copyright (C) 2003 Mark Spink -%% -%% 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, see <http://www.gnu.org/licenses/>. - function surf = nrbrevolve(curve,pnt,vec,theta) % @@ -22,7 +7,7 @@ % % srf = nrbrevolve(crv,pnt,vec[,ang]) % -% Parameters: +% INPUT: % % crv : NURBS curve to revolve, see nrbmak. % @@ -32,6 +17,10 @@ % vec : Vector defining the direction of the rotation axis. % % ang : Angle to revolve the curve, default 2*pi +% +% OUTPUT: +% +% srf : constructed surface % % Description: % @@ -66,14 +55,30 @@ % 6) rotate and vectrans the surface back into position % by reversing 1 and 2. % +% +% Copyright (C) 2000 Mark Spink +% +% 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. -% D.M. Spink -% Copyright (c) 2000. +% 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, see <http://www.gnu.org/licenses/>. if nargin < 3 error('Not enough arguments to construct revolved surface'); end +if size(curve.knots,2) == 3 + error('The function nrbrevolve is not yet ready for volumes') +end + if nargin < 4 theta = 2.0*pi; end @@ -116,6 +121,36 @@ RY = vecroty(angx); surf = nrbtform(surf,T*RY*RX); +end +%!demo +%! sphere = nrbrevolve(nrbcirc(1,[],0.0,pi),[0.0 0.0 0.0],[1.0 0.0 0.0]); +%! nrbplot(sphere,[40 40],'light','on'); +%! title('Ball and tori - surface construction by revolution'); +%! hold on; +%! torus = nrbrevolve(nrbcirc(0.2,[0.9 1.0]),[0.0 0.0 0.0],[1.0 0.0 0.0]); +%! nrbplot(torus,[40 40],'light','on'); +%! nrbplot(nrbtform(torus,vectrans([-1.8])),[20 10],'light','on'); +%! hold off; +%!demo +%! pnts = [3.0 5.5 5.5 1.5 1.5 4.0 4.5; +%! 0.0 0.0 0.0 0.0 0.0 0.0 0.0; +%! 0.5 1.5 4.5 3.0 7.5 6.0 8.5]; +%! crv = nrbmak(pnts,[0 0 0 1/4 1/2 3/4 3/4 1 1 1]); +%! +%! xx = vecrotz(deg2rad(25))*vecroty(deg2rad(15))*vecrotx(deg2rad(20)); +%! nrb = nrbtform(crv,vectrans([5 5])*xx); +%! +%! pnt = [5 5 0]'; +%! vec = xx*[0 0 1 1]'; +%! srf = nrbrevolve(nrb,pnt,vec(1:3)); +%! +%! p = nrbeval(srf,{linspace(0.0,1.0,100) linspace(0.0,1.0,100)}); +%! surfl(squeeze(p(1,:,:)),squeeze(p(2,:,:)),squeeze(p(3,:,:))); +%! title('Construct of a 3D surface by revolution of a curve.'); +%! shading interp; +%! colormap(copper); +%! axis equal; +%! hold off
--- a/extra/nurbs/inst/nrbruled.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/nrbruled.m Mon Mar 22 18:13:09 2010 +0000 @@ -1,18 +1,3 @@ -%% Copyright (C) 2003 Mark Spink -%% -%% 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, see <http://www.gnu.org/licenses/>. - function srf = nrbruled(crv1, crv2) % NRBRULED: Construct a ruled surface between two NURBS curves. @@ -21,11 +6,13 @@ % % srf = nrbruled(crv1, crv2) % -% Parameters: +% INPUT: % % crv1 : First NURBS curve, see nrbmak. % % crv2 : Second NURBS curve, see nrbmak. +% +% OUTPUT: % % srf : Ruled NURBS surface. % @@ -42,9 +29,21 @@ % line = nrbline([-1 0.5 1],[1 0.5 1]); % srf = nrbruled(cir,line); % nrbplot(srf,[20 20]); +% +% Copyright (C) 2000 Mark Spink +% +% 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. -% D.M. Spink -% Copyright (c) 2000. +% 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, see <http://www.gnu.org/licenses/>. if iscell(crv1.knots) | iscell(crv2.knots) error('Both NURBS must be curves'); @@ -76,3 +75,15 @@ coefs(:,:,2) = crv2.coefs; srf = nrbmak(coefs, {crv1.knots, [0 0 1 1]}); +end + +%!demo +%! pnts = [0.5 1.5 4.5 3.0 7.5 6.0 8.5; +%! 3.0 5.5 5.5 1.5 1.5 4.0 4.5; +%! 0.0 0.0 0.0 0.0 0.0 0.0 0.0]; +%! crv1 = nrbmak(pnts,[0 0 0 1/4 1/2 3/4 3/4 1 1 1]); +%! crv2 = nrbtform(nrbcirc(4,[4.5;0],pi,0.0),vectrans([0.0 4.0 -4.0])); +%! srf = nrbruled(crv1,crv2); +%! nrbplot(srf,[40 20]); +%! title('Ruled surface construction from two NURBS curves.'); +%! hold off \ No newline at end of file
--- a/extra/nurbs/inst/nrbsrfgradient.m Sun Mar 21 18:40:41 2010 +0000 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,149 +0,0 @@ -%% Copyright (C) 2009 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 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., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA - -%% -*- texinfo -*- -%% @deftypefn {Function File} {[@var{dzdx}, @var{dzdy}]=} nrbsrfgradient (@var{nrb}, @var{nrbder}, @var{u}, @var{v}) -%% NRBSRFGRADIENT: Compute the gradient of a NURBS surface. -%% @seealso{nrbderiv} -%% @end deftypefn - -%% Author: Carlo de Falco <cdf _AT_ users.sourceforge.net> -%% Created: 2009-03-17 - -function [dzdx, dzdy] = nrbsrfgradient (nrb, nrbder, u, v) - - if ((nargin ~= 4) || (nargout>2)) - print_usage(); - end - - [np, dp] = nrbdeval (nrb, nrbder, {u, v}); - - dxdu = squeeze (dp{1}(1,:,:)); - dydu = squeeze (dp{1}(2,:,:)); - dzdu = squeeze (dp{1}(3,:,:)); - dxdv = squeeze (dp{2}(1,:,:)); - dydv = squeeze (dp{2}(2,:,:)); - dzdv = squeeze (dp{2}(3,:,:)); - - detjac = dxdu.*dydv - dydu.*dxdv; - - dzdx = ( dydv .* dzdu - dydu .*dzdv)./detjac; - dzdy = (-dxdv .* dzdu + dxdu .*dzdv)./detjac; - -end - - -%!shared nrb, cntl, k, rec -%! -%!test -%! p = 5; n = 5; x = y = linspace(0,3,n+1); [cntl(1,:,:), cntl(2,:,:)] = meshgrid(x,y); -%! cntl(4,:,:) = 1; -%! k{1} = k{2} = [zeros(1,p) linspace(0,1,n-p+2) ones(1,p)]; -%! rec = nrbmak(cntl,k); -%! rec = nrbtform(rec, vecrotz(pi/3)); -%! rec.coefs(3,:,:) = rec.coefs(1,:,:); -%! recd = nrbderiv(rec); -%! [P, w] = nrbeval(rec, {linspace(0,1,20), -%! linspace(0,1,15)}); -%! [dzdx, dzdy] = nrbsrfgradient (rec, recd, linspace(0,1,20), -%! linspace(0,1,15)); -%! quiver(squeeze(P(1,:,:)), squeeze(P(2,:,:)), dzdx, dzdy) -%! assert(dzdx, ones(20,15), 10*eps); -%! assert(dzdy, zeros(20,15), 10*eps); -%! -%!test -%! rec = nrbtform(rec, vecrotz(pi/2)); -%! recd = nrbderiv(rec); -%! [P, w] = nrbeval(rec, {linspace(0,1,20), -%! linspace(0,1,15)}); -%! [dzdx, dzdy] = nrbsrfgradient (rec, recd, linspace(0,1,20), -%! linspace(0,1,15)); -%! quiver(squeeze(P(1,:,:)), squeeze(P(2,:,:)), dzdx, dzdy) -%! assert(dzdy, ones(20,15), 10*eps); -%! assert(dzdx, zeros(20,15), 10*eps); -%! -%!test -%! rec = nrbtform(rec, vecrotz(-pi/4)); -%! recd = nrbderiv(rec); -%! [P, w] = nrbeval(rec, {linspace(0,1,20), -%! linspace(0,1,15)}); -%! [dzdx, dzdy] = nrbsrfgradient (rec, recd, linspace(0,1,20), -%! linspace(0,1,15)); -%! quiver(squeeze(P(1,:,:)), squeeze(P(2,:,:)), dzdx, dzdy) -%! assert(dzdy, dzdx, 10*eps); - -%!shared isrf -%!test -%! mcp = 3; -%! ncp = 2; -%! p = 2; -%! q = 2; -%! cntl = ones(4,mcp+1,ncp+1); -%! cntl(3,:,:) = 0; -%! -%! sq2_1 = sqrt(2)-1; -%! sq2_2 = sqrt(2)/2; -%! -%! cntl(1,:,:) = [1 3/2 2; -%! 1 3/2 2; -%! sq2_1 sq2_1*3/2 sq2_1*2; -%! 0 0 0]; -%! -%! cntl(2,:,:) = [0 0 0; -%! sq2_1 sq2_1*3/2 sq2_1*2; -%! 1 3/2 2; -%! 1 3/2 2]; -%! -%! cntl(4,2:3,1:3) = (1+1/sqrt(2))/2; -%! -%! cntl(1,:,:) .*= cntl(4,:,:); -%! cntl(2,:,:) .*= cntl(4,:,:); -%! -%! u_knot=[0 0 0 1/2 1 1 1]; -%! v_knot=[0 0 0 1 1 1]; -%! knots = {u_knot, v_knot}; -%! -%! srf = nrbmak(cntl, knots); -%! isrf = nrbkntins(srf, {setdiff(linspace(0,1,15), u_knot), ... -%! setdiff(linspace(0,1,15), v_knot)}); -%! isrf.coefs(3,:,:) = isrf.coefs(1,:,:); -%! isrfd = nrbderiv(isrf); -%! [P, w] = nrbeval(isrf, {linspace(0,1,20), -%! linspace(0,1,15)}); -%! [dzdx, dzdy] = nrbsrfgradient (isrf, isrfd, linspace(0,1,20), -%! linspace(0,1,15)); -%! assert(dzdx, ones(size(dzdx)),10*eps); -%! assert(dzdy, zeros(size(dzdx)),10*eps); -%! -%!test -%! isrf = nrbtform(isrf, vecrotz(pi/4)); -%! isrfd = nrbderiv(isrf); -%! [P, w] = nrbeval(isrf, {linspace(0,1,20), -%! linspace(0,1,15)}); -%! [dzdx, dzdy] = nrbsrfgradient (isrf, isrfd, linspace(0,1,20), -%! linspace(0,1,15)); -%! assert(dzdx, dzdy,100*eps); -%! -%!test -%! isrf = nrbtform(isrf, vecrotz(-pi/4)); -%! isrf = nrbtform(isrf, vecroty(pi/4)); -%! isrfd = nrbderiv(isrf); -%! [P, w] = nrbeval(isrf, {linspace(0,1,20), -%! linspace(0,1,15)}); -%! [dzdx, dzdy] = nrbsrfgradient (isrf, isrfd, linspace(0,1,20), -%! linspace(0,1,15)); -%! assert(dzdx, zeros(size(dzdx)),100*eps); -%! assert(dzdy, zeros(size(dzdx)),100*eps);
--- a/extra/nurbs/inst/nrbsurfderiveval.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/nrbsurfderiveval.m Mon Mar 22 18:13:09 2010 +0000 @@ -1,65 +1,81 @@ -%% Copyright (C) 2009 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 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, see <http://www.gnu.org/licenses/>. +function skl = nrbsurfderiveval (srf, uv, d) -%% usage: skl = nrbsurfderiveval (srf, [u; v], d) -%% OUTPUT : skl (i, j, k, l) = i-th component derived j-1,k-1 times at the -%% l-th point. -%% -%% Adaptation of algorithm A4.4 from the NURBS book +% +% SURFDERIVEVAL: Compute the derivatives of a NURBS surface +% +% usage: skl = nrbsurfderiveval (srf, [u; v], d) +% +% INPUT: +% +% srf : NURBS surface structure, see nrbmak +% +% u, v : parametric coordinates of the point where we compute the +% derivatives +% +% d : number of partial derivatives to compute +% +% +% OUTPUT: +% +% skl (i, j, k, l) = i-th component derived j-1,k-1 times at the +% l-th point. +% +% Adaptation of algorithm A4.4 from the NURBS book, pg137 +% +% Copyright (C) 2009 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 2 of the License, or +% (at your option) any later version. -function skl = nrbsurfderiveval (srf, uv, d) +% 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, see <http://www.gnu.org/licenses/>. - for iu = 1:size(uv, 2); + for iu = 1:size(uv, 2); wders = squeeze (surfderiveval (srf.number(1)-1, srf.order(1)-1, ... srf.knots{1}, srf.number(2)-1, ... srf.order(2)-1, srf.knots{2}, ... squeeze (srf.coefs(4, :, :)), uv(1,iu), ... uv(2,iu), d)); - for idim = 1:3 + for idim = 1:3 Aders = squeeze (surfderiveval (srf.number(1)-1, srf.order(1)-1, ... srf.knots{1}, srf.number(2)-1, ... srf.order(2)-1, srf.knots{2}, ... - squeeze (srf.coefs(idim, :, :)), uv(1,iu), + squeeze (srf.coefs(idim, :, :)), uv(1,iu),... uv(2,iu), d)); - for k=0:d - for l=0:d-k + for k=0:d + for l=0:d-k v = Aders(k+1, l+1); - for j=1:l - v -= bincoeff(l,j)*wders(1,j+1)*skl(idim, k+1, l-j+1,iu); - endfor - for i=1:k - v -= bincoeff(k+1,i+1)*wders(i+1,1)*skl(idim, k-i+1, l+1, iu); + for j=1:l + v = v - nchoosek(l,j)*wders(1,j+1)*skl(idim, k+1, l-j+1,iu); + end + for i=1:k + v = v - nchoosek(k+1,i+1)*wders(i+1,1)*skl(idim, k-i+1, l+1, iu); v2 =0; - for j=1:l - v2 += bincoeff(l+1,j+1)*wders(i+1); - endfor - v -= bincoeff(k,i)*v2; - endfor + for j=1:l + v2 = v2 + nchoosek(l+1,j+1)*wders(i+1); + end + v = v - nchoosek(k,i)*v2; + end skl(idim, k+1, l+1, iu) = v/wders(1,1); - endfor - endfor - endfor + end + end + end - endfor + end -endfunction +end %!test -%! k = [0 0 1 1]; +%! k = [0 0 1 1]; %! c = [0 1]; %! [coef(2,:,:), coef(1,:,:)] = meshgrid (c, c); %! coef(3,:,:) = coef(1,:,:); @@ -67,11 +83,12 @@ %! [u, v] = meshgrid (linspace(0,1,11)); %! uv = [u(:)';v(:)']; %! skl = nrbsurfderiveval (srf, uv, 0); -%! assert (squeeze (skl (1:2,1,1,:)), nrbeval (srf, uv)(1:2,:), 1e3*eps) +%! aux = nrbeval(srf,uv); +%! assert (squeeze (skl (1:2,1,1,:)), aux(1:2,:), 1e3*eps) %!test -%! k = [0 0 1 1]; +%! k = [0 0 1 1]; %! c = [0 1]; %! [coef(2,:,:), coef(1,:,:)] = meshgrid (c, c); %! coef(3,:,:) = coef(1,:,:); @@ -80,7 +97,8 @@ %! [u, v] = meshgrid (linspace(0,1,11)); %! uv = [u(:)';v(:)']; %! skl = nrbsurfderiveval (srf, uv, 0); -%! assert (squeeze (skl (1:2,1,1,:)), nrbeval (srf, uv)(1:2,:), 1e3*eps) +%! aux = nrbeval(srf,uv); +%! assert (squeeze (skl (1:2,1,1,:)), aux(1:2,:), 1e3*eps) %!shared srf, uv %!test @@ -124,9 +142,9 @@ %! assert (squeeze (skl (1:2,1,2,:)), der{2}(1:2,:), 1e3*eps) %! %!test -%! p = q = 3; +%! p = 3; q = 3; %! mcp = 5; ncp = 5; -%! Lx = Ly = 10*rand(1); +%! Lx = 10*rand(1); Ly = Lx; %! srf = nrbdegelev (nrb4surf ([0 0], [Lx, 0], [0 Ly], [Lx Ly]), [p-1, q-1]); %! %%srf = nrbkntins (srf, {linspace(0,1,mcp-p+2)(2:end-1), linspace(0,1,ncp-q+2)(2:end-1)}); %! %%srf.coefs = permute (srf.coefs, [1 3 2]); @@ -157,14 +175,16 @@ %! P = squeeze(skl(:,1,1,:)); %! dPdx = squeeze(skl(:,2,1,:)); %! d2Pdx2 = squeeze(skl(:,3,1,:)); -%! assert (squeeze (skl (1:2,1,1,:)), nrbeval (srf, uv)(1:2,:), 1e3*eps) +%! aux = nrbeval(srf,uv); +%! assert (squeeze (skl (1:2,1,1,:)), aux(1:2,:), 1e3*eps) %!assert(P(3,:), 2*(P(1,:)-P(1,:).^2),100*eps) %!assert(dPdx(3,:), 2-4*P(1,:), 100*eps) %!assert(d2Pdx2(3,:), -4+0*P(1,:), 100*eps) %! %!test %! skl = nrbsurfderiveval (srf, uv, 0); -%! assert (squeeze (skl (1:2,1,1,:)), nrbeval (srf, uv)(1:2,:), 1e3*eps) +%! aux = nrbeval (srf, uv); +%! assert (squeeze (skl (1:2,1,1,:)), aux(1:2,:), 1e3*eps) %!shared dPdu, d2Pdu2, P, srf, uv %!test @@ -178,18 +198,20 @@ %! dPdu = squeeze(skl(:,2,1,:)); %! dPdv = squeeze(skl(:,1,2,:)); %! d2Pdu2 = squeeze(skl(:,3,1,:)); -%! assert (squeeze (skl (1:2,1,1,:)), nrbeval (srf, uv)(1:2,:), 1e3*eps) +%! aux = nrbeval(srf,uv); +%! assert (squeeze (skl (1:2,1,1,:)), aux(1:2,:), 1e3*eps) %!assert(dPdu(2,:), 3-4*P(1,:),100*eps) %!assert(d2Pdu2(2,:), -4+0*P(1,:),100*eps) %! %!test %! skl = nrbsurfderiveval (srf, uv, 0); -%! assert (squeeze (skl (1:2,1,1,:)), nrbeval (srf, uv)(1:2,:), 1e3*eps) +%! aux = nrbeval(srf,uv); +%! assert (squeeze (skl (1:2,1,1,:)), aux(1:2,:), 1e3*eps) %!test %! srf = nrb4surf([0 0], [1 0], [0 1], [1 1]); %! geo = nrbdegelev (srf, [3 3]); -%! geo.coefs (4, 2:end-1, 2:end-1) += .1 * rand (1, geo.number(1)-2, geo.number(2)-2); +%1 geo.coefs (4, 2:end-1, 2:end-1) = geo.coefs (4, 2:end-1, 2:end-1) + .1 * rand (1, geo.number(1)-2, geo.number(2)-2); %! geo = nrbkntins (geo, {[.1:.1:.9], [.2:.2:.8]}); %! [u, v] = meshgrid (linspace(0,1,10)); %! uv = [u(:)';v(:)'];
--- a/extra/nurbs/inst/nrbtestcrv.m Sun Mar 21 18:40:41 2010 +0000 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,14 +0,0 @@ -function crv = nrbtestcrv -% NRBTESTCRV: Constructs a simple test curve. - -% D.M. Spink -% Copyright (c) 2000 - -pnts = [0.5 1.5 4.5 3.0 7.5 6.0 8.5; - 3.0 5.5 5.5 1.5 1.5 4.0 4.5; - 0.0 0.0 0.0 0.0 0.0 0.0 0.0]; -crv = nrbmak(pnts,[0 0 0 1/4 1/2 3/4 3/4 1 1 1]); - - - -
--- a/extra/nurbs/inst/nrbtestsrf.m Sun Mar 21 18:40:41 2010 +0000 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,40 +0,0 @@ -function srf = nrbtestsrf -% NRBTESTSRF: Constructs a simple test surface. - -% D.M. Spink -% Copyright (c) 2000. - -% allocate multi-dimensional array of control points -pnts = zeros(3,5,5); - -% define a grid of control points -% in this case a regular grid of u,v points -% pnts(3,u,v) -% - -pnts(:,:,1) = [ 0.0 3.0 5.0 8.0 10.0; % w*x - 0.0 0.0 0.0 0.0 0.0; % w*y - 2.0 2.0 7.0 7.0 8.0]; % w*z - -pnts(:,:,2) = [ 0.0 3.0 5.0 8.0 10.0; - 3.0 3.0 3.0 3.0 3.0; - 0.0 0.0 5.0 5.0 7.0]; - -pnts(:,:,3) = [ 0.0 3.0 5.0 8.0 10.0; - 5.0 5.0 5.0 5.0 5.0; - 0.0 0.0 5.0 5.0 7.0]; - -pnts(:,:,4) = [ 0.0 3.0 5.0 8.0 10.0; - 8.0 8.0 8.0 8.0 8.0; - 5.0 5.0 8.0 8.0 10.0]; - -pnts(:,:,5) = [ 0.0 3.0 5.0 8.0 10.0; - 10.0 10.0 10.0 10.0 10.0; - 5.0 5.0 8.0 8.0 10.0]; - -% knots -knots{1} = [0 0 0 1/3 2/3 1 1 1]; % knots along u -knots{2} = [0 0 0 1/3 2/3 1 1 1]; % knots along v - -% make and draw nurbs surface -srf = nrbmak(pnts,knots);
--- a/extra/nurbs/inst/nrbtform.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/nrbtform.m Mon Mar 22 18:13:09 2010 +0000 @@ -1,5 +1,4 @@ function nurbs = nrbtform(nurbs,tmat) - % % NRBTFOR: Apply transformation matrix to the NURBS. % @@ -7,13 +6,15 @@ % % tnurbs = nrbtform(nurbs,tmatrix); % -% Parameters: +% INPUT: % % nurbs : NURBS data structure (see nrbmak for details). % % tmatrix : Transformation matrix, a matrix of size (4,4) defining % a single or multiple transformations. -% +% +% OUTPUT: +% % tnurbs : The return transformed NURBS data structure. % % Description: @@ -29,24 +30,51 @@ % Rotate a square by 45 degrees about the z axis. % % rsqr = nrbtform(nrbrect(), vecrotz(deg2rad(45))); -% nrbplot(rsqr, 10); +% nrbplot(rsqr, 1000); % -% See: +% See also: % % vecscale, vectrans, vecrotx, vecroty, vecrotz +% +% Copyright (C) 2000 Mark Spink +% +% 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. -% D.M. Spink -% Copyright (c) 2000 +% 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, see <http://www.gnu.org/licenses/>. if nargin < 2 error('Not enough input arguments!'); end; if iscell(nurbs.knots) + if size(nurbs.knots,2) == 2 % NURBS is a surface [dim,nu,nv] = size(nurbs.coefs); nurbs.coefs = reshape(tmat*reshape(nurbs.coefs,dim,nu*nv),[dim nu nv]); + elseif size(nurbs.knots,2) == 3 + % NURBS is a volume + error('The function nrbtform is not yet ready for volumes') + end else % NURBS is a curve nurbs.coefs = tmat*nurbs.coefs; end + +end + +%!demo +%! xx = vectrans([2.0 1.0])*vecroty(pi/8)*vecrotx(pi/4)*vecscale([1.0 2.0]); +%! c0 = nrbtform(nrbcirc, xx); +%! nrbplot(c0,50); +%! grid on +%! title('Construction of an ellipse by transforming a unit circle.'); +%! hold off \ No newline at end of file
--- a/extra/nurbs/inst/nrbtransp.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/nrbtransp.m Mon Mar 22 18:13:09 2010 +0000 @@ -6,9 +6,11 @@ % % tsrf = nrbtransp(srf) % -% Parameters: +% INPUT: % % srf : NURBS surface, see nrbmak. +% +% OUTPUT: % % tsrf : NURBS surface with U and V diretions transposed. % @@ -16,12 +18,38 @@ % % Utility function that transposes a NURBS surface, by swapping U and % V directions. NURBS curves cannot be transposed. +% +% Copyright (C) 2000 Mark Spink +% +% 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. -% D.M. Spink -% Copyright (c) 2000. +% 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, see <http://www.gnu.org/licenses/>. if ~iscell(srf.knots) - error(' A NURBSs curve cannot be transposed.'); + error(' A NURBS curve cannot be transposed.'); +elseif size(srf.knots,2) == 3 + error('The transposition of NURBS volumes has not been implemented.'); end tsrf = nrbmak(permute(srf.coefs,[1 3 2]), fliplr(srf.knots)); + +end + +%!demo +%! srf = nrb4surf([0 0 0], [1 0 1], [0 1 1], [1 1 2]); +%! nrbplot(srf,[20 5]); +%! title('Plane surface and its transposed (translated)') +%! hold on +%! srf.coefs(3,:,:) = srf.coefs(3,:,:) + 10; +%! srf = nrbtransp(srf); +%! nrbplot(srf,[20 5]); +%! hold off \ No newline at end of file
--- a/extra/nurbs/inst/numbasisfun.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/numbasisfun.m Mon Mar 22 18:13:09 2010 +0000 @@ -1,18 +1,3 @@ -%% Copyright (C) 2009 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 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, see <http://www.gnu.org/licenses/>. - function B = numbasisfun (iv, uv, p, U) % NUMBASISFUN: List non-zero Basis functions for B-Spline in a given knot-span @@ -32,17 +17,33 @@ % % N - Basis functions (numel(u)x(p+1)) % +% Copyright (C) 2009 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 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, see <http://www.gnu.org/licenses/>. B = bsxfun (@(a, b) a+b,iv-p, (0:p).').'; +end + %!test %! n = 3; %! U = [0 0 0 1/2 1 1 1]; %! p = 2; %! u = linspace (0, 1, 10); %! s = findspan (n, p, u, U); -%! Bref = [0 0 0 0 0 1 1 1 1 1 -%! 1 1 1 1 1 2 2 2 2 2 +%! Bref = [0 0 0 0 0 1 1 1 1 1; ... +%! 1 1 1 1 1 2 2 2 2 2; ... %! 2 2 2 2 2 3 3 3 3 3].'; %! B = numbasisfun (s, u, p, U); %! assert (B, Bref) \ No newline at end of file
--- a/extra/nurbs/inst/onebasisfunc.m Sun Mar 21 18:40:41 2010 +0000 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,77 +0,0 @@ -## Copyright (C) 2009 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 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 Octave; see the file COPYING. If not, see -## <http://www.gnu.org/licenses/>. - -## usage: -## N = onebasisfunc (n, u, i, p, U) - -## Author: Carlo de Falco <cdf _AT_ users.sourceforge.net> -## Created: 2009-06-25 - -function N = onebasisfunc (n, u, i, p, U) - - m = length (U)-1; - n = mod (n, m); - - spans = mod (n:n+p, m); - splen = U(spans+2) - U(spans+1); - - N = 0; - if (! any (i == spans)) - return; - elseif (p == 0) - N = 1; - return; - endif - - span = find (i == spans); - - if ((span != 1) || (u > U(spans(span)+1))) - - ld = sum (splen(1:p)); - if (ld != 0) - ln = sum (splen(1:span-1)) + max (u - U(spans(span)+1), 0); - Nl = onebasisfun (n, u, i, p - 1, U); - N += Nl * ln/ld; - endif - endif - - if ((span != p+1) || (u < U(spans(span)+2))) - rd = sum (splen(end-p+1:end)); - if (rd != 0) - rn = sum (splen(span+1:end)) + max (U(spans(span)+2)-u, 0); - Nr = onebasisfun (n+1, u, i, p - 1, U); - N += Nr * rn/rd; - endif - endif - -endfunction - -%!demo -%! p = 1; -%! U = sort([linspace(0, 1, 11) .5 .5 .5 ]); -%! uv= linspace (0, 1, 201); -%! for num = 0 %:length(U) -%! ii = 1; -%! for u = uv -%! s = findspanc (u, U); -%! N(ii++) = onebasisfunc (num, u, s, p, U); -%! endfor -%! plot (uv, N, [num2str(mod(num,6)) "-"]); -%! hold on; -%! endfor -%! -%! hold off -%! \ No newline at end of file
--- a/extra/nurbs/inst/onebasisfunderc.m Sun Mar 21 18:40:41 2010 +0000 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,72 +0,0 @@ -## Copyright (C) 2009 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 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 Octave; see the file COPYING. If not, see -## <http://www.gnu.org/licenses/>. - -## usage: -## N = onebasisfunderc (n, u, i, p, U) - -## Author: Carlo de Falco <cdf _AT_ users.sourceforge.net> -## Created: 2009-06-25 - -function N = onebasisfunderc (n, u, i, p, U) - m = length (U)-1; - n = mod (n, m); - - spans = mod (n:n+p, m); - splen = U(spans+2) - U(spans+1); - - N = 0; - if (! any (i == spans) || p == 0) - return; - endif - - span = find (i == spans); - - if ((span != 1) || (u > U(spans(span)+1))) - ld = sum (splen(1:p)); - if (ld != 0) - ln = p; - Nl = onebasisfunc (n, u, i, p - 1, U); - N += Nl * ln/ld; - endif - endif - - if ((span != p+1) || (u < U(spans(span)+2))) - rd = sum (splen(end-p+1:end)); - if (rd != 0) - rn = p; - Nr = onebasisfunc (n+1, u, i, p - 1, U); - N -= Nr * rn/rd; - endif - endif - -endfunction - -%!demo -%! p = 1; -%! U = sort([linspace(0, 1, 11) ]); -%! uv= linspace (0, 1, 201); -%! for num = 0 %:length(U) -%! ii = 1; -%! for u = uv -%! s = findspanc (u, U); -%! N(ii++) = onebasisfunderc (num, u, s, p, U); -%! endfor -%! plot (uv, N, [num2str(mod(num,6)) "-"]); -%! hold on; -%! endfor -%! -%! hold off -%! \ No newline at end of file
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/extra/nurbs/inst/private/nrb_crv_basisfun__.m Mon Mar 22 18:13:09 2010 +0000 @@ -0,0 +1,20 @@ + function [B, nbfu] = nrb_crv_basisfun__ (points, nrb); +% __NRB_CRV_BASISFUN__: Undocumented internal function +% +% Copyright (C) 2009 Carlo de Falco +% This software comes with ABSOLUTELY NO WARRANTY; see the file +% COPYING for details. This is free software, and you are welcome +% to distribute it under the conditions laid out in COPYING. + n = size (nrb.coefs, 2) -1; + p = nrb.order -1; + u = points; + U = nrb.knots; + w = nrb.coefs(4,:); + + spu = findspan (n, p, u, U); + nbfu = numbasisfun (spu, u, p, U); + + N = w(nbfu+1) .* basisfun (spu, u, p, U); + B = bsxfun (@(x,y) x./y, N, sum (N,2)); + + end
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/extra/nurbs/inst/private/nrb_crv_basisfun_der__.m Mon Mar 22 18:13:09 2010 +0000 @@ -0,0 +1,32 @@ + function [Bu, nbfu] = nrb_crv_basisfun_der__ (points, nrb); +% __NRB_CRV_BASISFUN_DER__: Undocumented internal function +% +% Copyright (C) 2009 Carlo de Falco +% This software comes with ABSOLUTELY NO WARRANTY; see the file +% COPYING for details. This is free software, and you are welcome +% to distribute it under the conditions laid out in COPYING. + n = size (nrb.coefs, 2) -1; + p = nrb.order -1; + u = points; + U = nrb.knots; + w = nrb.coefs(4,:); + + spu = findspan (n, p, u, U); + nbfu = numbasisfun (spu, u, p, U); + + N = basisfun (spu, u, p, U); + + Nprime = basisfunder (spu, p, u, U, 1); + Nprime = squeeze(Nprime(:,2,:)); + + + [Dpc, Dpk] = bspderiv (p, w, U); + D = bspeval (p, w, U, u); + Dprime = bspeval (p-1, Dpc, Dpk, u); + + + Bu1 = bsxfun (@(np, d) np/d , Nprime.', D); + Bu2 = bsxfun (@(n, dp) n*dp, N.', Dprime./D.^2); + Bu = w(nbfu+1) .* (Bu1 - Bu2).'; + + end
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/extra/nurbs/inst/private/nrb_srf_basisfun__.m Mon Mar 22 18:13:09 2010 +0000 @@ -0,0 +1,41 @@ +function [B, N] = nrb_srf_basisfun__ (points, nrb); + +% __NRB_SRF_BASISFUN__: Undocumented internal function +% +% Copyright (C) 2009 Carlo de Falco +% This software comes with ABSOLUTELY NO WARRANTY; see the file +% COPYING for details. This is free software, and you are welcome +% to distribute it under the conditions laid out in COPYING. + + m = size (nrb.coefs, 2) -1; + n = size (nrb.coefs, 3) -1; + + p = nrb.order(1) -1; + q = nrb.order(2) -1; + + u = points(1,:); + v = points(2,:); + npt = length(u); + + U = nrb.knots{1}; + V = nrb.knots{2}; + + w = squeeze(nrb.coefs(4,:,:)); + + spu = findspan (m, p, u, U); + spv = findspan (n, q, v, V); + NuIkuk = basisfun (spu, u, p, U); + NvJkvk = basisfun (spv, v, q, V); + + indIkJk = nrbnumbasisfun (points, nrb); + + for k=1:npt + wIkaJkb(1:p+1, 1:q+1) = reshape (w(indIkJk(k, :)), p+1, q+1); + NuIkukaNvJkvk(1:p+1, 1:q+1) = (NuIkuk(k, :).' * NvJkvk(k, :)); + RIkJk(k, :) = reshape((NuIkukaNvJkvk .* wIkaJkb ./ sum(sum(NuIkukaNvJkvk .* wIkaJkb))),1,[]); + end + + B = RIkJk; + N = indIkJk; + + end \ No newline at end of file
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/extra/nurbs/inst/private/nrb_srf_basisfun_der__.m Mon Mar 22 18:13:09 2010 +0000 @@ -0,0 +1,52 @@ +function [Bu, Bv, N] = nrb_srf_basisfun_der__ (points, nrb); + +% __NRB_SRF_BASISFUN_DER__: Undocumented internal function +% +% Copyright (C) 2009 Carlo de Falco +% This software comes with ABSOLUTELY NO WARRANTY; see the file +% COPYING for details. This is free software, and you are welcome +% to distribute it under the conditions laid out in COPYING. + + m = size (nrb.coefs, 2) -1; + n = size (nrb.coefs, 3) -1; + + p = nrb.order(1) -1; + q = nrb.order(2) -1; + + u = points(1,:); + v = points(2,:); + npt = length(u); + + U = nrb.knots{1}; + V = nrb.knots{2}; + + w = squeeze(nrb.coefs(4,:,:)); + + spu = findspan (m, p, u, U); + spv = findspan (n, q, v, V); + N = nrbnumbasisfun (points, nrb); + + NuIkuk = basisfun (spu, u, p, U); + NvJkvk = basisfun (spv, v, q, V); + + NuIkukprime = basisfunder (spu, p, u, U, 1); + NuIkukprime = reshape (NuIkukprime(:,2,:), npt, []); + + NvJkvkprime = basisfunder (spv, q, v, V, 1); + NvJkvkprime = reshape (NvJkvkprime(:,2,:), npt, []); + + for k=1:npt + wIkaJkb(1:p+1, 1:q+1) = reshape (w(N(k, :)), p+1, q+1); + + Num = (NuIkuk(k, :).' * NvJkvk(k, :)) .* wIkaJkb; + Num_du = (NuIkukprime(k, :).' * NvJkvk(k, :)) .* wIkaJkb; + Num_dv = (NuIkuk(k, :).' * NvJkvkprime(k, :)) .* wIkaJkb; + Denom = sum(sum(Num)); + Denom_du = sum(sum(Num_du)); + Denom_dv = sum(sum(Num_dv)); + + Bu(k, :) = reshape((Num_du/Denom - Denom_du.*Num/Denom.^2),1,[]); + Bv(k, :) = reshape((Num_dv/Denom - Denom_dv.*Num/Denom.^2),1,[]); + end + +end \ No newline at end of file
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/extra/nurbs/inst/private/nrb_srf_numbasisfun__.m Mon Mar 22 18:13:09 2010 +0000 @@ -0,0 +1,34 @@ +function idx = nrb_srf_numbasisfun__ (points, nrb) + +% __NRB_SRF_NUMBASISFUN__: Undocumented internal function +% +% Copyright (C) 2009 Carlo de Falco +% This software comes with ABSOLUTELY NO WARRANTY; see the file +% COPYING for details. This is free software, and you are welcome +% to distribute it under the conditions laid out in COPYING. + + m = nrb.number(1)-1; + n = nrb.number(2)-1; + + npt = size(points,2); + u = points(1,:); + v = points(2,:); + + U = nrb.knots{1}; + V = nrb.knots{2}; + + p = nrb.order(1)-1; + q = nrb.order(2)-1; + + spu = findspan (m, p, u, U); + Ik = numbasisfun (spu, u, p, U); + + spv = findspan (n, q, v, V); + Jk = numbasisfun (spv, v, q, V); + + for k=1:npt + [Jkb, Ika] = meshgrid(Jk(k, :), Ik(k, :)); + idx(k, :) = sub2ind([m+1, n+1], Ika(:)+1, Jkb(:)+1); + end + +end
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/extra/nurbs/inst/private/onebasisfun__.m Mon Mar 22 18:13:09 2010 +0000 @@ -0,0 +1,30 @@ +function N = onebasisfun__ (u, p, U) + +% __ONEBASISFUN__: Undocumented internal function +% +% Copyright (C) 2009 Carlo de Falco +% This software comes with ABSOLUTELY NO WARRANTY; see the file +% COPYING for details. This is free software, and you are welcome +% to distribute it under the conditions laid out in COPYING. + + N = 0; + if (~ any (U <= u)) || (~ any (U >= u)) + return; + elseif (p == 0) + N = 1; + return; + end + + ln = u - U(1); + ld = U(end-1) - U(1); + if (ld ~= 0) + N = N + ln * onebasisfun__ (u, p-1, U(1:end-1))/ ld; + end + + dn = U(end) - u; + dd = U(end) - U(2); + if (dd ~= 0) + N = N + dn * onebasisfun__ (u, p-1, U(2:end))/ dd; + end + +end
--- a/extra/nurbs/inst/rad2deg.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/rad2deg.m Mon Mar 22 18:13:09 2010 +0000 @@ -6,10 +6,12 @@ % % rad = rad2deg(deg); % -% Parameters: +% INPUT: % % rad : Angle in radians. % +% OUTPUT: +% % deg : Angle in degrees. % % Description: @@ -22,9 +24,22 @@ % Convert 0.3 radians to degrees % % rad = deg2rad(0.3); +% +% Copyright (C) 2000 Mark Spink +% +% 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. -% D.M. Spink -% Copyright (c) 2000. +% 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, see <http://www.gnu.org/licenses/>. deg = 180.0*rad/pi; +end
--- a/extra/nurbs/inst/surfderivcpts.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/surfderivcpts.m Mon Mar 22 18:13:09 2010 +0000 @@ -1,42 +1,42 @@ -%% Copyright (C) 2009 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 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, see <http://www.gnu.org/licenses/>. +function pkl = surfderivcpts (n, p, U, m, q, V, P, d, r1, r2, s1, s2) +% +% SURFDERIVCPTS: Compute control points of n-th derivatives of a NURBS surface. +% +% usage: pkl = surfderivcpts (n, p, U, m, q, V, P, d) +% +% INPUT: +% n+1, m+1 = number of control points +% p, q = spline order +% U, V = knots +% P = control points +% d = derivative order +% OUTPUT: +% pkl (k+1, l+1, i+1, j+1) = i,jth control point +% of the surface differentiated k +% times in the u direction and l +% times in the v direction +% +% Adaptation of algorithm A3.7 from the NURBS book, pg114 +% +% Copyright (C) 2009 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 2 of the License, or +% (at your option) any later version. -%% SURFDERIVCPTS: Compute control points of n-th derivatives of a NURBS surface. -%% -%% usage: pkl = surfderivcpts (n, p, U, m, q, V, P, d) -%% -%% INPUT: -%% n+1, m+1 = number of control points -%% p, q = spline order -%% U, V = knots -%% P = control points -%% d = derivative order -%% OUTPUT: -%% pkl (k+1, l+1, i+1, j+1) = i,jth control point -%% of the surface differentiated k -%% times in the u direction and l -%% times in the v direction -%% -%% Adaptation of algorithm A3.7 from the NURBS book - -function pkl = surfderivcpts (n, p, U, m, q, V, P, d, r1, r2, s1, s2) +% 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, see <http://www.gnu.org/licenses/>. if (nargin <= 8) r1 = 0; r2 = n; s1 = 0; s2 = m; - endif + end r = r2-r1; s = s2-s1; @@ -46,10 +46,10 @@ temp = curvederivcpts (n, p, U, P(:,j+1:end), du, r1, r2); for k=0:du for i=0:r-k - pkl (k+1, 1, i+1, j-s1+1) = temp (k+1, i+1); - endfor - endfor - endfor + pkl (k+1, 1, i+1, j-s1+1) = temp (k+1, i+1); + end + end + end for k=0:du for i=0:r-k @@ -57,20 +57,20 @@ temp = curvederivcpts (m, q, V(s1+1:end), pkl(k+1, 1, i+1, :), ... dd, 0, s); for l=1:dd - for j=0:s-l - pkl (k+1, l+1, i+1, j+1) = temp (l+1, j+1); - endfor - endfor - endfor - endfor + for j=0:s-l + pkl (k+1, l+1, i+1, j+1) = temp (l+1, j+1); + end + end + end + end -endfunction +end %!test %! coefs = cat(3,[0 0; 0 1],[1 1; 0 1]); %! knots = {[0 0 1 1] [0 0 1 1]}; %! plane = nrbmak(coefs,knots); -%! pkl = surfderivcpts (plane.number(1)-1, plane.order(1)-1, -%! plane.knots{1}, plane.number(2)-1, -%! plane.order(2)-1, plane.knots{2}, +%! pkl = surfderivcpts (plane.number(1)-1, plane.order(1)-1,... +%! plane.knots{1}, plane.number(2)-1,... +%! plane.order(2)-1, plane.knots{2}, ... %! squeeze (plane.coefs(1,:,:)), 1);
--- a/extra/nurbs/inst/surfderiveval.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/surfderiveval.m Mon Mar 22 18:13:09 2010 +0000 @@ -1,59 +1,61 @@ -%% Copyright (C) 2009 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 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, see <http://www.gnu.org/licenses/>. +function skl = surfderiveval (n, p, U, m, q, V, P, u, v, d) +% +% SURFDERIVEVAL: Compute the derivatives of a B-spline surface +% +% usage: pkl = surfderiveval (n, p, U, m, q, V, P, d) +% +% INPUT: +% n+1, m+1 = number of control points +% p, q = spline order +% U, V = knots +% P = control points +% u,v = evaluation points +% d = derivative order +% OUTPUT: +% skl (k+1, l+1) = surface differentiated k +% times in the u direction and l +% times in the v direction +% +% Adaptation of algorithm A3.8 from the NURBS book, pg115 +% +% Copyright (C) 2009 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 2 of the License, or +% (at your option) any later version. -%% usage: pkl = surfderivcpts (n, p, U, m, q, V, P, d) -%% -%% INPUT: -%% n+1, m+1 = number of control points -%% p, q = spline order -%% U, V = knots -%% P = control points -%% u,v = evaluation points -%% d = derivative order -%% OUTPUT: -%% skl (k+1, l+1) = surface differentiated k -%% times in the u direction and l -%% times in the v direction -%% -%% Adaptation of algorithm A3.8 from the NURBS book - -function skl = surfderiveval (n, p, U, m, q, V, P, u, v, d) +% 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, see <http://www.gnu.org/licenses/>. du = min (d, p); for k=p+1:d for l=0:d-k skl(k+1,l+1)=0; - endfor - endfor + end + end dv = min (d, q); for l=q+1:d for k=0:d-l skl(k+1,l+1)=0; - endfor - endfor + end + end uspan = findspan (n, p, u, U); for ip=0:p Nu(1:ip+1,ip+1) = basisfun (uspan, u, ip, U)'; - endfor + end vspan = findspan (m, q, v, V); for ip=0:q Nv(1:ip+1,ip+1) = basisfun (vspan, v, ip, V)'; - endfor + end pkl = surfderivcpts (n, p, U, m, q, V, P, d, uspan-p, uspan, ... vspan-q, vspan); @@ -63,16 +65,16 @@ for l = 0:dd skl(k+1,l+1) =0; for i=0:q-l - tmp = 0; - for j = 0:p-k - tmp += Nu(j+1,p-k+1) * pkl(k+1,l+1,j+1,i+1); - endfor - skl(k+1,l+1) += Nv(i+1,q-l+1)*tmp; - endfor - endfor - endfor + tmp = 0; + for j = 0:p-k + tmp = tmp + Nu(j+1,p-k+1) * pkl(k+1,l+1,j+1,i+1); + end + skl(k+1,l+1) = skl(k+1,l+1) + Nv(i+1,q-l+1)*tmp; + end + end + end -endfunction +end %!shared srf %!test @@ -80,21 +82,21 @@ %! c = [0 1/2 1]; %! [coef(2,:,:), coef(1,:,:)] = meshgrid (c, c); %! srf = nrbmak (coef, {k, k}); -%! skl = surfderiveval (srf.number(1)-1, -%! srf.order(1)-1, -%! srf.knots{1}, -%! srf.number(2)-1, -%! srf.order(2)-1, -%! srf.knots{2}, +%! skl = surfderiveval (srf.number(1)-1, ... +%! srf.order(1)-1, ... +%! srf.knots{1}, ... +%! srf.number(2)-1, ... +%! srf.order(2)-1, ... +%! srf.knots{2},... %! squeeze(srf.coefs(1,:,:)), .5, .5, 1) ; %! assert (skl, [.5 0; 1 0]) %!test %! srf = nrbkntins (srf, {[], rand(1,2)}); -%! skl = surfderiveval (srf.number(1)-1, -%! srf.order(1)-1, -%! srf.knots{1}, -%! srf.number(2)-1, -%! srf.order(2)-1, -%! srf.knots{2}, +%! skl = surfderiveval (srf.number(1)-1,... +%! srf.order(1)-1, ... +%! srf.knots{1},... +%! srf.number(2)-1,... +%! srf.order(2)-1, ... +%! srf.knots{2},... %! squeeze(srf.coefs(1,:,:)), .5, .5, 1) ; %! assert (skl, [.5 0; 1 0], 100*eps) \ No newline at end of file
--- a/extra/nurbs/inst/tbasisfun.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/tbasisfun.m Mon Mar 22 18:13:09 2010 +0000 @@ -1,88 +1,84 @@ -## Copyright (C) 2009 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 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 Octave; see the file COPYING. If not, see -## <http://www.gnu.org/licenses/>. +function N = tbasisfun (u, p, U) +% +% TBASISFUN: Compute a B- or T-Spline basis function from its local knot vector. +% +% usage: +% +% N = tbasisfun (u, p, U) +% N = tbasisfun ([u; v], [p q], {U, V}) +% N = tbasisfun ([u; v; w], [p q r], {U, V, W}) +% +% INPUT: +% u or [u; v] : points in parameter space where the basis function is to be +% evaluated +% +% U or {U, V} : local knot vector +% +% p or [p q] : polynomial order of the bais function +% +% OUTPUT: +% N : basis function evaluated at the given parametric points +% +% Copyright (C) 2009 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 2 of the License, or +% (at your option) any later version. -## TBASISFUN: Compute a B- or T-Spline basis function from its local knot vector. -## -## usage: -## -## N = tbasisfun (u, p, U) -## N = tbasisfun ([u; v], [p q], {U, V}) -## -## INPUT: -## u or [u; v] : points in parameter space where the basis function is to be -## evaluated -## -## U or {U, V} : local knot vector -## -## p or [p q] : polynomial order of the bais function -## -## OUTPUT: -## N : basis function evaluated at the given parametric points - -## Author: Carlo de Falco <cdf _AT_ users.sourceforge.net> -## Created: 2009-06-25 - -function N = tbasisfun (u, p, U) +% 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, see <http://www.gnu.org/licenses/>. - if (! iscell (U)) + if (~ iscell (U)) U = sort (U); assert (numel (U) == p+2) + N = zeros(1,numel(u)); for ii=1:numel(u) N(ii) = onebasisfun__ (u(ii), p, U); - endfor + end - else - + elseif size(U,2) == 2 + U{1} = sort(U{1}); U{2} = sort(U{2}); + assert (numel(U{1}) == p(1) && numel(U{2}) == p(2)) + + Nu = zeros(1,numel(u)); Nv = zeros(1,numel(u)); for ii=1:numel(u(1,:)) Nu(ii) = onebasisfun__ (u(1,ii), p(1), U{1}); - endfor + end for ii=1:numel(u(1,:)) Nv(ii) = onebasisfun__ (u(2,ii), p(2), U{2}); - endfor + end N = Nu.*Nv; - endif -endfunction - -function N = onebasisfun__ (u, p, U) - - N = 0; - if (! any (U <= u)) | (! any (U >= u)) - return; - elseif (p == 0) - N = 1; - return; - endif + elseif size(U,2) == 3 + U{1} = sort(U{1}); U{2} = sort(U{2}); U{3} = sort(U{3}); + assert (numel(U{1}) == p(1) && numel(U{2}) == p(2) && numel(U{3}) == p(3)) + + Nu = zeros(1,numel(u)); Nv = zeros(1,numel(u)); Nw = zeros(1,numel(u)); + for ii=1:numel(u(1,:)) + Nu(ii) = onebasisfun__ (u(1,ii), p(1), U{1}); + end - ln = u - U(1); - ld = U(end-1) - U(1); - if (ld != 0) - N += ln * onebasisfun__ (u, p-1, U(1:end-1))/ ld; - endif + for ii=1:numel(u(1,:)) + Nv(ii) = onebasisfun__ (u(2,ii), p(2), U{2}); + end - dn = U(end) - u; - dd = U(end) - U(2); - if (dd != 0) - N += dn * onebasisfun__ (u, p-1, U(2:end))/ dd; - endif - -endfunction + for ii=1:numel(u(1,:)) + Nw(ii) = onebasisfun__ (u(3,ii), p(3), U{3}); + end + + N = Nu.*Nv.*Nw; + end + +end %!demo %! U = {[0 0 1/2 1 1], [0 0 0 1 1]}; @@ -91,3 +87,5 @@ %! u = [X(:), Y(:)]'; %! N = tbasisfun (u, p, U); %! surf (X, Y, reshape (N, size(X))) +%! title('Basis function associated to a local knot vector') +%! hold off \ No newline at end of file
--- a/extra/nurbs/inst/vecangle.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/vecangle.m Mon Mar 22 18:13:09 2010 +0000 @@ -8,10 +8,12 @@ % % ang = vecmag2(num,dum) % -% Parameters: +% INPUT: % % num : Numerator, vector of size (1,nv). % dem : Denominator, vector of size (1,nv). +% +% OUTPUT: % ang : Arctangents, row vector of angles. % % Description: @@ -25,11 +27,24 @@ % Find the atan(1.2,2.0) and atan(1.5,3.4) using vecangle % % ang = vecangle([1.2 1.5], [2.0 3.4]); +% +% Copyright (C) 2000 Mark Spink +% +% 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. -% D.M. Spink -% Copyright (c) 2000. +% 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, see <http://www.gnu.org/licenses/>. ang = atan2(num,den); index = find(ang < 0.0); ang(index) = 2*pi+ang(index); +end
--- a/extra/nurbs/inst/veccross.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/veccross.m Mon Mar 22 18:13:09 2010 +0000 @@ -6,11 +6,13 @@ % % cross = veccross(vec1,vec2); % -% Parameters: +% INPUT: % % vec1 : An array of column vectors represented by a matrix of % vec2 size (dim,nv), where is the dimension of the vector and % nv the number of vectors. +% +% OUTPUT: % % cross : Array of column vectors, each element is corresponding % to the cross product of the respective components in vec1 @@ -27,9 +29,21 @@ % (5.1,0.0,2.3) and (2.5,3.2,4.0) % % cross = veccross([2.3 5.1; 3.4 0.0; 5.6 2.3],[1.2 2.5; 4.5 3.2; 1.2 4.0]); +% +% Copyright (C) 2000 Mark Spink +% +% 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. -% D.M. Spink -% Copyright (c) 2000. +% 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, see <http://www.gnu.org/licenses/>. if size(vec1,1) == 2 % 2D vector @@ -42,3 +56,4 @@ vec1(1,:).*vec2(2,:)-vec1(2,:).*vec2(1,:)]; end +end
--- a/extra/nurbs/inst/vecdot.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/vecdot.m Mon Mar 22 18:13:09 2010 +0000 @@ -6,12 +6,14 @@ % % dot = vecdot(vec1,vec2); % -% Parameters: +% INPUT: % % vec1 : An array of column vectors represented by a matrix of % vec2 size (dim,nv), where is the dimension of the vector and % nv the number of vectors. -% +% +% OUTPUT: +% % dot : Row vector of scalars, each element corresponding to % the dot product of the respective components in vec1 and % vec2. @@ -27,9 +29,22 @@ % (5.1,0.0,2.3) and (2.5,3.2,4.0) % % dot = vecdot([2.3 5.1; 3.4 0.0; 5.6 2.3],[1.2 2.5; 4.5 3.2; 1.2 4.0]); +% +% Copyright (C) 2000 Mark Spink +% +% 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. -% D.M. Spink -% Copyright (c) 2000 +% 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, see <http://www.gnu.org/licenses/>. dot = sum(vec1.*vec2); +end
--- a/extra/nurbs/inst/vecmag.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/vecmag.m Mon Mar 22 18:13:09 2010 +0000 @@ -6,12 +6,14 @@ % % mvec = vecmag(vec) % -% Parameters: +% INPUT: % % vec : An array of column vectors represented by a matrix of % size (dim,nv), where is the dimension of the vector and % nv the number of vectors. -% +% +% OUTPUT: +% % mvec : Magnitude of the vectors, vector of size (1,nv). % % Description: @@ -23,9 +25,22 @@ % Find the magnitude of the two vectors (0.0,2.0,1.3) and (1.5,3.4,2.3) % % mvec = vecmag([0.0 1.5; 2.0 3.4; 1.3 2.3]); +% +% Copyright (C) 2000 Mark Spink +% +% 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. -% D.M. Spink -% Copyright (c) 2000. +% 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, see <http://www.gnu.org/licenses/>. mag = sqrt(sum(vec.^2)); +end
--- a/extra/nurbs/inst/vecmag2.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/vecmag2.m Mon Mar 22 18:13:09 2010 +0000 @@ -6,11 +6,13 @@ % % mvec = vecmag2(vec) % -% Parameters: +% INPUT: % % vec : An array of column vectors represented by a matrix of % size (dim,nv), where dim is the dimension of the vector and % nv the number of vectors. +% +% OUTPUT: % % mvec : Squared magnitude of the vectors, vector of size (1,nv). % @@ -24,9 +26,22 @@ % and (1.5,3.4,2.3) % % mvec = vecmag2([0.0 1.5; 2.0 3.4; 1.3 2.3]); +% +% Copyright (C) 2000 Mark Spink +% +% 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. -% D.M. Spink -% Copyright (c) 2000. +% 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, see <http://www.gnu.org/licenses/>. mag = sum(vec.^2); +end
--- a/extra/nurbs/inst/vecnorm.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/vecnorm.m Mon Mar 22 18:13:09 2010 +0000 @@ -6,11 +6,13 @@ % % nvec = vecnorn(vec); % -% Parameters: +% INPUT: % % vec : An array of column vectors represented by a matrix of % size (dim,nv), where is the dimension of the vector and % nv the number of vectors. +% +% OUTPUT: % % nvec : Normalised vectors, matrix the smae size as vec. % @@ -23,9 +25,22 @@ % Normalise the two vectors (0.0,2.0,1.3) and (1.5,3.4,2.3) % % nvec = vecnorm([0.0 1.5; 2.0 3.4; 1.3 2.3]); +% +% Copyright (C) 2000 Mark Spink +% +% 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. -% D.M. Spink -% Copyright (c) 2000. +% 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, see <http://www.gnu.org/licenses/>. nvec = vec./repmat(sqrt(sum(vec.^2)),[size(vec,1) ones(1,ndims(vec)-1)]); +end
--- a/extra/nurbs/inst/vecrotx.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/vecrotx.m Mon Mar 22 18:13:09 2010 +0000 @@ -6,9 +6,11 @@ % % rx = vecrotx(angle); % -% Parameters: +% INPUT: % % angle : rotation angle defined in radians +% +% OUTPUT: % % rx : (4x4) Transformation matrix. % @@ -34,13 +36,27 @@ % trans = vecrotx(%pi/4); % rline = nrbtform(line, trans); % -% See: +% See also: % % nrbtform +% +% Copyright (C) 2000 Mark Spink +% +% 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. -% Dr D.M. Spink -% Copyright (c) 2000. +% 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, see <http://www.gnu.org/licenses/>. sn = sin(angle); cn = cos(angle); rx = [1 0 0 0; 0 cn -sn 0; 0 sn cn 0; 0 0 0 1]; + +end
--- a/extra/nurbs/inst/vecroty.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/vecroty.m Mon Mar 22 18:13:09 2010 +0000 @@ -6,10 +6,12 @@ % % ry = vecroty(angle); % -% Parameters: +% INPUT: % % angle : rotation angle defined in radians % +% OUTPUT: +% % ry : (4x4) Transformation matrix. % % @@ -34,13 +36,27 @@ % trans = vecroty(%pi/4); % rline = nrbtform(line, trans); % -% See: +% See also: % % nrbtform +% +% Copyright (C) 2000 Mark Spink +% +% 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. -% Dr D.M. Spink -% Copyright (c) 2000. +% 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, see <http://www.gnu.org/licenses/>. sn = sin(angle); cn = cos(angle); ry = [cn 0 sn 0; 0 1 0 0; -sn 0 cn 0; 0 0 0 1]; + +end
--- a/extra/nurbs/inst/vecrotz.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/vecrotz.m Mon Mar 22 18:13:09 2010 +0000 @@ -6,9 +6,11 @@ % % rz = vecrotz(angle); % -% Parameters: +% INPUT: % % angle : rotation angle defined in radians +% +% OUTPUT: % % rz : (4x4) Transformation matrix. % @@ -34,13 +36,27 @@ % trans = vecrotz(%pi/4); % rline = nrbtform(line, trans); % -% See: +% See also: % % nrbtform +% +% Copyright (C) 2000 Mark Spink +% +% 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. -% Dr D.M. Spink -% Copyright (c) 2000. +% 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, see <http://www.gnu.org/licenses/>. sn = sin(angle); cn = cos(angle); rz = [cn -sn 0 0; sn cn 0 0; 0 0 1 0; 0 0 0 1]; + +end
--- a/extra/nurbs/inst/vecscale.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/vecscale.m Mon Mar 22 18:13:09 2010 +0000 @@ -1,18 +1,3 @@ -%% Copyright (C) 2003 Mark Spink, 2007 Daniel Claxton -%% -%% 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, see <http://www.gnu.org/licenses/>. - function ss = vecscale(vector) % @@ -22,10 +7,12 @@ % % ss = vecscale(svec) % -% Parameters: +% INPUT: % % svec : A vectors defining the scaling along the x,y and z axes. % i.e. [sx, sy, sy] +% +% OUTPUT: % % ss : Scaling Transformation Matrix % @@ -49,12 +36,24 @@ % trans = vecscale([3.0 2.0 4.0]); % sline = nrbtform(line, trans); % -% See: +% See also: % % nrbtform +% +% Copyright (C) 2000 Mark Spink, 2007 Daniel Claxton +% +% 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. -% Dr D.M. Spink -% Copyright (c) 2000. +% 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, see <http://www.gnu.org/licenses/>. if nargin < 1 error('Scaling vector not specified'); @@ -62,3 +61,5 @@ s = [vector(:);0;0]; ss = [s(1) 0 0 0; 0 s(2) 0 0; 0 0 s(3) 0; 0 0 0 1]; + +end
--- a/extra/nurbs/inst/vectrans.m Sun Mar 21 18:40:41 2010 +0000 +++ b/extra/nurbs/inst/vectrans.m Mon Mar 22 18:13:09 2010 +0000 @@ -6,10 +6,12 @@ % % st = vectrans(tvec) % -% Parameters: +% INPUT: % % tvec : A vectors defining the translation along the x,y and % z axes. i.e. [tx, ty, ty] +% +% OUTPUT: % % st : Translation Transformation Matrix % @@ -33,12 +35,25 @@ % trans = vectrans([3.0 2.0 4.0]); % tline = nrbtform(line, trans); % -% See: +% See also: % % nrbtform +% +% Copyright (C) 2000 Mark Spink +% +% 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. -% Dr D.M. Spink -% Copyright (c) 2000. +% 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, see <http://www.gnu.org/licenses/>. + if nargin < 1 error('Translation vector required'); @@ -46,3 +61,5 @@ v = [vector(:);0;0]; dd = [1 0 0 v(1); 0 1 0 v(2); 0 0 1 v(3); 0 0 0 1]; + +end