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: s9clipit.c,v 1.1 1994-04-21 12:10:42 boh Exp $
45 *
46 */
47
48
49 #define S9CLIPIT
50
51 #include "sislP.h"
52
53 #if defined(SISLNEEDPROTOTYPES)
54 void
s9clipit(double epar11[],double epar12[],double epar21[],double epar22[],SISLSurf * psurf1,SISLSurf * psurf2,double euval[],double evval[],double esval[],double etval[],double aepsge,double gpnt1[],double gpnt2[],double gpar1[],double gpar2[],int * jstat)55 s9clipit(double epar11[],double epar12[],double epar21[],double epar22[],
56 SISLSurf *psurf1,SISLSurf *psurf2,double euval[],double evval[],double esval[],
57 double etval[],double aepsge,double gpnt1[],double gpnt2[],
58 double gpar1[],double gpar2[],int *jstat)
59 #else
60 void s9clipit(epar11,epar12,epar21,epar22,psurf1,psurf2,
61 euval,evval,esval,etval,
62 aepsge,gpnt1,gpnt2,gpar1,gpar2,jstat)
63 double epar11[];
64 double epar12[];
65 double epar21[];
66 double epar22[];
67 SISLSurf *psurf1;
68 SISLSurf *psurf2;
69 double euval[];
70 double evval[];
71 double esval[];
72 double etval[];
73 double aepsge;
74 double gpnt1[];
75 double gpnt2[];
76 double gpar1[];
77 double gpar2[];
78 int *jstat;
79 #endif
80 /*
81 *********************************************************************
82 *
83 * PURPOSE : To clip the intersection curve between epar1 and epar2
84 * against the parameter boundary of the patch 1 defined
85 * by euval and evval, and the parameter boundary of the
86 * patch 2 defined by esval and etval.
87 *
88 *
89 * INPUT : epar11 - Parameter pair of start point in first surface
90 * epar12 - Parameter pair of start point in second surface
91 * epar21 - Parameter pair of end point in first surface
92 * epar22 - Parameter pair of end point in second surface
93 * psurf1 - Description of B-spline surface 1
94 * psurf2 - Description of B-spline surface 2
95 * euval - Parameter interval in first parameter direction
96 * evval - Parameter interval in second parameter direction
97 * esval - Parameter interval in third parameter direction
98 * etval - Parameter interval in fourth parameter direction
99 * aepsge - Absolute tolerance
100 *
101 *
102 * OUTPUT : gpnt1 - 0-2 Derivatives + normal of result of iteration
103 * in B-spline surface 1
104 * gpnt2 - 0-2 Derivatives + normal of result of iteration
105 * in B-spline surface 2
106 * gpar1 - Parameter pair of result of iteration in B-spline
107 * surface 1
108 * gpar2 - Parameter pair of result of iteration in B-spline
109 * surface 2
110 * jstat - status messages
111 * = 2 : Iteration diverged or to many iterations
112 * = 1 : ok, The line cross the boundary, point
113 * found
114 * = 0 : ok, The line does not cross the boundary
115 * < 0 : error
116 *
117 *
118 * METHOD :
119 * USE : The function is only working i 3-D
120 *
121 * REFERENCES :
122 *
123 *-
124 * CALLS :
125 *
126 * WRITTEN BY : Tor Dokken, SI, Oslo, Norway, 4-July-1988
127 * Revised by : Tor Dokken, SI, Oslo, Norway, 4-APril-1989
128 * Correction of steps over two boundaries
129 *
130 *********************************************************************
131 */
132 {
133 int kpos=0; /* Position of error */
134 int klfu=0; /* Variable used as pointer in knot vector */
135 int klfv=0; /* Variable used as pointer in knot vector */
136 int klfs=0; /* Variable used as pointer in knot vector */
137 int klft=0; /* Variable used as pointer in knot vector */
138 int kder=2; /* Number of derivatives to be calculated */
139 int kstat; /* Local status variable */
140 int kdir; /* Parameter direction of tpar */
141 int kcross; /* Control variable in while loop */
142 int knbit; /* Number of iterations */
143 int krem; /* Remember status */
144 int kbound; /* Numbering of boundary */
145 double tpar; /* Constant parameter value */
146 double spnt1[21]; /* Coordinates of point */
147 double spnt2[21]; /* Coordinates of point */
148 double spar11[2],spar12[2]; /* Local parameter values */
149 double spar21[2],spar22[2]; /* Local parameter values */
150 double spar31[2],spar32[2]; /* Local parameter values */
151
152 /* Make local copy of parameters */
153
154 memcopy(spar11,epar11,2,DOUBLE);
155 memcopy(spar12,epar12,2,DOUBLE);
156 memcopy(spar21,epar21,2,DOUBLE);
157 memcopy(spar22,epar22,2,DOUBLE);
158
159 kcross = 1;
160 knbit = 0;
161
162 while(kcross && knbit<8)
163 {
164
165 /* Find intersection between boundary of parameter area and patch */
166
167 s1330(spar11,spar12,spar21,spar22,euval,evval,esval,etval,
168 &kbound,spar31,spar32,&kstat);
169 if (kstat<0) goto error;
170
171 /* Remember status so that the line can be updated properly */
172
173 krem = kstat;
174
175 if (kstat<2 || kbound == 0)
176 {
177 /* No intersection */
178
179 kcross = 0;
180 }
181 else
182 {
183
184 /* Calculate start point for iteration */
185
186 s1421(psurf1,kder,spar31,&klfu,&klfv,spnt1,spnt1+18,&kstat);
187 if (kstat<0) goto error;
188
189 s1421(psurf2,kder,spar32,&klfs,&klft,spnt2,spnt2+18,&kstat);
190 if (kstat<0) goto error;
191
192 if (kbound==1)
193 {
194 kdir = 1;
195 tpar = euval[0];
196 }
197 else if (kbound==2)
198 {
199 kdir = 2;
200 tpar = evval[1];
201 }
202 else if (kbound==3)
203 {
204 kdir = 1;
205 tpar = euval[1];
206 }
207 else if (kbound==4)
208 {
209 kdir = 2;
210 tpar = evval[0];
211 }
212 else if (kbound==5)
213 {
214 kdir = 3;
215 tpar = esval[0];
216 }
217 else if (kbound==6)
218 {
219 kdir = 4;
220 tpar = etval[1];
221 }
222 else if (kbound==7)
223 {
224 kdir = 3;
225 tpar = esval[1];
226 }
227 else if (kbound==8)
228 {
229 kdir = 4;
230 tpar = etval[0];
231 }
232
233
234 /* Iterate to boundary intersection */
235
236 s9boundit(spnt1,spnt2,spar31,spar32,psurf1,psurf2,tpar,kdir,aepsge,
237 gpnt1,gpnt2,gpar1,gpar2,&kstat);
238 if (kstat<0) goto error;
239 if (kstat==2) goto war02;
240
241 /* Iteration converged, copy output if new loop necessary */
242
243 if (krem == 2)
244 {
245 /* spar1 was inside, update spar2 */
246
247 memcopy(spar21,gpar1,2,double);
248 memcopy(spar22,gpar2,2,double);
249 }
250 else
251 {
252 /* spar2 was inside, update spar1 */
253
254 memcopy(spar11,gpar1,2,double);
255 memcopy(spar12,gpar2,2,double);
256 }
257
258 /* Update number of iterations */
259
260 knbit++;
261 }
262 }
263
264
265
266 /* Problem solved if kcross==0. In this case we might have two cases:
267 - iteration not used: Then knbit=0
268 - iteration used : Then knbit>0
269
270 if kcross==1, then we have stopped on the condition knbit>7, and we
271 have no success. */
272
273 if (kcross==0 && knbit ==0)
274 {
275 /* Iteration not used because boundary not crossed */
276
277 *jstat = 0;
278 }
279 else if (kcross==0 && knbit > 0)
280 {
281 /* Boundary crossed, point found, more than one intersection point
282 is possible, check which to return */
283
284 if (spar11[0] == epar11[0] && spar11[1] == epar11[1] &&
285 spar12[0] == epar12[0] && spar12[1] == epar12[1] )
286 {
287 memcopy(gpar1,spar21,2,double);
288 memcopy(gpar2,spar22,2,double);
289 }
290 else
291 {
292 memcopy(gpar1,spar11,2,double);
293 memcopy(gpar2,spar12,2,double);
294
295 /* Calculate crossing point, only necessary when we step into
296 the patch since we already could have stepped out and this
297 point is recorded in gpnt1 */
298
299 s1421(psurf1,kder,gpar1,&klfu,&klfv,gpnt1,gpnt1+18,&kstat);
300 if (kstat<0) goto error;
301
302 s1421(psurf2,kder,gpar2,&klfu,&klfv,gpnt2,gpnt2+18,&kstat);
303 if (kstat<0) goto error;
304 }
305
306 *jstat = 1;
307 }
308 else
309 {
310 /* To many tries */
311 goto war02;
312 }
313 goto out;
314
315 /* Iteration without success */
316
317 war02:
318 *jstat = 2;
319 goto out;
320
321 /* Error in lower level routine. */
322
323 error :
324 *jstat = kstat;
325 s6err("s9clipit",*jstat,kpos);
326 goto out;
327
328 out:
329 return;
330 }
331