diff options
Diffstat (limited to 'src/glu/mesa/nurbsutl.c')
-rw-r--r-- | src/glu/mesa/nurbsutl.c | 1309 |
1 files changed, 1309 insertions, 0 deletions
diff --git a/src/glu/mesa/nurbsutl.c b/src/glu/mesa/nurbsutl.c new file mode 100644 index 00000000000..b46348a50b2 --- /dev/null +++ b/src/glu/mesa/nurbsutl.c @@ -0,0 +1,1309 @@ + +/* + * Mesa 3-D graphics library + * Version: 3.3 + * Copyright (C) 1995-2000 Brian Paul + * + * This library is free software; you can redistribute it and/or + * modify it under the terms of the GNU Library General Public + * License as published by the Free Software Foundation; either + * version 2 of the License, or (at your option) any later version. + * + * This library is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + * Library General Public License for more details. + * + * You should have received a copy of the GNU Library General Public + * License along with this library; if not, write to the Free + * Software Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA. + */ + + +/* + * NURBS implementation written by Bogdan Sikorski ([email protected]) + * See README2 for more info. + */ + + +#ifdef PC_HEADER +#include "all.h" +#else +#include <math.h> +#include <stdlib.h> +#include "gluP.h" +#include "nurbs.h" +#endif + + +GLenum test_knot(GLint nknots, GLfloat * knot, GLint order) +{ + GLsizei i; + GLint knot_mult; + GLfloat tmp_knot; + + tmp_knot = knot[0]; + knot_mult = 1; + for (i = 1; i < nknots; i++) { + if (knot[i] < tmp_knot) + return GLU_NURBS_ERROR4; + if (fabs(tmp_knot - knot[i]) > EPSILON) { + if (knot_mult > order) + return GLU_NURBS_ERROR5; + knot_mult = 1; + tmp_knot = knot[i]; + } + else + ++knot_mult; + } + return GLU_NO_ERROR; +} + +static int +/* qsort function */ +#if defined(WIN32) && !defined(OPENSTEP) + __cdecl +#endif +knot_sort(const void *a, const void *b) +{ + GLfloat x, y; + + x = *((GLfloat *) a); + y = *((GLfloat *) b); + if (fabs(x - y) < EPSILON) + return 0; + if (x > y) + return 1; + return -1; +} + +/* insert into dest knot all values within the valid range from src knot */ +/* that do not appear in dest */ +void +collect_unified_knot(knot_str_type * dest, knot_str_type * src, + GLfloat maximal_min_knot, GLfloat minimal_max_knot) +{ + GLfloat *src_knot, *dest_knot; + GLint src_t_min, src_t_max, dest_t_min, dest_t_max; + GLint src_nknots, dest_nknots; + GLint i, j, k, new_cnt; + GLboolean not_found_flag; + + src_knot = src->unified_knot; + dest_knot = dest->unified_knot; + src_t_min = src->t_min; + src_t_max = src->t_max; + dest_t_min = dest->t_min; + dest_t_max = dest->t_max; + src_nknots = src->unified_nknots; + dest_nknots = dest->unified_nknots; + + k = new_cnt = dest_nknots; + for (i = src_t_min; i <= src_t_max; i++) + if (src_knot[i] - maximal_min_knot > -EPSILON && + src_knot[i] - minimal_max_knot < EPSILON) { + not_found_flag = GL_TRUE; + for (j = dest_t_min; j <= dest_t_max; j++) + if (fabs(dest_knot[j] - src_knot[i]) < EPSILON) { + not_found_flag = GL_FALSE; + break; + } + if (not_found_flag) { + /* knot from src is not in dest - add this knot to dest */ + dest_knot[k++] = src_knot[i]; + ++new_cnt; + ++(dest->t_max); /* the valid range widens */ + ++(dest->delta_nknots); /* increment the extra knot value counter */ + } + } + dest->unified_nknots = new_cnt; + qsort((void *) dest_knot, (size_t) new_cnt, (size_t) sizeof(GLfloat), + &knot_sort); +} + +/* basing on the new common knot range for all attributes set */ +/* t_min and t_max values for each knot - they will be used later on */ +/* by explode_knot() and calc_new_ctrl_pts */ +static void +set_new_t_min_t_max(knot_str_type * geom_knot, knot_str_type * color_knot, + knot_str_type * normal_knot, knot_str_type * texture_knot, + GLfloat maximal_min_knot, GLfloat minimal_max_knot) +{ + GLuint t_min = 0, t_max = 0, cnt = 0; + + if (minimal_max_knot - maximal_min_knot < EPSILON) { + /* knot common range empty */ + geom_knot->t_min = geom_knot->t_max = 0; + color_knot->t_min = color_knot->t_max = 0; + normal_knot->t_min = normal_knot->t_max = 0; + texture_knot->t_min = texture_knot->t_max = 0; + } + else { + if (geom_knot->unified_knot != NULL) { + cnt = geom_knot->unified_nknots; + for (t_min = 0; t_min < cnt; t_min++) + if (fabs((geom_knot->unified_knot)[t_min] - maximal_min_knot) < + EPSILON) break; + for (t_max = cnt - 1; t_max; t_max--) + if (fabs((geom_knot->unified_knot)[t_max] - minimal_max_knot) < + EPSILON) break; + } + else if (geom_knot->nknots) { + cnt = geom_knot->nknots; + for (t_min = 0; t_min < cnt; t_min++) + if (fabs((geom_knot->knot)[t_min] - maximal_min_knot) < EPSILON) + break; + for (t_max = cnt - 1; t_max; t_max--) + if (fabs((geom_knot->knot)[t_max] - minimal_max_knot) < EPSILON) + break; + } + geom_knot->t_min = t_min; + geom_knot->t_max = t_max; + if (color_knot->unified_knot != NULL) { + cnt = color_knot->unified_nknots; + for (t_min = 0; t_min < cnt; t_min++) + if (fabs((color_knot->unified_knot)[t_min] - maximal_min_knot) < + EPSILON) break; + for (t_max = cnt - 1; t_max; t_max--) + if (fabs((color_knot->unified_knot)[t_max] - minimal_max_knot) < + EPSILON) break; + color_knot->t_min = t_min; + color_knot->t_max = t_max; + } + if (normal_knot->unified_knot != NULL) { + cnt = normal_knot->unified_nknots; + for (t_min = 0; t_min < cnt; t_min++) + if (fabs((normal_knot->unified_knot)[t_min] - maximal_min_knot) < + EPSILON) break; + for (t_max = cnt - 1; t_max; t_max--) + if (fabs((normal_knot->unified_knot)[t_max] - minimal_max_knot) < + EPSILON) break; + normal_knot->t_min = t_min; + normal_knot->t_max = t_max; + } + if (texture_knot->unified_knot != NULL) { + cnt = texture_knot->unified_nknots; + for (t_min = 0; t_min < cnt; t_min++) + if (fabs((texture_knot->unified_knot)[t_min] - maximal_min_knot) + < EPSILON) + break; + for (t_max = cnt - 1; t_max; t_max--) + if (fabs((texture_knot->unified_knot)[t_max] - minimal_max_knot) + < EPSILON) + break; + texture_knot->t_min = t_min; + texture_knot->t_max = t_max; + } + } +} + +/* modify all knot valid ranges in such a way that all have the same */ +/* range, common to all knots */ +/* do this by knot insertion */ +GLenum +select_knot_working_range(GLUnurbsObj * nobj, knot_str_type * geom_knot, + knot_str_type * color_knot, + knot_str_type * normal_knot, + knot_str_type * texture_knot) +{ + GLint max_nknots; + GLfloat maximal_min_knot, minimal_max_knot; + GLint i; + + /* find the maximum modified knot length */ + max_nknots = geom_knot->nknots; + if (color_knot->unified_knot) + max_nknots += color_knot->nknots; + if (normal_knot->unified_knot) + max_nknots += normal_knot->nknots; + if (texture_knot->unified_knot) + max_nknots += texture_knot->nknots; + maximal_min_knot = (geom_knot->knot)[geom_knot->t_min]; + minimal_max_knot = (geom_knot->knot)[geom_knot->t_max]; + /* any attirb data ? */ + if (max_nknots != geom_knot->nknots) { + /* allocate space for the unified knots */ + if ((geom_knot->unified_knot = + (GLfloat *) malloc(sizeof(GLfloat) * max_nknots)) == NULL) { + call_user_error(nobj, GLU_OUT_OF_MEMORY); + return GLU_ERROR; + } + /* copy the original knot to the unified one */ + geom_knot->unified_nknots = geom_knot->nknots; + for (i = 0; i < geom_knot->nknots; i++) + (geom_knot->unified_knot)[i] = (geom_knot->knot)[i]; + if (color_knot->unified_knot) { + if ((color_knot->knot)[color_knot->t_min] - maximal_min_knot > + EPSILON) + maximal_min_knot = (color_knot->knot)[color_knot->t_min]; + if (minimal_max_knot - (color_knot->knot)[color_knot->t_max] > + EPSILON) + minimal_max_knot = (color_knot->knot)[color_knot->t_max]; + if ((color_knot->unified_knot = + (GLfloat *) malloc(sizeof(GLfloat) * max_nknots)) == NULL) { + free(geom_knot->unified_knot); + call_user_error(nobj, GLU_OUT_OF_MEMORY); + return GLU_ERROR; + } + /* copy the original knot to the unified one */ + color_knot->unified_nknots = color_knot->nknots; + for (i = 0; i < color_knot->nknots; i++) + (color_knot->unified_knot)[i] = (color_knot->knot)[i]; + } + if (normal_knot->unified_knot) { + if ((normal_knot->knot)[normal_knot->t_min] - maximal_min_knot > + EPSILON) + maximal_min_knot = (normal_knot->knot)[normal_knot->t_min]; + if (minimal_max_knot - (normal_knot->knot)[normal_knot->t_max] > + EPSILON) + minimal_max_knot = (normal_knot->knot)[normal_knot->t_max]; + if ((normal_knot->unified_knot = + (GLfloat *) malloc(sizeof(GLfloat) * max_nknots)) == NULL) { + free(geom_knot->unified_knot); + free(color_knot->unified_knot); + call_user_error(nobj, GLU_OUT_OF_MEMORY); + return GLU_ERROR; + } + /* copy the original knot to the unified one */ + normal_knot->unified_nknots = normal_knot->nknots; + for (i = 0; i < normal_knot->nknots; i++) + (normal_knot->unified_knot)[i] = (normal_knot->knot)[i]; + } + if (texture_knot->unified_knot) { + if ((texture_knot->knot)[texture_knot->t_min] - maximal_min_knot > + EPSILON) + maximal_min_knot = (texture_knot->knot)[texture_knot->t_min]; + if (minimal_max_knot - (texture_knot->knot)[texture_knot->t_max] > + EPSILON) + minimal_max_knot = (texture_knot->knot)[texture_knot->t_max]; + if ((texture_knot->unified_knot = + (GLfloat *) malloc(sizeof(GLfloat) * max_nknots)) == NULL) { + free(geom_knot->unified_knot); + free(color_knot->unified_knot); + free(normal_knot->unified_knot); + call_user_error(nobj, GLU_OUT_OF_MEMORY); + return GLU_ERROR; + } + /* copy the original knot to the unified one */ + texture_knot->unified_nknots = texture_knot->nknots; + for (i = 0; i < texture_knot->nknots; i++) + (texture_knot->unified_knot)[i] = (texture_knot->knot)[i]; + } + /* work on the geometry knot with all additional knot values */ + /* appearing in attirbutive knots */ + if (minimal_max_knot - maximal_min_knot < EPSILON) { + /* empty working range */ + geom_knot->unified_nknots = 0; + color_knot->unified_nknots = 0; + normal_knot->unified_nknots = 0; + texture_knot->unified_nknots = 0; + } + else { + if (color_knot->unified_knot) + collect_unified_knot(geom_knot, color_knot, maximal_min_knot, + minimal_max_knot); + if (normal_knot->unified_knot) + collect_unified_knot(geom_knot, normal_knot, maximal_min_knot, + minimal_max_knot); + if (texture_knot->unified_knot) + collect_unified_knot(geom_knot, texture_knot, maximal_min_knot, + minimal_max_knot); + /* since we have now built the "unified" geometry knot */ + /* add same knot values to all attributive knots */ + if (color_knot->unified_knot) + collect_unified_knot(color_knot, geom_knot, maximal_min_knot, + minimal_max_knot); + if (normal_knot->unified_knot) + collect_unified_knot(normal_knot, geom_knot, maximal_min_knot, + minimal_max_knot); + if (texture_knot->unified_knot) + collect_unified_knot(texture_knot, geom_knot, maximal_min_knot, + minimal_max_knot); + } + } + set_new_t_min_t_max(geom_knot, color_knot, normal_knot, texture_knot, + maximal_min_knot, minimal_max_knot); + return GLU_NO_ERROR; +} + +void +free_unified_knots(knot_str_type * geom_knot, knot_str_type * color_knot, + knot_str_type * normal_knot, knot_str_type * texture_knot) +{ + if (geom_knot->unified_knot) + free(geom_knot->unified_knot); + if (color_knot->unified_knot) + free(color_knot->unified_knot); + if (normal_knot->unified_knot) + free(normal_knot->unified_knot); + if (texture_knot->unified_knot) + free(texture_knot->unified_knot); +} + +GLenum explode_knot(knot_str_type * the_knot) +{ + GLfloat *knot, *new_knot; + GLint nknots, n_new_knots = 0; + GLint t_min, t_max; + GLint ord; + GLsizei i, j, k; + GLfloat tmp_float; + + if (the_knot->unified_knot) { + knot = the_knot->unified_knot; + nknots = the_knot->unified_nknots; + } + else { + knot = the_knot->knot; + nknots = the_knot->nknots; + } + ord = the_knot->order; + t_min = the_knot->t_min; + t_max = the_knot->t_max; + + for (i = t_min; i <= t_max;) { + tmp_float = knot[i]; + for (j = 0; j < ord && (i + j) <= t_max; j++) + if (fabs(tmp_float - knot[i + j]) > EPSILON) + break; + n_new_knots += ord - j; + i += j; + } + /* alloc space for new_knot */ + if ( + (new_knot = + (GLfloat *) malloc(sizeof(GLfloat) * (nknots + n_new_knots + 1))) == NULL) { + return GLU_OUT_OF_MEMORY; + } + /* fill in new knot */ + for (j = 0; j < t_min; j++) + new_knot[j] = knot[j]; + for (i = j; i <= t_max; i++) { + tmp_float = knot[i]; + for (k = 0; k < ord; k++) { + new_knot[j++] = knot[i]; + if (tmp_float == knot[i + 1]) + i++; + } + } + for (i = t_max + 1; i < (int) nknots; i++) + new_knot[j++] = knot[i]; + /* fill in the knot structure */ + the_knot->new_knot = new_knot; + the_knot->delta_nknots += n_new_knots; + the_knot->t_max += n_new_knots; + return GLU_NO_ERROR; +} + +GLenum calc_alphas(knot_str_type * the_knot) +{ + GLfloat tmp_float; + int i, j, k, m, n; + int order; + GLfloat *alpha, *alpha_new, *tmp_alpha; + GLfloat denom; + GLfloat *knot, *new_knot; + + + knot = the_knot->knot; + order = the_knot->order; + new_knot = the_knot->new_knot; + n = the_knot->nknots - the_knot->order; + m = n + the_knot->delta_nknots; + if ((alpha = (GLfloat *) malloc(sizeof(GLfloat) * n * m)) == NULL) { + return GLU_OUT_OF_MEMORY; + } + if ((alpha_new = (GLfloat *) malloc(sizeof(GLfloat) * n * m)) == NULL) { + free(alpha); + return GLU_OUT_OF_MEMORY; + } + for (j = 0; j < m; j++) { + for (i = 0; i < n; i++) { + if ((knot[i] <= new_knot[j]) && (new_knot[j] < knot[i + 1])) + tmp_float = 1.0; + else + tmp_float = 0.0; + alpha[i + j * n] = tmp_float; + } + } + for (k = 1; k < order; k++) { + for (j = 0; j < m; j++) + for (i = 0; i < n; i++) { + denom = knot[i + k] - knot[i]; + if (fabs(denom) < EPSILON) + tmp_float = 0.0; + else + tmp_float = (new_knot[j + k] - knot[i]) / denom * + alpha[i + j * n]; + denom = knot[i + k + 1] - knot[i + 1]; + if (fabs(denom) > EPSILON) + tmp_float += (knot[i + k + 1] - new_knot[j + k]) / denom * + alpha[(i + 1) + j * n]; + alpha_new[i + j * n] = tmp_float; + } + tmp_alpha = alpha_new; + alpha_new = alpha; + alpha = tmp_alpha; + } + the_knot->alpha = alpha; + free(alpha_new); + return GLU_NO_ERROR; +} + +GLenum +calc_new_ctrl_pts(GLfloat * ctrl, GLint stride, knot_str_type * the_knot, + GLint dim, GLfloat ** new_ctrl, GLint * ncontrol) +{ + GLsizei i, j, k, l, m, n; + GLsizei index1, index2; + GLfloat *alpha; + GLfloat *new_knot; + + new_knot = the_knot->new_knot; + n = the_knot->nknots - the_knot->order; + alpha = the_knot->alpha; + + m = the_knot->t_max + 1 - the_knot->t_min - the_knot->order; + k = the_knot->t_min; + /* allocate space for new control points */ + if ((*new_ctrl = (GLfloat *) malloc(sizeof(GLfloat) * dim * m)) == NULL) { + return GLU_OUT_OF_MEMORY; + } + for (j = 0; j < m; j++) { + for (l = 0; l < dim; l++) + (*new_ctrl)[j * dim + l] = 0.0; + for (i = 0; i < n; i++) { + index1 = i + (j + k) * n; + index2 = i * stride; + for (l = 0; l < dim; l++) + (*new_ctrl)[j * dim + l] += alpha[index1] * ctrl[index2 + l]; + } + } + *ncontrol = (GLint) m; + return GLU_NO_ERROR; +} + +static GLint +calc_factor(GLfloat * pts, GLint order, GLint indx, GLint stride, + GLfloat tolerance, GLint dim) +{ + GLdouble model[16], proj[16]; + GLint viewport[4]; + GLdouble x, y, z, w, winx1, winy1, winz, winx2, winy2; + GLint i; + GLdouble len, dx, dy; + + glGetDoublev(GL_MODELVIEW_MATRIX, model); + glGetDoublev(GL_PROJECTION_MATRIX, proj); + glGetIntegerv(GL_VIEWPORT, viewport); + if (dim == 4) { + w = (GLdouble) pts[indx + 3]; + x = (GLdouble) pts[indx] / w; + y = (GLdouble) pts[indx + 1] / w; + z = (GLdouble) pts[indx + 2] / w; + gluProject(x, y, z, model, proj, viewport, &winx1, &winy1, &winz); + len = 0.0; + for (i = 1; i < order; i++) { + w = (GLdouble) pts[indx + i * stride + 3]; + x = (GLdouble) pts[indx + i * stride] / w; + y = (GLdouble) pts[indx + i * stride + 1] / w; + z = (GLdouble) pts[indx + i * stride + 2] / w; + if (gluProject + (x, y, z, model, proj, viewport, &winx2, &winy2, &winz)) { + dx = winx2 - winx1; + dy = winy2 - winy1; + len += sqrt(dx * dx + dy * dy); + } + winx1 = winx2; + winy1 = winy2; + } + } + else { + x = (GLdouble) pts[indx]; + y = (GLdouble) pts[indx + 1]; + if (dim == 2) + z = 0.0; + else + z = (GLdouble) pts[indx + 2]; + gluProject(x, y, z, model, proj, viewport, &winx1, &winy1, &winz); + len = 0.0; + for (i = 1; i < order; i++) { + x = (GLdouble) pts[indx + i * stride]; + y = (GLdouble) pts[indx + i * stride + 1]; + if (dim == 2) + z = 0.0; + else + z = (GLdouble) pts[indx + i * stride + 2]; + if (gluProject + (x, y, z, model, proj, viewport, &winx2, &winy2, &winz)) { + dx = winx2 - winx1; + dy = winy2 - winy1; + len += sqrt(dx * dx + dy * dy); + } + winx1 = winx2; + winy1 = winy2; + } + } + len /= tolerance; + return ((GLint) len + 1); +} + +/* we can't use the Mesa evaluators - no way to get the point coords */ +/* so we use our own Bezier point calculus routines */ +/* because I'm lazy, I reuse the ones from eval.c */ + +static void +bezier_curve(GLfloat * cp, GLfloat * out, GLfloat t, + GLuint dim, GLuint order, GLint offset) +{ + GLfloat s, powert; + GLuint i, k, bincoeff; + + if (order >= 2) { + bincoeff = order - 1; + s = 1.0 - t; + + for (k = 0; k < dim; k++) + out[k] = s * cp[k] + bincoeff * t * cp[offset + k]; + + for (i = 2, cp += 2 * offset, powert = t * t; i < order; + i++, powert *= t, cp += offset) { + bincoeff *= order - i; + bincoeff /= i; + + for (k = 0; k < dim; k++) + out[k] = s * out[k] + bincoeff * powert * cp[k]; + } + } + else { /* order=1 -> constant curve */ + + for (k = 0; k < dim; k++) + out[k] = cp[k]; + } +} + +static GLint +calc_parametric_factor(GLfloat * pts, GLint order, GLint indx, GLint stride, + GLfloat tolerance, GLint dim) +{ + GLdouble model[16], proj[16]; + GLint viewport[4]; + GLdouble x, y, z, w, x1, y1, z1, x2, y2, z2, x3, y3, z3; + GLint i; + GLint P; + GLfloat bez_pt[4]; + GLdouble len = 0.0, tmp, z_med; + + P = 2 * (order + 2); + glGetDoublev(GL_MODELVIEW_MATRIX, model); + glGetDoublev(GL_PROJECTION_MATRIX, proj); + glGetIntegerv(GL_VIEWPORT, viewport); + z_med = (viewport[2] + viewport[3]) * 0.5; + switch (dim) { + case 4: + for (i = 1; i < P; i++) { + bezier_curve(pts + indx, bez_pt, (GLfloat) i / (GLfloat) P, 4, + order, stride); + w = (GLdouble) bez_pt[3]; + x = (GLdouble) bez_pt[0] / w; + y = (GLdouble) bez_pt[1] / w; + z = (GLdouble) bez_pt[2] / w; + gluProject(x, y, z, model, proj, viewport, &x3, &y3, &z3); + z3 *= z_med; + bezier_curve(pts + indx, bez_pt, (GLfloat) (i - 1) / (GLfloat) P, 4, + order, stride); + w = (GLdouble) bez_pt[3]; + x = (GLdouble) bez_pt[0] / w; + y = (GLdouble) bez_pt[1] / w; + z = (GLdouble) bez_pt[2] / w; + gluProject(x, y, z, model, proj, viewport, &x1, &y1, &z1); + z1 *= z_med; + bezier_curve(pts + indx, bez_pt, (GLfloat) (i + 1) / (GLfloat) P, 4, + order, stride); + w = (GLdouble) bez_pt[3]; + x = (GLdouble) bez_pt[0] / w; + y = (GLdouble) bez_pt[1] / w; + z = (GLdouble) bez_pt[2] / w; + gluProject(x, y, z, model, proj, viewport, &x2, &y2, &z2); + z2 *= z_med; + /* calc distance between point (x3,y3,z3) and line segment */ + /* <x1,y1,z1><x2,y2,z2> */ + x = x2 - x1; + y = y2 - y1; + z = z2 - z1; + tmp = sqrt(x * x + y * y + z * z); + x /= tmp; + y /= tmp; + z /= tmp; + tmp = x3 * x + y3 * y + z3 * z - x1 * x - y1 * y - z1 * z; + x = x1 + x * tmp - x3; + y = y1 + y * tmp - y3; + z = z1 + z * tmp - z3; + tmp = sqrt(x * x + y * y + z * z); + if (tmp > len) + len = tmp; + } + break; + case 3: + for (i = 1; i < P; i++) { + bezier_curve(pts + indx, bez_pt, (GLfloat) i / (GLfloat) P, 3, + order, stride); + x = (GLdouble) bez_pt[0]; + y = (GLdouble) bez_pt[1]; + z = (GLdouble) bez_pt[2]; + gluProject(x, y, z, model, proj, viewport, &x3, &y3, &z3); + z3 *= z_med; + bezier_curve(pts + indx, bez_pt, (GLfloat) (i - 1) / (GLfloat) P, 3, + order, stride); + x = (GLdouble) bez_pt[0]; + y = (GLdouble) bez_pt[1]; + z = (GLdouble) bez_pt[2]; + gluProject(x, y, z, model, proj, viewport, &x1, &y1, &z1); + z1 *= z_med; + bezier_curve(pts + indx, bez_pt, (GLfloat) (i + 1) / (GLfloat) P, 3, + order, stride); + x = (GLdouble) bez_pt[0]; + y = (GLdouble) bez_pt[1]; + z = (GLdouble) bez_pt[2]; + gluProject(x, y, z, model, proj, viewport, &x2, &y2, &z2); + z2 *= z_med; + /* calc distance between point (x3,y3,z3) and line segment */ + /* <x1,y1,z1><x2,y2,z2> */ + x = x2 - x1; + y = y2 - y1; + z = z2 - z1; + tmp = sqrt(x * x + y * y + z * z); + x /= tmp; + y /= tmp; + z /= tmp; + tmp = x3 * x + y3 * y + z3 * z - x1 * x - y1 * y - z1 * z; + x = x1 + x * tmp - x3; + y = y1 + y * tmp - y3; + z = z1 + z * tmp - z3; + tmp = sqrt(x * x + y * y + z * z); + if (tmp > len) + len = tmp; + } + break; + case 2: + for (i = 1; i < P; i++) { + bezier_curve(pts + indx, bez_pt, (GLfloat) i / (GLfloat) P, 2, + order, stride); + x = (GLdouble) bez_pt[0]; + y = (GLdouble) bez_pt[1]; + z = 0.0; + gluProject(x, y, z, model, proj, viewport, &x3, &y3, &z3); + z3 *= z_med; + bezier_curve(pts + indx, bez_pt, (GLfloat) (i - 1) / (GLfloat) P, 2, + order, stride); + x = (GLdouble) bez_pt[0]; + y = (GLdouble) bez_pt[1]; + z = 0.0; + gluProject(x, y, z, model, proj, viewport, &x1, &y1, &z1); + z1 *= z_med; + bezier_curve(pts + indx, bez_pt, (GLfloat) (i + 1) / (GLfloat) P, 2, + order, stride); + x = (GLdouble) bez_pt[0]; + y = (GLdouble) bez_pt[1]; + z = 0.0; + gluProject(x, y, z, model, proj, viewport, &x2, &y2, &z2); + z2 *= z_med; + /* calc distance between point (x3,y3,z3) and line segment */ + /* <x1,y1,z1><x2,y2,z2> */ + x = x2 - x1; + y = y2 - y1; + z = z2 - z1; + tmp = sqrt(x * x + y * y + z * z); + x /= tmp; + y /= tmp; + z /= tmp; + tmp = x3 * x + y3 * y + z3 * z - x1 * x - y1 * y - z1 * z; + x = x1 + x * tmp - x3; + y = y1 + y * tmp - y3; + z = z1 + z * tmp - z3; + tmp = sqrt(x * x + y * y + z * z); + if (tmp > len) + len = tmp; + } + break; + + } + if (len < tolerance) + return (order); + else + return (GLint) (sqrt(len / tolerance) * (order + 2) + 1); +} + +static GLenum +calc_sampling_3D(new_ctrl_type * new_ctrl, GLfloat tolerance, GLint dim, + GLint uorder, GLint vorder, GLint ** ufactors, + GLint ** vfactors) +{ + GLfloat *ctrl; + GLint tmp_factor1, tmp_factor2; + GLint ufactor_cnt, vfactor_cnt; + GLint offset1, offset2, offset3; + GLint i, j; + + ufactor_cnt = new_ctrl->s_bezier_cnt; + vfactor_cnt = new_ctrl->t_bezier_cnt; + if ((*ufactors = (GLint *) malloc(sizeof(GLint) * ufactor_cnt * 3)) + == NULL) { + return GLU_OUT_OF_MEMORY; + } + if ((*vfactors = (GLint *) malloc(sizeof(GLint) * vfactor_cnt * 3)) + == NULL) { + free(*ufactors); + return GLU_OUT_OF_MEMORY; + } + ctrl = new_ctrl->geom_ctrl; + offset1 = new_ctrl->geom_t_stride * vorder; + offset2 = new_ctrl->geom_s_stride * uorder; + for (j = 0; j < vfactor_cnt; j++) { + *(*vfactors + j * 3 + 1) = tmp_factor1 = calc_factor(ctrl, vorder, + j * offset1, dim, + tolerance, dim); + /* loop ufactor_cnt-1 times */ + for (i = 1; i < ufactor_cnt; i++) { + tmp_factor2 = calc_factor(ctrl, vorder, + j * offset1 + i * offset2, dim, tolerance, + dim); + if (tmp_factor2 > tmp_factor1) + tmp_factor1 = tmp_factor2; + } + /* last time for the opposite edge */ + *(*vfactors + j * 3 + 2) = tmp_factor2 = calc_factor(ctrl, vorder, + j * offset1 + + i * offset2 - + new_ctrl-> + geom_s_stride, dim, + tolerance, dim); + if (tmp_factor2 > tmp_factor1) + *(*vfactors + j * 3) = tmp_factor2; + else + *(*vfactors + j * 3) = tmp_factor1; + } + offset3 = new_ctrl->geom_s_stride; + offset2 = new_ctrl->geom_s_stride * uorder; + for (j = 0; j < ufactor_cnt; j++) { + *(*ufactors + j * 3 + 1) = tmp_factor1 = calc_factor(ctrl, uorder, + j * offset2, + offset3, tolerance, + dim); + /* loop vfactor_cnt-1 times */ + for (i = 1; i < vfactor_cnt; i++) { + tmp_factor2 = calc_factor(ctrl, uorder, + j * offset2 + i * offset1, offset3, + tolerance, dim); + if (tmp_factor2 > tmp_factor1) + tmp_factor1 = tmp_factor2; + } + /* last time for the opposite edge */ + *(*ufactors + j * 3 + 2) = tmp_factor2 = calc_factor(ctrl, uorder, + j * offset2 + + i * offset1 - + new_ctrl-> + geom_t_stride, + offset3, tolerance, + dim); + if (tmp_factor2 > tmp_factor1) + *(*ufactors + j * 3) = tmp_factor2; + else + *(*ufactors + j * 3) = tmp_factor1; + } + return GL_NO_ERROR; +} + +static GLenum +calc_sampling_param_3D(new_ctrl_type * new_ctrl, GLfloat tolerance, GLint dim, + GLint uorder, GLint vorder, GLint ** ufactors, + GLint ** vfactors) +{ + GLfloat *ctrl; + GLint tmp_factor1, tmp_factor2; + GLint ufactor_cnt, vfactor_cnt; + GLint offset1, offset2, offset3; + GLint i, j; + + ufactor_cnt = new_ctrl->s_bezier_cnt; + vfactor_cnt = new_ctrl->t_bezier_cnt; + if ((*ufactors = (GLint *) malloc(sizeof(GLint) * ufactor_cnt * 3)) + == NULL) { + return GLU_OUT_OF_MEMORY; + } + if ((*vfactors = (GLint *) malloc(sizeof(GLint) * vfactor_cnt * 3)) + == NULL) { + free(*ufactors); + return GLU_OUT_OF_MEMORY; + } + ctrl = new_ctrl->geom_ctrl; + offset1 = new_ctrl->geom_t_stride * vorder; + offset2 = new_ctrl->geom_s_stride * uorder; + for (j = 0; j < vfactor_cnt; j++) { + *(*vfactors + j * 3 + 1) = tmp_factor1 = + calc_parametric_factor(ctrl, vorder, j * offset1, dim, tolerance, + dim); + /* loop ufactor_cnt-1 times */ + for (i = 1; i < ufactor_cnt; i++) { + tmp_factor2 = calc_parametric_factor(ctrl, vorder, + j * offset1 + i * offset2, dim, + tolerance, dim); + if (tmp_factor2 > tmp_factor1) + tmp_factor1 = tmp_factor2; + } + /* last time for the opposite edge */ + *(*vfactors + j * 3 + 2) = tmp_factor2 = + calc_parametric_factor(ctrl, vorder, + j * offset1 + i * offset2 - + new_ctrl->geom_s_stride, dim, tolerance, dim); + if (tmp_factor2 > tmp_factor1) + *(*vfactors + j * 3) = tmp_factor2; + else + *(*vfactors + j * 3) = tmp_factor1; + } + offset3 = new_ctrl->geom_s_stride; + offset2 = new_ctrl->geom_s_stride * uorder; + for (j = 0; j < ufactor_cnt; j++) { + *(*ufactors + j * 3 + 1) = tmp_factor1 = + calc_parametric_factor(ctrl, uorder, j * offset2, offset3, tolerance, + dim); + /* loop vfactor_cnt-1 times */ + for (i = 1; i < vfactor_cnt; i++) { + tmp_factor2 = calc_parametric_factor(ctrl, uorder, + j * offset2 + i * offset1, + offset3, tolerance, dim); + if (tmp_factor2 > tmp_factor1) + tmp_factor1 = tmp_factor2; + } + /* last time for the opposite edge */ + *(*ufactors + j * 3 + 2) = tmp_factor2 = + calc_parametric_factor(ctrl, uorder, + j * offset2 + i * offset1 - + new_ctrl->geom_t_stride, offset3, tolerance, + dim); + if (tmp_factor2 > tmp_factor1) + *(*ufactors + j * 3) = tmp_factor2; + else + *(*ufactors + j * 3) = tmp_factor1; + } + return GL_NO_ERROR; +} + +static GLenum +calc_sampling_2D(GLfloat * ctrl, GLint cnt, GLint order, + GLfloat tolerance, GLint dim, GLint ** factors) +{ + GLint factor_cnt; + GLint tmp_factor; + GLint offset; + GLint i; + + factor_cnt = cnt / order; + if ((*factors = (GLint *) malloc(sizeof(GLint) * factor_cnt)) == NULL) { + return GLU_OUT_OF_MEMORY; + } + offset = order * dim; + for (i = 0; i < factor_cnt; i++) { + tmp_factor = calc_factor(ctrl, order, i * offset, dim, tolerance, dim); + if (tmp_factor == 0) + (*factors)[i] = 1; + else + (*factors)[i] = tmp_factor; + } + return GL_NO_ERROR; +} + +static void +set_sampling_and_culling(GLUnurbsObj * nobj) +{ + if (nobj->auto_load_matrix == GL_FALSE) { + GLint i; + GLfloat m[4]; + + glPushAttrib((GLbitfield) (GL_VIEWPORT_BIT | GL_TRANSFORM_BIT)); + for (i = 0; i < 4; i++) + m[i] = nobj->sampling_matrices.viewport[i]; + glViewport(m[0], m[1], m[2], m[3]); + glMatrixMode(GL_PROJECTION); + glPushMatrix(); + glLoadMatrixf(nobj->sampling_matrices.proj); + glMatrixMode(GL_MODELVIEW); + glPushMatrix(); + glLoadMatrixf(nobj->sampling_matrices.model); + } +} + +static void +revert_sampling_and_culling(GLUnurbsObj * nobj) +{ + if (nobj->auto_load_matrix == GL_FALSE) { + glMatrixMode(GL_MODELVIEW); + glPopMatrix(); + glMatrixMode(GL_PROJECTION); + glPopMatrix(); + glPopAttrib(); + } +} + +GLenum +glu_do_sampling_3D(GLUnurbsObj * nobj, new_ctrl_type * new_ctrl, + GLint ** sfactors, GLint ** tfactors) +{ + GLint dim; + GLenum err; + + *sfactors = NULL; + *tfactors = NULL; + dim = nobj->surface.geom.dim; + set_sampling_and_culling(nobj); + if ((err = calc_sampling_3D(new_ctrl, nobj->sampling_tolerance, dim, + nobj->surface.geom.sorder, + nobj->surface.geom.torder, sfactors, + tfactors)) == GLU_ERROR) { + revert_sampling_and_culling(nobj); + call_user_error(nobj, err); + return GLU_ERROR; + } + revert_sampling_and_culling(nobj); + return GLU_NO_ERROR; +} + +GLenum +glu_do_sampling_uv(GLUnurbsObj * nobj, new_ctrl_type * new_ctrl, + GLint ** sfactors, GLint ** tfactors) +{ + GLint s_cnt, t_cnt, i; + GLint u_steps, v_steps; + + s_cnt = new_ctrl->s_bezier_cnt; + t_cnt = new_ctrl->t_bezier_cnt; + *sfactors = NULL; + *tfactors = NULL; + if ((*sfactors = (GLint *) malloc(sizeof(GLint) * s_cnt * 3)) + == NULL) { + return GLU_OUT_OF_MEMORY; + } + if ((*tfactors = (GLint *) malloc(sizeof(GLint) * t_cnt * 3)) + == NULL) { + free(*sfactors); + return GLU_OUT_OF_MEMORY; + } + u_steps = nobj->u_step; + v_steps = nobj->v_step; + for (i = 0; i < s_cnt; i++) { + *(*sfactors + i * 3) = u_steps; + *(*sfactors + i * 3 + 1) = u_steps; + *(*sfactors + i * 3 + 2) = u_steps; + } + for (i = 0; i < t_cnt; i++) { + *(*tfactors + i * 3) = v_steps; + *(*tfactors + i * 3 + 1) = v_steps; + *(*tfactors + i * 3 + 2) = v_steps; + } + return GLU_NO_ERROR; +} + + +GLenum +glu_do_sampling_param_3D(GLUnurbsObj * nobj, new_ctrl_type * new_ctrl, + GLint ** sfactors, GLint ** tfactors) +{ + GLint dim; + GLenum err; + + *sfactors = NULL; + *tfactors = NULL; + dim = nobj->surface.geom.dim; + set_sampling_and_culling(nobj); + if ( + (err = + calc_sampling_param_3D(new_ctrl, nobj->parametric_tolerance, dim, + nobj->surface.geom.sorder, + nobj->surface.geom.torder, sfactors, + tfactors)) == GLU_ERROR) { + revert_sampling_and_culling(nobj); + call_user_error(nobj, err); + return GLU_ERROR; + } + revert_sampling_and_culling(nobj); + return GLU_NO_ERROR; +} + + +static GLenum +glu_do_sampling_2D(GLUnurbsObj * nobj, GLfloat * ctrl, GLint cnt, GLint order, + GLint dim, GLint ** factors) +{ + GLenum err; + + set_sampling_and_culling(nobj); + err = calc_sampling_2D(ctrl, cnt, order, nobj->sampling_tolerance, dim, + factors); + revert_sampling_and_culling(nobj); + return err; +} + + +static GLenum +glu_do_sampling_u(GLUnurbsObj * nobj, GLfloat * ctrl, GLint cnt, GLint order, + GLint dim, GLint ** factors) +{ + GLint i; + GLint u_steps; + + cnt /= order; + if ((*factors = (GLint *) malloc(sizeof(GLint) * cnt)) + == NULL) { + return GLU_OUT_OF_MEMORY; + } + u_steps = nobj->u_step; + for (i = 0; i < cnt; i++) + (*factors)[i] = u_steps; + return GLU_NO_ERROR; +} + + +static GLenum +glu_do_sampling_param_2D(GLUnurbsObj * nobj, GLfloat * ctrl, GLint cnt, + GLint order, GLint dim, GLint ** factors) +{ + GLint i; + GLint u_steps; + GLfloat tolerance; + + set_sampling_and_culling(nobj); + tolerance = nobj->parametric_tolerance; + cnt /= order; + if ((*factors = (GLint *) malloc(sizeof(GLint) * cnt)) + == NULL) { + revert_sampling_and_culling(nobj); + return GLU_OUT_OF_MEMORY; + } + u_steps = nobj->u_step; + for (i = 0; i < cnt; i++) { + (*factors)[i] = calc_parametric_factor(ctrl, order, 0, + dim, tolerance, dim); + + } + revert_sampling_and_culling(nobj); + return GLU_NO_ERROR; +} + +GLenum +glu_do_sampling_crv(GLUnurbsObj * nobj, GLfloat * ctrl, GLint cnt, + GLint order, GLint dim, GLint ** factors) +{ + GLenum err; + + *factors = NULL; + switch (nobj->sampling_method) { + case GLU_PATH_LENGTH: + if ((err = glu_do_sampling_2D(nobj, ctrl, cnt, order, dim, factors)) != + GLU_NO_ERROR) { + call_user_error(nobj, err); + return GLU_ERROR; + } + break; + case GLU_DOMAIN_DISTANCE: + if ((err = glu_do_sampling_u(nobj, ctrl, cnt, order, dim, factors)) != + GLU_NO_ERROR) { + call_user_error(nobj, err); + return GLU_ERROR; + } + break; + case GLU_PARAMETRIC_ERROR: + if ( + (err = + glu_do_sampling_param_2D(nobj, ctrl, cnt, order, dim, + factors)) != GLU_NO_ERROR) { + call_user_error(nobj, err); + return GLU_ERROR; + } + break; + default: + abort(); + } + + return GLU_NO_ERROR; +} + +/* TODO - i don't like this culling - this one just tests if at least one */ +/* ctrl point lies within the viewport . Also the point_in_viewport() */ +/* should be included in the fnctions for efficiency reasons */ + +static GLboolean +point_in_viewport(GLfloat * pt, GLint dim) +{ + GLdouble model[16], proj[16]; + GLint viewport[4]; + GLdouble x, y, z, w, winx, winy, winz; + + glGetDoublev(GL_MODELVIEW_MATRIX, model); + glGetDoublev(GL_PROJECTION_MATRIX, proj); + glGetIntegerv(GL_VIEWPORT, viewport); + if (dim == 3) { + x = (GLdouble) pt[0]; + y = (GLdouble) pt[1]; + z = (GLdouble) pt[2]; + gluProject(x, y, z, model, proj, viewport, &winx, &winy, &winz); + } + else { + w = (GLdouble) pt[3]; + x = (GLdouble) pt[0] / w; + y = (GLdouble) pt[1] / w; + z = (GLdouble) pt[2] / w; + gluProject(x, y, z, model, proj, viewport, &winx, &winy, &winz); + } + if ((GLint) winx >= viewport[0] && (GLint) winx < viewport[2] && + (GLint) winy >= viewport[1] && (GLint) winy < viewport[3]) + return GL_TRUE; + return GL_FALSE; +} + +GLboolean +fine_culling_test_3D(GLUnurbsObj * nobj, GLfloat * pts, GLint s_cnt, + GLint t_cnt, GLint s_stride, GLint t_stride, GLint dim) +{ + GLint i, j; + + if (nobj->culling == GL_FALSE) + return GL_FALSE; + set_sampling_and_culling(nobj); + + if (dim == 3) { + for (i = 0; i < s_cnt; i++) + for (j = 0; j < t_cnt; j++) + if (point_in_viewport(pts + i * s_stride + j * t_stride, dim)) { + revert_sampling_and_culling(nobj); + return GL_FALSE; + } + } + else { + for (i = 0; i < s_cnt; i++) + for (j = 0; j < t_cnt; j++) + if (point_in_viewport(pts + i * s_stride + j * t_stride, dim)) { + revert_sampling_and_culling(nobj); + return GL_FALSE; + } + } + revert_sampling_and_culling(nobj); + return GL_TRUE; +} + +/*GLboolean +fine_culling_test_3D(GLUnurbsObj *nobj,GLfloat *pts,GLint s_cnt,GLint t_cnt, + GLint s_stride,GLint t_stride, GLint dim) +{ + GLint visible_cnt; + GLfloat feedback_buffer[5]; + GLsizei buffer_size; + GLint i,j; + + if(nobj->culling==GL_FALSE) + return GL_FALSE; + buffer_size=5; + set_sampling_and_culling(nobj); + + glFeedbackBuffer(buffer_size,GL_2D,feedback_buffer); + glRenderMode(GL_FEEDBACK); + if(dim==3) + { + for(i=0;i<s_cnt;i++) + { + glBegin(GL_LINE_LOOP); + for(j=0;j<t_cnt;j++) + glVertex3fv(pts+i*s_stride+j*t_stride); + glEnd(); + } + for(j=0;j<t_cnt;j++) + { + glBegin(GL_LINE_LOOP); + for(i=0;i<s_cnt;i++) + glVertex3fv(pts+i*s_stride+j*t_stride); + glEnd(); + } + } + else + { + for(i=0;i<s_cnt;i++) + { + glBegin(GL_LINE_LOOP); + for(j=0;j<t_cnt;j++) + glVertex4fv(pts+i*s_stride+j*t_stride); + glEnd(); + } + for(j=0;j<t_cnt;j++) + { + glBegin(GL_LINE_LOOP); + for(i=0;i<s_cnt;i++) + glVertex4fv(pts+i*s_stride+j*t_stride); + glEnd(); + } + } + visible_cnt=glRenderMode(GL_RENDER); + + revert_sampling_and_culling(nobj); + return (GLboolean)(visible_cnt==0); +}*/ + +GLboolean +fine_culling_test_2D(GLUnurbsObj * nobj, GLfloat * pts, GLint cnt, + GLint stride, GLint dim) +{ + GLint i; + + if (nobj->culling == GL_FALSE) + return GL_FALSE; + set_sampling_and_culling(nobj); + + if (dim == 3) { + for (i = 0; i < cnt; i++) + if (point_in_viewport(pts + i * stride, dim)) { + revert_sampling_and_culling(nobj); + return GL_FALSE; + } + } + else { + for (i = 0; i < cnt; i++) + if (point_in_viewport(pts + i * stride, dim)) { + revert_sampling_and_culling(nobj); + return GL_FALSE; + } + } + revert_sampling_and_culling(nobj); + return GL_TRUE; +} + +/*GLboolean +fine_culling_test_2D(GLUnurbsObj *nobj,GLfloat *pts,GLint cnt, + GLint stride, GLint dim) +{ + GLint visible_cnt; + GLfloat feedback_buffer[5]; + GLsizei buffer_size; + GLint i; + + if(nobj->culling==GL_FALSE) + return GL_FALSE; + buffer_size=5; + set_sampling_and_culling(nobj); + + glFeedbackBuffer(buffer_size,GL_2D,feedback_buffer); + glRenderMode(GL_FEEDBACK); + glBegin(GL_LINE_LOOP); + if(dim==3) + { + for(i=0;i<cnt;i++) + glVertex3fv(pts+i*stride); + } + else + { + for(i=0;i<cnt;i++) + glVertex4fv(pts+i*stride); + } + glEnd(); + visible_cnt=glRenderMode(GL_RENDER); + + revert_sampling_and_culling(nobj); + return (GLboolean)(visible_cnt==0); +}*/ |