src/common/cone.c

/* Copyright (c) 1991 Regents of the University of California */

#ifndef lint
static char SCCSid[] = "$SunId$ LBL";
#endif

/*
 *  cone.c - routines for making cones
 *
 *     2/12/86
 */

#include  "standard.h"

#include  "object.h"

#include  "otypes.h"

#include  "cone.h"

/*
 *     In general, a cone may be any one of a cone, a cylinder, a ring,
 *  a cup (inverted cone), or a tube (inverted cylinder).
 *     Most cones are specified with a starting point and radius and
 *  an ending point and radius.  In the cases of a cylinder or tube,
 *  only one radius is needed.  In the case of a ring, a normal direction
 *  is specified instead of a second endpoint.
 *
 *      mtype (cone|cup) name
 *      0
 *      0
 *      8 P0x P0y P0z P1x P1y P1z R0 R1
 *
 *      mtype (cylinder|tube) name
 *      0
 *      0
 *      7 P0x P0y P0z P1x P1y P1z R
 *
 *      mtype ring name
 *      0
 *      0
 *      8 Px Py Pz Nx Ny Nz R0 R1
 */


CONE *
getcone(o, getxf)                       /* get cone structure */
register OBJREC  *o;
int  getxf;
{
        double  fabs(), sqrt();
        int  sgn0, sgn1;
        register CONE  *co;

        if ((co = (CONE *)o->os) == NULL) {

                co = (CONE *)malloc(sizeof(CONE));
                if (co == NULL)
                        error(SYSTEM, "out of memory in makecone");

                co->ca = o->oargs.farg;
                                                /* get radii */
                if (o->otype == OBJ_CYLINDER || o->otype == OBJ_TUBE) {
                        if (o->oargs.nfargs != 7)
                                goto argcerr;
                        if (co->ca[6] < -FTINY) {
                                objerror(o, WARNING, "negative radius");
                                o->otype = o->otype == OBJ_CYLINDER ?
                                                OBJ_TUBE : OBJ_CYLINDER;
                                co->ca[6] = -co->ca[6];
                        } else if (co->ca[6] <= FTINY)
                                goto raderr;
                        co->r0 = co->r1 = 6;
                } else {
                        if (o->oargs.nfargs != 8)
                                goto argcerr;
                        if (co->ca[6] < -FTINY) sgn0 = -1;
                        else if (co->ca[6] > FTINY) sgn0 = 1;
                        else sgn0 = 0;
                        if (co->ca[7] < -FTINY) sgn1 = -1;
                        else if (co->ca[7] > FTINY) sgn1 = 1;
                        else sgn1 = 0;
                        if (sgn0+sgn1 == 0)
                                goto raderr;
                        if (sgn0 < 0 || sgn1 < 0) {
                                objerror(o, o->otype==OBJ_RING?USER:WARNING,
                                        "negative radii");
                                o->otype = o->otype == OBJ_CONE ?
                                                OBJ_CUP : OBJ_CONE;
                        }
                        co->ca[6] = co->ca[6]*sgn0;
                        co->ca[7] = co->ca[7]*sgn1;
                        co->r0 = 6;
                        co->r1 = 7;
                        if (fabs(co->ca[7] - co->ca[6]) <= FTINY) {
                                if (o->otype == OBJ_RING)
                                        goto raderr;
                                o->otype = o->otype == OBJ_CONE ?
                                                OBJ_CYLINDER : OBJ_TUBE;
                                o->oargs.nfargs = 7;
                                co->r1 = 6;
                        }
                }
                                                /* get axis orientation */
                co->p0 = 0;
                if (o->otype == OBJ_RING) {
                        if (co->ca[6] > co->ca[7]) {    /* make r0 smaller */
                                co->r0 = 7;
                                co->r1 = 6;
                        }
                        co->p1 = 0;
                        VCOPY(co->ad, o->oargs.farg+3);
                } else {
                        co->p1 = 3;
                        co->ad[0] = co->ca[3] - co->ca[0];
                        co->ad[1] = co->ca[4] - co->ca[1];
                        co->ad[2] = co->ca[5] - co->ca[2];
                }
                co->al = normalize(co->ad);
                if (co->al == 0.0)
                        objerror(o, USER, "zero orientation");
                                        /* compute axis and side lengths */
                if (o->otype == OBJ_RING) {
                        co->al = 0.0;
                        co->sl = co->ca[co->r1] - co->ca[co->r0];
                } else if (o->otype == OBJ_CONE || o->otype == OBJ_CUP) {
                        co->sl = co->ca[7] - co->ca[6];
                        co->sl = sqrt(co->sl*co->sl + co->al*co->al);
                } else { /* OBJ_CYLINDER || OBJ_TUBE */
                        co->sl = co->al;
                }
                co->tm = NULL;
                o->os = (char *)co;
        }
        if (getxf && co->tm == NULL)
                conexform(co);
        return(co);

argcerr:
        objerror(o, USER, "bad # arguments");
raderr:
        objerror(o, USER, "illegal radii");
}


