1 /* 2 * Copyright (C) 1998, 2000-2007, 2010, 2011, 2012, 2013 SINTEF ICT, 3 * Applied Mathematics, Norway. 4 * 5 * Contact information: E-mail: tor.dokken@sintef.no 6 * SINTEF ICT, Department of Applied Mathematics, 7 * P.O. Box 124 Blindern, 8 * 0314 Oslo, Norway. 9 * 10 * This file is part of SISL. 11 * 12 * SISL is free software: you can redistribute it and/or modify 13 * it under the terms of the GNU Affero General Public License as 14 * published by the Free Software Foundation, either version 3 of the 15 * License, or (at your option) any later version. 16 * 17 * SISL is distributed in the hope that it will be useful, 18 * but WITHOUT ANY WARRANTY; without even the implied warranty of 19 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 20 * GNU Affero General Public License for more details. 21 * 22 * You should have received a copy of the GNU Affero General Public 23 * License along with SISL. If not, see 24 * <http://www.gnu.org/licenses/>. 25 * 26 * In accordance with Section 7(b) of the GNU Affero General Public 27 * License, a covered work must retain the producer line in every data 28 * file that is created or manipulated using SISL. 29 * 30 * Other Usage 31 * You can be released from the requirements of the license by purchasing 32 * a commercial license. Buying such a license is mandatory as soon as you 33 * develop commercial activities involving the SISL library without 34 * disclosing the source code of your own applications. 35 * 36 * This file may be used in accordance with the terms contained in a 37 * written agreement between you and SINTEF ICT. 38 */ 39 40 #include "sisl-copyright.h" 41 42 /* 43 * 44 * $Id: sh6degen.c,v 1.2 2001-03-19 15:59:07 afr Exp $ 45 * 46 */ 47 48 49 #define SH6DEGEN 50 51 #include "sislP.h" 52 53 #if defined(SISLNEEDPROTOTYPES) 54 static void sh6degen_geom(SISLObject *,SISLObject *,double [],double [], 55 int *); 56 #else 57 static void sh6degen_geom(); 58 #endif 59 60 #if defined(SISLNEEDPROTOTYPES) 61 void sh6degen(SISLObject * po1,SISLObject * po2,SISLIntdat ** pintdat,double aepsge,int * jstat)62 sh6degen(SISLObject * po1, SISLObject * po2, SISLIntdat ** pintdat, 63 double aepsge, int *jstat) 64 #else 65 void sh6degen(po1, po2, pintdat, aepsge, jstat) 66 SISLObject *po1; 67 SISLObject *po2; 68 SISLIntdat **pintdat; 69 double aepsge; 70 int *jstat; 71 #endif 72 /* 73 ********************************************************************* 74 * 75 ********************************************************************* 76 * 77 * PURPOSE : To replace any curves which are degenerate (to a 78 * point) by points and to remove any subsequent 79 * duplicate points. 80 * 81 * 82 * INPUT : pintdat - Pointer to a pointer to intersection data. 83 * 84 * 85 * OUTPUT : jstat - status messages 86 * = 0 : OK 87 * < 0 : error 88 * 89 * 90 * METHOD : Run through the Intlists. If an Intlist is degenerate, 91 * then it is a point. If this same point occurs elsewhere 92 * in pintdat (as an Intpt or an end point of another 93 * Intlist) then delete it. Otherwise add it to the 94 * array of Intpoints. 95 * At the end, the Intpoint array may have increased 96 * and the Intlist array may have decresed. 97 * 98 * 99 * REFERENCES : 100 * 101 *- 102 * CALLS : sh6getnext - Get next point in Intlist. 103 * sh6getother - Get next point in Intlist. 104 * sh6idkpt - Remove intersection point. 105 * s6length - Lenght of vector. 106 * s6dist - Distance between points. 107 * sh6degen_geom - Fetch geometry. 108 * newIntpt - Create new Intpt structure. 109 * freeIntlist - Free an intlist structure. 110 * 111 * WRITTEN BY :Mike Floater, 20/2/92. 112 * REWISED BY: Vibeke Skytt, SI, 06.92. 113 * 114 ********************************************************************* 115 */ 116 { 117 int kstat; /* Local status variable. */ 118 int kpos = 0; /* Position of error. */ 119 int ki,kj, kl; /* Loop variables. */ 120 int kpoint; /* Number of points in list. */ 121 int kstoremark; /* Save marker of point. */ 122 int knumbpt; /* Number of points in intersection list. */ 123 int knpar; /* Number of parameter directions. */ 124 int kexact; /* Indicates which parameter direction is exact. */ 125 SISLIntpt** vpoint; /* Array of Int points. */ 126 int ipoint; /* Number of Intpoints. */ 127 int ipmax; /* Size of vpoint array. */ 128 SISLIntlist** vlist; /* Array of Intlists. */ 129 int ilist; /* Number of Intlists. */ 130 int ilmax; /* Size of vlist array. */ 131 int isdegen; /* Flag. Is curve degenerate? */ 132 double spar[4]; /* Parameter value on intersection curve of 133 the objects. */ 134 SISLIntpt *pprev, *pcurr, *pnext; /* Pointers to Intpts. */ 135 SISLIntpt *pfirst, *plast; /* Pointers to Intpts. */ 136 SISLIntpt *qpt; 137 SISLObject *qo2=SISL_NULL; /* Either po2 or dummy point. */ 138 int newilist; /* New number of Intlists. */ 139 int idim; /* Dimension of space. */ 140 double geomfirst[3]; /* geometric info --- position etc. */ 141 double geomcurr[3]; /* geometric info --- position etc. */ 142 double geommid[3]; /* geometric info --- position etc. */ 143 double len; /* Distance between two vectors. */ 144 145 knpar = po1->iobj + po2->iobj; 146 *jstat = 0; 147 148 if(*pintdat == SISL_NULL) goto out; 149 150 /* Set up local variables. */ 151 152 vpoint = (*pintdat) -> vpoint; 153 ipoint = (*pintdat) -> ipoint; 154 ipmax = (*pintdat) -> ipmax; 155 vlist = (*pintdat) -> vlist; 156 ilist = (*pintdat) -> ilist; 157 ilmax = (*pintdat) -> ilmax; 158 159 if(ipoint == 0 && ilist == 0) goto out; 160 161 /* Make dimension of geometry space. */ 162 163 if (po1->iobj == SISLPOINT) idim = po1->p1->idim; 164 else if (po1->iobj == SISLCURVE) idim = po1->c1->idim; 165 else idim = po1->s1->idim; 166 167 /* UJK, feb 93, often po1 and po2 is the same geometry ie po2 168 is a dummy, so let create a dummy point in those cases. */ 169 if ((*pintdat) -> vpoint[0]->ipar == po1->iobj) 170 { 171 if ((qo2 = newObject(SISLPOINT)) == SISL_NULL) goto err101; 172 knpar = po1->iobj; 173 } 174 else qo2 = po2; 175 176 177 for(ki=0; ki<ilist; ki++) 178 { 179 /* Find out whether the curve vlist[ki] is degenerate. 180 We assume that if all the guide points in vlist[ki] 181 are equal then it is degenerate. This should be 182 correct if SISL has done its job properly before 183 reaching this stage of the intersection algorithm. 184 Otherwise it clearly is not degenerate. */ 185 186 187 isdegen = TRUE; 188 189 pfirst = vlist[ki]->pfirst; 190 plast = vlist[ki]->plast; 191 knumbpt = vlist[ki]->inumb; 192 193 /* VSK, 02-93. Check if the list describes a trimming area. 194 In that case no reduction is to be performed. */ 195 196 if (pfirst->iinter == SI_TRIM) continue; 197 198 kstat = 0; 199 sh6degen_geom(po1,qo2,pfirst->epar,geomfirst,&kstat); 200 if(kstat < 0) goto error; 201 202 pprev = pfirst; 203 204 pcurr = sh6getnext(pprev,vlist[ki]->ind_first); 205 206 /* Loop through the guide points and compare with pfirst. */ 207 208 while(pcurr != plast) 209 { 210 kstat = 0; 211 sh6degen_geom(po1,qo2,pcurr->epar,geomcurr,&kstat); 212 if(kstat < 0) goto error; 213 214 len = s6dist(geomcurr,geomfirst,idim); 215 if(len > aepsge) 216 { 217 isdegen = FALSE; 218 break; 219 } 220 221 sh6getother(pcurr,pprev,&pnext,&kstat); 222 if(kstat < 0) goto error; 223 224 pprev = pcurr; 225 pcurr = pnext; 226 } 227 228 /* Compare plast with pfirst. */ 229 230 if(isdegen == TRUE) 231 { 232 kstat = 0; 233 sh6degen_geom(po1,qo2,plast->epar,geomcurr,&kstat); 234 if(kstat < 0) goto error; 235 236 len = s6dist(geomcurr,geomfirst,idim); 237 if(len > aepsge) 238 { 239 isdegen = FALSE; 240 } 241 } 242 243 /* Test if it is only 2 points in the intersection list. In that 244 case it may be a cyclic constant parameter curve, and it is 245 necessary to check in an internal point. */ 246 247 if (isdegen && knumbpt == 2) 248 { 249 for (kexact=0, kj=0; kj<knpar; kj++) 250 { 251 if (pfirst->curve_dir[vlist[ki]->ind_first] & 252 (1 << (kj + 1))) 253 kexact = kj + 1; 254 spar[kj] = (double)0.5*(pfirst->epar[kj]+plast->epar[kj]); 255 256 } 257 258 /* UJK, 930811: 259 kstat = (kexact == 0 || kexact > po1->iobj) ? 1 : 0; */ 260 kstat = (kexact > po1->iobj) ? 1 : 0; 261 sh6degen_geom(po1,qo2,spar,geommid,&kstat); 262 263 len = s6dist(geommid,geomfirst,idim); 264 if(len > aepsge) 265 { 266 isdegen = FALSE; 267 } 268 else 269 { 270 len = s6dist(geommid,geomcurr,idim); 271 if(len > aepsge) 272 { 273 isdegen = FALSE; 274 } 275 } 276 277 } 278 279 if(isdegen == TRUE) 280 { 281 /* Curve vlist[ki] is degenerate, i.e. it's just a point. 282 Free the list in pintdat, but make sure that there is one 283 point left. If there is a point that is an endpoint of an 284 other list, and that are equal to this degenerate list, 285 this endpoint should represent these point. 286 287 First pass through the list and mark all points as removed. */ 288 289 pcurr = pfirst; 290 kstoremark = pfirst->marker; 291 pnext = sh6getnext(pcurr,vlist[ki]->ind_first); 292 kpoint = vlist[ki]->inumb; 293 294 for (kl=0; kl<kpoint; kl++) 295 { 296 /* Mark the point as removed. */ 297 298 pcurr->marker = -99; 299 300 pprev = pcurr; 301 pcurr = pnext; 302 303 if (kl < kpoint-1) 304 { 305 sh6getother(pcurr,pprev,&pnext,&kstat); 306 if (kstat < 0) goto error; 307 } 308 } 309 310 /* Free the list. */ 311 312 freeIntlist(vlist[ki]); 313 vlist[ki] = SISL_NULL; 314 315 /* Check the intersection points. If no point represent the 316 degenerated curve,restore the first point of the curve. */ 317 318 for(kj=0; kj<ipoint; kj++) 319 { 320 qpt = vpoint[kj]; 321 if (qpt->marker == -99) continue; 322 kstat = 0; 323 sh6degen_geom(po1,qo2,qpt->epar,geomcurr,&kstat); 324 if(kstat < 0) goto error; 325 326 len = s6dist(geomfirst,geomcurr,idim); 327 if(len <= aepsge) break; 328 } 329 if (kj >= ipoint) pfirst->marker = kstoremark; 330 331 } 332 } 333 334 /* Tidy up the vlist array afterwards since some of the 335 lists have been freed. */ 336 337 newilist = 0; 338 339 for(ki=0; ki<ilist; ki++) 340 { 341 if(vlist[ki] != SISL_NULL) 342 { 343 newilist++; 344 if(newilist < ki+1) 345 { 346 vlist[newilist-1] = vlist[ki]; 347 vlist[ki] = SISL_NULL; 348 } 349 } 350 } 351 352 ilist = newilist; 353 354 /* Pass through the points and remove all marked points. */ 355 356 for (kj=0; kj<(*pintdat) -> ipoint; ) 357 { 358 qpt = (*pintdat)->vpoint[kj]; 359 if (qpt->marker == -99) 360 { 361 /* Remove point. */ 362 363 sh6idkpt (pintdat, &qpt, 0, &kstat); 364 if (kstat < 0) goto error; 365 } 366 else kj++; 367 } 368 369 /* Reset the integer pintdat variables with the local 370 ones in case they have changed. */ 371 372 (*pintdat)->ipmax = ipmax; 373 (*pintdat)->ilist = ilist; 374 (*pintdat)->ilmax = ilmax; 375 376 377 /* Degeneracies removed. Finish. */ 378 379 *jstat = 0; 380 goto out; 381 382 /* Error in alloc. */ 383 err101: 384 *jstat = -101; 385 s6err ("sh6degen", *jstat, kpos); 386 goto out; 387 388 /* Error in lower level routine. */ 389 error: 390 *jstat = kstat; 391 s6err ("sh6degen", *jstat, kpos); 392 goto out; 393 394 out: 395 if (qo2 && qo2 != po2) freeObject(qo2); 396 397 } 398 399 400 #if defined(SISLNEEDPROTOTYPES) 401 static void sh6degen_geom(SISLObject * po1,SISLObject * po2,double epar[],double geom[],int * jstat)402 sh6degen_geom(SISLObject *po1,SISLObject *po2,double epar[], 403 double geom[],int *jstat) 404 #else 405 static void sh6degen_geom(po1,po2,epar,geom,jstat) 406 SISLObject *po1; 407 SISLObject *po2; 408 double epar[]; 409 double geom[]; 410 int *jstat; 411 #endif 412 /* 413 ********************************************************************* 414 * 415 * PURPOSE : Evaluate object in an intersection point between the 416 * two objects, po1 and po2. 417 * 418 * 419 * INPUT : po1 - First object in intersection. 420 * po2 - Second object, in the case of a surface/curve - analytic 421 * intersection, po2 = po1 = object pointing at 422 * the surface or curve. 423 * epar - Parameter values of intersection point inherited 424 * from the actual intersection point. 425 * jstat - = 0 : Evaluate first object. 426 * = 1 : Evaluate second object. 427 * 428 * 429 * OUTPUT : geom - Position of intersection point in geometry space. 430 * jstat - status messages 431 * > 0 : warning 432 * = 0 : ok 433 * < 0 : error 434 * 435 * 436 ********************************************************************* 437 */ 438 { 439 int kstat = 0; /* Local status variable. */ 440 int kdim; /* Dimension. */ 441 int kder = 0; /* Evaluate position. */ 442 int kleft1 = 0; /* Parameter to the evaluator. */ 443 int kleft2 = 0; /* Parameter to the evaluator. */ 444 int keval = *jstat; /* Indicates which object to evaluate. */ 445 double sder[3]; /* Position and derivatives. */ 446 double snorm[3]; /* Dummy normal in surf. evalutation. */ 447 448 *jstat = 0; 449 if (keval == 0) 450 { 451 if (po1->iobj == SISLSURFACE) 452 { 453 /* Evaluate the surface. */ 454 455 kdim = po1->s1->idim; 456 s1421(po1->s1,kder,epar,&kleft1,&kleft2,sder,snorm,&kstat); 457 if (kstat < 0) goto error; 458 } 459 460 else if (po1->iobj == SISLCURVE) 461 { 462 463 /* Evaluate the curve. */ 464 465 kdim = po1->c1->idim; 466 s1221(po1->c1,kder,epar[0],&kleft1,sder,&kstat); 467 if (kstat < 0) goto error; 468 } 469 470 else if (po1->iobj == SISLPOINT) 471 { 472 473 /* Copy the value of the point. */ 474 475 kdim = po1->p1->idim; 476 memcopy(sder,po1->p1->ecoef,kdim,DOUBLE); 477 } 478 } 479 else if (keval == 1) 480 { 481 if (po2->iobj == SISLSURFACE) 482 { 483 /* Evaluate the surface. */ 484 485 kdim = po2->s1->idim; 486 s1421(po2->s1,kder,epar+po1->iobj,&kleft1,&kleft2,sder,snorm,&kstat); 487 if (kstat < 0) goto error; 488 } 489 490 else if (po2->iobj == SISLCURVE) 491 { 492 493 /* Evaluate the curve. */ 494 495 kdim = po2->c1->idim; 496 s1221(po2->c1,kder,epar[po1->iobj],&kleft1,sder,&kstat); 497 if (kstat < 0) goto error; 498 } 499 500 else if (po2->iobj == SISLPOINT) 501 { 502 503 /* Copy the value of the point. */ 504 505 kdim = po2->p1->idim; 506 memcopy(sder,po1->p1->ecoef,kdim,DOUBLE); 507 } 508 } 509 510 511 /* Copy result to output. */ 512 513 memcopy(geom,sder,kdim,DOUBLE); 514 515 *jstat = 0; 516 goto out; 517 518 /* Error in lower level function. */ 519 520 error : 521 *jstat = kstat; 522 goto out; 523 524 out : 525 return; 526 } 527