1function [p,w] = nrbeval(nurbs,tt) 2% 3% NRBEVAL: Evaluate a NURBS at parametric points. 4% 5% Calling Sequences: 6% 7% [p,w] = nrbeval(crv,ut) 8% [p,w] = nrbeval(srf,{ut,vt}) 9% [p,w] = nrbeval(vol,{ut,vt,wt}) 10% [p,w] = nrbeval(srf,pts) 11% 12% INPUT: 13% 14% crv : NURBS curve, see nrbmak. 15% 16% srf : NURBS surface, see nrbmak. 17% 18% vol : NURBS volume, see nrbmak. 19% 20% ut : Parametric evaluation points along U direction. 21% 22% vt : Parametric evaluation points along V direction. 23% 24% wt : Parametric evaluation points along W direction. 25% 26% pts : Array of scattered points in parametric domain 27% 28% OUTPUT: 29% 30% p : Evaluated points on the NURBS curve, surface or volume as 31% Cartesian coordinates (x,y,z). If w is included on the lhs argument 32% list the points are returned as homogeneous coordinates (wx,wy,wz). 33% 34% w : Weights of the homogeneous coordinates of the evaluated 35% points. Note inclusion of this argument changes the type 36% of coordinates returned in p (see above). 37% 38% Description: 39% 40% Evaluation of NURBS curves, surfaces or volume at parametric points along 41% the U, V and W directions. Either homogeneous coordinates are returned 42% if the weights are requested in the lhs arguments, or as Cartesian coordinates. 43% This function utilises the 'C' interface bspeval. 44% 45% Examples: 46% 47% Evaluate the NURBS circle at twenty points from 0.0 to 1.0 48% 49% nrb = nrbcirc; 50% ut = linspace(0.0,1.0,20); 51% p = nrbeval(nrb,ut); 52% 53% See also: 54% 55% bspeval 56% 57% Copyright (C) 2000 Mark Spink 58% Copyright (C) 2010 Carlo de Falco 59% Copyright (C) 2010, 2011, 2015 Rafael Vazquez 60% 61% This program is free software: you can redistribute it and/or modify 62% it under the terms of the GNU General Public License as published by 63% the Free Software Foundation, either version 3 of the License, or 64% (at your option) any later version. 65 66% This program is distributed in the hope that it will be useful, 67% but WITHOUT ANY WARRANTY; without even the implied warranty of 68% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 69% GNU General Public License for more details. 70% 71% You should have received a copy of the GNU General Public License 72% along with this program. If not, see <http://www.gnu.org/licenses/>. 73 74if (nargin < 2) 75 error('Not enough input arguments'); 76end 77 78foption = 1; % output format 3D cartesian coordinates 79if (nargout == 2) 80 foption = 0; % output format 4D homogenous coordinates 81end 82 83if (~isstruct(nurbs)) 84 error('NURBS representation is not structure!'); 85end 86 87if (~strcmp(nurbs.form,'B-NURBS')) 88 error('Not a recognised NURBS representation'); 89end 90 91if (iscell(nurbs.knots)) 92 if (size(nurbs.knots,2) == 3) 93 %% NURBS structure represents a volume 94 95 num1 = nurbs.number(1); 96 num2 = nurbs.number(2); 97 num3 = nurbs.number(3); 98 degree = nurbs.order-1; 99 100 if (iscell(tt)) 101 nt1 = numel (tt{1}); 102 nt2 = numel (tt{2}); 103 nt3 = numel (tt{3}); 104 105 %% evaluate along the w direction 106 val = reshape (nurbs.coefs, 4*num1*num2, num3); 107 val = bspeval (degree(3), val, nurbs.knots{3}, tt{3}); 108 val = reshape (val, [4 num1 num2 nt3]); 109 110 %% Evaluate along the v direction 111 val = permute (val, [1 2 4 3]); 112 val = reshape (val, 4*num1*nt3, num2); 113 val = bspeval (degree(2), val, nurbs.knots{2}, tt{2}); 114 val = reshape (val, [4 num1 nt3 nt2]); 115 val = permute (val, [1 2 4 3]); 116 117 %% Evaluate along the u direction 118 val = permute (val, [1 3 4 2]); 119 val = reshape (val, 4*nt2*nt3, num1); 120 val = bspeval (degree(1), val, nurbs.knots{1}, tt{1}); 121 val = reshape (val, [4 nt2 nt3 nt1]); 122 val = permute (val, [1 4 2 3]); 123 pnts = val; 124 125 p = pnts(1:3,:,:,:); 126 w = pnts(4,:,:,:); 127 if (foption) 128 p = p./repmat(w,[3 1 1 1]); 129 end 130 131 else 132 133 %% Evaluate at scattered points 134 %% tt(1,:) represents the u direction 135 %% tt(2,:) represents the v direction 136 %% tt(3,:) represents the w direction 137 138 st = size(tt); 139 if (st(1) ~= 3 && st(2) == 3 && numel(st) == 2) 140 tt = tt'; 141 st = size (tt); 142 end 143 nt = prod(st(2:end)); 144 145 tt = reshape (tt, [3, nt]); 146 147 %% evaluate along the w direction 148 val = reshape(nurbs.coefs,4*num1*num2,num3); 149 val = bspeval(degree(3),val,nurbs.knots{3},tt(3,:)); 150 val = reshape(val,[4 num1 num2 nt]); 151 152 %% evaluate along the v direction 153 val2 = zeros(4*num1,nt); 154 for v = 1:nt 155 coefs = reshape(val(:,:,:,v),4*num1,num2); 156 val2(:,v) = bspeval(degree(2),coefs,nurbs.knots{2},tt(2,v)); 157 end 158 val2 = reshape(val2,[4 num1 nt]); 159 160 %% evaluate along the u direction 161 pnts = zeros(4,nt); 162 for v = 1:nt 163 coefs = reshape (val2(:,:,v), [4 num1]); 164 pnts(:,v) = bspeval(degree(1),coefs,nurbs.knots{1},tt(1,v)); 165 end 166 167 w = pnts(4,:); 168 p = pnts(1:3,:); 169 if (foption) 170 p = p./repmat(w,[3, 1]); 171 end 172 173 if (numel(st) ~= 2) 174 w = reshape (w, [st(2:end)]); 175 p = reshape (p, [3, st(2:end)]); 176 end 177 end 178 179 elseif (size(nurbs.knots,2) == 2) 180 %% NURBS structure represents a surface 181 182 num1 = nurbs.number(1); 183 num2 = nurbs.number(2); 184 degree = nurbs.order-1; 185 186 if (iscell(tt)) 187 %% Evaluate over a [u,v] grid 188 %% tt{1} represents the u direction 189 %% tt{2} represents the v direction 190 191 nt1 = length(tt{1}); 192 nt2 = length(tt{2}); 193 194 %% Evaluate along the v direction 195 val = reshape(nurbs.coefs,4*num1,num2); 196 val = bspeval(degree(2),val,nurbs.knots{2},tt{2}); 197 val = reshape(val,[4 num1 nt2]); 198 199 %% Evaluate along the u direction 200 val = permute(val,[1 3 2]); 201 val = reshape(val,4*nt2,num1); 202 val = bspeval(degree(1),val,nurbs.knots{1},tt{1}); 203 val = reshape(val,[4 nt2 nt1]); 204 val = permute(val,[1 3 2]); 205 206 w = val(4,:,:); 207 p = val(1:3,:,:); 208 if (foption) 209 p = p./repmat(w,[3 1 1]); 210 end 211 212 else 213 214 %% Evaluate at scattered points 215 %% tt(1,:) represents the u direction 216 %% tt(2,:) represents the v direction 217 218 st = size(tt); 219 if (st(1) ~= 2 && st(2) == 2 && numel(st) == 2) 220 tt = tt'; 221 st = size (tt); 222 end 223 nt = prod(st(2:end)); 224 225 tt = reshape (tt, [2, nt]); 226 227 val = reshape(nurbs.coefs,4*num1,num2); 228 val = bspeval(degree(2),val,nurbs.knots{2},tt(2,:)); 229 val = reshape(val,[4 num1 nt]); 230 231 232 %% evaluate along the u direction 233 pnts = zeros(4,nt); 234 for v = 1:nt 235 coefs = reshape (val(:,:,v), [4 num1]); 236 pnts(:,v) = bspeval(degree(1),coefs,nurbs.knots{1},tt(1,v)); 237 end 238 239 w = pnts(4,:); 240 p = pnts(1:3,:); 241 if (foption) 242 p = p./repmat(w,[3, 1]); 243 end 244 245 if (numel(st) ~= 2) 246 w = reshape (w, [st(2:end)]); 247 p = reshape (p, [3, st(2:end)]); 248 end 249 250 end 251 252 end 253else 254 255 %% NURBS structure represents a curve 256 %% tt represent a vector of parametric points in the u direction 257 258 if (iscell (tt) && numel (tt) == 1) 259 tt = cell2mat (tt); 260 end 261 262 st = size (tt); 263 264 val = bspeval(nurbs.order-1,nurbs.coefs,nurbs.knots,tt(:)'); 265 266 w = val(4,:); 267 p = val(1:3,:); 268 if foption 269 p = p./repmat(w,3,1); 270 end 271 272 if (st(1) ~= 1 || numel(st) ~= 2) 273 w = reshape (w, st); 274 p = reshape (p, [3, st]); 275 end 276 277end 278 279end 280 281%!demo 282%! srf = nrbtestsrf; 283%! p = nrbeval(srf,{linspace(0.0,1.0,20) linspace(0.0,1.0,20)}); 284%! h = surf(squeeze(p(1,:,:)),squeeze(p(2,:,:)),squeeze(p(3,:,:))); 285%! title('Test surface.'); 286%! hold off 287 288%!test 289%! knots{1} = [0 0 0 1 1 1]; 290%! knots{2} = [0 0 0 .5 1 1 1]; 291%! knots{3} = [0 0 0 0 1 1 1 1]; 292%! cx = [0 0.5 1]; nx = length(cx); 293%! cy = [0 0.25 0.75 1]; ny = length(cy); 294%! cz = [0 1/3 2/3 1]; nz = length(cz); 295%! coefs(1,:,:,:) = repmat(reshape(cx,nx,1,1),[1 ny nz]); 296%! coefs(2,:,:,:) = repmat(reshape(cy,1,ny,1),[nx 1 nz]); 297%! coefs(3,:,:,:) = repmat(reshape(cz,1,1,nz),[nx ny 1]); 298%! coefs(4,:,:,:) = 1; 299%! nurbs = nrbmak(coefs, knots); 300%! x = rand(5,1); y = rand(5,1); z = rand(5,1); 301%! tt = [x y z]'; 302%! points = nrbeval(nurbs,tt); 303%! 304%! assert(points,tt,1e-10) 305%! 306%!test 307%! knots{1} = [0 0 0 1 1 1]; 308%! knots{2} = [0 0 0 0 1 1 1 1]; 309%! knots{3} = [0 0 1 1]; 310%! cx = [0 0 1]; nx = length(cx); 311%! cy = [0 0 0 1]; ny = length(cy); 312%! cz = [0 1]; nz = length(cz); 313%! coefs(1,:,:,:) = repmat(reshape(cx,nx,1,1),[1 ny nz]); 314%! coefs(2,:,:,:) = repmat(reshape(cy,1,ny,1),[nx 1 nz]); 315%! coefs(3,:,:,:) = repmat(reshape(cz,1,1,nz),[nx ny 1]); 316%! coefs(4,:,:,:) = 1; 317%! nurbs = nrbmak(coefs, knots); 318%! x = rand(5,1); y = rand(5,1); z = rand(5,1); 319%! tt = [x y z]'; 320%! points = nrbeval(nurbs,tt); 321%! assert(points,[x.^2 y.^3 z]',1e-10); 322%! 323%!test 324%! knots{1} = [0 0 0 1 1 1]; 325%! knots{2} = [0 0 0 0 1 1 1 1]; 326%! knots{3} = [0 0 1 1]; 327%! cx = [0 0 1]; nx = length(cx); 328%! cy = [0 0 0 1]; ny = length(cy); 329%! cz = [0 1]; nz = length(cz); 330%! coefs(1,:,:,:) = repmat(reshape(cx,nx,1,1),[1 ny nz]); 331%! coefs(2,:,:,:) = repmat(reshape(cy,1,ny,1),[nx 1 nz]); 332%! coefs(3,:,:,:) = repmat(reshape(cz,1,1,nz),[nx ny 1]); 333%! coefs(4,:,:,:) = 1; 334%! coefs = coefs([2 1 3 4],:,:,:); 335%! nurbs = nrbmak(coefs, knots); 336%! x = rand(5,1); y = rand(5,1); z = rand(5,1); 337%! tt = [x y z]'; 338%! points = nrbeval(nurbs,tt); 339%! [y.^3 x.^2 z]'; 340%! assert(points,[y.^3 x.^2 z]',1e-10); 341