|
|
1.1 root 1: #include "ideal.h"
2: #include "y.tab.h"
3:
4: boolean llinter (x1,y1, x2,y2, z1,w1, z2,w2, alpha, beta, collinear)
5: /* if this function returns TRUE,
6: /* then alpha[(x1,y1),(x2,y2)] = beta[(z1,w1),(z2,w2)] */
7: float x1,y1, x2,y2, z1,w1, z2,w2;
8: float *alpha;
9: float *beta;
10: boolean *collinear;
11: {
12: float A1, B1, C1, A2, B2, C2, D, x, y;
13: dprintf "(%f,%f) -- (%f,%f)\n", x1,y1, x2,y2);
14: dprintf "(%f,%f) -- (%f,%f)\n", z1,w1, z2,w2);
15: A1 = y1 - y2;
16: B1 = x2 - x1;
17: C1 = -B1*y1 - A1*x1;
18: A2 = w1 - w2;
19: B2 = z2 - z1;
20: C2 = -B2*w1 - A2*z1;
21: D = A1*B2 - A2*B1;
22: if (fabs(D) < EPSILON) {
23: *collinear = arecollinear(x1,y1,x2,y2,z1,w1);
24: dprintf "%s\n", (*collinear)?"coincident":"disjoint");
25: return (FALSE);
26: }
27: *collinear = FALSE;
28: x = (B1*C2 - B2*C1)/D;
29: y = (A2*C1 - A1*C2)/D;
30: if (fabs(x2 - x1) > EPSILON) {
31: *alpha = (x - x1)/(x2 - x1);
32: } else if (fabs(y2 - y1) > EPSILON) {
33: *alpha = (y - y1)/(y2 - y1);
34: } else fprintf (stderr, "ideal: llinter: can't happen\n");
35: if (fabs(z2 - z1) > EPSILON) {
36: *beta = (x - z1)/(z2 - z1);
37: } else if (fabs(w2 - w1) > EPSILON) {
38: *beta = (y - w1)/(w2 - w1);
39: } else fprintf (stderr, "ideal: llinter: can't happen\n");
40: dprintf "intersection alpha = %f; beta = %f\n", *alpha, *beta);
41: return (TRUE);
42: } /* llinter */
43:
44: boolean lcinter (x1,y1, x2,y2, x0,y0, r, alpha1,theta1, alpha2,theta2)
45: /* if this function returns TRUE,
46: /* then alpha1[(x1,y1),(x2,y2)] = (x0,y0) + r*cis(theta1)
47: /* and alpha2[(x1,y1),(x2,y2)] = (x0,y0) + r*cis(theta2) */
48: float x1,y1, x2,y2, x0,y0, r;
49: float *alpha1;
50: float *theta1;
51: float *alpha2;
52: float *theta2;
53: {
54: float dx1, dx2, dy1, dy2;
55: float A, B, C, D;
56:
57: dprintf "intersection parameters:\n");
58: dprintf "%f, %f -- %f, %f\n", x1, y1, x2, y2);
59: dprintf "%f, %f (%f)\n", x0, y0, r);
60: r = fabs(r);
61: dx1 = x1 - x0;
62: dx2 = x2 - x1;
63: dy1 = y1 - y0;
64: dy2 = y2 - y1;
65: A = dx2*dx2 + dy2*dy2;
66: dprintf "A=%f\n", A);
67: if (A < EPSILON) {
68: if (fabs (hypot (dx1, dy1) - r) < EPSILON) {
69: *alpha1 = *alpha2 = 0.0;
70: *theta1 = atan2 (dy1, dx1);
71: *theta2 = *theta1 = rprin (*theta1);
72: dprintf "alpha1 = alpha2 = %f theta1 = theta2 = %f\n",
73: *alpha1, *theta1);
74: return (TRUE);
75: }
76: else
77: return (FALSE);
78: }
79: B = 2*(dx1*dx2 + dy1*dy2);
80: C = dx1*dx1 + dy1*dy1 - r*r;
81: D = B*B - 4*A*C;
82: dprintf "B=%f C=%f D=%f\n", B, C, D);
83: if (D < -EPSILON)
84: return (FALSE);
85: if (fabs(D) < EPSILON)
86: D = 0.0;
87: D = sqrt(D);
88: *alpha1 = (-B + D)/(2.0*A);
89: *theta1 = rprin (atan2 (dy1 + *alpha1*dy2, dx1 + *alpha1*dx2));
90: *alpha2 = (-B - D)/(2.0*A);
91: *theta2 = rprin (atan2 (dy1 + *alpha2*dy2, dx1 + *alpha2*dx2));
92: dprintf "intersection alpha1 = %f, theta1 = %f\n", *alpha1, *theta1);
93: dprintf "intersection alpha2 = %f, theta2 = %f\n", *alpha2, *theta2);
94: return (TRUE);
95: }
96:
97: boolean ccinter (x0,y0,r0, x1,y1,r1, theta1,phi1, theta2,phi2)
98: /* if this function returns TRUE,
99: /* then (x0,y0) + r0*cis(theta1) = (x1,y1) + r1*cis(phi1)
100: /* and (x0,y0) + r0*cis(theta2) = (x1,y1) + r1*cis(phi2) */
101: float x0,y0,r0;
102: float x1,y1,r1;
103: float *theta1;
104: float *phi1;
105: float *theta2;
106: float *phi2;
107: {
108: float xcoeff, ycoeff, const;
109: float u1, v1, u2, v2;
110: boolean lncrc;
111:
112: dprintf "intersection parameters\n");
113: dprintf "%f %f (%f)\n", x0, y0, r0);
114: dprintf "%f %f (%f)\n", x1, y1, r1);
115:
116: r0 = fabs(r0);
117: r1 = fabs(r1);
118: xcoeff = 2*(x1 - x0);
119: ycoeff = 2*(y1 - y0);
120: const = r0*r0 - x0*x0 - y0*y0 - r1*r1 + x1*x1 + y1*y1;
121: if (fabs(xcoeff) < EPSILON && fabs(ycoeff) < EPSILON)
122: return (FALSE);
123: if (fabs(xcoeff) < EPSILON) {
124: u1 = 0.0;
125: u2 = 1.0;
126: v1 = v2 = const/ycoeff;
127: } else if (fabs(ycoeff) < EPSILON) {
128: v1 = 0.0;
129: v2 = 1.0;
130: u1 = u2 = const/xcoeff;
131: } else if (fabs(const) < EPSILON) {
132: u1 = 0.0;
133: v1 = 0.0;
134: u2 = 1.0;
135: v2 = (const - 1.0/xcoeff)/ycoeff;
136: } else {
137: u1 = 0.0;
138: v1 = const/ycoeff;
139: u2 = const/xcoeff;
140: v2 = 0.0;
141: }
142: lncrc = lcinter (u1,v1, u2,v2, x1,y1,r1, theta1,phi1, theta2,phi2);
143: if (lncrc) {
144: *phi1 = rprin (*phi1);
145: *phi2 = rprin (*phi2);
146: *theta1 = atan2 (y1 + r1*sin(*phi1) - y0, x1 + r1*cos(*phi1) - x0);
147: *theta2 = atan2 (y1 + r1*sin(*phi2) - y0, x1 + r1*cos(*phi2) - x0);
148: *theta1 = rprin (*theta1);
149: *theta2 = rprin (*theta2);
150: dprintf "intersection theta1 = %f phi1 = %f\n", *theta1, *phi1);
151: dprintf "intersection theta2 = %f phi2 = %f\n", *theta2, *phi2);
152: return (TRUE);
153: } else
154: return (FALSE);
155: }
This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.