Back to index

radiance  4R0+20100331
cone.c
Go to the documentation of this file.
00001 #ifndef lint
00002 static const char    RCSid[] = "$Id: cone.c,v 2.9 2003/07/27 22:12:01 schorsch Exp $";
00003 #endif
00004 /*
00005  *  cone.c - routines for making cones
00006  */
00007 
00008 #include "copyright.h"
00009 
00010 #include  "standard.h"
00011 
00012 #include  "object.h"
00013 
00014 #include  "otypes.h"
00015 
00016 #include  "cone.h"
00017 
00018 /*
00019  *     In general, a cone may be any one of a cone, a cylinder, a ring,
00020  *  a cup (inverted cone), or a tube (inverted cylinder).
00021  *     Most cones are specified with a starting point and radius and
00022  *  an ending point and radius.  In the cases of a cylinder or tube,
00023  *  only one radius is needed.  In the case of a ring, a normal direction
00024  *  is specified instead of a second endpoint.
00025  *
00026  *     mtype (cone|cup) name
00027  *     0
00028  *     0
00029  *     8 P0x P0y P0z P1x P1y P1z R0 R1
00030  *
00031  *     mtype (cylinder|tube) name
00032  *     0
00033  *     0
00034  *     7 P0x P0y P0z P1x P1y P1z R
00035  *
00036  *     mtype ring name
00037  *     0
00038  *     0
00039  *     8 Px Py Pz Nx Ny Nz R0 R1
00040  */
00041 
00042 
00043 CONE *
00044 getcone(o, getxf)                  /* get cone structure */
00045 register OBJREC  *o;
00046 int  getxf;
00047 {
00048        int  sgn0, sgn1;
00049        register CONE  *co;
00050 
00051        if ((co = (CONE *)o->os) == NULL) {
00052 
00053               co = (CONE *)malloc(sizeof(CONE));
00054               if (co == NULL)
00055                      error(SYSTEM, "out of memory in makecone");
00056 
00057               co->ca = o->oargs.farg;
00058                                           /* get radii */
00059               if ((o->otype == OBJ_CYLINDER) | (o->otype == OBJ_TUBE)) {
00060                      if (o->oargs.nfargs != 7)
00061                             goto argcerr;
00062                      if (co->ca[6] < -FTINY) {
00063                             objerror(o, WARNING, "negative radius");
00064                             o->otype = o->otype == OBJ_CYLINDER ?
00065                                           OBJ_TUBE : OBJ_CYLINDER;
00066                             co->ca[6] = -co->ca[6];
00067                      } else if (co->ca[6] <= FTINY)
00068                             goto raderr;
00069                      co->p0 = 0; co->p1 = 3;
00070                      co->r0 = co->r1 = 6;
00071               } else {
00072                      if (o->oargs.nfargs != 8)
00073                             goto argcerr;
00074                      if (co->ca[6] < -FTINY) sgn0 = -1;
00075                      else if (co->ca[6] > FTINY) sgn0 = 1;
00076                      else sgn0 = 0;
00077                      if (co->ca[7] < -FTINY) sgn1 = -1;
00078                      else if (co->ca[7] > FTINY) sgn1 = 1;
00079                      else sgn1 = 0;
00080                      if (sgn0+sgn1 == 0)
00081                             goto raderr;
00082                      if ((sgn0 < 0) | (sgn1 < 0)) {
00083                             objerror(o, o->otype==OBJ_RING?USER:WARNING,
00084                                    "negative radii");
00085                             o->otype = o->otype == OBJ_CONE ?
00086                                           OBJ_CUP : OBJ_CONE;
00087                      }
00088                      co->ca[6] = co->ca[6]*sgn0;
00089                      co->ca[7] = co->ca[7]*sgn1;
00090                      if (co->ca[7] - co->ca[6] > FTINY) {
00091                             if (o->otype == OBJ_RING)
00092                                    co->p0 = co->p1 = 0;
00093                             else {
00094                                    co->p0 = 0; co->p1 = 3;
00095                             }
00096                             co->r0 = 6; co->r1 = 7;
00097                      } else if (co->ca[6] - co->ca[7] > FTINY) {
00098                             if (o->otype == OBJ_RING)
00099                                    co->p0 = co->p1 = 0;
00100                             else {
00101                                    co->p0 = 3; co->p1 = 0;
00102                             }
00103                             co->r0 = 7; co->r1 = 6;
00104                      } else {
00105                             if (o->otype == OBJ_RING)
00106                                    goto raderr;
00107                             o->otype = o->otype == OBJ_CONE ?
00108                                           OBJ_CYLINDER : OBJ_TUBE;
00109                             o->oargs.nfargs = 7;
00110                             co->p0 = 0; co->p1 = 3;
00111                             co->r0 = co->r1 = 6;
00112                      }
00113               }
00114                                           /* get axis orientation */
00115               if (o->otype == OBJ_RING)
00116                      VCOPY(co->ad, o->oargs.farg+3);
00117               else {
00118                      co->ad[0] = CO_P1(co)[0] - CO_P0(co)[0];
00119                      co->ad[1] = CO_P1(co)[1] - CO_P0(co)[1];
00120                      co->ad[2] = CO_P1(co)[2] - CO_P0(co)[2];
00121               }
00122               co->al = normalize(co->ad);
00123               if (co->al == 0.0)
00124                      objerror(o, USER, "zero orientation");
00125                                    /* compute axis and side lengths */
00126               if (o->otype == OBJ_RING) {
00127                      co->al = 0.0;
00128                      co->sl = CO_R1(co) - CO_R0(co);
00129               } else if ((o->otype == OBJ_CONE) | (o->otype == OBJ_CUP)) {
00130                      co->sl = co->ca[7] - co->ca[6];
00131                      co->sl = sqrt(co->sl*co->sl + co->al*co->al);
00132               } else { /* OBJ_CYLINDER or OBJ_TUBE */
00133                      co->sl = co->al;
00134               }
00135               co->tm = NULL;
00136               o->os = (char *)co;
00137        }
00138        if (getxf && co->tm == NULL)
00139               conexform(co);
00140        return(co);
00141 
00142 argcerr:
00143        objerror(o, USER, "bad # arguments");
00144 raderr:
00145        objerror(o, USER, "illegal radii");
00146        return NULL; /* pro forma return */
00147 }
00148 
00149 
00150 void
00151 freecone(o)                 /* free memory associated with cone */
00152 OBJREC  *o;
00153 {
00154        register CONE  *co = (CONE *)o->os;
00155 
00156        if (co == NULL)
00157               return;
00158        if (co->tm != NULL)
00159               free((void *)co->tm);
00160        free((void *)co);
00161        o->os = NULL;
00162 }
00163 
00164 
00165 void
00166 conexform(co)               /* get cone transformation matrix */
00167 register CONE  *co;
00168 {
00169        MAT4  m4;
00170        register double  d;
00171        register int  i;
00172 
00173        co->tm = (RREAL (*)[4])malloc(sizeof(MAT4));
00174        if (co->tm == NULL)
00175               error(SYSTEM, "out of memory in conexform");
00176 
00177                             /* translate to origin */
00178        setident4(co->tm);
00179        if (co->r0 == co->r1)
00180               d = 0.0;
00181        else
00182               d = CO_R0(co) / (CO_R1(co) - CO_R0(co));
00183        for (i = 0; i < 3; i++)
00184               co->tm[3][i] = d*(CO_P1(co)[i] - CO_P0(co)[i])
00185                             - CO_P0(co)[i];
00186        
00187                             /* rotate to positive z-axis */
00188        setident4(m4);
00189        d = co->ad[1]*co->ad[1] + co->ad[2]*co->ad[2];
00190        if (d <= FTINY*FTINY) {
00191               m4[0][0] = 0.0;
00192               m4[0][2] = co->ad[0];
00193               m4[2][0] = -co->ad[0];
00194               m4[2][2] = 0.0;
00195        } else {
00196               d = sqrt(d);
00197               m4[0][0] = d;
00198               m4[1][0] = -co->ad[0]*co->ad[1]/d;
00199               m4[2][0] = -co->ad[0]*co->ad[2]/d;
00200               m4[1][1] = co->ad[2]/d;
00201               m4[2][1] = -co->ad[1]/d;
00202               m4[0][2] = co->ad[0];
00203               m4[1][2] = co->ad[1];
00204               m4[2][2] = co->ad[2];
00205        }
00206        multmat4(co->tm, co->tm, m4);
00207 
00208                             /* scale z-axis */
00209        if ((co->p0 != co->p1) & (co->r0 != co->r1)) {
00210               setident4(m4);
00211               m4[2][2] = (CO_R1(co) - CO_R0(co)) / co->al;
00212               multmat4(co->tm, co->tm, m4);
00213        }
00214 }