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:
45 *
46 */
47 #define S1541
48
49 #include "sislP.h"
50
51 #if defined(SISLNEEDPROTOTYPES)
52 void
s1541(SISLCurve * pc1,int npol,double ebder[],int ileft[],double eder[],int * jstat)53 s1541 ( SISLCurve *pc1,
54 int npol,
55 double ebder[],
56 int ileft[],
57 double eder[],
58 int *jstat )
59
60 #else
61 void
62 s1541 ( pc1, npol, ebder, ileft, eder, jstat )
63 SISLCurve* pc1;
64 int npol;
65 double ebder[];
66 int ileft[];
67 double eder[];
68 int *jstat;
69
70 #endif
71 /*
72 *********************************************************************
73 *
74 * PURPOSE: Given a (polynomial) spline curve pc1 and
75 * preevaluated basis functions (using s1540()) on an
76 * npol polyline, calculate the 3D positions on
77 * that polyline (eder).
78 *
79 * INPUT: pc1 - the spline curve
80 * npol - the number of polyline points
81 * ebder - array containing the basis values
82 * B(ax[0 ],i0-k+1),...,B(ax[0 ],i0)
83 * B(ax[1 ],i1-k+1),...,B(ax[1 ],i1)
84 * : :
85 * B(ax[m1-1],im1-1-k+1),...,B(ax[m1-1],im1-1)
86 *
87 * ileft - ileft[i] <= ti < ileft[i] + 1
88 * (exception for ti == in : ileft[ti] = n-1)
89 *
90 * OUTPUT: eder - contains the 3D polyline points
91 * jstat - Status messages
92 *
93 * = 0 : Ok.
94 * < 0 : Error.
95 *
96 *
97 * METHOD:
98 *
99 * WRITTEN BY: Geir Westgaard, SINTEF, Oslo, November 1999
100 * REVISED BY: Vibeke Skytt, SINTEF, Dec. 2006. Allow dimension different from 3
101 * and rational curves
102 *
103 *********************************************************************
104 */
105 {
106 int m = 0, m1 = 0;
107 int my, my1;
108 int np, i, j;
109 int ik;
110 double bas;
111 double* ecoef = SISL_NULL;
112 int kdim = pc1->idim;
113 double scratch[4];
114 double *xyz = NULL;
115 int krat = (pc1->ikind == 2 || pc1->ikind == 4);
116
117
118 /* Check the input. */
119
120 //if ( pc1->idim != 3 ) goto err104;
121
122 if (krat)
123 kdim++;
124
125 if (kdim > 4)
126 {
127 if ((xyz = newarray(kdim, DOUBLE)) == SISL_NULL)
128 goto err101;
129 }
130 else xyz = scratch;
131
132
133 /* Set input to local variables. */
134
135 ik = pc1 -> ik;
136 ecoef = (krat) ? pc1->rcoef : pc1 -> ecoef;
137
138
139 for ( np = 0; np < npol; np++ )
140 {
141 my = ileft[ np ] - ik;
142
143 for (j=0; j<kdim; j++)
144 xyz[j] = 0.0;
145
146 for ( i = 0; i < ik; i++ )
147 {
148 my++;
149 my1 = kdim*my;
150 bas = ebder[ m1++ ];
151
152 for (j=0; j<kdim; j++)
153 xyz[j] += ecoef[my1+j]*bas;
154 }
155
156 if (krat)
157 {
158 for (j=0; j<pc1->idim; j++)
159 xyz[j] /= xyz[pc1->idim];
160 }
161
162 for (j=0; j<pc1->idim; j++)
163 eder[m++] = xyz[j];
164 }
165
166
167 /* Successful computations. */
168
169 *jstat = 0;
170 goto out;
171
172 /* Error in input, crv->idim != 3 */
173 //err104: *jstat = -104;
174 // s6err( "s1541", *jstat, 0 );
175 // goto out;
176
177 err101: *jstat = -101;
178 s6err( "s1541", *jstat, 0 );
179 goto out;
180
181 out:
182 if (xyz != SISL_NULL && xyz != scratch)
183 freearray(xyz);
184 return;
185
186 }
187