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 * tobezier.c++
37 *
38 */
39
40#include "glimports.h"
41#include "myassert.h"
42#include "mystdio.h"
43#include "mystring.h"
44#include "quilt.h"
45#include "knotvector.h"
46
47/* local type definitions */
48struct Breakpt {		/* breakpoints	*/
49    Knot		value;		/* value	*/
50    int			multi;		/* multiplicity	*/
51    int			def;		/* deficit */
52};
53
54struct Knotspec {		/* knotvector format */
55    long		order;		/* order of spline  */
56    Knot_ptr		inkbegin;	/* input knot sequence */
57    Knot_ptr		inkend;		/* location after last knot */
58    Knot_ptr		outkbegin;	/* in-process knot subsequence */
59    Knot_ptr		outkend;	/* location after last knot */
60    Knot_ptr		kleft;		/* */
61    Knot_ptr		kright;		/* */
62    Knot_ptr		kfirst;		/* */
63    Knot_ptr		klast;		/* */
64    Knot_ptr		sbegin;		/* conversion factor values */
65    Breakpt *		bbegin;		/* in-process breakpoints */
66    Breakpt *		bend;		/* last breakpoint */
67    int			ncoords;	/* coordinates per control point */
68    int			prestride;	/* stride between input points */
69    int			poststride;	/* stride between output points	*/
70    int 		preoffset;	/* scaled point offset	*/
71    int 		postoffset;	/* scaled point offset	*/
72    int 		prewidth;	/* width of dimension	*/
73    int 		postwidth;	/* width of dimension	*/
74    int 		istransformed;	/* was dimension transformed */
75    Knotspec *		next;   	/* next knotspec */
76    Knotspec *		kspectotrans;   /* knotspec in transformation direction */
77
78			Knotspec( void );
79			~Knotspec( void );
80    void		factors( void );
81    void		insert( REAL * );
82    void		preselect();
83    void		select( void );
84    void		copy( INREAL *, REAL * );
85    void		breakpoints( void );
86    void		knots( void );
87    void		transform( REAL * );
88    void		showpts( REAL * );
89
90    void		pt_io_copy( REAL *, INREAL * );
91    void		pt_oo_copy( REAL *, REAL * );
92    void		pt_oo_sum( REAL*, REAL*, REAL*, Knot, Knot );
93};
94
95struct Splinespec {		/* a non-uniform tensor element */
96			Splinespec( int );
97                        ~Splinespec(void);
98    Knotspec		*kspec;	/* format of each param. dir. */
99    int			dim;		/* domain dimension */
100    REAL *		outcpts;	/* Bezier control points */
101
102    void		kspecinit( Knotvector & );
103    void		kspecinit( Knotvector &, Knotvector & );
104    void		select( void );
105    void		layout( long );
106    void		setupquilt( Quilt_ptr );
107    void		copy( INREAL * );
108    void		transform( void );
109};
110
111/*-----------------------------------------------------------------------------
112 * Quilt::toBezier - convert from NURBS to rational Bezier
113 *-----------------------------------------------------------------------------
114 */
115
116void
117Quilt::toBezier(
118    Knotvector& knotvector,	/* a knot vector */
119    INREAL *ctlpts,		/* input contol points */
120    long ncoords )		/* number of coordinates per control point */
121{
122    Splinespec spline( 1 );
123    spline.kspecinit( knotvector );
124    spline.select();
125    spline.layout( ncoords );
126    spline.setupquilt( this );
127    spline.copy( ctlpts );
128    spline.transform();
129}
130
131void
132Quilt::toBezier(
133    Knotvector& sknotvector,	/* a knot vector */
134    Knotvector& tknotvector,	/* a knot vector */
135    INREAL *ctlpts,		/* input contol points */
136    long ncoords )		/* number of coordinates per control point */
137{
138    Splinespec spline( 2 );
139    spline.kspecinit( sknotvector, tknotvector );
140    spline.select();
141    spline.layout( ncoords );
142    spline.setupquilt( this );
143    spline.copy( ctlpts );
144    spline.transform();
145}
146Splinespec::Splinespec( int dimen )
147{
148    dim = dimen;
149}
150
151Splinespec::~Splinespec( void )
152{
153    /* Note: do NOT delete 'outcpts' here since its address (not contents)
154     * is copied in 'cpts' in this file in function Splinespec::setupquilt().
155     * This block of memory will eventually be deleted in file quilt.c++ in
156     * function Quilt::deleteMe() through 'cpts' so do NOT delete it here!
157     */
158    Knotspec *ktrav= kspec;         //start at beginning of list
159    while (ktrav != 0) {            //any items to delete?
160       Knotspec *deleteThis= ktrav; //remember to delete this
161       ktrav= ktrav->next;          //go to next item if any
162       delete deleteThis;           //delete it
163    }
164} /* ~Splinespec() */
165
166/*-----------------------------------------------------------------------------
167 * Splinespec::kspecinit - initialize Splinespec structure
168 *
169 * Client: Quilt::toBezier
170 *-----------------------------------------------------------------------------
171 */
172
173void
174Splinespec::kspecinit( Knotvector& knotvector )
175{
176    kspec = new Knotspec;
177    kspec->inkbegin = knotvector.knotlist;
178    kspec->inkend = knotvector.knotlist + knotvector.knotcount;
179    kspec->prestride = (int) knotvector.stride;
180    kspec->order = knotvector.order;
181    kspec->next = NULL;
182}
183
184void
185Splinespec::kspecinit( Knotvector& sknotvector, Knotvector& tknotvector )
186{
187    kspec = new Knotspec;
188    Knotspec *tkspec = new Knotspec;
189
190    kspec->inkbegin = sknotvector.knotlist;
191    kspec->inkend = sknotvector.knotlist + sknotvector.knotcount;
192    kspec->prestride = (int) sknotvector.stride;
193    kspec->order = sknotvector.order;
194    kspec->next = tkspec;
195
196    tkspec->inkbegin = tknotvector.knotlist;
197    tkspec->inkend = tknotvector.knotlist + tknotvector.knotcount;
198    tkspec->prestride = (int) tknotvector.stride;
199    tkspec->order = tknotvector.order;
200    tkspec->next = NULL;
201}
202
203
204/*-----------------------------------------------------------------------------
205 * Splinespec::select - select the subsegments to copy
206 *
207 * Client: gl_quilt_to_bezier
208 *-----------------------------------------------------------------------------
209 */
210
211void
212Splinespec::select( )
213{
214    for( Knotspec *knotspec = kspec; knotspec; knotspec = knotspec->next ) {
215	knotspec->preselect();
216	knotspec->select();
217    }
218}
219
220/*-----------------------------------------------------------------------------
221 * Splinespec::layout -
222 *
223 * Client: gl_quilt_to_bezier
224 *-----------------------------------------------------------------------------
225 */
226
227void
228Splinespec::layout( long ncoords )
229{
230
231    long stride = ncoords;
232    for( Knotspec *knotspec = kspec; knotspec; knotspec=knotspec->next ) {
233	knotspec->poststride = (int) stride;
234	stride *= ((knotspec->bend-knotspec->bbegin)*knotspec->order + knotspec->postoffset);
235        knotspec->preoffset  *= knotspec->prestride;
236	knotspec->prewidth  *= knotspec->poststride;
237	knotspec->postwidth *= knotspec->poststride;
238	knotspec->postoffset *= knotspec->poststride;
239        knotspec->ncoords = (int) ncoords;
240    }
241    outcpts = new REAL[stride];
242    assert( outcpts != 0 );
243}
244
245/*-----------------------------------------------------------------------------
246 * Splinespec::copy - copy the control points of current subobject
247 *
248 * Client: gl_quilt_to_bezier
249 *-----------------------------------------------------------------------------
250 */
251
252void
253Splinespec::copy( INREAL *incpts )
254{
255    kspec->copy( incpts, outcpts );
256}
257
258/*-----------------------------------------------------------------------------
259 * Splinespec::setupquilt - assign all quilt variables from knotspec
260 *
261 * Client: gl_quilt_to_bezier
262 *-----------------------------------------------------------------------------
263 */
264
265void
266Splinespec::setupquilt( Quilt_ptr quilt )
267{
268    Quiltspec_ptr qspec = quilt->qspec;
269    quilt->eqspec = qspec + dim;
270    for( Knotspec *knotspec = kspec; knotspec; knotspec=knotspec->next, qspec++ ) {
271	qspec->stride	= knotspec->poststride;
272	qspec->width	= knotspec->bend - knotspec->bbegin;
273	qspec->order	= (int) knotspec->order;
274	qspec->offset	= knotspec->postoffset;
275	qspec->index	= 0;
276	qspec->bdry[0]	= (knotspec->kleft == knotspec->kfirst) ? 1 : 0;
277	qspec->bdry[1]	= (knotspec->kright == knotspec->klast) ? 1 : 0;
278	qspec->breakpoints = new Knot[qspec->width+1];
279	Knot_ptr k =  qspec->breakpoints;
280	for( Breakpt *bk = knotspec->bbegin; bk <= knotspec->bend; bk++ )
281	    *(k++) = bk->value;
282    }
283    quilt->cpts = outcpts;
284    quilt->next = 0;
285}
286
287/*-----------------------------------------------------------------------------
288 * Splinespec::transform - convert a spline to Bezier format
289 *
290 * Client: gl_quilt_to_bezier
291 *-----------------------------------------------------------------------------
292 */
293
294void
295Splinespec::transform( void )
296{
297    Knotspec *knotspec;
298    for( knotspec = kspec; knotspec; knotspec=knotspec->next )
299        knotspec->istransformed = 0;
300
301    for( knotspec = kspec; knotspec; knotspec=knotspec->next ) {
302	for( Knotspec *kspec2 = kspec; kspec2; kspec2=kspec2->next )
303	    kspec2->kspectotrans = knotspec;
304	kspec->transform( outcpts );
305	knotspec->istransformed = 1;
306    }
307}
308
309
310/*-----------------------------------------------------------------------------
311 * Knotspec::Knotspec -  constuct a knot spec
312 *-----------------------------------------------------------------------------
313 */
314
315Knotspec::Knotspec( void )
316{
317    bbegin = 0;
318    sbegin = 0;
319    outkbegin = 0;
320}
321
322/*-----------------------------------------------------------------------------
323 * Knotspec::copy -  copy the control points along minor direction
324 *
325 * Client: Splinespec::copy
326 *-----------------------------------------------------------------------------
327 */
328
329void
330Knotspec::copy( INREAL *inpt, REAL *outpt )
331{
332    inpt = (INREAL *) (((char *) inpt) + preoffset);
333
334    if( next ) {
335        for( REAL *lpt=outpt+prewidth; outpt != lpt; outpt += poststride ) {
336	    next->copy( inpt, outpt );
337	    inpt = (INREAL *) (((char *) inpt) + prestride);
338	}
339   } else {
340        for( REAL *lpt=outpt+prewidth; outpt != lpt; outpt += poststride ) {
341	    pt_io_copy( outpt, inpt );
342	    inpt = (INREAL *) (((char *) inpt) + prestride);
343	}
344     }
345}
346
347/*-----------------------------------------------------------------------------
348 * Knotspec::showpts - print out points before transformation
349 *
350 * Client: Knotspec::select
351 *-----------------------------------------------------------------------------
352 */
353void
354Knotspec::showpts( REAL *outpt )
355{
356    if( next ) {
357        for( REAL *lpt=outpt+prewidth; outpt != lpt; outpt += poststride )
358	    next->showpts( outpt );
359    } else {
360        for( REAL *lpt=outpt+prewidth; outpt != lpt; outpt += poststride )
361	    _glu_dprintf(  "show %g %g %g\n", outpt[0], outpt[1], outpt[2] );
362    }
363}
364
365/*-----------------------------------------------------------------------------
366 * Knotspec::factors - precompute scale factors
367 *	   - overwrites knot vector, actual new knot vector is NOT produced
368 *
369 * Client: Knotspec::select
370 *-----------------------------------------------------------------------------
371 */
372
373void
374Knotspec::factors( void )
375{
376    Knot *mid = (outkend - 1) - order + bend->multi;
377    Knot_ptr fptr = sbegin;
378
379    for( Breakpt *bpt = bend; bpt >= bbegin; bpt-- ) {
380    	mid -= bpt->multi;		// last knot less than knot to insert
381	int def = bpt->def - 1;		// number of knots to insert
382	if( def <= 0 ) continue;
383	Knot kv = bpt->value;		// knot to insert
384
385	Knot *kf = (mid-def) + (order-1);
386	for( Knot *kl = kf + def; kl != kf; kl-- ) {
387	    Knot *kh, *kt;
388	    for( kt=kl, kh=mid; kt != kf; kh--, kt-- )
389		*(fptr++) = (kv - *kh) / (*kt - *kh);
390	    *kl = kv;
391	}
392    }
393}
394
395/*-----------------------------------------------------------------------------
396 * Knotspec::insert - convert subobject in direction of kspec into Bezier
397 *
398 * Client: Knotspec::transform
399 *-----------------------------------------------------------------------------
400 */
401
402void
403Knotspec::insert( REAL *p )
404{
405    Knot_ptr fptr = sbegin;
406    REAL *srcpt = p + prewidth - poststride;
407    REAL *dstpt = p + postwidth + postoffset - poststride;
408    Breakpt *bpt = bend;
409
410   for( REAL *pend = srcpt - poststride*bpt->def; srcpt != pend; pend +=poststride ) {
411	REAL *p1 = srcpt;
412	for( REAL *p2 = srcpt-poststride; p2 != pend; p1 = p2, p2 -= poststride ) {
413	    pt_oo_sum( p1, p1, p2, *fptr, 1.0-*fptr );
414	    fptr++;
415	}
416    }
417
418    for( --bpt; bpt >= bbegin; bpt-- ) {
419
420	for( int multi = bpt->multi; multi > 0; multi-- ) {
421	    pt_oo_copy( dstpt, srcpt );
422	    dstpt -= poststride;
423	    srcpt -= poststride;
424	}
425
426	for( REAL *pend = srcpt - poststride*bpt->def; srcpt != pend; pend +=poststride, dstpt-=poststride ) {
427	    pt_oo_copy( dstpt, srcpt );
428	    REAL *p1 = srcpt;
429
430	    for( REAL *p2 = srcpt-poststride; p2 != pend; p1=p2, p2 -= poststride ) {
431		pt_oo_sum( p1, p1, p2, *fptr, 1.0-*fptr );
432		fptr++;
433	    }
434	}
435    }
436}
437
438/*-----------------------------------------------------------------------------
439 * Knotspec::preselect - initialize kspec for processing
440 *
441 * Client: Splinespec::select
442 *-----------------------------------------------------------------------------
443 */
444
445void
446Knotspec::preselect( void )
447{
448    Knot kval;
449
450    /* position klast after last knot of "last" breakpoint */
451    for( klast = inkend - order, kval = *klast; klast != inkend; klast++ )
452	if( ! identical( *klast, kval ) ) break;
453
454    /* position kfirst after last knot of "first" breakpoint */
455    for( kfirst = inkbegin+order-1, kval= *kfirst;  kfirst != inkend; kfirst++ )
456	if( ! identical( *kfirst, kval ) ) break;
457
458    /* compute multiplicity of first breakpoint */
459    Knot_ptr k;
460    for( k  = kfirst - 1; k >= inkbegin; k-- )
461	if( ! identical( kval, *k ) ) break;
462    k++;
463
464    /* allocate space for breakpoints -
465       use worst case estimate on number of breakpoints */
466
467    bbegin = new Breakpt[(klast - kfirst)+1];
468    /* record multiplicity and value of first breakpoint */
469    bbegin->multi = kfirst - k;
470    bbegin->value = kval;
471    bend = bbegin;
472
473    kleft = kright = kfirst;
474}
475
476
477/*-----------------------------------------------------------------------------
478 * Knotspec::select - Knotspec::select segments and precompute scale factors
479 *
480 * Client: Splinespec::select
481 *-----------------------------------------------------------------------------
482 */
483
484void
485Knotspec::select( void )
486{
487    breakpoints();
488    knots();
489    factors();
490
491    preoffset	= kleft - (inkbegin + order);
492    postwidth	= (int)((bend - bbegin) * order);
493    prewidth 	= (int)((outkend - outkbegin) - order);
494    postoffset  = (bbegin->def > 1) ? (bbegin->def-1) : 0;
495}
496
497/*-----------------------------------------------------------------------------
498 * Knotspec::breakpoints - compute breakpoints for knotspec
499 *
500 * Client: Knotspec::select
501 *-----------------------------------------------------------------------------
502 */
503
504void
505Knotspec::breakpoints( void )
506{
507    Breakpt *ubpt	= bbegin;
508    Breakpt *ubend	= bend;
509    long    nfactors  	= 0;
510
511    ubpt->value	= ubend->value;
512    ubpt->multi	= ubend->multi;
513
514    kleft = kright;
515
516    for( ; kright != klast; kright++ ) {
517        if ( identical(*kright,ubpt->value) ) {
518	    (ubpt->multi)++;
519	} else {
520    	    ubpt->def = (int) (order - ubpt->multi);
521    	    nfactors += (ubpt->def * (ubpt->def - 1)) / 2;
522	    (++ubpt)->value = *kright;
523	    ubpt->multi = 1;
524	}
525    }
526    ubpt->def = (int) (order - ubpt->multi);
527    nfactors += (ubpt->def * (ubpt->def - 1)) / 2;
528
529    bend = ubpt;
530
531    if( nfactors ) {
532        sbegin = new Knot[nfactors];
533    } else {
534	sbegin = NULL;
535    }
536}
537
538
539/*-----------------------------------------------------------------------------
540 * Knotspec::knots - copy relevant subsequence of knots into temporary area
541 *
542 * Client: Knotspec::select
543 *-----------------------------------------------------------------------------
544 */
545
546void
547Knotspec::knots( void )
548{
549    Knot_ptr inkpt = kleft - order;
550    Knot_ptr inkend = kright  + bend->def;
551
552    /* allocate space for knots and factors */
553    outkbegin = new Knot[inkend-inkpt];
554    Knot_ptr outkpt;
555    for( outkpt = outkbegin; inkpt != inkend; inkpt++, outkpt++ )
556	*outkpt = *inkpt;
557
558    outkend = outkpt;
559}
560
561
562/*-----------------------------------------------------------------------------
563 * Knotspec::transform -	convert a spline along a given direction
564 *
565 * Client: Splienspec::transform
566 *-----------------------------------------------------------------------------
567 */
568
569void
570Knotspec::transform( REAL *p )
571{
572   if( next ) {
573	if( this == kspectotrans ) {
574	    next->transform( p );
575	} else {
576	    if( istransformed ) {
577		p += postoffset;
578		for( REAL *pend = p + postwidth; p != pend; p += poststride )
579		    next->transform( p );
580	    } else {
581		REAL *pend = p + prewidth;
582		for( ; p != pend; p += poststride )
583		    next->transform( p );
584	    }
585	}
586   } else {
587	if( this == kspectotrans ) {
588	    insert( p );
589	} else {
590	    if( istransformed ) {
591		p += postoffset;
592		for( REAL *pend = p + postwidth; p != pend; p += poststride )
593		    kspectotrans->insert( p );
594	    } else {
595		REAL *pend = p + prewidth;
596		for( ; p != pend; p += poststride )
597		    kspectotrans->insert( p );
598	    }
599	}
600   }
601}
602
603/*-----------------------------------------------------------------------------
604 * Knotspec::~Knotspec - free space alocated for knotspec
605 *-----------------------------------------------------------------------------
606 */
607
608Knotspec::~Knotspec( void )
609{
610    if( bbegin ) delete[] bbegin;
611    if( sbegin ) delete[] sbegin;
612    if( outkbegin ) delete[] outkbegin;
613}
614
615
616/*-----------------------------------------------------------------------------
617 * pt_io_copy - make internal copy of input cntrl pt. of x coords
618 *-----------------------------------------------------------------------------
619 */
620
621void
622Knotspec::pt_io_copy( REAL *topt, INREAL *frompt )
623{
624    switch( ncoords ) {
625    case 4:
626        topt[3] = (REAL) frompt[3];
627    case 3:
628        topt[2] = (REAL) frompt[2];
629    case 2:
630        topt[1] = (REAL) frompt[1];
631    case 1:
632        topt[0] = (REAL) frompt[0];
633	break;
634    default: {
635	    for( int i = 0; i < ncoords; i++ )
636		*topt++ = (REAL) *frompt++;
637	}
638    }
639}
640
641/*-----------------------------------------------------------------------------
642 * pt_oo_copy - make internal copy of internal cntrl pt. of x coords
643 *-----------------------------------------------------------------------------
644 */
645
646void
647Knotspec::pt_oo_copy( REAL *topt, REAL *frompt )
648{
649    switch( ncoords ) {
650    case 4:
651        topt[3] = frompt[3];
652    case 3:
653        topt[2] = frompt[2];
654    case 2:
655        topt[1] = frompt[1];
656    case 1:
657        topt[0] = frompt[0];
658	break;
659    default:
660	memcpy( topt, frompt, ncoords * sizeof( REAL ) );
661    }
662}
663
664/*-----------------------------------------------------------------------------
665 * pt_oo_sum - compute affine combination of internal cntrl pts
666 *-----------------------------------------------------------------------------
667 */
668
669void
670Knotspec::pt_oo_sum( REAL *x, REAL *y, REAL *z, Knot a, Knot b )
671{
672    switch( ncoords ) {
673    case 4:
674        x[3] = a * y[3]  +  b * z[3];
675    case 3:
676        x[2] = a * y[2]  +  b * z[2];
677    case 2:
678        x[1] = a * y[1]  +  b * z[1];
679    case 1:
680        x[0] = a * y[0]  +  b * z[0];
681	break;
682    default: {
683          for( int i = 0; i < ncoords; i++ )
684              *x++ = a * *y++   +   b * *z++;
685    }
686    }
687}
688