nurbs_曲面拟合_NURBS曲面拟合_拟合_

function [ic,ik] = bspdegelev(d,c,k,t) % % Function Name: % % bspdegevel - Degree elevate a univariate B-Spline. % % Calling Sequence: % % [ic,ik] = bspdegelev(d,c,k,t) % % Parameters: % % 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. % % ic : Control points of the new B-Spline. % % ik : Knot vector of the new B-Spline. % % Description: % % Degree elevate a univariate B-Spline. This function provides an % interface to a toolbox 'C' routine. [mc,nc] = size(c); % % int bspdegelev(int d, double *c, int mc, int nc, double *k, int nk, % int t, int *nh, double *ic, double *ik) % { % int row,col; % % int ierr = 0; % int i, j, q, s, m, ph, ph2, mpi, mh, r, a, b, cind, oldr, mul; % int n, lbz, rbz, save, tr, kj, first, kind, last, bet, ii; % double inv, ua, ub, numer, den, alf, gam; % double **bezalfs, **bpts, **ebpts, **Nextbpts, *alfs; % % double **ctrl = vec2mat(c, mc, nc); % ic = zeros(mc,nc*(t)); % double **ictrl = vec2mat(ic, mc, nc*(t+1)); % n = nc - 1; % n = nc - 1; % bezalfs = zeros(d+1,d+t+1); % bezalfs = matrix(d+1,d+t+1); bpts = zeros(mc,d+1); % bpts = matrix(mc,d+1); ebpts = zeros(mc,d+t+1); % ebpts = matrix(mc,d+t+1); Nextbpts = zeros(mc,d+1); % Nextbpts = matrix(mc,d+1); alfs = zeros(d,1); % alfs = (double *) mxMalloc(d*sizeof(double)); % m = n + d + 1; % m = n + d + 1; ph = d + t; % ph = d + t; ph2 = floor(ph / 2); % ph2 = ph / 2; % % // compute bezier degree elevation coefficeients bezalfs(1,1) = 1; % bezalfs[0][0] = bezalfs[ph][d] = 1.0; bezalfs(d+1,ph+1) = 1; % for i=1:ph2 % for (i = 1; i <= ph2; i++) { inv = 1/bincoeff(ph,i); % inv = 1.0 / bincoeff(ph,i); mpi = min(d,i); % mpi = min(d,i); % for j=max(0,i-t):mpi % for (j = max(0,i-t); j <= mpi; j++) bezalfs(j+1,i+1) = inv*bincoeff(d,j)*bincoeff(t,i-j); % bezalfs[i][j] = inv * bincoeff(d,j) * bincoeff(t,i-j); end end % } % for i=ph2+1:ph-1 % for (i = ph2+1; i <= ph-1; i++) { mpi = min(d,i); % mpi = min(d, i); for j=max(0,i-t):mpi % for (j = max(0,i-t); j <= mpi; j++) bezalfs(j+1,i+1) = bezalfs(d-j+1,ph-i+1); % bezalfs[i][j] = bezalfs[ph-i][d-j]; end end % } % mh = ph; % mh = ph; kind = ph+1; % kind = ph+1; r = -1; % r = -1; a = d; % a = d; b = d+1; % b = d+1; cind = 1; % cind = 1; ua = k(1); % ua = k[0]; % for ii=0:mc-1 % for (ii = 0; ii < mc; ii++) ic(ii+1,1) = c(ii+1,1); % ictrl[0][ii] = ctrl[0][ii]; end % for i=0:ph % for (i = 0; i <= ph; i++) ik(i+1) = ua; % ik[i] = ua; end % % // initialise first bezier seg for i=0:d % for (i = 0; i <= d; i++) for ii=0:mc-1 % for (ii = 0; ii < mc; ii++) bpts(ii+1,i+1) = c(ii+1,i+1); % bpts[i][ii] = ctrl[i][ii]; end end % % // big loop thru knot vector while b < m % while (b < m) { i = b; % i = b; while b < m && k(b+1) == k(b+2) % while (b < m && k[b] == k[b+1]) b = b + 1; % b++; end % mul = b - i + 1; % mul = b - i + 1; mh = mh + mul + t; % mh += mul + t; ub = k(b+1); % ub = k[b]; oldr = r; % oldr = r; r = d - mul; % r = d - mul; % % // insert knot u(b) r times if oldr > 0 % if (oldr > 0) lbz = floor((oldr+2)/2); % lbz = (oldr+2) / 2; else % else lbz = 1; % lbz = 1; end % if r > 0 % if (r > 0) rbz = ph - floor((r+1)/2); % rbz = ph - (r+1)/2; else % else rbz = ph; % rbz = ph; end