1 #define INKSCAPE_HELPER_GEOM_CPP
3 /**
4 * Specific geometry functions for Inkscape, not provided my lib2geom.
5 *
6 * Author:
7 * Johan Engelen <goejendaagh@zonnet.nl>
8 *
9 * Copyright (C) 2008 Johan Engelen
10 *
11 * Released under GNU GPL
12 */
14 #include "helper/geom.h"
15 #include <typeinfo>
16 #include <2geom/pathvector.h>
17 #include <2geom/path.h>
18 #include <2geom/transforms.h>
19 #include <2geom/rect.h>
20 #include <2geom/coord.h>
21 #include <2geom/sbasis-to-bezier.h>
22 #include <libnr/nr-convert2geom.h>
23 #include <glibmm.h>
25 using Geom::X;
26 using Geom::Y;
28 #define NR_HUGE 1e18
30 //#################################################################################
31 // BOUNDING BOX CALCULATIONS
33 /* Fast bbox calculation */
34 /* Thanks to Nathan Hurst for suggesting it */
35 static void
36 cubic_bbox (Geom::Coord x000, Geom::Coord y000, Geom::Coord x001, Geom::Coord y001, Geom::Coord x011, Geom::Coord y011, Geom::Coord x111, Geom::Coord y111, Geom::Rect &bbox)
37 {
38 Geom::Coord a, b, c, D;
40 bbox[0].extendTo(x111);
41 bbox[1].extendTo(y111);
43 /*
44 * xttt = s * (s * (s * x000 + t * x001) + t * (s * x001 + t * x011)) + t * (s * (s * x001 + t * x011) + t * (s * x011 + t * x111))
45 * xttt = s * (s2 * x000 + s * t * x001 + t * s * x001 + t2 * x011) + t * (s2 * x001 + s * t * x011 + t * s * x011 + t2 * x111)
46 * xttt = s * (s2 * x000 + 2 * st * x001 + t2 * x011) + t * (s2 * x001 + 2 * st * x011 + t2 * x111)
47 * xttt = s3 * x000 + 2 * s2t * x001 + st2 * x011 + s2t * x001 + 2st2 * x011 + t3 * x111
48 * xttt = s3 * x000 + 3s2t * x001 + 3st2 * x011 + t3 * x111
49 * xttt = s3 * x000 + (1 - s) 3s2 * x001 + (1 - s) * (1 - s) * 3s * x011 + (1 - s) * (1 - s) * (1 - s) * x111
50 * xttt = s3 * x000 + (3s2 - 3s3) * x001 + (3s - 6s2 + 3s3) * x011 + (1 - 2s + s2 - s + 2s2 - s3) * x111
51 * xttt = (x000 - 3 * x001 + 3 * x011 - x111) * s3 +
52 * ( 3 * x001 - 6 * x011 + 3 * x111) * s2 +
53 * ( 3 * x011 - 3 * x111) * s +
54 * ( x111)
55 * xttt' = (3 * x000 - 9 * x001 + 9 * x011 - 3 * x111) * s2 +
56 * ( 6 * x001 - 12 * x011 + 6 * x111) * s +
57 * ( 3 * x011 - 3 * x111)
58 */
60 a = 3 * x000 - 9 * x001 + 9 * x011 - 3 * x111;
61 b = 6 * x001 - 12 * x011 + 6 * x111;
62 c = 3 * x011 - 3 * x111;
64 /*
65 * s = (-b +/- sqrt (b * b - 4 * a * c)) / 2 * a;
66 */
67 if (fabs (a) < Geom::EPSILON) {
68 /* s = -c / b */
69 if (fabs (b) > Geom::EPSILON) {
70 double s, t, xttt;
71 s = -c / b;
72 if ((s > 0.0) && (s < 1.0)) {
73 t = 1.0 - s;
74 xttt = s * s * s * x000 + 3 * s * s * t * x001 + 3 * s * t * t * x011 + t * t * t * x111;
75 bbox[0].extendTo(xttt);
76 }
77 }
78 } else {
79 /* s = (-b +/- sqrt (b * b - 4 * a * c)) / 2 * a; */
80 D = b * b - 4 * a * c;
81 if (D >= 0.0) {
82 Geom::Coord d, s, t, xttt;
83 /* Have solution */
84 d = sqrt (D);
85 s = (-b + d) / (2 * a);
86 if ((s > 0.0) && (s < 1.0)) {
87 t = 1.0 - s;
88 xttt = s * s * s * x000 + 3 * s * s * t * x001 + 3 * s * t * t * x011 + t * t * t * x111;
89 bbox[0].extendTo(xttt);
90 }
91 s = (-b - d) / (2 * a);
92 if ((s > 0.0) && (s < 1.0)) {
93 t = 1.0 - s;
94 xttt = s * s * s * x000 + 3 * s * s * t * x001 + 3 * s * t * t * x011 + t * t * t * x111;
95 bbox[0].extendTo(xttt);
96 }
97 }
98 }
100 a = 3 * y000 - 9 * y001 + 9 * y011 - 3 * y111;
101 b = 6 * y001 - 12 * y011 + 6 * y111;
102 c = 3 * y011 - 3 * y111;
104 if (fabs (a) < Geom::EPSILON) {
105 /* s = -c / b */
106 if (fabs (b) > Geom::EPSILON) {
107 double s, t, yttt;
108 s = -c / b;
109 if ((s > 0.0) && (s < 1.0)) {
110 t = 1.0 - s;
111 yttt = s * s * s * y000 + 3 * s * s * t * y001 + 3 * s * t * t * y011 + t * t * t * y111;
112 bbox[1].extendTo(yttt);
113 }
114 }
115 } else {
116 /* s = (-b +/- sqrt (b * b - 4 * a * c)) / 2 * a; */
117 D = b * b - 4 * a * c;
118 if (D >= 0.0) {
119 Geom::Coord d, s, t, yttt;
120 /* Have solution */
121 d = sqrt (D);
122 s = (-b + d) / (2 * a);
123 if ((s > 0.0) && (s < 1.0)) {
124 t = 1.0 - s;
125 yttt = s * s * s * y000 + 3 * s * s * t * y001 + 3 * s * t * t * y011 + t * t * t * y111;
126 bbox[1].extendTo(yttt);
127 }
128 s = (-b - d) / (2 * a);
129 if ((s > 0.0) && (s < 1.0)) {
130 t = 1.0 - s;
131 yttt = s * s * s * y000 + 3 * s * s * t * y001 + 3 * s * t * t * y011 + t * t * t * y111;
132 bbox[1].extendTo(yttt);
133 }
134 }
135 }
136 }
138 Geom::Rect
139 bounds_fast_transformed(Geom::PathVector const & pv, Geom::Matrix const & t)
140 {
141 return bounds_exact_transformed(pv, t); //use this as it is faster for now! :)
142 // return Geom::bounds_fast(pv * t);
143 }
145 Geom::Rect
146 bounds_exact_transformed(Geom::PathVector const & pv, Geom::Matrix const & t)
147 {
148 Geom::Rect bbox;
150 if (pv.empty())
151 return bbox;
153 Geom::Point initial = pv.front().initialPoint() * t;
154 bbox = Geom::Rect(initial, initial); // obtain well defined bbox as starting point to unionWith
156 for (Geom::PathVector::const_iterator it = pv.begin(); it != pv.end(); ++it) {
157 bbox.expandTo(it->initialPoint() * t);
159 // don't loop including closing segment, since that segment can never increase the bbox
160 for (Geom::Path::const_iterator cit = it->begin(); cit != it->end_open(); ++cit) {
161 Geom::Curve const &c = *cit;
163 if( dynamic_cast<Geom::LineSegment const*>(&c) ||
164 dynamic_cast<Geom::HLineSegment const*>(&c) ||
165 dynamic_cast<Geom::VLineSegment const*>(&c) )
166 {
167 bbox.expandTo( c.finalPoint() * t );
168 }
169 else if(Geom::CubicBezier const *cubic_bezier = dynamic_cast<Geom::CubicBezier const *>(&c))
170 {
171 Geom::Point c0 = (*cubic_bezier)[0] * t;
172 Geom::Point c1 = (*cubic_bezier)[1] * t;
173 Geom::Point c2 = (*cubic_bezier)[2] * t;
174 Geom::Point c3 = (*cubic_bezier)[3] * t;
175 cubic_bbox( c0[0], c0[1],
176 c1[0], c1[1],
177 c2[0], c2[1],
178 c3[0], c3[1],
179 bbox );
180 }
181 else
182 {
183 // should handle all not-so-easy curves:
184 Geom::Curve *ctemp = cit->transformed(t);
185 bbox.unionWith( ctemp->boundsExact());
186 delete ctemp;
187 }
188 }
189 }
190 //return Geom::bounds_exact(pv * t);
191 return bbox;
192 }
196 static void
197 geom_line_wind_distance (Geom::Coord x0, Geom::Coord y0, Geom::Coord x1, Geom::Coord y1, Geom::Point const &pt, int *wind, Geom::Coord *best)
198 {
199 Geom::Coord Ax, Ay, Bx, By, Dx, Dy, s;
200 Geom::Coord dist2;
202 /* Find distance */
203 Ax = x0;
204 Ay = y0;
205 Bx = x1;
206 By = y1;
207 Dx = x1 - x0;
208 Dy = y1 - y0;
209 const Geom::Coord Px = pt[X];
210 const Geom::Coord Py = pt[Y];
212 if (best) {
213 s = ((Px - Ax) * Dx + (Py - Ay) * Dy) / (Dx * Dx + Dy * Dy);
214 if (s <= 0.0) {
215 dist2 = (Px - Ax) * (Px - Ax) + (Py - Ay) * (Py - Ay);
216 } else if (s >= 1.0) {
217 dist2 = (Px - Bx) * (Px - Bx) + (Py - By) * (Py - By);
218 } else {
219 Geom::Coord Qx, Qy;
220 Qx = Ax + s * Dx;
221 Qy = Ay + s * Dy;
222 dist2 = (Px - Qx) * (Px - Qx) + (Py - Qy) * (Py - Qy);
223 }
225 if (dist2 < (*best * *best)) *best = sqrt (dist2);
226 }
228 if (wind) {
229 /* Find wind */
230 if ((Ax >= Px) && (Bx >= Px)) return;
231 if ((Ay >= Py) && (By >= Py)) return;
232 if ((Ay < Py) && (By < Py)) return;
233 if (Ay == By) return;
234 /* Ctach upper y bound */
235 if (Ay == Py) {
236 if (Ax < Px) *wind -= 1;
237 return;
238 } else if (By == Py) {
239 if (Bx < Px) *wind += 1;
240 return;
241 } else {
242 Geom::Coord Qx;
243 /* Have to calculate intersection */
244 Qx = Ax + Dx * (Py - Ay) / Dy;
245 if (Qx < Px) {
246 *wind += (Dy > 0.0) ? 1 : -1;
247 }
248 }
249 }
250 }
252 static void
253 geom_cubic_bbox_wind_distance (Geom::Coord x000, Geom::Coord y000,
254 Geom::Coord x001, Geom::Coord y001,
255 Geom::Coord x011, Geom::Coord y011,
256 Geom::Coord x111, Geom::Coord y111,
257 Geom::Point const &pt,
258 Geom::Rect *bbox, int *wind, Geom::Coord *best,
259 Geom::Coord tolerance)
260 {
261 Geom::Coord x0, y0, x1, y1, len2;
262 int needdist, needwind, needline;
264 const Geom::Coord Px = pt[X];
265 const Geom::Coord Py = pt[Y];
267 needdist = 0;
268 needwind = 0;
269 needline = 0;
271 if (bbox) cubic_bbox (x000, y000, x001, y001, x011, y011, x111, y111, *bbox);
273 x0 = MIN (x000, x001);
274 x0 = MIN (x0, x011);
275 x0 = MIN (x0, x111);
276 y0 = MIN (y000, y001);
277 y0 = MIN (y0, y011);
278 y0 = MIN (y0, y111);
279 x1 = MAX (x000, x001);
280 x1 = MAX (x1, x011);
281 x1 = MAX (x1, x111);
282 y1 = MAX (y000, y001);
283 y1 = MAX (y1, y011);
284 y1 = MAX (y1, y111);
286 if (best) {
287 /* Quickly adjust to endpoints */
288 len2 = (x000 - Px) * (x000 - Px) + (y000 - Py) * (y000 - Py);
289 if (len2 < (*best * *best)) *best = (Geom::Coord) sqrt (len2);
290 len2 = (x111 - Px) * (x111 - Px) + (y111 - Py) * (y111 - Py);
291 if (len2 < (*best * *best)) *best = (Geom::Coord) sqrt (len2);
293 if (((x0 - Px) < *best) && ((y0 - Py) < *best) && ((Px - x1) < *best) && ((Py - y1) < *best)) {
294 /* Point is inside sloppy bbox */
295 /* Now we have to decide, whether subdivide */
296 /* fixme: (Lauris) */
297 if (((y1 - y0) > 5.0) || ((x1 - x0) > 5.0)) {
298 needdist = 1;
299 } else {
300 needline = 1;
301 }
302 }
303 }
304 if (!needdist && wind) {
305 if ((y1 >= Py) && (y0 < Py) && (x0 < Px)) {
306 /* Possible intersection at the left */
307 /* Now we have to decide, whether subdivide */
308 /* fixme: (Lauris) */
309 if (((y1 - y0) > 5.0) || ((x1 - x0) > 5.0)) {
310 needwind = 1;
311 } else {
312 needline = 1;
313 }
314 }
315 }
317 if (needdist || needwind) {
318 Geom::Coord x00t, x0tt, xttt, x1tt, x11t, x01t;
319 Geom::Coord y00t, y0tt, yttt, y1tt, y11t, y01t;
320 Geom::Coord s, t;
322 t = 0.5;
323 s = 1 - t;
325 x00t = s * x000 + t * x001;
326 x01t = s * x001 + t * x011;
327 x11t = s * x011 + t * x111;
328 x0tt = s * x00t + t * x01t;
329 x1tt = s * x01t + t * x11t;
330 xttt = s * x0tt + t * x1tt;
332 y00t = s * y000 + t * y001;
333 y01t = s * y001 + t * y011;
334 y11t = s * y011 + t * y111;
335 y0tt = s * y00t + t * y01t;
336 y1tt = s * y01t + t * y11t;
337 yttt = s * y0tt + t * y1tt;
339 geom_cubic_bbox_wind_distance (x000, y000, x00t, y00t, x0tt, y0tt, xttt, yttt, pt, NULL, wind, best, tolerance);
340 geom_cubic_bbox_wind_distance (xttt, yttt, x1tt, y1tt, x11t, y11t, x111, y111, pt, NULL, wind, best, tolerance);
341 } else if (1 || needline) {
342 geom_line_wind_distance (x000, y000, x111, y111, pt, wind, best);
343 }
344 }
346 static void
347 geom_curve_bbox_wind_distance(Geom::Curve const & c, Geom::Matrix const &m,
348 Geom::Point const &pt,
349 Geom::Rect *bbox, int *wind, Geom::Coord *dist,
350 Geom::Coord tolerance, Geom::Rect const *viewbox,
351 Geom::Point &p0) // pass p0 through as it represents the last endpoint added (the finalPoint of last curve)
352 {
353 if( dynamic_cast<Geom::LineSegment const*>(&c) ||
354 dynamic_cast<Geom::HLineSegment const*>(&c) ||
355 dynamic_cast<Geom::VLineSegment const*>(&c) )
356 {
357 Geom::Point pe = c.finalPoint() * m;
358 if (bbox) {
359 bbox->expandTo(pe);
360 }
361 if (dist || wind) {
362 if (wind) { // we need to pick fill, so do what we're told
363 geom_line_wind_distance (p0[X], p0[Y], pe[X], pe[Y], pt, wind, dist);
364 } else { // only stroke is being picked; skip this segment if it's totally outside the viewbox
365 Geom::Rect swept(p0, pe);
366 if (!viewbox || swept.intersects(*viewbox))
367 geom_line_wind_distance (p0[X], p0[Y], pe[X], pe[Y], pt, wind, dist);
368 }
369 }
370 p0 = pe;
371 }
372 else if(Geom::CubicBezier const *cubic_bezier = dynamic_cast<Geom::CubicBezier const *>(&c)) {
373 Geom::Point p1 = (*cubic_bezier)[1] * m;
374 Geom::Point p2 = (*cubic_bezier)[2] * m;
375 Geom::Point p3 = (*cubic_bezier)[3] * m;
377 // get approximate bbox from handles (convex hull property of beziers):
378 Geom::Rect swept(p0, p3);
379 swept.expandTo(p1);
380 swept.expandTo(p2);
382 if (!viewbox || swept.intersects(*viewbox)) { // we see this segment, so do full processing
383 geom_cubic_bbox_wind_distance ( p0[X], p0[Y],
384 p1[X], p1[Y],
385 p2[X], p2[Y],
386 p3[X], p3[Y],
387 pt,
388 bbox, wind, dist, tolerance);
389 } else {
390 if (wind) { // if we need fill, we can just pretend it's a straight line
391 geom_line_wind_distance (p0[X], p0[Y], p3[X], p3[Y], pt, wind, dist);
392 } else { // otherwise, skip it completely
393 }
394 }
395 p0 = p3;
396 } else {
397 //this case handles sbasis as well as all other curve types
398 Geom::Path sbasis_path = Geom::cubicbezierpath_from_sbasis(c.toSBasis(), 0.1);
400 //recurse to convert the new path resulting from the sbasis to svgd
401 for(Geom::Path::iterator iter = sbasis_path.begin(); iter != sbasis_path.end(); ++iter) {
402 geom_curve_bbox_wind_distance(*iter, m, pt, bbox, wind, dist, tolerance, viewbox, p0);
403 }
404 }
405 }
407 /* Calculates...
408 and returns ... in *wind and the distance to ... in *dist.
409 Returns bounding box in *bbox if bbox!=NULL.
410 */
411 void
412 pathv_matrix_point_bbox_wind_distance (Geom::PathVector const & pathv, Geom::Matrix const &m, Geom::Point const &pt,
413 Geom::Rect *bbox, int *wind, Geom::Coord *dist,
414 Geom::Coord tolerance, Geom::Rect const *viewbox)
415 {
416 if (pathv.empty()) {
417 if (wind) *wind = 0;
418 if (dist) *dist = NR_HUGE;
419 return;
420 }
422 // remember last point of last curve
423 Geom::Point p0(0,0);
425 // remembering the start of subpath
426 Geom::Point p_start(0,0);
427 bool start_set = false;
429 for (Geom::PathVector::const_iterator it = pathv.begin(); it != pathv.end(); ++it) {
431 if (start_set) { // this is a new subpath
432 if (wind && (p0 != p_start)) // for correct fill picking, each subpath must be closed
433 geom_line_wind_distance (p0[X], p0[Y], p_start[X], p_start[Y], pt, wind, dist);
434 }
435 p0 = it->initialPoint() * m;
436 p_start = p0;
437 start_set = true;
438 if (bbox) {
439 bbox->expandTo(p0);
440 }
442 // loop including closing segment if path is closed
443 for (Geom::Path::const_iterator cit = it->begin(); cit != it->end_default(); ++cit) {
444 geom_curve_bbox_wind_distance(*cit, m, pt, bbox, wind, dist, tolerance, viewbox, p0);
445 }
446 }
448 if (start_set) {
449 if (wind && (p0 != p_start)) // for correct picking, each subpath must be closed
450 geom_line_wind_distance (p0[X], p0[Y], p_start[X], p_start[Y], pt, wind, dist);
451 }
452 }
454 // temporary wrapper
455 void
456 pathv_matrix_point_bbox_wind_distance (Geom::PathVector const & pathv, NR::Matrix const &m, NR::Point const &pt,
457 NR::Rect *bbox, int *wind, NR::Coord *dist,
458 NR::Coord tolerance, NR::Rect const *viewbox)
459 {
460 Geom::Rect _bbox;
461 if (bbox)
462 _bbox = to_2geom(*bbox);
463 Geom::Coord _dist;
464 if (dist)
465 _dist = *dist;
466 Geom::Rect _viewbox;
467 if (viewbox)
468 _viewbox = to_2geom(*viewbox);
470 pathv_matrix_point_bbox_wind_distance( pathv, to_2geom(m), to_2geom(pt),
471 bbox ? &_bbox : NULL,
472 wind,
473 dist ? &_dist : NULL,
474 tolerance,
475 viewbox ? &_viewbox : NULL );
477 if (bbox)
478 *bbox = from_2geom(_bbox);
479 if (dist)
480 *dist = _dist;
481 }
482 //#################################################################################
484 /*
485 * Converts all segments in all paths to Geom::LineSegment or Geom::HLineSegment or
486 * Geom::VLineSegment or Geom::CubicBezier.
487 */
488 Geom::PathVector
489 pathv_to_linear_and_cubic_beziers( Geom::PathVector const &pathv )
490 {
491 Geom::PathVector output;
493 for (Geom::PathVector::const_iterator pit = pathv.begin(); pit != pathv.end(); ++pit) {
494 output.push_back( Geom::Path() );
495 output.back().start( pit->initialPoint() );
496 output.back().close( pit->closed() );
498 for (Geom::Path::const_iterator cit = pit->begin(); cit != pit->end_open(); ++cit) {
499 if( dynamic_cast<Geom::LineSegment const*>(&*cit) ||
500 dynamic_cast<Geom::HLineSegment const*>(&*cit) ||
501 dynamic_cast<Geom::VLineSegment const*>(&*cit) )
502 {
503 output.back().append(*cit);
504 }
505 else if(Geom::CubicBezier const *cubic_bezier = dynamic_cast<Geom::CubicBezier const*>(&*cit)) {
506 (void)cubic_bezier;
507 output.back().append(*cit);
508 }
509 else {
510 // convert all other curve types to cubicbeziers
511 Geom::Path cubicbezier_path = Geom::cubicbezierpath_from_sbasis(cit->toSBasis(), 0.1);
513 for(Geom::Path::iterator iter = cubicbezier_path.begin(); iter != cubicbezier_path.end(); ++iter) {
514 output.back().append(*iter);
515 }
516 }
517 }
518 }
520 return output;
521 }
526 /*
527 Local Variables:
528 mode:c++
529 c-file-style:"stroustrup"
530 c-file-offsets:((innamespace . 0)(inline-open . 0)(case-label . +))
531 indent-tabs-mode:nil
532 fill-column:99
533 End:
534 */
535 // vim: filetype=cpp:expandtab:shiftwidth=4:tabstop=8:softtabstop=4:encoding=utf-8:textwidth=99 :