Bug Summary

File:egads2cgt.c
Warning:line 163, column 12
Although the value stored to 'last' is used in the enclosing expression, the value is never actually read from 'last'

Annotated Source Code

Press '?' to see keyboard shortcuts

clang -cc1 -triple x86_64-pc-linux-gnu -analyze -disable-free -disable-llvm-verifier -discard-value-names -main-file-name egads2cgt.c -analyzer-store=region -analyzer-opt-analyze-nested-blocks -analyzer-checker=core -analyzer-checker=apiModeling -analyzer-checker=unix -analyzer-checker=deadcode -analyzer-checker=security.insecureAPI.UncheckedReturn -analyzer-checker=security.insecureAPI.getpw -analyzer-checker=security.insecureAPI.gets -analyzer-checker=security.insecureAPI.mktemp -analyzer-checker=security.insecureAPI.mkstemp -analyzer-checker=security.insecureAPI.vfork -analyzer-checker=nullability.NullPassedToNonnull -analyzer-checker=nullability.NullReturnedFromNonnull -analyzer-output plist -w -setup-static-analyzer -mrelocation-model pic -pic-level 2 -mthread-model posix -mframe-pointer=none -fmath-errno -fno-rounding-math -masm-verbose -mconstructor-aliases -munwind-tables -target-cpu x86-64 -dwarf-column-info -fno-split-dwarf-inlining -debugger-tuning=gdb -resource-dir /usr/lib/llvm-10/lib/clang/10.0.0 -D REVISION=7.8 -I /home/jenkins/workspace/ESP_Stanalizer/LINUX64/CAPS/scan-build/ESP/LINUX64/include -internal-isystem /usr/local/include -internal-isystem /usr/lib/llvm-10/lib/clang/10.0.0/include -internal-externc-isystem /usr/include/x86_64-linux-gnu -internal-externc-isystem /include -internal-externc-isystem /usr/include -O2 -Wno-dangling-else -Wno-parentheses -Wno-unused-result -Wno-format-truncation -fdebug-compilation-dir /home/jenkins/workspace/ESP_Stanalizer/LINUX64/CAPS/scan-build/CAPS/aim/cgt -ferror-limit 19 -fmessage-length 0 -fgnuc-version=4.2.1 -fobjc-runtime=gcc -fdiagnostics-show-option -vectorize-loops -vectorize-slp -analyzer-output=html -faddrsig -o /home/jenkins/workspace/ESP_Stanalizer/LINUX64/CAPS/scan-build/scanCAPS/2024-04-01-115124-71240-1 -x c egads2cgt.c
1/*
2 * EGADS: Electronic Geometry Aircraft Design System
3 *
4 * Tesselate EGADS bodies and write:
5 * CART3D unstructured surface triangulation tri file
6 * with FaceID in Component field
7 * TESS information about vertex/geometry ownership
8 * PLOT3D structured surface grid file from quadding algorithm
9 * PLOT3D structured surface grid file from direct uv parameters space
10 * evaluation (untrimmed)
11 *
12 * Copyright 2014-2024, Massachusetts Institute of Technology
13 * Licensed under The GNU Lesser General Public License, version 2.1
14 * See http://www.opensource.org/licenses/lgpl-2.1.php
15 *
16 */
17
18#include "egads.h"
19#include <math.h>
20#include <string.h>
21
22//#define INFACEOCC
23
24#ifdef INFACEOCC
25extern int EG_inFaceOCC( const egObject *face, double tol, const double *uv );
26#endif
27extern int aflr4egads( ego model, ego *tesses );
28
29#ifdef WIN32
30#define ULONGunsigned long unsigned long long
31#define snprintf _snprintf
32#else
33#define ULONGunsigned long unsigned long
34#endif
35
36//#define DEBUG
37
38#define FUZZ1.e-14 1.e-14
39#define PI3.1415926535897931159979635 3.1415926535897931159979635
40#define MAXUVS1024 1024 /* max Us or Vs for untrimmed data */
41#define MAXSEG1023 1023
42#define CROSS(a,b,c)a[0] = (b[1]*c[2]) - (b[2]*c[1]); a[1] = (b[2]*c[0]) - (b[0]*
c[2]); a[2] = (b[0]*c[1]) - (b[1]*c[0])
a[0] = (b[1]*c[2]) - (b[2]*c[1]);\
43 a[1] = (b[2]*c[0]) - (b[0]*c[2]);\
44 a[2] = (b[0]*c[1]) - (b[1]*c[0])
45
46
47typedef struct {
48 int nu;
49 int nv;
50 double us[MAXUVS1024];
51 double vs[MAXUVS1024];
52} UVmap;
53
54
55typedef struct parmSeg {
56 struct parmSeg *prev; /* for front -- seg connection or NULL */
57 double parms[2]; /* for front -- dir is second value (-1,1) */
58 double alen; /* front only */
59 double size;
60 struct parmSeg *next;
61} parmSeg;
62
63
64/* Global vars */
65
66static char version[5] = "1.17";
67static int aflr4 = 0;
68static int wrtqud = 0;
69static int wrtuv = 0;
70static double mxang = 0.0;
71static double mxedg = 0.0;
72static double mxchd = 0.0;
73static double ggf = 1.2;
74
75
76/* front handling routines for geometric growth smoothing */
77/*@-immediatetrans@*/
78static void addSeg(int *nSeg, parmSeg *segs, parmSeg seg)
79{
80 parmSeg *link;
81
82 if (*nSeg >= MAXSEG1023) {
83 printf(" ERROR: No more room for Segments!\n");
84 exit(EXIT_FAILURE1);
85 }
86 segs[*nSeg] = seg;
87 link = seg.prev;
88 if (link != NULL((void*)0)) {
89 if (link->next != NULL((void*)0)) {
90 printf(" ERROR: Double hit for prev/next!\n");
91 exit(EXIT_FAILURE1);
92 }
93 link->next = &segs[*nSeg];
94 }
95 link = seg.next;
96 if (link != NULL((void*)0)) {
97 if (link->prev != NULL((void*)0)) {
98 printf(" ERROR: Double hit for next/prev!\n");
99 exit(EXIT_FAILURE1);
100 }
101 link->prev = &segs[*nSeg];
102 }
103 *nSeg += 1;
104}
105/*@+immediatetrans@*/
106
107static void addFront(parmSeg **frontsx, parmSeg **poolx, double params[2],
108 double alen, double size, /*@null@*/ parmSeg *connect)
109{
110 parmSeg *seg, *fronts, *pool, *link;
111
112 fronts = *frontsx;
113 pool = *poolx;
114
115 /* get a free segment */
116 if (pool != NULL((void*)0)) {
117 seg = pool;
118 pool = seg->next;
119 } else {
120 seg = (parmSeg *) EG_alloc(sizeof(parmSeg));
121 if (seg == NULL((void*)0)) {
122 printf(" ERROR: allocating Front Segment!\n");
123 exit(EXIT_FAILURE1);
124 }
125 }
126 seg->prev = connect;
127 seg->parms[0] = params[0];
128 seg->parms[1] = params[1];
129 seg->alen = alen;
130 seg->size = size;
131 seg->next = NULL((void*)0);
132 if (fronts == NULL((void*)0)) {
133 fronts = seg;
134 } else {
135 link = fronts;
136 while (link->next != NULL((void*)0)) link = link->next;
137 link->next = seg;
138 }
139
140 *frontsx = fronts;
141 *poolx = pool;
142}
143
144
145static void delFront(parmSeg **frontsx, parmSeg **poolx, parmSeg *front)
146{
147 int hit = 0;
148 parmSeg *fronts, *pool, *link, *last;
149
150 fronts = *frontsx;
151 pool = *poolx;
152
153 if (front == NULL((void*)0)) {
154 printf(" ERROR: NULL Front to delete!\n");
155 exit(EXIT_FAILURE1);
156 }
157
158 /* delete a segment */
159 if (front == fronts) {
160 fronts = front->next;
161 hit++;
162 } else {
163 link = last = fronts;
Although the value stored to 'last' is used in the enclosing expression, the value is never actually read from 'last'
164 while (link->next != NULL((void*)0)) {
165 last = link;
166 link = link->next;
167 if (link == front) {
168 last->next = front->next;
169 hit++;
170 break;
171 }
172 }
173 }
174 if (hit == 0) {
175 printf(" ERROR: Front not found to delete!\n");
176 exit(EXIT_FAILURE1);
177 }
178
179 front->next = NULL((void*)0);
180 if (pool == NULL((void*)0)) {
181 pool = front;
182 } else {
183 link = pool;
184 while (link->next != NULL((void*)0)) link = link->next;
185 link->next = front;
186 }
187
188 *frontsx = fronts;
189 *poolx = pool;
190}
191
192
193static void freeFront(/*@null@*/ parmSeg *fronts, /*@null@*/ parmSeg *pool)
194{
195 parmSeg *link, *last;
196
197 if (fronts != NULL((void*)0)) {
198 link = fronts;
199 while (link != NULL((void*)0)) {
200 last = link;
201 link = link->next;
202 EG_free(last);
203 }
204 }
205
206 if (pool != NULL((void*)0)) {
207 link = pool;
208 while (link != NULL((void*)0)) {
209 last = link;
210 link = link->next;
211 EG_free(last);
212 }
213 }
214}
215
216
217/*
218 * calculates and returns a complete Body tessellation
219 * note that this is consistent with the EGADS global numbers
220 *
221 * where: body - ego of a body tessellation
222 * nface - number of faces in the body
223 * face - the EGADS Face objects
224 * nvert - Number of vertices (returned)
225 * verts - coordinates (returned) 3*nverts in len -- freeable
226 * ntriang - number of triangles (returned)
227 * triang - triangle indices (returned) 3*ntriang in len
228 * freeable
229 * comp - Cart3D component ID per triangle -- freeable
230 */
231
232static int
233bodyTessellation(ego tess, int nface, int *nvert, double **verts,
234 int *ntriang, int **triang, int **comp)
235
236{
237 int status, i, j, k, base, npts, ntri, *tri, *table, *compID;
238 int plen, tlen;
239 const int *tris, *tric, *ptype, *pindex;
240 double *xyzs;
241 const double *points, *uv;
242
243 *nvert = *ntriang = 0;
244 *verts = NULL((void*)0);
245 *triang = NULL((void*)0);
246 *comp = NULL((void*)0);
247 npts = ntri = 0;
248
249 for (i = 1; i <= nface; i++) {
250 status = EG_getTessFace(tess, i, &plen, &points, &uv, &ptype, &pindex,
251 &tlen, &tris, &tric);
252 if (status != EGADS_SUCCESS0) {
253 printf(" Face %d: EG_getTessFace status = %d (bodyTessellation)!\n",
254 i, status);
255 } else {
256 npts += plen;
257 ntri += tlen;
258 }
259 }
260
261 /* get the memory associated with the points */
262
263 table = (int *) EG_alloc(2*npts*sizeof(int));
264 if (table == NULL((void*)0)) {
265 printf(" Error: Can not allocate node table (bodyTessellation)!\n");
266 return EGADS_MALLOC-4;
267 }
268
269 xyzs = (double *) EG_alloc(3*npts*sizeof(double));
270 if (xyzs == NULL((void*)0)) {
271 printf(" Error: Can not allocate XYZs (bodyTessellation)!\n");
272 EG_free(table);
273 return EGADS_MALLOC-4;
274 }
275
276 /* zipper up the edges -- a Face at a time */
277
278 npts = 0;
279 for (j = 1; j <= nface; j++) {
280 status = EG_getTessFace(tess, j, &plen, &points, &uv, &ptype, &pindex,
281 &tlen, &tris, &tric);
282 if (status != EGADS_SUCCESS0) continue;
283
284 for (i = 0; i < plen; i++) {
285 table[2*npts ] = ptype[i];
286 table[2*npts+1] = pindex[i];
287 xyzs[3*npts ] = points[3*i ];
288 xyzs[3*npts+1] = points[3*i+1];
289 xyzs[3*npts+2] = points[3*i+2];
290
291 /* for non-interior pts -- try to match with others */
292
293 if (ptype[i] != -1)
294 for (k = 0; k < npts; k++)
295 if ((table[2*k]==table[2*npts]) && (table[2*k+1]==table[2*npts+1])) {
296 table[2*npts ] = k;
297 table[2*npts+1] = 0;
298 break;
299 }
300
301 npts++;
302 }
303
304 }
305
306 /* fill up the whole triangle list -- a Face at a time */
307
308 tri = (int *) EG_alloc(3*ntri*sizeof(int));
309 if (tri == NULL((void*)0)) {
310 printf(" Error: Can not allocate triangles (bodyTessellation)!\n");
311 EG_free(xyzs);
312 EG_free(table);
313 return EGADS_MALLOC-4;
314 }
315 compID = (int *) EG_alloc(ntri*sizeof(int));
316 if (compID == NULL((void*)0)) {
317 printf(" Error: Can not allocate components (bodyTessellation)!\n");
318 EG_free(tri);
319 EG_free(xyzs);
320 EG_free(table);
321 return EGADS_MALLOC-4;
322 }
323
324 ntri = base = 0;
325 for (j = 1; j <= nface; j++) {
326
327 /* get the face tessellation and store it away */
328 status = EG_getTessFace(tess, j, &plen, &points, &uv, &ptype, &pindex,
329 &tlen, &tris, &tric);
330 if (status != EGADS_SUCCESS0) continue;
331
332 for (i = 0; i < tlen; i++, ntri++) {
333
334 k = tris[3*i ] + base;
335
336 if (table[2*k-1] == 0) {
337 tri[3*ntri ] = table[2*k-2] + 1;
338 } else {
339 tri[3*ntri ] = k;
340 }
341
342 k = tris[3*i+1] + base;
343
344 if (table[2*k-1] == 0) {
345 tri[3*ntri+1] = table[2*k-2] + 1;
346 } else {
347 tri[3*ntri+1] = k;
348 }
349
350 k = tris[3*i+2] + base;
351
352 if (table[2*k-1] == 0) {
353 tri[3*ntri+2] = table[2*k-2] + 1;
354 } else {
355 tri[3*ntri+2] = k;
356 }
357
358 compID[ntri] = j;
359 }
360 base += plen;
361 }
362
363
364 /* remove the unused points -- crunch the point list
365 * NOTE: the returned pointer verts has the full length (not realloc'ed)
366 */
367 for (i = 0; i < npts; i++) table[i] = 0;
368 for (i = 0; i < 3*ntri; i++) table[tri[i]-1]++;
369 for (plen = i = 0; i < npts; i++) {
370 if (table[i] == 0) continue;
371 xyzs[3*plen ] = xyzs[3*i ];
372 xyzs[3*plen+1] = xyzs[3*i+1];
373 xyzs[3*plen+2] = xyzs[3*i+2];
374 plen++;
375 table[i] = plen;
376 }
377
378 /* reset the triangle indices */
379 for (i = 0; i < 3*ntri; i++) {
380 k = tri[i]-1;
381 tri[i] = table[k];
382 }
383
384 EG_free(table);
385
386 *nvert = plen;
387 *verts = xyzs;
388 *ntriang = ntri;
389 *triang = tri;
390 *comp = compID;
391
392 return EGADS_SUCCESS0;
393
394}
395
396
397static int strroot(char *s)
398{
399
400 /* Return string index location in string s by locating last . */
401
402 int i, iroot, len;
403
404 len = strlen(s);
405 iroot = len-1;
406
407 i = len-1;
408 while ( s[i] != '.' && i > 0 ) {
409 i--;
410 }
411 if ( i >= 0 ) iroot = i;
412
413 return iroot;
414}
415
416
417static void writeAttr(FILE *fp, ego obj, /*@null@*/ char *filter)
418{
419 int i, j, n, stat, nattr, atype, alen, flen;
420 const int *ints;
421 const double *reals;
422 const char *name, *str;
423
424 stat = EG_attributeNum(obj, &nattr);
425 if (stat != EGADS_SUCCESS0) return;
426
427 for (n = i = 0; i < nattr; i++) {
428 stat = EG_attributeGet(obj, i+1, &name, &atype, &alen, &ints, &reals, &str);
429 if (stat != EGADS_SUCCESS0) continue;
430 if ((atype != ATTRINT1) && (atype != ATTRREAL2) && (atype != ATTRSTRING3))
431 continue;
432 if (filter != NULL((void*)0)) {
433 flen = strlen(filter);
434 if (flen > strlen(name)) continue;
435 for (j = 0; j < flen; j++)
436 if (filter[j] != name[j]) break;
437 if (j != flen) continue;
438 }
439 n++;
440 }
441 fprintf(fp, " %6d\n", n);
442 if (n == 0) return;
443
444 for (i = 0; i < nattr; i++) {
445 stat = EG_attributeGet(obj, i+1, &name, &atype, &alen, &ints, &reals, &str);
446 if (stat != EGADS_SUCCESS0) continue;
447 if ((atype != ATTRINT1) && (atype != ATTRREAL2) && (atype != ATTRSTRING3))
448 continue;
449 if (filter != NULL((void*)0)) {
450 flen = strlen(filter);
451 if (flen > strlen(name)) continue;
452 for (j = 0; j < flen; j++)
453 if (filter[j] != name[j]) break;
454 if (j != flen) continue;
455 }
456 if (atype == ATTRSTRING3) alen = strlen(str);
457 fprintf(fp, " %6d %6d %s\n", atype, alen, name);
458 if (atype == ATTRSTRING3) {
459 fprintf(fp, " %s\n", str);
460 } else if (atype == ATTRREAL2) {
461 for (j = 0; j < alen; j++) {
462 fprintf(fp, " %20.13le", reals[j]);
463 if ((j+1)%4 == 0) fprintf(fp,"\n");
464 }
465 if (j%4 != 0) fprintf(fp, "\n");
466 } else {
467 for (j = 0; j < alen; j++) {
468 fprintf(fp, " %10d", ints[j]);
469 if ((j+1)%8 == 0) fprintf(fp,"\n");
470 }
471 if (j%8 != 0) fprintf(fp, "\n");
472 }
473 }
474
475}
476
477
478static void getNodeSpacing(ego tess, int nedge, ego *edges, double *spacing)
479{
480 int i, j, m, n, nvert, oclass, mtype, status, *senses;
481 double dist, range[2];
482 const double *xyzs, *ts;
483 ego body, geom, *objs;
484
485 /* compute the smallest segment spacing touching a Node from the Edge
486 tessellations */
487
488 status = EG_statusTessBody(tess, &body, &n, &nvert);
489 if (status != EGADS_SUCCESS0) {
490 printf(" EG_statusTessBody = %d (getNodeSpacing)!\n", status);
491 return;
492 }
493
494 for (i = 0; i < nedge; i++) {
495 status = EG_getTopology(edges[i], &geom, &oclass, &mtype, range, &n, &objs,
496 &senses);
497 if (status != EGADS_SUCCESS0) {
498 printf(" %d: EG_getTopology = %d (getNodeSpacing)!\n", i+1, status);
499 continue;
500 }
501 if (mtype == DEGENERATE5) continue;
502
503 status = EG_getTessEdge(tess, i+1, &m, &xyzs, &ts);
504 if (status != EGADS_SUCCESS0) {
505 printf(" %d: EG_getTessEdge = %d (getNodeSpacing)!\n", i+1, status);
506 continue;
507 }
508 if (m == 0) continue;
509 j = EG_indexBodyTopo(body, objs[0]);
510 if (j <= EGADS_SUCCESS0) {
511 printf(" %d: EG_indexBodyTopo 0 = %d (getNodeSpacing)!\n", i+1, status);
512 continue;
513 }
514 dist = sqrt((xyzs[0]-xyzs[3])*(xyzs[0]-xyzs[3]) +
515 (xyzs[1]-xyzs[4])*(xyzs[1]-xyzs[4]) +
516 (xyzs[2]-xyzs[5])*(xyzs[2]-xyzs[5]));
517 if (spacing[j-1] == 0.0) {
518 spacing[j-1] = dist;
519 } else {
520 if (dist < spacing[j-1]) spacing[j-1] = dist;
521 }
522
523 if (mtype == TWONODE2) j = EG_indexBodyTopo(body, objs[1]);
524 if (j <= EGADS_SUCCESS0) {
525 printf(" %d: EG_indexBodyTopo 1 = %d (getNodeSpacing)!\n", i+1, status);
526 continue;
527 }
528 dist = sqrt((xyzs[3*m-3]-xyzs[3*m-6])*(xyzs[3*m-3]-xyzs[3*m-6]) +
529 (xyzs[3*m-2]-xyzs[3*m-5])*(xyzs[3*m-2]-xyzs[3*m-5]) +
530 (xyzs[3*m-1]-xyzs[3*m-4])*(xyzs[3*m-1]-xyzs[3*m-4]));
531 if (spacing[j-1] == 0.0) {
532 spacing[j-1] = dist;
533 } else {
534 if (dist < spacing[j-1]) spacing[j-1] = dist;
535 }
536 }
537}
538
539
540static void fillIn(double parm, double delta, int *n, double *parms)
541{
542 int i, ii, hit;
543
544 /* adjust the parameter sequence by inserting Node positions and spacing
545 before, after or both */
546
547 for (i = 0; i < *n; i++)
548 if (parms[i] > parm) break;
549
550 ii = i-1;
551 if (ii == -1) ii = 0;
552 if (i == *n) i = *n-1;
553#ifdef DEBUG
554 printf(" %lf %lf %lf", parms[ii], parm, parms[i]);
555#endif
556 /* do we insert the value? */
557 if ((fabs(parm-parms[ii]) <= delta) || (fabs(parm-parms[i]) <= delta)) {
558 hit = i;
559 if (fabs(parm-parms[ii]) <= delta) hit = ii;
560#ifdef DEBUG
561 printf(" Node @ %d", hit);
562#endif
563 } else {
564 if (*n == MAXUVS1024) return;
565 hit = i;
566 for (i = *n; i >= hit; i--) parms[i] = parms[i-1];
567 parms[hit] = parm;
568 *n += 1;
569#ifdef DEBUG
570 printf(" Add @ %d", hit);
571#endif
572 }
573
574 /* look forward */
575 if (hit != *n-1)
576 if (parms[hit+1]-parms[hit] > 2.0*delta) {
577 if (*n == MAXUVS1024) return;
578 for (i = *n; i >= hit+1; i--) parms[i] = parms[i-1];
579 parms[hit+1] = parms[hit] + delta;
580 *n += 1;
581#ifdef DEBUG
582 printf(" Add+");
583#endif
584 }
585
586 /* look back */
587 if (hit != 0)
588 if (parms[hit]-parms[hit-1] > 2.0*delta) {
589 if (*n == MAXUVS1024) return;
590 for (i = *n; i >= hit+1; i--) parms[i] = parms[i-1];
591 parms[hit] -= delta;
592 *n += 1;
593#ifdef DEBUG
594 printf(" Add-");
595#endif
596 }
597#ifdef DEBUG
598 printf("\n");
599#endif
600}
601
602
603static void insertNodeSpacing(ego tess, ego face, int iface, double *spacing,
604 UVmap *map)
605{
606 int i, m, n, status;
607 double udist, vdist, result[18], *u1, *v1;
608 const int *pindex, *ptype, *tris, *tric;
609 const double *xyzs, *uvs;
610
611 u1 = &result[3];
612 v1 = &result[6];
613
614 status = EG_getTessFace(tess, iface, &m, &xyzs, &uvs, &ptype, &pindex,
615 &n, &tris, &tric);
616 if (status != EGADS_SUCCESS0) {
617 printf(" %d: EG_getTessFace = %d (insertNodeSpacing)!\n", iface, status);
618 return;
619 }
620 for (i = 0; i < m; i++) {
621 if (ptype[i] != 0) continue;
622#ifdef DEBUG
623 printf(" Face %d: Node = %d, spacing = %lf uv = %lf %lf\n",
624 iface, pindex[i], spacing[pindex[i]-1], uvs[2*i ], uvs[2*i+1]);
625#endif
626 status = EG_evaluate(face, &uvs[2*i], result);
627 if (status != EGADS_SUCCESS0) {
628 printf(" %d: EG_evaluate = %d (insertNodeSpacing)!\n", iface, status);
629 return;
630 }
631 udist = sqrt(u1[0]*u1[0] + u1[1]*u1[1] + u1[2]*u1[2]);
632 if (udist != 0.0) {
633 udist = spacing[pindex[i]-1]/udist;
634 fillIn(uvs[2*i ], udist, &map->nu, map->us);
635 if (map->nu == MAXUVS1024) break;
636 }
637 vdist = sqrt(v1[0]*v1[0] + v1[1]*v1[1] + v1[2]*v1[2]);
638 if (vdist != 0.0) {
639 vdist = spacing[pindex[i]-1]/vdist;
640 fillIn(uvs[2*i+1], vdist, &map->nv, map->vs);
641 }
642#ifdef DEBUG
643 printf(" udelta = %lf vdelta = %lf\n", udist, vdist);
644#endif
645 }
646}
647
648
649static void dumpFront(int nSeg, parmSeg *segs, parmSeg *front)
650{
651 int i;
652 parmSeg *link;
653
654 link = front;
655 while (link != NULL((void*)0)) {
656 printf(" F %p parms = %lf %lf, len = %lf\n size = %lf connect = %p\n",
657 link, link->parms[0], link->parms[1], link->alen, link->size,
658 link->prev);
659 link = link->next;
660 }
661
662 for (i = 0; i < nSeg; i++) {
663 link = &segs[i];
664 printf(" S %p parms = %lf %lf, size = %lf\n prev = %p next = %p\n",
665 link, link->parms[0], link->parms[1], link->size,
666 link->prev, link->next);
667 }
668}
669
670/*@-immediatetrans@*/
671static void patchSeg(int nSeg, parmSeg *segs, parmSeg **frontx,
672 parmSeg **poolx, parmSeg *fill)
673{
674 int i;
675 parmSeg *front, *pool, *last, *link;
676
677 /* look in our current segments */
678 for (i = 0; i < nSeg; i++) {
679 if (fabs(fill->parms[0] - segs[i].parms[1]) < FUZZ1.e-14) fill->prev = &segs[i];
680 if (fabs(fill->parms[1] - segs[i].parms[0]) < FUZZ1.e-14) fill->next = &segs[i];
681 }
682
683 front = *frontx;
684 pool = *poolx;
685
686 link = front;
687 while (link != NULL((void*)0)) {
688 last = link;
689 link = link->next;
690 if (last->parms[1] < 0) {
691 if (fabs(fill->parms[1] - last->parms[0]) < FUZZ1.e-14)
692 delFront(&front, &pool, last);
693 } else {
694 if (fabs(fill->parms[0] - last->parms[0]) < FUZZ1.e-14)
695 delFront(&front, &pool, last);
696 }
697 }
698
699 *frontx = front;
700 *poolx = pool;
701}
702/*@+immediatetrans@*/
703
704static void smoothParm(int *np, double *parms, double *r, double *q, double fact)
705{
706 int i, j, k, m, n, mark[MAXUVS1024], nSeg = 0;
707 double smals, ave, alen, dist, qi, pi, size, fra, params[2], al[MAXUVS1024];
708 parmSeg segs[MAXSEG1023], seg, *link, *smallest, *other;
709 parmSeg *first = NULL((void*)0), *front = NULL((void*)0), *pool = NULL((void*)0);
710
711 n = *np;
712 if (n <= 2) return;
713 al[0] = 0.0;
714 for (i = 0; i < n-1; i++) {
715 mark[i] = 1; /* available */
716 al[i+1] = al[i] + r[i]; /* compute segment arclength */
717 }
718
719 /* get the smallest & average segment */
720 smals = ave = r[0];
721 j = 0;
722 for (i = 1; i < n-1; i++) {
723 ave += r[i];
724 if (r[i] < smals) {
725 j = i;
726 smals = r[i];
727 }
728 }
729 ave /= (n-1);
730 if (fabs(smals-ave)/ave < 0.1) return;
731
732 /* add the first segment */
733 seg.prev = NULL((void*)0);
734 seg.parms[0] = parms[j];
735 seg.parms[1] = parms[j+1];
736 seg.alen = 0.0;
737 seg.size = r[j];
738 seg.next = NULL((void*)0);
739 addSeg(&nSeg, segs, seg);
740 mark[j] = -1; /* used */
741
742 /* initialize the front */
743 if (j == 0) {
744 params[0] = parms[j+1];
745 params[1] = +1.0;
746 addFront(&front, &pool, params, al[1], r[0], &segs[0]);
747 first = &segs[0];
748 params[0] = parms[n-1];
749 params[1] = -1.0;
750 addFront(&front, &pool, params, al[n-1], r[n-2], NULL((void*)0));
751 } else if (j == n-2) {
752 params[0] = parms[0];
753 params[1] = +1.0;
754 addFront(&front, &pool, params, al[0], r[0], NULL((void*)0));
755 params[0] = parms[n-2];
756 params[1] = -1.0;
757 addFront(&front, &pool, params, al[j], r[j], &segs[0]);
758 } else {
759 params[0] = parms[0];
760 params[1] = +1.0;
761 addFront(&front, &pool, params, al[0], r[0], NULL((void*)0));
762 params[0] = parms[j];
763 params[1] = -1.0;
764 addFront(&front, &pool, params, al[j], r[j], &segs[0]);
765 params[0] = parms[j+1];
766 params[1] = +1.0;
767 addFront(&front, &pool, params, al[j+1], r[j], &segs[0]);
768 params[0] = parms[n-1];
769 params[1] = -1.0;
770 addFront(&front, &pool, params, al[n-1], r[n-2], NULL((void*)0));
771 }
772
773 do {
774 smals *= fact;
775 /* first look at available original segments */
776 j = -1;
777 for (i = 0; i < n-1; i++)
778 if (mark[i] == 1)
779 if (r[i] < smals)
780 if (j == -1) {
781 ave = r[i];
782 j = i;
783 } else {
784 if (ave > r[i]) {
785 ave = r[i];
786 j = i;
787 }
788 }
789#ifdef DEBUG
790 printf(" next = %d ave = %lf small = %lf nSeg = %d\n",
791 j, ave, smals, nSeg);
792#endif
793
794 if (j != -1) {
795 /* add the next segment */
796 seg.prev = NULL((void*)0);
797 seg.parms[0] = parms[j];
798 seg.parms[1] = parms[j+1];
799 seg.alen = 0.0;
800 seg.size = r[j];
801 seg.next = NULL((void*)0);
802 patchSeg(nSeg, segs, &front, &pool, &seg);
803 addSeg(&nSeg, segs, seg);
804 if (seg.parms[0] == parms[0]) first = &segs[nSeg-1];
805 mark[j] = -1; /* used */
806 smals = r[j];
807 if ((seg.prev == NULL((void*)0)) && (j != 0)) {
808 params[0] = parms[j];
809 params[1] = -1.0;
810 addFront(&front, &pool, params, al[j], r[j], &segs[nSeg-1]);
811 }
812 if ((seg.next == NULL((void*)0)) && (j != n-2)) {
813 params[0] = parms[j+1];
814 params[1] = +1.0;
815 addFront(&front, &pool, params, al[j+1], r[j], &segs[nSeg-1]);
816 }
817 continue;
818 }
819 if (front == NULL((void*)0)) continue;
820
821 /* now find the smallest front segment */
822 smallest = link = front;
823 while (link != NULL((void*)0)) {
824 if (smallest->size > link->size) smallest = link;
825 link = link->next;
826 }
827
828 /* look ahead for fronts colliding */
829 link = front;
830 ave = al[n-1];
831 other = NULL((void*)0);
832 while (link != NULL((void*)0)) {
833 if (smallest->parms[1]*link->parms[1] < 0.0) {
834 if (smallest->parms[1] < 0) {
835 if (link->parms[0] < smallest->parms[0]) {
836 alen = fabs(smallest->alen - link->alen);
837 if (alen < ave) {
838 other = link;
839 ave = alen;
840 }
841 }
842 } else {
843 if (link->parms[0] > smallest->parms[0]) {
844 alen = fabs(smallest->alen - link->alen);
845 if (alen < ave) {
846 other = link;
847 ave = alen;
848 }
849 }
850 }
851 }
852 link = link->next;
853 }
854 if (other != NULL((void*)0)) {
855 dist = 0.5*(other->size + smallest->size);
856 m = ave/dist + .49;
857 alen = smallest->size/other->size;
858 if (alen < 1.0) alen = 1.0/alen;
859#ifdef DEBUG
860 printf(" alen = %lf delta = %lf divs = %d ratio = %lf\n",
861 ave, dist, m, alen);
862#endif
863 if (m == 0) m = 1;
864 if (alen > 10.0) m = 0;
865 if ((m > 0) && (m < 10)) { /* 10 is arbitrary */
866 if (smallest->parms[0] > other->parms[0]) {
867 link = smallest;
868 smallest = other;
869 other = link;
870 }
871
872 /* mark original segments */
873 alen = smallest->alen;
874 for (j = 0; j < n-1; j++)
875 if ((alen >= al[j]) && (alen <= al[j+1])) break;
876 if (j == n-1) {
877 printf(" ERROR: Cannot Interpolate 0 %lf %lf -- collide!\n",
878 alen, al[n-1]);
879 exit(EXIT_FAILURE1);
880 }
881 dist = (alen - al[j])/(al[j+1] - al[j]);
882 if ((fabs(dist-1.0) < FUZZ1.e-14) && (mark[j] == -1) && (j != n-2)) {
883 j++;
884 dist = 0.0;
885 }
886 alen = other->alen;
887 for (k = 0; k < n-1; k++)
888 if ((alen >= al[k]) && (alen <= al[k+1])) break;
889 if (k == n-1) {
890 if (fabs(alen-al[n-1]) > FUZZ1.e-14) {
891 printf(" ERROR: Cannot Interpolate 1 %lf %lf -- collide!\n",
892 alen, al[n-1]);
893 exit(EXIT_FAILURE1);
894 } else {
895 k--;
896 }
897 }
898 dist = (alen - al[k])/(al[k+1] - al[k]);
899 if ((fabs(dist) < FUZZ1.e-14) && (mark[k] == -1) && (k != 0)) {
900 k--;
901 dist = 1.0;
902 }
903 for (i = j; i <= k; i++) {
904 if (mark[i] == -1) {
905 printf(" ERROR: Used Segment %d -- %d %d -- collide!\n", i, j, k);
906 exit(EXIT_FAILURE1);
907 }
908 mark[i] = 0;
909 }
910
911 /* make new segments -- scale parameter by arcLength */
912 alen = 0.0;
913 qi = smallest->size;
914 for (i = 0; i < m; i++) {
915 ave = i+1;
916 ave /= m;
917 size = (1.0-ave)*smallest->size + ave*other->size;
918 alen += 0.5*(qi + size);
919 qi = size;
920 }
921 qi = smallest->size;
922 pi = smallest->parms[0];
923 seg.alen = dist = 0.0;
924 for (i = 0; i < m; i++) {
925 ave = i+1;
926 ave /= m;
927 size = (1.0-ave)*smallest->size + ave*other->size;
928 dist += 0.5*(qi + size);
929 fra = dist/alen;
930 seg.size = 0.5*(qi + size);
931 seg.parms[0] = pi;
932 seg.parms[1] = (1.0-fra)*smallest->parms[0] + fra*other->parms[0];
933 seg.prev = NULL((void*)0);
934 seg.next = NULL((void*)0);
935 qi = size;
936 pi = seg.parms[1];
937 patchSeg(nSeg, segs, &front, &pool, &seg);
938 addSeg(&nSeg, segs, seg);
939 if (seg.parms[0] == parms[0]) first = &segs[nSeg-1];
940 }
941 smals /= fact;
942 continue;
943 }
944 }
945
946 /* add a single segment from the front's smallest */
947 smals = fact*smallest->size;
948 alen = smallest->parms[1]*smals + smallest->alen;
949 for (j = 0; j < n-1; j++)
950 if ((alen >= al[j]) && (alen <= al[j+1])) break;
951 if (j == n-1) {
952 printf(" ERROR: Cannot Interpolate alen = %lf [%lf %lf] %p!\n",
953 alen, al[0], al[n-1], smallest);
954 dumpFront(nSeg, segs, front);
955 exit(EXIT_FAILURE1);
956 }
957 dist = (alen - al[j])/(al[j+1] - al[j]);
958 if (mark[j] == -1) {
959 printf(" ERROR: Used Segment -- dist = %lf!\n", dist);
960 exit(EXIT_FAILURE1);
961 }
962 mark[j] = 0;
963 qi = q[j] + dist*(q[j+1] - q[j]);
964 pi = smals/qi;
965 seg.prev = NULL((void*)0);
966 seg.alen = 0.0;
967 seg.size = smals;
968 seg.next = NULL((void*)0);
969 params[1] = smallest->parms[1];
970 if (params[1] > 0) {
971 seg.parms[0] = smallest->parms[0];
972 seg.parms[1] = params[0] = smallest->parms[0] + pi;
973 } else {
974 seg.parms[0] = params[0] = smallest->parms[0] - pi;
975 seg.parms[1] = smallest->parms[0];
976 }
977 patchSeg(nSeg, segs, &front, &pool, &seg);
978 addSeg(&nSeg, segs, seg);
979 if (seg.parms[0] == parms[0]) first = &segs[nSeg-1];
980 addFront(&front, &pool, params, alen, smals, &segs[nSeg-1]);
981
982 } while (front != NULL((void*)0));
983
984 /* reset our parameter sequence from the linked list */
985 if (first == NULL((void*)0)) {
986 printf(" ERROR: No First Segment!\n");
987 dumpFront(nSeg, segs, front);
988 exit(EXIT_FAILURE1);
989 }
990 link = first;
991 i = 0;
992 pi = parms[n-1];
993 while (link != NULL((void*)0)) {
994 if (i == 0) parms[0] = link->parms[0];
995 parms[i+1] = link->parms[1];
996 i++;
997 link = link->next;
998 }
999 if ((i != nSeg) || (parms[i] != pi)) {
1000 printf(" ERROR: Finialization -- %d %d %lf %lf!\n", i, nSeg, parms[i], pi);
1001 dumpFront(nSeg, segs, front);
1002 exit(EXIT_FAILURE1);
1003 }
1004 *np = nSeg+1;
1005
1006 freeFront(front, pool);
1007}
1008
1009
1010static void smoothMap(ego face, double factor, UVmap *map)
1011{
1012 int i, j, stat;
1013 double uv[2], *u1, *v1, r[MAXUVS1024], q[MAXUVS1024], result[18];
1014 double dp[3*MAXUVS1024], xyzs[3*MAXUVS1024];
1015
1016 if (factor < 1.0) return;
1017 u1 = &result[3];
1018 v1 = &result[6];
1019
1020 /* look in the U direction */
1021 for (i = 0; i < map->nu; i++) r[i] = q[i] = 0.0;
1022 for (j = 0; j < map->nv; j++) {
1023 uv[1] = map->vs[j];
1024 for (i = 0; i < map->nu; i++) {
1025 uv[0] = map->us[i];
1026 stat = EG_evaluate(face, uv, result);
1027 if (stat != EGADS_SUCCESS0) {
1028 printf(" smoothMap: Fill U EG_evaluate = %d!\n", stat);
1029 return;
1030 }
1031 xyzs[3*i ] = result[0];
1032 xyzs[3*i+1] = result[1];
1033 xyzs[3*i+2] = result[2];
1034 dp[3*i ] = u1[0];
1035 dp[3*i+1] = u1[1];
1036 dp[3*i+2] = u1[2];
1037 }
1038 q[0] += sqrt(dp[0]*dp[0] + dp[1]*dp[1] + dp[2]*dp[2]);
1039 for (i = 1; i < map->nu; i++) {
1040 q[i ] += sqrt(dp[3*i ]*dp[3*i ] + dp[3*i+1]*dp[3*i+1] +
1041 dp[3*i+2]*dp[3*i+2]);
1042 r[i-1] += sqrt((xyzs[3*i-3]-xyzs[3*i ])*(xyzs[3*i-3]-xyzs[3*i ]) +
1043 (xyzs[3*i-2]-xyzs[3*i+1])*(xyzs[3*i-2]-xyzs[3*i+1]) +
1044 (xyzs[3*i-1]-xyzs[3*i+2])*(xyzs[3*i-1]-xyzs[3*i+2]));
1045 }
1046 }
1047 q[0] /= map->nv;
1048 for (i = 0; i < map->nu-1; i++) {
1049 q[i+1] /= map->nv;
1050 r[i] /= map->nv;
1051 }
1052#ifdef DEBUG
1053 printf(" nUs = %d\n %lf %lf\n", map->nu, map->us[0], q[0]);
1054 for (i = 1; i < map->nu; i++)
1055 printf(" %lf %lf %lf\n", map->us[i], r[i-1], q[i]);
1056#endif
1057 smoothParm(&map->nu, map->us, r, q, factor);
1058#ifdef DEBUG
1059 for (i = 0; i < map->nu; i++) r[i] = q[i] = 0.0;
1060 for (j = 0; j < map->nv; j++) {
1061 uv[1] = map->vs[j];
1062 for (i = 0; i < map->nu; i++) {
1063 uv[0] = map->us[i];
1064 stat = EG_evaluate(face, uv, result);
1065 if (stat != EGADS_SUCCESS0) {
1066 printf(" smoothMap: Fill U EG_evaluate = %d!\n", stat);
1067 return;
1068 }
1069 xyzs[3*i ] = result[0];
1070 xyzs[3*i+1] = result[1];
1071 xyzs[3*i+2] = result[2];
1072 dp[3*i ] = u1[0];
1073 dp[3*i+1] = u1[1];
1074 dp[3*i+2] = u1[2];
1075 }
1076 q[0] += sqrt(dp[0]*dp[0] + dp[1]*dp[1] + dp[2]*dp[2]);
1077 for (i = 1; i < map->nu; i++) {
1078 q[i ] += sqrt(dp[3*i ]*dp[3*i ] + dp[3*i+1]*dp[3*i+1] +
1079 dp[3*i+2]*dp[3*i+2]);
1080 r[i-1] += sqrt((xyzs[3*i-3]-xyzs[3*i ])*(xyzs[3*i-3]-xyzs[3*i ]) +
1081 (xyzs[3*i-2]-xyzs[3*i+1])*(xyzs[3*i-2]-xyzs[3*i+1]) +
1082 (xyzs[3*i-1]-xyzs[3*i+2])*(xyzs[3*i-1]-xyzs[3*i+2]));
1083 }
1084 }
1085 q[0] /= map->nv;
1086 for (i = 0; i < map->nu-1; i++) {
1087 q[i+1] /= map->nv;
1088 r[i] /= map->nv;
1089 }
1090 printf(" nUs = %d\n %lf %lf\n", map->nu, map->us[0], q[0]);
1091 r[map->nu-1] = r[map->nu-2];
1092 for (i = 1; i < map->nu; i++)
1093 printf(" %lf %lf %lf %lf\n", map->us[i], r[i-1], q[i], r[i]/r[i-1]);
1094#endif
1095
1096 /* look in the V direction */
1097 for (j = 0; j < map->nv; j++) r[j] = q[j] = 0.0;
1098 for (i = 0; i < map->nu; i++) {
1099 uv[0] = map->us[i];
1100 for (j = 0; j < map->nv; j++) {
1101 uv[1] = map->vs[j];
1102 stat = EG_evaluate(face, uv, result);
1103 if (stat != EGADS_SUCCESS0) {
1104 printf(" smoothMap: Fill V EG_evaluate = %d!\n", stat);
1105 return;
1106 }
1107 xyzs[3*j ] = result[0];
1108 xyzs[3*j+1] = result[1];
1109 xyzs[3*j+2] = result[2];
1110 dp[3*j ] = v1[0];
1111 dp[3*j+1] = v1[1];
1112 dp[3*j+2] = v1[2];
1113 }
1114 q[0] += sqrt(dp[0]*dp[0] + dp[1]*dp[1] + dp[2]*dp[2]);
1115 for (j = 1; j < map->nv; j++) {
1116 q[j ] += sqrt(dp[3*j ]*dp[3*j ] + dp[3*j+1]*dp[3*j+1] +
1117 dp[3*j+2]*dp[3*j+2]);
1118 r[j-1] += sqrt((xyzs[3*j-3]-xyzs[3*j ])*(xyzs[3*j-3]-xyzs[3*j ]) +
1119 (xyzs[3*j-2]-xyzs[3*j+1])*(xyzs[3*j-2]-xyzs[3*j+1]) +
1120 (xyzs[3*j-1]-xyzs[3*j+2])*(xyzs[3*j-1]-xyzs[3*j+2]));
1121 }
1122 }
1123 q[0] /= map->nu;
1124 for (j = 0; j < map->nv-1; j++) {
1125 q[j+1] /= map->nu;
1126 r[j] /= map->nu;
1127 }
1128#ifdef DEBUG
1129 printf(" nVs = %d\n %lf %lf\n", map->nv, map->vs[0], q[0]);
1130 for (j = 1; j < map->nv; j++)
1131 printf(" %lf %lf %lf\n", map->vs[j], r[j-1], q[j]);
1132#endif
1133 smoothParm(&map->nv, map->vs, r, q, factor);
1134#ifdef DEBUG
1135 for (j = 0; j < map->nv; j++) r[j] = q[j] = 0.0;
1136 for (i = 0; i < map->nu; i++) {
1137 uv[0] = map->us[i];
1138 for (j = 0; j < map->nv; j++) {
1139 uv[1] = map->vs[j];
1140 stat = EG_evaluate(face, uv, result);
1141 if (stat != EGADS_SUCCESS0) {
1142 printf(" smoothMap: Fill V EG_evaluate = %d!\n", stat);
1143 return;
1144 }
1145 xyzs[3*j ] = result[0];
1146 xyzs[3*j+1] = result[1];
1147 xyzs[3*j+2] = result[2];
1148 dp[3*j ] = v1[0];
1149 dp[3*j+1] = v1[1];
1150 dp[3*j+2] = v1[2];
1151 }
1152 q[0] += sqrt(dp[0]*dp[0] + dp[1]*dp[1] + dp[2]*dp[2]);
1153 for (j = 1; j < map->nv; j++) {
1154 q[j ] += sqrt(dp[3*j ]*dp[3*j ] + dp[3*j+1]*dp[3*j+1] +
1155 dp[3*j+2]*dp[3*j+2]);
1156 r[j-1] += sqrt((xyzs[3*j-3]-xyzs[3*j ])*(xyzs[3*j-3]-xyzs[3*j ]) +
1157 (xyzs[3*j-2]-xyzs[3*j+1])*(xyzs[3*j-2]-xyzs[3*j+1]) +
1158 (xyzs[3*j-1]-xyzs[3*j+2])*(xyzs[3*j-1]-xyzs[3*j+2]));
1159 }
1160 }
1161 q[0] /= map->nu;
1162 for (j = 0; j < map->nv-1; j++) {
1163 q[j+1] /= map->nu;
1164 r[j] /= map->nu;
1165 }
1166 printf(" nVs = %d\n %lf %lf\n", map->nv, map->vs[0], q[0]);
1167 r[map->nv-1] = r[map->nv-2];
1168 for (j = 1; j < map->nv; j++)
1169 printf(" %lf %lf %lf %lf\n", map->vs[j], r[j-1], q[j], r[j]/r[j-1]);
1170#endif
1171}
1172
1173
1174static void updateMap(ego face, double mxedg, double sag, double angle, UVmap *map)
1175{
1176 int i, j, max, stat, cnt;
1177 double d, dist, dot, last, uv[2], xyz[3], *u1, *v1, result[18], xyzs[3*MAXUVS1024];
1178
1179 dot = cos(PI3.1415926535897931159979635*angle/180.0);
1180 u1 = &result[3];
1181 v1 = &result[6];
1182
1183 /* sag -- look in U direction */
1184 do {
1185 max = -1;
1186 dist = 0.0;
1187 for (j = 0; j < map->nv; j++) {
1188 uv[1] = map->vs[j];
1189 for (i = 0; i < map->nu; i++) {
1190 uv[0] = map->us[i];
1191 stat = EG_evaluate(face, uv, result);
1192 if (stat != EGADS_SUCCESS0) {
1193 printf(" updateMap: Fill U EG_evaluate = %d!\n", stat);
1194 return;
1195 }
1196 xyzs[3*i ] = result[0];
1197 xyzs[3*i+1] = result[1];
1198 xyzs[3*i+2] = result[2];
1199 }
1200 for (i = 1; i < map->nu; i++) {
1201 xyz[0] = 0.5*(xyzs[3*i-3] + xyzs[3*i ]);
1202 xyz[1] = 0.5*(xyzs[3*i-2] + xyzs[3*i+1]);
1203 xyz[2] = 0.5*(xyzs[3*i-1] + xyzs[3*i+2]);
1204 uv[0] = 0.5*(map->us[i-1] + map->us[i]);
1205 stat = EG_evaluate(face, uv, result);
1206 if (stat != EGADS_SUCCESS0) {
1207 printf(" updateMap: Half U EG_evaluate = %d!\n", stat);
1208 return;
1209 }
1210 d = sqrt((xyz[0]-result[0])*(xyz[0]-result[0]) +
1211 (xyz[1]-result[1])*(xyz[1]-result[1]) +
1212 (xyz[2]-result[2])*(xyz[2]-result[2]));
1213 if (d <= sag) continue;
1214 if (d > dist) {
1215 max = i;
1216 dist = d;
1217 }
1218 }
1219 }
1220 if (map->nu == MAXUVS1024) break;
1221 if (max != -1) {
1222 for (i = map->nu; i >= max; i--) map->us[i] = map->us[i-1];
1223 map->us[max] = 0.5*(map->us[max+1] + map->us[max-1]);
1224 map->nu++;
1225 }
1226 } while (max != -1);
1227
1228 /* sag -- look in V direction */
1229 do {
1230 max = -1;
1231 dist = 0.0;
1232 for (i = 0; i < map->nu; i++) {
1233 uv[0] = map->us[i];
1234 for (j = 0; j < map->nv; j++) {
1235 uv[1] = map->vs[j];
1236 stat = EG_evaluate(face, uv, result);
1237 if (stat != EGADS_SUCCESS0) {
1238 printf(" updateMap: Fill V EG_evaluate = %d!\n", stat);
1239 return;
1240 }
1241 xyzs[3*j ] = result[0];
1242 xyzs[3*j+1] = result[1];
1243 xyzs[3*j+2] = result[2];
1244 }
1245 for (j = 1; j < map->nv; j++) {
1246 xyz[0] = 0.5*(xyzs[3*j-3] + xyzs[3*j ]);
1247 xyz[1] = 0.5*(xyzs[3*j-2] + xyzs[3*j+1]);
1248 xyz[2] = 0.5*(xyzs[3*j-1] + xyzs[3*j+2]);
1249 uv[1] = 0.5*(map->vs[j-1] + map->vs[j]);
1250 stat = EG_evaluate(face, uv, result);
1251 if (stat != EGADS_SUCCESS0) {
1252 printf(" updateMap: Half V EG_evaluate = %d!\n", stat);
1253 return;
1254 }
1255 d = sqrt((xyz[0]-result[0])*(xyz[0]-result[0]) +
1256 (xyz[1]-result[1])*(xyz[1]-result[1]) +
1257 (xyz[2]-result[2])*(xyz[2]-result[2]));
1258 if (d <= sag) continue;
1259 if (d > dist) {
1260 max = j;
1261 dist = d;
1262 }
1263 }
1264 }
1265 if (map->nv == MAXUVS1024) break;
1266 if (max != -1) {
1267 for (j = map->nv; j >= max; j--) map->vs[j] = map->vs[j-1];
1268 map->vs[max] = 0.5*(map->vs[max+1] + map->vs[max-1]);
1269 map->nv++;
1270 }
1271 } while (max != -1);
1272
1273 /* angle -- look in U direction */
1274 last = -1.0;
1275 cnt = 0;
1276 do {
1277 max = -1;
1278 dist = 1.0;
1279 for (j = 0; j < map->nv; j++) {
1280 uv[1] = map->vs[j];
1281 for (i = 0; i < map->nu; i++) {
1282 uv[0] = map->us[i];
1283 stat = EG_evaluate(face, uv, result);
1284 if (stat != EGADS_SUCCESS0) {
1285 printf(" updateMap: Fill Ua EG_evaluate = %d!\n", stat);
1286 return;
1287 }
1288 CROSS(xyz, u1, v1)xyz[0] = (u1[1]*v1[2]) - (u1[2]*v1[1]); xyz[1] = (u1[2]*v1[0]
) - (u1[0]*v1[2]); xyz[2] = (u1[0]*v1[1]) - (u1[1]*v1[0])
;
1289 d = sqrt(xyz[0]*xyz[0] + xyz[1]*xyz[1] + xyz[2]*xyz[2]);
1290 if (d == 0) d = 1.0;
1291 xyzs[3*i ] = xyz[0]/d;
1292 xyzs[3*i+1] = xyz[1]/d;
1293 xyzs[3*i+2] = xyz[2]/d;
1294 }
1295 for (i = 1; i < map->nu; i++) {
1296 if (sqrt(xyzs[3*i-3]*xyzs[3*i-3] + xyzs[3*i-2]*xyzs[3*i-2] +
1297 xyzs[3*i-1]*xyzs[3*i-1]) == 0.0) continue;
1298 if (sqrt(xyzs[3*i ]*xyzs[3*i ] + xyzs[3*i+1]*xyzs[3*i+1] +
1299 xyzs[3*i+2]*xyzs[3*i+2]) == 0.0) continue;
1300 d = xyzs[3*i ]*xyzs[3*i-3] + xyzs[3*i+1]*xyzs[3*i-2] +
1301 xyzs[3*i+2]*xyzs[3*i-1];
1302 if (d >= dot) continue;
1303 if (d < -0.2) continue; /* there is some problem - C0? */
1304 if (d < dist) {
1305 max = i;
1306 dist = d;
1307 }
1308 }
1309 }
1310 if (map->nu == MAXUVS1024) break;
1311 if (max != -1) {
1312/* printf(" Umax = %d (%d) %lf (%lf)\n", max, map->nu, dist, dot); */
1313 for (i = map->nu; i >= max; i--) map->us[i] = map->us[i-1];
1314 map->us[max] = 0.5*(map->us[max+1] + map->us[max-1]);
1315 map->nu++;
1316 }
1317 if (dist >= last) cnt++;
1318 last = dist;
1319 } while ((max != -1) && (cnt < 10)) ;
1320
1321 /* angle -- look in V direction */
1322 last = -1.0;
1323 cnt = 0;
1324 do {
1325 max = -1;
1326 dist = 1.0;
1327 for (i = 0; i < map->nu; i++) {
1328 uv[0] = map->us[i];
1329 for (j = 0; j < map->nv; j++) {
1330 uv[1] = map->vs[j];
1331 stat = EG_evaluate(face, uv, result);
1332 if (stat != EGADS_SUCCESS0) {
1333 printf(" updateMap: Fill Va EG_evaluate = %d!\n", stat);
1334 return;
1335 }
1336 CROSS(xyz, u1, v1)xyz[0] = (u1[1]*v1[2]) - (u1[2]*v1[1]); xyz[1] = (u1[2]*v1[0]
) - (u1[0]*v1[2]); xyz[2] = (u1[0]*v1[1]) - (u1[1]*v1[0])
;
1337 d = sqrt(xyz[0]*xyz[0] + xyz[1]*xyz[1] + xyz[2]*xyz[2]);
1338 if (d == 0) d = 1.0;
1339 xyzs[3*j ] = xyz[0]/d;
1340 xyzs[3*j+1] = xyz[1]/d;
1341 xyzs[3*j+2] = xyz[2]/d;
1342 }
1343 for (j = 1; j < map->nv; j++) {
1344 if (sqrt(xyzs[3*j-3]*xyzs[3*j-3] + xyzs[3*j-2]*xyzs[3*j-2] +
1345 xyzs[3*j-1]*xyzs[3*j-1]) == 0.0) continue;
1346 if (sqrt(xyzs[3*j ]*xyzs[3*j ] + xyzs[3*j+1]*xyzs[3*j+1] +
1347 xyzs[3*j+2]*xyzs[3*j+2]) == 0.0) continue;
1348 d = xyzs[3*j ]*xyzs[3*j-3] + xyzs[3*j+1]*xyzs[3*j-2] +
1349 xyzs[3*j+2]*xyzs[3*j-1];
1350 if (d >= dot) continue;
1351 if (d < -0.2) continue; /* there is some problem - C0? */
1352 if (d < dist) {
1353 max = j;
1354 dist = d;
1355 }
1356 }
1357 }
1358 if (map->nv == MAXUVS1024) break;
1359 if (max != -1) {
1360/* printf(" Vmax = %d (%d) %lf (%lf)\n", max, map->nv, dist, dot); */
1361 for (j = map->nv; j >= max; j--) map->vs[j] = map->vs[j-1];
1362 map->vs[max] = 0.5*(map->vs[max+1] + map->vs[max-1]);
1363 map->nv++;
1364 }
1365 if (dist >= last) cnt++;
1366 last = dist;
1367 } while ((max != -1) && (cnt < 10)) ;
1368
1369 /* maxedge -- look in U direction */
1370 do {
1371 max = -1;
1372 dist = 0.0;
1373 for (j = 0; j < map->nv; j++) {
1374 uv[1] = map->vs[j];
1375 for (i = 0; i < map->nu; i++) {
1376 uv[0] = map->us[i];
1377 stat = EG_evaluate(face, uv, result);
1378 if (stat != EGADS_SUCCESS0) {
1379 printf(" updateMap: Fill U EG_evaluate = %d!\n", stat);
1380 return;
1381 }
1382 xyzs[3*i ] = result[0];
1383 xyzs[3*i+1] = result[1];
1384 xyzs[3*i+2] = result[2];
1385 }
1386 for (i = 1; i < map->nu; i++) {
1387 xyz[0] = xyzs[3*i-3] - xyzs[3*i ];
1388 xyz[1] = xyzs[3*i-2] - xyzs[3*i+1];
1389 xyz[2] = xyzs[3*i-1] - xyzs[3*i+2];
1390 d = sqrt(xyz[0]*xyz[0] + xyz[1]*xyz[1] + xyz[2]*xyz[2]);
1391 if (d <= mxedg) continue;
1392 if (d > dist) {
1393 max = i;
1394 dist = d;
1395 }
1396 }
1397 }
1398 if (map->nu == MAXUVS1024) break;
1399 if (max != -1) {
1400 for (i = map->nu; i >= max; i--) map->us[i] = map->us[i-1];
1401 map->us[max] = 0.5*(map->us[max+1] + map->us[max-1]);
1402 map->nu++;
1403 }
1404 } while (max != -1);
1405
1406 /* sag -- look in V direction */
1407 do {
1408 max = -1;
1409 dist = 0.0;
1410 for (i = 0; i < map->nu; i++) {
1411 uv[0] = map->us[i];
1412 for (j = 0; j < map->nv; j++) {
1413 uv[1] = map->vs[j];
1414 stat = EG_evaluate(face, uv, result);
1415 if (stat != EGADS_SUCCESS0) {
1416 printf(" updateMap: Fill V EG_evaluate = %d!\n", stat);
1417 return;
1418 }
1419 xyzs[3*j ] = result[0];
1420 xyzs[3*j+1] = result[1];
1421 xyzs[3*j+2] = result[2];
1422 }
1423 for (j = 1; j < map->nv; j++) {
1424 xyz[0] = xyzs[3*j-3] - xyzs[3*j ];
1425 xyz[1] = xyzs[3*j-2] - xyzs[3*j+1];
1426 xyz[2] = xyzs[3*j-1] - xyzs[3*j+2];
1427 d = sqrt(xyz[0]*xyz[0] + xyz[1]*xyz[1] + xyz[2]*xyz[2]);
1428 if (d <= mxedg) continue;
1429 if (d > dist) {
1430 max = j;
1431 dist = d;
1432 }
1433 }
1434 }
1435 if (map->nv == MAXUVS1024) break;
1436 if (max != -1) {
1437 for (j = map->nv; j >= max; j--) map->vs[j] = map->vs[j-1];
1438 map->vs[max] = 0.5*(map->vs[max+1] + map->vs[max-1]);
1439 map->nv++;
1440 }
1441 } while (max != -1);
1442
1443}
1444
1445
1446static int parse_args(const int argc, char *argv[], char *ifile)
1447{
1448
1449 int i;
1450
1451 printf("********** egads2cgt version %s **********\n", version);
1452
1453 for (i = 1; i <= argc-1; i++) {
1454
1455 if (strcmp("-i", argv[i]) == 0) {
1456 i++;
1457 strcpy(ifile, argv[i]);
1458 } else if (strcmp("-aflr4", argv[i]) == 0) {
1459 aflr4 = 1;
1460 } else if (strcmp("-q", argv[i]) == 0) {
1461 wrtqud = 1;
1462 } else if (strcmp("-uv", argv[i]) == 0) {
1463 wrtuv = 1;
1464 } else if (strcmp("-maxa", argv[i]) == 0) {
1465 i++;
1466 sscanf(argv[i], "%lf", &mxang);
1467 } else if (strcmp("-maxe", argv[i]) == 0) {
1468 i++;
1469 sscanf(argv[i], "%lf", &mxedg);
1470 } else if (strcmp("-maxc", argv[i]) == 0) {
1471 i++;
1472 sscanf(argv[i], "%lf", &mxchd);
1473 } else if (strcmp("-ggf", argv[i]) == 0) {
1474 i++;
1475 sscanf(argv[i], "%lf", &ggf);
1476 } else {
1477 printf("Usage: egads2cgt [argument list] (defaults are in parenthesis)\n");
1478 printf(" -i input geometry filename <*.egads, *.stp, *.igs> ()\n");
1479 printf(" -aflr4 use AFLR4 for surface triangulation\n");
1480 printf(" -q write structured patches from quadding scheme to plot3d surface grid file\n");
1481 printf(" -uv write structured patches from uv evaluation to plot3d surface grid file\n");
1482 printf(" -ggf factr <geometric growth factor with -uv for isocline smoothing> (1.2)\n");
1483 printf(" -maxa mxang <Max allow dihedral angle (deg)> (15.0)\n");
1484 printf(" -maxe mxedg <Max allow edge length> (0.025 * size)\n");
1485 printf(" -maxc mxchd <Max allow chord-height tolerance> (0.001 * size)\n");
1486 return 1;
1487 }
1488 }
1489
1490 return 0;
1491}
1492
1493
1494int main(int argc, char *argv[])
1495
1496{
1497 int i, j, k, m, mm, n, status, oclass, mtype, nbody, nvert, ntriang;
1498 int nnode, nedge, nloop, nface, dim, ptype, pindex;
1499 int qnverts, npatch, npatchtot, nn, ip, ii, jj, iroot, iper, count;
1500 const int *qptype, *qpindex, *pvindex, *pvbounds;
1501 int *triang, *comp, *nipatch, *njpatch, *senses;
1502 char ifile[120], filename[120], rootname[120], trifilename[120];
1503 char p3dfilename[120], uvfilename[120];
1504 const char *OCCrev;
1505 double params[3], box[6], size, *verts, range[4], uv[2], result[18];
1506 double qparams[3], alen, tol, *ndist;
1507 const double *qxyz, *quv;
1508 FILE *fp;
1509 ego context, model, geom, *bodies, tess, body;
1510 ego *dum, *nodes, *edges, *loops, *faces, *tesses;
1511 UVmap *map;
1512
1513 /* look at EGADS revision */
1514 EG_revision(&i, &j, &OCCrev);
1515 printf("\n Using EGADS %2d.%02d %s\n\n", i, j, OCCrev);
1516
1517 /* get arguments */
1518 status = parse_args(argc, argv, ifile);
1519 if (status != 0) exit(status);
1520
1521 printf("mxang: %lf\n", mxang);
1522 printf("mxedg: %lf\n", mxedg);
1523 printf("mxchd: %lf\n", mxchd);
1524 if (wrtuv == 1) printf("ggf: %lf\n", ggf);
1525
1526 /* Build filenames from rootname */
1527 iroot = strroot(ifile);
1528 printf("iroot: %d\n", iroot);
1529 strncpy(rootname, ifile, iroot);
1530 rootname[iroot] = 0;
1531 printf("rootname: %s\n", rootname);
1532
1533 /* Initialize */
1534 status = EG_open(&context);
1535 if (status != EGADS_SUCCESS0) {
1536 printf(" EG_open = %d!\n\n", status);
1537 return 1;
1538 }
1539 status = EG_loadModel(context, 0, ifile, &model);
1540 if (status != EGADS_SUCCESS0) {
1541 printf(" EG_loadModel = %d\n\n", status);
1542 EG_close(context);
1543 return 1;
1544 }
1545 status = EG_getBoundingBox(model, box);
1546 if (status != EGADS_SUCCESS0) {
1547 printf(" EG_getBoundingBox = %d\n\n", status);
1548 EG_deleteObject(model);
1549 EG_close(context);
1550 return 1;
1551 }
1552
1553 size = sqrt((box[0]-box[3])*(box[0]-box[3]) + (box[1]-box[4])*(box[1]-box[4]) +
1554 (box[2]-box[5])*(box[2]-box[5]));
1555
1556 /* Get all Bodies */
1557 status = EG_getTopology(model, &geom, &oclass, &mtype, NULL((void*)0), &nbody, &bodies,
1558 &triang);
1559 if (status != EGADS_SUCCESS0) {
1560 printf(" EG_getTopology = %d\n\n", status);
1561 EG_deleteObject(model);
1562 EG_close(context);
1563 return 1;
1564 }
1565 printf(" Number of Bodies = %d\n\n", nbody);
1566
1567 /* Set tesselation parameters */
1568 if (mxedg == 0.0) {
1569 params[0] = 0.025*size;
1570 } else if (mxedg > 0.0) {
1571 params[0] = mxedg;
1572 } else if (mxedg < 0.0) {
1573 params[0] = fabs(mxedg*size);
1574 }
1575
1576 if (mxchd == 0.0) {
1577 params[1] = 0.001*size;
1578 } else if (mxchd > 0.0) {
1579 params[1] = mxchd;
1580 } else if (mxchd < 0.0) {
1581 params[1] = fabs(mxchd*size);
1582 }
1583 if (mxang == 0.0) {
1584 params[2] = 15.0;
1585 } else {
1586 params[2] = fabs(mxang);
1587 }
1588 printf(" Tess params: %lf %lf %lf\n\n", params[0], params[1], params[2]);
1589 tesses = (ego *) EG_alloc(nbody*sizeof(ego));
1590 if (tesses == NULL((void*)0)) {
1591 printf(" Alloc error on %d EGOs\n\n", nbody);
1592 EG_deleteObject(model);
1593 EG_close(context);
1594 return 1;
1595 }
1596 for (i = 0; i < nbody; i++) tesses[i] = NULL((void*)0);
1597
1598 /* use AFLR4? */
1599 if (aflr4 == 1) {
1600 status = aflr4egads(model, tesses);
1601 if (status != EGADS_SUCCESS0) {
1602 printf(" aflr4egads = %d\n\n", status);
1603 for (i = 0; i < nbody; i++)
1604 if (tesses[i] != NULL((void*)0)) EG_deleteObject(tesses[i]);
1605 EG_free(tesses);
1606 EG_deleteObject(model);
1607 EG_close(context);
1608 return 1;
1609 }
1610 printf(" Tessellations completed!\n\n");
1611 }
1612
1613 /* ---------------- Tesselate each body ---------------- */
1614
1615 for (i = 0; i < nbody; i++) {
1616 snprintf(filename, 120, "%s.%3.3d.p3d", rootname, i+1);
1617
1618 /* Add attribute */
1619 qparams[0] = 0.15;
1620 qparams[1] = 10.0;
1621 qparams[2] = 0.0;
1622 status = EG_attributeAdd(bodies[i], ".qParams", ATTRREAL2, 3, NULL((void*)0), qparams,
1623 NULL((void*)0));
1624 if (status != EGADS_SUCCESS0) {
1625 printf(" Body %d: EG_attributeAdd = %d\n", i+1, status);
1626 }
1627
1628 mtype = 0;
1629 status = EG_getTopology(bodies[i], &geom, &oclass, &mtype, NULL((void*)0), &j,
1630 &dum, &triang);
1631 if (status != EGADS_SUCCESS0) {
1632 printf(" Body %d: EG_getTopology = %d\n", i+1, status);
1633 continue;
1634 }
1635
1636 status = EG_tolerance(bodies[i], &tol);
1637 if (status != EGADS_SUCCESS0) {
1638 printf(" Body %d: EG_tolerance = %d\n", i+1, status);
1639 continue;
1640 }
1641 if (mtype == WIREBODY6) {
1642 printf(" Body %2d: Type = WireBody -- tolerance = %le\n", i+1, tol);
1643 } else if (mtype == FACEBODY7) {
1644 printf(" Body %2d: Type = FaceBody -- tolerance = %le\n", i+1, tol);
1645 } else if (mtype == SHEETBODY8) {
1646 printf(" Body %2d: Type = SheetBody -- tolerance = %le\n", i+1, tol);
1647 } else {
1648 printf(" Body %2d: Type = SolidBody -- tolerance = %le\n", i+1, tol);
1649 }
1650
1651 if (aflr4 == 0) {
1652 status = EG_makeTessBody(bodies[i], params, &tess);
1653 if (status != EGADS_SUCCESS0) {
1654 printf(" EG_makeTessBody %d = %d\n", i, status);
1655 continue;
1656 }
1657 printf(" Tessellation completed!\n");
1658 tesses[i] = tess;
1659 } else {
1660 tess = tesses[i];
1661 }
1662
1663 status = EG_getBodyTopos(bodies[i], NULL((void*)0), NODE20, &nnode, &nodes);
1664 if (status != EGADS_SUCCESS0) {
1665 printf(" EG_getBodyTopos N %d = %d\n", i, status);
1666 continue;
1667 }
1668 status = EG_getBodyTopos(bodies[i], NULL((void*)0), EDGE21, &nedge, &edges);
1669 if (status != EGADS_SUCCESS0) {
1670 printf(" EG_getBodyTopos E %d = %d\n", i, status);
1671 EG_free(nodes);
1672 continue;
1673 }
1674 status = EG_getBodyTopos(bodies[i], NULL((void*)0), LOOP22, &nloop, &loops);
1675 if (status != EGADS_SUCCESS0) {
1676 printf(" EG_getBodyTopos L %d = %d\n", i, status);
1677 EG_free(edges);
1678 EG_free(nodes);
1679 continue;
1680 }
1681 status = EG_getBodyTopos(bodies[i], NULL((void*)0), FACE23, &nface, &faces);
1682 if (status != EGADS_SUCCESS0) {
1683 printf(" EG_getBodyTopos F %d = %d\n", i, status);
1684 EG_free(loops);
1685 EG_free(edges);
1686 EG_free(nodes);
1687 continue;
1688 }
1689
1690 /* ------------------- Write out tess owner file ----------------------- */
1691
1692 status = EG_statusTessBody(tess, &body, &n, &nvert);
1693 if (status != EGADS_SUCCESS0) {
1694 printf(" EG_statusTessBody = %d!\n", status);
1695 EG_free(faces);
1696 EG_free(loops);
1697 EG_free(edges);
1698 EG_free(nodes);
1699 continue;
1700 }
1701 comp = (int *) EG_alloc(nnode*sizeof(int));
1702 if (comp == NULL((void*)0)) {
1703 printf(" Cannot Allocate %d ints!\n", nnode);
1704 EG_free(faces);
1705 EG_free(loops);
1706 EG_free(edges);
1707 EG_free(nodes);
1708 continue;
1709 }
1710 for (n = 1; n <= nvert; n++) {
1711 status = EG_getGlobal(tess, n, &ptype, &pindex, NULL((void*)0));
1712 if (status != EGADS_SUCCESS0) {
1713 printf(" Vert %d: EG_getGlobal = %d\n", n, status);
1714 continue;
1715 }
1716 if (ptype != 0) continue;
1717 comp[pindex-1] = n;
1718 }
1719
1720 /* write it out */
1721 snprintf(trifilename, 120, "%s.%3.3d.tess", rootname, i+1);
1722 fp = fopen(trifilename, "w");
1723 if (fp == NULL((void*)0)) {
1724 printf(" Cannot Open file %s -- NO FILE WRITTEN\n", trifilename);
1725 EG_free(comp);
1726 EG_free(faces);
1727 EG_free(loops);
1728 EG_free(edges);
1729 EG_free(nodes);
1730 continue;
1731 }
1732 printf(" Writing EGADS tess file: %s\n", trifilename);
1733
1734 /* header */
1735 fprintf(fp, " %6d %6d %6d %6d %6d\n", mtype, nnode, nedge, nloop, nface);
1736 writeAttr(fp, bodies[i], NULL((void*)0));
1737
1738 /* Nodes */
1739 for (n = 1; n <= nnode; n++) {
1740 fprintf(fp, " %6d %6d\n", n, comp[n-1]);
1741 writeAttr(fp, nodes[n-1], NULL((void*)0));
1742 }
1743 EG_free(comp);
1744 EG_free(nodes);
1745
1746 /* Edges */
1747 for (n = 1; n <= nedge; n++) {
1748 status = EG_getTopology(edges[n-1], &geom, &oclass, &mtype, range, &j,
1749 &dum, &triang);
1750 if (status != EGADS_SUCCESS0) {
1751 printf(" %d: EG_getTopology = %d!\n", n, status);
1752 fprintf(fp, " %6d %6d\n", n, 0);
1753 fprintf(fp, " %6d\n", 0);
1754 continue;
1755 }
1756 if (mtype == DEGENERATE5) {
1757 fprintf(fp, " %6d %6d\n", n, 0);
1758 writeAttr(fp, edges[n-1], NULL((void*)0));
1759 continue;
1760 }
1761
1762 status = EG_getTessEdge(tess, n, &m, &qxyz, &quv);
1763 if (status != EGADS_SUCCESS0) {
1764 printf(" %d: EG_getTessEdge = %d!\n", n, status);
1765 fprintf(fp, " %6d %6d\n", n, 0);
1766 writeAttr(fp, edges[n-1], NULL((void*)0));
1767 continue;
1768 }
1769
1770 status = EG_arcLength(edges[n-1], range[0], range[1], &alen);
1771 if (status != EGADS_SUCCESS0) {
1772 printf(" %d: EG_arcLength = %d!\n", n, status);
1773 } else {
1774 if (alen < 0.01)
1775 printf(" Edge %2d: arc length = %le, number of pts = %d\n",
1776 n, alen, m);
1777 }
1778
1779 fprintf(fp, " %6d %6d\n", n, m);
1780 for (j = 1; j <= m; j++) {
1781 dim = 0;
1782 if (j == 1) {
1783 status = EG_indexBodyTopo(bodies[i], dum[0]);
1784 if (status > EGADS_SUCCESS0) dim = status;
1785 }
1786 if (j == m) {
1787 if (mtype == TWONODE2) {
1788 status = EG_indexBodyTopo(bodies[i], dum[1]);
1789 } else {
1790 status = EG_indexBodyTopo(bodies[i], dum[0]);
1791 }
1792 if (status > EGADS_SUCCESS0) dim = status;
1793 }
1794 status = EG_localToGlobal(tess, -n, j, &k);
1795 if (status != EGADS_SUCCESS0)
1796 printf(" %d/%d: EG_localToGlobal Edge = %d!\n", n, j, status);
1797 fprintf(fp, " %6d %20.13le %6d ", k, quv[j-1], dim);
1798 if (j%2 == 0) fprintf(fp,"\n");
1799 }
1800 if ((j-1)%2 != 0) fprintf(fp,"\n");
1801 writeAttr(fp, edges[n-1], NULL((void*)0));
1802 }
1803
1804 /* Loops */
1805 for (j = 1; j <= nloop; j++) {
1806 status = EG_getTopology(loops[j-1], &geom, &oclass, &mtype, NULL((void*)0), &m,
1807 &dum, &senses);
1808 if (status != EGADS_SUCCESS0) {
1809 printf(" Body %d: EG_getTopology L %d = %d\n", i+1, j, status);
1810 continue;
1811 }
1812 fprintf(fp, " %6d %6d\n", j, m);
1813 for (k = 0; k < m; k++)
1814 fprintf(fp, " %6d %6d\n",
1815 EG_indexBodyTopo(bodies[i], dum[k]), senses[k]);
1816 writeAttr(fp, loops[j-1], NULL((void*)0));
1817 }
1818 EG_free(loops);
1819
1820 /* Faces */
1821 for (j = 1; j <= nface; j++) {
1822 status = EG_getTopology(faces[j-1], &geom, &oclass, &mtype, NULL((void*)0), &mm,
1823 &dum, &senses);
1824 if (status != EGADS_SUCCESS0) {
1825 printf(" Body %d: EG_getTopology F %d = %d\n", i+1, j, status);
1826 continue;
1827 }
1828 printf(" Face %2d: surface type = %d\n", j, geom->mtype);
1829 status = EG_getTessFace(tess, j, &m, &qxyz, &quv, &qptype, &qpindex,
1830 &n, &pvindex, &pvbounds);
1831 if (status != EGADS_SUCCESS0) {
1832 printf(" %d: EG_getTessFace = %d!\n", j, status);
1833 continue;
1834 }
1835 fprintf(fp, " %6d %6d %6d %6d\n", j, mm, mtype, m);
1836 for (k = 0; k < mm; k++)
1837 fprintf(fp, " %6d", EG_indexBodyTopo(bodies[i], dum[k]));
1838 fprintf(fp, "\n");
1839 for (k = 1; k <= m; k++) {
1840 status = EG_localToGlobal(tess, j, k, &n);
1841 if (status != EGADS_SUCCESS0)
1842 printf(" %d/%d: EG_localToGlobal Face = %d!\n", j, k, status);
1843 fprintf(fp, " %6d %20.13le %20.13le %6d %6d\n",
1844 n, quv[2*k-2], quv[2*k-1], qptype[k-1], qpindex[k-1]);
1845 }
1846 writeAttr(fp, faces[j-1], NULL((void*)0));
1847 }
1848 fclose(fp);
1849
1850 /* ---------- Write quads file (plot3d format) for each body ----------- */
1851
1852 if (wrtqud == 1) {
1853
1854 qparams[0] = qparams[1] = qparams[2] = 0.0;
1855
1856 /* First count total number of quad patches to write and write to file */
1857 npatchtot = 0;
1858 EG_setOutLevel(context, 0);
1859 for (n = 0; n < nface; n++) {
1860 status = EG_makeQuads(tess, qparams, n+1);
1861 if (EGADS_SUCCESS0 == status) {
1862 status = EG_getQuads(tess, n+1, &qnverts, &qxyz, &quv, &qptype,
1863 &qpindex, &npatch);
1864 npatchtot = npatchtot + npatch;
1865 }
1866 }
1867 EG_setOutLevel(context, 1);
1868 printf(" Total number of quad patches %d \n", npatchtot);
1869
1870 if (npatchtot > 0) {
1871 snprintf(p3dfilename, 120, "%s.%3.3d.p3d", rootname, i+1);
1872 printf(" Writing PLOT3D Quadded file: %s\n", p3dfilename);
1873
1874 nipatch = njpatch = NULL((void*)0);
1875 fp = fopen(p3dfilename, "w");
1876 if (fp == NULL((void*)0)) {
1877 printf(" Error Opening %s!\n", p3dfilename);
1878 goto bailQuad;
1879 }
1880 fprintf(fp, "%d\n", npatchtot);
1881
1882 /* Now get dimensions of all patches and write to file */
1883 nipatch = (int *) EG_alloc(npatchtot*sizeof(int));
1884 njpatch = (int *) EG_alloc(npatchtot*sizeof(int));
1885 if ((nipatch == NULL((void*)0)) || (njpatch == NULL((void*)0))) {
1886 printf(" Malloc ERROR on %d Patches!\n", npatchtot);
1887 goto bailQuad;
1888 }
1889
1890 EG_setOutLevel(context, 0);
1891 for (nn = n = 0; n < nface; n++) {
1892 status = EG_makeQuads(tess, qparams, n+1) ;
1893 if (EGADS_SUCCESS0 == status) {
1894 status = EG_getQuads(tess, n+1, &qnverts, &qxyz, &quv, &qptype,
1895 &qpindex, &npatch);
1896 for (ip = 0; ip < npatch; ip++, nn++) {
1897 status = EG_getPatch(tess, n+1, ip+1, &nipatch[nn], &njpatch[nn],
1898 &pvindex, &pvbounds);
1899 }
1900 }
1901 }
1902
1903 for (ip = 0; ip < npatchtot; ip++)
1904 fprintf(fp, "%d %d %d\n", nipatch[ip], njpatch[ip],1);
1905
1906 /* Write x,y,z of each quad patch */
1907 for (nn = n = 0; n < nface; n++) {
1908
1909 status = EG_makeQuads(tess, qparams, n+1);
1910 if (status == EGADS_SUCCESS0) {
1911 status = EG_getQuads(tess, n+1, &qnverts, &qxyz, &quv, &qptype,
1912 &qpindex, &npatch);
1913 for (ip = 0; ip < npatch; ip++, nn++) {
1914 status = EG_getPatch(tess, n+1, ip+1, &nipatch[nn], &njpatch[nn],
1915 &pvindex, &pvbounds);
1916 for (count = jj = 0; jj < njpatch[nn]; jj++) {
1917 for (ii = 0; ii < nipatch[nn]; ii++) {
1918 m = pvindex[jj*nipatch[nn] + ii] - 1;
1919 fprintf(fp, "%20.13le ", qxyz[3*m]);
1920 count++;
1921 if (count%5 == 0) fprintf(fp, "\n");
1922 }
1923 }
1924 fprintf(fp, "\n");
1925
1926 for (count = jj = 0; jj < njpatch[nn]; jj++) {
1927 for (ii = 0; ii < nipatch[nn]; ii++) {
1928 m = pvindex[jj*nipatch[nn] + ii] - 1;
1929 fprintf(fp, "%20.13le ", qxyz[3*m+1]);
1930 count++;
1931 if (count%5 == 0) fprintf(fp, "\n");
1932 }
1933 }
1934 fprintf(fp, "\n");
1935
1936 for (count = jj = 0; jj < njpatch[nn]; jj++) {
1937 for (ii = 0; ii < nipatch[nn]; ii++) {
1938 m = pvindex[jj*nipatch[nn] + ii] - 1;
1939 fprintf(fp, "%20.13le ", qxyz[3*m+2]);
1940 count++;
1941 if (count%5 == 0) fprintf(fp, "\n");
1942 }
1943 }
1944 fprintf(fp, "\n");
1945
1946 }
1947 }
1948 }
1949 EG_setOutLevel(context, 1);
1950 bailQuad:
1951 if (nipatch != NULL((void*)0)) EG_free(nipatch);
1952 if (njpatch != NULL((void*)0)) EG_free(njpatch);
1953 if (fp != NULL((void*)0)) fclose(fp);
1954 }
1955 }
1956
1957 /* ----------- Write uv file (plot3d format) for each Body ------------ */
1958
1959 if (wrtuv == 1) {
1960
1961 ndist = (double *) EG_alloc(nnode*sizeof(double));
1962 map = (UVmap *) EG_alloc(nface*sizeof(UVmap));
1963 if ((map == NULL((void*)0)) || (ndist == NULL((void*)0))) {
1964 printf(" ERROR Allocating %d Maps!\n", nface);
1965 if (map != NULL((void*)0)) EG_free(map);
1966 if (ndist != NULL((void*)0)) EG_free(ndist);
1967 } else {
1968 for (n = 0; n < nnode; n++) ndist[n] = 0.0;
1969 getNodeSpacing(tess, nedge, edges, ndist);
1970 for (n = 0; n < nnode; n++)
1971 printf(" Node %2d: spacing = %lf\n", n+1, ndist[n]);
1972 snprintf(uvfilename, 120, "%s.%3.3d.uv", rootname, i+1);
1973 printf(" Writing PLOT3D Untrimmed file: %s\n", uvfilename);
1974
1975 fp = fopen(uvfilename, "w");
1976 if (fp == NULL((void*)0)) {
1977 printf(" Error Opening %s!\n", uvfilename);
1978 } else {
1979
1980 /* write out the uv file a Face at a time */
1981 for (n = 0; n < nface; n++) {
1982 status = EG_getRange(faces[n], range, &iper);
1983 if (status != EGADS_SUCCESS0)
1984 printf(" %d: EG_getRange = %d\n", n+1, status);
1985 /* start out 9x9 regardless & squeeze in a little */
1986 map[n].nu = 9;
1987 map[n].nv = 9;
1988 for (jj = 0; jj < map[n].nv; jj++)
1989 map[n].vs[jj] = range[2] +
1990 (range[3] - range[2])*(double)jj/(map[n].nv-1);
1991 map[n].vs[0] += 1.e-5*(range[3] - range[2]);
1992 map[n].vs[8] -= 1.e-5*(range[3] - range[2]);
1993 for (ii = 0; ii < map[n].nu; ii++)
1994 map[n].us[ii] = range[0] +
1995 (range[1] - range[0])*(double)ii/(map[n].nu-1);
1996 map[n].us[0] += 1.e-5*(range[1] - range[0]);
1997 map[n].us[8] -= 1.e-5*(range[1] - range[0]);
1998 /* enhance based on usual tessellation parameters */
1999 updateMap(faces[n], params[0], params[1], params[2], &map[n]);
2000 /* insert Node spacings into both U and V */
2001 insertNodeSpacing(tess, faces[n], n+1, ndist, &map[n]);
2002 /* smooth the Us and Vs based on the geometric growth factor */
2003 smoothMap(faces[n], ggf, &map[n]);
2004 }
2005 fprintf(fp, "%d\n", nface);
2006 for (n = 0; n < nface; n++)
2007 fprintf(fp,"%d %d %d \n", map[n].nu, map[n].nv, 1);
2008
2009 for (n = 0; n < nface; n++) {
2010
2011 status = EG_getTopology(faces[n], &geom, &oclass, &mtype, NULL((void*)0), &mm,
2012 &dum, &senses);
2013 if (status != EGADS_SUCCESS0) {
2014 printf(" Body %d: EG_getTopology F %d = %d\n", i+1, n+1, status);
2015 continue;
2016 }
2017
2018 for (dim = 0; dim < 3; dim++) {
2019 for (count = jj = 0; jj < map[n].nv; jj++) {
2020 uv[1] = map[n].vs[jj];
2021 for (ii = 0; ii < map[n].nu; ii++) {
2022 if ( mtype == SFORWARD1 ) {
2023 uv[0] = map[n].us[ii];
2024 } else {
2025 uv[0] = map[n].us[map[n].nu-ii-1];
2026 }
2027 status = EG_evaluate(faces[n], uv, result);
2028 if (status != EGADS_SUCCESS0)
2029 printf(" %d: EG_evaluate = %d\n", n+1, status);
2030 fprintf(fp, "%20.13le ", result[dim]);
2031 count++;
2032 if (count%5 == 0) fprintf(fp, "\n");
2033 }
2034
2035 if (count%5 != 0) fprintf(fp, "\n");
2036 }
2037 }
2038
2039 for (count = jj = 0; jj < map[n].nv; jj++) {
2040 uv[1] = map[n].vs[jj];
2041 for (ii = 0; ii < map[n].nu; ii++) {
2042 if (mtype == SFORWARD1) {
2043 uv[0] = map[n].us[ii];
2044 } else {
2045 uv[0] = map[n].us[map[n].nu-ii-1];
2046 }
2047#ifdef INFACEOCC
2048 status = EG_inFaceOCC(faces[n], 1.5*tol, uv);
2049#else
2050 status = EG_inFace(faces[n], uv);
2051#endif
2052 if (status < 0)
2053 printf(" face,ii,jj %d %d %d: inFace = %d\n",
2054 n+1, ii, jj, status);
2055 fprintf(fp, "%d ", 1-status);
2056 count++;
2057 if (count%15 == 0) fprintf(fp, "\n");
2058 }
2059 if (count%15 != 0) fprintf(fp, "\n");
2060 }
2061
2062 }
2063 fclose(fp);
2064 }
2065 EG_free(map);
2066 EG_free(ndist);
2067 }
2068 }
2069 EG_free(edges);
2070 EG_free(faces);
2071
2072 /* --------- Zip up the tessellation and write in CART3D format --------- */
2073
2074 status = bodyTessellation(tess, nface, &nvert, &verts,
2075 &ntriang, &triang, &comp);
2076 if (status != EGADS_SUCCESS0) continue;
2077
2078 /* write it out */
2079 snprintf(trifilename, 120, "%s.%3.3d.tri", rootname, i+1);
2080
2081 fp = fopen(trifilename, "w");
2082 if (fp == NULL((void*)0)) {
2083 printf(" Can not Open file %s -- NO FILE WRITTEN\n", trifilename);
2084 EG_free(verts);
2085 EG_free(triang);
2086 continue;
2087 }
2088 printf(" Writing CART3D tri file: %s\n", trifilename);
2089
2090 /* header */
2091 fprintf(fp, "%d %d\n", nvert, ntriang);
2092
2093 /* ...vertList */
2094 for (j = 0; j < nvert; j++)
2095 fprintf(fp, " %20.13le %20.13le %20.13le\n",
2096 verts[3*j ], verts[3*j+1], verts[3*j+2]);
2097
2098 /* ...Connectivity */
2099 for (j = 0; j < ntriang; j++)
2100 fprintf(fp, "%6d %6d %6d\n", triang[3*j], triang[3*j+1], triang[3*j+2]);
2101
2102 /* ...Component list*/
2103 for (j = 0; j < ntriang; j++)
2104 fprintf(fp, "%6d\n", comp[j]);
2105
2106 fclose(fp);
2107
2108 printf(" # verts = %d, # tris = %d\n\n", nvert, ntriang);
2109
2110 EG_free(verts);
2111 EG_free(triang);
2112 EG_free(comp);
2113 }
2114
2115 for (i = 0; i < nbody; i++)
2116 if (tesses[i] != NULL((void*)0)) EG_deleteObject(tesses[i]);
2117 EG_free(tesses);
2118 status = EG_deleteObject(model);
2119 if (status != EGADS_SUCCESS0) printf(" EG_deleteObject = %d\n", status);
2120
2121 EG_close(context);
2122
2123 return 0;
2124}