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: s1609.c,v 1.2 2001-03-19 15:58:51 afr Exp $ 45 * 46 */ 47 48 49 #define S1609 50 51 #include "sislP.h" 52 53 #if defined(SISLNEEDPROTOTYPES) 54 void s1609(SISLCurve * pc1,SISLCurve * pc2,double aepsge,double eps1[],double epf[],double eps2[],double aradius,double enorm[],int itype,int idim,int ik,SISLCurve ** rc,double * ct11,double * ctf1,double * ct21,double * ctf2,int * jstat)55 s1609(SISLCurve *pc1,SISLCurve *pc2,double aepsge, 56 double eps1[],double epf[],double eps2[],double aradius, 57 double enorm[],int itype,int idim,int ik,SISLCurve **rc, 58 double *ct11,double *ctf1,double *ct21,double *ctf2, 59 int *jstat) 60 #else 61 void s1609(pc1,pc2,aepsge,eps1,epf,eps2,aradius,enorm,itype,idim,ik,rc, 62 ct11,ctf1,ct21,ctf2,jstat) 63 SISLCurve *pc1; 64 SISLCurve *pc2; 65 double aepsge; 66 double eps1[]; 67 double epf[]; 68 double eps2[]; 69 double aradius; 70 double enorm[]; 71 int itype; 72 int idim; 73 int ik; 74 SISLCurve **rc; 75 double *ct11; 76 double *ctf1; 77 double *ct21; 78 double *ctf2; 79 int *jstat; 80 #endif 81 /* 82 ********************************************************************* 83 * 84 * PURPOSE : To calculate a fillet curve given radius between two curves. 85 * 86 * INPUT : pc1 - The first input curve. 87 * pc2 - The second input curve. 88 * aepsge - Geometry resolution. 89 * eps1 - SISLPoint telling that the fillet should be put on 90 * the side of curve 1 where eps1 is situated. 91 * epf - SISLPoint indicating where the fillet curve should go. 92 * eps1 together with epf indicates the direction of 93 * the start tangent of the curve, while epf together 94 * with eps2 indicates the direction of the end tangent 95 * of the curve. If more than one position of the fillet 96 * curve is possible, the closest curve to epf is chosen. 97 * eps2 - SISLPoint telling that the fillet should be put on the 98 * side of curve 2 where eps2 is situated. 99 * aradius - The radius to be used on the fillet if a circular 100 * fillet possible, else a conic or a quadratic 101 * polynomial curve is used, approximating the circular 102 * fillet. 103 * enorm - Direction normal to the plane the fillet curve 104 * should lie close to. (only used in 3-d fillet 105 * calculations). 106 * itype - Indicator of type of fillet. 107 * = 1 - Circle, interpolating tangent on first curve, 108 * not on curve 2. 109 * = 2 - Conic if possible 110 * else - Polynomial segment 111 * idim - Dimension of space. 112 * ik - Order of fillet curve. 113 * 114 * OUTPUT : rc - Fillet curve produced 115 * ct11 - Parameter value of the end of curve 1 not affected by 116 * the fillet. 117 * ctf1 - Parameter value of the point on curve 1 where the 118 * fillet starts. 119 * ct21 - Parameter value of the end of curve 2 not affected by 120 * the fillet. 121 * ctf2 - Parameter value of the point on curve 2 where the 122 * fillet ends. 123 * > 0 : warning 124 * = 0 : ok 125 * < 0 : error 126 * 127 * METHOD : 128 * 129 * USE : 130 * 131 * REFERENCES : 132 * 133 *- 134 * CALLS : 135 * 136 * 137 * WRITTEN BY : Qyvind Hjelle, SI, Oslo, Norway. 28. Nov 1988 138 * Reviced by : Tor Dokken, SI, Oslo, Norway, August 1989 139 * 140 ********************************************************************* 141 */ 142 { 143 SISLIntcurve **qic1=SISL_NULL; /* SISLObject containing intervals if any */ 144 SISLIntcurve **qic2=SISL_NULL; /* SISLObject containing intervals if any */ 145 SISLIntcurve **qic3=SISL_NULL; /* SISLObject containing intervals if any */ 146 SISLCurve *qc1=SISL_NULL; /* Offset of first curve */ 147 SISLCurve *qc2=SISL_NULL; /* Offset of second curve */ 148 149 int kstat; /* Status variable */ 150 int kpos=0; /* Position of error */ 151 int kleft=0; /* Pointer in knot vector */ 152 int kn; /* The number of B-splines vertices */ 153 int kk; /* The polynomial order of the curve. */ 154 double *st; /* Pointer to the first element of the knot vector 155 of the curve. The knot vector has [kn+kk] 156 elements. */ 157 158 int kcrv1,kcrv2; 159 int kcrv3 ; /* Number of intervals */ 160 int kpt1,kpt2,kpt3; /* Number of points */ 161 int ki; /* Loop variable */ 162 double tpar1,tpar2; /* Parameter values */ 163 double spnt1[6]; /* SISLPoint and derivate */ 164 double spnt2[6]; /* SISLPoint and derivate */ 165 double *spar1=SISL_NULL; /* Pointer to parameter values */ 166 double *spar2=SISL_NULL; /* Pointer to parameter values */ 167 double *spar3=SISL_NULL; /* Pointer to parameter values */ 168 double *spar4=SISL_NULL; /* Pointer to parameter values */ 169 double snorm[3]; /* Normal vector */ 170 double sdiff[3]; /* Difference vector */ 171 double tprod; /* Scalar product */ 172 double tdir1,tdir2; /* Direction of curves */ 173 double trad1,trad2; /* Offsets with correct direction */ 174 double tmin; /* Minimum distance */ 175 double tdist; /* Distance */ 176 177 /* Check dimensions */ 178 179 if (idim != 2 && idim != 3) goto err105; 180 if (pc1->idim != pc2->idim) goto err106; 181 182 /* Check if curves are correct */ 183 184 s1707(pc1,&kstat); 185 if (kstat < 0) goto error; 186 187 s1707(pc2,&kstat); 188 if (kstat < 0) goto error; 189 190 /* Calculate closest point to eps1 */ 191 192 s1953(pc1,eps1,idim,REL_COMP_RES,aepsge,&kpt1,&spar1,&kcrv1,&qic1,&kstat); 193 if (kstat < 0) goto error; 194 195 /* Remember closest point */ 196 197 if (kpt1 > 0) 198 tpar1 = spar1[0]; 199 else if (kcrv1>0) 200 { 201 SISLIntcurve *q1= *qic1; 202 if (q1->ipar1 ==1) 203 tpar1 = q1 -> epar1[0]; 204 else if (q1->ipar2 ==1) 205 tpar1 = q1 -> epar2[0]; 206 else 207 goto errxxx; 208 } 209 210 /* Calculate point and tangent on curve at tpar1 */ 211 212 s1221(pc1,1,tpar1,&kleft,spnt1,&kstat); 213 if (kstat < 0) goto error; 214 215 216 /* Calculate closest point to eps2 */ 217 218 s1953(pc2,eps2,idim,REL_COMP_RES,aepsge,&kpt2,&spar2,&kcrv2,&qic2,&kstat); 219 if (kstat < 0) goto error; 220 221 222 /* Remember closest point */ 223 224 if (kpt2 > 0) 225 tpar2 = spar2[0]; 226 else if (kcrv1>0) 227 { 228 SISLIntcurve *q1= *qic2; 229 if (q1->ipar1 ==1) 230 tpar2 = q1 -> epar1[0]; 231 else if (q1->ipar2 ==1) 232 tpar2 = q1 -> epar2[0]; 233 else 234 goto errxxx; 235 } 236 237 /* Calculate point and tangent on curve at tpar1 */ 238 239 s1221(pc2,1,tpar2,&kleft,spnt2,&kstat); 240 if (kstat<0) goto error; 241 242 243 /* Determine if offset to be used on curve 1 is positive or negative */ 244 245 if (idim==2) 246 { 247 snorm[0] = -spnt1[3]; 248 snorm[1] = spnt1[2]; 249 s6diff(eps1,spnt1,idim,sdiff); 250 tprod = s6scpr(snorm,sdiff,idim); 251 s6diff(epf,spnt1,idim,sdiff); 252 tdir1 = s6scpr(sdiff,spnt1+idim,idim); 253 } 254 else 255 { 256 s6diff(epf,spnt1,idim,sdiff); 257 tdir1 = s6scpr(spnt1+idim,sdiff,idim); 258 s6crss(sdiff,spnt1+idim,snorm); 259 tprod = s6scpr(snorm,enorm,idim); 260 } 261 262 /* Set direction of offset of curve 1 */ 263 264 if (tprod >= DZERO) 265 trad1 = aradius; 266 else 267 trad1 = -aradius; 268 269 /* Determine if offset to be used on curve 1 is positive or negative */ 270 271 if (idim==2) 272 { 273 snorm[0] = -spnt2[3]; 274 snorm[1] = spnt2[2]; 275 s6diff(eps2,spnt2,idim,sdiff); 276 tprod = s6scpr(snorm,sdiff,idim); 277 s6diff(epf,spnt2,idim,sdiff); 278 tdir2 = s6scpr(sdiff,spnt2+idim,idim); 279 } 280 else 281 { 282 s6diff(epf,spnt2,idim,sdiff); 283 tdir2 = s6scpr(spnt2+idim,sdiff,idim); 284 s6crss(sdiff,spnt2+idim,snorm); 285 tprod = s6scpr(snorm,enorm,idim); 286 } 287 288 289 /* Determine if offset to be used on curve 1 is positive or negative */ 290 291 if (tprod >= DZERO) 292 trad2 = aradius; 293 else 294 trad2 = -aradius; 295 296 /* Make curve offset from curve 1 */ 297 298 s1360(pc1,trad1,aepsge,enorm,(double)0.0,idim,&qc1,&kstat); 299 if (kstat<0) goto error; 300 301 302 /* Make curve offset from curve 2 */ 303 304 s1360(pc2,trad2,aepsge,enorm,(double)0.0,idim,&qc2,&kstat); 305 if (kstat<0) goto error; 306 307 /* Find closest point between the two curves, 308 if intersection does not succeed 309 use closest point between curves */ 310 311 s1857(qc1,qc2,REL_COMP_RES,aepsge,&kpt3,&spar3,&spar4,&kcrv3,&qic3,&kstat); 312 if (kstat<0) goto error; 313 314 if (kpt3 == 0 && kcrv3 == 0) 315 { 316 s1955(qc1,qc2,REL_COMP_RES,aepsge, 317 &kpt3,&spar3,&spar4,&kcrv3,&qic3,&kstat); 318 if (kstat<0) goto error; 319 } 320 321 /* Find point closest to ept */ 322 323 tmin=HUGE; 324 325 if (kpt3>0) 326 { 327 /* Find intersection point closest to epf */ 328 329 for (ki=0;ki<kpt3;ki++) 330 { 331 s1221(pc1,0,spar3[ki],&kleft,spnt1,&kstat); 332 if (kstat<0) goto error; 333 334 tdist = s6dist(spnt1,epf,idim); 335 if (tdist<=tmin) 336 { 337 *ctf1 = spar3[ki]; 338 *ctf2 = spar4[ki]; 339 tmin = tdist; 340 } 341 } 342 } 343 else if (kcrv3>0) 344 { 345 SISLIntcurve *q1= *qic3; 346 for (ki=0; ki < q1->ipoint ; ki++) 347 { 348 s1221(pc1,0,q1->epar1[ki],&kleft,spnt1,&kstat); 349 if (kstat<0) goto error; 350 tdist = s6dist(spnt1,epf,idim); 351 if (tdist<=tmin) 352 { 353 *ctf1 = q1->epar1[ki]; 354 *ctf2 = q1->epar2[ki]; 355 tmin = tdist; 356 } 357 } 358 } 359 360 if (tdist == HUGE) goto errxxx; 361 362 /* Set parameter values indicating which part of curves remains after 363 the fillet */ 364 365 st = pc1->et; 366 kk = pc1->ik; 367 kn = pc1->in; 368 369 if (tdir1>=0) 370 *ct11 = st[kk-1]; 371 else 372 *ct11 = st[kn]; 373 374 /* Set parameter values indicating which part of curve two remains after 375 the fillet */ 376 377 st = pc2->et; 378 kk = pc2->ik; 379 kn = pc2->in; 380 381 if (tdir2>=0) 382 *ct21 = st[kk-1]; 383 else 384 *ct21 = st[kn]; 385 386 s1607(pc1,pc2,aepsge,*ct11,*ctf1,*ct21,*ctf2,itype,idim,ik,rc,&kstat); 387 if (kstat<0) goto error; 388 389 *jstat = 0; 390 391 goto out; 392 393 394 /* Error in input, conflicting dimensions */ 395 396 err106: 397 *jstat = -106; 398 s6err("s1609",*jstat,kpos); 399 goto out; 400 401 /* Dimension nmot equal to 2 or 3 */ 402 403 err105: 404 *jstat = -105; 405 s6err("s1609",*jstat,kpos); 406 goto out; 407 408 409 /* No fillet produced */ 410 411 errxxx: 412 *jstat = -1; 413 goto out; 414 415 /* Error in lower level function */ 416 417 error: 418 *jstat = kstat; 419 s6err("s1609",*jstat,kpos); 420 goto out; 421 422 out: 423 if (qc1 != SISL_NULL) freeCurve(qc1); 424 if (qc2 != SISL_NULL) freeCurve(qc2); 425 if (qic1 != SISL_NULL) freeIntcrvlist(qic1,kcrv1); 426 if (qic2 != SISL_NULL) freeIntcrvlist(qic2,kcrv2); 427 if (qic3 != SISL_NULL) freeIntcrvlist(qic3,kcrv3); 428 if (spar1 != SISL_NULL) freearray(spar1); 429 if (spar2 != SISL_NULL) freearray(spar2); 430 if (spar3 != SISL_NULL) freearray(spar3); 431 if (spar4 != SISL_NULL) freearray(spar4); 432 433 return; 434 } 435