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: s6addcurve.c,v 1.2 2001-03-19 15:59:00 afr Exp $
45  *
46  */
47 
48 
49 #define S6ADDCURVE
50 
51 #include "sislP.h"
52 
53 #if defined(SISLNEEDPROTOTYPES)
54 static void s6addcurve_s9moveknots(double [],int,double,double,
55 				   double [],int *);
56 #else
57 static void s6addcurve_s9moveknots();
58 #endif
59 
60 #if defined(SISLNEEDPROTOTYPES)
s6addcurve(SISLCurve * pc1,SISLCurve * pc2,int isign,SISLCurve ** rcurve,int * jstat)61 void s6addcurve(SISLCurve *pc1,SISLCurve *pc2,int isign,
62 		SISLCurve **rcurve,int *jstat)
63 #else
64 void s6addcurve(pc1,pc2,isign,rcurve,jstat)
65      SISLCurve *pc1;
66      SISLCurve *pc2;
67      int isign;
68      SISLCurve **rcurve;
69      int *jstat;
70 #endif
71 /*
72 *********************************************************************
73 *
74 * PURPOSE    : Compute the sum or difference between two B-spline
75 *              curves depending on the constang isign. That is find
76 *              the curve pc1 + isign*pc2.
77 *
78 *
79 *
80 * INPUT      : pc1      - First curve in expression.
81 *              pc2      - Second curve in expression.
82 *              isign    - Sign of second curve.
83 *
84 *
85 * OUTPUT     : rcurve   - Sum/difference between the input curves.
86 *              jstat    - status messages
87 *                         > 0 : warning
88 *                         = 0 : ok
89 *                         < 0 : error
90 *
91 *
92 * METHOD     : Express the curves on the same knot vector. Compute the
93 *              sum difference between corresponding vertices.
94 *
95 *
96 * REFERENCES :
97 *
98 *
99 * USE        :
100 *
101 *-
102 * CALLS      : s1932  - Express curve-set in given basis.
103 *              s1363  - Pick parametrisation of curve.
104 *              s6takeunion - Union of two ordered vectors.
105 *              make_cv_kreg - Make curve k-regular.
106 *              s1333_count - Count continuity in ends of closed curves.
107 *              make_cv_cyclic - Represent curve in a periodic basis.
108 *              newCurve  -   Create new curve-object.
109 *              s6addcurve_s9moveknots() - Move knotvector into new parameter interval.
110 *
111 *
112 * WRITTEN BY : Vibeke Skytt, SI, 08.90.
113 * REWISED BY : Vibeke Skytt, SI, 05.92. To treat periodic curves.
114 *
115 *********************************************************************
116 */
117 {
118   int kstat = 0;         /* Status variable.  */
119   int ki;                /* Counter.          */
120   int kdim = pc1->idim;  /* Dimension of geometry space.              */
121   int knbcrv = 2;        /* Number of input curves.                   */
122   int kknot1 = pc1->in + pc1->ik;  /* Number of knots of first curve. */
123   int kknot2 = pc2->in + pc2->ik;  /* Number of knots of second curve. */
124   int korder;            /* Order of output curve.     */
125   int ktau;              /* Number of knots / number of vertices of
126 			    output curve.                             */
127   int kcont;             /* Continuity at ends in periodic case.      */
128   int kopen;             /* Open/closed/periodic parameter.           */
129   int kkind = pc1->ikind;  /* Kind of curves.                         */
130   int kcopy = 1;         /* Copy arrays when creating new curve.      */
131   double tmin1,tmax1;    /* Parameter interval of first input curve.  */
132   double tmin2,tmax2;    /* Parameter interval of second input curve. */
133   double *st = SISL_NULL;     /* Knot vector of second curve.              */
134   double *stau = SISL_NULL;   /* Knot vector of output curve.              */
135   double *scoef = SISL_NULL;  /* Vertices of input curves represented in
136 			    the same knot vector, then of output curve. */
137   SISLCurve *ucurves[2];     /* Local pointer array to input curves.      */
138 
139   ucurves[0] = ucurves[1] = SISL_NULL;
140 
141   /* Test input.  */
142 
143   if (kdim != pc2->idim) goto err106;
144 
145   /* Make curves k-regular.  */
146 
147   make_cv_kreg(pc1,ucurves,&kstat);
148   if (kstat < 0) goto error;
149 
150   make_cv_kreg(pc2,ucurves+1,&kstat);
151   if (kstat < 0) goto error;
152 
153   /* Allocate scratch for local knot vector.  */
154 
155   if ((st = newarray(kknot2,DOUBLE)) == SISL_NULL) goto err101;
156 
157   /* Fetch endparameter values of the input curves.  */
158 
159   s1363(ucurves[0],&tmin1,&tmax1,&kstat);
160   if (kstat < 0) goto error;
161 
162   s1363(ucurves[1],&tmin2,&tmax2,&kstat);
163   if (kstat < 0) goto error;
164 
165   if (DNEQUAL(tmin1,tmin2) || DNEQUAL(tmax1,tmax2))
166     {
167       /* Move the knot vector of the second curve into the parameter
168 	 interval of the first curve.  */
169 
170        s6addcurve_s9moveknots(ucurves[1]->et,kknot2,tmin1,tmax1,st,&kstat);
171       if (kstat < 0) goto error;
172     }
173   else memcopy(st,ucurves[1]->et,kknot2,DOUBLE);
174 
175   /* Find the union of the knot vectors of the two curves.  */
176 
177   s6takeunion(ucurves[0]->et,kknot1,st,kknot2,&stau,&ktau,&kstat);
178   if (kstat < 0) goto error;
179 
180   korder = MAX(ucurves[0]->ik,ucurves[1]->ik);
181   ktau -= korder;
182 
183   /* Express the curves in the basis.  */
184 
185   s1932(2,ucurves,tmin1,tmax1,stau,ktau,korder,&scoef,&kstat);
186   if (kstat < 0) goto error;
187 
188   /* Find the coefficients of the sum/difference curve.  */
189 
190   for (ki=0; ki<ktau*kdim; ki++)
191     scoef[ki] += (double)isign*scoef[ktau*kdim+ki];
192 
193   /* Present the sum/difference as a curve.  */
194 
195   if ((*rcurve = newCurve(ktau,korder,stau,scoef,kkind,kdim,kcopy)) == SISL_NULL)
196     goto err101;
197 
198   (*rcurve)->cuopen = kopen = MAX(pc1->cuopen,pc2->cuopen);
199   if (kopen == SISL_CRV_PERIODIC)
200   {
201      /* Represent the output curve cyclic. First find continuity. */
202 
203      s1333_count(knbcrv,ucurves,&kcont,&kstat);
204      if (kstat < 0) goto error;
205 
206      make_cv_cyclic(*rcurve,kcont,&kstat);
207      if (kstat < 0) goto error;
208   }
209 
210   /* Task performed.  */
211 
212   *jstat = 0;
213   goto out;
214 
215   /* Error in scratch allocation.  */
216 
217   err101 :
218     *jstat = -101;
219   goto out;
220 
221   /* Error in input. Conflicting dimensions.  */
222 
223   err106 :
224     *jstat = -106;
225   goto out;
226 
227   /* Error in lower order routine.  */
228 
229   error :
230     *jstat = kstat;
231   goto out;
232 
233   out :
234 
235     /* Free scratch occupied by local arrays.  */
236 
237     if (st != SISL_NULL) freearray(st);
238   if (stau != SISL_NULL) freearray(stau);
239   if (scoef != SISL_NULL) freearray(scoef);
240   if (ucurves[0] != SISL_NULL) freeCurve(ucurves[0]);
241   if (ucurves[1] != SISL_NULL) freeCurve(ucurves[1]);
242 
243   return;
244 }
245 
246 #if defined(SISLNEEDPROTOTYPES)
247 static void
s6addcurve_s9moveknots(double et[],int inmbknot,double astart,double aend,double etmoved[],int * jstat)248   s6addcurve_s9moveknots(double et[],int inmbknot,double astart,
249 			 double aend,double etmoved[],int *jstat)
250 #else
251 static void s6addcurve_s9moveknots(et,inmbknot,astart,aend,etmoved,jstat)
252      double et[];
253      int inmbknot;
254      double astart;
255      double aend;
256      double etmoved[];
257      int *jstat;
258 #endif
259 /*
260 *********************************************************************
261 *
262 * PURPOSE    : Move knot vector into the parameter interval
263 *              [astart,aend].
264 *
265 *
266 *
267 * INPUT      : et       - Knot vector.
268 *              inmbknot - Number of knots in the vector et.
269 *              astart   - Start parameter value of new parameter interval.
270 *              aend     - End parameter value of new parameter interval.
271 *
272 *
273 * OUTPUT     : etmoved  - Knot vector on new parameter interval.
274 *              jstat    - status messages
275 *                         = 1 : Parameter interval changed. Thus,
276 *                               reparametrization is performed.
277 *                         = 0 : ok
278 *                         < 0 : error
279 *
280 *
281 * METHOD     :
282 *
283 *
284 * REFERENCES :
285 *
286 *
287 * USE        :
288 *
289 *-
290 * CALLS      :
291 *
292 *
293 * WRITTEN BY : Vibeke Skytt, SI, 07.90.
294 *
295 *********************************************************************
296 */
297 {
298   int kwarn = 0;             /* Indicates if reparametrization is performed. */
299   double tstart = et[0];     /* Start parameter of original knot vector.     */
300   double tend = et[inmbknot-1];   /* End parameter of original knot vector.  */
301   double tfac = (aend - astart)/(tend - tstart);  /* Factor used when moving
302 						     knot vector.            */
303   double *s1,*s2,*s3;        /* Pointers used to traverse knot vectors.      */
304 
305   /* Test if reparametrization is performed.  */
306 
307   if (DNEQUAL(aend-astart,tend-tstart)) kwarn = 1;
308 
309   /* Put knot vector into new parameter interval.  */
310 
311   for (s1=et,s2=et+inmbknot,s3=etmoved; s1<s2; s1++,s3++)
312     *s3 = astart + (*s1 - tstart)*tfac;
313 
314   *jstat = kwarn;
315   return;
316 }
317 
318 
319