1 /* 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: 9 ** 10 ** http://oss.sgi.com/projects/FreeB 11 ** 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. 17 ** 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. 23 ** 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. 33 */ 34 35 /* 36 * arctessellator.c++ 37 * 38 */ 39 40 //#include "glimports.h" 41 //#include "mystdio.h" 42 //#include "myassert.h" 43 #include "arctess.h" 44 //#include "bufpool.h" 45 #include "simplemath.h" 46 #include "bezierarc.h" 47 //#include "trimvertex.h" 48 #include "trimvertpool.h" 49 50 #define NOELIMINATION 51 52 #define steps_function(large, small, rate) (max(1, 1+ (int) ((large-small)/rate))); 53 54 /*----------------------------------------------------------------------------- 55 * ArcTessellator - construct an ArcTessellator 56 *----------------------------------------------------------------------------- 57 */ 58 59 ArcTessellator::ArcTessellator( TrimVertexPool& t, Pool& p ) 60 : pwlarcpool(p), trimvertexpool(t) 61 { 62 } 63 64 /*----------------------------------------------------------------------------- 65 * ~ArcTessellator - destroy an ArcTessellator 66 *----------------------------------------------------------------------------- 67 */ 68 69 ArcTessellator::~ArcTessellator( void ) 70 { 71 } 72 73 /*----------------------------------------------------------------------------- 74 * bezier - construct a bezier arc and attach it to an Arc 75 *----------------------------------------------------------------------------- 76 */ 77 78 void 79 ArcTessellator::bezier( Arc *arc, REAL s1, REAL s2, REAL t1, REAL t2 ) 80 { 81 assert( arc != 0 ); 82 assert( ! arc->isTessellated() ); 83 84 #ifndef NDEBUG 85 switch( arc->getside() ) { 86 case arc_left: 87 assert( s1 == s2 ); 88 assert( t2 < t1 ); 89 break; 90 case arc_right: 91 assert( s1 == s2 ); 92 assert( t1 < t2 ); 93 break; 94 case arc_top: 95 assert( t1 == t2 ); 96 assert( s2 < s1 ); 97 break; 98 case arc_bottom: 99 assert( t1 == t2 ); 100 assert( s1 < s2 ); 101 break; 102 case arc_none: 103 (void) abort(); 104 break; 105 } 106 #endif 107 108 TrimVertex *p = trimvertexpool.get(2); 109 arc->pwlArc = new(pwlarcpool) PwlArc( 2, p ); 110 p[0].param[0] = s1; 111 p[0].param[1] = t1; 112 p[1].param[0] = s2; 113 p[1].param[1] = t2; 114 assert( (s1 == s2) || (t1 == t2) ); 115 arc->setbezier(); 116 } 117 118 119 /*----------------------------------------------------------------------------- 120 * pwl_left - construct a left boundary pwl arc and attach it to an arc 121 *----------------------------------------------------------------------------- 122 */ 123 124 void 125 ArcTessellator::pwl_left( Arc *arc, REAL s, REAL t1, REAL t2, REAL rate ) 126 { 127 assert( t2 < t1 ); 128 129 /* if(rate <= 0.06) rate = 0.06;*/ 130 /* int nsteps = 1 + (int) ((t1 - t2) / rate ); */ 131 int nsteps = steps_function(t1, t2, rate); 132 133 134 REAL stepsize = (t1 - t2) / (REAL) nsteps; 135 136 TrimVertex *newvert = trimvertexpool.get( nsteps+1 ); 137 int i; 138 for( i = nsteps; i > 0; i-- ) { 139 newvert[i].param[0] = s; 140 newvert[i].param[1] = t2; 141 t2 += stepsize; 142 } 143 newvert[i].param[0] = s; 144 newvert[i].param[1] = t1; 145 146 arc->makeSide( new(pwlarcpool) PwlArc( nsteps+1, newvert ), arc_left ); 147 } 148 149 /*----------------------------------------------------------------------------- 150 * pwl_right - construct a right boundary pwl arc and attach it to an arc 151 *----------------------------------------------------------------------------- 152 */ 153 154 void 155 ArcTessellator::pwl_right( Arc *arc, REAL s, REAL t1, REAL t2, REAL rate ) 156 { 157 assert( t1 < t2 ); 158 159 /* if(rate <= 0.06) rate = 0.06;*/ 160 161 /* int nsteps = 1 + (int) ((t2 - t1) / rate ); */ 162 int nsteps = steps_function(t2,t1,rate); 163 REAL stepsize = (t2 - t1) / (REAL) nsteps; 164 165 TrimVertex *newvert = trimvertexpool.get( nsteps+1 ); 166 int i; 167 for( i = 0; i < nsteps; i++ ) { 168 newvert[i].param[0] = s; 169 newvert[i].param[1] = t1; 170 t1 += stepsize; 171 } 172 newvert[i].param[0] = s; 173 newvert[i].param[1] = t2; 174 175 arc->makeSide( new(pwlarcpool) PwlArc( nsteps+1, newvert ), arc_right ); 176 } 177 178 179 /*----------------------------------------------------------------------------- 180 * pwl_top - construct a top boundary pwl arc and attach it to an arc 181 *----------------------------------------------------------------------------- 182 */ 183 184 void 185 ArcTessellator::pwl_top( Arc *arc, REAL t, REAL s1, REAL s2, REAL rate ) 186 { 187 assert( s2 < s1 ); 188 189 /* if(rate <= 0.06) rate = 0.06;*/ 190 191 /* int nsteps = 1 + (int) ((s1 - s2) / rate ); */ 192 int nsteps = steps_function(s1,s2,rate); 193 REAL stepsize = (s1 - s2) / (REAL) nsteps; 194 195 TrimVertex *newvert = trimvertexpool.get( nsteps+1 ); 196 int i; 197 for( i = nsteps; i > 0; i-- ) { 198 newvert[i].param[0] = s2; 199 newvert[i].param[1] = t; 200 s2 += stepsize; 201 } 202 newvert[i].param[0] = s1; 203 newvert[i].param[1] = t; 204 205 arc->makeSide( new(pwlarcpool) PwlArc( nsteps+1, newvert ), arc_top ); 206 } 207 208 /*----------------------------------------------------------------------------- 209 * pwl_bottom - construct a bottom boundary pwl arc and attach it to an arc 210 *----------------------------------------------------------------------------- 211 */ 212 213 void 214 ArcTessellator::pwl_bottom( Arc *arc, REAL t, REAL s1, REAL s2, REAL rate ) 215 { 216 assert( s1 < s2 ); 217 218 /* if(rate <= 0.06) rate = 0.06;*/ 219 220 /* int nsteps = 1 + (int) ((s2 - s1) / rate ); */ 221 int nsteps = steps_function(s2,s1,rate); 222 REAL stepsize = (s2 - s1) / (REAL) nsteps; 223 224 TrimVertex *newvert = trimvertexpool.get( nsteps+1 ); 225 int i; 226 for( i = 0; i < nsteps; i++ ) { 227 newvert[i].param[0] = s1; 228 newvert[i].param[1] = t; 229 s1 += stepsize; 230 } 231 newvert[i].param[0] = s2; 232 newvert[i].param[1] = t; 233 234 arc->makeSide( new(pwlarcpool) PwlArc( nsteps+1, newvert ), arc_bottom ); 235 } 236 237 /*----------------------------------------------------------------------------- 238 * pwl - construct a pwl arc and attach it to an arc 239 *----------------------------------------------------------------------------- 240 */ 241 242 void 243 ArcTessellator::pwl( Arc *arc, REAL s1, REAL s2, REAL t1, REAL t2, REAL rate ) 244 { 245 246 /* if(rate <= 0.06) rate = 0.06;*/ 247 248 int snsteps = 1 + (int) (glu_abs(s2 - s1) / rate ); 249 int tnsteps = 1 + (int) (glu_abs(t2 - t1) / rate ); 250 int nsteps = max(1,max( snsteps, tnsteps )); 251 252 REAL sstepsize = (s2 - s1) / (REAL) nsteps; 253 REAL tstepsize = (t2 - t1) / (REAL) nsteps; 254 TrimVertex *newvert = trimvertexpool.get( nsteps+1 ); 255 long i; 256 for( i = 0; i < nsteps; i++ ) { 257 newvert[i].param[0] = s1; 258 newvert[i].param[1] = t1; 259 s1 += sstepsize; 260 t1 += tstepsize; 261 } 262 newvert[i].param[0] = s2; 263 newvert[i].param[1] = t2; 264 265 /* arc->makeSide( new(pwlarcpool) PwlArc( nsteps+1, newvert ), arc_bottom ); */ 266 arc->pwlArc = new(pwlarcpool) PwlArc( nsteps+1, newvert ); 267 268 arc->clearbezier(); 269 arc->clearside( ); 270 } 271 272 273 /*----------------------------------------------------------------------------- 274 * tessellateLinear - constuct a linear pwl arc and attach it to an Arc 275 *----------------------------------------------------------------------------- 276 */ 277 278 void 279 ArcTessellator::tessellateLinear( Arc *arc, REAL geo_stepsize, REAL arc_stepsize, int isrational ) 280 { 281 assert( arc->pwlArc == NULL ); 282 REAL s1, s2, t1, t2; 283 284 //we don't need to scale by arc_stepsize if the trim curve 285 //is piecewise linear. Reason: In pwl_right, pwl_left, pwl_top, pwl_left, 286 //and pwl, the nsteps is computed by deltaU (or V) /stepsize. 287 //The quantity deltaU/arc_stepsize doesn't have any meaning. And 288 //it causes problems: see bug 517641 289 REAL stepsize = geo_stepsize; /* * arc_stepsize*/; 290 291 BezierArc *b = arc->bezierArc; 292 293 if( isrational ) { 294 s1 = b->cpts[0] / b->cpts[2]; 295 t1 = b->cpts[1] / b->cpts[2]; 296 s2 = b->cpts[b->stride+0] / b->cpts[b->stride+2]; 297 t2 = b->cpts[b->stride+1] / b->cpts[b->stride+2]; 298 } else { 299 s1 = b->cpts[0]; 300 t1 = b->cpts[1]; 301 s2 = b->cpts[b->stride+0]; 302 t2 = b->cpts[b->stride+1]; 303 } 304 if( s1 == s2 ) 305 if( t1 < t2 ) 306 pwl_right( arc, s1, t1, t2, stepsize ); 307 else 308 pwl_left( arc, s1, t1, t2, stepsize ); 309 else if( t1 == t2 ) 310 if( s1 < s2 ) 311 pwl_bottom( arc, t1, s1, s2, stepsize ); 312 else 313 pwl_top( arc, t1, s1, s2, stepsize ); 314 else 315 pwl( arc, s1, s2, t1, t2, stepsize ); 316 } 317 318 /*----------------------------------------------------------------------------- 319 * tessellateNonlinear - constuct a nonlinear pwl arc and attach it to an Arc 320 *----------------------------------------------------------------------------- 321 */ 322 323 void 324 ArcTessellator::tessellateNonlinear( Arc *arc, REAL geo_stepsize, REAL arc_stepsize, int isrational ) 325 { 326 assert( arc->pwlArc == NULL ); 327 328 REAL stepsize = geo_stepsize * arc_stepsize; 329 330 BezierArc *bezierArc = arc->bezierArc; 331 332 REAL size; //bounding box size of the curve in UV 333 { 334 int i,j; 335 REAL min_u, min_v, max_u,max_v; 336 min_u = max_u = bezierArc->cpts[0]; 337 min_v = max_v = bezierArc->cpts[1]; 338 for(i=1, j=bezierArc->stride; i<bezierArc->order; i++, j+= bezierArc->stride) 339 { 340 if(bezierArc->cpts[j] < min_u) 341 min_u = bezierArc->cpts[j]; 342 if(bezierArc->cpts[j] > max_u) 343 max_u = bezierArc->cpts[j]; 344 if(bezierArc->cpts[j+1] < min_v) 345 min_v = bezierArc->cpts[j+1]; 346 if(bezierArc->cpts[j+1] > max_v) 347 max_v = bezierArc->cpts[j+1]; 348 } 349 350 size = max_u - min_u; 351 if(size < max_v - min_v) 352 size = max_v - min_v; 353 } 354 355 /*int nsteps = 1 + (int) (1.0/stepsize);*/ 356 357 int nsteps = (int) (size/stepsize); 358 if(nsteps <=0) 359 nsteps=1; 360 361 TrimVertex *vert = trimvertexpool.get( nsteps+1 ); 362 REAL dp = 1.0/nsteps; 363 364 365 arc->pwlArc = new(pwlarcpool) PwlArc(); 366 arc->pwlArc->pts = vert; 367 368 if( isrational ) { 369 REAL pow_u[MAXORDER], pow_v[MAXORDER], pow_w[MAXORDER]; 370 trim_power_coeffs( bezierArc, pow_u, 0 ); 371 trim_power_coeffs( bezierArc, pow_v, 1 ); 372 trim_power_coeffs( bezierArc, pow_w, 2 ); 373 374 /* compute first point exactly */ 375 REAL *b = bezierArc->cpts; 376 vert->param[0] = b[0]/b[2]; 377 vert->param[1] = b[1]/b[2]; 378 379 /* strength reduction on p = dp * step would introduce error */ 380 int step; 381 #ifndef NOELIMINATION 382 int ocanremove = 0; 383 #endif 384 long order = bezierArc->order; 385 for( step=1, ++vert; step<nsteps; step++, vert++ ) { 386 REAL p = dp * step; 387 REAL u = pow_u[0]; 388 REAL v = pow_v[0]; 389 REAL w = pow_w[0]; 390 for( int i = 1; i < order; i++ ) { 391 u = u * p + pow_u[i]; 392 v = v * p + pow_v[i]; 393 w = w * p + pow_w[i]; 394 } 395 vert->param[0] = u/w; 396 vert->param[1] = v/w; 397 #ifndef NOELIMINATION 398 REAL ds = glu_abs(vert[0].param[0] - vert[-1].param[0]); 399 REAL dt = glu_abs(vert[0].param[1] - vert[-1].param[1]); 400 int canremove = (ds<geo_stepsize && dt<geo_stepsize) ? 1 : 0; 401 REAL ods=0.0, odt=0.0; 402 403 if( ocanremove && canremove ) { 404 REAL nds = ds + ods; 405 REAL ndt = dt + odt; 406 if( nds<geo_stepsize && ndt<geo_stepsize ) { 407 // remove previous point 408 --vert; 409 vert[0].param[0] = vert[1].param[0]; 410 vert[0].param[1] = vert[1].param[1]; 411 ods = nds; 412 odt = ndt; 413 ocanremove = 1; 414 } else { 415 ocanremove = canremove; 416 ods = ds; 417 odt = dt; 418 } 419 } else { 420 ocanremove = canremove; 421 ods = ds; 422 odt = dt; 423 } 424 #endif 425 } 426 427 /* compute last point exactly */ 428 b += (order - 1) * bezierArc->stride; 429 vert->param[0] = b[0]/b[2]; 430 vert->param[1] = b[1]/b[2]; 431 432 } else { 433 REAL pow_u[MAXORDER], pow_v[MAXORDER]; 434 trim_power_coeffs( bezierArc, pow_u, 0 ); 435 trim_power_coeffs( bezierArc, pow_v, 1 ); 436 437 /* compute first point exactly */ 438 REAL *b = bezierArc->cpts; 439 vert->param[0] = b[0]; 440 vert->param[1] = b[1]; 441 442 /* strength reduction on p = dp * step would introduce error */ 443 int step; 444 #ifndef NOELIMINATION 445 int ocanremove = 0; 446 #endif 447 long order = bezierArc->order; 448 for( step=1, ++vert; step<nsteps; step++, vert++ ) { 449 REAL p = dp * step; 450 REAL u = pow_u[0]; 451 REAL v = pow_v[0]; 452 for( int i = 1; i < bezierArc->order; i++ ) { 453 u = u * p + pow_u[i]; 454 v = v * p + pow_v[i]; 455 } 456 vert->param[0] = u; 457 vert->param[1] = v; 458 #ifndef NOELIMINATION 459 REAL ds = glu_abs(vert[0].param[0] - vert[-1].param[0]); 460 REAL dt = glu_abs(vert[0].param[1] - vert[-1].param[1]); 461 int canremove = (ds<geo_stepsize && dt<geo_stepsize) ? 1 : 0; 462 REAL ods=0.0, odt=0.0; 463 464 if( ocanremove && canremove ) { 465 REAL nds = ds + ods; 466 REAL ndt = dt + odt; 467 if( nds<geo_stepsize && ndt<geo_stepsize ) { 468 // remove previous point 469 --vert; 470 vert[0].param[0] = vert[1].param[0]; 471 vert[0].param[1] = vert[1].param[1]; 472 ods = nds; 473 odt = ndt; 474 ocanremove = 1; 475 } else { 476 ocanremove = canremove; 477 ods = ds; 478 odt = dt; 479 } 480 } else { 481 ocanremove = canremove; 482 ods = ds; 483 odt = dt; 484 } 485 #endif 486 } 487 488 /* compute last point exactly */ 489 b += (order - 1) * bezierArc->stride; 490 vert->param[0] = b[0]; 491 vert->param[1] = b[1]; 492 } 493 arc->pwlArc->npts = vert - arc->pwlArc->pts + 1; 494 /* 495 for( TrimVertex *vt=pwlArc->pts; vt != vert-1; vt++ ) { 496 if( tooclose( vt[0].param[0], vt[1].param[0] ) ) 497 vt[1].param[0] = vt[0].param[0]; 498 if( tooclose( vt[0].param[1], vt[1].param[1] ) ) 499 vt[1].param[1] = vt[0].param[1]; 500 } 501 */ 502 } 503 504 const REAL ArcTessellator::gl_Bernstein[][MAXORDER][MAXORDER] = { 505 { 506 {1, 0, 0, 0, 0, 0, 0, 0 }, 507 {0, 0, 0, 0, 0, 0, 0, 0 }, 508 {0, 0, 0, 0, 0, 0, 0, 0 }, 509 {0, 0, 0, 0, 0, 0, 0, 0 }, 510 {0, 0, 0, 0, 0, 0, 0, 0 }, 511 {0, 0, 0, 0, 0, 0, 0, 0 }, 512 {0, 0, 0, 0, 0, 0, 0, 0 }, 513 {0, 0, 0, 0, 0, 0, 0, 0 } 514 }, 515 { 516 {-1, 1, 0, 0, 0, 0, 0, 0 }, 517 {1, 0, 0, 0, 0, 0, 0, 0 }, 518 {0, 0, 0, 0, 0, 0, 0, 0 }, 519 {0, 0, 0, 0, 0, 0, 0, 0 }, 520 {0, 0, 0, 0, 0, 0, 0, 0 }, 521 {0, 0, 0, 0, 0, 0, 0, 0 }, 522 {0, 0, 0, 0, 0, 0, 0, 0 }, 523 {0, 0, 0, 0, 0, 0, 0, 0 } 524 }, 525 { 526 {1, -2, 1, 0, 0, 0, 0, 0 }, 527 {-2, 2, 0, 0, 0, 0, 0, 0 }, 528 {1, 0, 0, 0, 0, 0, 0, 0 }, 529 {0, 0, 0, 0, 0, 0, 0, 0 }, 530 {0, 0, 0, 0, 0, 0, 0, 0 }, 531 {0, 0, 0, 0, 0, 0, 0, 0 }, 532 {0, 0, 0, 0, 0, 0, 0, 0 }, 533 {0, 0, 0, 0, 0, 0, 0, 0 } 534 }, 535 { 536 {-1, 3, -3, 1, 0, 0, 0, 0 }, 537 {3, -6, 3, 0, 0, 0, 0, 0 }, 538 {-3, 3, 0, 0, 0, 0, 0, 0 }, 539 {1, 0, 0, 0, 0, 0, 0, 0 }, 540 {0, 0, 0, 0, 0, 0, 0, 0 }, 541 {0, 0, 0, 0, 0, 0, 0, 0 }, 542 {0, 0, 0, 0, 0, 0, 0, 0 }, 543 {0, 0, 0, 0, 0, 0, 0, 0 } 544 }, 545 { 546 {1, -4, 6, -4, 1, 0, 0, 0 }, 547 {-4, 12, -12, 4, 0, 0, 0, 0 }, 548 {6, -12, 6, 0, 0, 0, 0, 0 }, 549 {-4, 4, 0, 0, 0, 0, 0, 0 }, 550 {1, 0, 0, 0, 0, 0, 0, 0 }, 551 {0, 0, 0, 0, 0, 0, 0, 0 }, 552 {0, 0, 0, 0, 0, 0, 0, 0 }, 553 {0, 0, 0, 0, 0, 0, 0, 0 } 554 }, 555 { 556 {-1, 5, -10, 10, -5, 1, 0, 0 }, 557 {5, -20, 30, -20, 5, 0, 0, 0 }, 558 {-10, 30, -30, 10, 0, 0, 0, 0 }, 559 {10, -20, 10, 0, 0, 0, 0, 0 }, 560 {-5, 5, 0, 0, 0, 0, 0, 0 }, 561 {1, 0, 0, 0, 0, 0, 0, 0 }, 562 {0, 0, 0, 0, 0, 0, 0, 0 }, 563 {0, 0, 0, 0, 0, 0, 0, 0 } 564 }, 565 { 566 {1, -6, 15, -20, 15, -6, 1, 0 }, 567 {-6, 30, -60, 60, -30, 6, 0, 0 }, 568 {15, -60, 90, -60, 15, 0, 0, 0 }, 569 {-20, 60, -60, 20, 0, 0, 0, 0 }, 570 {15, -30, 15, 0, 0, 0, 0, 0 }, 571 {-6, 6, 0, 0, 0, 0, 0, 0 }, 572 {1, 0, 0, 0, 0, 0, 0, 0 }, 573 {0, 0, 0, 0, 0, 0, 0, 0 } 574 }, 575 { 576 {-1, 7, -21, 35, -35, 21, -7, 1 }, 577 {7, -42, 105, -140, 105, -42, 7, 0 }, 578 {-21, 105, -210, 210, -105, 21, 0, 0 }, 579 {35, -140, 210, -140, 35, 0, 0, 0 }, 580 {-35, 105, -105, 35, 0, 0, 0, 0 }, 581 {21, -42, 21, 0, 0, 0, 0, 0 }, 582 {-7, 7, 0, 0, 0, 0, 0, 0 }, 583 {1, 0, 0, 0, 0, 0, 0, 0 } 584 }}; 585 586 587 /*----------------------------------------------------------------------------- 588 * trim_power_coeffs - compute power basis coefficients from bezier coeffients 589 *----------------------------------------------------------------------------- 590 */ 591 void 592 ArcTessellator::trim_power_coeffs( BezierArc *bez_arc, REAL *p, int coord ) 593 { 594 int stride = bez_arc->stride; 595 int order = bez_arc->order; 596 REAL *base = bez_arc->cpts + coord; 597 598 REAL const (*mat)[MAXORDER][MAXORDER] = &gl_Bernstein[order-1]; 599 REAL const (*lrow)[MAXORDER] = &(*mat)[order]; 600 601 /* WIN32 didn't like the following line within the for-loop */ 602 REAL const (*row)[MAXORDER] = &(*mat)[0]; 603 for( ; row != lrow; row++ ) { 604 REAL s = 0.0; 605 REAL *point = base; 606 REAL const *mlast = *row + order; 607 for( REAL const *m = *row; m != mlast; m++, point += stride ) 608 s += *(m) * (*point); 609 *(p++) = s; 610 } 611 } 612