1 /* Dia -- an diagram creation/manipulation program
2 * Copyright (C) 1998 Alexander Larsson
3 *
4 * Matrix code from GIMP:app/transform_tool.c
5 * Copyright (C) 1995 Spencer Kimball and Peter Mattis
6 *
7 * This program is free software; you can redistribute it and/or modify
8 * it under the terms of the GNU General Public License as published by
9 * the Free Software Foundation; either version 2 of the License, or
10 * (at your option) any later version.
11 *
12 * This program is distributed in the hope that it will be useful,
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 * GNU General Public License for more details.
16 *
17 * You should have received a copy of the GNU General Public License
18 * along with this program; if not, write to the Free Software
19 * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
20 */
21
22 /* include glib.h first, so we don't pick up its inline functions as well */
23 #include <config.h>
24
25 #include <glib.h>
26 /* include normal versions of the inlined functions here ... */
27 #undef G_INLINE_FUNC
28 #define G_INLINE_FUNC extern
29 #define G_CAN_INLINE 1
30 #include "geometry.h"
31
32 #include "object.h"
33 #include "units.h"
34
35 void
rectangle_union(Rectangle * r1,const Rectangle * r2)36 rectangle_union(Rectangle *r1, const Rectangle *r2)
37 {
38 r1->top = MIN( r1->top, r2->top );
39 r1->bottom = MAX( r1->bottom, r2->bottom );
40 r1->left = MIN( r1->left, r2->left );
41 r1->right = MAX( r1->right, r2->right );
42 }
43
44 void
int_rectangle_union(IntRectangle * r1,const IntRectangle * r2)45 int_rectangle_union(IntRectangle *r1, const IntRectangle *r2)
46 {
47 r1->top = MIN( r1->top, r2->top );
48 r1->bottom = MAX( r1->bottom, r2->bottom );
49 r1->left = MIN( r1->left, r2->left );
50 r1->right = MAX( r1->right, r2->right );
51 }
52
53 void
rectangle_intersection(Rectangle * r1,const Rectangle * r2)54 rectangle_intersection(Rectangle *r1, const Rectangle *r2)
55 {
56 r1->top = MAX( r1->top, r2->top );
57 r1->bottom = MIN( r1->bottom, r2->bottom );
58 r1->left = MAX( r1->left, r2->left );
59 r1->right = MIN( r1->right, r2->right );
60
61 /* Is rectangle empty? */
62 if ( (r1->top >= r1->bottom) ||
63 (r1->left >= r1->right) ) {
64 r1->top = r1->bottom = r1->left = r1->right = 0.0;
65 }
66 }
67
68 int
rectangle_intersects(const Rectangle * r1,const Rectangle * r2)69 rectangle_intersects(const Rectangle *r1, const Rectangle *r2)
70 {
71 if ( (r1->right < r2->left) ||
72 (r1->left > r2->right) ||
73 (r1->top > r2->bottom) ||
74 (r1->bottom < r2->top) )
75 return FALSE;
76
77 return TRUE;
78 }
79
80 int
point_in_rectangle(const Rectangle * r,const Point * p)81 point_in_rectangle(const Rectangle* r, const Point *p)
82 {
83 if ( (p->x < r->left) ||
84 (p->x > r->right) ||
85 (p->y > r->bottom) ||
86 (p->y < r->top) )
87 return FALSE;
88
89 return TRUE;
90 }
91
92 int
rectangle_in_rectangle(const Rectangle * outer,const Rectangle * inner)93 rectangle_in_rectangle(const Rectangle* outer, const Rectangle *inner)
94 {
95 if ( (inner->left < outer->left) ||
96 (inner->right > outer->right) ||
97 (inner->top < outer->top) ||
98 (inner->bottom > outer->bottom))
99 return FALSE;
100
101 return TRUE;
102 }
103
104 void
rectangle_add_point(Rectangle * r,const Point * p)105 rectangle_add_point(Rectangle *r, const Point *p)
106 {
107 if (p->x < r->left)
108 r->left = p->x;
109 else if (p->x > r->right)
110 r->right = p->x;
111
112 if (p->y < r->top)
113 r->top = p->y;
114 else if (p->y > r->bottom)
115 r->bottom = p->y;
116 }
117
118
119 /*
120 * This function estimates the distance from a point to a rectangle.
121 * If the point is in the rectangle, 0.0 is returned. Otherwise the
122 * distance in a manhattan metric from the point to the nearest point
123 * on the rectangle is returned.
124 */
125 real
distance_rectangle_point(const Rectangle * rect,const Point * point)126 distance_rectangle_point(const Rectangle *rect, const Point *point)
127 {
128 real dx = 0.0;
129 real dy = 0.0;
130
131 if (point->x < rect->left ) {
132 dx = rect->left - point->x;
133 } else if (point->x > rect->right ) {
134 dx = point->x - rect->right;
135 }
136
137 if (point->y < rect->top ) {
138 dy = rect->top - point->y;
139 } else if (point->y > rect->bottom ) {
140 dy = point->y - rect->bottom;
141 }
142
143 return dx + dy;
144 }
145
146 /*
147 * This function estimates the distance from a point to a line segment
148 * specified by two endpoints.
149 * If the point is on the line segment, 0.0 is returned. Otherwise the
150 * distance in the R^2 metric from the point to the nearest point
151 * on the line segment is returned. Does one sqrt per call.
152 * Philosophical bug: line_width is ignored iff point is beyond
153 * end of line segment.
154 */
155 real
distance_line_point(const Point * line_start,const Point * line_end,real line_width,const Point * point)156 distance_line_point(const Point *line_start, const Point *line_end,
157 real line_width, const Point *point)
158 {
159 Point v1, v2;
160 real v1_lensq;
161 real projlen;
162 real perp_dist;
163
164 v1 = *line_end;
165 point_sub(&v1,line_start);
166
167 v2 = *point;
168 point_sub(&v2, line_start);
169
170 v1_lensq = point_dot(&v1,&v1);
171
172 if (v1_lensq<0.000001) {
173 return sqrt(point_dot(&v2,&v2));
174 }
175
176 projlen = point_dot(&v1,&v2) / v1_lensq;
177
178 if (projlen<0.0) {
179 return sqrt(point_dot(&v2,&v2));
180 }
181
182 if (projlen>1.0) {
183 Point v3 = *point;
184 point_sub(&v3,line_end);
185 return sqrt(point_dot(&v3,&v3));
186 }
187
188 point_scale(&v1, projlen);
189 point_sub(&v1, &v2);
190
191 perp_dist = sqrt(point_dot(&v1,&v1));
192
193 perp_dist -= line_width / 2.0;
194 if (perp_dist < 0.0) {
195 perp_dist = 0.0;
196 }
197
198 return perp_dist;
199 }
200
201 /* returns 1 if the line crosses the ray from (-Inf, rayend.y) to rayend */
202 static int
line_crosses_ray(const Point * line_start,const Point * line_end,const Point * rayend)203 line_crosses_ray(const Point *line_start,
204 const Point *line_end, const Point *rayend)
205 {
206 if ((line_start->y <= rayend->y && line_end->y > rayend->y) || /* upward crossing */
207 (line_start->y > rayend->y && line_end->y <= rayend->y)) { /* downward crossing */
208 real vt = (rayend->y - line_start->y) / (line_end->y - line_start->y);
209 if (rayend->x < line_start->x + vt * (line_end->x - line_start->x)) /* intersect */
210 return 1;
211 }
212 return 0;
213 }
214
215 real
distance_polygon_point(const Point * poly,guint npoints,real line_width,const Point * point)216 distance_polygon_point(const Point *poly, guint npoints, real line_width,
217 const Point *point)
218 {
219 guint i, last = npoints - 1;
220 real line_dist = G_MAXFLOAT;
221 guint crossings = 0;
222
223 /* calculate ray crossings and line distances */
224 for (i = 0; i < npoints; i++) {
225 real dist;
226
227 crossings += line_crosses_ray(&poly[last], &poly[i], point);
228 dist = distance_line_point(&poly[last], &poly[i], line_width, point);
229 line_dist = MIN(line_dist, dist);
230 last = i;
231 }
232 /* If there is an odd number of ray crossings, we are inside the polygon.
233 * Otherwise, return the minium distance from a line segment */
234 if (crossings % 2 == 1)
235 return 0.0;
236 else
237 return line_dist;
238 }
239
240 /* number of segments to use in bezier curve approximation */
241 #define NBEZ_SEGS 10
242
243 /* if cross is not NULL, it will be incremented for each ray crossing */
244 static real
bez_point_distance_and_ray_crosses(const Point * b1,const Point * b2,const Point * b3,const Point * b4,real line_width,const Point * point,guint * cross)245 bez_point_distance_and_ray_crosses(const Point *b1,
246 const Point *b2, const Point *b3,
247 const Point *b4,
248 real line_width, const Point *point,
249 guint *cross)
250 {
251 static gboolean calculated_coeff = FALSE;
252 static real coeff[NBEZ_SEGS+1][4];
253 int i;
254 real line_dist = G_MAXFLOAT;
255 Point prev, pt;
256
257 if (!calculated_coeff) {
258 for (i = 0; i <= NBEZ_SEGS; i++) {
259 real t1 = ((real)i)/NBEZ_SEGS, t2 = t1*t1, t3 = t1*t2;
260 real it1 = 1-t1, it2 = it1*it1, it3 = it1*it2;
261
262 coeff[i][0] = it3;
263 coeff[i][1] = 3 * t1 * it2;
264 coeff[i][2] = 3 * t2 * it1;
265 coeff[i][3] = t3;
266 }
267 }
268 calculated_coeff = TRUE;
269
270 prev.x = coeff[0][0] * b1->x + coeff[0][1] * b2->x +
271 coeff[0][2] * b3->x + coeff[0][3] * b4->x;
272 prev.y = coeff[0][0] * b1->y + coeff[0][1] * b2->y +
273 coeff[0][2] * b3->y + coeff[0][3] * b4->y;
274 for (i = 1; i <= NBEZ_SEGS; i++) {
275 real dist;
276
277 pt.x = coeff[i][0] * b1->x + coeff[i][1] * b2->x +
278 coeff[i][2] * b3->x + coeff[i][3] * b4->x;
279 pt.y = coeff[i][0] * b1->y + coeff[i][1] * b2->y +
280 coeff[i][2] * b3->y + coeff[i][3] * b4->y;
281
282 dist = distance_line_point(&prev, &pt, line_width, point);
283 line_dist = MIN(line_dist, dist);
284 if (cross)
285 *cross += line_crosses_ray(&prev, &pt, point);
286
287 prev = pt;
288 }
289 return line_dist;
290 }
291
292 real
distance_bez_seg_point(const Point * b1,const Point * b2,const Point * b3,const Point * b4,real line_width,const Point * point)293 distance_bez_seg_point(const Point *b1, const Point *b2,
294 const Point *b3, const Point *b4,
295 real line_width, const Point *point)
296 {
297 return bez_point_distance_and_ray_crosses(b1, b2, b3, b4,
298 line_width, point, NULL);
299 }
300
301 real
distance_bez_line_point(const BezPoint * b,guint npoints,real line_width,const Point * point)302 distance_bez_line_point(const BezPoint *b, guint npoints,
303 real line_width, const Point *point)
304 {
305 Point last;
306 guint i;
307 real line_dist = G_MAXFLOAT;
308
309 g_return_val_if_fail(b[0].type == BEZ_MOVE_TO, -1);
310
311 last = b[0].p1;
312
313 for (i = 1; i < npoints; i++) {
314 real dist;
315
316 switch (b[i].type) {
317 case BEZ_MOVE_TO:
318 /*g_assert_not_reached();*/
319 g_warning("BEZ_MOVE_TO found half way through a bezier line");
320 break;
321 case BEZ_LINE_TO:
322 dist = distance_line_point(&last, &b[i].p1, line_width, point);
323 line_dist = MIN(line_dist, dist);
324 last = b[i].p1;
325 break;
326 case BEZ_CURVE_TO:
327 dist = bez_point_distance_and_ray_crosses(&last, &b[i].p1, &b[i].p2,
328 &b[i].p3, line_width, point,
329 NULL);
330 line_dist = MIN(line_dist, dist);
331 last = b[i].p3;
332 break;
333 }
334 }
335 return line_dist;
336 }
337
338 real
distance_bez_shape_point(const BezPoint * b,guint npoints,real line_width,const Point * point)339 distance_bez_shape_point(const BezPoint *b, guint npoints,
340 real line_width, const Point *point)
341 {
342 Point last;
343 guint i;
344 real line_dist = G_MAXFLOAT;
345 guint crossings = 0;
346
347 g_return_val_if_fail(b[0].type == BEZ_MOVE_TO, -1);
348
349 last = b[0].p1;
350
351 for (i = 1; i < npoints; i++) {
352 real dist;
353
354 switch (b[i].type) {
355 case BEZ_MOVE_TO:
356 /*g_assert_not_reached();*/
357 g_warning("BEZ_MOVE_TO found half way through a bezier shape");
358 break;
359 case BEZ_LINE_TO:
360 dist = distance_line_point(&last, &b[i].p1, line_width, point);
361 crossings += line_crosses_ray(&last, &b[i].p1, point);
362 line_dist = MIN(line_dist, dist);
363 last = b[i].p1;
364 break;
365 case BEZ_CURVE_TO:
366 dist = bez_point_distance_and_ray_crosses(&last, &b[i].p1, &b[i].p2,
367 &b[i].p3, line_width, point,
368 &crossings);
369 line_dist = MIN(line_dist, dist);
370 last = b[i].p3;
371 break;
372 }
373 }
374 /* If there is an odd number of ray crossings, we are inside the polygon.
375 * Otherwise, return the minium distance from a line segment */
376 if (crossings % 2 == 1)
377 return 0.0;
378 else
379 return line_dist;
380 }
381
382 real
distance_ellipse_point(const Point * centre,real width,real height,real line_width,const Point * point)383 distance_ellipse_point(const Point *centre, real width, real height,
384 real line_width, const Point *point)
385 {
386 /* A faster intersection method would be to scaled the ellipse and the
387 * point to where the ellipse is a circle, intersect with the circle,
388 * then scale back.
389 */
390 real w2 = width*width, h2 = height*height;
391 real scale, rad, dist;
392 Point pt;
393
394 /* find the radius of the ellipse in the appropriate direction */
395
396 /* find the point of intersection between line (x=cx+(px-cx)t; y=cy+(py-cy)t)
397 * and ellipse ((x-cx)^2)/(w/2)^2 + ((y-cy)^2)/(h/2)^2 = 1 */
398 /* radius along ray is sqrt((px-cx)^2 * t^2 + (py-cy)^2 * t^2) */
399
400 pt = *point;
401 point_sub(&pt, centre);
402 pt.x *= pt.x;
403 pt.y *= pt.y; /* pt = (point - centre).^2 */
404
405 scale = w2 * h2 / (4*h2*pt.x + 4*w2*pt.y);
406 rad = sqrt((pt.x + pt.y)*scale) + line_width/2;
407
408 dist = sqrt(pt.x + pt.y);
409
410 if (dist <= rad)
411 return 0.0;
412 else
413 return dist - rad;
414 }
415
416
417 void
transform_point(Matrix m,Point * src,Point * dest)418 transform_point (Matrix m, Point *src, Point *dest)
419 {
420 real xx, yy, ww;
421
422 xx = m[0][0] * src->x + m[0][1] * src->y + m[0][2];
423 yy = m[1][0] * src->x + m[1][1] * src->y + m[1][2];
424 ww = m[2][0] * src->x + m[2][1] * src->y + m[2][2];
425
426 if (!ww)
427 ww = 1.0;
428
429 dest->x = xx / ww;
430 dest->y = yy / ww;
431 }
432
433 void
mult_matrix(Matrix m1,Matrix m2)434 mult_matrix (Matrix m1, Matrix m2)
435 {
436 Matrix result;
437 int i, j, k;
438
439 for (i = 0; i < 3; i++)
440 for (j = 0; j < 3; j++)
441 {
442 result [i][j] = 0.0;
443 for (k = 0; k < 3; k++)
444 result [i][j] += m1 [i][k] * m2[k][j];
445 }
446
447 /* copy the result into matrix 2 */
448 for (i = 0; i < 3; i++)
449 for (j = 0; j < 3; j++)
450 m2 [i][j] = result [i][j];
451 }
452
453 void
identity_matrix(Matrix m)454 identity_matrix (Matrix m)
455 {
456 int i, j;
457
458 for (i = 0; i < 3; i++)
459 for (j = 0; j < 3; j++)
460 m[i][j] = (i == j) ? 1 : 0;
461
462 }
463
464 void
translate_matrix(Matrix m,real x,real y)465 translate_matrix (Matrix m, real x, real y)
466 {
467 Matrix trans;
468
469 identity_matrix (trans);
470 trans[0][2] = x;
471 trans[1][2] = y;
472 mult_matrix (trans, m);
473 }
474
475 void
scale_matrix(Matrix m,real x,real y)476 scale_matrix (Matrix m, real x, real y)
477 {
478 Matrix scale;
479
480 identity_matrix (scale);
481 scale[0][0] = x;
482 scale[1][1] = y;
483 mult_matrix (scale, m);
484 }
485
486 void
rotate_matrix(Matrix m,real theta)487 rotate_matrix (Matrix m, real theta)
488 {
489 Matrix rotate;
490 real cos_theta, sin_theta;
491
492 cos_theta = cos (theta);
493 sin_theta = sin (theta);
494
495 identity_matrix (rotate);
496 rotate[0][0] = cos_theta;
497 rotate[0][1] = -sin_theta;
498 rotate[1][0] = sin_theta;
499 rotate[1][1] = cos_theta;
500 mult_matrix (rotate, m);
501 }
502
503 void
xshear_matrix(Matrix m,real shear)504 xshear_matrix (Matrix m, real shear)
505 {
506 Matrix shear_m;
507
508 identity_matrix (shear_m);
509 shear_m[0][1] = shear;
510 mult_matrix (shear_m, m);
511 }
512
513 void
yshear_matrix(Matrix m,real shear)514 yshear_matrix (Matrix m, real shear)
515 {
516 Matrix shear_m;
517
518 identity_matrix (shear_m);
519 shear_m[1][0] = shear;
520 mult_matrix (shear_m, m);
521 }
522
523 /* find the determinate for a 3x3 matrix */
524 static real
determinate(Matrix m)525 determinate (Matrix m)
526 {
527 int i;
528 double det = 0;
529
530 for (i = 0; i < 3; i ++)
531 {
532 det += m[0][i] * m[1][(i+1)%3] * m[2][(i+2)%3];
533 det -= m[2][i] * m[1][(i+1)%3] * m[0][(i+2)%3];
534 }
535
536 return det;
537 }
538
539
540
541 void
point_convex(Point * dst,const Point * src1,const Point * src2,real alpha)542 point_convex(Point *dst, const Point *src1, const Point *src2, real alpha)
543 {
544 /* Make convex combination of src1 and src2:
545 dst = alpha * src1 + (1-alpha) * src2;
546 */
547 point_copy(dst,src1);
548 point_scale(dst,alpha);
549 point_add_scaled(dst,src2,1.0 - alpha);
550 }
551
552
553
554 /* The following functions were originally taken from Graphics Gems III,
555 pg 193 and corresponding code pgs 496-499
556
557 They were modified to work within the dia coordinate system
558 */
559 /* return angle subtended by two vectors */
dot2(Point * p1,Point * p2)560 real dot2(Point *p1, Point *p2)
561 {
562 real d, t;
563 d = sqrt(((p1->x*p1->x)+(p1->y*p1->y)) *
564 ((p2->x*p2->x)+(p2->y*p2->y)));
565 if ( d != 0.0 )
566 {
567 t = (p1->x*p2->x+p1->y*p2->y)/d;
568 return (dia_acos(t));
569 }
570 else
571 {
572 return 0.0;
573 }
574 }
575
576 /* Find a,b,c in ax + by + c = 0 for line p1, p2 */
line_coef(real * a,real * b,real * c,Point * p1,Point * p2)577 void line_coef(real *a, real *b, real *c, Point *p1, Point *p2)
578 {
579 *c = ( p2->x * p1->y ) - ( p1->x * p2->y );
580 *a = p2->y - p1->y;
581 *b = p1->x - p2->x;
582 }
583
584 /* Return signed distance from line
585 ax + by + c = 0 to point p. */
line_to_point(real a,real b,real c,Point * p)586 real line_to_point(real a, real b , real c, Point *p) {
587 real d, lp;
588 d = sqrt((a*a)+(b*b));
589 if ( d == 0.0 ) {
590 lp = 0.0;
591 } else {
592 lp = (a*p->x+b*p->y+c)/d;
593 }
594 return (lp);
595 }
596
597 /* Given line l = ax + by + c = 0 and point p,
598 compute x,y so point perp is perpendicular to line l. */
point_perp(Point * p,real a,real b,real c,Point * perp)599 void point_perp(Point *p, real a, real b, real c, Point *perp) {
600 real d, cp;
601 perp->x = 0.0;
602 perp->y = 0.0;
603 d = a*a + b*b;
604 cp = a*p->y-b*p->x;
605 if ( d != 0.0 ) {
606 perp->x = (-a*c-b*cp)/d;
607 perp->y = (a*cp-b*c)/d;
608 }
609 return;
610 }
611
612 /* Compute a circular arc fillet between lines L1 (p1 to p2)
613 and L2 (p3 to p4) with radius r.
614 The circle center is c.
615 The start angle is pa
616 The end angle is aa,
617 The points p1-p4 will be modified as necessary */
fillet(Point * p1,Point * p2,Point * p3,Point * p4,real r,Point * c,real * pa,real * aa)618 void fillet(Point *p1, Point *p2, Point *p3, Point *p4,
619 real r, Point *c, real *pa, real *aa) {
620 real a1, b1, c1; /* Coefficients for L1 */
621 real a2, b2, c2; /* Coefficients for L2 */
622 real d, d1, d2;
623 real c1p, c2p;
624 real rr;
625 real start_angle, stop_angle;
626 Point mp,gv1,gv2;
627 gboolean righthand = FALSE;
628
629 line_coef(&a1,&b1,&c1,p1,p2);
630 line_coef(&a2,&b2,&c2,p3,p4);
631
632 if ( (a1*b2) == (a2*b1) ) /* Parallel or coincident lines */
633 {
634 return;
635 }
636
637 mp.x = (p3->x + p4->x) / 2.0; /* Find midpoint of p3 p4 */
638 mp.y = (p3->y + p4->y) / 2.0;
639 d1 = line_to_point(a1, b1, c1, &mp); /* Find distance p1 p2 to
640 midpoint p3 p4 */
641 if ( d1 == 0.0 ) return; /* p1p2 to p3 */
642
643 mp.x = (p1->x + p2->x) / 2.0; /* Find midpoint of p1 p2 */
644 mp.y = (p1->y + p2->y) / 2.0;
645 d2 = line_to_point(a2, b2, c2, &mp); /* Find distance p3 p4 to
646 midpoint p1 p2 */
647 if ( d2 == 0.0 ) return;
648
649 rr = r;
650 if ( d1 <= 0.0 ) rr = -rr;
651
652 c1p = c1-rr*sqrt((a1*a1)+(b1*b1)); /* Line parallell l1 at d */
653 rr = r;
654
655 if ( d2 <= 0.0 ) rr = -rr;
656 c2p = c2-rr*sqrt((a2*a2)+(b2*b2)); /* Line parallell l2 at d */
657 d = a1*b2-a2*b1;
658
659 c->x = ( c2p*b1-c1p*b2 ) / d; /* Intersect constructed lines */
660 c->y = ( c1p*a2-c2p*a1 ) / d; /* to find center of arc */
661
662 point_perp(c,a1,b1,c1,p2); /* Clip or extend lines as required */
663 point_perp(c,a2,b2,c2,p3);
664
665 /* need to negate the y coords to calculate angles correctly */
666 gv1.x = p2->x-c->x; gv1.y = -(p2->y-c->y);
667 gv2.x = p3->x-c->x; gv2.y = -(p3->y-c->y);
668
669 start_angle = atan2(gv1.y,gv1.x); /* Beginning angle for arc */
670 stop_angle = dot2(&gv1,&gv2);
671 if ( point_cross(&gv1,&gv2) < 0.0 ) {
672 stop_angle = -stop_angle; /* Angle subtended by arc */
673 /* also this means we need to swap our angles later on */
674 righthand = TRUE;
675 }
676 /* now calculate the actual angles in a form that the draw arc function
677 of the renderer can use */
678 start_angle = start_angle*180.0/G_PI;
679 stop_angle = start_angle + stop_angle*180.0/G_PI;
680 while (start_angle < 0.0) start_angle += 360.0;
681 while (stop_angle < 0.0) stop_angle += 360.0;
682 /* swap the start and stop if we had to negate the cross product */
683 if ( righthand ) {
684 real tmp = start_angle;
685 start_angle = stop_angle;
686 stop_angle = tmp;
687 }
688 *pa = start_angle;
689 *aa = stop_angle;
690 }
691
692 int
three_point_circle(const Point * p1,const Point * p2,const Point * p3,Point * center,real * radius)693 three_point_circle (const Point *p1, const Point *p2, const Point *p3,
694 Point* center, real* radius)
695 {
696 const real epsilon = 1e-4;
697 real x1 = p1->x;
698 real y1 = p1->y;
699 real x2 = p2->x;
700 real y2 = p2->y;
701 real x3 = p3->x;
702 real y3 = p3->y;
703 real ma, mb;
704 if (fabs(x2 - x1) < epsilon)
705 return 0;
706 if (fabs(x3 - x2) < epsilon)
707 return 0;
708 ma = (y2 - y1) / (x2 - x1);
709 mb = (y3 - y2) / (x3 - x2);
710 if (fabs (mb - ma) < epsilon)
711 return 0;
712 center->x = (ma*mb*(y1-y3)+mb*(x1+x2)-ma*(x2+x3))/(2*(mb-ma));
713 if (fabs(ma)>epsilon)
714 center->y = -1/ma*(center->x - (x1+x2)/2.0) + (y1+y2)/2.0;
715 else if (fabs(mb)>epsilon)
716 center->y = -1/mb*(center->x - (x2+x3)/2.0) + (y2+y3)/2.0;
717 else
718 return 0;
719 *radius = distance_point_point(center, p1);
720 return 1;
721 }
722
723
724 /* moved this here since it is being reused by rounded polyline code*/
725 real
point_cross(Point * p1,Point * p2)726 point_cross(Point *p1, Point *p2)
727 {
728 return p1->x * p2->y - p2->x * p1->y;
729 }
730
731 /* moved this here since it is used by line and bezier */
732 /* Computes one point between end and objmid which
733 * touches the edge of the object */
734 Point
calculate_object_edge(Point * objmid,Point * end,DiaObject * obj)735 calculate_object_edge(Point *objmid, Point *end, DiaObject *obj)
736 {
737 #define MAXITER 25
738 #ifdef TRACE_DIST
739 Point trace[MAXITER];
740 real disttrace[MAXITER];
741 #endif
742 Point mid1, mid2, mid3;
743 real dist;
744 int i = 0;
745
746 mid1 = *objmid;
747 mid2.x = (objmid->x+end->x)/2;
748 mid2.y = (objmid->y+end->y)/2;
749 mid3 = *end;
750
751 /* If the other end is inside the object */
752 dist = obj->ops->distance_from(obj, &mid3);
753 if (dist < 0.001) return mid1;
754
755
756 do {
757 dist = obj->ops->distance_from(obj, &mid2);
758 if (dist < 0.0000001) {
759 mid1 = mid2;
760 } else {
761 mid3 = mid2;
762 }
763 mid2.x = (mid1.x + mid3.x)/2;
764 mid2.y = (mid1.y + mid3.y)/2;
765 #ifdef TRACE_DIST
766 trace[i] = mid2;
767 disttrace[i] = dist;
768 #endif
769 i++;
770 } while (i < MAXITER && (dist < 0.0000001 || dist > 0.001));
771
772 #ifdef TRACE_DIST
773 if (i == MAXITER) {
774 for (i = 0; i < MAXITER; i++) {
775 printf("%d: %f, %f: %f\n", i, trace[i].x, trace[i].y, disttrace[i]);
776 }
777 printf("i = %d, dist = %f\n", i, dist);
778 }
779 #endif
780
781 return mid2;
782 }
783
784 /*!
785 * \brief asin wrapped to limit to valid result range
786 *
787 * Although clipping the value might hide some miscalculation
788 * elsewhere this should still be used in all of Dia. Continuing
789 * calculation with bogus values might final fall on our foots
790 * because rendering libraries might not be graceful.
791 * See https://bugzilla.gnome.org/show_bug.cgi?id=710818
792 */
793 real
dia_asin(real x)794 dia_asin (real x)
795 {
796 /* clamp to valid range */
797 if (x <= -1.0)
798 return -M_PI/2;
799 if (x >= 1.0)
800 return M_PI/2;
801 return asin (x);
802 }
803
804 /*!
805 * \brief acos wrapped to limit to valid result range
806 */
807 real
dia_acos(real x)808 dia_acos (real x)
809 {
810 /* clamp to valid range */
811 if (x <= -1.0)
812 return M_PI;
813 if (x >= 1.0)
814 return 0.0;
815 return acos (x);
816 }
817