FreeWRL/FreeX3D  3.0.0
Component_NURBS.c
1 
2 /****************************************************************************
3  This file is part of the FreeWRL/FreeX3D Distribution.
4 
5  Copyright 2009 CRC Canada. (http://www.crc.gc.ca)
6 
7  FreeWRL/FreeX3D is free software: you can redistribute it and/or modify
8  it under the terms of the GNU Lesser Public License as published by
9  the Free Software Foundation, either version 3 of the License, or
10  (at your option) any later version.
11 
12  FreeWRL/FreeX3D is distributed in the hope that it will be useful,
13  but WITHOUT ANY WARRANTY; without even the implied warranty of
14  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15  GNU General Public License for more details.
16 
17  You should have received a copy of the GNU General Public License
18  along with FreeWRL/FreeX3D. If not, see <http://www.gnu.org/licenses/>.
19 ****************************************************************************/
20 
21 
22 /*******************************************************************
23 
24  X3D NURBS Component
25 
26 *********************************************************************/
27 
28 #include <config.h>
29 #include <system.h>
30 #include <display.h>
31 #include <internal.h>
32 
33 #include <libFreeWRL.h>
34 
35 #include "../vrml_parser/Structs.h"
36 #include "../main/headers.h"
37 
38 #include "Collision.h"
39 #include "LinearAlgebra.h"
40 #include "../opengl/Frustum.h"
41 #include "../opengl/Material.h"
42 #include "Component_Geometry3D.h"
43 #include "../opengl/OpenGL_Utils.h"
44 #include "../opengl/Textures.h"
45 
46 #include "Component_Shape.h"
47 #include "../scenegraph/RenderFuncs.h"
48 #include "../vrml_parser/CRoutes.h"
49 #include "Polyrep.h"
50 #include <float.h>
51 #if defined(_MSC_VER) && _MSC_VER < 1500
52 #define cosf cos
53 #define sinf sin
54 #endif
55 
56 typedef struct pComponent_NURBS{
57  void *nada;// = 0;
58 
60 void *Component_NURBS_constructor(){
61  void *v = MALLOCV(sizeof(struct pComponent_NURBS));
62  memset(v,0,sizeof(struct pComponent_NURBS));
63  return v;
64 }
65 void Component_NURBS_init(struct tComponent_NURBS *t){
66  //public
67  //private
68  t->prv = Component_NURBS_constructor();
69  {
71  p->nada = NULL;
72  }
73 }
74 //ppComponent_NURBS p = (ppComponent_NURBS)gglobal()->Component_NURBS.prv;
75 
76 /*
77 Specs:
78 http://www.web3d.org/documents/specifications/19775-1/V3.3/Part01/components/nurbs.html
79 Examples:
80 http://www.web3d.org/x3d/content/examples/Basic/NURBS/
81 
82 Book: The NURBS Book
83 Piegl, Les and Tiller, Wayne; The NURBS Book, 2nd Edition, Springer-Verlag (Berlin), 1997, ISBN: 3-540-61545-8.
84 - about $120 new, softcover or kindle
85 - some university libraries have it, you don't absolutely need it, the swept and swung ideas are in it
86 
87 Nov 2014 - dug9 did some nurbs using libnurbs2, and following The Nurbs Book
88  - didn't finish the component
89  - personal notes say Nurbs Level 1 done
90 July 2016
91 Not sure what's done and what's not, nothing indicaded on Conformace page:
92 http://freewrl.sourceforge.net/conformance.html
93  code
94 Level 1 Conformance page review
95 CoordinateDouble (not mentioned) DONE
96 NurbsCurve Not Implemented DONE
97 NurbsOrientationInterpolatorNot Implemented stub
98 NurbsPatchSurface Not Implemented DONE
99 NurbsPositionInterpolator Not Implemented stub
100 NurbsSurfaceInterpolator Not Implemented stub
101 NurbsTextureCoordinate Not Implemented -
102 Level 2
103 NurbsSet Not Implemented -
104 Level 3
105 NurbsCurve2D Not Implemented DONE
106 ContourPolyline2D Not Implemented DONE
107 NurbsSweptSurface Not Implemented -
108 NurbsSwungSurface Not Implemented -
109 Level 4
110 Contour2D Not Implemented DONE
111 NurbsTrimmedSurface Not Implemented DONE?
112 
113 Perl: all the above including CoordinateDouble are in perl
114 
115 Analysis
116 dug9, 2016:
117 Chapter 12 of Redbook has lots of glu nurbs functions.
118 freewrl doesn't use opengl glu.
119 There is some source for it:
120  http://mesa3d.org/
121  ftp://ftp.freedesktop.org/pub/mesa/glu/
122  http://oss.sgi.com/projects/ogl-sample/
123 So we're left to re-implement the hard way.
124 I find nurbs libs are always disappointing in documentation.
125 I think that's because there's not a lot to nurbs, mostly plumbing and little meat:
126 - Each tesselation point Ci is computed as a weighted blend of control points: Ci = f(Pn,wn)
127 - There are a few blending functions: linear, quadratic, cubic - in theory could go higher
128 - there's a trick for doing sharp corners: knots (repeating a control point)
129 - hard part: adaptive triangulation, might need callback plumbing, libnurbs2 has some of that
130 Did I make a mistake trying to use libnurbs in 2014? Is it too much lib plumbing for too little meat?
131 Here's a smaller MIT .c
132 https://github.com/retuxx/tinyspline
133 Or you can just read a bit more and manually implement nurbs.
134 
135 
136 HARD PARTS:
137 1. nurbsSet:
138  Q1.how smooth between nurbs 'tiles'? C0, C1, C2 continuity?
139  Q1a. how to find matching/joining/adjacent edges from 2+ patches?
140  Q1b. how to use adjoinging patch info to modify/blend on/near/in-overalap of the join?
141 
142 2. trimmed surfaces - procedure?
143  Hint: however we triangulate font glyphs in Component_Text would be a good example,
144  - might be able to borrow code from it
145  Terminology:
146  Tesselate: Nurbs term for interpolate/blend a given target interpolation point
147  Triangulate: take cloud of points, and join them with edges to make surface triangles
148 
149  Fuzzy guess algo:
150  a) tesselate 2D curve in uv space of surface => A2dtess
151  b) tesselate 3D surface tesselation grid by itself => Btess
152  [option skip points inside A2dtess here instead of c)]
153  c) use point-in-poly test to remove Btess surface tesselation points
154  inside A2dtess tesselated curve: Btess - inside(A2dtess) => Ctess
155  d) go over A2dtess 2D curve point by point, interpolating/blending/tesselating in 3D like they
156  were surface tesselation points A2dtess => Dtess
157  e) add Dtess points to Ctess points => Etess
158  f) triangulate Etess => Etris
159  g) follow A2dtess points in Etris, swapping triangles so triangle edges are along A2dtess
160  h) do point-in-poly of Etris 2D triangle centroids, using A2dtess polygon,
161  to remove triangles inside A2dtess => Ftris
162  Ftris should be a trimmed surface
163 
164 3. Adaptive triangulation
165  a regular grid in uv space does not transform into a uniform grid in xyz space
166  with nurbs. And so there's a lot of plumbing in glu and other libs for
167  refining the xyz grid a) in xyz space or b) in screenspace (pixels), using some
168  measure of error and/or chordlength.
169 
170  http://people.scs.carleton.ca/~c_shu/Publications/stmalo99.pdf
171  - Triangulating Trimmed NURBS Surfaces
172  http://www.saccade.com/writing/graphics/RE-PARAM.PDF
173  - Arc Length Parameterization of Spline Curves
174  - uses yet-another nurbs curve to interpolate XY space regular length chords
175 
176 
177 */
178 void free_polyrep(struct X3D_PolyRep *rep){
179  //see also delete_polyrep - did dug9 duplicate the function or is it different?
180  if(rep){
181  rep->ntri = 0;
182  rep->transparency = 0;
183  //Q. are any of these added to GC tables? If not..
184  glDeleteBuffers(VBO_COUNT, rep->VBO_buffers);
185  FREE_IF_NZ(rep->actualCoord);
186  FREE_IF_NZ(rep->cindex);
187  FREE_IF_NZ(rep->colindex);
188  FREE_IF_NZ(rep->GeneratedTexCoords[0]);
189  FREE_IF_NZ(rep->norindex);
190  FREE_IF_NZ(rep->normal);
191  FREE_IF_NZ(rep->tcindex);
192  FREE_IF_NZ(rep);
193  }
194 }
195 struct X3D_PolyRep * create_polyrep(){
196  int i;
197  struct X3D_PolyRep *polyrep;
198 
199  polyrep = MALLOC(struct X3D_PolyRep *, sizeof(struct X3D_PolyRep));
200  memset(polyrep,0,sizeof(struct X3D_PolyRep));
201  polyrep->ntri = -1;
202  //polyrep->cindex = 0; polyrep->actualCoord = 0; polyrep->colindex = 0; polyrep->color = 0;
203  //polyrep->norindex = 0; polyrep->normal = 0; polyrep->flat_normal = 0; polyrep->GeneratedTexCoords = 0;
204  //polyrep->tri_indices = 0; polyrep->wire_indices = 0; polyrep->actualFog = 0;
205  //polyrep->tcindex = 0;
206  //polyrep->tcoordtype = 0;
207  //polyrep->last_index_type = 0; polyrep->last_normal_type = 0;
208  polyrep->streamed = FALSE;
209 
210  /* for Collision, default texture generation */
211  polyrep->minVals[0] = 999999.9f;
212  polyrep->minVals[1] = 999999.9f;
213  polyrep->minVals[2] = 999999.9f;
214  polyrep->maxVals[0] = -999999.9f;
215  polyrep->maxVals[1] = -999999.9f;
216  polyrep->maxVals[2] = -999999.9f;
217 
218  for (i=0; i<VBO_COUNT; i++)
219  polyrep->VBO_buffers[i] = 0;
220 
221  /* printf ("generating buffers for node %p, type %s\n",p,stringNodeType(p->_nodeType)); */
222  glGenBuffers(1,&polyrep->VBO_buffers[VERTEX_VBO]);
223  glGenBuffers(1,&polyrep->VBO_buffers[INDEX_VBO]);
224  //glGenBuffers(1,&polyrep->VBO_buffers[NORMAL_VBO]);
225  //glGenBuffers(1,&polyrep->VBO_buffers[TEXTURE_VBO0+0]);
226 
227 
228 
229  /* printf ("they are %u %u %u %u\n",polyrep->VBO_buffers[0],polyrep->VBO_buffers[1],polyrep->VBO_buffers[2],polyrep->VBO_buffers[3]); */
230  return polyrep;
231 }
232 
233 
234 
235 #ifdef NURBS_LIB
236 //START MIT LIC >>>>>>>>
237 //some algorithms from "The Nurbs Book", Les Piegl et al
238 int uniformKnot(int n, int p, float *U){
239  int j, k, m, mm;
240  float uniform;
241  m = n + p + 1;
242  k = 0;
243  uniform = 1.0f/(float)(n-p);
244  for(j=0;j<p;k++,j++){
245  U[k] = 0.0f;
246  }
247  mm = n - p + 1;
248  for(j=0;j<mm;j++,k++){
249  U[k] = j*uniform;
250  }
251  for(j=0;j<p;j++,k++){
252  U[k] = 1.0f;
253  }
254  U[8] = 1.0f;
255  printf("U= ");
256  for(j=0;j<m+1;j++){
257  printf(" U[%d]=%f",j,U[j]);
258  }
259  return 1;
260 }
261 
262 //ALGORITHM A2.1 p.68 Piegl
263 int FindSpan(int n, int p, float u, float *U)
264 {
265  /* Determine the knot span index: where u is in U[i]
266  Input:
267  n - # of control points == m - p - 1
268  p - degree of curve = power + 1 ie linear 2, quadratic 3, cubic 4
269  U - knot vector [0 ... m-1]
270  u - scalar curve parameter in range u0 - um
271  Return:
272  knot span index ie if u is between U[i] and U[i+1] return i
273  Internal:
274  order = p + 1
275  m = number of knots = n + order
276  Algorithm:
277  limit the search range between p and m - p - 1 (2 and 4 for this example)
278  assume clamped/pinned ends
279  Example:
280  U = { 0,0,0,1,2,3,4,4,5,5,5 } m = 11
281  spnidx 0 1 2 3 4 5 6 7 8 9
282  u = 2.5 ^ span index == 4
283  u = .0001 ^ span index == 2
284  u = 0 ^ span index == 2
285  u = .4999 ^ span index = 4
286  u = 5 ^ span index = 4
287 
288  */
289  if(1){
290  //dug9 algo, simpler linear search
291  int i, span, m, order;
292  order = p + 1;
293  m = n + order;
294  span = p;
295  for(i=p;i<n;i++){
296  span = i;
297  if(u >= U[i] && u < U[i+1])
298  break;
299  }
300  return span;
301  }else{
302  int low, high, mid;
303  //if(u == U[n+1]) return n;
304  if(u == U[n]) return n-1; //this prevents divide by zero when u = 1
305  low = p; high = n+1; mid = (low+high)/2;
306  while(u < U[mid] || u >= U[mid+1]){
307  if(u < U[mid]) high = mid;
308  else low = mid;
309  mid = (low + high)/2;
310  }
311  return mid;
312  }
313 }
314 //ALGORITHM A2.2 p.70 Piegl
315 int BasisFuns(int span, float u, int p, float *U, float *N){
316  /* Compute the non-vanishing Basis functions
317  Input:
318  span = knot span: which knots is this u in between: if between U[i] and U[i+1], span == i
319  u - scalar curve parameter in range u0 - um
320  p - degree of curve = power + 1 ie linear 2, quadratic 3, cubic 4
321  U - knot vector [0 ... m-1]
322  Output:
323  N - precomputed rational bernstein basis functions for a given span
324  - these are blending weights that say how much of each surrounding
325  control point is used in a given span
326  */
327  int j, r;
328  float left[5], right[5], saved, temp;
329  //float testzero;
330  N[0] =1.0f;
331  for(j=1;j<=p;j++){
332  left[j] = u - U[span+1 - j];
333  right[j] = U[span+j] - u;
334  saved = 0.0f;
335  for(r=0;r<j;r++){
336  //testzero = right[r+1]+left[j-r];
337  //if(fabs(testzero) < .00001)
338  // printf("ouch divide by zero\n");
339  temp = N[r]/(right[r+1]+left[j-r]);
340  N[r] = saved + right[r+1]*temp;
341  saved = left[j-r]*temp;
342  }
343  N[j] = saved;
344  }
345  return 1;
346 }
347 
348 //ALGORITHM A4.1 p.124 Piegl
349 int CurvePoint(int n, int p, float* U, float *Pw, float u, float *C )
350 {
351  /* Compute point on rational B-spline curve
352  Input:
353  n - # of control points == m - p - 1
354  p - degree of curve linear 1, quadratic 2, cubic 3
355  U[] - knot vector [0 ... m], m = n + p + 1
356  Pw[] - control point vector
357  where w means rational/homogenous: Pw[i] = {wi*xi,wi*yi,wi*zi,wi}
358  u - scalar curve parameter in range u0 - um
359  Output:
360  C - 3D point = Cw/w
361  Internal:
362  span = knot span: which knots is this u in between: if between U[i] and U[i+1], span == i
363  N[] - precomputed rational bernstein basis functions for a given span
364  - these are blending weights that say how much of each surrounding control point is used in a given span
365  w - weight, assuming it's uniform
366  */
367  int span,i,j;
368  float N[100], w;
369  float Cw[4];
370  span = FindSpan(n,p,u,U);
371  BasisFuns(span,u,p,U,N);
372  w = 1.0f;
373  for(i=0;i<4;i++) Cw[i] = 0.0f;
374  //Cw[3] = w;
375  for(j=0;j<=p;j++){
376  for(i=0;i<4;i++){
377  Cw[i] += N[j]*Pw[(span-p+j)*4 + i];
378  }
379  }
380  for(i=0;i<3;i++)
381  C[i] = Cw[i]/Cw[3];
382 
383  return 1;
384 }
385 
386 
387 
388 
389 
390 //ALGORITHM A4.3 p.134 Piegl
391 /* example call:
392 ok = SurfacePoint( node->uDimension,node->uOrder-1,node->uKnot.p,
393  node->vDimension,node->vOrder-1,node->vKnot.p,
394  node->controlPoint.p,uv[0],uv[1],xyz);
395 */
396 int SurfacePoint(int n,int p,float *U,
397  int m, int q,float *V,
398  float *Pw,float u,float v,float *S)
399 {
400  /* Compute point on rational B-Spline surface S(u,v)
401  Input:
402  u direction:
403  n - # of control points
404  p - degree of curve linear 1, quadratic 2, cubic 3
405  U[] - knot vector [0 ... n + p + 1]
406  u - scalar curve parameter
407  v direction:
408  m - # of control points
409  q - degree of curve linear 1, quadratic 2, cubic 3
410  V[] - knot vector [0 ... m + q + 1]
411  v - scalar curve parameter
412  Pw[] - control point vector
413  where w means rational/homogenous: Pw[i] = {wi*xi,wi*yi,wi*zi,wi}
414  Output:
415  S - output 3D point = Sw/w
416  */
417  int uspan, vspan, i, l, k;
418  float Nu[100], Nv[100], temp[6][4], Sw[4];
419 
420  uspan = FindSpan(n,p,u,U);
421  BasisFuns(uspan,u,p,U,Nu);
422  vspan = FindSpan(m,q,v,V);
423  BasisFuns(vspan,v,q,V,Nv);
424  for(l=0;l<=q;l++){
425  for(i=0;i<4;i++)
426  temp[l][i] = 0.0f;
427  for(k=0;k<=p;k++){
428  //temp[l] += Nu[k]*Pw[uspan-p+k][vspan-q+l];
429  for(i=0;i<4;i++)
430  temp[l][i] += Nu[k]*Pw[((uspan-p+k)*n + (vspan-q+l))*4 + i];
431 
432  }
433  }
434  for(i=0;i<4;i++) Sw[i] = 0.0f;
435  for(l=0;l<=q;l++){
436  for(i=0;i<4;i++)
437  Sw[i] += Nv[l]*temp[l][i];
438  }
439  for(i=0;i<3;i++)
440  S[i] = Sw[i]/Sw[3];
441  return 1;
442 }
443 // <<<<< END MIT LIC
444 
445 #include <libnurbs2.h>
446 static int DEBG = 0; //glu nurbs surface and trim calls
447 static int DEBGC = 0; //curve calls
448 
449 //defined in Component_RigidBodyPhysics
450 int NNC0(struct X3D_Node* node);
451 void MNC0(struct X3D_Node* node);
452 void MNX0(struct X3D_Node* node);
453 #define NNC(A) NNC0(X3D_NODE(A)) //node needs compiling
454 #define MNC(A) MNC0(X3D_NODE(A)) //mark node compiled
455 #define MNX(A) MNX0(X3D_NODE(A)) //mark node changed
456 #define PPX(A) getTypeNode(X3D_NODE(A)) //possible proto expansion
457 
458 void CALLBACK nurbsError(GLenum errorCode)
459 {
460  //const GLubyte *estring;
461 
462  //estring = gluErrorString(errorCode);
463  //fprintf (stderr, "Nurbs Error: %s\n", estring);
464  printf("ouch from nurbsError\n");
465  // exit (0);
466 }
467 
468 //curve
469 void CALLBACK nurbscurveBegincb(GLenum type, void *ud)
470 {
471  struct X3D_NurbsCurve *node = (struct X3D_NurbsCurve *)ud;
472  if(DEBGC) printf("nurbscurveBegin\n");
473 }
474 void CALLBACK nurbscurveVertexcb(GLfloat *vertex, void *ud)
475 {
476  int i, np,ns;
477  struct SFVec3f *pp;
478  struct X3D_NurbsCurve *node = (struct X3D_NurbsCurve *)ud;
479  ns = node->__points.n;
480  np = node->__numPoints;
481  if(np+1 > ns) {
482  ns = np *2;
483  node->__points.p = REALLOC(node->__points.p,ns * sizeof(struct SFVec3f));
484  node->__points.n = ns;
485  }
486  pp = &node->__points.p[np];
487  for(i=0;i<3;i++)
488  pp->c[i] = vertex[i];
489  node->__numPoints ++;
490  //node->__points.n++;
491  if(DEBGC) printf("nurbscurveVertex\n");
492 }
493 void CALLBACK nurbscurveNormalcb(GLfloat *nml, void *ud)
494 {
495  struct X3D_NurbsCurve *node = (struct X3D_NurbsCurve *)ud;
496  if(DEBGC) printf("nurbscurveNormal\n");
497 }
498 void CALLBACK nurbscurveEndcb(void *ud)
499 {
500  struct X3D_NurbsCurve *node = (struct X3D_NurbsCurve *)ud;
501  //node->__numPoints = node->__points.n;
502  if(DEBGC) printf("nurbscurveEnd\n");
503 }
504 
505 
506 
507 int generateUniformKnotVector(int order, int ncontrol, float *knots){
508  //produced pinned uniform knot vector
509  //caller: please malloc knots = malloc( (ncontrol + order ) * sizeof(float))
510  // http://www.saccade.com/writing/graphics/KnotVectors.pdf
511  //maximum nuber of equalvalue consecutive knots:
512  // a) in middle of knot vector: <= order-1
513  // b) at start and end of knot vector: <= order (for pinned uniform)
514  // exmple order = 4 + ncontrol = 6 => 10 knots
515  // 0 0 0 0 .33 .66 1 1 1 1
516  // example order = 3 + ncontrol = 3 => 6 knots
517  // 0 0 0 1 1 1
518  // example order = 2 + ncontrol = 2 => 4 knots
519  // 0 0 1 1
520  //number of knots == ncontrol + order
521  int j,k,m;
522  float uniform;
523  m = ncontrol - order;
524  k = 0;
525  uniform = 1.0f/(float)(m + 1);
526  for(j=0;j<order;k++,j++){
527  knots[k] = 0.0f;
528  }
529  for(j=0;j<m;j++,k++){
530  knots[k] =uniform*(float)(j+1);
531  }
532  for(j=0;j<order;j++,k++){
533  knots[k] = 1.0f;
534  }
535  return m;
536 }
537 int knotsOK(int order, int ncontrol, int nknots, double *knots){
538  int ok = TRUE;
539 
540  if(nknots < 2 || nknots != ncontrol + order )
541  ok = FALSE;
542  if(ok){
543  int i,nconsec = 1;
544  double lastval = knots[0];
545  for(i=1;i<nknots;i++){
546  if(lastval == knots[i]) nconsec++;
547  else nconsec = 1;
548  if(nconsec > order)
549  ok = false;
550  if(knots[i] < lastval)
551  ok = false;
552  if(!ok) break;
553  lastval = knots[i];
554  }
555  }
556  return ok;
557 }
558 int knotsOKf(int order, int ncontrol, int nknots, float *knots){
559  int ok = TRUE;
560 
561  if(nknots < 2 || nknots != ncontrol + order )
562  ok = FALSE;
563  if(ok){
564  int i, nconsec = 1;
565  double lastval = knots[0];
566  for(i=1;i<nknots;i++){
567  if(lastval == knots[i]) nconsec++;
568  else nconsec = 1;
569  if(nconsec > order)
570  ok = false;
571  if(knots[i] < lastval)
572  ok = false;
573  if(!ok) break;
574  lastval = knots[i];
575  }
576  }
577  return ok;
578 }
579 void compile_ContourPolyline2D(struct X3D_ContourPolyline2D *node){
580  MARK_NODE_COMPILED;
581  if(node->point.n && !node->controlPoint.n){
582  int i;
583  //version v3.0 had a mfvec2f point field, version 3.1+ changed to mfvec2d controlPoint field
584  node->controlPoint.p = MALLOC(struct SFVec2d*,node->point.n * sizeof(struct SFVec2d));
585  for(i=0;i<node->point.n;i++)
586  float2double(node->controlPoint.p[i].c,node->point.p[i].c,2);
587  }
588 }
589 
590 void compile_NurbsCurve(struct X3D_NurbsCurve *node){
591  MARK_NODE_COMPILED
592  {
593  int i,j, n, nk;
594  GLfloat *xyzw, *knots;
595  nk = n = 0;
596  xyzw = knots = NULL;
597  if(node->controlPoint){
598  if(node->controlPoint->_nodeType == NODE_CoordinateDouble){
599  struct Multi_Vec3d *mfd;
600  mfd = &((struct X3D_CoordinateDouble *)(node->controlPoint))->point;
601  n = mfd->n;
602  xyzw = MALLOC(void *, n * 4 * sizeof(GLfloat));
603  for(i=0;i<mfd->n;i++){
604  for(j=0;j<3;j++){
605  xyzw[i*4 + j] = (float) mfd->p[i].c[j];
606  }
607  }
608  }else if(node->controlPoint->_nodeType == NODE_Coordinate){
609  struct Multi_Vec3f *mff;
610  mff = &((struct X3D_Coordinate *)(node->controlPoint))->point;
611  n = mff->n;
612  xyzw = MALLOC(void *, n * 4 * sizeof(GLfloat));
613  for(i=0;i<mff->n;i++){
614  for(j=0;j<3;j++){
615  xyzw[i*4 + j] = mff->p[i].c[j];
616  }
617  }
618  }
619  }else{
620  n = 0;
621  }
622  if(node->weight.n && node->weight.n == n){
623  double w;
624  int m,im;
625  m = min(node->weight.n, n);
626  for(i=0;i<n;i++){
627  im = i < m ? i : m-1;
628  w = node->weight.p[im];
629  xyzw[i*4 + 3] = (float)w;
630  }
631  }else{
632  for(i=0;i<n;i++) xyzw[i*4 + 3] = 1.0;
633  }
634  //if(node->knot.n && node->knot.n == n + node->order ){
635  if(knotsOK(node->order,n,node->knot.n,node->knot.p)){
636 
637  nk = node->knot.n;
638  knots = MALLOC(void *, nk * sizeof(GLfloat));
639  for(i=0;i<nk;i++){
640  knots[i] = (GLfloat)node->knot.p[i];
641  }
642  //printf("good knot nk=%d\n",nk);
643  //for(int ii=0;ii<nk;ii++)
644  // printf("[%d]=%f \n",ii,knots[ii]);
645 
646  }else{
647  static int once = 0;
648  //generate uniform knot vector
649  nk = n + node->order ;
650  //caller: please malloc knots = malloc( (ncontrol + order ) * sizeof(float))
651  knots = MALLOC(void *, nk *sizeof(GLfloat));
652  generateUniformKnotVector(node->order,n, knots);
653  if(!once){
654  int ii;
655  printf("bad knot vector, replacing with:\n");
656  for(ii=0;ii<nk;ii++)
657  printf("[%d]=%f \n",ii,knots[ii]);
658  once = 1;
659  }
660  //nk = 0;
661  }
662 
663  if(n && nk && nk >= n){
664  GLUnurbsObj *theNurb;
665  int ntess, mtess;
666  mtess = node->order + 1;
667  ntess = node->tessellation;
668  theNurb = gluNewNurbsRenderer();
669  gluNurbsProperty(theNurb, GLU_NURBS_MODE, GLU_NURBS_TESSELLATOR);
670  gluNurbsCallbackData(theNurb,(GLvoid*)node);
671  if(0){
672  //chord length or automatic - not implemented properly nor tested thoroughly
673  //if you do chord length in pixels, you need to manually pass in sampling matrices
674  //and somehow you need to trigger a recompile: another call to this compile_
675  // as avatar/viewer moves closer (father) from the nurbs node
676  double model[16], proj[16];
677  float modelf[16], projf[16];
678  int viewPort[10];
679  if(ntess > 0)
680  mtess = ntess;
681  else if(ntess < 0)
682  mtess = -ntess;
683  node->__points.p = MALLOC(void *, sizeof(struct SFVec3f)*n*10); // just a guess to get started
684  node->__points.n = n*10; //.n will be used for realloc test in callbacks
685 
686  gluNurbsProperty(theNurb, GLU_SAMPLING_TOLERANCE, (float)(mtess)); //25.0);
687  if(ntess < 0)
688  gluNurbsProperty(theNurb,GLU_SAMPLING_METHOD,GLU_PATH_LENGTH); //pixels, the default
689  else
690  gluNurbsProperty(theNurb,GLU_SAMPLING_METHOD,GLU_PARAMETRIC_TOLERANCE);
691  gluNurbsProperty(theNurb, GLU_AUTO_LOAD_MATRIX,GL_FALSE);
692 
693  FW_GL_GETDOUBLEV(GL_MODELVIEW_MATRIX, model);
694  FW_GL_GETDOUBLEV(GL_PROJECTION_MATRIX, proj);
695  FW_GL_GETINTEGERV(GL_VIEWPORT, viewPort);
696  for(i=0;i<16;i++){
697  modelf[i] = (float)model[i];
698  projf[i] = (float)proj[i];
699  }
700  gluLoadSamplingMatrices(theNurb,modelf,projf,viewPort);
701  }
702  if(1){
703  //uniform spacing of sampling points in u (or uv) parameter space - works
704  //node must specify tesselation value. see specs for interpretation
705  // http://www.web3d.org/documents/specifications/19775-1/V3.3/Part01/components/nurbs.html#NurbsCurve
706  if(ntess > 0)
707  mtess = max(mtess,ntess+1);
708  else if(ntess < 0)
709  mtess = max(mtess,(-ntess * n) + 1);
710  else
711  mtess = max(mtess,2*n + 1);
712  mtess = (int)((float)mtess * node->_tscale);
713  node->__points.p = MALLOC(void *, sizeof(struct SFVec3f)*mtess+1);
714  node->__points.n = mtess; //.n will be used for realloc test in callbacks
715  gluNurbsProperty(theNurb,GLU_SAMPLING_METHOD,GLU_DOMAIN_DISTANCE);
716  gluNurbsProperty(theNurb,GLU_U_STEP,(GLfloat)mtess);
717  }
718  gluNurbsProperty(theNurb, GLU_DISPLAY_MODE, GLU_FILL);
719  gluNurbsCallback(theNurb, GLU_ERROR, nurbsError);
720  gluNurbsCallback(theNurb, GLU_NURBS_BEGIN_DATA, nurbscurveBegincb);
721  gluNurbsCallback(theNurb, GLU_NURBS_VERTEX_DATA, nurbscurveVertexcb);
722  gluNurbsCallback(theNurb, GLU_NURBS_NORMAL_DATA, nurbscurveNormalcb);
723  gluNurbsCallback(theNurb, GLU_NURBS_END_DATA, nurbscurveEndcb);
724  gluBeginCurve(theNurb);
725  node->__numPoints = 0;
726  gluNurbsCurve(theNurb,nk,knots,4,xyzw,node->order,GL_MAP1_VERTEX_4);
727  gluEndCurve(theNurb);
728  gluDeleteNurbsRenderer(theNurb);
729  node->__points.n = node->__numPoints;
730  }
731  }
732 }
733 
734 void render_NurbsCurve(struct X3D_NurbsCurve *node){
735  ttglobal tg = gglobal();
736  COMPILE_IF_REQUIRED
737  // from Arc2D
738  if (node->__numPoints>0) {
739  // for BoundingBox calculations
740  setExtent( node->EXTENT_MAX_X, node->EXTENT_MIN_X,
741  node->EXTENT_MAX_Y, node->EXTENT_MIN_Y, 0.0f,0.0f,X3D_NODE(node));
742 
743  //OLDCODE GET_COLOUR_POINTER
744  LIGHTING_OFF
745  DISABLE_CULL_FACE
746  //OLDCODE DO_COLOUR_POINTER
747 
748  FW_GL_VERTEX_POINTER (3,GL_FLOAT,0,(GLfloat *)node->__points.p);
749  sendArraysToGPU (GL_LINE_STRIP, 0, node->__numPoints);
750  tg->Mainloop.trisThisLoop += node->__numPoints;
751  }
752 }
753 
754 void compile_NurbsTextureCoordinate(struct X3D_NurbsTextureCoordinate *node){
755  //get knots from double to float and QC the knots
756  int nc, nu, nku, nkv, nv, i,j;
757  float *knotsu, *knotsv, *xyzw;
758 
759  struct Multi_Vec2f *mff;
760  mff = &node->controlPoint;
761  nc = mff->n;
762  xyzw = MALLOC(void *, nc * 4 * sizeof(GLfloat));
763  for(i=0;i<mff->n;i++){
764  for(j=0;j<2;j++){
765  xyzw[i*4 + j] = mff->p[i].c[j];
766  }
767  xyzw[i*4 + 2] = 0.0f; //z == 0 for 2D
768  xyzw[i*4 + 3] = 1.0f; //homogenous 1
769  }
770  nu = node->uDimension;
771  nv = node->vDimension;
772  if(node->weight.n && node->weight.n == nc){
773  double w;
774  int m,im;
775  m = min(node->weight.n, nc);
776  for(i=0;i<nc;i++){
777  im = i < m ? i : m-1;
778  w = node->weight.p[im];
779  xyzw[i*4 + 3] = (float) w;
780  }
781  }else{
782  for(i=0;i<nc;i++) xyzw[i*4 + 3] = 1.0;
783  }
784  nu = node->uDimension;
785  nv = node->vDimension;
786  //int knotsOK(int order, int ncontrol, int nknots, double *knots)
787  //if(node->uKnot.n && node->uKnot.n == nu + node->uOrder ){
788  if(knotsOK(node->uOrder,nu,node->uKnot.n,node->uKnot.p)){
789  //could do another check: max number of consecutive equal value knots == order
790  //could do another check: knot values == or ascending
791  nku = node->uKnot.n;
792  knotsu = MALLOC(void *, nku * sizeof(GLfloat));
793  for(i=0;i<nku;i++){
794  knotsu[i] = (GLfloat)node->uKnot.p[i];
795  }
796  if(DEBG){
797  int ii;
798  printf("good u knot vector nk=%d\n",nku);
799  for(ii=0;ii<nku;ii++)
800  printf("[%d]=%f \n",ii,knotsu[ii]);
801  }
802 
803  }else{
804  //generate uniform knot vector
805  static int once = 0;
806  nku = nu + node->uOrder ;
807  //caller: please malloc knots = MALLOC(void *, (ncontrol + order ) * sizeof(float))
808  knotsu = MALLOC(void *, nku *sizeof(GLfloat));
809  generateUniformKnotVector(node->uOrder,nu, knotsu);
810  if(!once){
811  int ii;
812  printf("bad u knot vector given, replacing with:\n");
813  for(ii=0;ii<nku;ii++)
814  printf("[%d]=%f \n",ii,knotsu[ii]);
815  once = 1;
816  }
817  //nk = 0;
818  }
819 
820  if(knotsOK(node->vOrder,nv,node->vKnot.n,node->vKnot.p)){
821  //if(node->vKnot.n && node->vKnot.n == nv + node->vOrder ){
822  nkv = node->vKnot.n;
823  knotsv = MALLOC(void *, nkv * sizeof(GLfloat));
824  for(i=0;i<nkv;i++){
825  knotsv[i] = (GLfloat)node->vKnot.p[i];
826  }
827  if(DEBG){
828  int ii;
829  printf("good v knot vector nk=%d\n",nkv);
830  for(ii=0;ii<nkv;ii++)
831  printf("[%d]=%f \n",ii,knotsv[ii]);
832  }
833 
834  }else{
835  static int once = 0;
836  //generate uniform knot vector
837  nkv = nv + node->vOrder ;
838  //caller: please malloc knots = MALLOC(void *, (ncontrol + order ) * sizeof(float))
839  knotsv = MALLOC(void *, nkv *sizeof(GLfloat));
840  generateUniformKnotVector(node->vOrder,nv, knotsv);
841  if(!once){
842  int ii;
843  printf("bad v knot vector given, replacing with:\n");
844  for(ii=0;ii<nkv;ii++)
845  printf("[%d]=%f \n",ii,knotsv[ii]);
846  once = 1;
847  }
848  if(!knotsOKf(node->vOrder,nv,nkv,knotsv))
849  printf("ouch still not right knot vector\n");
850  //nk = 0;
851  }
852  node->_uKnot.p = knotsu;
853  node->_uKnot.n = nku;
854  node->_vKnot.p = knotsv;
855  node->_vKnot.n = nkv;
856  node->_controlPoint.p = (struct SFVec4f*)xyzw;
857  node->_controlPoint.n = nc;
858  MNC(node);
859 
860 }
861 int getNurbsSurfacePoint(struct X3D_Node *nurbsSurfaceNode, float *uv, float *xyz){
862  int ret = 0;
863  if(nurbsSurfaceNode){
864  switch(nurbsSurfaceNode->_nodeType){
865  case NODE_NurbsTextureCoordinate:
866  {
867  struct X3D_NurbsTextureCoordinate *node = (struct X3D_NurbsTextureCoordinate *)PPX(nurbsSurfaceNode);
868  if(NNC(node)) compile_NurbsTextureCoordinate(node);
869  ret = SurfacePoint( node->uDimension,node->uOrder-1,node->_uKnot.p,
870  node->vDimension,node->vOrder-1,node->_vKnot.p,
871  (float *)node->_controlPoint.p,uv[0],uv[1],xyz);
872  }
873  break;
874  case NODE_NurbsPatchSurface:
875  break;
876  case NODE_NurbsTrimmedSurface:
877  break;
878  default:
879  break;
880  }
881  }
882 
883  return ret;
884 }
885 
886 /* GenPolyrep functions assume a node inherits from X3DGeometryNode.
887  NurbsPatchSurface inherits from X3DParametricGeometryNode.
888  So we can't use all the genpolyrep stuff, until we have our node compiled into one.
889  Then we can delegate back to generic polyrep functions
890 render - can delegate to polyrep
891 rendray - can delegate to polyrep
892 (x make - we will do a compile instead)
893 collide - can delegate to polyrep
894 compile - custom: we will convert our parametric surface to a polyrep in here
895 
896  */
897 
898  //stripstate - used for capturing the callback data when using gluNurbs in TESSELATOR mode
899  //we use TESSELATOR mode instead of RENDER mode because we use GLES2, not desktop GL.
900  //the libnurbs we have has the GL rendering stuff ifdefed out.
901  //so we capture the data in the callbacks, and then we can do what we normally do with
902  //mesh data in GLES2
903 
904  struct stripState{
905  int type;
906  struct Vector pv; //vector of vertex points
907  struct Vector nv; //vector of normals
908  struct Vector tv; //vector of texcoords
909 };
910 
911 //surface
912 // strategy: on begin, we'll store the gl_ type, and zero __points
913 // then accumulate the __points
914 // and on end, convert the __points to ->_intern->polyrep points and triangle indices
915 // using the stored type as a guide
916 void CALLBACK nurbssurfBegincb(GLenum type, void *ud)
917 {
918  struct stripState ss;
919  struct Vector * strips = (struct Vector *)ud;
920  if(0) if(DEBG) printf("callback nurbsSurfaceBegin\n");
921  if(0){
922  printf("nurbssurfBegin type = ");
923  switch(type){
924  case GL_QUAD_STRIP: printf("QUAD_STRIP");break;
925  case GL_TRIANGLE_STRIP: printf("TRIANGLE_STRIP");break;
926  case GL_TRIANGLE_FAN: printf("TRIANGLE_FAN");break;
927  case GL_TRIANGLES: printf("TRIANGLES");break;
928  default:
929  printf("not sure %x %d",type,type);
930  }
931  printf("\n");
932  }
933  ss.nv.n = 0;
934  ss.nv.allocn = 0;
935  ss.nv.data = NULL;
936  ss.tv.n = 0;
937  ss.tv.allocn = 0;
938  ss.tv.data = NULL;
939  ss.pv.n = 0;
940  ss.pv.allocn = 0;
941  ss.pv.data = NULL;
942  ss.type = type;
943  vector_pushBack(struct stripState,strips,ss);
944 }
945 void CALLBACK nurbssurfVertexcb(GLfloat *vertex, void *ud)
946 {
947  struct stripState *ss;
948  struct SFVec3f pp;
949  struct Vector * strips = (struct Vector *)ud;
950  ss = vector_get_ptr(struct stripState,strips,strips->n -1);
951  memcpy(&pp,vertex,sizeof(struct SFVec3f));
952  vector_pushBack(struct SFVec3f,&ss->pv,pp);
953  //vector_set(struct stripState,strips,strips->n-1,ss);
954 
955  if(0) printf("callback nurbssurfVertex %f %f %f\n",vertex[0],vertex[1],vertex[2]);
956 }
957 void CALLBACK nurbssurfNormalcb(GLfloat *nml, void *ud)
958 {
959  struct stripState *ss;
960  struct SFVec3f pp;
961  struct Vector * strips = (struct Vector *)ud;
962  ss = vector_get_ptr(struct stripState,strips,strips->n -1);
963  memcpy(&pp,nml,sizeof(struct SFVec3f));
964  vector_pushBack(struct SFVec3f,&ss->nv,pp);
965  //vector_set(struct stripState,strips,strips->n-1,ss);
966 
967  if(0) printf("callback nurbssurfNormal\n");
968 }
969 void CALLBACK nurbssurfEndcb(void *ud)
970 {
971  struct stripState *ss;
972  struct Vector * strips = (struct Vector *)ud;
973  ss = vector_get_ptr(struct stripState,strips,strips->n -1);
974  if(0){
975  int i;
976  printf("nurbssurfEnd #p %d #n %d\n",ss->pv.n, ss->nv.n);
977  for(i=0;i<ss->pv.n;i++){
978  struct SFVec3f pp = vector_get(struct SFVec3f,&ss->pv,i);
979  printf("%f %f %f\n",pp.c[0],pp.c[1],pp.c[2]);
980  }
981  }
982  if(0) if(DEBG) printf("callback nurbsSurfaceEnd\n");
983 
984 }
985 void CALLBACK nurbssurfTexcoordcb(GLfloat *tCrd, void *ud){
986  static int count = 0;
987  struct stripState *ss;
988  struct SFVec2f tp;
989  struct Vector * strips = (struct Vector *)ud;
990  ss = vector_get_ptr(struct stripState,strips,strips->n -1);
991  memcpy(&tp,tCrd,sizeof(struct SFVec2f));
992  vector_pushBack(struct SFVec2f,&ss->tv,tp);
993  //vector_set(struct stripState,strips,strips->n-1,ss);
994  //printf("%f %f %f\n",tCrd[0],tCrd[1],0.0f);
995  //count++;
996  //if(count % 50 == 0)
997  // printf("\n");
998  if(0) if(DEBG)
999  printf("callback nurbssufTexcoordcb\n");
1000 }
1001 
1002 
1003 
1004 #define GL_QUAD_STRIP 0x0008
1005 
1006 static int USETXCOORD = 1;
1007 void convert_strips_to_polyrep(struct Vector * strips,struct X3D_NurbsTrimmedSurface *node){
1008  //this is a bit like compile_polyrep, except the virt_make is below
1009 
1010  int i, j, npoints, np, ni, ntri, nindex, ntc;
1011  struct stripState *ss;
1012  struct X3D_PolyRep *rep_, *polyrep;
1013  struct X3D_TextureCoordinate * tcnode, tcnode0;
1014  GLuint *cindex, *norindex, *tcindex;
1015  float *tcoord = NULL;
1016 
1017  //from compile_polyrep:
1018  //node = X3D_NODE(innode);
1019  //virt = virtTable[node->_nodeType];
1020 
1021  /* first time through; make the intern structure for this polyrep node */
1022  if(node->_intern){
1023  polyrep = node->_intern;
1024  FREE_IF_NZ(polyrep->cindex);
1025  FREE_IF_NZ(polyrep->actualCoord);
1026  FREE_IF_NZ(polyrep->GeneratedTexCoords[0]);
1027  FREE_IF_NZ(polyrep->colindex);
1028  FREE_IF_NZ(polyrep->color);
1029  FREE_IF_NZ(polyrep->norindex);
1030  FREE_IF_NZ(polyrep->normal);
1031  FREE_IF_NZ(polyrep->flat_normal);
1032  FREE_IF_NZ(polyrep->tcindex);
1033  FREE_IF_NZ(polyrep->wire_indices);
1034  //glDeleteBuffers(VBO_COUNT,polyrep->VBO_buffers); //streampoly checks if 0 before doing a new one
1035  }
1036  if(!node->_intern)
1037  node->_intern = create_polyrep();
1038 
1039  rep_ = polyrep = node->_intern;
1040 
1041 
1042  /* if multithreading, tell the rendering loop that we are regenning this one */
1043  /* if singlethreading, this'll be set to TRUE before it is tested */
1044 
1045  //<< END FROM Compile_polyrep
1046 
1047  // Start Virt_make_polyrep section >>>
1048  //texcoord
1049  if(USETXCOORD){
1050  rep_->ntexdim[0] = 2;
1051  rep_->tcoordtype = NODE_TextureCoordinate; //??
1052  rep_->ntcoord = 1;
1053  }
1054  tcnode = &tcnode0; //createNewX3DNode(NODE_TextureCoordinate);
1055 
1056  npoints = nindex = ntc = 0;
1057  for(i=0;i<strips->n;i++){
1058  ss = vector_get_ptr(struct stripState,strips,i);
1059  npoints += ss->pv.n;
1060  ntc += ss->tv.n;
1061  switch(ss->type){
1062  case GL_QUAD_STRIP: nindex += (ss->pv.n -2)/2 * 5;break;
1063  case GL_TRIANGLE_STRIP: nindex += (ss->pv.n -2);break;
1064  case GL_TRIANGLE_FAN: nindex += (ss->pv.n -2);break;
1065  case GL_TRIANGLES: nindex += (ss->pv.n -2);break;
1066  default:
1067  nindex += (ss->pv.n -2);
1068  }
1069  }
1070  if (npoints > 0)
1071  {
1072  //printf("npoints %d ntc %d\n",npoints,ntc);
1073  rep_->actualCoord = MALLOC(void *, npoints * 3 * sizeof(float));
1074  rep_->normal = MALLOC(void *, npoints * 3 * sizeof(float));
1075  //if(USETXCOORD) rep->GeneratedTexCoords[0] = MALLOC(void *, npoints * 2 * sizeof(float));
1076  //if(USETXCOORD){
1077  // tcnode->point.p = MALLOC(void*, npoints * 2 * sizeof(float));
1078  // tcnode->point.n = npoints;
1079  //}
1080  //rep->t
1081  }
1082  rep_->ntri = ntri = nindex; //we'll over-malloc
1083 
1084  if (rep_->ntri > 0)
1085  {
1086  cindex = rep_->cindex = MALLOC(GLuint *, sizeof(GLuint)*3*(ntri));
1087  norindex = rep_->norindex = MALLOC(GLuint *,sizeof(GLuint)*3*ntri);
1088  //if(USETXCOORD) rep_->tcindex = MALLOC(void *, ntri * 4 * sizeof(GLuint));
1089  tcindex = rep_->tcindex = MALLOC(GLuint*, sizeof(GLuint)*3*(ntri));
1090  // colindex = rep_->colindex = MALLOC(GLuint *, sizeof(*(rep_->colindex))*3*(ntri));
1091 
1092  //FREE_IF_NZ(rep_->GeneratedTexCoords[0]);
1093  // we'll pass a X3D_TexCoordinate node //rep_->GeneratedTexCoords[0]
1094  tcoord = MALLOC (float *, sizeof (float) * ntri * 2 * 3);
1095 
1096  }
1097 
1098  np = 0;
1099  ni = 0;
1100  ntri = 0;
1101  for(i=0;i<strips->n;i++){
1102  ss = vector_get_ptr(struct stripState,strips,i);
1103 //printf("ss.pv.n=%d nv.n=%d tv.n=%d\n",ss.pv.n,ss.nv.n,ss.tv.n);
1104  memcpy(&rep_->actualCoord[np*3],ss->pv.data,ss->pv.n * 3 * sizeof(float));
1105  memcpy(&rep_->normal[np*3],ss->nv.data,ss->nv.n * 3 * sizeof(float));
1106  if(USETXCOORD && tcoord) memcpy(&tcoord[np*2],ss->tv.data,ss->tv.n * 2 * sizeof(float));
1107  switch(ss->type){
1108  case GL_QUAD_STRIP:
1109  for(j=0;j<ss->pv.n -2;j+=2){
1110  rep_->cindex[ni++] = np+j;
1111  rep_->cindex[ni++] = np+j+1;
1112  rep_->cindex[ni++] = np+j+3;
1113  //rep->cindex[ni++] = -1;
1114  rep_->cindex[ni++] = np+j+3;
1115  rep_->cindex[ni++] = np+j+2;
1116  rep_->cindex[ni++] = np+j;
1117  //rep->cindex[ni++] = -1;
1118  //memcpy(&rep->norindex[ntri*4],&rep->cindex[ntri*4],2*4*sizeof(int));
1119  memcpy(&rep_->norindex[ntri*3],&rep_->cindex[ntri*3],2*3*sizeof(int));
1120  if(USETXCOORD) memcpy(&rep_->tcindex[ntri*3],&rep_->cindex[ntri*3],2*3*sizeof(int));
1121  ntri += 2;
1122  }
1123  break;
1124  case GL_TRIANGLE_STRIP:
1125  nindex += (ss->pv.n -2);
1126  break;
1127  case GL_TRIANGLE_FAN:
1128  //nindex += (ss.pv.n -2);
1129  for(j=0;j<ss->pv.n -2;j+=1){
1130  rep_->cindex[ni++] = np;
1131  rep_->cindex[ni++] = np+j+1;
1132  rep_->cindex[ni++] = np+j+2;
1133  memcpy(&rep_->norindex[ntri*3],&rep_->cindex[ntri*3],3*sizeof(int));
1134  if(USETXCOORD) memcpy(&rep_->tcindex[ntri*3],&rep_->cindex[ntri*3],3*sizeof(int));
1135  ntri += 1;
1136  }
1137  break;
1138  case GL_TRIANGLES:
1139  nindex += (ss->pv.n -2);
1140  break;
1141  default:
1142  nindex += (ss->pv.n -2);
1143  }
1144  np += ss->pv.n;
1145  }
1146  rep_->ntri = ntri;
1147  if(node->texCoord && node->texCoord->_nodeType == NODE_NurbsTextureCoordinate){
1148  static FILE *fp = NULL;
1149  for(i=0;i<np;i++){
1150  float stru[4];
1151  float xyz[4];
1152  //treat above callback texturecoord as uv,
1153  //and do lookup on NurbsTextureCoordinate surface to get new texture coord st
1154  xyz[0] = tcoord[i*2 + 1];
1155  xyz[1] = tcoord[i*2 + 0]; //don't know why but I have to swap xy to match octaga
1156  xyz[2] = 0.0f;
1157  xyz[3] = tcoord[i*2 + 3];
1158  getNurbsSurfacePoint(node->texCoord, xyz, stru);
1159  //memcpy(&tcoord[i*2],stru,2*sizeof(float));
1160  tcoord[i*2 + 0] = stru[0];
1161  tcoord[i*2 + 1] = stru[1];
1162  //if(1){
1163  // static int once = 0;
1164  // if(!once) fp= fopen("nurbstexturecoord.txt","w+");
1165  // once = 1;
1166  // fprintf(fp,"%d uv %f %f st %f %f\n",i,tcoord[i*2 +0],tcoord[i*2 + 1],stru[0],stru[1]);
1167  //}
1168  }
1169  //if(fp) fclose(fp);
1170  }
1171  tcnode->point.p = (struct SFVec2f*)tcoord;
1172  tcnode->point.n = np;
1173  if(0) for(i=0;i<tcnode->point.n;i++){
1174  printf("%d %f %f\n",i,tcnode->point.p[i].c[0],tcnode->point.p[i].c[1]);
1175  if(i % 50 == 0)
1176  printf("\n");
1177  }
1178 
1179  //END virt_make_polyrep section <<<<<<
1180 
1181  //FROM Compile_polyrep
1182  if (polyrep->ntri != 0) {
1183  //float *fogCoord = NULL;
1184  stream_polyrep(node, NULL,NULL,NULL,NULL, tcnode);
1185  /* and, tell the rendering process that this shape is now compiled */
1186  }
1187  FREE_IF_NZ(tcnode->point.p);
1188  //else wait for set_coordIndex to be converted to coordIndex
1189  //MARK POLYREP COMPILED
1190  polyrep->irep_change = node->_change;
1191 
1192  /*
1193  // dump then can copy and paste to x3d or wrl IndexedFaceSet.coordIndex and Coordinate.point fields
1194  FILE * fp = fopen("IFS_DUMP.txt","w+");
1195  fprintf(fp,"#vertices %d\n",np);
1196  for(i=0;i<np;i++){
1197  fprintf(fp,"%f %f %f\n",rep->actualCoord[i*3 +0],rep->actualCoord[i*3 +1],rep->actualCoord[i*3 +2]);
1198  }
1199  fprintf(fp,"#face indices %d\n",ni);
1200  for(i=0;i<ni;i++){
1201  fprintf(fp,"%d ",rep->cindex[i]);
1202  if((ni+1) % 3 == 0)
1203  fprintf(fp,"%d ",-1);
1204  }
1205  fprintf(fp,"\n");
1206  fclose(fp);
1207  */
1208 
1209 }
1210 
1211 void compile_NurbsSurface(struct X3D_NurbsPatchSurface *node, struct Multi_Node *trim){
1212  MARK_NODE_COMPILED
1213 
1214  {
1215  int i,j, n, nu, nv, nku, nkv;
1216  GLfloat *xyzw, *knotsu, *knotsv;
1217  ppComponent_NURBS p = (ppComponent_NURBS)gglobal()->Component_NURBS.prv;
1218 
1219  nku = nkv = nu = nv = n = 0;
1220  xyzw = knotsu = knotsv = NULL;
1221  // I should call something like:
1222  // struct Multi_Vec3f *getCoordinate (struct X3D_Node *innode, char *str);
1223  // to get the control points - it will do proto expansion, compile the coordinate node if needed
1224  // (should do something similar with texcoord when implemented,
1225  // as I think it needs conversion from controlpoint spacing to sampling/tesselation spacing for use in polyrep)
1226  // here's an amature shortcut that returns doubles
1227  if(node->controlPoint){
1228  if(node->controlPoint->_nodeType == NODE_CoordinateDouble){
1229  struct Multi_Vec3d *mfd;
1230  mfd = &((struct X3D_CoordinateDouble *)(node->controlPoint))->point;
1231  n = mfd->n;
1232  xyzw = MALLOC(void *, n * 4 * sizeof(GLfloat));
1233  for(i=0;i<mfd->n;i++){
1234  for(j=0;j<3;j++){
1235  xyzw[i*4 + j] = (float)mfd->p[i].c[j];
1236  }
1237  }
1238  }else if(node->controlPoint->_nodeType == NODE_Coordinate){
1239  struct Multi_Vec3f *mff;
1240  mff = &((struct X3D_Coordinate *)(node->controlPoint))->point;
1241  n = mff->n;
1242  xyzw = MALLOC(void *, n * 4 * sizeof(GLfloat));
1243  for(i=0;i<mff->n;i++){
1244  for(j=0;j<3;j++){
1245  xyzw[i*4 + j] = mff->p[i].c[j];
1246  }
1247  }
1248  }
1249  }else{
1250  n = 0;
1251  }
1252  if(node->weight.n && node->weight.n == n){
1253  double w;
1254  int m,im;
1255  m = min(node->weight.n, n);
1256  for(i=0;i<n;i++){
1257  im = i < m ? i : m-1;
1258  w = node->weight.p[im];
1259  xyzw[i*4 + 3] = (float)w;
1260  }
1261  }else{
1262  for(i=0;i<n;i++) xyzw[i*4 + 3] = 1.0;
1263  }
1264  nu = node->uDimension;
1265  nv = node->vDimension;
1266  //int knotsOK(int order, int ncontrol, int nknots, double *knots)
1267  //if(node->uKnot.n && node->uKnot.n == nu + node->uOrder ){
1268  if(knotsOK(node->uOrder,nu,node->uKnot.n,node->uKnot.p)){
1269  //could do another check: max number of consecutive equal value knots == order
1270  //could do another check: knot values == or ascending
1271  nku = node->uKnot.n;
1272  knotsu = MALLOC(void *, nku * sizeof(GLfloat));
1273  for(i=0;i<nku;i++){
1274  knotsu[i] = (GLfloat)node->uKnot.p[i];
1275  }
1276  if(DEBG){
1277  int ii;
1278  printf("good u knot vector nk=%d\n",nku);
1279  for(ii=0;ii<nku;ii++)
1280  printf("[%d]=%f \n",ii,knotsu[ii]);
1281  }
1282 
1283  }else{
1284  //generate uniform knot vector
1285  static int once = 0;
1286  nku = nu + node->uOrder ;
1287  //caller: please malloc knots = MALLOC(void *, (ncontrol + order ) * sizeof(float))
1288  knotsu = MALLOC(void *, nku *sizeof(GLfloat));
1289  generateUniformKnotVector(node->uOrder,nu, knotsu);
1290  if(!once){
1291  int ii;
1292  printf("bad u knot vector given, replacing with:\n");
1293  for(ii=0;ii<nku;ii++)
1294  printf("[%d]=%f \n",ii,knotsu[ii]);
1295  once = 1;
1296  }
1297  //nk = 0;
1298  }
1299 
1300  if(knotsOK(node->vOrder,nv,node->vKnot.n,node->vKnot.p)){
1301  //if(node->vKnot.n && node->vKnot.n == nv + node->vOrder ){
1302  nkv = node->vKnot.n;
1303  knotsv = MALLOC(void *, nkv * sizeof(GLfloat));
1304  for(i=0;i<nkv;i++){
1305  knotsv[i] = (GLfloat)node->vKnot.p[i];
1306  }
1307  if(DEBG){
1308  int ii;
1309  printf("good v knot vector nk=%d\n",nkv);
1310  for(ii=0;ii<nkv;ii++)
1311  printf("[%d]=%f \n",ii,knotsv[ii]);
1312  }
1313 
1314  }else{
1315  static int once = 0;
1316  //generate uniform knot vector
1317  nkv = nv + node->vOrder ;
1318  //caller: please malloc knots = MALLOC(void *, (ncontrol + order ) * sizeof(float))
1319  knotsv = MALLOC(void *, nkv *sizeof(GLfloat));
1320  generateUniformKnotVector(node->vOrder,nv, knotsv);
1321  if(!once){
1322  int ii;
1323  printf("bad v knot vector given, replacing with:\n");
1324  for(ii=0;ii<nkv;ii++)
1325  printf("[%d]=%f \n",ii,knotsv[ii]);
1326  once = 1;
1327  }
1328  //nk = 0;
1329  }
1330 
1331  if(n && nku && nkv){
1332  static GLUnurbsObj *theNurb = NULL;
1333  int ntessu, ntessv, mtessu, mtessv;
1334  struct Vector * strips;
1335  struct X3D_Node * texCoordNode;
1336  int texcoordnodeIsGenerated;
1337 
1338  mtessu = node->uOrder + 1;
1339  ntessu = node->uTessellation;
1340  mtessv = node->vOrder + 1;
1341  ntessv = node->vTessellation;
1342 
1343  if(DEBG) printf("gluNewNurbsRenderer\n");
1344  if(!theNurb)
1345  theNurb = gluNewNurbsRenderer();
1346  gluNurbsProperty(theNurb, GLU_NURBS_MODE, GLU_NURBS_TESSELLATOR);
1347  if(0){
1348  //chord length or automatic - not implemented properly nor tested thoroughly
1349  //if you do chord length in pixels, you need to manually pass in sampling matrices
1350  //and somehow you need to trigger a recompile: another call to this compile_
1351  // as avatar/viewer moves closer (father) from the nurbs node
1352  double model[16], proj[16];
1353  float modelf[16], projf[16];
1354  int viewPort[10];
1355  if(ntessu > 0)
1356  mtessu = ntessu;
1357  else if(ntessu < 0)
1358  mtessu = -ntessu;
1359 
1360  gluNurbsProperty(theNurb, GLU_SAMPLING_TOLERANCE, (float)(mtessu)); //25.0);
1361  if(ntessu < 0)
1362  gluNurbsProperty(theNurb,GLU_SAMPLING_METHOD,GLU_PATH_LENGTH); //pixels, the default
1363  else
1364  gluNurbsProperty(theNurb,GLU_SAMPLING_METHOD,GLU_PARAMETRIC_TOLERANCE);
1365  gluNurbsProperty(theNurb, GLU_AUTO_LOAD_MATRIX,GL_FALSE);
1366 
1367  FW_GL_GETDOUBLEV(GL_MODELVIEW_MATRIX, model);
1368  FW_GL_GETDOUBLEV(GL_PROJECTION_MATRIX, proj);
1369  FW_GL_GETINTEGERV(GL_VIEWPORT, viewPort);
1370  for(i=0;i<16;i++){
1371  modelf[i] = (float)model[i];
1372  projf[i] = (float)proj[i];
1373  }
1374  gluLoadSamplingMatrices(theNurb,modelf,projf,viewPort);
1375  }
1376  if(1){
1377  //uniform spacing of sampling points in u (or uv) parameter space - works
1378  //node must specify tesselation value. see specs for interpretation
1379  // http://www.web3d.org/documents/specifications/19775-1/V3.3/Part01/components/nurbs.html#NurbsCurve
1380  if(ntessu > 0)
1381  mtessu = max(mtessu,ntessu+1);
1382  else if(ntessu < 0)
1383  mtessu = max(mtessu,(-ntessu * nu) + 1);
1384  else
1385  mtessu = max(mtessu,2*nu + 1);
1386 
1387  if(ntessv > 0)
1388  mtessv = max(mtessv,ntessv+1);
1389  else if(ntessv < 0)
1390  mtessv = max(mtessv,(-ntessv * nv) + 1);
1391  else
1392  mtessv = max(mtessv,2*nv + 1);
1393  mtessu = (int)((float)mtessu * node->_tscale);
1394  mtessv = (int)((float)mtessv * node->_tscale);
1395 
1396  gluNurbsProperty(theNurb,GLU_SAMPLING_METHOD,GLU_DOMAIN_DISTANCE);
1397  gluNurbsProperty(theNurb,GLU_U_STEP,(GLfloat)mtessu);
1398  gluNurbsProperty(theNurb,GLU_V_STEP,(GLfloat)mtessv);
1399  }
1400  gluNurbsProperty(theNurb, GLU_DISPLAY_MODE, GLU_FILL);
1401 
1402  gluNurbsCallback(theNurb, GLU_ERROR, nurbsError);
1403  gluNurbsCallback(theNurb, GLU_NURBS_BEGIN_DATA, nurbssurfBegincb);
1404  gluNurbsCallback(theNurb, GLU_NURBS_VERTEX_DATA, nurbssurfVertexcb);
1405  gluNurbsCallback(theNurb, GLU_NURBS_NORMAL_DATA, nurbssurfNormalcb);
1406  gluNurbsCallback(theNurb, GLU_NURBS_END_DATA, nurbssurfEndcb);
1407  gluNurbsCallback(theNurb, GLU_NURBS_TEXTURE_COORD_DATA, nurbssurfTexcoordcb);
1408 
1409  strips = newVector(struct stripState,20);
1410  gluNurbsCallbackData(theNurb,(GLvoid*)strips);
1411 
1412  if(DEBG) printf("gluBeginSurface \n");
1413  gluBeginSurface(theNurb);
1414  gluNurbsSurface(theNurb,nku,knotsu,nkv,knotsv,4,4*nu,xyzw,node->uOrder,node->vOrder,GL_MAP2_VERTEX_4);
1415  /*
1416  TextureCoordinate handling
1417  https://www.opengl.org/discussion_boards/showthread.php/127668-Texture-mapping-for-NURBS
1418  https://www.cs.drexel.edu/~david/Classes/ICG/Lectures/Lecture7.pdf p.56 of slides
1419 
1420  texture coordinate hypotheses:
1421  1. if regular texture coordinate node supplied,
1422  - H1a: use same order and knot vector as control, or
1423  - H1b: use linear order=2 and knot vector
1424  - specify the texturecoordinate points as control
1425  2. if no texture coordinate node node supplied,
1426  H2a:
1427  - compute defaults using relative spatial distance between xyz control
1428  along u (row) and v (column) directions
1429  - apply #1
1430  H2b:
1431  - compute defaults using equal spacing
1432  along u (row) and v (column) directions
1433  - apply #1
1434  H2c:
1435  - leave texcoords blank and let stream_polyrep supply defaults
1436  3. if nurbstexturecoordinate node supplied,
1437  H3a:
1438  - set texture surface control to 0 0 1 1
1439  - use linear order 2
1440  H3b:
1441  - do H2a or H2b
1442  Both:
1443  - interpret the texturecallback points as uv
1444  - use uv to lookup st using piegl surface interpolator
1445  on separate surface in nurbstexturecoordinate node
1446  options: a) in texture callback b) when converting to polyrep
1447  */
1448  texcoordnodeIsGenerated = FALSE; //for possible garbage collection
1449  texCoordNode = node->texCoord;
1450  if(!texCoordNode){
1451  //2 no texcoord node supplied
1452  switch('b'){
1453  case 'a':
1454  //H2a compute uniform defaults using relative spatial distance of xyz in u and v
1455  break;
1456  case 'b':
1457  // H2b: - compute defaults using equal spacing
1458  {
1459  float du, dv, uu, vv;
1460  int jj,j,k;
1461  struct X3D_TextureCoordinate *texCoord = createNewX3DNode0(NODE_TextureCoordinate);
1462  texcoordnodeIsGenerated = TRUE;
1463  texCoord->point.p = MALLOC(struct SFVec2f*,nu * nv * sizeof(struct SFVec2f));
1464  du = 1.0f / (float)max(1,(nu -1));
1465  dv = 1.0f / (float)max(1,(nv -1));
1466  vv = 0.0f;
1467  jj = 0;
1468  for(k=0;k<nv;k++){
1469  if(k == nv-1) vv = 1.0f; //they like end exact on 1.0f
1470  uu = 0.0f;
1471  for(j=0;j<nu;j++){
1472  if(j == nu-1) uu = 1.0f; //they like end exact on 1.0f
1473  texCoord->point.p[jj].c[0] = uu;
1474  texCoord->point.p[jj].c[1] = vv;
1475  uu += du;
1476  jj++;
1477  }
1478  vv += dv;
1479  }
1480  texCoord->point.n = jj;
1481  texCoordNode = X3D_NODE(texCoord);
1482  }
1483  break;
1484  default:
1485  //H2c skip - nada
1486  break;
1487  }
1488  }
1489  if(texCoordNode){
1490  //USETEXCOORD = TRUE
1491  if(texCoordNode->_nodeType == NODE_TextureCoordinate){
1492  //TEXCOORDTYPE = 1 //interpret nurbsurftexcoordcb coords as texture coords
1493  struct X3D_TextureCoordinate *texCoord = (struct X3D_TextureCoordinate *)texCoordNode;
1494  float *control2D = (float*)texCoord->point.p;
1495  int nctrl = texCoord->point.n;
1496  if(1){
1497  //H1a: use same order and knot vector as controlPoints
1498  // except using texcoord.point as nurbscontrol
1499  //CONFIRMED when using H2b 'b' above to generate missing TexCoord node
1500  gluNurbsSurface(theNurb,nku,knotsu,nkv,knotsv,2,2*nu,control2D,node->uOrder,node->vOrder,GL_MAP2_TEXTURE_COORD_2);
1501  }else{
1502  //H1b: use linear order=2 and knot vectors, order from main surface
1503  //DISCONFIRMED: the texturecoord callback is never called
1504  float *tknotsu, *tknotsv;
1505  int k;
1506 
1507  tknotsu = MALLOC(float *, (nu+2) *sizeof(GLfloat));
1508  tknotsv = MALLOC(float *, (nv+2) *sizeof(GLfloat));
1509  generateUniformKnotVector(2,nu,tknotsu);
1510  generateUniformKnotVector(2,nv,tknotsv);
1511  printf("tknotsu = [");
1512  for(k=0;k<nu+2;k++) printf("%f ",tknotsu[k]);
1513  printf("]\ntknotsv = [");
1514  for(k=0;k<nv+2;k++) printf("%f ",tknotsv[k]);
1515  printf("]\n");
1516  gluNurbsSurface(theNurb,nu+2,knotsu,nv+2,knotsv,4,4*nu,control2D,2,2,GL_MAP2_TEXTURE_COORD_2);
1517  }
1518  } else if(texCoordNode->_nodeType == NODE_NurbsTextureCoordinate){
1519  //TEXCOORDTYPE = 2 /interpret nurbsurftexcoordcb coords as uv coords,
1520  // use to lookup st via piegl interpolation of NurbsTextureCoordinate nurbs surface
1521  if(0){
1522  //H3a:
1523  //- set texture surface control to 0 0 1 1
1524  //- use linear order 2
1525  float tknots[4] = {0.0f, 0.0f, 1.0f, 1.0f};
1526  float unit_control2D [8] = {0.0f, 0.0f, 1.0f, 0.0f, 1.0f, 1.0f, 0.0f, 1.0f};
1527  struct X3D_NurbsTextureCoordinate *texCoord = (struct X3D_NurbsTextureCoordinate *)texCoordNode;
1528  gluNurbsSurface(theNurb,4,tknots,4,tknots,2,2*2,unit_control2D,2,2,GL_MAP2_TEXTURE_COORD_2);
1529  }else{
1530  //H3b same as H2a or H2b
1531  // H2b: - compute defaults using equal spacing
1532  {
1533  float du, dv, uu, vv;
1534  int jj,j,k;
1535  float *control2D;
1536  struct X3D_TextureCoordinate *texCoord = createNewX3DNode0(NODE_TextureCoordinate);
1537  texcoordnodeIsGenerated = TRUE;
1538  FREE_IF_NZ(texCoord->point.p);
1539  texCoord->point.p = MALLOC(struct SFVec2f*,nu * nv * sizeof(struct SFVec2f));
1540  du = 1.0f / (float)max(1,(nu -1));
1541  dv = 1.0f / (float)max(1,(nv -1));
1542  vv = 0.0f;
1543  jj = 0;
1544  for(k=0;k<nv;k++){
1545  if(k == nv-1) vv = 1.0f; //they like end exact on 1.0f
1546  uu = 0.0f;
1547  for(j=0;j<nu;j++){
1548  if(j == nu-1) uu = 1.0f; //they like end exact on 1.0f
1549  texCoord->point.p[jj].c[0] = uu;
1550  texCoord->point.p[jj].c[1] = vv;
1551  uu += du;
1552  jj++;
1553  }
1554  vv += dv;
1555  }
1556  texCoord->point.n = jj;
1557  texCoordNode = X3D_NODE(texCoord);
1558  //H1a: use same order and knot vector as controlPoints
1559  // except using texcoord.point as nurbscontrol
1560  //CONFIRMED when using H2b 'b' above to generate missing TexCoord node
1561  control2D = (float*)texCoord->point.p;
1562  gluNurbsSurface(theNurb,nku,knotsu,nkv,knotsv,2,2*nu,control2D,node->uOrder,node->vOrder,GL_MAP2_TEXTURE_COORD_2);
1563 
1564  }
1565 
1566  }
1567 
1568  }
1569  }
1570  if(trim){
1571  int i;
1572 
1573  if(0){
1574  if(DEBG) printf("gluBeginTrim \n");
1575  gluBeginTrim (theNurb);
1576  //outside border H: scene author is responsible
1577  if(1){
1578  // counter clockwise, simple 4 corner uv from redbook sample
1579  GLfloat edgePt[5][2] = {{0.0, 0.0}, {1.0, 0.0}, {1.0, 1.0}, {0.0, 1.0}, {0.0, 0.0}};
1580  if(DEBG) printf("gluPwlCurve 0\n");
1581  gluPwlCurve (theNurb, 5, &edgePt[0][0], 2, GLU_MAP1_TRIM_2);
1582  }else{
1583  // 2 x (node.utessselation u + node.vtesselation v) + 1 edge values
1584  // except I get a nurbs error with the following
1585  GLfloat *edges = MALLOC(void *, (2 * ntessu + 2 * ntessv) *2*sizeof(GLfloat));
1586  GLfloat uspan, vspan;
1587  uspan = 1.0f/(float)(ntessu -1);
1588  vspan = 1.0f/(float)(ntessv -1);
1589  for(i=0;i<ntessu-1;i++){
1590  edges[i*2 +0] = (float)(i)*uspan;
1591  edges[i*2 +1] = 0.0;
1592  edges[(ntessu+ntessv+i)*2 + 0] = (float)(ntessu - 1 - i)*uspan;
1593  edges[(ntessu+ntessv+i)*2 + 1] = 1.0;
1594  }
1595  for(i=0;i<ntessv;i++){
1596  edges[(ntessu+i)*2 + 0] = 1.0;
1597  edges[(ntessu+i)*2 + 1] = (float)(i)*vspan;
1598  edges[(ntessu+ntessv+ntessu+i)*2 + 0] = 0.0;
1599  edges[(ntessu+ntessv+ntessu+i)*2 + 1] = (float)(ntessv - 1 - i)*vspan;
1600  }
1601  //close curve
1602  edges[((ntessu -1)*2 + (ntessv -1)*2)*2 + 0] = 0.0;
1603  edges[((ntessu -1)*2 + (ntessv -1)*2)*2 + 1] = 0.0;
1604  if(DEBG) printf("gluPwlCurve 1\n");
1605  gluPwlCurve (theNurb, 2*(ntessu -1 + ntessv -1) +1, edges, 2, GLU_MAP1_TRIM_2);
1606  }
1607  if(DEBG) printf("gluEndTrim\n");
1608  gluEndTrim (theNurb);
1609  }
1610 
1611  //interior cutouts
1612  if(0){
1613  //redbook example trim curves - these work
1614  GLfloat curvePt[4][2] = /* clockwise */
1615  {{0.25, 0.5}, {0.25, 0.75}, {0.75, 0.75}, {0.75, 0.5}};
1616  GLfloat curveKnots[8] =
1617  {0.0, 0.0, 0.0, 0.0, 1.0, 1.0, 1.0, 1.0};
1618  GLfloat pwlPt[4][2] = /* clockwise */
1619  {{0.75, 0.5}, {0.5, 0.25}, {0.25, 0.5}};
1620 
1621  if(DEBG) printf("gluBeginTrim A\n");
1622  gluBeginTrim (theNurb);
1623  if(DEBG) printf("gluNurbsCurve A\n");
1624  gluNurbsCurve (theNurb, 8, curveKnots, 2,
1625  &curvePt[0][0], 4, GLU_MAP1_TRIM_2);
1626  if(DEBG) printf("gluPwlCurve A\n");
1627  gluPwlCurve (theNurb, 3, &pwlPt[0][0], 2, GLU_MAP1_TRIM_2);
1628  if(DEBG) printf("gluEndTrim A\n");
1629  gluEndTrim (theNurb);
1630  }
1631  if(1)
1632  for(i=0;i<trim->n;i++){
1633  int m;
1634  struct X3D_Contour2D * tc = (struct X3D_Contour2D *)trim->p[i];
1635  if(DEBG) printf("gluBeginTrim B\n");
1636  gluBeginTrim (theNurb);
1637  for(m=0;m<tc->children.n;m++)
1638  {
1639  int j,k,nk,dim;
1640  struct X3D_ContourPolyline2D *cp2d;
1641  struct X3D_NurbsCurve2D *nc2d;
1642  struct X3D_Node *ctr = tc->children.p[m]; //trim->p[i];
1643  GLfloat *cknot, *ctrl, *cweight;
1644  cknot = ctrl = cweight = NULL;
1645  switch(ctr->_nodeType){
1646  case NODE_ContourPolyline2D:
1647  cp2d = (struct X3D_ContourPolyline2D *)ctr;
1648  ctrl = MALLOC(void *, cp2d->controlPoint.n * 2*sizeof(GLfloat));
1649  for(j=0;j<cp2d->controlPoint.n;j++) {
1650  for(k=0;k<2;k++)
1651  ctrl[j*2 + k] = (float)cp2d->controlPoint.p[j].c[k];
1652  }
1653  if(DEBG) printf("gluPwlCurve B\n");
1654  gluPwlCurve (theNurb, cp2d->controlPoint.n, ctrl, 2, GLU_MAP1_TRIM_2);
1655  FREE_IF_NZ(ctrl);
1656  break;
1657  case NODE_NurbsCurve2D:
1658  nc2d = (struct X3D_NurbsCurve2D *)ctr;
1659  dim = 2;
1660  nk = nc2d->controlPoint.n + nc2d->order;
1661  if(nk == nc2d->knot.n)
1662 
1663  cknot = MALLOC(void *, nk * sizeof(GLfloat));
1664  if(nc2d->weight.n){ // == 3){
1665  dim = 3;
1666  cweight = MALLOC(void *, nc2d->controlPoint.n * sizeof(GLfloat));
1667  if(nc2d->weight.n == nc2d->controlPoint.n){
1668  for(j=0;j<nc2d->weight.n;j++) cweight[j] = (float)nc2d->weight.p[j];
1669  }else{
1670  for(j=0;j<nc2d->controlPoint.n;j++) cweight[j] = 1.0f;
1671  }
1672  }
1673  ctrl = MALLOC(void *, nc2d->controlPoint.n * dim*sizeof(GLfloat));
1674  for(j=0;j<nc2d->controlPoint.n;j++) {
1675  for(k=0;k<2;k++){
1676  ctrl[j*dim + k] = (float)nc2d->controlPoint.p[j].c[k];
1677  if(dim == 3) ctrl[j*dim + k] *= cweight[j];
1678  }
1679  if(dim == 3) ctrl[j*dim + dim-1] = (float)cweight[j];
1680  }
1681  if(knotsOK(nc2d->order,nc2d->controlPoint.n,nc2d->knot.n,nc2d->knot.p)){
1682  //if(nk == nc2d->knot.n){
1683  for(j=0;j<nc2d->knot.n;j++)
1684  cknot[j] = (float)nc2d->knot.p[j];
1685  }else{
1686  generateUniformKnotVector(nc2d->order,nc2d->controlPoint.n,cknot);
1687  printf("replacing nurbscurve2D knotvector with:\n");
1688  for(j=0;j<nk;j++){
1689  printf("%f ",cknot[j]);
1690  }
1691  printf("\n");
1692 
1693  }
1694  if(DEBGC) {
1695  printf("knot %d ={",nc2d->knot.n);
1696  for(j=0;j<nc2d->knot.n;j++){
1697  printf("%f ",cknot[j]);
1698  }
1699  printf("}\n");
1700  printf("control %d = {\n",nc2d->controlPoint.n);
1701  for(j=0;j<nc2d->controlPoint.n;j++) {
1702  for(k=0;k<dim;k++) printf("%f \n",ctrl[j*dim +k]);
1703  printf("\n");
1704  }
1705  printf("}\n");
1706  }
1707  if(DEBG) printf("gluNurbsCurve B\n");
1708  if(1){
1709  int mtess, ntess;
1710  mtess = nc2d->order + 1;
1711  ntess = nc2d->tessellation;
1712  if(0){
1713  //chord length or automatic - not implemented properly nor tested thoroughly
1714  //if you do chord length in pixels, you need to manually pass in sampling matrices
1715  //and somehow you need to trigger a recompile: another call to this compile_
1716  // as avatar/viewer moves closer (father) from the nurbs node
1717  double model[16], proj[16];
1718  float modelf[16], projf[16];
1719  int viewPort[10];
1720  if(ntess > 0)
1721  mtess = ntess;
1722  else if(ntess < 0)
1723  mtess = -ntess;
1724 
1725  gluNurbsProperty(theNurb, GLU_SAMPLING_TOLERANCE, (float)(mtess)); //25.0);
1726  if(ntess < 0)
1727  gluNurbsProperty(theNurb,GLU_SAMPLING_METHOD,GLU_PATH_LENGTH); //pixels, the default
1728  else
1729  gluNurbsProperty(theNurb,GLU_SAMPLING_METHOD,GLU_PARAMETRIC_TOLERANCE);
1730  gluNurbsProperty(theNurb, GLU_AUTO_LOAD_MATRIX,GL_FALSE);
1731 
1732  FW_GL_GETDOUBLEV(GL_MODELVIEW_MATRIX, model);
1733  FW_GL_GETDOUBLEV(GL_PROJECTION_MATRIX, proj);
1734  FW_GL_GETINTEGERV(GL_VIEWPORT, viewPort);
1735  for(i=0;i<16;i++){
1736  modelf[i] = (float)model[i];
1737  projf[i] = (float)proj[i];
1738  }
1739  gluLoadSamplingMatrices(theNurb,modelf,projf,viewPort);
1740  }
1741  if(1){
1742  //uniform spacing of sampling points in u (or uv) parameter space - works
1743  //node must specify tesselation value. see specs for interpretation
1744  // http://www.web3d.org/documents/specifications/19775-1/V3.3/Part01/components/nurbs.html#NurbsCurve
1745  if(ntess > 0)
1746  mtess = max(mtess,ntess+1);
1747  else if(ntess < 0)
1748  mtess = max(mtess,(-ntess * n) + 1);
1749  else
1750  mtess = max(mtess,2*n + 1);
1751  gluNurbsProperty(theNurb,GLU_SAMPLING_METHOD,GLU_DOMAIN_DISTANCE);
1752  gluNurbsProperty(theNurb,GLU_U_STEP,(GLfloat)mtess);
1753  }
1754  }
1755  gluNurbsCurve (theNurb, nk, cknot, dim, ctrl, nc2d->order, GLU_MAP1_TRIM_2);
1756  break;
1757  default:
1758  ConsoleMessage("%s %d","unknown trimming contour node",ctr->_nodeType);
1759  }
1760  FREE_IF_NZ(ctrl);
1761  FREE_IF_NZ(cknot);
1762  FREE_IF_NZ(cweight);
1763  }
1764  if(DEBG) printf("gluEndTrim B\n");
1765  gluEndTrim (theNurb);
1766  }
1767  }
1768  if(DEBG) printf("gluEndSurface \n");
1769  gluEndSurface(theNurb);
1770  if(DEBG) printf("gluDeleteNurbsRenderer \n");
1771  if(0) gluDeleteNurbsRenderer(theNurb);
1772 
1773  //convert points to polyrep
1774  convert_strips_to_polyrep(strips,(struct X3D_NurbsTrimmedSurface*)node);
1775 
1776  if(texcoordnodeIsGenerated){
1777  struct X3D_TextureCoordinate *texCoord = (struct X3D_TextureCoordinate *)texCoordNode;
1778  FREE_IF_NZ(texCoord->point.p);
1779  FREE_IF_NZ(texCoordNode);
1780  }
1781  //vector_releaseData(,strips);
1782  if(strips){
1783  for(i=0;i<vectorSize(strips);i++){
1784  struct stripState *ss = vector_get_ptr(struct stripState,strips,i);
1785  FREE_IF_NZ(ss->pv.data); ss->pv.allocn = 0; ss->pv.n = 0;
1786  FREE_IF_NZ(ss->nv.data); ss->nv.allocn = 0; ss->nv.n = 0;
1787  FREE_IF_NZ(ss->tv.data); ss->tv.allocn = 0; ss->tv.n = 0;
1788  }
1789  FREE_IF_NZ(strips->data); strips->allocn = 0; strips->n = 0;
1790  FREE_IF_NZ(strips);
1791  }
1792  FREE_IF_NZ(xyzw);
1793  }
1794  FREE_IF_NZ(knotsv);
1795  FREE_IF_NZ(knotsu);
1796  }
1797 
1798 }
1799 
1800 void render_ray_polyrep(void *node);
1801 void collide_genericfaceset(void *node);
1802 void render_polyrep(void *node);
1803 
1804 void compile_NurbsPatchSurface(struct X3D_NurbsPatchSurface *node){
1805  MARK_NODE_COMPILED
1806  compile_NurbsSurface(node, NULL);
1807 }
1808 void rendray_NurbsPatchSurface (struct X3D_NurbsPatchSurface *node) {
1809  COMPILE_IF_REQUIRED
1810  if (!node->_intern) return;
1811  render_ray_polyrep(node);
1812 }
1813 
1814 void collide_NurbsPatchSurface (struct X3D_NurbsPatchSurface *node) {
1815  COMPILE_IF_REQUIRED
1816  if (!node->_intern) return;
1817  collide_genericfaceset(node);
1818 }
1819 
1820 void render_NurbsPatchSurface (struct X3D_NurbsPatchSurface *node) {
1821  COMPILE_IF_REQUIRED
1822  if (!node->_intern) return;
1823  CULL_FACE(node->solid)
1824  render_polyrep(node);
1825 }
1826 
1827 void compile_NurbsTrimmedSurface(struct X3D_NurbsTrimmedSurface *node){
1828  MARK_NODE_COMPILED
1829  compile_NurbsSurface((struct X3D_NurbsPatchSurface *)node, &node->trimmingContour);
1830 }
1831 
1832 void rendray_NurbsTrimmedSurface (struct X3D_NurbsTrimmedSurface *node) {
1833  COMPILE_IF_REQUIRED
1834  if (!node->_intern) return;
1835  render_ray_polyrep(node);
1836 }
1837 
1838 void collide_NurbsTrimmedSurface (struct X3D_NurbsTrimmedSurface *node) {
1839  COMPILE_IF_REQUIRED
1840  if (!node->_intern) return;
1841  collide_genericfaceset(node);
1842 }
1843 
1844 void render_NurbsTrimmedSurface (struct X3D_NurbsTrimmedSurface *node) {
1845  COMPILE_IF_REQUIRED
1846  if (!node->_intern) return;
1847  CULL_FACE(node->solid)
1848  render_polyrep(node);
1849 }
1850 
1851 void do_NurbsPositionInterpolator (void *node) {
1852  struct X3D_NurbsPositionInterpolator *px;
1853  struct X3D_Coordinate *control;
1854  float fraction; //knotrange[2],
1855  double *weight;
1856  float cw[4]; //*xyzw,
1857 
1858  if (!node) return;
1859  px = (struct X3D_NurbsPositionInterpolator *) PPX(node);
1860  if(NNC(px)){
1861  float *knots;
1862  int i,nk, n;
1863 
1864  MNC(px);
1865  if(!px->_OK == TRUE){
1866  px->_OK = FALSE;
1867  if(!(px->knot.n > 1))
1868  return;
1869  if(!px->controlPoint )
1870  return;
1871  if(px->controlPoint->_nodeType != NODE_Coordinate)
1872  return;
1873  control = (struct X3D_Coordinate *) px->controlPoint;
1874  if(control->point.n < 2)
1875  return;
1876  n = control->point.n;
1877  weight = NULL;
1878  if(px->weight.n == n)
1879  weight = px->weight.p;
1880  px->_xyzw.p = MALLOC(struct SFVec4f*,control->point.n * sizeof(struct SFVec4f));
1881  px->_xyzw.n = n;
1882  for(i=0;i<control->point.n;i++){
1883  float xyzw[4], wt;
1884  wt = weight ? (float)weight[i] : 1.0f;
1885  veccopy3f(xyzw,control->point.p[i].c);
1886  vecscale3f(xyzw,xyzw,wt);
1887  xyzw[3] = wt;
1888  veccopy4f(px->_xyzw.p[i].c,xyzw);
1889  }
1890  if(knotsOK(px->order,px->_xyzw.n,px->knot.n,px->knot.p)){
1891 
1892  nk = px->knot.n;
1893  knots = MALLOC(void *, nk * sizeof(GLfloat));
1894  for(i=0;i<nk;i++){
1895  knots[i] = (GLfloat)px->knot.p[i];
1896  }
1897  //printf("good knot nk=%d\n",nk);
1898  //for(int ii=0;ii<nk;ii++)
1899  // printf("[%d]=%f \n",ii,knots[ii]);
1900 
1901  }else{
1902  static int once = 0;
1903  //generate uniform knot vector
1904  nk = n + px->order ;
1905  //caller: please malloc knots = malloc( (ncontrol + order ) * sizeof(float))
1906  knots = MALLOC(void *, nk *sizeof(GLfloat));
1907  generateUniformKnotVector(px->order,n, knots);
1908  if(!once){
1909  int ii;
1910  printf("bad knot vector, replacing with:\n");
1911  for(ii=0;ii<nk;ii++)
1912  printf("[%d]=%f \n",ii,knots[ii]);
1913  once = 1;
1914  }
1915  //nk = 0;
1916  }
1917  px->_knot.p = knots;
1918  px->_knot.n = nk;
1919  px->_knotrange.c[0] = px->_knot.p[0];
1920  px->_knotrange.c[1] = px->_knot.p[px->_knot.n-1];
1921 
1922  px->_OK = TRUE;
1923  }
1924  }
1925  if(!px->_OK)
1926  return;
1927 
1928  fraction = max(px->_knotrange.c[0],px->set_fraction);
1929  fraction = min(px->_knotrange.c[1],px->set_fraction);
1930  CurvePoint(px->_xyzw.n, px->order-1, px->_knot.p, (float*)px->_xyzw.p, fraction, cw );
1931  veccopy3f(px->value_changed.c,cw);
1932  MARK_EVENT (node, offsetof (struct X3D_NurbsPositionInterpolator, value_changed));
1933 
1934  #ifdef SEVERBOSE
1935  printf("do_PositionInt: Position/Vec3f interp, node %u kin %d kvin %d set_fraction %f\n",
1936  node, kin, kvin, px->set_fraction);
1937  #endif
1938 
1939 
1940 
1941 }
1942 
1943 /* NurbsOrientationInterpolator
1944  Called during the "events_processed" section of the event loop,
1945  so this is called ONLY when there is something required to do, thus
1946  there is no need to look at whether it is active or not
1947  */
1948 void do_NurbsOrientationInterpolator (void *node) {
1950  struct X3D_Coordinate *control;
1951  float fraction; //knotrange[2],
1952  double *weight;
1953  //float *xyzw; //, cw[4];
1954 
1955  if (!node) return;
1956  px = (struct X3D_NurbsOrientationInterpolator *) PPX(node);
1957  if(NNC(px)){
1958  float *knots;
1959  int i,nk, n;
1960 
1961  MNC(px);
1962  if(!px->_OK == TRUE){
1963  px->_OK = FALSE;
1964  if(!(px->knot.n > 1))
1965  return;
1966  if(!px->controlPoint )
1967  return;
1968  if(px->controlPoint->_nodeType != NODE_Coordinate)
1969  return;
1970  control = (struct X3D_Coordinate *) px->controlPoint;
1971  if(control->point.n < 2)
1972  return;
1973  n = control->point.n;
1974  weight = NULL;
1975  if(px->weight.n == n)
1976  weight = px->weight.p;
1977  px->_xyzw.p = MALLOC(struct SFVec4f*,control->point.n * sizeof(struct SFVec4f));
1978  px->_xyzw.n = n;
1979  for(i=0;i<control->point.n;i++){
1980  float xyzw[4], wt;
1981  wt = weight ? (float)weight[i] : 1.0f;
1982  veccopy3f(xyzw,control->point.p[i].c);
1983  vecscale3f(xyzw,xyzw,wt);
1984  xyzw[3] = wt;
1985  veccopy4f(px->_xyzw.p[i].c,xyzw);
1986  }
1987  if(knotsOK(px->order,px->_xyzw.n,px->knot.n,px->knot.p)){
1988 
1989  nk = px->knot.n;
1990  knots = MALLOC(void *, nk * sizeof(GLfloat));
1991  for(i=0;i<nk;i++){
1992  knots[i] = (GLfloat)px->knot.p[i];
1993  }
1994  //printf("good knot nk=%d\n",nk);
1995  //for(int ii=0;ii<nk;ii++)
1996  // printf("[%d]=%f \n",ii,knots[ii]);
1997 
1998  }else{
1999  static int once = 0;
2000  //generate uniform knot vector
2001  nk = n + px->order ;
2002  //caller: please malloc knots = malloc( (ncontrol + order ) * sizeof(float))
2003  knots = MALLOC(void *, nk *sizeof(GLfloat));
2004  generateUniformKnotVector(px->order,n, knots);
2005  if(!once){
2006  int ii;
2007  printf("bad knot vector, replacing with:\n");
2008  for(ii=0;ii<nk;ii++)
2009  printf("[%d]=%f \n",ii,knots[ii]);
2010  once = 1;
2011  }
2012  //nk = 0;
2013  }
2014  px->_knot.p = knots;
2015  px->_knot.n = nk;
2016  px->_knotrange.c[0] = px->_knot.p[0];
2017  px->_knotrange.c[1] = px->_knot.p[px->_knot.n-1];
2018 
2019  px->_OK = TRUE;
2020  }
2021  }
2022  if(!px->_OK)
2023  return;
2024 
2025  fraction = max(px->_knotrange.c[0],px->set_fraction);
2026  fraction = min(px->_knotrange.c[1],px->set_fraction);
2027  if(1){
2028  //DELTA METHOD: instead of using a piegl formula for curve derivitive,
2029  //we sample 2 points near the u of interest, and take the difference
2030  //to get a slope vector
2031  float f1, f2, cw1[4], cw2[4], dir[3], rot[4];
2032 
2033 
2034  f1 = fraction;
2035  f2 = fraction + .01f;
2036  if(f2 > 1.0f){
2037  f1 = fraction - .01f;
2038  f2 = fraction;
2039  }
2040 
2041  CurvePoint(px->_xyzw.n, px->order-1, px->_knot.p, (float*)px->_xyzw.p, f1, cw1 );
2042  CurvePoint(px->_xyzw.n, px->order-1, px->_knot.p, (float*)px->_xyzw.p, f2, cw2 );
2043  vecdif3f(dir,cw2,cw1);
2044  vecnormalize3f(dir,dir);
2045 
2046  if(1){
2047  //1-direction vector relative to x-axis
2048  //now have 2 direction vectors - one at start, one at end
2049  //do a difference to get relative rotation from start
2050  float perp[3], dirx[3], sine;
2051  memset(dirx,0,3*sizeof(float));
2052  dirx[0] = 1.0f;
2053  veccross3f(perp,dirx,dir);
2054  sine = veclength3f(perp);
2055  if(sine == 0.0f){
2056  //no net rotation from start
2057  float default_direction [4] = {1.0f, 0.0f, 0.0f, 0.0f};
2058  veccopy4f(rot,default_direction);
2059  }else{
2060  float cosine, angle;
2061  vecnormalize3f(perp,perp);
2062  cosine = vecdot3f(dir,dirx);
2063  angle = atan2f(sine,cosine);
2064  veccopy3f(rot,perp);
2065  rot[3] = angle;
2066  }
2067  }else if(0){
2068  //2-direction vector to axis angle method
2069  //now have 2 direction vectors - one at start, one at end
2070  //do a difference to get relative rotation from start
2071  float perp[3],dir0[3];
2072 
2073  CurvePoint(px->_xyzw.n, px->order-1, px->_knot.p, (float*)px->_xyzw.p, 0.0f, cw1 );
2074  CurvePoint(px->_xyzw.n, px->order-1, px->_knot.p, (float*)px->_xyzw.p, .001f, cw2 );
2075  vecdif3f(dir0,cw2,cw1);
2076  vecnormalize3f(dir0,dir0);
2077 
2078  veccross3f(perp,dir,dir0);
2079  if(veclength3f(perp) == 0.0f){
2080  //no net rotation from start
2081  float default_direction [4] = {1.0f, 0.0f, 0.0f, 0.0f};
2082  veccopy4f(rot,default_direction);
2083  printf(".");
2084  }else{
2085  float cosine,angle;
2086  vecnormalize3f(perp,perp);
2087  cosine = vecdot3f(dir0,dir);
2088  angle = acosf(cosine);
2089  veccopy3f(rot,perp);
2090  rot[3] = -angle;
2091  }
2092  }
2093  veccopy4f(px->value_changed.c,rot);
2094  }else{
2095  float default_direction [4] = {1.0f, 0.0f, 0.0f, 0.0f};
2096  veccopy4f(px->value_changed.c,default_direction);
2097  }
2098  MARK_EVENT (node, offsetof (struct X3D_NurbsOrientationInterpolator, value_changed));
2099 
2100  #ifdef SEVERBOSE
2101  printf("do_NurbsOrientInt: set_fraction %f value_changed %f %f %f %f\n",fraction,
2102  px->value_changed.c[0],px->value_changed.c[1],px->value_changed.c[2],px->value_changed.c[3] );
2103  #endif
2104 
2105 }
2106 
2107 void do_NurbsSurfaceInterpolator (void *_node) {
2108  float uv[2], xyzw[4];
2109  int ok;
2110  struct X3D_NurbsSurfaceInterpolator *node;
2111 
2112  if (!_node) return;
2113  node = (struct X3D_NurbsSurfaceInterpolator *) _node;
2114  if(node->_nodeType != NODE_NurbsSurfaceInterpolator) return;
2115 
2116  if(NNC(node)){
2117  MNC(node);
2118  if(node->_OK != TRUE){
2119  int i,j, n, nu, nv, nku, nkv;
2120  GLfloat *xyzw, *knotsu, *knotsv;
2121 
2122  nku = nkv = nu = nv = n = 0;
2123  xyzw = knotsu = knotsv = NULL;
2124  if(node->controlPoint){
2125  if(node->controlPoint->_nodeType == NODE_CoordinateDouble){
2126  struct Multi_Vec3d *mfd;
2127  mfd = &((struct X3D_CoordinateDouble *)(node->controlPoint))->point;
2128  n = mfd->n;
2129  xyzw = MALLOC(void *, n * 4 * sizeof(GLfloat));
2130  for(i=0;i<mfd->n;i++){
2131  for(j=0;j<3;j++){
2132  xyzw[i*4 + j] = (float)mfd->p[i].c[j];
2133  }
2134  }
2135  }else if(node->controlPoint->_nodeType == NODE_Coordinate){
2136  struct Multi_Vec3f *mff;
2137  mff = &((struct X3D_Coordinate *)(node->controlPoint))->point;
2138  n = mff->n;
2139  xyzw = MALLOC(void *, n * 4 * sizeof(GLfloat));
2140  for(i=0;i<mff->n;i++){
2141  for(j=0;j<3;j++){
2142  xyzw[i*4 + j] = mff->p[i].c[j];
2143  }
2144  }
2145  }
2146  }else{
2147  n = 0;
2148  }
2149 
2150  if(node->weight.n && node->weight.n == n){
2151  double w;
2152  int m,im;
2153  m = min(node->weight.n, n);
2154  for(i=0;i<n;i++){
2155  im = i < m ? i : m-1;
2156  w = node->weight.p[im];
2157  xyzw[i*4 + 3] = (float)w;
2158  }
2159  }else{
2160  for(i=0;i<n;i++) xyzw[i*4 + 3] = 1.0;
2161  }
2162  nu = node->uDimension;
2163  nv = node->vDimension;
2164  if(nu * nv != n) return;
2165  if(nu < node->uOrder) return;
2166  if(nv < node->vOrder) return;
2167  //int knotsOK(int order, int ncontrol, int nknots, double *knots)
2168  //if(node->uKnot.n && node->uKnot.n == nu + node->uOrder ){
2169  if(knotsOK(node->uOrder,nu,node->uKnot.n,node->uKnot.p)){
2170  //could do another check: max number of consecutive equal value knots == order
2171  //could do another check: knot values == or ascending
2172  nku = node->uKnot.n;
2173  knotsu = MALLOC(void *, nku * sizeof(GLfloat));
2174  for(i=0;i<nku;i++){
2175  knotsu[i] = (GLfloat)node->uKnot.p[i];
2176  }
2177  if(DEBG){
2178  int ii;
2179  printf("good u knot vector nk=%d\n",nku);
2180  for(ii=0;ii<nku;ii++)
2181  printf("[%d]=%f \n",ii,knotsu[ii]);
2182  }
2183 
2184  }else{
2185  //generate uniform knot vector
2186  static int once = 0;
2187  nku = nu + node->uOrder ;
2188  //caller: please malloc knots = MALLOC(void *, (ncontrol + order ) * sizeof(float))
2189  knotsu = MALLOC(void *, nku *sizeof(GLfloat));
2190  generateUniformKnotVector(node->uOrder,nu, knotsu);
2191  if(!once){
2192  int ii;
2193  printf("bad u knot vector given, replacing with:\n");
2194  for(ii=0;ii<nku;ii++)
2195  printf("[%d]=%f \n",ii,knotsu[ii]);
2196  once = 1;
2197  }
2198  //nk = 0;
2199  }
2200 
2201  if(knotsOK(node->vOrder,nv,node->vKnot.n,node->vKnot.p)){
2202  //if(node->vKnot.n && node->vKnot.n == nv + node->vOrder ){
2203  nkv = node->vKnot.n;
2204  knotsv = MALLOC(void *, nkv * sizeof(GLfloat));
2205  for(i=0;i<nkv;i++){
2206  knotsv[i] = (GLfloat)node->vKnot.p[i];
2207  }
2208  if(DEBG){
2209  int ii;
2210  printf("good v knot vector nk=%d\n",nkv);
2211  for(ii=0;ii<nkv;ii++)
2212  printf("[%d]=%f \n",ii,knotsv[ii]);
2213  }
2214 
2215  }else{
2216  static int once = 0;
2217  //generate uniform knot vector
2218  nkv = nv + node->vOrder ;
2219  //caller: please malloc knots = MALLOC(void *, (ncontrol + order ) * sizeof(float))
2220  knotsv = MALLOC(void *, nkv *sizeof(GLfloat));
2221  generateUniformKnotVector(node->vOrder,nv, knotsv);
2222  if(!once){
2223  int ii;
2224  printf("bad v knot vector given, replacing with:\n");
2225  for(ii=0;ii<nkv;ii++)
2226  printf("[%d]=%f \n",ii,knotsv[ii]);
2227  once = 1;
2228  }
2229  //nk = 0;
2230  }
2231  node->_controlPoint.p = (struct SFVec4f*)xyzw;
2232  node->_controlPoint.n = n;
2233  node->_uKnot.p = knotsu;
2234  node->_uKnot.n = nku;
2235  node->_vKnot.p = knotsv;
2236  node->_vKnot.n = nkv;
2237  node->_OK = TRUE;
2238  }
2239 
2240  }
2241 
2242  if(!node->_OK) return;
2243 
2244  veccopy2f(uv,node->set_fraction.c);
2245 
2246  ok = SurfacePoint(node->uDimension,node->uOrder-1,node->_uKnot.p,
2247  node->vDimension,node->vOrder-1,node->_vKnot.p,
2248  (float *)node->_controlPoint.p,uv[1],uv[0],xyzw);
2249 
2250  veccopy3f(node->position_changed.c,xyzw);
2251  if(1){
2252  //DELTA method to get normal
2253  //u direction
2254  float udir[3], vdir[3], normal[3], xyz1[4];
2255  ok = SurfacePoint(node->uDimension,node->uOrder-1,node->_uKnot.p,
2256  node->vDimension,node->vOrder-1,node->_vKnot.p,
2257  (float *)node->_controlPoint.p,uv[1],uv[0]+.01f,xyz1);
2258  vecdif3f(udir,xyz1,xyzw);
2259  ok = SurfacePoint(node->uDimension,node->uOrder-1,node->_uKnot.p,
2260  node->vDimension,node->vOrder-1,node->_vKnot.p,
2261  (float *)node->_controlPoint.p,uv[1]+.01f,uv[0],xyz1);
2262  vecdif3f(vdir,xyz1,xyzw);
2263  veccross3f(normal,udir,vdir);
2264  vecnormalize3f(normal,normal);
2265  veccopy3f(node->normal_changed.c,normal);
2266 
2267 
2268  }
2269 
2270  MARK_EVENT (_node, offsetof (struct X3D_NurbsSurfaceInterpolator, position_changed));
2271  MARK_EVENT (_node, offsetof (struct X3D_NurbsSurfaceInterpolator, normal_changed));
2272 
2273  #ifdef SEVERBOSE
2274  printf ("Pos/Col, new value (%f %f %f)\n",
2275  px->value_changed.c[0],px->value_changed.c[1],px->value_changed.c[2]);
2276  #endif
2277 
2278 }
2279 
2280 //SWUNG
2281 
2282 void compile_NurbsSwungSurface(struct X3D_NurbsSwungSurface *node){
2283  struct X3D_NurbsPatchSurface *patch;
2284  struct X3D_Coordinate *controlPoint;
2285  struct X3D_NurbsCurve2D *trajectoryxz;
2286  struct X3D_NurbsCurve2D *profileyz;
2287  int nt, np,ic,j,i,k;
2288  double *xyzp, *xyzt;
2289  float *xyz;
2290 
2291  MARK_NODE_COMPILED
2292  //strategy: generate 3D control net from curves,
2293  // then delegate to NurbsPatchSurface
2294  //Swung:
2295  patch = (struct X3D_NurbsPatchSurface*) node->_patch;
2296  controlPoint = (struct X3D_Coordinate*)patch->controlPoint;
2297  if(!patch){
2298  patch = (struct X3D_NurbsPatchSurface*)createNewX3DNode(NODE_NurbsPatchSurface);
2299  controlPoint = (struct X3D_Coordinate*)createNewX3DNode(NODE_Coordinate);
2300  node->_patch = X3D_NODE(patch);
2301  patch->controlPoint = X3D_NODE(controlPoint);
2302  }
2303 
2304  trajectoryxz = (struct X3D_NurbsCurve2D *)node->trajectoryCurve;
2305  profileyz = (struct X3D_NurbsCurve2D *)node->profileCurve;
2306 
2307  nt = trajectoryxz->controlPoint.n;
2308  np = profileyz->controlPoint.n;
2309  xyzp = (double*)profileyz->controlPoint.p;
2310  xyzt = (double*)trajectoryxz->controlPoint.p;
2311  xyz = MALLOC(float*,nt * np * 3 * sizeof(float));
2312  controlPoint->point.p = (struct SFVec3f*)xyz;
2313  controlPoint->point.n = nt * np;
2314  ic = 0;
2315  for(j=0;j<nt;j++){
2316  float pt[3];
2317  double2float(pt,&xyzt[j*2],2);
2318  for(i=0;i<np;i++){
2319  float pp[3];
2320  float cosine, sine, swingangle;
2321  double2float(pp,&xyzp[2*i],2);
2322  swingangle = atan2f(pt[1],pt[0]);
2323  cosine = cosf(swingangle);
2324  sine = sinf(swingangle);
2325  xyz[ic*3 + 0] = pt[0] + cosine * pp[0];
2326  xyz[ic*3 + 1] = pp[1];
2327  xyz[ic*3 + 2] = pt[1] + sine * pp[0];
2328  ic++;
2329  }
2330  }
2331  patch->solid = node->solid;
2332  //u will be profile,
2333  patch->uDimension = np;
2334  patch->uKnot.p = malloc(profileyz->knot.n * sizeof(double));
2335  memcpy(patch->uKnot.p,profileyz->knot.p,profileyz->knot.n * sizeof(double));
2336  patch->uKnot.n = profileyz->knot.n;
2337  patch->uOrder = profileyz->order;
2338  patch->uTessellation = (int)((float)profileyz->tessellation * profileyz->_tscale);
2339  //v will be trajectory
2340  patch->vDimension = nt;
2341  patch->vKnot.p = malloc(trajectoryxz->knot.n * sizeof(double));
2342  memcpy(patch->vKnot.p,trajectoryxz->knot.p,trajectoryxz->knot.n * sizeof(double));
2343  patch->vKnot.n = trajectoryxz->knot.n;
2344  patch->vOrder = trajectoryxz->order;
2345  patch->vTessellation = (int)((float)trajectoryxz->tessellation * profileyz->_tscale);
2346  if(0){
2347  int ic = 0;
2348  for(j=0;j<nt;j++){
2349  for(k=0;k<np;k++){
2350  printf("%f %f %f,",xyz[ic*3 + 0], xyz[ic*3 +1], xyz[ic*3 +2]);
2351  ic++;
2352  }
2353  printf("\n");
2354  }
2355  printf("uDimension=%d vDimension=%d nc=%d\n",np,nt,ic);
2356  }
2357  compile_NurbsPatchSurface((struct X3D_NurbsPatchSurface*)node->_patch);
2358 }
2359 void rendray_NurbsSwungSurface (struct X3D_NurbsSwungSurface *node) {
2360  COMPILE_IF_REQUIRED
2361  if (!node->_intern) return;
2362  render_ray_polyrep(node->_patch);
2363 }
2364 
2365 void collide_NurbsSwungSurface (struct X3D_NurbsSwungSurface *node) {
2366  COMPILE_IF_REQUIRED
2367  if (!node->_intern) return;
2368  collide_genericfaceset(node->_patch);
2369 }
2370 
2371 void render_NurbsSwungSurface (struct X3D_NurbsSwungSurface *node) {
2372  struct X3D_NurbsPatchSurface *patch;
2373  COMPILE_IF_REQUIRED
2374  if (!node->_patch->_intern)
2375  return;
2376  patch = (struct X3D_NurbsPatchSurface*) node->_patch;
2377  CULL_FACE(patch->solid)
2378  render_polyrep(X3D_NODE(patch));
2379 }
2380 
2381 
2382 //SWEPT
2383 int compute_tessellation(int tessellation, int order, int ncontrol ){
2384  /* for a given axis on a nurbs surface or nurbs curve,
2385  and given the node->tesselation value, order and nDimension or ncontrol for the axis
2386  return a computed tesselation value
2387  */
2388  int mtessu, ntessu; //mtessv, ntessv,
2389  mtessu = order + 1;
2390  ntessu = tessellation;
2391 
2392  if(ntessu > 0)
2393  mtessu = max(mtessu,ntessu+1);
2394  else if(ntessu < 0)
2395  mtessu = max(mtessu,(-ntessu * ncontrol) + 1);
2396  else
2397  mtessu = max(mtessu,2*ncontrol + 1);
2398  return mtessu;
2399 }
2400 void compute_knotvector(int order, int ncontrol, int nknots, double *knots, int *newnknots, float **newknots, float *range){
2401  //for a given axis, QCs the knot vector, and replaces if necessary
2402  //will malloc the necessary newknot vector - caller responsible
2403  //range: pass in float[2] already allocated
2404  int nku,i,ii;
2405  float *knotsu;
2406  if(knotsOK(order,ncontrol,nknots,knots)){
2407  //could do another check: max number of consecutive equal value knots == order
2408  //could do another check: knot values == or ascending
2409  nku = nknots;
2410  knotsu = MALLOC(void *, nku * sizeof(GLfloat));
2411  for(i=0;i<nku;i++){
2412  knotsu[i] = (GLfloat)knots[i];
2413  }
2414  if(DEBG){
2415  printf("good u knot vector nk=%d\n",nku);
2416  for(ii=0;ii<nku;ii++)
2417  printf("[%d]=%f \n",ii,knotsu[ii]);
2418  }
2419  *newknots = knotsu;
2420  *newnknots = nku;
2421  }else{
2422  //generate uniform knot vector
2423  static int once = 0;
2424  nku = ncontrol + order ;
2425  //caller: please malloc knots = MALLOC(void *, (ncontrol + order ) * sizeof(float))
2426  knotsu = MALLOC(void *, nku *sizeof(GLfloat));
2427  generateUniformKnotVector(order,ncontrol, knotsu);
2428  if(!once){
2429  printf("bad u knot vector given, replacing with:\n");
2430  for(ii=0;ii<nku;ii++)
2431  printf("[%d]=%f \n",ii,knotsu[ii]);
2432  once = 1;
2433  }
2434  *newknots = knotsu;
2435  *newnknots= nku;
2436  //nk = 0;
2437  }
2438  range[0] = knotsu[0];
2439  range[1] = knotsu[nku-1];
2440 }
2441 void compute_weightedcontrol(double *xyz, int dim, int nc, int nweight, double *weights, float **cxyzw){
2442 /* will malloc xyzw the right size, and return, caller responsible
2443  the opengl and piegl functions take floats
2444  dim = 3 for incoming xyz control
2445  dim = 2 for incoming xy (2D) control
2446  outgoing is [w*x,w*y,w*z,w]
2447 */
2448  int i,j;
2449  float *xyzw;
2450 
2451  xyzw = MALLOC(float *, nc * 4 * sizeof(GLfloat));
2452  for(i=0;i<nc;i++)
2453  xyzw[i*4 +0] = xyzw[i*4 +1] =xyzw[i*4 +2] =xyzw[i*4 +3] = 0.0f;
2454  for(i=0;i<nc;i++){
2455  for(j=0;j<dim;j++){
2456  xyzw[i*4 + j] = (float)xyz[i*dim + j];
2457  }
2458  xyzw[i*4 + 3] = 1.0f;
2459  }
2460  if(nweight && nweight == nc ){
2461  for(i=0;i<nc;i++){
2462  float wt = (float)weights[i];
2463  xyzw[i*4 + 3] = wt;
2464  vecscale3f(&xyzw[i*4],&xyzw[i*4],wt);
2465  }
2466  }
2467  *cxyzw = xyzw;
2468 }
2469 void compute_doublecontrol(struct X3D_Node *controlPoint, int *nc, double** xyz ){
2470  /* rather than switch-casing elsewhere in the code, we'll get both types into double here
2471 
2472  */
2473  int n,i;
2474  double *xyzd = NULL;
2475  n = 0;
2476  switch(controlPoint->_nodeType){
2477  case NODE_Coordinate:
2478  {
2479  struct X3D_Coordinate *tcoord = (struct X3D_Coordinate*) controlPoint;
2480  n = tcoord->point.n;
2481  xyzd = MALLOC(double *,n * 3 * sizeof(double));
2482  for(i=0;i<tcoord->point.n;i++)
2483  float2double(&xyzd[i*3],tcoord->point.p[i].c,3);
2484  }
2485  break;
2486  case NODE_CoordinateDouble:
2487  {
2488  struct X3D_CoordinateDouble *tcoord = (struct X3D_CoordinateDouble *)controlPoint;
2489  n = tcoord->point.n;
2490  xyzd = MALLOC(double *,n * 3 * sizeof(double));
2491  for(i=0;i<tcoord->point.n;i++)
2492  veccopyd(&xyzd[i*3],tcoord->point.p[i].c);
2493 
2494  }
2495  break;
2496  default:
2497  break;
2498  }
2499  *nc = n;
2500  *xyz = xyzd;
2501 
2502 }
2503 float *vecabs3f(float *res, float *p){
2504  int i;
2505  for(i=0;i<3;i++)
2506  res[i] = fabsf(p[i]);
2507  return res;
2508 }
2509 int ivecdominantdirection3f(int *irank, float *p){
2510  float rmax, rmin, vabs[3];
2511  int i,iret = 0;
2512  vecabs3f(vabs,p);
2513  rmax = max(vabs[0],max(vabs[1],vabs[2]));
2514  rmin = min(vabs[0],min(vabs[1],vabs[2]));
2515  for(i=0;i<3;i++){
2516  irank[i] = 1;
2517  if(vabs[i] == rmax) {
2518  irank[i] = 2;
2519  iret = i;
2520  }
2521  if(vabs[i] == rmin) irank[i] = 0;
2522  }
2523  return iret;
2524 }
2525 void convert_mesh_to_polyrep(float *xyz, int npts, float *nxyz, int* tindex, int ntri, struct X3D_Node *node){
2526  //this is a bit like compile_polyrep, except the virt_make is below
2527 
2528  struct X3D_PolyRep *rep_, *polyrep;
2529  GLuint *norindex; //*cindex,
2530 
2531  //from compile_polyrep:
2532  //node = X3D_NODE(innode);
2533  //virt = virtTable[node->_nodeType];
2534 
2535  /* first time through; make the intern structure for this polyrep node */
2536  if(node->_intern){
2537  polyrep = node->_intern;
2538  FREE_IF_NZ(polyrep->cindex);
2539  FREE_IF_NZ(polyrep->actualCoord);
2540  FREE_IF_NZ(polyrep->GeneratedTexCoords[0]);
2541  FREE_IF_NZ(polyrep->colindex);
2542  FREE_IF_NZ(polyrep->color);
2543  FREE_IF_NZ(polyrep->norindex);
2544  FREE_IF_NZ(polyrep->normal);
2545  FREE_IF_NZ(polyrep->flat_normal);
2546  FREE_IF_NZ(polyrep->tcindex);
2547  }
2548  if(!node->_intern)
2549  node->_intern = create_polyrep();
2550 
2551  rep_ = polyrep = node->_intern;
2552 
2553 
2554  /* if multithreading, tell the rendering loop that we are regenning this one */
2555  /* if singlethreading, this'll be set to TRUE before it is tested */
2556 
2557  //<< END FROM Compile_polyrep
2558 
2559  // Start Virt_make_polyrep section >>>
2560 
2561  if (npts > 0)
2562  {
2563  //printf("npoints %d ntc %d\n",npoints,ntc);
2564  rep_->actualCoord = xyz;
2565  rep_->normal = nxyz;
2566  }
2567  rep_->ntri = ntri;
2568 
2569  if (rep_->ntri > 0)
2570  {
2571  rep_->cindex = tindex;
2572  norindex = rep_->norindex = MALLOC(GLuint *,sizeof(GLuint)*3*ntri);
2573  memcpy(norindex,tindex,sizeof(GLuint)*3*ntri);
2574  }
2575 
2576 
2577  //END virt_make_polyrep section <<<<<<
2578 
2579  //FROM Compile_polyrep
2580  if (polyrep->ntri != 0) {
2581  //float *fogCoord = NULL;
2582  stream_polyrep(node, NULL,NULL,NULL,NULL, NULL);
2583  /* and, tell the rendering process that this shape is now compiled */
2584  }
2585  //else wait for set_coordIndex to be converted to coordIndex
2586  polyrep->irep_change = node->_change;
2587 
2588  /*
2589  // dump then can copy and paste to x3d or wrl IndexedFaceSet.coordIndex and Coordinate.point fields
2590  FILE * fp = fopen("IFS_DUMP.txt","w+");
2591  fprintf(fp,"#vertices %d\n",np);
2592  for(i=0;i<np;i++){
2593  fprintf(fp,"%f %f %f\n",rep->actualCoord[i*3 +0],rep->actualCoord[i*3 +1],rep->actualCoord[i*3 +2]);
2594  }
2595  fprintf(fp,"#face indices %d\n",ni);
2596  for(i=0;i<ni;i++){
2597  fprintf(fp,"%d ",rep->cindex[i]);
2598  if((ni+1) % 3 == 0)
2599  fprintf(fp,"%d ",-1);
2600  }
2601  fprintf(fp,"\n");
2602  fclose(fp);
2603  */
2604 
2605 }
2606 void compile_NurbsSweptSurface(struct X3D_NurbsSweptSurface *node){
2607  int nt, np;
2608  double *xyzx;
2609  double *xyzt;
2610  struct X3D_NurbsCurve *trajectory;
2611  struct X3D_NurbsCurve2D *xsection;
2612  MARK_NODE_COMPILED
2613  //Swept:
2614 
2615  trajectory = (struct X3D_NurbsCurve *)node->trajectoryCurve;
2616  xyzt = NULL;
2617  nt = 0;
2618  /*
2619  switch(trajectory->controlPoint->_nodeType){
2620  case NODE_Coordinate:
2621  {
2622  struct X3D_Coordinate *tcoord = (struct X3D_Coordinate*) trajectory->controlPoint;
2623  nt = tcoord->point.n;
2624  xyzt = MALLOC(double *,nt * sizeof(double));
2625  for(int i=0;i<tcoord->point.n;i++)
2626  float2double(&xyzt[i*3],tcoord->point.p[i].c,3);
2627  }
2628  break;
2629  case NODE_CoordinateDouble:
2630  {
2631  struct X3D_CoordinateDouble *tcoord = (struct X3D_CoordinateDouble *)trajectory->controlPoint;
2632  nt = tcoord->point.n;
2633  xyzt = MALLOC(double *,nt * sizeof(double));
2634  for(int i=0;i<tcoord->point.n;i++)
2635  veccopyd(&xyzt[i*3],tcoord->point.p[i].c);
2636 
2637  }
2638  break;
2639  default:
2640  break;
2641  }
2642  */
2643  compute_doublecontrol(trajectory->controlPoint,&nt,&xyzt);
2644  xsection = (struct X3D_NurbsCurve2D *)node->crossSectionCurve;
2645 
2646  np = xsection->controlPoint.n;
2647  xyzx = (double*)xsection->controlPoint.p;
2648 
2649  // "The Nurbs Book", Les Piegl, Wayne Tiller, 2nd, 1997, Springer
2650  // piegl p.472 10.4 Swept Surfaces
2651  //ALGO Method 1.
2652  // method 1. S(u,v) = T(v) + C(u)
2653  // - has a precise NURBS definition
2654  // - (no need for spine, B up vector, planes, skinning)
2655  // - just set up the Suv control net and weights, and delegate to Patch
2656  // P(i,j) = Tj + Qi (control points)
2657  // w(i,j) = wT(j) x wC(i) (weights)
2658  // - but piegl Figure 10.11 p.474 shows the results of this sweep:
2659  // * its good for nearly linear trajectories
2660  // x but not good if you turn a 90 corner
2661  // - the profile isn't rotated with the trajectory
2662  // - so the 'tube' will flatten
2663  //ALGO 2 Method 2.
2664  // if you want the tube to stay open ie profile rotates with trajectory curve,
2665  // then you need to implement method 2. and for that
2666  // 1. compute tesselation points along trajectory curve (use piegl CurvePoint)
2667  // - get direction vector aka Tangent T of curve using Delta or Derivs, as xsection normal
2668  // - project up vector from last profile (1st profile up is arbitrary)
2669  // Piegl p.483 formula 10.27:
2670  // B0 - arbitrary unit vector perpendicular/orthogonal to trajectory tangent vector at v0
2671  // Ti = T'(vi)/|T'(vi)| //get the tangent vector along the trajectory curve, can use delta or Derivs
2672  // bi = Bi-1 - (Bi-1 * Ti)Ti
2673  // Bi = bi/|bi|
2674  // 2. comupte tesselated cross section aka xsection aka profile (use piegl CurvePoint)
2675  // 3. for each trajectory tesselation point:
2676  // a) insert up- and tangent- oriented xsection points
2677  // b) skin: join current xsection points with last with triangles
2678  //
2679  if(!strcmp(node->method->strptr,"FULL"))
2680  node->_method = 2;
2681  if(!strcmp(node->method->strptr,"TRANSLATE"))
2682  node->_method = 1;
2683  if(node->_method == 1){
2684  //ALGO 1 Suv = T(v) + C(u)
2685  struct X3D_NurbsPatchSurface *patch;
2686  struct X3D_Coordinate *controlPoint;
2687  float *xyz;
2688  int ic,j,i;
2689  double *weight;
2690 
2691  patch = (struct X3D_NurbsPatchSurface*)node->_patch;
2692  controlPoint = (struct X3D_Coordinate*)patch->controlPoint;
2693  if(!patch){
2694  patch = (struct X3D_NurbsPatchSurface*) createNewX3DNode(NODE_NurbsPatchSurface);
2695  controlPoint = (struct X3D_Coordinate*) createNewX3DNode(NODE_Coordinate);
2696  node->_patch = X3D_NODE(patch);
2697  patch->controlPoint = X3D_NODE(controlPoint);
2698  }
2699  xyz = MALLOC(float*,nt * np * 3 * sizeof(float));
2700  controlPoint->point.p = (struct SFVec3f*) xyz;
2701  controlPoint->point.n = nt * np;
2702 
2703  ic = 0;
2704  for(j=0;j<nt;j++){
2705  float pt[3];
2706  double2float(pt,&xyzt[j*3],3);
2707  for(i=0;i<np;i++){
2708  float pp[3];
2709  double2float(pp,&xyzx[2*i],2);
2710  xyz[ic*3 + 0] = pt[0] + pp[0];
2711  xyz[ic*3 + 1] = pt[1] + pp[1];
2712  xyz[ic*3 + 2] = pt[2];
2713  ic++;
2714  }
2715  }
2716 
2717  weight = NULL;
2718  if((trajectory->weight.n && trajectory->weight.n == nt) ||
2719  (xsection->weight.n && xsection->weight.n == np)){
2720  //we have some non-default weights, so apply the piegl formula p.473
2721  // w(i,j) = wT(j) x wC(i) (weights)
2722  weight = MALLOC(double*,nt * np * sizeof(double));
2723  for(j=0;j<nt;j++){
2724  double wtTj = trajectory->weight.p[j];
2725  for(i=0;i<np;i++)
2726  weight[j*np + i] = wtTj * xsection->weight.p[i];
2727  }
2728  }
2729 
2730 
2731  if(weight){
2732  patch->weight.p = weight;
2733  patch->weight.n = nt * np;
2734  }
2735  patch->solid = node->solid;
2736  //u will be profile,
2737  patch->uDimension = np;
2738  patch->uKnot.p = malloc(xsection->knot.n * sizeof(double));
2739  memcpy(patch->uKnot.p,xsection->knot.p,xsection->knot.n * sizeof(double));
2740  patch->uKnot.n = xsection->knot.n;
2741  patch->uOrder = xsection->order;
2742  patch->uTessellation = (int)((float)xsection->tessellation * xsection->_tscale);
2743  //v will be trajectory
2744  patch->vDimension = nt;
2745  patch->vKnot.p = malloc(xsection->knot.n * sizeof(double));
2746  memcpy(patch->vKnot.p,trajectory->knot.p,trajectory->knot.n * sizeof(double));
2747  patch->vKnot.n = trajectory->knot.n;
2748  patch->vOrder = trajectory->order;
2749  patch->vTessellation = (int)((float)trajectory->tessellation * trajectory->_tscale);
2750  if(0){
2751  int j,k,ic = 0;
2752  for(j=0;j<nt;j++){
2753  for(k=0;k<np;k++){
2754  printf("%f %f %f,",xyz[ic*3 + 0], xyz[ic*3 +1], xyz[ic*3 +2]);
2755  ic++;
2756  }
2757  printf("\n");
2758  }
2759  printf("uDimension=%d vDimension=%d nc=%d\n",np,nt,ic);
2760  }
2761  compile_NurbsPatchSurface((struct X3D_NurbsPatchSurface*)node->_patch);
2762  } //end method == 1
2763  if(node->_method == 2){
2764  //ALGO 2 skinning like extrusion
2765  int mtessv, mtessu, nku, nkv;
2766  int i,DBGSW;
2767  float *knotsu,*knotsv,*xyzwu,*xyzwv, urange[2],vrange[2];
2768  float *Tv, *Tangentv, *Bup, *pts;
2769  float *Qu, *Nu;
2770  float *normals;
2771  int *idx;
2772  int ic, it;
2773  float matB0[9];
2774  int ntri;
2775 
2776 
2777  int mtessu1, mtessv1;
2778 
2779  mtessu = compute_tessellation(xsection->tessellation,xsection->order,np);
2780  mtessu = (int)((float)mtessu * xsection->_tscale);
2781  mtessv = compute_tessellation(trajectory->tessellation,trajectory->order,nt);
2782  mtessv = (int)((float)mtessv * trajectory->_tscale);
2783  compute_knotvector(xsection->order,np,xsection->knot.n,xsection->knot.p,&nku,&knotsu,urange);
2784  compute_knotvector(trajectory->order,nt,trajectory->knot.n,trajectory->knot.p,&nkv,&knotsv,vrange);
2785  compute_weightedcontrol(xyzt,3,nt, trajectory->weight.n, trajectory->weight.p, &xyzwv);
2786  compute_weightedcontrol(xyzx,2,np, xsection->weight.n, xsection->weight.p, &xyzwu);
2787  DBGSW = FALSE;
2788  if(DBGSW){
2789  int i;
2790  printf("np %d mtessu %d nku %d, nt %d mtessv %d, nkv %d",np,mtessu,nku,nt,mtessv,nkv);
2791  printf("trajectory nt %d points:\n",nt);
2792  for(i=0;i<nt;i++)
2793  printf("%d %f %f %f %f\n",i,xyzwv[i*4 + 0],xyzwv[i*4 + 1],xyzwv[i*4 + 2],xyzwv[i*4 + 3]);
2794  printf("xsection np %d points:\n",np);
2795  for(i=0;i<np;i++)
2796  printf("%d %f %f %f %f\n",i,xyzwu[i*4 + 0],xyzwu[i*4 + 1],xyzwu[i*4 + 2],xyzwu[i*4 + 3]);
2797  }
2798  // 1. compute tesselation points along trajectory curve (use piegl CurvePoint)
2799  // - get direction vector aka Tangent T of curve using Delta or Derivs, as xsection normal
2800  // - project up-vector from last profile (1st profile up is arbitrary)
2801  // Piegl p.483 formula 10.27:
2802  // B0 - arbitrary unit vector perpendicular/orthogonal to trajectory tangent vector at v0
2803  // Ti = T'(vi)/|T'(vi)| //get the tangent vector Tangentv along the trajectory curve, can use delta or Derivs
2804  // bi = Bi-1 - (Bi-1 dot Ti)Ti //where dot is the dot product
2805  // Bi = bi/|bi|
2806  mtessu1 = mtessu + 1;
2807  mtessv1 = mtessv + 1;
2808  Tv = MALLOC(float*,(mtessv1)*3*sizeof(float));
2809  Tangentv = MALLOC(float*,(mtessv1)*3*sizeof(float));
2810  Bup = MALLOC(float*,(mtessv1)*3*sizeof(float));
2811  for(i=0;i<mtessv1;i++){
2812  float cw[4], cw1[4], delta[3], v;
2813  v = (float)i*(vrange[1]-vrange[0])/(float)mtessv;
2814  CurvePoint(nt, trajectory->order-1, knotsv, xyzwv, v, cw );
2815  veccopy3f(&Tv[i*3],cw);
2816  //- get direction vector aka Tangent T of curve using Delta or Derivs, as xsection normal
2817  CurvePoint(nt, trajectory->order-1, knotsv, xyzwv, v+.01f, cw1 );
2818  vecdif3f(delta,cw1,cw);
2819  // Ti = T'(vi)/|T'(vi)| //get the tangent vector Tangentv along the trajectory curve, can use delta or Derivs
2820  vecnormalize3f(delta,delta);
2821  veccopy3f(&Tangentv[i*3],delta);
2822  //- project up-vector from last profile (1st profile up is arbitrary)
2823  if(i==0){
2824  // B0 - arbitrary unit vector perpendicular/orthogonal to trajectory tangent vector at v0
2825  int k,irank[3],idom, inondom;
2826  float perp[3], perp2[3];
2827  idom = ivecdominantdirection3f(irank,delta);
2828  inondom = idom + 1 > 2 ? 0 : idom + 1;
2829  for(k=0;k<3;k++) if(irank[k] == 0) inondom = k;
2830  memset(perp,0,3*sizeof(float));
2831  perp[inondom] = 1.0; //close to perpendicular, but not quite
2832  veccross3f(perp2,delta,perp); //perp2 perpendicular to perp and delta
2833  veccross3f(perp,perp2,delta); //another cross to get perp exactly perpendicular
2834  //perp should be perpendicular to delta[0]
2835  veccopy3f(&Bup[i*3],perp);
2836  }else{
2837  // bi = Bi-1 - (Bi-1 dot Ti)Ti
2838  // Bi = bi/|bi|
2839  float bi[3], bi1dotti, tiscaled[3];
2840  bi1dotti = vecdot3f(&Bup[(i-1)*3],&Tangentv[i*3]);
2841  vecdif3f(bi,&Bup[(i-1)*3],vecscale3f(tiscaled,&Tangentv[i*3],bi1dotti));
2842  vecnormalize3f(&Bup[i*3],bi);
2843  }
2844  }
2845  if(DBGSW){
2846  printf("trajectory T:\n");
2847  for(i=0;i<mtessv1;i++){
2848  printf("%d [%f %f %f] \n",i,
2849  Tv[i*3 +0],Tv[i*3 +1],Tv[i*3 +2]);
2850  }
2851  printf("trajectory T', B:\n");
2852  for(i=0;i<mtessv1;i++){
2853  printf("%d [%f %f %f] [%f %f %f] \n",i,
2854  Tangentv[i*3 +0],Tangentv[i*3 +1],Tangentv[i*3 +2],
2855  Bup[i*3 +0],Bup[i*3 +1],Bup[i*3 +2]);
2856  }
2857 
2858  }
2859  //if(trajectory->closed){
2860  // //compute backward Bup, and average fore and aft Bup
2861  //}
2862  // 2. comupte tesselated cross section aka xsection aka profile (use piegl CurvePoint)
2863  Qu = MALLOC(float*,(mtessu+1)*3*sizeof(float)); //cross section points
2864  Nu = MALLOC(float*,(mtessu+1)*3*sizeof(float)); //cross section normals
2865  if(DBGSW) printf("Xsection tess pts:\n");
2866  for(i=0;i<mtessu1;i++){
2867  float u, cw[4], cw1[4], delta[3], normal[3];
2868  float zzz[3] = {0.0f,0.0f,1.0f};
2869  u = (float)i*(urange[1]-urange[0])/(float)mtessu;
2870  CurvePoint(np, xsection->order-1, knotsu, xyzwu, u, cw );
2871  veccopy3f(&Qu[i*3],cw);
2872  CurvePoint(np, xsection->order-1, knotsu, xyzwu, u+.01f, cw1 );
2873  vecdif3f(delta,cw1,cw);
2874  vecnormalize3f(delta,delta);
2875  veccross3f(normal,zzz,delta);
2876  veccopy3f(&Nu[i*3],normal);
2877  if(DBGSW) printf("%d %f %f %f\n",i,Qu[i*3 +0],Qu[i*3 +1],Qu[i*3 +2]);
2878 
2879  }
2880  // 3. for each trajectory tesselation point:
2881  // a) insert up- and tangent- oriented xsection points
2882  // b) skin: join current xsection points with last with triangles
2883  pts = MALLOC(float*,mtessu1 * mtessv1 * 3 * sizeof(float));
2884  normals = MALLOC(float*,mtessu1 * mtessv1 * 3 * sizeof(float));
2885  idx = MALLOC(int *, mtessu * mtessv * 2 * 3 * sizeof(int));
2886  ic = 0;
2887  it = 0;
2888  for(i=0;i<mtessv1;i++){
2889  //insert oriented xsection at T(v)
2890  float mat [9], matt[9];
2891  int j;
2892  //set up 3x3 rotation by using 3 perpendicular local unit vectors as rot mat rows
2893  //http://renderdan.blogspot.ca/2006/05/rotation-matrix-from-axis-vectors.html
2894 
2895  veccross3f(&mat[0],&Bup[i*3],&Tangentv[i*3]);
2896  veccopy3f(&mat[3],&Bup[i*3]);
2897  veccopy3f(&mat[6],&Tangentv[i*3]);
2898  mattranspose3f(matt,mat); //seems like I need columns, not rows, by experimentation
2899  if(i==0){
2900  //not sure but I think we should take off the
2901  //arbitrary rotation of the first crossection from all subsequent
2902  //so first xsection not rotated (you need to design with
2903  //your crosssection plane perpendicular to the start of your trajectory)
2904  //and subsequent are rotated with respect to first
2905  //Looks good
2906  memcpy(matB0,mat,9*sizeof(float));
2907  }
2908  matmultiply3f(mat,matt,matB0);
2909  for(j=0;j<mtessu1;j++){
2910  float pp[3], norm[3];
2911  matmultvec3f(pp, mat, &Qu[j*3] ); //orient profile point
2912  vecadd3f(pp,pp,&Tv[i*3]); //add on trajectory point
2913  veccopy3f(&pts[ic*3],pp);
2914  matmultvec3f(norm,mat,&Nu[j*3]);
2915  veccopy3f(&normals[ic*3],norm);
2916  ic++;
2917  }
2918  //connect to last xsection with triangles
2919  if(i > 0){
2920  int j,kk,mm;
2921  kk = (i-1)*mtessu1;
2922  mm = i*mtessu1;
2923  for(j=0;j<mtessu;j++){
2924  // 1 1 3
2925  // 0 2 2
2926  //first triangle
2927  idx[it++] = kk+j;
2928  idx[it++] = mm+j;
2929  idx[it++] = kk+j+1;
2930  //second triangle
2931  idx[it++] = kk+j+1;
2932  idx[it++] = mm+j;
2933  idx[it++] = mm+j+1;
2934  }
2935  }
2936  }
2937  //assign to something
2938  ntri = it/3;
2939  if(DBGSW){
2940  printf("ntri %d triangle indexes:\n",ntri);
2941  for(i=0;i<ntri;i++){
2942  int j;
2943  printf("%d [",i);
2944  for(j=0;j<3;j++){
2945  printf("%d ",idx[i*3 +j]);
2946  }
2947  printf("\n");
2948  }
2949  printf("triangle vertices:\n");
2950  for(i=0;i<ntri;i++){
2951  int j;
2952  for(j=0;j<3;j++){
2953  float pt[3];
2954  int ix = idx[i*3 +j];
2955  veccopy3f(pt,&pts[ix*3]);
2956  printf("%d %d %f %f %f\n",i,ix,pt[0],pt[1],pt[2]);
2957  }
2958  }
2959  }
2960 
2961  //send to polyrep
2962  convert_mesh_to_polyrep(pts,ic,normals,idx,ntri,X3D_NODE(node));
2963  }
2964 
2965 }
2966 void rendray_NurbsSweptSurface (struct X3D_NurbsSweptSurface *node) {
2967  COMPILE_IF_REQUIRED
2968  if(node->_method == 1){
2969  if (!node->_patch) return;
2970  render_ray_polyrep(node->_patch);
2971  }
2972  if(node->_method == 2){
2973  if(!node->_intern) return;
2974  render_ray_polyrep(node);
2975  }
2976 }
2977 
2978 void collide_NurbsSweptSurface (struct X3D_NurbsSweptSurface *node) {
2979  COMPILE_IF_REQUIRED
2980  if(node->_method == 1){
2981  if (!node->_patch) return;
2982  collide_genericfaceset(node->_patch);
2983  }
2984  if(node->_method == 2){
2985  if (!node->_intern) return;
2986  collide_genericfaceset(node);
2987  }
2988 }
2989 
2990 void render_NurbsSweptSurface (struct X3D_NurbsSweptSurface *node) {
2991  COMPILE_IF_REQUIRED
2992  if(node->_method == 1){
2993  struct X3D_NurbsPatchSurface *patch;
2994  if (!node->_patch->_intern)
2995  return;
2996  patch = (struct X3D_NurbsPatchSurface *)node->_patch;
2997  CULL_FACE(patch->solid)
2998  render_polyrep(patch);
2999  }
3000  if(node->_method == 2){
3001  if (!node->_intern)
3002  return;
3003  //CULL_FACE(node->solid)
3004  render_polyrep(node);
3005 
3006  }
3007 }
3008 void compile_NurbsSet(struct X3D_NurbsSet *node){
3009  int i;
3010  MARK_NODE_COMPILED
3011  //tessellationScale
3012  for(i=0;i<node->geometry.n;i++){
3013  struct X3D_Node *gn = (struct X3D_Node*)node->geometry.p[i];
3014  switch(gn->_nodeType){
3015  case NODE_NurbsCurve:
3016  {
3017  struct X3D_NurbsCurve* g = (struct X3D_NurbsCurve*)gn;
3018  g->_tscale = node->tessellationScale;
3019  MNX(g);
3020  }
3021  break;
3022  case NODE_NurbsCurve2D:
3023  {
3024  struct X3D_NurbsCurve2D* g = (struct X3D_NurbsCurve2D*)gn;
3025  g->_tscale = node->tessellationScale;
3026  MNX(g);
3027  }
3028  break;
3029  case NODE_NurbsPatchSurface:
3030  {
3031  struct X3D_NurbsPatchSurface* g = (struct X3D_NurbsPatchSurface*)gn;
3032  g->_tscale = node->tessellationScale;
3033  MNX(g);
3034  }
3035  break;
3036  case NODE_NurbsTrimmedSurface:
3037  {
3038  struct X3D_NurbsTrimmedSurface* g = (struct X3D_NurbsTrimmedSurface*)gn;
3039  g->_tscale = node->tessellationScale;
3040  MNX(g);
3041  }
3042  break;
3043  case NODE_NurbsSweptSurface:
3044  {
3045  struct X3D_NurbsSweptSurface* g = (struct X3D_NurbsSweptSurface*)gn;
3046  if(g->_method == 1){
3047  if(g->_patch){
3048  struct X3D_NurbsPatchSurface *patch = (struct X3D_NurbsPatchSurface *)g->_patch;
3049  patch->_tscale = node->tessellationScale;
3050  MNX(g);
3051  }
3052  }
3053  if(g->_method == 2){
3054  struct X3D_NurbsCurve2D *curve = (struct X3D_NurbsCurve2D *)g->crossSectionCurve;
3055  curve->_tscale = node->tessellationScale;
3056  MNX(g->crossSectionCurve);
3057  curve = (struct X3D_NurbsCurve2D *)g->trajectoryCurve;
3058  curve->_tscale = node->tessellationScale;
3059  MNX(g->trajectoryCurve);
3060  }
3061  }
3062  break;
3063  case NODE_NurbsSwungSurface:
3064  {
3065  struct X3D_NurbsSwungSurface* g = (struct X3D_NurbsSwungSurface*)gn;
3066  if(g->_patch){
3067  struct X3D_NurbsPatchSurface *patch = (struct X3D_NurbsPatchSurface *)g->_patch;
3068  patch->_tscale = node->tessellationScale;
3069  MNX(g);
3070  }
3071  }
3072  break;
3073  }
3074  }
3075 }
3076 void render_NurbsSet(struct X3D_NurbsSet *node){
3077  COMPILE_IF_REQUIRED
3078 }
3079 #else //LIB_NURBS
3080 void compile_ContourPolyline2D(struct X3D_ContourPolyline2D *node){}
3081 void compile_NurbsCurve(struct X3D_NurbsCurve *node){}
3082 void render_NurbsCurve(struct X3D_NurbsCurve *node){}
3083 void compile_NurbsPatchSurface(struct X3D_NurbsPatchSurface *node){}
3084 void rendray_NurbsPatchSurface (struct X3D_NurbsPatchSurface *node) {}
3085 void collide_NurbsPatchSurface (struct X3D_NurbsPatchSurface *node) {}
3086 void render_NurbsPatchSurface (struct X3D_NurbsPatchSurface *node) {}
3087 void compile_NurbsTrimmedSurface(struct X3D_NurbsTrimmedSurface *node){}
3088 void rendray_NurbsTrimmedSurface (struct X3D_NurbsTrimmedSurface *node) {}
3089 void collide_NurbsTrimmedSurface (struct X3D_NurbsTrimmedSurface *node) {}
3090 void render_NurbsTrimmedSurface (struct X3D_NurbsTrimmedSurface *node) {}
3091 void do_NurbsPositionInterpolator (void *node) {}
3092 void do_NurbsOrientationInterpolator (void *node){}
3093 void do_NurbsSurfaceInterpolator (void *_node){}
3094 void compile_NurbsSwungSurface(struct X3D_NurbsSwungSurface *node){}
3095 void rendray_NurbsSwungSurface (struct X3D_NurbsSwungSurface *node) {}
3096 void collide_NurbsSwungSurface (struct X3D_NurbsSwungSurface *node) {}
3097 void render_NurbsSwungSurface (struct X3D_NurbsSwungSurface *node) {}
3098 void compile_NurbsSweptSurface(struct X3D_NurbsSweptSurface *node){}
3099 void rendray_NurbsSweptSurface (struct X3D_NurbsSweptSurface *node) {}
3100 void collide_NurbsSweptSurface (struct X3D_NurbsSweptSurface *node) {}
3101 void render_NurbsSweptSurface (struct X3D_NurbsSweptSurface *node){}
3102 void compile_NurbsSet(struct X3D_NurbsSet *node){}
3103 void render_NurbsSet(struct X3D_NurbsSet *node){}
3104 #endif //LIB_NURBS
Definition: Vector.h:36