1cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 2cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell/* 3cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * Mesa 3-D graphics library 4cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * Version: 3.5 5cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * 622144ab7552f0799bcfca506bf4ffa7f70a06649Gareth Hughes * Copyright (C) 1999-2001 Brian Paul All Rights Reserved. 7cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * 8cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * Permission is hereby granted, free of charge, to any person obtaining a 9cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * copy of this software and associated documentation files (the "Software"), 10cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * to deal in the Software without restriction, including without limitation 11cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * the rights to use, copy, modify, merge, publish, distribute, sublicense, 12cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * and/or sell copies of the Software, and to permit persons to whom the 13cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * Software is furnished to do so, subject to the following conditions: 14cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * 15cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * The above copyright notice and this permission notice shall be included 16cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * in all copies or substantial portions of the Software. 17cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * 18cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS 19cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, 20cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL 21cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * BRIAN PAUL BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN 22cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN 23cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. 24cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell */ 25cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 26cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 27cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell/* 28cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * eval.c was written by 29cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * Bernd Barsuhn (bdbarsuh@cip.informatik.uni-erlangen.de) and 30cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * Volker Weiss (vrweiss@cip.informatik.uni-erlangen.de). 31cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * 32cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * My original implementation of evaluators was simplistic and didn't 33cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * compute surface normal vectors properly. Bernd and Volker applied 34cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * used more sophisticated methods to get better results. 35cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * 36cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * Thanks guys! 37cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell */ 38cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 39cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 40c223c6b663cd5db39ba19c2be74b88cc3b8f53f3Brian#include "main/glheader.h" 41c223c6b663cd5db39ba19c2be74b88cc3b8f53f3Brian#include "main/config.h" 42cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell#include "m_eval.h" 43cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 44cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwellstatic GLfloat inv_tab[MAX_EVAL_ORDER]; 45cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 46cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 47cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 48cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell/* 49cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * Horner scheme for Bezier curves 50cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * 51cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * Bezier curves can be computed via a Horner scheme. 52cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * Horner is numerically less stable than the de Casteljau 53cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * algorithm, but it is faster. For curves of degree n 54cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * the complexity of Horner is O(n) and de Casteljau is O(n^2). 55cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * Since stability is not important for displaying curve 56cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * points I decided to use the Horner scheme. 57cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * 58cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * A cubic Bezier curve with control points b0, b1, b2, b3 can be 59cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * written as 60cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * 61cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * (([3] [3] ) [3] ) [3] 62cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * c(t) = (([0]*s*b0 + [1]*t*b1)*s + [2]*t^2*b2)*s + [3]*t^2*b3 63cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * 64cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * [n] 65cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * where s=1-t and the binomial coefficients [i]. These can 66cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * be computed iteratively using the identity: 67cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * 68cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * [n] [n ] [n] 69cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * [i] = (n-i+1)/i * [i-1] and [0] = 1 70cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell */ 71cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 72cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 73cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwellvoid 74896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul_math_horner_bezier_curve(const GLfloat * cp, GLfloat * out, GLfloat t, 75cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell GLuint dim, GLuint order) 76cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell{ 77417ed16a88bd6c695e9792c2023e3f1737ee1e64Brian Paul GLfloat s, powert, bincoeff; 78417ed16a88bd6c695e9792c2023e3f1737ee1e64Brian Paul GLuint i, k; 79cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 80896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul if (order >= 2) { 81417ed16a88bd6c695e9792c2023e3f1737ee1e64Brian Paul bincoeff = (GLfloat) (order - 1); 827b9fe820a3fba3849864682fbb1cb512362934abKarl Schultz s = 1.0F - t; 83cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 84896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (k = 0; k < dim; k++) 85896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul out[k] = s * cp[k] + bincoeff * t * cp[dim + k]; 86cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 87896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (i = 2, cp += 2 * dim, powert = t * t; i < order; 88896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul i++, powert *= t, cp += dim) { 89417ed16a88bd6c695e9792c2023e3f1737ee1e64Brian Paul bincoeff *= (GLfloat) (order - i); 90417ed16a88bd6c695e9792c2023e3f1737ee1e64Brian Paul bincoeff *= inv_tab[i]; 91cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 92896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (k = 0; k < dim; k++) 93896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul out[k] = s * out[k] + bincoeff * powert * cp[k]; 94cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 95cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 96896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul else { /* order=1 -> constant curve */ 97896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul 98896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (k = 0; k < dim; k++) 99cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell out[k] = cp[k]; 100cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 101cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell} 102cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 103cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell/* 104cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * Tensor product Bezier surfaces 105cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * 106cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * Again the Horner scheme is used to compute a point on a 107cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * TP Bezier surface. First a control polygon for a curve 108cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * on the surface in one parameter direction is computed, 109cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * then the point on the curve for the other parameter 110cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * direction is evaluated. 111cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * 112cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * To store the curve control polygon additional storage 113cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * for max(uorder,vorder) points is needed in the 114cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * control net cn. 115cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell */ 116cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 117cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwellvoid 118896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul_math_horner_bezier_surf(GLfloat * cn, GLfloat * out, GLfloat u, GLfloat v, 119cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell GLuint dim, GLuint uorder, GLuint vorder) 120cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell{ 121896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul GLfloat *cp = cn + uorder * vorder * dim; 122896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul GLuint i, uinc = vorder * dim; 123cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 124896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul if (vorder > uorder) { 125896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul if (uorder >= 2) { 126417ed16a88bd6c695e9792c2023e3f1737ee1e64Brian Paul GLfloat s, poweru, bincoeff; 127417ed16a88bd6c695e9792c2023e3f1737ee1e64Brian Paul GLuint j, k; 128cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 129cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* Compute the control polygon for the surface-curve in u-direction */ 130896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (j = 0; j < vorder; j++) { 131896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul GLfloat *ucp = &cn[j * dim]; 132cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 133cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* Each control point is the point for parameter u on a */ 134cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* curve defined by the control polygons in u-direction */ 135417ed16a88bd6c695e9792c2023e3f1737ee1e64Brian Paul bincoeff = (GLfloat) (uorder - 1); 1367b9fe820a3fba3849864682fbb1cb512362934abKarl Schultz s = 1.0F - u; 137cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 138896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (k = 0; k < dim; k++) 139896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul cp[j * dim + k] = s * ucp[k] + bincoeff * u * ucp[uinc + k]; 140cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 141896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (i = 2, ucp += 2 * uinc, poweru = u * u; i < uorder; 142896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul i++, poweru *= u, ucp += uinc) { 143417ed16a88bd6c695e9792c2023e3f1737ee1e64Brian Paul bincoeff *= (GLfloat) (uorder - i); 144417ed16a88bd6c695e9792c2023e3f1737ee1e64Brian Paul bincoeff *= inv_tab[i]; 145cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 146896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (k = 0; k < dim; k++) 147896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul cp[j * dim + k] = 148896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul s * cp[j * dim + k] + bincoeff * poweru * ucp[k]; 149cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 150cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 151cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 152cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* Evaluate curve point in v */ 153cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell _math_horner_bezier_curve(cp, out, v, dim, vorder); 154cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 155896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul else /* uorder=1 -> cn defines a curve in v */ 156cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell _math_horner_bezier_curve(cn, out, v, dim, vorder); 157cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 158896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul else { /* vorder <= uorder */ 159896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul 160896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul if (vorder > 1) { 161cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell GLuint i; 162cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 163cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* Compute the control polygon for the surface-curve in u-direction */ 164896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (i = 0; i < uorder; i++, cn += uinc) { 165cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* For constant i all cn[i][j] (j=0..vorder) are located */ 166cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* on consecutive memory locations, so we can use */ 167cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* horner_bezier_curve to compute the control points */ 168cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 169896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul _math_horner_bezier_curve(cn, &cp[i * dim], v, dim, vorder); 170cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 171cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 172cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* Evaluate curve point in u */ 173cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell _math_horner_bezier_curve(cp, out, u, dim, uorder); 174cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 175896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul else /* vorder=1 -> cn defines a curve in u */ 176cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell _math_horner_bezier_curve(cn, out, u, dim, uorder); 177cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 178cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell} 179cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 180cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell/* 181cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * The direct de Casteljau algorithm is used when a point on the 182cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * surface and the tangent directions spanning the tangent plane 183cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * should be computed (this is needed to compute normals to the 184cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * surface). In this case the de Casteljau algorithm approach is 185cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * nicer because a point and the partial derivatives can be computed 186cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * at the same time. To get the correct tangent length du and dv 187cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * must be multiplied with the (u2-u1)/uorder-1 and (v2-v1)/vorder-1. 188cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * Since only the directions are needed, this scaling step is omitted. 189cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * 190cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * De Casteljau needs additional storage for uorder*vorder 191cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * values in the control net cn. 192cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell */ 193cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 194cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwellvoid 195896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul_math_de_casteljau_surf(GLfloat * cn, GLfloat * out, GLfloat * du, 196896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul GLfloat * dv, GLfloat u, GLfloat v, GLuint dim, 197cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell GLuint uorder, GLuint vorder) 198cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell{ 199896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul GLfloat *dcn = cn + uorder * vorder * dim; 2007b9fe820a3fba3849864682fbb1cb512362934abKarl Schultz GLfloat us = 1.0F - u, vs = 1.0F - v; 201cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell GLuint h, i, j, k; 202cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell GLuint minorder = uorder < vorder ? uorder : vorder; 203896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul GLuint uinc = vorder * dim; 204cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell GLuint dcuinc = vorder; 205cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 206cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* Each component is evaluated separately to save buffer space */ 207cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* This does not drasticaly decrease the performance of the */ 208cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* algorithm. If additional storage for (uorder-1)*(vorder-1) */ 209cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* points would be available, the components could be accessed */ 210cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* in the innermost loop which could lead to less cache misses. */ 211cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 212cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell#define CN(I,J,K) cn[(I)*uinc+(J)*dim+(K)] 213cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell#define DCN(I, J) dcn[(I)*dcuinc+(J)] 214896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul if (minorder < 3) { 215896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul if (uorder == vorder) { 216896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (k = 0; k < dim; k++) { 217cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* Derivative direction in u */ 218896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul du[k] = vs * (CN(1, 0, k) - CN(0, 0, k)) + 219896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul v * (CN(1, 1, k) - CN(0, 1, k)); 220cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 221cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* Derivative direction in v */ 222896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul dv[k] = us * (CN(0, 1, k) - CN(0, 0, k)) + 223896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul u * (CN(1, 1, k) - CN(1, 0, k)); 224cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 225cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* bilinear de Casteljau step */ 226896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul out[k] = us * (vs * CN(0, 0, k) + v * CN(0, 1, k)) + 227896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul u * (vs * CN(1, 0, k) + v * CN(1, 1, k)); 228cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 229cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 230896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul else if (minorder == uorder) { 231896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (k = 0; k < dim; k++) { 232cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* bilinear de Casteljau step */ 233896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(1, 0) = CN(1, 0, k) - CN(0, 0, k); 234896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(0, 0) = us * CN(0, 0, k) + u * CN(1, 0, k); 235cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 236896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (j = 0; j < vorder - 1; j++) { 237cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* for the derivative in u */ 238896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(1, j + 1) = CN(1, j + 1, k) - CN(0, j + 1, k); 239896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(1, j) = vs * DCN(1, j) + v * DCN(1, j + 1); 240cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 241cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* for the `point' */ 242896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(0, j + 1) = us * CN(0, j + 1, k) + u * CN(1, j + 1, k); 243896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(0, j) = vs * DCN(0, j) + v * DCN(0, j + 1); 244cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 245cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 246cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* remaining linear de Casteljau steps until the second last step */ 247896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (h = minorder; h < vorder - 1; h++) 248896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (j = 0; j < vorder - h; j++) { 249cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* for the derivative in u */ 250896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(1, j) = vs * DCN(1, j) + v * DCN(1, j + 1); 251cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 252cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* for the `point' */ 253896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(0, j) = vs * DCN(0, j) + v * DCN(0, j + 1); 254cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 255cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 256cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* derivative direction in v */ 257896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul dv[k] = DCN(0, 1) - DCN(0, 0); 258cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 259cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* derivative direction in u */ 260896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul du[k] = vs * DCN(1, 0) + v * DCN(1, 1); 261cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 262cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* last linear de Casteljau step */ 263896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul out[k] = vs * DCN(0, 0) + v * DCN(0, 1); 264cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 265cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 266896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul else { /* minorder == vorder */ 267896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul 268896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (k = 0; k < dim; k++) { 269cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* bilinear de Casteljau step */ 270896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(0, 1) = CN(0, 1, k) - CN(0, 0, k); 271896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(0, 0) = vs * CN(0, 0, k) + v * CN(0, 1, k); 272896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (i = 0; i < uorder - 1; i++) { 273cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* for the derivative in v */ 274896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(i + 1, 1) = CN(i + 1, 1, k) - CN(i + 1, 0, k); 275896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(i, 1) = us * DCN(i, 1) + u * DCN(i + 1, 1); 276cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 277cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* for the `point' */ 278896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(i + 1, 0) = vs * CN(i + 1, 0, k) + v * CN(i + 1, 1, k); 279896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(i, 0) = us * DCN(i, 0) + u * DCN(i + 1, 0); 280cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 281cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 282cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* remaining linear de Casteljau steps until the second last step */ 283896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (h = minorder; h < uorder - 1; h++) 284896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (i = 0; i < uorder - h; i++) { 285cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* for the derivative in v */ 286896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(i, 1) = us * DCN(i, 1) + u * DCN(i + 1, 1); 287cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 288cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* for the `point' */ 289896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(i, 0) = us * DCN(i, 0) + u * DCN(i + 1, 0); 290cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 291cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 292cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* derivative direction in u */ 293896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul du[k] = DCN(1, 0) - DCN(0, 0); 294cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 295cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* derivative direction in v */ 296896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul dv[k] = us * DCN(0, 1) + u * DCN(1, 1); 297cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 298cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* last linear de Casteljau step */ 299896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul out[k] = us * DCN(0, 0) + u * DCN(1, 0); 300cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 301cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 302cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 303896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul else if (uorder == vorder) { 304896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (k = 0; k < dim; k++) { 305cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* first bilinear de Casteljau step */ 306896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (i = 0; i < uorder - 1; i++) { 307896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(i, 0) = us * CN(i, 0, k) + u * CN(i + 1, 0, k); 308896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (j = 0; j < vorder - 1; j++) { 309896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(i, j + 1) = us * CN(i, j + 1, k) + u * CN(i + 1, j + 1, k); 310896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(i, j) = vs * DCN(i, j) + v * DCN(i, j + 1); 311cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 312cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 313cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 314cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* remaining bilinear de Casteljau steps until the second last step */ 315896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (h = 2; h < minorder - 1; h++) 316896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (i = 0; i < uorder - h; i++) { 317896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(i, 0) = us * DCN(i, 0) + u * DCN(i + 1, 0); 318896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (j = 0; j < vorder - h; j++) { 319896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(i, j + 1) = us * DCN(i, j + 1) + u * DCN(i + 1, j + 1); 320896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(i, j) = vs * DCN(i, j) + v * DCN(i, j + 1); 321cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 322cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 323cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 324cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* derivative direction in u */ 325896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul du[k] = vs * (DCN(1, 0) - DCN(0, 0)) + v * (DCN(1, 1) - DCN(0, 1)); 326cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 327cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* derivative direction in v */ 328896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul dv[k] = us * (DCN(0, 1) - DCN(0, 0)) + u * (DCN(1, 1) - DCN(1, 0)); 329cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 330cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* last bilinear de Casteljau step */ 331896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul out[k] = us * (vs * DCN(0, 0) + v * DCN(0, 1)) + 332896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul u * (vs * DCN(1, 0) + v * DCN(1, 1)); 333cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 334cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 335896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul else if (minorder == uorder) { 336896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (k = 0; k < dim; k++) { 337cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* first bilinear de Casteljau step */ 338896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (i = 0; i < uorder - 1; i++) { 339896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(i, 0) = us * CN(i, 0, k) + u * CN(i + 1, 0, k); 340896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (j = 0; j < vorder - 1; j++) { 341896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(i, j + 1) = us * CN(i, j + 1, k) + u * CN(i + 1, j + 1, k); 342896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(i, j) = vs * DCN(i, j) + v * DCN(i, j + 1); 343cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 344cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 345cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 346cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* remaining bilinear de Casteljau steps until the second last step */ 347896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (h = 2; h < minorder - 1; h++) 348896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (i = 0; i < uorder - h; i++) { 349896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(i, 0) = us * DCN(i, 0) + u * DCN(i + 1, 0); 350896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (j = 0; j < vorder - h; j++) { 351896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(i, j + 1) = us * DCN(i, j + 1) + u * DCN(i + 1, j + 1); 352896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(i, j) = vs * DCN(i, j) + v * DCN(i, j + 1); 353cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 354cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 355cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 356cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* last bilinear de Casteljau step */ 357896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(2, 0) = DCN(1, 0) - DCN(0, 0); 358896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(0, 0) = us * DCN(0, 0) + u * DCN(1, 0); 359896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (j = 0; j < vorder - 1; j++) { 360cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* for the derivative in u */ 361896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(2, j + 1) = DCN(1, j + 1) - DCN(0, j + 1); 362896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(2, j) = vs * DCN(2, j) + v * DCN(2, j + 1); 363896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul 364cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* for the `point' */ 365896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(0, j + 1) = us * DCN(0, j + 1) + u * DCN(1, j + 1); 366896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(0, j) = vs * DCN(0, j) + v * DCN(0, j + 1); 367cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 368cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 369cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* remaining linear de Casteljau steps until the second last step */ 370896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (h = minorder; h < vorder - 1; h++) 371896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (j = 0; j < vorder - h; j++) { 372cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* for the derivative in u */ 373896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(2, j) = vs * DCN(2, j) + v * DCN(2, j + 1); 374896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul 375cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* for the `point' */ 376896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(0, j) = vs * DCN(0, j) + v * DCN(0, j + 1); 377cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 378cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 379cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* derivative direction in v */ 380896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul dv[k] = DCN(0, 1) - DCN(0, 0); 381cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 382cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* derivative direction in u */ 383896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul du[k] = vs * DCN(2, 0) + v * DCN(2, 1); 384cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 385cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* last linear de Casteljau step */ 386896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul out[k] = vs * DCN(0, 0) + v * DCN(0, 1); 387cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 388cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 389896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul else { /* minorder == vorder */ 390896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul 391896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (k = 0; k < dim; k++) { 392cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* first bilinear de Casteljau step */ 393896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (i = 0; i < uorder - 1; i++) { 394896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(i, 0) = us * CN(i, 0, k) + u * CN(i + 1, 0, k); 395896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (j = 0; j < vorder - 1; j++) { 396896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(i, j + 1) = us * CN(i, j + 1, k) + u * CN(i + 1, j + 1, k); 397896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(i, j) = vs * DCN(i, j) + v * DCN(i, j + 1); 398cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 399cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 400cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 401cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* remaining bilinear de Casteljau steps until the second last step */ 402896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (h = 2; h < minorder - 1; h++) 403896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (i = 0; i < uorder - h; i++) { 404896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(i, 0) = us * DCN(i, 0) + u * DCN(i + 1, 0); 405896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (j = 0; j < vorder - h; j++) { 406896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(i, j + 1) = us * DCN(i, j + 1) + u * DCN(i + 1, j + 1); 407896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(i, j) = vs * DCN(i, j) + v * DCN(i, j + 1); 408cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 409cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 410cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 411cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* last bilinear de Casteljau step */ 412896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(0, 2) = DCN(0, 1) - DCN(0, 0); 413896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(0, 0) = vs * DCN(0, 0) + v * DCN(0, 1); 414896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (i = 0; i < uorder - 1; i++) { 415cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* for the derivative in v */ 416896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(i + 1, 2) = DCN(i + 1, 1) - DCN(i + 1, 0); 417896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(i, 2) = us * DCN(i, 2) + u * DCN(i + 1, 2); 418896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul 419cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* for the `point' */ 420896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(i + 1, 0) = vs * DCN(i + 1, 0) + v * DCN(i + 1, 1); 421896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(i, 0) = us * DCN(i, 0) + u * DCN(i + 1, 0); 422cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 423cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 424cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* remaining linear de Casteljau steps until the second last step */ 425896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (h = minorder; h < uorder - 1; h++) 426896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (i = 0; i < uorder - h; i++) { 427cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* for the derivative in v */ 428896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(i, 2) = us * DCN(i, 2) + u * DCN(i + 1, 2); 429896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul 430cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* for the `point' */ 431896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul DCN(i, 0) = us * DCN(i, 0) + u * DCN(i + 1, 0); 432cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 433cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 434cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* derivative direction in u */ 435896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul du[k] = DCN(1, 0) - DCN(0, 0); 436cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 437cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* derivative direction in v */ 438896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul dv[k] = us * DCN(0, 2) + u * DCN(1, 2); 439cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 440cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* last linear de Casteljau step */ 441896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul out[k] = us * DCN(0, 0) + u * DCN(1, 0); 442cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 443cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell } 444cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell#undef DCN 445cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell#undef CN 446cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell} 447cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 448cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 449cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell/* 450cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell * Do one-time initialization for evaluators. 451cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell */ 452896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paulvoid 453896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul_math_init_eval(void) 454cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell{ 455cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell GLuint i; 456cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell 457cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell /* KW: precompute 1/x for useful x. 458cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell */ 459896e8bd2d7eb1385ca89e71b7eac146577320e00Brian Paul for (i = 1; i < MAX_EVAL_ORDER; i++) 4607b9fe820a3fba3849864682fbb1cb512362934abKarl Schultz inv_tab[i] = 1.0F / i; 461cab974cf6c2dbfbf5dd5d291e1aae0f8eeb34290Keith Whitwell} 462