home *** CD-ROM | disk | FTP | other *** search
/ AmigActive 6 / AACD06.ISO / AACD / System / Mesa-3.1 / src-glu / nurbsutl.c < prev    next >
C/C++ Source or Header  |  1999-08-18  |  35KB  |  1,404 lines
  1. /* $Id: nurbsutl.c,v 1.1.1.1 1999/08/19 00:55:42 jtg Exp $ */
  2.  
  3. /*
  4.  * Mesa 3-D graphics library
  5.  * Version:  2.4
  6.  * Copyright (C) 1995-1997  Brian Paul
  7.  *
  8.  * This library is free software; you can redistribute it and/or
  9.  * modify it under the terms of the GNU Library General Public
  10.  * License as published by the Free Software Foundation; either
  11.  * version 2 of the License, or (at your option) any later version.
  12.  *
  13.  * This library is distributed in the hope that it will be useful,
  14.  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  15.  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  16.  * Library General Public License for more details.
  17.  *
  18.  * You should have received a copy of the GNU Library General Public
  19.  * License along with this library; if not, write to the Free
  20.  * Software Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
  21.  */
  22.  
  23.  
  24. /*
  25.  * $Log: nurbsutl.c,v $
  26.  * Revision 1.1.1.1  1999/08/19 00:55:42  jtg
  27.  * Imported sources
  28.  *
  29.  * Revision 1.8  1999/06/08 00:44:51  brianp
  30.  * OpenStep updates (pete@ohm.york.ac.uk)
  31.  *
  32.  * Revision 1.7  1998/07/26 02:07:59  brianp
  33.  * updated for Windows compilation per Ted Jump
  34.  *
  35.  * Revision 1.6  1997/10/29 02:02:20  brianp
  36.  * various MS Windows compiler changes (David Bucciarelli, v20 3dfx driver)
  37.  *
  38.  * Revision 1.5  1997/07/24 01:28:44  brianp
  39.  * changed precompiled header symbol from PCH to PC_HEADER
  40.  *
  41.  * Revision 1.4  1997/05/28 02:29:38  brianp
  42.  * added support for precompiled headers (PCH), inserted APIENTRY keyword
  43.  *
  44.  * Revision 1.3  1997/05/27 03:19:54  brianp
  45.  * minor clean-up
  46.  *
  47.  * Revision 1.2  1997/05/27 03:00:16  brianp
  48.  * incorporated Bogdan's new NURBS code
  49.  *
  50.  * Revision 1.1  1996/09/27 01:19:39  brianp
  51.  * Initial revision
  52.  *
  53.  */
  54.  
  55.  
  56. /*
  57.  * NURBS implementation written by Bogdan Sikorski (bogdan@cira.it)
  58.  * See README2 for more info.
  59.  */
  60.  
  61.  
  62. #ifdef PC_HEADER
  63. #include "all.h"
  64. #else
  65. #include <math.h>
  66. #include <stdlib.h>
  67. #include "gluP.h"
  68. #include "nurbs.h"
  69. #endif
  70.  
  71.  
  72. GLenum
  73. test_knot(GLint nknots, GLfloat *knot, GLint order)
  74. {
  75.     GLsizei i;
  76.     GLint knot_mult;
  77.     GLfloat tmp_knot;
  78.  
  79.     tmp_knot=knot[0];
  80.     knot_mult=1;
  81.     for(i=1;i<nknots;i++)
  82.     {
  83.         if(knot[i] < tmp_knot)
  84.             return GLU_NURBS_ERROR4;
  85.         if(fabs(tmp_knot-knot[i]) > EPSILON)
  86.         {
  87.             if(knot_mult>order)
  88.                 return GLU_NURBS_ERROR5;
  89.             knot_mult=1;
  90.             tmp_knot=knot[i];
  91.         }
  92.         else
  93.             ++knot_mult;
  94.     }
  95.     return GLU_NO_ERROR;
  96. }
  97.  
  98. static int
  99. /* qsort function */
  100. #if defined(WIN32) && !defined(OPENSTEP)
  101. __cdecl
  102. #endif
  103. knot_sort(const void *a, const void *b)
  104. {
  105.     GLfloat x,y;
  106.  
  107.     x=*((GLfloat *)a);
  108.     y=*((GLfloat *)b);
  109.     if(fabs(x-y) < EPSILON)
  110.         return 0;
  111.     if(x > y)
  112.         return 1;
  113.     return -1;
  114. }
  115.  
  116. /* insert into dest knot all values within the valid range from src knot */
  117. /* that do not appear in dest */
  118. void
  119. collect_unified_knot(knot_str_type *dest, knot_str_type *src,
  120.     GLfloat maximal_min_knot, GLfloat minimal_max_knot)
  121. {
  122.     GLfloat *src_knot,*dest_knot;
  123.     GLint src_t_min,src_t_max,dest_t_min,dest_t_max;
  124.     GLint src_nknots,dest_nknots;
  125.     GLint i,j,k,new_cnt;
  126.     GLboolean    not_found_flag;
  127.  
  128.     src_knot=src->unified_knot;
  129.     dest_knot=dest->unified_knot;
  130.     src_t_min=src->t_min;
  131.     src_t_max=src->t_max;
  132.     dest_t_min=dest->t_min;
  133.     dest_t_max=dest->t_max;
  134.     src_nknots=src->unified_nknots;
  135.     dest_nknots=dest->unified_nknots;
  136.  
  137.     k=new_cnt=dest_nknots;
  138.     for(i=src_t_min;i<=src_t_max;i++)
  139.         if(src_knot[i] - maximal_min_knot > -EPSILON &&
  140.         src_knot[i] - minimal_max_knot < EPSILON)
  141.     {
  142.         not_found_flag=GL_TRUE;
  143.         for(j=dest_t_min;j<=dest_t_max;j++)
  144.             if(fabs(dest_knot[j]-src_knot[i]) < EPSILON)
  145.             {
  146.                 not_found_flag=GL_FALSE;
  147.                 break;
  148.             }
  149.         if(not_found_flag)
  150.         {
  151.             /* knot from src is not in dest - add this knot to dest */
  152.             dest_knot[k++]=src_knot[i];
  153.             ++new_cnt;
  154.             ++(dest->t_max); /* the valid range widens */
  155.             ++(dest->delta_nknots); /* increment the extra knot value counter */
  156.         }
  157.     }
  158.     dest->unified_nknots=new_cnt;
  159.     qsort((void *)dest_knot,(size_t)new_cnt,(size_t)sizeof(GLfloat),
  160.         &knot_sort);
  161. }
  162.  
  163. /* basing on the new common knot range for all attributes set */
  164. /* t_min and t_max values for each knot - they will be used later on */
  165. /* by explode_knot() and calc_new_ctrl_pts */
  166. static void
  167. set_new_t_min_t_max(knot_str_type *geom_knot, knot_str_type *color_knot,
  168.     knot_str_type *normal_knot, knot_str_type *texture_knot,
  169.     GLfloat maximal_min_knot, GLfloat minimal_max_knot)
  170. {
  171.     GLuint    t_min,t_max,cnt;
  172.  
  173.     if(minimal_max_knot-maximal_min_knot < EPSILON)
  174.     {
  175.         /* knot common range empty */
  176.         geom_knot->t_min=geom_knot->t_max=0;
  177.         color_knot->t_min=color_knot->t_max=0;
  178.         normal_knot->t_min=normal_knot->t_max=0;
  179.         texture_knot->t_min=texture_knot->t_max=0;
  180.     }
  181.     else
  182.     {
  183.         if(geom_knot->unified_knot!=NULL)
  184.         {
  185.             cnt=geom_knot->unified_nknots;
  186.             for(t_min=0;t_min<cnt;t_min++)
  187.                 if(fabs((geom_knot->unified_knot)[t_min] - maximal_min_knot) <
  188.                         EPSILON)
  189.                     break;
  190.             for(t_max=cnt-1;t_max;t_max--)
  191.                 if(fabs((geom_knot->unified_knot)[t_max] - minimal_max_knot) < 
  192.                         EPSILON)
  193.                     break;
  194.         }
  195.         else
  196.         if(geom_knot->nknots)
  197.         {
  198.             cnt=geom_knot->nknots;
  199.             for(t_min=0;t_min<cnt;t_min++)
  200.                 if(fabs((geom_knot->knot)[t_min] - maximal_min_knot) < EPSILON)
  201.                     break;
  202.             for(t_max=cnt-1;t_max;t_max--)
  203.                 if(fabs((geom_knot->knot)[t_max] - minimal_max_knot) < EPSILON)
  204.                     break;
  205.         }
  206.         geom_knot->t_min=t_min;
  207.         geom_knot->t_max=t_max;
  208.         if(color_knot->unified_knot!=NULL)
  209.         {
  210.             cnt=color_knot->unified_nknots;
  211.             for(t_min=0;t_min<cnt;t_min++)
  212.                 if(fabs((color_knot->unified_knot)[t_min] - maximal_min_knot) <
  213.                         EPSILON)
  214.                     break;
  215.             for(t_max=cnt-1;t_max;t_max--)
  216.                 if(fabs((color_knot->unified_knot)[t_max] - minimal_max_knot) < 
  217.                         EPSILON)
  218.                     break;
  219.             color_knot->t_min=t_min;
  220.             color_knot->t_max=t_max;
  221.         }
  222.         if(normal_knot->unified_knot!=NULL)
  223.         {
  224.             cnt=normal_knot->unified_nknots;
  225.             for(t_min=0;t_min<cnt;t_min++)
  226.                 if(fabs((normal_knot->unified_knot)[t_min] - maximal_min_knot) <
  227.                         EPSILON)
  228.                     break;
  229.             for(t_max=cnt-1;t_max;t_max--)
  230.                 if(fabs((normal_knot->unified_knot)[t_max] - minimal_max_knot) < 
  231.                         EPSILON)
  232.                     break;
  233.             normal_knot->t_min=t_min;
  234.             normal_knot->t_max=t_max;
  235.         }
  236.         if(texture_knot->unified_knot!=NULL)
  237.         {
  238.             cnt=texture_knot->unified_nknots;
  239.             for(t_min=0;t_min<cnt;t_min++)
  240.                 if(fabs((texture_knot->unified_knot)[t_min] - maximal_min_knot) 
  241.                         < EPSILON)
  242.                     break;
  243.             for(t_max=cnt-1;t_max;t_max--)
  244.                 if(fabs((texture_knot->unified_knot)[t_max] - minimal_max_knot) 
  245.                         < EPSILON)
  246.                     break;
  247.             texture_knot->t_min=t_min;
  248.             texture_knot->t_max=t_max;
  249.         }
  250.     }
  251. }
  252.  
  253. /* modify all knot valid ranges in such a way that all have the same */
  254. /* range, common to all knots */
  255. /* do this by knot insertion */
  256. GLenum
  257. select_knot_working_range(GLUnurbsObj *nobj,knot_str_type *geom_knot,
  258.     knot_str_type *color_knot, knot_str_type *normal_knot,
  259.     knot_str_type *texture_knot)
  260. {
  261.     GLint max_nknots;
  262.     GLfloat maximal_min_knot,minimal_max_knot;
  263.     GLint i;
  264.  
  265.     /* find the maximum modified knot length */
  266.     max_nknots=geom_knot->nknots;
  267.     if(color_knot->unified_knot)
  268.         max_nknots+=color_knot->nknots;
  269.     if(normal_knot->unified_knot)
  270.         max_nknots+=normal_knot->nknots;
  271.     if(texture_knot->unified_knot)
  272.         max_nknots+=texture_knot->nknots;
  273.     maximal_min_knot=(geom_knot->knot)[geom_knot->t_min];
  274.     minimal_max_knot=(geom_knot->knot)[geom_knot->t_max];
  275.     /* any attirb data ? */
  276.     if(max_nknots!=geom_knot->nknots)
  277.     {
  278.         /* allocate space for the unified knots */
  279.         if((geom_knot->unified_knot=
  280.                 (GLfloat *)malloc(sizeof(GLfloat)*max_nknots))==NULL)
  281.         {
  282.             call_user_error(nobj,GLU_OUT_OF_MEMORY);
  283.             return GLU_ERROR;
  284.         }
  285.         /* copy the original knot to the unified one */
  286.         geom_knot->unified_nknots=geom_knot->nknots;
  287.         for(i=0;i<geom_knot->nknots;i++)
  288.             (geom_knot->unified_knot)[i]=(geom_knot->knot)[i];
  289.         if(color_knot->unified_knot)
  290.         {
  291.             if((color_knot->knot)[color_knot->t_min] - maximal_min_knot >
  292.                     EPSILON)
  293.                 maximal_min_knot=(color_knot->knot)[color_knot->t_min];
  294.             if(minimal_max_knot - (color_knot->knot)[color_knot->t_max] >
  295.                     EPSILON)
  296.                 minimal_max_knot=(color_knot->knot)[color_knot->t_max];
  297.             if((color_knot->unified_knot=
  298.                     (GLfloat *)malloc(sizeof(GLfloat)*max_nknots))==NULL)
  299.             {
  300.                 free(geom_knot->unified_knot);
  301.                 call_user_error(nobj,GLU_OUT_OF_MEMORY);
  302.                 return GLU_ERROR;
  303.             }
  304.             /* copy the original knot to the unified one */
  305.             color_knot->unified_nknots=color_knot->nknots;
  306.             for(i=0;i<color_knot->nknots;i++)
  307.                 (color_knot->unified_knot)[i]=(color_knot->knot)[i];
  308.         }
  309.         if(normal_knot->unified_knot)
  310.         {
  311.             if((normal_knot->knot)[normal_knot->t_min] - maximal_min_knot >
  312.                     EPSILON)
  313.                 maximal_min_knot=(normal_knot->knot)[normal_knot->t_min];
  314.             if(minimal_max_knot - (normal_knot->knot)[normal_knot->t_max] >
  315.                     EPSILON)
  316.                 minimal_max_knot=(normal_knot->knot)[normal_knot->t_max];
  317.             if((normal_knot->unified_knot=
  318.                     (GLfloat *)malloc(sizeof(GLfloat)*max_nknots))==NULL)
  319.             {
  320.                 free(geom_knot->unified_knot);
  321.                 free(color_knot->unified_knot);
  322.                 call_user_error(nobj,GLU_OUT_OF_MEMORY);
  323.                 return GLU_ERROR;
  324.             }
  325.             /* copy the original knot to the unified one */
  326.             normal_knot->unified_nknots=normal_knot->nknots;
  327.             for(i=0;i<normal_knot->nknots;i++)
  328.                 (normal_knot->unified_knot)[i]=(normal_knot->knot)[i];
  329.         }
  330.         if(texture_knot->unified_knot)
  331.         {
  332.             if((texture_knot->knot)[texture_knot->t_min] - maximal_min_knot >
  333.                     EPSILON)
  334.                 maximal_min_knot=(texture_knot->knot)[texture_knot->t_min];
  335.             if(minimal_max_knot - (texture_knot->knot)[texture_knot->t_max] >
  336.                     EPSILON)
  337.                 minimal_max_knot=(texture_knot->knot)[texture_knot->t_max];
  338.             if((texture_knot->unified_knot=
  339.                     (GLfloat *)malloc(sizeof(GLfloat)*max_nknots))==NULL)
  340.             {
  341.                 free(geom_knot->unified_knot);
  342.                 free(color_knot->unified_knot);
  343.                 free(normal_knot->unified_knot);
  344.                 call_user_error(nobj,GLU_OUT_OF_MEMORY);
  345.                 return GLU_ERROR;
  346.             }
  347.             /* copy the original knot to the unified one */
  348.             texture_knot->unified_nknots=texture_knot->nknots;
  349.             for(i=0;i<texture_knot->nknots;i++)
  350.                 (texture_knot->unified_knot)[i]=(texture_knot->knot)[i];
  351.         }
  352.         /* work on the geometry knot with all additional knot values */
  353.         /* appearing in attirbutive knots */
  354.         if(minimal_max_knot-maximal_min_knot < EPSILON)
  355.         {
  356.             /* empty working range */
  357.             geom_knot->unified_nknots=0;
  358.             color_knot->unified_nknots=0;
  359.             normal_knot->unified_nknots=0;
  360.             texture_knot->unified_nknots=0;
  361.         }
  362.         else
  363.         {
  364.             if(color_knot->unified_knot)
  365.                 collect_unified_knot(geom_knot,color_knot,maximal_min_knot,
  366.                     minimal_max_knot);
  367.             if(normal_knot->unified_knot)
  368.                 collect_unified_knot(geom_knot,normal_knot,maximal_min_knot,
  369.                     minimal_max_knot);
  370.             if(texture_knot->unified_knot)
  371.                 collect_unified_knot(geom_knot,texture_knot,maximal_min_knot,
  372.                     minimal_max_knot);
  373.             /* since we have now built the "unified" geometry knot */
  374.             /* add same knot values to all attributive knots */
  375.             if(color_knot->unified_knot)
  376.                 collect_unified_knot(color_knot,geom_knot,maximal_min_knot,
  377.                     minimal_max_knot);
  378.             if(normal_knot->unified_knot)
  379.                 collect_unified_knot(normal_knot,geom_knot,maximal_min_knot,
  380.                     minimal_max_knot);
  381.             if(texture_knot->unified_knot)
  382.                 collect_unified_knot(texture_knot,geom_knot,maximal_min_knot,
  383.                     minimal_max_knot);
  384.         }
  385.     }
  386.     set_new_t_min_t_max(geom_knot,color_knot,normal_knot,texture_knot,
  387.         maximal_min_knot,minimal_max_knot);
  388.     return GLU_NO_ERROR;
  389. }
  390.  
  391. void
  392. free_unified_knots(knot_str_type *geom_knot, knot_str_type *color_knot,
  393.     knot_str_type *normal_knot, knot_str_type *texture_knot)
  394. {
  395.     if(geom_knot->unified_knot)
  396.         free(geom_knot->unified_knot);
  397.     if(color_knot->unified_knot)
  398.         free(color_knot->unified_knot);
  399.     if(normal_knot->unified_knot)
  400.         free(normal_knot->unified_knot);
  401.     if(texture_knot->unified_knot)
  402.         free(texture_knot->unified_knot);
  403. }
  404.  
  405. GLenum
  406. explode_knot(knot_str_type *the_knot)
  407. {
  408.     GLfloat *knot,*new_knot;
  409.     GLint nknots,n_new_knots=0;
  410.     GLint t_min,t_max;
  411.     GLint ord;
  412.     GLsizei i,j,k;
  413.     GLfloat tmp_float;
  414.  
  415.     if(the_knot->unified_knot)
  416.     {
  417.         knot=the_knot->unified_knot;
  418.         nknots=the_knot->unified_nknots;
  419.     }
  420.     else
  421.     {
  422.         knot=the_knot->knot;
  423.         nknots=the_knot->nknots;
  424.     }
  425.     ord=the_knot->order;
  426.     t_min=the_knot->t_min;
  427.     t_max=the_knot->t_max;
  428.  
  429.     for(i=t_min;i<=t_max;)
  430.     {
  431.         tmp_float=knot[i];
  432.         for(j=0;j<ord && (i+j)<=t_max;j++)
  433.             if(fabs(tmp_float-knot[i+j])>EPSILON)
  434.                 break;
  435.         n_new_knots+=ord-j;
  436.         i+=j;
  437.     }
  438.     /* alloc space for new_knot */
  439.     if((new_knot=(GLfloat *)malloc(sizeof(GLfloat)*(nknots+n_new_knots)))==NULL)
  440.     {
  441.         return GLU_OUT_OF_MEMORY;
  442.     }
  443.     /* fill in new knot */
  444.     for(j=0;j<t_min;j++)
  445.         new_knot[j]=knot[j];
  446.     for(i=j;i<=t_max;i++)
  447.     {
  448.         tmp_float=knot[i];
  449.         for(k=0;k<ord;k++)
  450.         {
  451.             new_knot[j++]=knot[i];
  452.             if(tmp_float==knot[i+1])
  453.                 i++;
  454.         }
  455.     }
  456.     for(i=t_max+1;i<(int)nknots;i++)
  457.         new_knot[j++]=knot[i];
  458.     /* fill in the knot structure */
  459.     the_knot->new_knot=new_knot;
  460.     the_knot->delta_nknots+=n_new_knots;
  461.     the_knot->t_max+=n_new_knots;
  462.     return GLU_NO_ERROR;
  463. }
  464.  
  465. GLenum
  466. calc_alphas(knot_str_type *the_knot)
  467. {
  468.     GLfloat tmp_float;
  469.     int i,j,k,m,n;
  470.     int order;
  471.     GLfloat *alpha,*alpha_new,*tmp_alpha;
  472.     GLfloat denom;
  473.     GLfloat *knot,*new_knot;
  474.  
  475.  
  476.     knot=the_knot->knot;
  477.     order=the_knot->order;
  478.     new_knot=the_knot->new_knot;
  479.     n=the_knot->nknots-the_knot->order;
  480.     m=n+the_knot->delta_nknots;
  481.     if((alpha=(GLfloat *)malloc(sizeof(GLfloat)*n*m))==NULL)
  482.     {
  483.         return GLU_OUT_OF_MEMORY;
  484.     }
  485.     if((alpha_new=(GLfloat *)malloc(sizeof(GLfloat)*n*m))==NULL)
  486.     {
  487.         free(alpha);
  488.         return GLU_OUT_OF_MEMORY;
  489.     }
  490.     for(j=0;j<m;j++)
  491.     {
  492.         for(i=0;i<n;i++)
  493.         {
  494.             if((knot[i] <= new_knot[j]) && (new_knot[j] < knot[i+1]))
  495.                 tmp_float=1.0;
  496.             else
  497.                 tmp_float=0.0;
  498.             alpha[i+j*n]=tmp_float;
  499.         }
  500.     }
  501.     for(k=1;k<order;k++)
  502.     {
  503.         for(j=0;j<m;j++)
  504.             for(i=0;i<n;i++)
  505.             {
  506.                 denom=knot[i+k]-knot[i];
  507.                 if(fabs(denom)<EPSILON)
  508.                     tmp_float=0.0;
  509.                 else
  510.                     tmp_float=(new_knot[j+k]-knot[i])/denom*
  511.                         alpha[i+j*n];
  512.                 denom=knot[i+k+1]-knot[i+1];
  513.                 if(fabs(denom)>EPSILON)
  514.                     tmp_float+=(knot[i+k+1]-new_knot[j+k])/denom*
  515.                         alpha[(i+1)+j*n];
  516.                 alpha_new[i+j*n]=tmp_float;
  517.             }
  518.         tmp_alpha=alpha_new;
  519.         alpha_new=alpha;
  520.         alpha=tmp_alpha;
  521.     }
  522.     the_knot->alpha=alpha;
  523.     free(alpha_new);
  524.     return GLU_NO_ERROR;
  525. }
  526.  
  527. GLenum
  528. calc_new_ctrl_pts(GLfloat *ctrl,GLint stride,knot_str_type *the_knot,
  529.     GLint dim,GLfloat **new_ctrl,GLint *ncontrol)
  530. {
  531.     GLsizei i,j,k,l,m,n;
  532.     GLsizei index1,index2;
  533.     GLfloat *alpha;
  534.     GLfloat *new_knot;
  535.  
  536.     new_knot=the_knot->new_knot;
  537.     n=the_knot->nknots-the_knot->order;
  538.     alpha=the_knot->alpha;
  539.  
  540.     m=the_knot->t_max+1-the_knot->t_min-the_knot->order;
  541.     k=the_knot->t_min;
  542.     /* allocate space for new control points */
  543.     if((*new_ctrl=(GLfloat *)malloc(sizeof(GLfloat)*dim*m))==NULL)
  544.     {
  545.         return GLU_OUT_OF_MEMORY;
  546.     }
  547.     for(j=0;j<m;j++)
  548.     {
  549.         for(l=0;l<dim;l++)
  550.             (*new_ctrl)[j*dim+l]=0.0;
  551.         for(i=0;i<n;i++)
  552.         {
  553.             index1=i+(j+k)*n;
  554.             index2=i*stride;
  555.             for(l=0;l<dim;l++)
  556.                 (*new_ctrl)[j*dim+l]+=alpha[index1]*ctrl[index2+l];
  557.         }
  558.     }
  559.     *ncontrol=(GLint)m;
  560.     return GLU_NO_ERROR;
  561. }
  562.  
  563. static GLint
  564. calc_factor(GLfloat *pts,GLint order,GLint indx,GLint stride,GLfloat tolerance,
  565.     GLint dim)
  566. {
  567.     GLdouble model[16],proj[16];
  568.     GLint viewport[4];
  569.     GLdouble x,y,z,w,winx1,winy1,winz,winx2,winy2;
  570.     GLint i;
  571.     GLdouble len,dx,dy;
  572.  
  573.     glGetDoublev(GL_MODELVIEW_MATRIX,model);
  574.     glGetDoublev(GL_PROJECTION_MATRIX,proj);
  575.     glGetIntegerv(GL_VIEWPORT,viewport);
  576.     if(dim==4)
  577.     {
  578.         w=(GLdouble)pts[indx+3];
  579.         x=(GLdouble)pts[indx]/w;
  580.         y=(GLdouble)pts[indx+1]/w;
  581.         z=(GLdouble)pts[indx+2]/w;
  582.         gluProject(x,y,z,model,proj,viewport,&winx1,&winy1,&winz);
  583.         len=0.0;
  584.         for(i=1;i<order;i++)
  585.         {
  586.             w=(GLdouble)pts[indx+i*stride+3];
  587.             x=(GLdouble)pts[indx+i*stride]/w;
  588.             y=(GLdouble)pts[indx+i*stride+1]/w;
  589.             z=(GLdouble)pts[indx+i*stride+2]/w;
  590.             if(gluProject(x,y,z,model,proj,viewport,&winx2,&winy2,&winz))
  591.             {
  592.                 dx=winx2-winx1;
  593.                 dy=winy2-winy1;
  594.                 len+=sqrt(dx*dx+dy*dy);
  595.             }
  596.             winx1=winx2; winy1=winy2;
  597.         }
  598.     }
  599.     else
  600.     {
  601.         x=(GLdouble)pts[indx];
  602.         y=(GLdouble)pts[indx+1];
  603.         if(dim==2)
  604.             z=0.0;
  605.         else
  606.             z=(GLdouble)pts[indx+2];
  607.         gluProject(x,y,z,model,proj,viewport,&winx1,&winy1,&winz);
  608.         len=0.0;
  609.         for(i=1;i<order;i++)
  610.         {
  611.             x=(GLdouble)pts[indx+i*stride];
  612.             y=(GLdouble)pts[indx+i*stride+1];
  613.             if(dim==2)
  614.                 z=0.0;
  615.             else
  616.                 z=(GLdouble)pts[indx+i*stride+2];
  617.             if(gluProject(x,y,z,model,proj,viewport,&winx2,&winy2,&winz))
  618.             {
  619.                 dx=winx2-winx1;
  620.                 dy=winy2-winy1;
  621.                 len+=sqrt(dx*dx+dy*dy);
  622.             }
  623.             winx1=winx2; winy1=winy2;
  624.         }
  625.     }
  626.     len /= tolerance;
  627.     return ((GLint)len+1);
  628. }
  629.  
  630. /* we can't use the Mesa evaluators - no way to get the point coords */
  631. /* so we use our own Bezier point calculus routines */
  632. /* because I'm lazy, I reuse the ones from eval.c */
  633.  
  634. static void
  635. bezier_curve(GLfloat *cp, GLfloat *out, GLfloat t,
  636.     GLuint dim, GLuint order, GLint offset)
  637. {
  638.     GLfloat s, powert;
  639.     GLuint i, k, bincoeff;
  640.  
  641.     if(order >= 2)
  642.     { 
  643.         bincoeff = order-1;
  644.         s = 1.0-t;
  645.  
  646.         for(k=0; k<dim; k++)
  647.             out[k] = s*cp[k] + bincoeff*t*cp[offset+k];
  648.  
  649.         for(i=2, cp+=2*offset, powert=t*t; i<order; i++, powert*=t, cp +=offset)
  650.         {
  651.             bincoeff *= order-i;
  652.             bincoeff /= i;
  653.  
  654.             for(k=0; k<dim; k++)
  655.                 out[k] = s*out[k] + bincoeff*powert*cp[k];
  656.         }
  657.     }
  658.     else /* order=1 -> constant curve */
  659.     { 
  660.         for(k=0; k<dim; k++)
  661.             out[k] = cp[k];
  662.     } 
  663. }
  664.  
  665. static GLint
  666. calc_parametric_factor(GLfloat *pts,GLint order,GLint indx,GLint stride,
  667.     GLfloat tolerance,GLint dim)
  668. {
  669.     GLdouble model[16],proj[16];
  670.     GLint viewport[4];
  671.     GLdouble x,y,z,w,x1,y1,z1,x2,y2,z2,x3,y3,z3;
  672.     GLint i;
  673.     GLint P;
  674.     GLfloat bez_pt[4];
  675.     GLdouble len=0.0,tmp,z_med;
  676.  
  677.     P = 2*(order+2);
  678.     glGetDoublev(GL_MODELVIEW_MATRIX,model);
  679.     glGetDoublev(GL_PROJECTION_MATRIX,proj);
  680.     glGetIntegerv(GL_VIEWPORT,viewport);
  681.     z_med = (viewport[2] + viewport[3]) * 0.5;
  682.     switch(dim)
  683.     {
  684.         case 4:
  685.             for(i=1;i<P;i++)
  686.             {
  687.                 bezier_curve(pts+indx, bez_pt, (GLfloat)i/(GLfloat)P, 4,
  688.                     order,stride);
  689.                 w = (GLdouble)bez_pt[3];
  690.                 x = (GLdouble)bez_pt[0] / w;
  691.                 y = (GLdouble)bez_pt[1] / w;
  692.                 z = (GLdouble)bez_pt[2] / w;
  693.                 gluProject(x,y,z,model,proj,viewport,&x3,&y3,&z3);
  694.                 z3 *= z_med;
  695.                 bezier_curve(pts+indx, bez_pt, (GLfloat)(i-1)/(GLfloat)P, 4,
  696.                     order,stride);
  697.                 w = (GLdouble)bez_pt[3];
  698.                 x = (GLdouble)bez_pt[0] / w;
  699.                 y = (GLdouble)bez_pt[1] / w;
  700.                 z = (GLdouble)bez_pt[2] / w;
  701.                 gluProject(x,y,z,model,proj,viewport,&x1,&y1,&z1);
  702.                 z1 *= z_med;
  703.                 bezier_curve(pts+indx, bez_pt, (GLfloat)(i+1)/(GLfloat)P, 4,
  704.                     order,stride);
  705.                 w = (GLdouble)bez_pt[3];
  706.                 x = (GLdouble)bez_pt[0] / w;
  707.                 y = (GLdouble)bez_pt[1] / w;
  708.                 z = (GLdouble)bez_pt[2] / w;
  709.                 gluProject(x,y,z,model,proj,viewport,&x2,&y2,&z2);
  710.                 z2 *= z_med;
  711.                 /* calc distance between point (x3,y3,z3) and line segment */
  712.                 /* <x1,y1,z1><x2,y2,z2> */
  713.                 x = x2-x1;
  714.                 y = y2-y1;
  715.                 z = z2-z1;
  716.                 tmp = sqrt(x*x+y*y+z*z);
  717.                 x /= tmp;
  718.                 y /= tmp;
  719.                 z /= tmp;
  720.                 tmp = x3*x+y3*y+z3*z-x1*x-y1*y-z1*z;
  721.                 x = x1+x*tmp-x3;
  722.                 y = y1+y*tmp-y3;
  723.                 z = z1+z*tmp-z3;
  724.                 tmp = sqrt(x*x+y*y+z*z);
  725.                 if(tmp > len)
  726.                     len = tmp;
  727.             }
  728.             break;
  729.         case 3:
  730.             for(i=1;i<P;i++)
  731.             {
  732.                 bezier_curve(pts+indx, bez_pt, (GLfloat)i/(GLfloat)P, 3,
  733.                     order,stride);
  734.                 x = (GLdouble)bez_pt[0];
  735.                 y = (GLdouble)bez_pt[1];
  736.                 z = (GLdouble)bez_pt[2];
  737.                 gluProject(x,y,z,model,proj,viewport,&x3,&y3,&z3);
  738.                 z3 *= z_med;
  739.                 bezier_curve(pts+indx, bez_pt, (GLfloat)(i-1)/(GLfloat)P, 3,
  740.                     order,stride);
  741.                 x = (GLdouble)bez_pt[0];
  742.                 y = (GLdouble)bez_pt[1];
  743.                 z = (GLdouble)bez_pt[2];
  744.                 gluProject(x,y,z,model,proj,viewport,&x1,&y1,&z1);
  745.                 z1 *= z_med;
  746.                 bezier_curve(pts+indx, bez_pt, (GLfloat)(i+1)/(GLfloat)P, 3,
  747.                     order,stride);
  748.                 x = (GLdouble)bez_pt[0];
  749.                 y = (GLdouble)bez_pt[1];
  750.                 z = (GLdouble)bez_pt[2];
  751.                 gluProject(x,y,z,model,proj,viewport,&x2,&y2,&z2);
  752.                 z2 *= z_med;
  753.                 /* calc distance between point (x3,y3,z3) and line segment */
  754.                 /* <x1,y1,z1><x2,y2,z2> */
  755.                 x = x2-x1;
  756.                 y = y2-y1;
  757.                 z = z2-z1;
  758.                 tmp = sqrt(x*x+y*y+z*z);
  759.                 x /= tmp;
  760.                 y /= tmp;
  761.                 z /= tmp;
  762.                 tmp = x3*x+y3*y+z3*z-x1*x-y1*y-z1*z;
  763.                 x = x1+x*tmp-x3;
  764.                 y = y1+y*tmp-y3;
  765.                 z = z1+z*tmp-z3;
  766.                 tmp = sqrt(x*x+y*y+z*z);
  767.                 if(tmp > len)
  768.                     len = tmp;
  769.             }
  770.             break;
  771.         case 2:
  772.             for(i=1;i<P;i++)
  773.             {
  774.                 bezier_curve(pts+indx, bez_pt, (GLfloat)i/(GLfloat)P, 2,
  775.                     order,stride);
  776.                 x = (GLdouble)bez_pt[0];
  777.                 y = (GLdouble)bez_pt[1];
  778.                 z = 0.0;
  779.                 gluProject(x,y,z,model,proj,viewport,&x3,&y3,&z3);
  780.                 z3 *= z_med;
  781.                 bezier_curve(pts+indx, bez_pt, (GLfloat)(i-1)/(GLfloat)P, 2,
  782.                     order,stride);
  783.                 x = (GLdouble)bez_pt[0];
  784.                 y = (GLdouble)bez_pt[1];
  785.                 z = 0.0;
  786.                 gluProject(x,y,z,model,proj,viewport,&x1,&y1,&z1);
  787.                 z1 *= z_med;
  788.                 bezier_curve(pts+indx, bez_pt, (GLfloat)(i+1)/(GLfloat)P, 2,
  789.                     order,stride);
  790.                 x = (GLdouble)bez_pt[0];
  791.                 y = (GLdouble)bez_pt[1];
  792.                 z = 0.0;
  793.                 gluProject(x,y,z,model,proj,viewport,&x2,&y2,&z2);
  794.                 z2 *= z_med;
  795.                 /* calc distance between point (x3,y3,z3) and line segment */
  796.                 /* <x1,y1,z1><x2,y2,z2> */
  797.                 x = x2-x1;
  798.                 y = y2-y1;
  799.                 z = z2-z1;
  800.                 tmp = sqrt(x*x+y*y+z*z);
  801.                 x /= tmp;
  802.                 y /= tmp;
  803.                 z /= tmp;
  804.                 tmp = x3*x+y3*y+z3*z-x1*x-y1*y-z1*z;
  805.                 x = x1+x*tmp-x3;
  806.                 y = y1+y*tmp-y3;
  807.                 z = z1+z*tmp-z3;
  808.                 tmp = sqrt(x*x+y*y+z*z);
  809.                 if(tmp > len)
  810.                     len = tmp;
  811.             }
  812.             break;
  813.  
  814.     }
  815.     if(len < tolerance)
  816.         return (order);
  817.     else
  818.         return (GLint)(sqrt(len/tolerance)*(order+2)+1);
  819. }
  820.  
  821. static GLenum
  822. calc_sampling_3D(new_ctrl_type *new_ctrl, GLfloat tolerance, GLint dim,
  823.     GLint uorder, GLint vorder, GLint **ufactors, GLint **vfactors)
  824. {
  825.     GLfloat        *ctrl;
  826.     GLint        tmp_factor1,tmp_factor2;
  827.     GLint        ufactor_cnt,vfactor_cnt;
  828.     GLint        offset1,offset2,offset3;
  829.     GLint        i,j;
  830.  
  831.     ufactor_cnt=new_ctrl->s_bezier_cnt;
  832.     vfactor_cnt=new_ctrl->t_bezier_cnt;
  833.     if((*ufactors=(GLint *)malloc(sizeof(GLint)*ufactor_cnt*3))
  834.             ==NULL)
  835.     {
  836.         return GLU_OUT_OF_MEMORY;
  837.     }
  838.     if((*vfactors=(GLint *)malloc(sizeof(GLint)*vfactor_cnt*3))
  839.             ==NULL)
  840.     {
  841.         free(*ufactors);
  842.         return GLU_OUT_OF_MEMORY;
  843.     }
  844.     ctrl=new_ctrl->geom_ctrl;
  845.     offset1=new_ctrl->geom_t_stride*vorder;
  846.     offset2=new_ctrl->geom_s_stride*uorder;
  847.     for(j=0;j<vfactor_cnt;j++)
  848.     {
  849.         *(*vfactors+j*3+1)=tmp_factor1=calc_factor(ctrl,vorder,
  850.             j*offset1,dim,tolerance,dim);
  851.         /* loop ufactor_cnt-1 times */
  852.         for(i=1;i<ufactor_cnt;i++)
  853.         {
  854.             tmp_factor2=calc_factor(ctrl,vorder,
  855.                 j*offset1+i*offset2,dim,tolerance,dim);
  856.             if(tmp_factor2>tmp_factor1)
  857.                 tmp_factor1=tmp_factor2;
  858.         }
  859.         /* last time for the opposite edge */
  860.         *(*vfactors+j*3+2)=tmp_factor2=calc_factor(ctrl,vorder,
  861.             j*offset1+i*offset2-new_ctrl->geom_s_stride,
  862.             dim,tolerance,dim);
  863.         if(tmp_factor2>tmp_factor1)
  864.             *(*vfactors+j*3)=tmp_factor2;
  865.         else
  866.             *(*vfactors+j*3)=tmp_factor1;
  867.     }
  868.     offset3=new_ctrl->geom_s_stride;
  869.     offset2=new_ctrl->geom_s_stride*uorder;
  870.     for(j=0;j<ufactor_cnt;j++)
  871.     {
  872.         *(*ufactors+j*3+1)=tmp_factor1=calc_factor(ctrl,uorder,
  873.             j*offset2,offset3,tolerance,dim);
  874.         /* loop vfactor_cnt-1 times */
  875.         for(i=1;i<vfactor_cnt;i++)
  876.         {
  877.             tmp_factor2=calc_factor(ctrl,uorder,
  878.                 j*offset2+i*offset1,offset3,tolerance,dim);
  879.             if(tmp_factor2>tmp_factor1)
  880.                 tmp_factor1=tmp_factor2;
  881.         }
  882.         /* last time for the opposite edge */
  883.         *(*ufactors+j*3+2)=tmp_factor2=calc_factor(ctrl,uorder,
  884.             j*offset2+i*offset1-new_ctrl->geom_t_stride,
  885.             offset3,tolerance,dim);
  886.         if(tmp_factor2>tmp_factor1)
  887.             *(*ufactors+j*3)=tmp_factor2;
  888.         else
  889.             *(*ufactors+j*3)=tmp_factor1;
  890.     }
  891.     return GL_NO_ERROR;
  892. }
  893.  
  894. static GLenum
  895. calc_sampling_param_3D(new_ctrl_type *new_ctrl, GLfloat tolerance, GLint dim,
  896.     GLint uorder, GLint vorder, GLint **ufactors, GLint **vfactors)
  897. {
  898.     GLfloat        *ctrl;
  899.     GLint        tmp_factor1,tmp_factor2;
  900.     GLint        ufactor_cnt,vfactor_cnt;
  901.     GLint        offset1,offset2,offset3;
  902.     GLint        i,j;
  903.  
  904.     ufactor_cnt=new_ctrl->s_bezier_cnt;
  905.     vfactor_cnt=new_ctrl->t_bezier_cnt;
  906.     if((*ufactors=(GLint *)malloc(sizeof(GLint)*ufactor_cnt*3))
  907.             ==NULL)
  908.     {
  909.         return GLU_OUT_OF_MEMORY;
  910.     }
  911.     if((*vfactors=(GLint *)malloc(sizeof(GLint)*vfactor_cnt*3))
  912.             ==NULL)
  913.     {
  914.         free(*ufactors);
  915.         return GLU_OUT_OF_MEMORY;
  916.     }
  917.     ctrl=new_ctrl->geom_ctrl;
  918.     offset1=new_ctrl->geom_t_stride*vorder;
  919.     offset2=new_ctrl->geom_s_stride*uorder;
  920.     for(j=0;j<vfactor_cnt;j++)
  921.     {
  922.         *(*vfactors+j*3+1)=tmp_factor1=calc_parametric_factor(ctrl,vorder,
  923.             j*offset1,dim,tolerance,dim);
  924.         /* loop ufactor_cnt-1 times */
  925.         for(i=1;i<ufactor_cnt;i++)
  926.         {
  927.             tmp_factor2=calc_parametric_factor(ctrl,vorder,
  928.                 j*offset1+i*offset2,dim,tolerance,dim);
  929.             if(tmp_factor2>tmp_factor1)
  930.                 tmp_factor1=tmp_factor2;
  931.         }
  932.         /* last time for the opposite edge */
  933.         *(*vfactors+j*3+2)=tmp_factor2=calc_parametric_factor(ctrl,vorder,
  934.             j*offset1+i*offset2-new_ctrl->geom_s_stride,
  935.             dim,tolerance,dim);
  936.         if(tmp_factor2>tmp_factor1)
  937.             *(*vfactors+j*3)=tmp_factor2;
  938.         else
  939.             *(*vfactors+j*3)=tmp_factor1;
  940.     }
  941.     offset3=new_ctrl->geom_s_stride;
  942.     offset2=new_ctrl->geom_s_stride*uorder;
  943.     for(j=0;j<ufactor_cnt;j++)
  944.     {
  945.         *(*ufactors+j*3+1)=tmp_factor1=calc_parametric_factor(ctrl,uorder,
  946.             j*offset2,offset3,tolerance,dim);
  947.         /* loop vfactor_cnt-1 times */
  948.         for(i=1;i<vfactor_cnt;i++)
  949.         {
  950.             tmp_factor2=calc_parametric_factor(ctrl,uorder,
  951.                 j*offset2+i*offset1,offset3,tolerance,dim);
  952.             if(tmp_factor2>tmp_factor1)
  953.                 tmp_factor1=tmp_factor2;
  954.         }
  955.         /* last time for the opposite edge */
  956.         *(*ufactors+j*3+2)=tmp_factor2=calc_parametric_factor(ctrl,uorder,
  957.             j*offset2+i*offset1-new_ctrl->geom_t_stride,
  958.             offset3,tolerance,dim);
  959.         if(tmp_factor2>tmp_factor1)
  960.             *(*ufactors+j*3)=tmp_factor2;
  961.         else
  962.             *(*ufactors+j*3)=tmp_factor1;
  963.     }
  964.     return GL_NO_ERROR;
  965. }
  966.  
  967. static GLenum
  968. calc_sampling_2D(GLfloat *ctrl, GLint cnt, GLint order,
  969.     GLfloat tolerance, GLint dim, GLint **factors)
  970. {
  971.     GLint        factor_cnt;
  972.     GLint        tmp_factor;
  973.     GLint        offset;
  974.     GLint        i;
  975.  
  976.     factor_cnt=cnt/order;
  977.     if((*factors=(GLint *)malloc(sizeof(GLint)*factor_cnt))==NULL)
  978.     {
  979.         return GLU_OUT_OF_MEMORY;
  980.     }
  981.     offset=order*dim;
  982.     for(i=0;i<factor_cnt;i++)
  983.     {
  984.         tmp_factor=calc_factor(ctrl,order,i*offset,dim,tolerance,dim);
  985.         if(tmp_factor == 0)
  986.             (*factors)[i]=1;
  987.         else
  988.             (*factors)[i]=tmp_factor;
  989.     }
  990.     return GL_NO_ERROR;
  991. }
  992.  
  993. static void
  994. set_sampling_and_culling( GLUnurbsObj *nobj )
  995. {
  996.     if(nobj->auto_load_matrix==GL_FALSE)
  997.     {
  998.         GLint i;
  999.         GLfloat m[4];
  1000.  
  1001.         glPushAttrib( (GLbitfield) (GL_VIEWPORT_BIT | GL_TRANSFORM_BIT));
  1002.         for(i=0;i<4;i++)
  1003.             m[i]=nobj->sampling_matrices.viewport[i];
  1004.         glViewport(m[0],m[1],m[2],m[3]);
  1005.         glMatrixMode(GL_PROJECTION);
  1006.         glPushMatrix();
  1007.         glLoadMatrixf(nobj->sampling_matrices.proj);
  1008.         glMatrixMode(GL_MODELVIEW);
  1009.         glPushMatrix();
  1010.         glLoadMatrixf(nobj->sampling_matrices.model);
  1011.     }
  1012. }
  1013.  
  1014. static void
  1015. revert_sampling_and_culling( GLUnurbsObj *nobj )
  1016. {
  1017.     if(nobj->auto_load_matrix==GL_FALSE)
  1018.     {
  1019.         glMatrixMode(GL_MODELVIEW);
  1020.         glPopMatrix();
  1021.         glMatrixMode(GL_PROJECTION);
  1022.         glPopMatrix();
  1023.         glPopAttrib();
  1024.     }
  1025. }
  1026.  
  1027. GLenum
  1028. glu_do_sampling_3D( GLUnurbsObj *nobj, new_ctrl_type *new_ctrl,
  1029.     GLint **sfactors, GLint **tfactors)
  1030. {
  1031.     GLint    dim;
  1032.     GLenum    err;
  1033.  
  1034.     *sfactors=NULL;
  1035.     *tfactors=NULL;
  1036.     dim=nobj->surface.geom.dim;
  1037.     set_sampling_and_culling(nobj);
  1038.     if((err=calc_sampling_3D(new_ctrl,nobj->sampling_tolerance,dim,
  1039.         nobj->surface.geom.sorder,nobj->surface.geom.torder,
  1040.         sfactors,tfactors))==GLU_ERROR)
  1041.     {
  1042.         revert_sampling_and_culling(nobj);
  1043.         call_user_error(nobj,err);
  1044.         return GLU_ERROR;
  1045.     }
  1046.     revert_sampling_and_culling(nobj);
  1047.     return GLU_NO_ERROR;
  1048. }
  1049.  
  1050. GLenum
  1051. glu_do_sampling_uv( GLUnurbsObj *nobj, new_ctrl_type *new_ctrl,
  1052.     GLint **sfactors, GLint **tfactors)
  1053. {
  1054.     GLint    s_cnt, t_cnt, i;
  1055.     GLint    u_steps, v_steps;
  1056.  
  1057.     s_cnt = new_ctrl->s_bezier_cnt;
  1058.     t_cnt = new_ctrl->t_bezier_cnt;
  1059.     *sfactors=NULL;
  1060.     *tfactors=NULL;
  1061.     if((*sfactors=(GLint *)malloc(sizeof(GLint)*s_cnt*3))
  1062.             ==NULL)
  1063.     {
  1064.         return GLU_OUT_OF_MEMORY;
  1065.     }
  1066.     if((*tfactors=(GLint *)malloc(sizeof(GLint)*t_cnt*3))
  1067.             ==NULL)
  1068.     {
  1069.         free(*sfactors);
  1070.         return GLU_OUT_OF_MEMORY;
  1071.     }
  1072.     u_steps = nobj->u_step;
  1073.     v_steps = nobj->v_step;
  1074.     for(i=0; i<s_cnt; i++)
  1075.     {
  1076.         *(*sfactors+i*3) = u_steps;
  1077.         *(*sfactors+i*3+1) = u_steps;
  1078.         *(*sfactors+i*3+2) = u_steps;
  1079.     }
  1080.     for(i=0; i<t_cnt; i++)
  1081.     {
  1082.         *(*tfactors+i*3) = v_steps;
  1083.         *(*tfactors+i*3+1) = v_steps;
  1084.         *(*tfactors+i*3+2) = v_steps;
  1085.     }
  1086.     return GLU_NO_ERROR;
  1087. }
  1088.  
  1089. GLenum
  1090. glu_do_sampling_param_3D( GLUnurbsObj *nobj, new_ctrl_type *new_ctrl,
  1091.     GLint **sfactors, GLint **tfactors)
  1092. {
  1093.     GLint    dim;
  1094.     GLenum    err;
  1095.  
  1096.     *sfactors=NULL;
  1097.     *tfactors=NULL;
  1098.     dim=nobj->surface.geom.dim;
  1099.     set_sampling_and_culling(nobj);
  1100.     if((err=calc_sampling_param_3D(new_ctrl,nobj->parametric_tolerance,dim,
  1101.         nobj->surface.geom.sorder,nobj->surface.geom.torder,
  1102.         sfactors,tfactors))==GLU_ERROR)
  1103.     {
  1104.         revert_sampling_and_culling(nobj);
  1105.         call_user_error(nobj,err);
  1106.         return GLU_ERROR;
  1107.     }
  1108.     revert_sampling_and_culling(nobj);
  1109.     return GLU_NO_ERROR;
  1110. }
  1111.  
  1112. GLenum
  1113. glu_do_sampling_2D( GLUnurbsObj *nobj, GLfloat *ctrl, GLint cnt, GLint order,
  1114.     GLint dim, GLint **factors)
  1115. {
  1116.     GLenum err;
  1117.  
  1118.     set_sampling_and_culling(nobj);
  1119.     err=calc_sampling_2D(ctrl,cnt,order,nobj->sampling_tolerance,dim,
  1120.         factors);
  1121.     revert_sampling_and_culling(nobj);
  1122.     return err;
  1123. }
  1124.  
  1125.  
  1126. GLenum
  1127. glu_do_sampling_u( GLUnurbsObj *nobj, GLfloat *ctrl, GLint cnt, GLint order,
  1128.     GLint dim, GLint **factors)
  1129. {
  1130.     GLint    i;
  1131.     GLint    u_steps;
  1132.  
  1133.     cnt /= order;
  1134.     if((*factors=(GLint *)malloc(sizeof(GLint)*cnt))
  1135.             ==NULL)
  1136.     {
  1137.         return GLU_OUT_OF_MEMORY;
  1138.     }
  1139.     u_steps = nobj->u_step;
  1140.     for(i=0; i<cnt; i++)
  1141.         (*factors)[i] = u_steps;
  1142.     return GLU_NO_ERROR;
  1143. }
  1144.  
  1145. GLenum
  1146. glu_do_sampling_param_2D( GLUnurbsObj *nobj, GLfloat *ctrl, GLint cnt,
  1147.     GLint order, GLint dim, GLint **factors)
  1148. {
  1149.     GLint    i;
  1150.     GLint    u_steps;
  1151.     GLfloat    tolerance;
  1152.  
  1153.     set_sampling_and_culling(nobj);
  1154.     tolerance = nobj->parametric_tolerance;
  1155.     cnt /= order;
  1156.     if((*factors=(GLint *)malloc(sizeof(GLint)*cnt))
  1157.             ==NULL)
  1158.     {
  1159.         revert_sampling_and_culling(nobj);
  1160.         return GLU_OUT_OF_MEMORY;
  1161.     }
  1162.     u_steps = nobj->u_step;
  1163.     for(i=0; i<cnt; i++)
  1164.     {
  1165.         (*factors)[i] = calc_parametric_factor(ctrl,order,0,
  1166.             dim,tolerance,dim);
  1167.  
  1168.     }
  1169.     revert_sampling_and_culling(nobj);
  1170.     return GLU_NO_ERROR;
  1171. }
  1172.  
  1173. GLenum
  1174. glu_do_sampling_crv( GLUnurbsObj *nobj, GLfloat *ctrl, GLint cnt, GLint order,
  1175.     GLint dim, GLint **factors)
  1176. {
  1177.     GLenum err;
  1178.  
  1179.     *factors=NULL;
  1180.     switch(nobj->sampling_method)
  1181.     {
  1182.         case GLU_PATH_LENGTH:
  1183.             if((err=glu_do_sampling_2D(nobj,ctrl,cnt,order,dim,factors))!=
  1184.                 GLU_NO_ERROR)
  1185.             {
  1186.                 call_user_error(nobj,err);
  1187.                 return GLU_ERROR;
  1188.             }
  1189.             break;
  1190.         case GLU_DOMAIN_DISTANCE:
  1191.             if((err=glu_do_sampling_u(nobj,ctrl,cnt,order,dim,factors))!=
  1192.                 GLU_NO_ERROR)
  1193.             {
  1194.                 call_user_error(nobj,err);
  1195.                 return GLU_ERROR;
  1196.             }
  1197.             break;
  1198.         case GLU_PARAMETRIC_ERROR:
  1199.             if((err=glu_do_sampling_param_2D(nobj,ctrl,cnt,order,dim,factors))!=
  1200.                 GLU_NO_ERROR)
  1201.             {
  1202.                 call_user_error(nobj,err);
  1203.                 return GLU_ERROR;
  1204.             }
  1205.             break;
  1206.         default:
  1207.             abort();
  1208.     }
  1209.  
  1210.     return GLU_NO_ERROR;
  1211. }
  1212.  
  1213. /* TODO - i don't like this culling - this one just tests if at least one */
  1214. /* ctrl point lies within the viewport . Also the point_in_viewport() */
  1215. /* should be included in the fnctions for efficiency reasons */
  1216.  
  1217. static GLboolean
  1218. point_in_viewport(GLfloat *pt, GLint dim)
  1219. {
  1220.     GLdouble model[16],proj[16];
  1221.     GLint viewport[4];
  1222.     GLdouble x,y,z,w,winx,winy,winz;
  1223.  
  1224.     glGetDoublev(GL_MODELVIEW_MATRIX,model);
  1225.     glGetDoublev(GL_PROJECTION_MATRIX,proj);
  1226.     glGetIntegerv(GL_VIEWPORT,viewport);
  1227.     if(dim==3)
  1228.     {
  1229.         x=(GLdouble)pt[0];
  1230.         y=(GLdouble)pt[1];
  1231.         z=(GLdouble)pt[2];
  1232.         gluProject(x,y,z,model,proj,viewport,&winx,&winy,&winz);
  1233.     }
  1234.     else
  1235.     {
  1236.         w=(GLdouble)pt[3];
  1237.         x=(GLdouble)pt[0]/w;
  1238.         y=(GLdouble)pt[1]/w;
  1239.         z=(GLdouble)pt[2]/w;
  1240.         gluProject(x,y,z,model,proj,viewport,&winx,&winy,&winz);
  1241.     }
  1242.     if((GLint)winx >= viewport[0] && (GLint)winx < viewport[2] &&
  1243.             (GLint)winy >= viewport[1] && (GLint)winy < viewport[3])
  1244.         return GL_TRUE;
  1245.     return GL_FALSE;
  1246. }
  1247.  
  1248. GLboolean
  1249. fine_culling_test_3D(GLUnurbsObj *nobj,GLfloat *pts,GLint s_cnt,GLint t_cnt,
  1250.     GLint s_stride,GLint t_stride, GLint dim)
  1251. {
  1252.     GLint     i,j;
  1253.  
  1254.     if(nobj->culling==GL_FALSE)
  1255.         return GL_FALSE;
  1256.     set_sampling_and_culling(nobj);
  1257.     
  1258.     if(dim==3)
  1259.     {
  1260.         for(i=0;i<s_cnt;i++)
  1261.             for(j=0;j<t_cnt;j++)
  1262.                 if(point_in_viewport(pts+i*s_stride+j*t_stride,dim))
  1263.                 {
  1264.                     revert_sampling_and_culling(nobj);
  1265.                     return GL_FALSE;
  1266.                 }
  1267.     }
  1268.     else
  1269.     {
  1270.         for(i=0;i<s_cnt;i++)
  1271.             for(j=0;j<t_cnt;j++)
  1272.                 if(point_in_viewport(pts+i*s_stride+j*t_stride,dim))
  1273.                 {
  1274.                     revert_sampling_and_culling(nobj);
  1275.                     return GL_FALSE;
  1276.                 }
  1277.     }
  1278.     revert_sampling_and_culling(nobj);
  1279.     return GL_TRUE;
  1280. }
  1281.  
  1282. /*GLboolean
  1283. fine_culling_test_3D(GLUnurbsObj *nobj,GLfloat *pts,GLint s_cnt,GLint t_cnt,
  1284.     GLint s_stride,GLint t_stride, GLint dim)
  1285. {
  1286.     GLint        visible_cnt;
  1287.     GLfloat        feedback_buffer[5];
  1288.     GLsizei        buffer_size;
  1289.     GLint         i,j;
  1290.  
  1291.     if(nobj->culling==GL_FALSE)
  1292.         return GL_FALSE;
  1293.     buffer_size=5;
  1294.     set_sampling_and_culling(nobj);
  1295.     
  1296.     glFeedbackBuffer(buffer_size,GL_2D,feedback_buffer);
  1297.     glRenderMode(GL_FEEDBACK);
  1298.     if(dim==3)
  1299.     {
  1300.         for(i=0;i<s_cnt;i++)
  1301.         {
  1302.             glBegin(GL_LINE_LOOP);
  1303.             for(j=0;j<t_cnt;j++)
  1304.                 glVertex3fv(pts+i*s_stride+j*t_stride);
  1305.             glEnd();
  1306.         }
  1307.         for(j=0;j<t_cnt;j++)
  1308.         {
  1309.             glBegin(GL_LINE_LOOP);
  1310.             for(i=0;i<s_cnt;i++)
  1311.                 glVertex3fv(pts+i*s_stride+j*t_stride);
  1312.             glEnd();
  1313.         }
  1314.     }
  1315.     else
  1316.     {
  1317.         for(i=0;i<s_cnt;i++)
  1318.         {
  1319.             glBegin(GL_LINE_LOOP);
  1320.             for(j=0;j<t_cnt;j++)
  1321.                 glVertex4fv(pts+i*s_stride+j*t_stride);
  1322.             glEnd();
  1323.         }
  1324.         for(j=0;j<t_cnt;j++)
  1325.         {
  1326.             glBegin(GL_LINE_LOOP);
  1327.             for(i=0;i<s_cnt;i++)
  1328.                 glVertex4fv(pts+i*s_stride+j*t_stride);
  1329.             glEnd();
  1330.         }
  1331.     }
  1332.     visible_cnt=glRenderMode(GL_RENDER);
  1333.  
  1334.     revert_sampling_and_culling(nobj);
  1335.     return (GLboolean)(visible_cnt==0);
  1336. }*/
  1337.  
  1338. GLboolean
  1339. fine_culling_test_2D(GLUnurbsObj *nobj,GLfloat *pts,GLint cnt,
  1340.     GLint stride, GLint dim)
  1341. {
  1342.     GLint     i;
  1343.  
  1344.     if(nobj->culling==GL_FALSE)
  1345.         return GL_FALSE;
  1346.     set_sampling_and_culling(nobj);
  1347.     
  1348.     if(dim==3)
  1349.     {
  1350.         for(i=0;i<cnt;i++)
  1351.             if(point_in_viewport(pts+i*stride,dim))
  1352.             {
  1353.                 revert_sampling_and_culling(nobj);
  1354.                 return GL_FALSE;
  1355.             }
  1356.     }
  1357.     else
  1358.     {
  1359.         for(i=0;i<cnt;i++)
  1360.             if(point_in_viewport(pts+i*stride,dim))
  1361.             {
  1362.                 revert_sampling_and_culling(nobj);
  1363.                 return GL_FALSE;
  1364.             }
  1365.     }
  1366.     revert_sampling_and_culling(nobj);
  1367.     return GL_TRUE;
  1368. }
  1369.  
  1370. /*GLboolean
  1371. fine_culling_test_2D(GLUnurbsObj *nobj,GLfloat *pts,GLint cnt,
  1372.     GLint stride, GLint dim)
  1373. {
  1374.     GLint        visible_cnt;
  1375.     GLfloat        feedback_buffer[5];
  1376.     GLsizei        buffer_size;
  1377.     GLint         i;
  1378.  
  1379.     if(nobj->culling==GL_FALSE)
  1380.         return GL_FALSE;
  1381.     buffer_size=5;
  1382.     set_sampling_and_culling(nobj);
  1383.     
  1384.     glFeedbackBuffer(buffer_size,GL_2D,feedback_buffer);
  1385.     glRenderMode(GL_FEEDBACK);
  1386.     glBegin(GL_LINE_LOOP);
  1387.     if(dim==3)
  1388.     {
  1389.         for(i=0;i<cnt;i++)
  1390.             glVertex3fv(pts+i*stride);
  1391.     }
  1392.     else
  1393.     {
  1394.         for(i=0;i<cnt;i++)
  1395.             glVertex4fv(pts+i*stride);
  1396.     }
  1397.     glEnd();
  1398.     visible_cnt=glRenderMode(GL_RENDER);
  1399.  
  1400.     revert_sampling_and_culling(nobj);
  1401.     return (GLboolean)(visible_cnt==0);
  1402. }*/
  1403.  
  1404.