Edit
IABSD.fr/xenocara/lib/libGLU/src/libnurbs/internals/mapdescv.cc
jsg
- lib/libGLU/src/libnurbs/internals/mapdescv.cc
-
/* ** License Applicability. Except to the extent portions of this file are ** made subject to an alternative license as permitted in the SGI Free ** Software License B, Version 1.1 (the "License"), the contents of this ** file are subject only to the provisions of the License. You may not use ** this file except in compliance with the License. You may obtain a copy ** of the License at Silicon Graphics, Inc., attn: Legal Services, 1600 ** Amphitheatre Parkway, Mountain View, CA 94043-1351, or at: ** ** http://oss.sgi.com/projects/FreeB ** ** Note that, as provided in the License, the Software is distributed on an ** "AS IS" basis, with ALL EXPRESS AND IMPLIED WARRANTIES AND CONDITIONS ** DISCLAIMED, INCLUDING, WITHOUT LIMITATION, ANY IMPLIED WARRANTIES AND ** CONDITIONS OF MERCHANTABILITY, SATISFACTORY QUALITY, FITNESS FOR A ** PARTICULAR PURPOSE, AND NON-INFRINGEMENT. ** ** Original Code. The Original Code is: OpenGL Sample Implementation, ** Version 1.2.1, released January 26, 2000, developed by Silicon Graphics, ** Inc. The Original Code is Copyright (c) 1991-2000 Silicon Graphics, Inc. ** Copyright in any portions created by third parties is as indicated ** elsewhere herein. All Rights Reserved. ** ** Additional Notice Provisions: The application programming interfaces ** established by SGI in conjunction with the Original Code are The ** OpenGL(R) Graphics System: A Specification (Version 1.2.1), released ** April 1, 1999; The OpenGL(R) Graphics System Utility Library (Version ** 1.3), released November 4, 1998; and OpenGL(R) Graphics with the X ** Window System(R) (Version 1.3), released October 19, 1998. This software ** was created using the OpenGL(R) version 1.2.1 Sample Implementation ** published by SGI, but has not been independently verified as being ** compliant with the OpenGL(R) version 1.2.1 Specification. */ /* * mapdescv.c++ * */ #include "glimports.h" #include "mystdio.h" #include "myassert.h" #include "mystring.h" #include "mymath.h" #include "nurbsconsts.h" #include "mapdesc.h" /*-------------------------------------------------------------------------- * calcPartialVelocity - calculate maximum magnitude of a given partial * derivative *-------------------------------------------------------------------------- */ REAL Mapdesc::calcPartialVelocity ( REAL *p, int stride, int ncols, int partial, REAL range ) { REAL tmp[MAXORDER][MAXCOORDS]; REAL mag[MAXORDER]; assert( ncols <= MAXORDER ); int j, k, t; // copy inhomogeneous control points into temporary array for( j=0; j != ncols; j++ ) for( k=0; k != inhcoords; k++ ) tmp[j][k] = p[j*stride + k]; for( t=0; t != partial; t++ ) for( j=0; j != ncols-t-1; j++ ) for( k=0; k != inhcoords; k++ ) tmp[j][k] = tmp[j+1][k] - tmp[j][k]; // compute magnitude and store in mag array for( j=0; j != ncols-partial; j++ ) { mag[j] = 0.0; for( k=0; k != inhcoords; k++ ) mag[j] += tmp[j][k] * tmp[j][k]; } // compute scale factor REAL fac = 1; REAL invt = 1.0 / range; for( t = ncols-1; t != ncols-1-partial; t-- ) fac *= t * invt; // compute max magnitude of all entries in array REAL max = 0.0; for( j=0; j != ncols-partial; j++ ) if( mag[j] > max ) max = mag[j]; max = fac * sqrtf( (float) max ); return max; } /*-------------------------------------------------------------------------- * calcPartialVelocity - calculate maximum magnitude of a given partial * derivative *-------------------------------------------------------------------------- */ REAL Mapdesc::calcPartialVelocity ( REAL *dist, REAL *p, int rstride, int cstride, int nrows, int ncols, int spartial, int tpartial, REAL srange, REAL trange, int side ) { REAL tmp[MAXORDER][MAXORDER][MAXCOORDS]; REAL mag[MAXORDER][MAXORDER]; assert( nrows <= MAXORDER ); assert( ncols <= MAXORDER ); REAL *tp = &tmp[0][0][0]; REAL *mp = &mag[0][0]; const int istride = sizeof( tmp[0]) / sizeof( tmp[0][0][0] ); const int jstride = sizeof( tmp[0][0]) / sizeof( tmp[0][0][0] ); /* const int kstride = sizeof( tmp[0][0][0]) / sizeof( tmp[0][0][0] ); */ const int mistride = sizeof( mag[0]) / sizeof( mag[0][0] ); const int mjstride = sizeof( mag[0][0]) / sizeof( mag[0][0] ); const int idist = nrows * istride; const int jdist = ncols * jstride; /* const int kdist = inhcoords * kstride; */ const int id = idist - spartial * istride; const int jd = jdist - tpartial * jstride; { // copy control points REAL *ti = tp; REAL *qi = p; REAL *til = tp + idist; for( ; ti != til; ) { REAL *tj = ti; REAL *qj = qi; REAL *tjl = ti + jdist; for( ; tj != tjl; ) { for( int k=0; k != inhcoords; k++ ) { tj[k] = qj[k]; } tj += jstride; qj += cstride; } ti += istride; qi += rstride; } } { // compute (s)-partial derivative control points REAL *til = tp + idist - istride; const REAL *till = til - ( spartial * istride ); for( ; til != till; til -= istride ) for( REAL *ti = tp; ti != til; ti += istride ) for( REAL *tj = ti, *tjl = tj + jdist; tj != tjl; tj += jstride ) for( int k=0; k != inhcoords; k++ ) tj[k] = tj[k+istride] - tj[k]; } { // compute (s,t)-partial derivative control points REAL *tjl = tp + jdist - jstride; const REAL *tjll = tjl - ( tpartial * jstride ); for( ; tjl != tjll; tjl -= jstride ) for( REAL *tj = tp; tj != tjl; tj += jstride ) for( REAL *ti = tj, *til = ti + id; ti != til; ti += istride ) for( int k=0; k != inhcoords; k++ ) ti[k] = ti[k+jstride] - ti[k]; } REAL max = 0.0; { // compute magnitude and store in mag array memset( (void *) mp, 0, sizeof( mag ) ); for( REAL *ti = tp, *mi = mp, *til = tp + id; ti != til; ti += istride, mi += mistride ) for( REAL *tj = ti, *mj = mi, *tjl = ti + jd; tj != tjl; tj += jstride, mj += mjstride ) { for( int k=0; k != inhcoords; k++ ) *mj += tj[k] * tj[k]; if( *mj > max ) max = *mj; } } int i, j; // compute scale factor REAL fac = 1.0; { REAL invs = 1.0 / srange; REAL invt = 1.0 / trange; for( int s = nrows-1, slast = s-spartial; s != slast; s-- ) fac *= s * invs; for( int t = ncols-1, tlast = t-tpartial; t != tlast; t-- ) fac *= t * invt; } if( side == 0 ) { // compute max magnitude of first and last column dist[0] = 0.0; dist[1] = 0.0; for( i=0; i != nrows-spartial; i++ ) { j = 0; if( mag[i][j] > dist[0] ) dist[0] = mag[i][j]; j = ncols-tpartial-1; if( mag[i][j] > dist[1] ) dist[1] = mag[i][j]; } dist[0] = fac * sqrtf( dist[0] ); dist[1] = fac * sqrtf( dist[1] ); } else if( side == 1 ) { // compute max magnitude of first and last row dist[0] = 0.0; dist[1] = 0.0; for( j=0; j != ncols-tpartial; j++ ) { i = 0; if( mag[i][j] > dist[0] ) dist[0] = mag[i][j]; i = nrows-spartial-1; if( mag[i][j] > dist[1] ) dist[1] = mag[i][j]; } dist[0] = fac * sqrtf( dist[0] ); dist[1] = fac * sqrtf( dist[1] ); } max = fac * sqrtf( (float) max ); return max; }