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: s1387.c,v 1.2 2001-03-19 15:58:48 afr Exp $
45 *
46 */
47
48
49 #define S1387
50
51 #include "sislP.h"
52
53 #if defined(SISLNEEDPROTOTYPES)
54 void
s1387(SISLSurf * ps,int ik1,int ik2,SISLSurf ** rsnew,int * jstat)55 s1387(SISLSurf *ps,int ik1,int ik2,SISLSurf **rsnew,int *jstat)
56 #else
57 void s1387(ps,ik1,ik2,rsnew,jstat)
58 SISLSurf *ps;
59 int ik1;
60 int ik2;
61 SISLSurf **rsnew;
62 int *jstat;
63 #endif
64 /*
65 *********************************************************************
66 *
67 *********************************************************************
68 *
69 * PURPOSE : To express the B-spline surface as a B-spline surface
70 * of higher orders.
71 *
72 *
73 *
74 * INPUT : ps - Surface
75 * ik1 - New order in first direction
76 * ik2 - New order in second direction
77 *
78 *
79 * OUTPUT : rsnew - The resulting order elevated surface
80 * jstat - status messages
81 * = 1 : Input orders equal
82 * to surface orders.
83 * Pointer set to input
84 * surface.
85 * = 0 : ok
86 * < 0 : error
87 *
88 *
89 * METHOD : We treat the surface as a curve in both parameter directions
90 * and utilize the curve order elevation method twice
91 *
92 *
93 * REFERENCES : Larry L. Schumaker, Spline Functions:Basic Theory.
94 * Carl de Boor, A Practial Guide to Spline.
95 *
96 *-
97 * CALLS : newCurve - Allocate space for a new curve-object.
98 * freeCurve - Free space occupied by given curve-object.
99 *
100 * WRITTEN BY : Tor Dokken, SI 1988-11
101 * REVISED BY : Johannes Kaasa, SI, May 1992 (Introduced NURBS)
102 * REVISED BY : Christophe Birkeland, SI, July 1992 (Parameters in call s1750)
103 *
104 **********************************************************************/
105 {
106 SISLCurve *qc1 = SISL_NULL; /* Temporary curve */
107 SISLCurve *qc2 = SISL_NULL; /* Temporary curve */
108 SISLCurve *qc3 = SISL_NULL; /* Temporary curve */
109 SISLCurve *qc4 = SISL_NULL; /* Temporary curve */
110 int kk1; /* Order in first parameter direction */
111 int kk2; /* Order in second parameter direction */
112 int kn1; /* NumberOrder in first parameter direction */
113 int kn2; /* Order in second parameter direction */
114 int kdim; /* Dimension used in temporary calc */
115 int kstat = 0; /* Local parameter value */
116 int kpos = 0; /* Position of error */
117 double *st1 = SISL_NULL, *st2 = SISL_NULL; /* Pointers to knot vectors */
118 double *scoef = SISL_NULL; /* Pointer to coefficients */
119 int rdim; /* Potential rational dimension. */
120 double *rcoef; /* Potential rational vertices. */
121
122 *jstat = 0;
123
124
125 kk1 = ps->ik1;
126 kk2 = ps->ik2;
127 kn1 = ps->in1;
128 kn2 = ps->in2;
129
130 if (ps->ikind == 2 || ps->ikind == 4)
131 {
132 rdim = ps->idim + 1;
133 rcoef = ps->rcoef;
134 }
135 else
136 {
137 rdim = ps->idim;
138 rcoef = ps->ecoef;
139 }
140
141 if (ik1 < kk1 || ik2 < kk2)
142 goto err183;
143
144 if (ik1 == kk1 && ik2 == kk2)
145 goto war01;
146
147
148 /* Create curve representing the surface a a curve in the second parameter
149 direction, copy input arrays. */
150
151 kdim = (ps->in1) * rdim;
152
153 qc1 = newCurve (ps->in2, ps->ik2, ps->et2, rcoef, 1, kdim, 1);
154 if (qc1 == SISL_NULL)
155 goto err171;
156
157
158 /* Make the order elevation in second parameter direction. */
159
160 s1750 (qc1, ik2, &qc2, &kstat);
161 if (kstat < 0)
162 goto error;
163
164
165 /* Remember new knot vector in second parameter direction. */
166
167 kk2 = qc2->ik;
168 kn2 = qc2->in;
169 st2 = newarray (kk2 + kn2, DOUBLE);
170 if (st2 == SISL_NULL)
171 goto err101;
172
173 memcopy (st2, qc2->et, kk2 + kn2, DOUBLE);
174
175 /* Allocate space for turned parameter directions. */
176
177 scoef = newarray (kn1 * kn2 * rdim, DOUBLE);
178 if (scoef == SISL_NULL)
179 goto err101;
180
181
182 /* Turn parameter directions. */
183
184 s6chpar (qc2->ecoef, kn1, kn2, rdim, scoef);
185
186
187 /* Represent the surface a curve using the first knot vector. */
188
189 kdim = kn2 * rdim;
190
191 qc3 = newCurve (ps->in1, ps->ik1, ps->et1, scoef, 1, kdim, 1);
192 if (qc3 == SISL_NULL)
193 goto err101;
194
195
196 /* Make the order elevation in the first parameter direction. */
197
198 s1750 (qc3, ik1, &qc4, &kstat);
199 if (kstat < 0)
200 goto error;
201
202
203 /* Remember new knot vector in first parameter direction. */
204
205 kk1 = qc4->ik;
206 kn1 = qc4->in;
207 st1 = newarray (kk1 + kn1, DOUBLE);
208 if (st1 == SISL_NULL)
209 goto err101;
210
211 memcopy (st1, qc4->et, kk1 + kn1, DOUBLE);
212
213 /* Turn parameter directions of coefficients to match surface. */
214
215
216 /* Allocate space for turned parameter directions. */
217
218 scoef = increasearray (scoef, kn1 * kn2 * rdim, DOUBLE);
219 if (scoef == SISL_NULL)
220 goto err101;
221
222 s6chpar (qc4->ecoef, kn2, kn1, rdim, scoef);
223
224
225 /* Create surface object containing the order elevated surface. */
226
227 *rsnew = newSurf (kn1, kn2, kk1, kk2, st1, st2, scoef, (ps->ikind), (ps->idim), 1);
228 if (*rsnew == SISL_NULL)
229 goto err171;
230
231 /* Set periodicity flag according to that of the input surface. */
232
233 (*rsnew)->cuopen_1 = ps->cuopen_1;
234 (*rsnew)->cuopen_2 = ps->cuopen_2;
235
236 goto out;
237
238
239 /* Input orders equal to surface orders */
240
241 war01:
242 *jstat = 1;
243 *rsnew = ps;
244 goto out;
245
246 /* Error in space allocation */
247
248 err101:
249 *jstat = -101;
250 s6err ("s1387", *jstat, kpos);
251 goto out;
252
253 /* Could not create curve or surface. */
254
255 err171:
256 *jstat = -171;
257 s6err ("s1387", *jstat, kpos);
258 goto out;
259
260 /* Order(s) specified too low */
261
262 err183:
263 *jstat = -183;
264 s6err ("s1387", *jstat, kpos);
265 goto out;
266
267 /* Error in lower level function. */
268
269 error:
270 *jstat = kstat;
271 s6err ("s1387", *jstat, kpos);
272 goto out;
273
274 /* Free local used memory. */
275
276 out:
277 if (qc1 != SISL_NULL)
278 freeCurve (qc1);
279 if (qc2 != SISL_NULL)
280 freeCurve (qc2);
281 if (qc3 != SISL_NULL)
282 freeCurve (qc3);
283 if (qc4 != SISL_NULL)
284 freeCurve (qc4);
285 if (st1 != SISL_NULL)
286 freearray (st1);
287 if (st2 != SISL_NULL)
288 freearray (st2);
289 if (scoef != SISL_NULL)
290 freearray (scoef);
291
292 return;
293 }
294