freecone(o)                     /* free memory associated with cone */
OBJREC  *o;
{
        register CONE  *co = (CONE *)o->os;

        if (o->os == NULL)
                return;
        if (co->tm != NULL)
                free((char *)co->tm);
        free(o->os);
        o->os = NULL;
}


conexform(co)                   /* get cone transformation matrix */
register CONE  *co;
{
        double  sqrt(), fabs();
        double  m4[4][4];
        register double  d;
        register int  i;

        co->tm = (double (*)[4])malloc(sizeof(m4));
        if (co->tm == NULL)
                error(SYSTEM, "out of memory in conexform");

                                /* translate to origin */
        setident4(co->tm);
        if (co->r0 == co->r1)
                d = 0.0;
        else
                d = co->ca[co->r0] / (co->ca[co->r1] - co->ca[co->r0]);
        for (i = 0; i < 3; i++)
                co->tm[3][i] = d*(co->ca[co->p1+i] - co->ca[co->p0+i])
                                - co->ca[co->p0+i];
        
                                /* rotate to positive z-axis */
        setident4(m4);
        d = co->ad[1]*co->ad[1] + co->ad[2]*co->ad[2];
        if (d <= FTINY*FTINY) {
                m4[0][0] = 0.0;
                m4[0][2] = co->ad[0];
                m4[2][0] = -co->ad[0];
                m4[2][2] = 0.0;
        } else {
                d = sqrt(d);
                m4[0][0] = d;
                m4[1][0] = -co->ad[0]*co->ad[1]/d;
                m4[2][0] = -co->ad[0]*co->ad[2]/d;
                m4[1][1] = co->ad[2]/d;
                m4[2][1] = -co->ad[1]/d;
                m4[0][2] = co->ad[0];
                m4[1][2] = co->ad[1];
                m4[2][2] = co->ad[2];
        }
        multmat4(co->tm, co->tm, m4);

                                /* scale z-axis */
        setident4(m4);
        if (co->p0 != co->p1 && co->r0 != co->r1) {
                d = fabs(co->ca[co->r1] - co->ca[co->r0]);
                m4[2][2] = d/co->al;
        }
        multmat4(co->tm, co->tm, m4);
}
Revision:	1.5
Committed:	Wed Jul 17 14:10:16 1991 UTC (32 years, 9 months ago) by greg
Content type:	text/plain
Branch:	MAIN
Changes since 1.4:	+2 -0 lines
Log Message:	made call to free structures more robust
#	User	Rev	Content
1	greg	1.3	/* Copyright (c) 1991 Regents of the University of California */
2	greg	1.1
3			#ifndef lint
4			static char SCCSid[] = "$SunId$ LBL";
5			#endif
6
7			/*
8			* cone.c - routines for making cones
9			*
10			* 2/12/86
11			*/
12
13			#include "standard.h"
14
15			#include "object.h"
16
17			#include "otypes.h"
18
19			#include "cone.h"
20
21			/*
22			* In general, a cone may be any one of a cone, a cylinder, a ring,
23			* a cup (inverted cone), or a tube (inverted cylinder).
24			* Most cones are specified with a starting point and radius and
25			* an ending point and radius. In the cases of a cylinder or tube,
26			* only one radius is needed. In the case of a ring, a normal direction
27			* is specified instead of a second endpoint.
28			*
29			* mtype (cone\|cup) name
30			* 0
31			* 0
32			* 8 P0x P0y P0z P1x P1y P1z R0 R1
33			*
34			* mtype (cylinder\|tube) name
35			* 0
36			* 0
37			* 7 P0x P0y P0z P1x P1y P1z R
38			*
39			* mtype ring name
40			* 0
41			* 0
42			* 8 Px Py Pz Nx Ny Nz R0 R1
43			*/
44
45
46			CONE *
47			getcone(o, getxf) /* get cone structure */
48			register OBJREC *o;
49			int getxf;
50			{
51			double fabs(), sqrt();
52	greg	1.4	int sgn0, sgn1;
53	greg	1.1	register CONE *co;
54
55			if ((co = (CONE *)o->os) == NULL) {
56
57			co = (CONE *)malloc(sizeof(CONE));
58			if (co == NULL)
59			error(SYSTEM, "out of memory in makecone");
60
61			co->ca = o->oargs.farg;
62			/* get radii */
63			if (o->otype == OBJ_CYLINDER \|\| o->otype == OBJ_TUBE) {
64			if (o->oargs.nfargs != 7)
65			goto argcerr;
66	greg	1.4	if (co->ca[6] < -FTINY) {
67			objerror(o, WARNING, "negative radius");
68			o->otype = o->otype == OBJ_CYLINDER ?
69			OBJ_TUBE : OBJ_CYLINDER;
70			co->ca[6] = -co->ca[6];
71			} else if (co->ca[6] <= FTINY)
72	greg	1.1	goto raderr;
73			co->r0 = co->r1 = 6;
74			} else {
75			if (o->oargs.nfargs != 8)
76			goto argcerr;
77	greg	1.4	if (co->ca[6] < -FTINY) sgn0 = -1;
78			else if (co->ca[6] > FTINY) sgn0 = 1;
79			else sgn0 = 0;
80			if (co->ca[7] < -FTINY) sgn1 = -1;
81			else if (co->ca[7] > FTINY) sgn1 = 1;
82			else sgn1 = 0;
83			if (sgn0+sgn1 == 0)
84	greg	1.1	goto raderr;
85	greg	1.4	if (sgn0 < 0 \|\| sgn1 < 0) {
86			objerror(o, o->otype==OBJ_RING?USER:WARNING,
87			"negative radii");
88			o->otype = o->otype == OBJ_CONE ?
89			OBJ_CUP : OBJ_CONE;
90			}
91			co->ca[6] = co->ca[6]*sgn0;
92			co->ca[7] = co->ca[7]*sgn1;
93	greg	1.1	co->r0 = 6;
94			co->r1 = 7;
95	greg	1.4	if (fabs(co->ca[7] - co->ca[6]) <= FTINY) {
96			if (o->otype == OBJ_RING)
97			goto raderr;
98			o->otype = o->otype == OBJ_CONE ?
99			OBJ_CYLINDER : OBJ_TUBE;
100			o->oargs.nfargs = 7;
101			co->r1 = 6;
102			}
103	greg	1.1	}
104			/* get axis orientation */
105			co->p0 = 0;
106			if (o->otype == OBJ_RING) {
107			if (co->ca[6] > co->ca[7]) { /* make r0 smaller */
108			co->r0 = 7;
109			co->r1 = 6;
110			}
111			co->p1 = 0;
112			VCOPY(co->ad, o->oargs.farg+3);
113			} else {
114			co->p1 = 3;
115			co->ad[0] = co->ca[3] - co->ca[0];
116			co->ad[1] = co->ca[4] - co->ca[1];
117			co->ad[2] = co->ca[5] - co->ca[2];
118			}
119			co->al = normalize(co->ad);
120			if (co->al == 0.0)
121			objerror(o, USER, "zero orientation");
122			/* compute axis and side lengths */
123			if (o->otype == OBJ_RING) {
124			co->al = 0.0;
125			co->sl = co->ca[co->r1] - co->ca[co->r0];
126			} else if (o->otype == OBJ_CONE \|\| o->otype == OBJ_CUP) {
127			co->sl = co->ca[7] - co->ca[6];
128			co->sl = sqrt(co->slco->sl + co->alco->al);
129			} else { /* OBJ_CYLINDER \|\| OBJ_TUBE */
130			co->sl = co->al;
131			}
132			co->tm = NULL;
133	greg	1.2	o->os = (char *)co;
134	greg	1.1	}
135			if (getxf && co->tm == NULL)
136			conexform(co);
137			return(co);
138
139			argcerr:
140			objerror(o, USER, "bad # arguments");
141			raderr:
142			objerror(o, USER, "illegal radii");
143			}
144
145
146			freecone(o) /* free memory associated with cone */
147			OBJREC *o;
148			{
149			register CONE co = (CONE )o->os;
150
151	greg	1.5	if (o->os == NULL)
152			return;
153	greg	1.1	if (co->tm != NULL)
154			free((char *)co->tm);
155			free(o->os);
156			o->os = NULL;
157			}
158
159
160			conexform(co) /* get cone transformation matrix */
161			register CONE *co;
162			{
163			double sqrt(), fabs();
164			double m4[4][4];
165			register double d;
166			register int i;
167
168			co->tm = (double (*)[4])malloc(sizeof(m4));
169			if (co->tm == NULL)
170			error(SYSTEM, "out of memory in conexform");
171
172			/* translate to origin */
173			setident4(co->tm);
174			if (co->r0 == co->r1)
175			d = 0.0;
176			else
177			d = co->ca[co->r0] / (co->ca[co->r1] - co->ca[co->r0]);
178			for (i = 0; i < 3; i++)
179			co->tm[3][i] = d*(co->ca[co->p1+i] - co->ca[co->p0+i])
180			- co->ca[co->p0+i];
181
182			/* rotate to positive z-axis */
183			setident4(m4);
184			d = co->ad[1]co->ad[1] + co->ad[2]co->ad[2];
185			if (d <= FTINY*FTINY) {
186			m4[0][0] = 0.0;
187			m4[0][2] = co->ad[0];
188			m4[2][0] = -co->ad[0];
189			m4[2][2] = 0.0;
190			} else {
191			d = sqrt(d);
192			m4[0][0] = d;
193			m4[1][0] = -co->ad[0]*co->ad[1]/d;
194			m4[2][0] = -co->ad[0]*co->ad[2]/d;
195			m4[1][1] = co->ad[2]/d;
196			m4[2][1] = -co->ad[1]/d;
197			m4[0][2] = co->ad[0];
198			m4[1][2] = co->ad[1];
199			m4[2][2] = co->ad[2];
200			}
201			multmat4(co->tm, co->tm, m4);
202
203			/* scale z-axis */
204			setident4(m4);
205			if (co->p0 != co->p1 && co->r0 != co->r1) {
206			d = fabs(co->ca[co->r1] - co->ca[co->r0]);
207			m4[2][2] = d/co->al;
208			}
209			multmat4(co->tm, co->tm, m4);
210			}