2 ** License Applicability. Except to the extent portions of this file are
3 ** made subject to an alternative license as permitted in the SGI Free
4 ** Software License B, Version 1.1 (the "License"), the contents of this
5 ** file are subject only to the provisions of the License. You may not use
6 ** this file except in compliance with the License. You may obtain a copy
7 ** of the License at Silicon Graphics, Inc., attn: Legal Services, 1600
8 ** Amphitheatre Parkway, Mountain View, CA 94043-1351, or at:
10 ** http://oss.sgi.com/projects/FreeB
12 ** Note that, as provided in the License, the Software is distributed on an
13 ** "AS IS" basis, with ALL EXPRESS AND IMPLIED WARRANTIES AND CONDITIONS
14 ** DISCLAIMED, INCLUDING, WITHOUT LIMITATION, ANY IMPLIED WARRANTIES AND
15 ** CONDITIONS OF MERCHANTABILITY, SATISFACTORY QUALITY, FITNESS FOR A
16 ** PARTICULAR PURPOSE, AND NON-INFRINGEMENT.
18 ** Original Code. The Original Code is: OpenGL Sample Implementation,
19 ** Version 1.2.1, released January 26, 2000, developed by Silicon Graphics,
20 ** Inc. The Original Code is Copyright (c) 1991-2000 Silicon Graphics, Inc.
21 ** Copyright in any portions created by third parties is as indicated
22 ** elsewhere herein. All Rights Reserved.
24 ** Additional Notice Provisions: The application programming interfaces
25 ** established by SGI in conjunction with the Original Code are The
26 ** OpenGL(R) Graphics System: A Specification (Version 1.2.1), released
27 ** April 1, 1999; The OpenGL(R) Graphics System Utility Library (Version
28 ** 1.3), released November 4, 1998; and OpenGL(R) Graphics with the X
29 ** Window System(R) (Version 1.3), released October 19, 1998. This software
30 ** was created using the OpenGL(R) version 1.2.1 Sample Implementation
31 ** published by SGI, but has not been independently verified as being
32 ** compliant with the OpenGL(R) version 1.2.1 Specification.
42 #include "bezierEval.h"
48 #define TOLERANCE 0.0001
55 #define MAX_DIMENSION 4
58 static void normalize(float vec
[3]);
59 static void crossProduct(float x
[3], float y
[3], float ret
[3]);
61 static void bezierCurveEvalfast(float u0
, float u1
, int order
, float *ctlpoints
, int stride
, int dimension
, float u
, float retpoint
[]);
64 static float binomialCoefficients
[8][8] = {
75 void bezierCurveEval(float u0
, float u1
, int order
, float *ctlpoints
, int stride
, int dimension
, float u
, float retpoint
[])
77 float uprime
= (u
-u0
)/(u1
-u0
);
78 float *ctlptr
= ctlpoints
;
79 float oneMinusX
= 1.0f
-uprime
;
83 for(k
=0; k
<dimension
; k
++)
84 retpoint
[k
] = (*(ctlptr
+ k
));
86 for(i
=1; i
<order
; i
++){
89 for(k
=0; k
<dimension
; k
++) {
90 retpoint
[k
] = retpoint
[k
]*oneMinusX
+ ctlptr
[k
]* binomialCoefficients
[order
-1][i
] * XPower
;
97 /*order = degree +1 >=1.
99 void bezierCurveEvalfast(float u0
, float u1
, int order
, float *ctlpoints
, int stride
, int dimension
, float u
, float retpoint
[])
101 float uprime
= (u
-u0
)/(u1
-u0
);
102 float buf
[MAX_ORDER
][MAX_ORDER
][MAX_DIMENSION
];
103 float* ctlptr
= ctlpoints
;
105 for(i
=0; i
<order
; i
++) {
106 for(j
=0; j
<dimension
; j
++)
107 buf
[0][i
][j
] = ctlptr
[j
];
110 for(r
=1; r
<order
; r
++){
111 for(i
=0; i
<order
-r
; i
++) {
112 for(j
=0; j
<dimension
; j
++)
113 buf
[r
][i
][j
] = (1-uprime
)*buf
[r
-1][i
][j
] + uprime
*buf
[r
-1][i
+1][j
];
117 for(j
=0; j
<dimension
; j
++)
118 retpoint
[j
] = buf
[order
-1][0][j
];
123 /*order = degree +1 >=1.
125 void bezierCurveEvalDer(float u0
, float u1
, int order
, float *ctlpoints
, int stride
, int dimension
, float u
, float retDer
[])
129 float *ctlptr
= ctlpoints
;
131 float buf
[MAX_ORDER
][MAX_DIMENSION
];
133 for(k
=0; k
<dimension
; k
++)
136 for(i
=0; i
<order
-1; i
++){
137 for(k
=0; k
<dimension
; k
++) {
138 buf
[i
][k
] = (ctlptr
[stride
+k
] - ctlptr
[k
])*(order
-1)/width
;
143 bezierCurveEval(u0
, u1
, order
-1, (float*) buf
, MAX_DIMENSION
, dimension
, u
, retDer
);
146 void bezierCurveEvalDerGen(int der
, float u0
, float u1
, int order
, float *ctlpoints
, int stride
, int dimension
, float u
, float retDer
[])
149 float *ctlptr
= ctlpoints
;
151 float buf
[MAX_ORDER
][MAX_ORDER
][MAX_DIMENSION
];
153 for(i
=0; i
<order
; i
++){
154 for(k
=0; k
<dimension
; k
++){
155 buf
[0][i
][k
] = ctlptr
[k
];
161 for(r
=1; r
<=der
; r
++){
162 for(i
=0; i
<order
-r
; i
++){
163 for(k
=0; k
<dimension
; k
++){
164 buf
[r
][i
][k
] = (buf
[r
-1][i
+1][k
] - buf
[r
-1][i
][k
])*(order
-r
)/width
;
169 bezierCurveEval(u0
, u1
, order
-der
, (float *) (buf
[der
]), MAX_DIMENSION
, dimension
, u
, retDer
);
172 /*the Bezier bivarite polynomial is:
173 * sum[i:0,uorder-1][j:0,vorder-1] { ctlpoints[i*ustride+j*vstride] * B(i)*B(j)
174 * where B(i) and B(j) are basis functions
176 void bezierSurfEvalDerGen(int uder
, int vder
, float u0
, float u1
, int uorder
, float v0
, float v1
, int vorder
, int dimension
, float *ctlpoints
, int ustride
, int vstride
, float u
, float v
, float ret
[])
179 float newPoints
[MAX_ORDER
][MAX_DIMENSION
];
181 for(i
=0; i
<uorder
; i
++){
183 bezierCurveEvalDerGen(vder
, v0
, v1
, vorder
, ctlpoints
+ustride
*i
, vstride
, dimension
, v
, newPoints
[i
]);
187 bezierCurveEvalDerGen(uder
, u0
, u1
, uorder
, (float *) newPoints
, MAX_DIMENSION
, dimension
, u
, ret
);
191 /*division by w is performed*/
192 void bezierSurfEval(float u0
, float u1
, int uorder
, float v0
, float v1
, int vorder
, int dimension
, float *ctlpoints
, int ustride
, int vstride
, float u
, float v
, float ret
[])
194 bezierSurfEvalDerGen(0, 0, u0
, u1
, uorder
, v0
, v1
, vorder
, dimension
, ctlpoints
, ustride
, vstride
, u
, v
, ret
);
195 if(dimension
== 4) /*homogeneous*/{
202 void bezierSurfEvalNormal(float u0
, float u1
, int uorder
, float v0
, float v1
, int vorder
, int dimension
, float *ctlpoints
, int ustride
, int vstride
, float u
, float v
, float retNormal
[])
206 assert(dimension
>=3 && dimension
<=4);
207 bezierSurfEvalDerGen(1,0, u0
, u1
, uorder
, v0
, v1
, vorder
, dimension
, ctlpoints
, ustride
, vstride
, u
, v
, partialU
);
208 bezierSurfEvalDerGen(0,1, u0
, u1
, uorder
, v0
, v1
, vorder
, dimension
, ctlpoints
, ustride
, vstride
, u
, v
, partialV
);
210 if(dimension
== 3){/*inhomogeneous*/
211 crossProduct(partialU
, partialV
, retNormal
);
213 normalize(retNormal
);
217 else { /*homogeneous*/
218 float val
[4]; /*the point coordinates (without derivative)*/
219 float newPartialU
[MAX_DIMENSION
];
220 float newPartialV
[MAX_DIMENSION
];
222 bezierSurfEvalDerGen(0,0, u0
, u1
, uorder
, v0
, v1
, vorder
, dimension
, ctlpoints
, ustride
, vstride
, u
, v
, val
);
225 newPartialU
[i
] = partialU
[i
] * val
[3] - val
[i
] * partialU
[3];
226 newPartialV
[i
] = partialV
[i
] * val
[3] - val
[i
] * partialV
[3];
228 crossProduct(newPartialU
, newPartialV
, retNormal
);
229 normalize(retNormal
);
233 /*if size is 0, then nothing is done*/
234 static void normalize(float vec
[3])
236 float size
= (float)sqrt(vec
[0]*vec
[0] + vec
[1]*vec
[1] + vec
[2]*vec
[2]);
241 fprintf(stderr
, "Warning: in oglBSpline.c normal is 0\n");
246 vec
[0] = vec
[0]/size
;
247 vec
[1] = vec
[1]/size
;
248 vec
[2] = vec
[2]/size
;
253 static void crossProduct(float x
[3], float y
[3], float ret
[3])
255 ret
[0] = x
[1]*y
[2] - y
[1]*x
[2];
256 ret
[1] = x
[2]*y
[0] - y
[2]*x
[0];
257 ret
[2] = x
[0]*y
[1] - y
[0]*x
[1];