From 60aeb869a2fe7980724b933484cb871c72078ca6 Mon Sep 17 00:00:00 2001 From: "Johan B. C. Engelen" Date: Sat, 11 Oct 2014 19:05:39 +0200 Subject: Fix "Argument with 'nonnull' attribute passed null" API bug in extension color parameter construction. Clang static analyzer has found an "Argument with 'nonnull' attribute passed null". The applied fix here checks if the defaulthex argument is nullptr before making the call. A default initialization value is added to the constructor in case defaulthex == nullptr. (with additional whitespace improvements) (bzr r13587) --- src/extension/param/color.cpp | 12 +++++++----- 1 file changed, 7 insertions(+), 5 deletions(-) (limited to 'src') diff --git a/src/extension/param/color.cpp b/src/extension/param/color.cpp index 0a2598c56..5bd70359f 100644 --- a/src/extension/param/color.cpp +++ b/src/extension/param/color.cpp @@ -59,9 +59,10 @@ guint32 ParamColor::set( guint32 in, SPDocument * /*doc*/, Inkscape::XML::Node * return _value; } -ParamColor::ParamColor (const gchar * name, const gchar * guitext, const gchar * desc, const Parameter::_scope_t scope, bool gui_hidden, const gchar * gui_tip, Inkscape::Extension::Extension * ext, Inkscape::XML::Node * xml) : - Parameter(name, guitext, desc, scope, gui_hidden, gui_tip, ext), - _changeSignal(0) +ParamColor::ParamColor(const gchar *name, const gchar *guitext, const gchar *desc, const Parameter::_scope_t scope, + bool gui_hidden, const gchar *gui_tip, Inkscape::Extension::Extension *ext, + Inkscape::XML::Node *xml) + : Parameter(name, guitext, desc, scope, gui_hidden, gui_tip, ext), _value(0), _changeSignal(0) { const char * defaulthex = NULL; if (xml->firstChild() != NULL) @@ -75,7 +76,8 @@ ParamColor::ParamColor (const gchar * name, const gchar * guitext, const gchar * if (!paramval.empty()) defaulthex = paramval.data(); - _value = atoi(defaulthex); + if (defaulthex) + _value = atoi(defaulthex); } void ParamColor::string(std::string &string) const @@ -87,7 +89,7 @@ void ParamColor::string(std::string &string) const Gtk::Widget *ParamColor::get_widget( SPDocument * /*doc*/, Inkscape::XML::Node * /*node*/, sigc::signal * changeSignal ) { - if (_gui_hidden) return NULL; + if (_gui_hidden) return NULL; _changeSignal = new sigc::signal(*changeSignal); Gtk::HBox * hbox = Gtk::manage(new Gtk::HBox(false, 4)); -- cgit v1.2.3 From a38add027810b9725f81cebb099282a8656e3aa3 Mon Sep 17 00:00:00 2001 From: "Johan B. C. Engelen" Date: Sat, 11 Oct 2014 19:15:03 +0200 Subject: Fix uninitialized variable use in svg-color.cpp. (additionally reduce scope of variable) Bug found using clang static analyzer. The bug happens with color_out not being initialized upon definition; then 1. if(icc) == true 2. if(prof) == true 3. if(trans) == false 4. *r = color_out[0]; <-- un-init use Fixed by initializing color_out (with black). (bzr r13588) --- src/svg/svg-color.cpp | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) (limited to 'src') diff --git a/src/svg/svg-color.cpp b/src/svg/svg-color.cpp index 5108b7702..c9f22f8a4 100644 --- a/src/svg/svg-color.cpp +++ b/src/svg/svg-color.cpp @@ -511,18 +511,18 @@ sp_svg_create_color_hash() void icc_color_to_sRGB(SVGICCColor* icc, guchar* r, guchar* g, guchar* b) { - guchar color_out[4]; - guchar color_in[4]; if (icc) { g_message("profile name: %s", icc->colorProfile.c_str()); Inkscape::ColorProfile* prof = SP_ACTIVE_DOCUMENT->profileManager->find(icc->colorProfile.c_str()); if ( prof ) { + guchar color_out[4] = {0,0,0,0}; cmsHTRANSFORM trans = prof->getTransfToSRGB8(); if ( trans ) { std::vector comps = colorspace::getColorSpaceInfo( prof ); size_t count = CMSSystem::getChannelCount( prof ); size_t cap = std::min(count, comps.size()); + guchar color_in[4]; for (size_t i = 0; i < cap; i++) { color_in[i] = static_cast((((gdouble)icc->colors[i]) * 256.0) * (gdouble)comps[i].scale); g_message("input[%d]: %d", (int)i, (int)color_in[i]); -- cgit v1.2.3 From ccb75de1cadbe4c2e8ad833cbbcd1bc8f79cd548 Mon Sep 17 00:00:00 2001 From: "Johan B. C. Engelen" Date: Sat, 11 Oct 2014 19:21:54 +0200 Subject: LPEVonKoch : remove redundant assignment (is taken care of by the following if-statement) Found by clang static analyzer (bzr r13589) --- src/live_effects/lpe-vonkoch.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) (limited to 'src') diff --git a/src/live_effects/lpe-vonkoch.cpp b/src/live_effects/lpe-vonkoch.cpp index 8b0b716fe..7e396e84a 100644 --- a/src/live_effects/lpe-vonkoch.cpp +++ b/src/live_effects/lpe-vonkoch.cpp @@ -133,7 +133,7 @@ LPEVonKoch::doEffect_path (std::vector const & path_in) for (unsigned k = 0; k < path_in.size(); k++){ path_in_complexity+=path_in[k].size(); } - double complexity = std::pow(transforms.size(),nbgenerations)*path_in_complexity; + double complexity; if (drawall.get_value()){ int k = transforms.size(); if(k>1){ -- cgit v1.2.3 From 82565dc1a86e232f1a809624b14daa2c831b5282 Mon Sep 17 00:00:00 2001 From: Alvin Penner Date: Sat, 11 Oct 2014 17:55:07 -0400 Subject: maintain page size units when resizing page to drawing (Bug 1310787) Fixed bugs: - https://launchpad.net/bugs/1310787 (bzr r13590) --- src/document.cpp | 33 ++++++++++++++++----------------- 1 file changed, 16 insertions(+), 17 deletions(-) (limited to 'src') diff --git a/src/document.cpp b/src/document.cpp index f94a9f04f..f7edf5ed3 100644 --- a/src/document.cpp +++ b/src/document.cpp @@ -706,7 +706,12 @@ void SPDocument::fitToRect(Geom::Rect const &rect, bool with_margins) double const h = rect.height(); double const old_height = getHeight().value("px"); - Inkscape::Util::Unit const *px = unit_table.getUnit("px"); + Inkscape::Util::Unit const *nv_units = unit_table.getUnit("px"); + SPNamedView *nv = sp_document_namedview(this, NULL); + if (nv != NULL) { + if (nv->getAttribute("units")) + nv_units = unit_table.getUnit(nv->getAttribute("units")); + } /* in px */ double margin_top = 0.0; @@ -714,22 +719,16 @@ void SPDocument::fitToRect(Geom::Rect const &rect, bool with_margins) double margin_right = 0.0; double margin_bottom = 0.0; - SPNamedView *nv = sp_document_namedview(this, NULL); - if (with_margins && nv) { if (nv != NULL) { - gchar const * const units_abbr = nv->getAttribute("units"); - Inkscape::Util::Unit const *margin_units = NULL; - if (units_abbr) { - margin_units = unit_table.getUnit(units_abbr); - } - if (!margin_units) { - margin_units = px; - } - margin_top = nv->getMarginLength("fit-margin-top",margin_units, px, w, h, false); - margin_left = nv->getMarginLength("fit-margin-left",margin_units, px, w, h, true); - margin_right = nv->getMarginLength("fit-margin-right",margin_units, px, w, h, true); - margin_bottom = nv->getMarginLength("fit-margin-bottom",margin_units, px, w, h, false); + margin_top = nv->getMarginLength("fit-margin-top", nv_units, unit_table.getUnit("px"), w, h, false); + margin_left = nv->getMarginLength("fit-margin-left", nv_units, unit_table.getUnit("px"), w, h, true); + margin_right = nv->getMarginLength("fit-margin-right", nv_units, unit_table.getUnit("px"), w, h, true); + margin_bottom = nv->getMarginLength("fit-margin-bottom", nv_units, unit_table.getUnit("px"), w, h, false); + margin_top = Inkscape::Util::Quantity::convert(margin_top, nv_units, "px"); + margin_left = Inkscape::Util::Quantity::convert(margin_left, nv_units, "px"); + margin_right = Inkscape::Util::Quantity::convert(margin_right, nv_units, "px"); + margin_bottom = Inkscape::Util::Quantity::convert(margin_bottom, nv_units, "px"); } } @@ -738,8 +737,8 @@ void SPDocument::fitToRect(Geom::Rect const &rect, bool with_margins) rect.max() + Geom::Point(margin_right, margin_top)); - setWidth(Inkscape::Util::Quantity(rect_with_margins.width(), px)); - setHeight(Inkscape::Util::Quantity(rect_with_margins.height(), px)); + setWidth(Inkscape::Util::Quantity(Inkscape::Util::Quantity::convert(rect_with_margins.width(), "px", nv_units), nv_units)); + setHeight(Inkscape::Util::Quantity(Inkscape::Util::Quantity::convert(rect_with_margins.height(), "px", nv_units), nv_units)); Geom::Translate const tr( Geom::Point(0, old_height - rect_with_margins.height()) -- cgit v1.2.3 From f369e0ed4157c33e09c426b40ec1e39c7591d2c6 Mon Sep 17 00:00:00 2001 From: "Johan B. C. Engelen" Date: Sun, 12 Oct 2014 12:04:17 +0200 Subject: LPEVonKoch: improve rev. 13589 by removing the else branch for complexity calculation for better human readability (bzr r13591) --- src/live_effects/lpe-vonkoch.cpp | 6 ++---- 1 file changed, 2 insertions(+), 4 deletions(-) (limited to 'src') diff --git a/src/live_effects/lpe-vonkoch.cpp b/src/live_effects/lpe-vonkoch.cpp index 7e396e84a..709c9e56e 100644 --- a/src/live_effects/lpe-vonkoch.cpp +++ b/src/live_effects/lpe-vonkoch.cpp @@ -132,8 +132,8 @@ LPEVonKoch::doEffect_path (std::vector const & path_in) int path_in_complexity = 0; for (unsigned k = 0; k < path_in.size(); k++){ path_in_complexity+=path_in[k].size(); - } - double complexity; + } + double complexity = std::pow(transforms.size(), nbgenerations) * path_in_complexity; if (drawall.get_value()){ int k = transforms.size(); if(k>1){ @@ -141,8 +141,6 @@ LPEVonKoch::doEffect_path (std::vector const & path_in) }else{ complexity = nbgenerations*k*path_in_complexity; } - }else{ - complexity = std::pow(transforms.size(),nbgenerations)*path_in_complexity; } if (complexity > double(maxComplexity)){ g_warning("VonKoch lpe's output too complex. Effect bypassed."); -- cgit v1.2.3 From 6fda35fe16784f75a6a98b4d44141177616a06ef Mon Sep 17 00:00:00 2001 From: "Johan B. C. Engelen" Date: Sun, 12 Oct 2014 12:36:25 +0200 Subject: Mass fix whitespace (clang-format) before working on bug (bzr r13593) --- src/trace/potrace/trace.cpp | 2059 ++++++++++++++++++++++--------------------- 1 file changed, 1040 insertions(+), 1019 deletions(-) (limited to 'src') diff --git a/src/trace/potrace/trace.cpp b/src/trace/potrace/trace.cpp index f1e88a908..2bdefb90a 100644 --- a/src/trace/potrace/trace.cpp +++ b/src/trace/potrace/trace.cpp @@ -16,124 +16,129 @@ #include "trace.h" #include "progress.h" -#define INFTY 10000000 /* it suffices that this is longer than any - path; it need not be really infinite */ -#define COS179 -0.999847695156 /* the cosine of 179 degrees */ +#define INFTY 10000000 // it suffices that this is longer than any + // path; it need not be really infinite +#define COS179 -0.999847695156 // the cosine of 179 degrees /* ---------------------------------------------------------------------- */ #define SAFE_MALLOC(var, n, typ) \ - if ((var = (typ *)malloc((n)*sizeof(typ))) == NULL) goto malloc_error + if ((var = (typ *)malloc((n)*sizeof(typ))) == NULL) goto malloc_error /* ---------------------------------------------------------------------- */ /* auxiliary functions */ /* return a direction that is 90 degrees counterclockwise from p2-p0, but then restricted to one of the major wind directions (n, nw, w, etc) */ -static inline point_t dorth_infty(dpoint_t p0, dpoint_t p2) { - point_t r; - - r.y = sign(p2.x-p0.x); - r.x = -sign(p2.y-p0.y); +static inline point_t dorth_infty(dpoint_t p0, dpoint_t p2) +{ + point_t r; + + r.y = sign(p2.x - p0.x); + r.x = -sign(p2.y - p0.y); - return r; + return r; } /* return (p1-p0)x(p2-p0), the area of the parallelogram */ -static inline double dpara(dpoint_t p0, dpoint_t p1, dpoint_t p2) { - double x1, y1, x2, y2; +static inline double dpara(dpoint_t p0, dpoint_t p1, dpoint_t p2) +{ + double x1, y1, x2, y2; - x1 = p1.x-p0.x; - y1 = p1.y-p0.y; - x2 = p2.x-p0.x; - y2 = p2.y-p0.y; + x1 = p1.x - p0.x; + y1 = p1.y - p0.y; + x2 = p2.x - p0.x; + y2 = p2.y - p0.y; - return x1*y2 - x2*y1; + return x1 * y2 - x2 * y1; } /* ddenom/dpara have the property that the square of radius 1 centered at p1 intersects the line p0p2 iff |dpara(p0,p1,p2)| <= ddenom(p0,p2) */ -static inline double ddenom(dpoint_t p0, dpoint_t p2) { - point_t r = dorth_infty(p0, p2); +static inline double ddenom(dpoint_t p0, dpoint_t p2) +{ + point_t r = dorth_infty(p0, p2); - return r.y*(p2.x-p0.x) - r.x*(p2.y-p0.y); + return r.y * (p2.x - p0.x) - r.x * (p2.y - p0.y); } /* return 1 if a <= b < c < a, in a cyclic sense (mod n) */ -static inline int cyclic(int a, int b, int c) { - if (a<=c) { - return (a<=b && blen; - sums_t *sums = pp->sums; - - double x, y, x2, xy, y2; - double k; - double a, b, c, lambda2, l; - int r=0; /* rotations from i to j */ - - while (j>=n) { - j-=n; - r+=1; - } - while (i>=n) { - i-=n; - r-=1; - } - while (j<0) { - j+=n; - r-=1; - } - while (i<0) { - i+=n; - r+=1; - } - - x = sums[j+1].x-sums[i].x+r*sums[n].x; - y = sums[j+1].y-sums[i].y+r*sums[n].y; - x2 = sums[j+1].x2-sums[i].x2+r*sums[n].x2; - xy = sums[j+1].xy-sums[i].xy+r*sums[n].xy; - y2 = sums[j+1].y2-sums[i].y2+r*sums[n].y2; - k = j+1-i+r*n; - - ctr->x = x/k; - ctr->y = y/k; - - a = (x2-(double)x*x/k)/k; - b = (xy-(double)x*y/k)/k; - c = (y2-(double)y*y/k)/k; - - lambda2 = (a+c+sqrt((a-c)*(a-c)+4*b*b))/2; /* larger e.value */ - - /* now find e.vector for lambda2 */ - a -= lambda2; - c -= lambda2; - - if (fabs(a) >= fabs(c)) { - l = sqrt(a*a+b*b); - if (l!=0) { - dir->x = -b/l; - dir->y = a/l; +static void pointslope(privpath_t *pp, int i, int j, dpoint_t *ctr, dpoint_t *dir) +{ + /* assume ilen; + sums_t *sums = pp->sums; + + double x, y, x2, xy, y2; + double k; + double a, b, c, lambda2, l; + int r = 0; /* rotations from i to j */ + + while (j >= n) { + j -= n; + r += 1; + } + while (i >= n) { + i -= n; + r -= 1; + } + while (j < 0) { + j += n; + r -= 1; + } + while (i < 0) { + i += n; + r += 1; + } + + x = sums[j + 1].x - sums[i].x + r * sums[n].x; + y = sums[j + 1].y - sums[i].y + r * sums[n].y; + x2 = sums[j + 1].x2 - sums[i].x2 + r * sums[n].x2; + xy = sums[j + 1].xy - sums[i].xy + r * sums[n].xy; + y2 = sums[j + 1].y2 - sums[i].y2 + r * sums[n].y2; + k = j + 1 - i + r * n; + + ctr->x = x / k; + ctr->y = y / k; + + a = (x2 - (double)x * x / k) / k; + b = (xy - (double)x * y / k) / k; + c = (y2 - (double)y * y / k) / k; + + lambda2 = (a + c + sqrt((a - c) * (a - c) + 4 * b * b)) / 2; /* larger e.value */ + + /* now find e.vector for lambda2 */ + a -= lambda2; + c -= lambda2; + + if (fabs(a) >= fabs(c)) { + l = sqrt(a * a + b * b); + if (l != 0) { + dir->x = -b / l; + dir->y = a / l; + } + } else { + l = sqrt(c * c + b * b); + if (l != 0) { + dir->x = -c / l; + dir->y = b / l; + } } - } else { - l = sqrt(c*c+b*b); - if (l!=0) { - dir->x = -c/l; - dir->y = b/l; + if (l == 0) { + dir->x = dir->y = 0; /* sometimes this can happen when k=4: + the two eigenvalues coincide */ } - } - if (l==0) { - dir->x = dir->y = 0; /* sometimes this can happen when k=4: - the two eigenvalues coincide */ - } } /* the type of (affine) quadratic forms, represented as symmetric 3x3 @@ -142,155 +147,158 @@ static void pointslope(privpath_t *pp, int i, int j, dpoint_t *ctr, dpoint_t *di typedef double quadform_t[3][3]; /* Apply quadratic form Q to vector w = (w.x,w.y) */ -static inline double quadform(quadform_t Q, dpoint_t w) { - double v[3]; - int i, j; - double sum; - - v[0] = w.x; - v[1] = w.y; - v[2] = 1; - sum = 0.0; - - for (i=0; i<3; i++) { - for (j=0; j<3; j++) { - sum += v[i] * Q[i][j] * v[j]; +static inline double quadform(quadform_t Q, dpoint_t w) +{ + double v[3]; + int i, j; + double sum; + + v[0] = w.x; + v[1] = w.y; + v[2] = 1; + sum = 0.0; + + for (i = 0; i < 3; i++) { + for (j = 0; j < 3; j++) { + sum += v[i] * Q[i][j] * v[j]; + } } - } - return sum; + return sum; } /* calculate p1 x p2 */ -static inline int xprod(point_t p1, point_t p2) { - return p1.x*p2.y - p1.y*p2.x; -} +static inline int xprod(point_t p1, point_t p2) { return p1.x * p2.y - p1.y * p2.x; } /* calculate (p1-p0)x(p3-p2) */ -static inline double cprod(dpoint_t p0, dpoint_t p1, dpoint_t p2, dpoint_t p3) { - double x1, y1, x2, y2; +static inline double cprod(dpoint_t p0, dpoint_t p1, dpoint_t p2, dpoint_t p3) +{ + double x1, y1, x2, y2; - x1 = p1.x - p0.x; - y1 = p1.y - p0.y; - x2 = p3.x - p2.x; - y2 = p3.y - p2.y; + x1 = p1.x - p0.x; + y1 = p1.y - p0.y; + x2 = p3.x - p2.x; + y2 = p3.y - p2.y; - return x1*y2 - x2*y1; + return x1 * y2 - x2 * y1; } /* calculate (p1-p0)*(p2-p0) */ -static inline double iprod(dpoint_t p0, dpoint_t p1, dpoint_t p2) { - double x1, y1, x2, y2; +static inline double iprod(dpoint_t p0, dpoint_t p1, dpoint_t p2) +{ + double x1, y1, x2, y2; - x1 = p1.x - p0.x; - y1 = p1.y - p0.y; - x2 = p2.x - p0.x; - y2 = p2.y - p0.y; + x1 = p1.x - p0.x; + y1 = p1.y - p0.y; + x2 = p2.x - p0.x; + y2 = p2.y - p0.y; - return x1*x2 + y1*y2; + return x1 * x2 + y1 * y2; } /* calculate (p1-p0)*(p3-p2) */ -static inline double iprod1(dpoint_t p0, dpoint_t p1, dpoint_t p2, dpoint_t p3) { - double x1, y1, x2, y2; +static inline double iprod1(dpoint_t p0, dpoint_t p1, dpoint_t p2, dpoint_t p3) +{ + double x1, y1, x2, y2; - x1 = p1.x - p0.x; - y1 = p1.y - p0.y; - x2 = p3.x - p2.x; - y2 = p3.y - p2.y; + x1 = p1.x - p0.x; + y1 = p1.y - p0.y; + x2 = p3.x - p2.x; + y2 = p3.y - p2.y; - return x1*x2 + y1*y2; + return x1 * x2 + y1 * y2; } /* calculate distance between two points */ -static inline double ddist(dpoint_t p, dpoint_t q) { - return sqrt(sq(p.x-q.x)+sq(p.y-q.y)); -} +static inline double ddist(dpoint_t p, dpoint_t q) { return sqrt(sq(p.x - q.x) + sq(p.y - q.y)); } /* calculate point of a bezier curve */ -static inline dpoint_t bezier(double t, dpoint_t p0, dpoint_t p1, dpoint_t p2, dpoint_t p3) { - double s = 1-t; - dpoint_t res; +static inline dpoint_t bezier(double t, dpoint_t p0, dpoint_t p1, dpoint_t p2, dpoint_t p3) +{ + double s = 1 - t; + dpoint_t res; - /* Note: a good optimizing compiler (such as gcc-3) reduces the - following to 16 multiplications, using common subexpression - elimination. */ + /* Note: a good optimizing compiler (such as gcc-3) reduces the + following to 16 multiplications, using common subexpression + elimination. */ - res.x = s*s*s*p0.x + 3*(s*s*t)*p1.x + 3*(t*t*s)*p2.x + t*t*t*p3.x; - res.y = s*s*s*p0.y + 3*(s*s*t)*p1.y + 3*(t*t*s)*p2.y + t*t*t*p3.y; + res.x = s * s * s * p0.x + 3 * (s * s * t) * p1.x + 3 * (t * t * s) * p2.x + t * t * t * p3.x; + res.y = s * s * s * p0.y + 3 * (s * s * t) * p1.y + 3 * (t * t * s) * p2.y + t * t * t * p3.y; - return res; + return res; } /* calculate the point t in [0..1] on the (convex) bezier curve (p0,p1,p2,p3) which is tangent to q1-q0. Return -1.0 if there is no solution in [0..1]. */ -static double tangent(dpoint_t p0, dpoint_t p1, dpoint_t p2, dpoint_t p3, dpoint_t q0, dpoint_t q1) { - double A, B, C; /* (1-t)^2 A + 2(1-t)t B + t^2 C = 0 */ - double a, b, c; /* a t^2 + b t + c = 0 */ - double d, s, r1, r2; - - A = cprod(p0, p1, q0, q1); - B = cprod(p1, p2, q0, q1); - C = cprod(p2, p3, q0, q1); - - a = A - 2*B + C; - b = -2*A + 2*B; - c = A; - - d = b*b - 4*a*c; - - if (a==0 || d<0) { - return -1.0; - } - - s = sqrt(d); - - r1 = (-b + s) / (2 * a); - r2 = (-b - s) / (2 * a); - - if (r1 >= 0 && r1 <= 1) { - return r1; - } else if (r2 >= 0 && r2 <= 1) { - return r2; - } else { - return -1.0; - } +static double tangent(dpoint_t p0, dpoint_t p1, dpoint_t p2, dpoint_t p3, dpoint_t q0, dpoint_t q1) +{ + double A, B, C; /* (1-t)^2 A + 2(1-t)t B + t^2 C = 0 */ + double a, b, c; /* a t^2 + b t + c = 0 */ + double d, s, r1, r2; + + A = cprod(p0, p1, q0, q1); + B = cprod(p1, p2, q0, q1); + C = cprod(p2, p3, q0, q1); + + a = A - 2 * B + C; + b = -2 * A + 2 * B; + c = A; + + d = b * b - 4 * a * c; + + if (a == 0 || d < 0) { + return -1.0; + } + + s = sqrt(d); + + r1 = (-b + s) / (2 * a); + r2 = (-b - s) / (2 * a); + + if (r1 >= 0 && r1 <= 1) { + return r1; + } else if (r2 >= 0 && r2 <= 1) { + return r2; + } else { + return -1.0; + } } /* ---------------------------------------------------------------------- */ /* Preparation: fill in the sum* fields of a path (used for later rapid summing). Return 0 on success, 1 with errno set on failure. */ -static int calc_sums(privpath_t *pp) { - int i, x, y; - int n = pp->len; - - SAFE_MALLOC(pp->sums, pp->len+1, sums_t); - - /* origin */ - pp->x0 = pp->pt[0].x; - pp->y0 = pp->pt[0].y; - - /* preparatory computation for later fast summing */ - pp->sums[0].x2 = pp->sums[0].xy = pp->sums[0].y2 = pp->sums[0].x = pp->sums[0].y = 0; - for (i=0; ipt[i].x - pp->x0; - y = pp->pt[i].y - pp->y0; - pp->sums[i+1].x = pp->sums[i].x + x; - pp->sums[i+1].y = pp->sums[i].y + y; - pp->sums[i+1].x2 = pp->sums[i].x2 + x*x; - pp->sums[i+1].xy = pp->sums[i].xy + x*y; - pp->sums[i+1].y2 = pp->sums[i].y2 + y*y; - } - return 0; - - malloc_error: - return 1; +static int calc_sums(privpath_t *pp) +{ + int i, x, y; + int n = pp->len; + + SAFE_MALLOC(pp->sums, pp->len + 1, sums_t); + + /* origin */ + pp->x0 = pp->pt[0].x; + pp->y0 = pp->pt[0].y; + + /* preparatory computation for later fast summing */ + pp->sums[0].x2 = pp->sums[0].xy = pp->sums[0].y2 = pp->sums[0].x = pp->sums[0].y = 0; + for (i = 0; i < n; i++) { + x = pp->pt[i].x - pp->x0; + y = pp->pt[i].y - pp->y0; + pp->sums[i + 1].x = pp->sums[i].x + x; + pp->sums[i + 1].y = pp->sums[i].y + y; + pp->sums[i + 1].x2 = pp->sums[i].x2 + x * x; + pp->sums[i + 1].xy = pp->sums[i].xy + x * y; + pp->sums[i + 1].y2 = pp->sums[i].y2 + y * y; + } + return 0; + +malloc_error: + return 1; } /* ---------------------------------------------------------------------- */ /* Stage 1: determine the straight subpaths (Sec. 2.2.1). Fill in the - "lon" component of a path object (based on pt/len). For each i, + "lon" component of a path object (based on pt/len). For each i, lon[i] is the furthest index such that a straight line can be drawn from i to lon[i]. Return 1 on error with errno set, else 0. */ @@ -318,206 +326,208 @@ static int calc_sums(privpath_t *pp) { substantial. */ /* returns 0 on success, 1 on error with errno set */ -static int calc_lon(privpath_t *pp) { - point_t *pt = pp->pt; - int n = pp->len; - int i, j, k, k1; - int ct[4], dir; - point_t constraint[2]; - point_t cur; - point_t off; - int *pivk = NULL; /* pivk[n] */ - int *nc = NULL; /* nc[n]: next corner */ - point_t dk; /* direction of k-k1 */ - int a, b, c, d; - - SAFE_MALLOC(pivk, n, int); - SAFE_MALLOC(nc, n, int); - - /* initialize the nc data structure. Point from each point to the - furthest future point to which it is connected by a vertical or - horizontal segment. We take advantage of the fact that there is - always a direction change at 0 (due to the path decomposition - algorithm). But even if this were not so, there is no harm, as - in practice, correctness does not depend on the word "furthest" - above. */ - k = 0; - for (i=n-1; i>=0; i--) { - if (pt[i].x != pt[k].x && pt[i].y != pt[k].y) { - k = i+1; /* necessarily ilon, n, int); - - /* determine pivot points: for each i, let pivk[i] be the furthest k - such that all j with i=0; i--) { - ct[0] = ct[1] = ct[2] = ct[3] = 0; - - /* keep track of "directions" that have occurred */ - dir = (3+3*(pt[mod(i+1,n)].x-pt[i].x)+(pt[mod(i+1,n)].y-pt[i].y))/2; - ct[dir]++; - - constraint[0].x = 0; - constraint[0].y = 0; - constraint[1].x = 0; - constraint[1].y = 0; - - /* find the next k such that no straight line from i to k */ - k = nc[i]; - k1 = i; - while (1) { - - dir = (3+3*sign(pt[k].x-pt[k1].x)+sign(pt[k].y-pt[k1].y))/2; - ct[dir]++; - - /* if all four "directions" have occurred, cut this path */ - if (ct[0] && ct[1] && ct[2] && ct[3]) { - pivk[i] = k1; - goto foundk; - } - - cur.x = pt[k].x - pt[i].x; - cur.y = pt[k].y - pt[i].y; - - /* see if current constraint is violated */ - if (xprod(constraint[0], cur) < 0 || xprod(constraint[1], cur) > 0) { - goto constraint_viol; - } - - /* else, update constraint */ - if (abs(cur.x) <= 1 && abs(cur.y) <= 1) { - /* no constraint */ - } else { - off.x = cur.x + ((cur.y>=0 && (cur.y>0 || cur.x<0)) ? 1 : -1); - off.y = cur.y + ((cur.x<=0 && (cur.x<0 || cur.y<0)) ? 1 : -1); - if (xprod(constraint[0], off) >= 0) { - constraint[0] = off; - } - off.x = cur.x + ((cur.y<=0 && (cur.y<0 || cur.x<0)) ? 1 : -1); - off.y = cur.y + ((cur.x>=0 && (cur.x>0 || cur.y<0)) ? 1 : -1); - if (xprod(constraint[1], off) <= 0) { - constraint[1] = off; - } - } - k1 = k; - k = nc[k1]; - if (!cyclic(k,i,k1)) { - break; - } - } - constraint_viol: - /* k1 was the last "corner" satisfying the current constraint, and - k is the first one violating it. We now need to find the last - point along k1..k which satisfied the constraint. */ - dk.x = sign(pt[k].x-pt[k1].x); - dk.y = sign(pt[k].y-pt[k1].y); - cur.x = pt[k1].x - pt[i].x; - cur.y = pt[k1].y - pt[i].y; - /* find largest integer j such that xprod(constraint[0], cur+j*dk) - >= 0 and xprod(constraint[1], cur+j*dk) <= 0. Use bilinearity - of xprod. */ - a = xprod(constraint[0], cur); - b = xprod(constraint[0], dk); - c = xprod(constraint[1], cur); - d = xprod(constraint[1], dk); - /* find largest integer j such that a+j*b>=0 and c+j*d<=0. This - can be solved with integer arithmetic. */ - j = INFTY; - if (b<0) { - j = floordiv(a,-b); - } - if (d>0) { - j = min(j, floordiv(-c,d)); +static int calc_lon(privpath_t *pp) +{ + point_t *pt = pp->pt; + int n = pp->len; + int i, j, k, k1; + int ct[4], dir; + point_t constraint[2]; + point_t cur; + point_t off; + int *pivk = NULL; /* pivk[n] */ + int *nc = NULL; /* nc[n]: next corner */ + point_t dk; /* direction of k-k1 */ + int a, b, c, d; + + SAFE_MALLOC(pivk, n, int); + SAFE_MALLOC(nc, n, int); + + /* initialize the nc data structure. Point from each point to the + furthest future point to which it is connected by a vertical or + horizontal segment. We take advantage of the fact that there is + always a direction change at 0 (due to the path decomposition + algorithm). But even if this were not so, there is no harm, as + in practice, correctness does not depend on the word "furthest" + above. */ + k = 0; + for (i = n - 1; i >= 0; i--) { + if (pt[i].x != pt[k].x && pt[i].y != pt[k].y) { + k = i + 1; /* necessarily ilon[n-1]=j; - for (i=n-2; i>=0; i--) { - if (cyclic(i+1,pivk[i],j)) { - j=pivk[i]; + + SAFE_MALLOC(pp->lon, n, int); + + /* determine pivot points: for each i, let pivk[i] be the furthest k + such that all j with i= 0; i--) { + ct[0] = ct[1] = ct[2] = ct[3] = 0; + + /* keep track of "directions" that have occurred */ + dir = (3 + 3 * (pt[mod(i + 1, n)].x - pt[i].x) + (pt[mod(i + 1, n)].y - pt[i].y)) / 2; + ct[dir]++; + + constraint[0].x = 0; + constraint[0].y = 0; + constraint[1].x = 0; + constraint[1].y = 0; + + /* find the next k such that no straight line from i to k */ + k = nc[i]; + k1 = i; + while (1) { + + dir = (3 + 3 * sign(pt[k].x - pt[k1].x) + sign(pt[k].y - pt[k1].y)) / 2; + ct[dir]++; + + /* if all four "directions" have occurred, cut this path */ + if (ct[0] && ct[1] && ct[2] && ct[3]) { + pivk[i] = k1; + goto foundk; + } + + cur.x = pt[k].x - pt[i].x; + cur.y = pt[k].y - pt[i].y; + + /* see if current constraint is violated */ + if (xprod(constraint[0], cur) < 0 || xprod(constraint[1], cur) > 0) { + goto constraint_viol; + } + + /* else, update constraint */ + if (abs(cur.x) <= 1 && abs(cur.y) <= 1) { + /* no constraint */ + } else { + off.x = cur.x + ((cur.y >= 0 && (cur.y > 0 || cur.x < 0)) ? 1 : -1); + off.y = cur.y + ((cur.x <= 0 && (cur.x < 0 || cur.y < 0)) ? 1 : -1); + if (xprod(constraint[0], off) >= 0) { + constraint[0] = off; + } + off.x = cur.x + ((cur.y <= 0 && (cur.y < 0 || cur.x < 0)) ? 1 : -1); + off.y = cur.y + ((cur.x >= 0 && (cur.x > 0 || cur.y < 0)) ? 1 : -1); + if (xprod(constraint[1], off) <= 0) { + constraint[1] = off; + } + } + k1 = k; + k = nc[k1]; + if (!cyclic(k, i, k1)) { + break; + } + } + constraint_viol: + /* k1 was the last "corner" satisfying the current constraint, and + k is the first one violating it. We now need to find the last + point along k1..k which satisfied the constraint. */ + dk.x = sign(pt[k].x - pt[k1].x); + dk.y = sign(pt[k].y - pt[k1].y); + cur.x = pt[k1].x - pt[i].x; + cur.y = pt[k1].y - pt[i].y; + /* find largest integer j such that xprod(constraint[0], cur+j*dk) + >= 0 and xprod(constraint[1], cur+j*dk) <= 0. Use bilinearity + of xprod. */ + a = xprod(constraint[0], cur); + b = xprod(constraint[0], dk); + c = xprod(constraint[1], cur); + d = xprod(constraint[1], dk); + /* find largest integer j such that a+j*b>=0 and c+j*d<=0. This + can be solved with integer arithmetic. */ + j = INFTY; + if (b < 0) { + j = floordiv(a, -b); + } + if (d > 0) { + j = min(j, floordiv(-c, d)); + } + pivk[i] = mod(k1 + j, n); + foundk: + ; + } /* for i */ + + /* clean up: for each i, let lon[i] be the largest k such that for + all i' with i<=i'lon[n - 1] = j; + for (i = n - 2; i >= 0; i--) { + if (cyclic(i + 1, pivk[i], j)) { + j = pivk[i]; + } + pp->lon[i] = j; } - pp->lon[i]=j; - } - for (i=n-1; cyclic(mod(i+1,n),j,pp->lon[i]); i--) { - pp->lon[i] = j; - } + for (i = n - 1; cyclic(mod(i + 1, n), j, pp->lon[i]); i--) { + pp->lon[i] = j; + } - free(pivk); - free(nc); - return 0; + free(pivk); + free(nc); + return 0; - malloc_error: - free(pivk); - free(nc); - return 1; +malloc_error: + free(pivk); + free(nc); + return 1; } /* ---------------------------------------------------------------------- */ -/* Stage 2: calculate the optimal polygon (Sec. 2.2.2-2.2.4). */ +/* Stage 2: calculate the optimal polygon (Sec. 2.2.2-2.2.4). */ /* Auxiliary function: calculate the penalty of an edge from i to j in the given path. This needs the "lon" and "sum*" data. */ -static double penalty3(privpath_t *pp, int i, int j) { - int n = pp->len; - point_t *pt = pp->pt; - sums_t *sums = pp->sums; - - /* assume 0<=i=n) { - j -= n; - r = 1; - } - - /* critical inner loop: the "if" gives a 4.6 percent speedup */ - if (r == 0) { - x = sums[j+1].x - sums[i].x; - y = sums[j+1].y - sums[i].y; - x2 = sums[j+1].x2 - sums[i].x2; - xy = sums[j+1].xy - sums[i].xy; - y2 = sums[j+1].y2 - sums[i].y2; - k = j+1 - i; - } else { - x = sums[j+1].x - sums[i].x + sums[n].x; - y = sums[j+1].y - sums[i].y + sums[n].y; - x2 = sums[j+1].x2 - sums[i].x2 + sums[n].x2; - xy = sums[j+1].xy - sums[i].xy + sums[n].xy; - y2 = sums[j+1].y2 - sums[i].y2 + sums[n].y2; - k = j+1 - i + n; - } - - px = (pt[i].x + pt[j].x) / 2.0 - pt[0].x; - py = (pt[i].y + pt[j].y) / 2.0 - pt[0].y; - ey = (pt[j].x - pt[i].x); - ex = -(pt[j].y - pt[i].y); - - a = ((x2 - 2*x*px) / k + px*px); - b = ((xy - x*py - y*px) / k + px*py); - c = ((y2 - 2*y*py) / k + py*py); - - s = ex*ex*a + 2*ex*ey*b + ey*ey*c; - - return sqrt(s); +static double penalty3(privpath_t *pp, int i, int j) +{ + int n = pp->len; + point_t *pt = pp->pt; + sums_t *sums = pp->sums; + + /* assume 0<=i= n) { + j -= n; + r = 1; + } + + /* critical inner loop: the "if" gives a 4.6 percent speedup */ + if (r == 0) { + x = sums[j + 1].x - sums[i].x; + y = sums[j + 1].y - sums[i].y; + x2 = sums[j + 1].x2 - sums[i].x2; + xy = sums[j + 1].xy - sums[i].xy; + y2 = sums[j + 1].y2 - sums[i].y2; + k = j + 1 - i; + } else { + x = sums[j + 1].x - sums[i].x + sums[n].x; + y = sums[j + 1].y - sums[i].y + sums[n].y; + x2 = sums[j + 1].x2 - sums[i].x2 + sums[n].x2; + xy = sums[j + 1].xy - sums[i].xy + sums[n].xy; + y2 = sums[j + 1].y2 - sums[i].y2 + sums[n].y2; + k = j + 1 - i + n; + } + + px = (pt[i].x + pt[j].x) / 2.0 - pt[0].x; + py = (pt[i].y + pt[j].y) / 2.0 - pt[0].y; + ey = (pt[j].x - pt[i].x); + ex = -(pt[j].y - pt[i].y); + + a = ((x2 - 2 * x * px) / k + px * px); + b = ((xy - x * py - y * px) / k + px * py); + c = ((y2 - 2 * y * py) / k + py * py); + + s = ex * ex * a + 2 * ex * ey * b + ey * ey * c; + + return sqrt(s); } /* find the optimal polygon. Fill in the m and po components. Return 1 @@ -525,109 +535,109 @@ static double penalty3(privpath_t *pp, int i, int j) { is in the polygon. Fixme: implement cyclic version. */ static int bestpolygon(privpath_t *pp) { - int i, j, m, k; - int n = pp->len; - double *pen = NULL; /* pen[n+1]: penalty vector */ - int *prev = NULL; /* prev[n+1]: best path pointer vector */ - int *clip0 = NULL; /* clip0[n]: longest segment pointer, non-cyclic */ - int *clip1 = NULL; /* clip1[n+1]: backwards segment pointer, non-cyclic */ - int *seg0 = NULL; /* seg0[m+1]: forward segment bounds, m<=n */ - int *seg1 = NULL; /* seg1[m+1]: backward segment bounds, m<=n */ - double thispen; - double best; - int c; - - SAFE_MALLOC(pen, n+1, double); - SAFE_MALLOC(prev, n+1, int); - SAFE_MALLOC(clip0, n, int); - SAFE_MALLOC(clip1, n+1, int); - SAFE_MALLOC(seg0, n+1, int); - SAFE_MALLOC(seg1, n+1, int); - - /* calculate clipped paths */ - for (i=0; ilon[mod(i-1,n)]-1,n); - if (c == i) { - c = mod(i+1,n); + int i, j, m, k; + int n = pp->len; + double *pen = NULL; /* pen[n+1]: penalty vector */ + int *prev = NULL; /* prev[n+1]: best path pointer vector */ + int *clip0 = NULL; /* clip0[n]: longest segment pointer, non-cyclic */ + int *clip1 = NULL; /* clip1[n+1]: backwards segment pointer, non-cyclic */ + int *seg0 = NULL; /* seg0[m+1]: forward segment bounds, m<=n */ + int *seg1 = NULL; /* seg1[m+1]: backward segment bounds, m<=n */ + double thispen; + double best; + int c; + + SAFE_MALLOC(pen, n + 1, double); + SAFE_MALLOC(prev, n + 1, int); + SAFE_MALLOC(clip0, n, int); + SAFE_MALLOC(clip1, n + 1, int); + SAFE_MALLOC(seg0, n + 1, int); + SAFE_MALLOC(seg1, n + 1, int); + + /* calculate clipped paths */ + for (i = 0; i < n; i++) { + c = mod(pp->lon[mod(i - 1, n)] - 1, n); + if (c == i) { + c = mod(i + 1, n); + } + if (c < i) { + clip0[i] = n; + } else { + clip0[i] = c; + } } - if (c < i) { - clip0[i] = n; - } else { - clip0[i] = c; + + /* calculate backwards path clipping, non-cyclic. j <= clip0[i] iff + clip1[j] <= i, for i,j=0..n. */ + j = 1; + for (i = 0; i < n; i++) { + while (j <= clip0[i]) { + clip1[j] = i; + j++; + } + } + + /* calculate seg0[j] = longest path from 0 with j segments */ + i = 0; + for (j = 0; i < n; j++) { + seg0[j] = i; + i = clip0[i]; + } + seg0[j] = n; + m = j; + + /* calculate seg1[j] = longest path to n with m-j segments */ + i = n; + for (j = m; j > 0; j--) { + seg1[j] = i; + i = clip1[i]; } - } - - /* calculate backwards path clipping, non-cyclic. j <= clip0[i] iff - clip1[j] <= i, for i,j=0..n. */ - j = 1; - for (i=0; i= clip1[i]; k--) { + thispen = penalty3(pp, k, i) + pen[k]; + if (best < 0 || thispen < best) { + prev[i] = k; + best = thispen; + } + } + pen[i] = best; + } } - } - - /* calculate seg0[j] = longest path from 0 with j segments */ - i = 0; - for (j=0; i0; j--) { - seg1[j] = i; - i = clip1[i]; - } - seg1[0] = 0; - - /* now find the shortest path with m segments, based on penalty3 */ - /* note: the outer 2 loops jointly have at most n iterations, thus - the worst-case behavior here is quadratic. In practice, it is - close to linear since the inner loop tends to be short. */ - pen[0]=0; - for (j=1; j<=m; j++) { - for (i=seg1[j]; i<=seg0[j]; i++) { - best = -1; - for (k=seg0[j-1]; k>=clip1[i]; k--) { - thispen = penalty3(pp, k, i) + pen[k]; - if (best < 0 || thispen < best) { - prev[i] = k; - best = thispen; - } - } - pen[i] = best; + + pp->m = m; + SAFE_MALLOC(pp->po, m, int); + + /* read off shortest path */ + for (i = n, j = m - 1; i > 0; j--) { + i = prev[i]; + pp->po[j] = i; } - } - - pp->m = m; - SAFE_MALLOC(pp->po, m, int); - - /* read off shortest path */ - for (i=n, j=m-1; i>0; j--) { - i = prev[i]; - pp->po[j] = i; - } - - free(pen); - free(prev); - free(clip0); - free(clip1); - free(seg0); - free(seg1); - return 0; - - malloc_error: - free(pen); - free(prev); - free(clip0); - free(clip1); - free(seg0); - free(seg1); - return 1; + + free(pen); + free(prev); + free(clip0); + free(clip1); + free(seg0); + free(seg1); + return 0; + +malloc_error: + free(pen); + free(prev); + free(clip0); + free(clip1); + free(seg0); + free(seg1); + return 1; } /* ---------------------------------------------------------------------- */ @@ -638,266 +648,269 @@ static int bestpolygon(privpath_t *pp) if it lies outside. Return 1 with errno set on error; 0 on success. */ -static int adjust_vertices(privpath_t *pp) { - int m = pp->m; - int *po = pp->po; - int n = pp->len; - point_t *pt = pp->pt; - int x0 = pp->x0; - int y0 = pp->y0; - - dpoint_t *ctr = NULL; /* ctr[m] */ - dpoint_t *dir = NULL; /* dir[m] */ - quadform_t *q = NULL; /* q[m] */ - double v[3]; - double d; - int i, j, k, l; - dpoint_t s; - int r; - - SAFE_MALLOC(ctr, m, dpoint_t); - SAFE_MALLOC(dir, m, dpoint_t); - SAFE_MALLOC(q, m, quadform_t); - - r = privcurve_init(&pp->curve, m); - if (r) { - goto malloc_error; - } - - /* calculate "optimal" point-slope representation for each line - segment */ - for (i=0; im; + int *po = pp->po; + int n = pp->len; + point_t *pt = pp->pt; + int x0 = pp->x0; + int y0 = pp->y0; + + dpoint_t *ctr = NULL; /* ctr[m] */ + dpoint_t *dir = NULL; /* dir[m] */ + quadform_t *q = NULL; /* q[m] */ + double v[3]; + double d; + int i, j, k, l; + dpoint_t s; + int r; + + SAFE_MALLOC(ctr, m, dpoint_t); + SAFE_MALLOC(dir, m, dpoint_t); + SAFE_MALLOC(q, m, quadform_t); + + r = privcurve_init(&pp->curve, m); + if (r) { + goto malloc_error; } - - while(1) { - /* minimize the quadratic form Q on the unit square */ - /* find intersection */ -#ifdef HAVE_GCC_LOOP_BUG - /* work around gcc bug #12243 */ - free(NULL); -#endif - - det = Q[0][0]*Q[1][1] - Q[0][1]*Q[1][0]; - if (det != 0.0) { - w.x = (-Q[0][2]*Q[1][1] + Q[1][2]*Q[0][1]) / det; - w.y = ( Q[0][2]*Q[1][0] - Q[1][2]*Q[0][0]) / det; - break; - } - - /* matrix is singular - lines are parallel. Add another, - orthogonal axis, through the center of the unit square */ - if (Q[0][0]>Q[1][1]) { - v[0] = -Q[0][1]; - v[1] = Q[0][0]; - } else if (Q[1][1]) { - v[0] = -Q[1][1]; - v[1] = Q[1][0]; - } else { - v[0] = 1; - v[1] = 0; - } - d = sq(v[0]) + sq(v[1]); - v[2] = - v[1] * s.y - v[0] * s.x; - for (l=0; l<3; l++) { - for (k=0; k<3; k++) { - Q[l][k] += v[l] * v[k] / d; - } - } + /* calculate "optimal" point-slope representation for each line + segment */ + for (i = 0; i < m; i++) { + j = po[mod(i + 1, m)]; + j = mod(j - po[i], n) + po[i]; + pointslope(pp, po[i], j, &ctr[i], &dir[i]); } - dx = fabs(w.x-s.x); - dy = fabs(w.y-s.y); - if (dx <= .5 && dy <= .5) { - pp->curve.vertex[i].x = w.x+x0; - pp->curve.vertex[i].y = w.y+y0; - continue; + + /* represent each line segment as a singular quadratic form; the + distance of a point (x,y) from the line segment will be + (x,y,1)Q(x,y,1)^t, where Q=q[i]. */ + for (i = 0; i < m; i++) { + d = sq(dir[i].x) + sq(dir[i].y); + if (d == 0.0) { + for (j = 0; j < 3; j++) { + for (k = 0; k < 3; k++) { + q[i][j][k] = 0; + } + } + } else { + v[0] = dir[i].y; + v[1] = -dir[i].x; + v[2] = -v[1] * ctr[i].y - v[0] * ctr[i].x; + for (l = 0; l < 3; l++) { + for (k = 0; k < 3; k++) { + q[i][l][k] = v[l] * v[k] / d; + } + } + } } - /* the minimum was not in the unit square; now minimize quadratic - on boundary of square */ - min = quadform(Q, s); - xmin = s.x; - ymin = s.y; + /* now calculate the "intersections" of consecutive segments. + Instead of using the actual intersection, we find the point + within a given unit square which minimizes the square distance to + the two lines. */ + for (i = 0; i < m; i++) { + quadform_t Q; + dpoint_t w; + double dx, dy; + double det; + double min, cand; /* minimum and candidate for minimum of quad. form */ + double xmin, ymin; /* coordinates of minimum */ + int z; + + /* let s be the vertex, in coordinates relative to x0/y0 */ + s.x = pt[po[i]].x - x0; + s.y = pt[po[i]].y - y0; + + /* intersect segments i-1 and i */ + + j = mod(i - 1, m); + + /* add quadratic forms */ + for (l = 0; l < 3; l++) { + for (k = 0; k < 3; k++) { + Q[l][k] = q[j][l][k] + q[i][l][k]; + } + } + + while (1) { +/* minimize the quadratic form Q on the unit square */ +/* find intersection */ - if (Q[0][0] == 0.0) { - goto fixx; - } - for (z=0; z<2; z++) { /* value of the y-coordinate */ - w.y = s.y-0.5+z; - w.x = - (Q[0][1] * w.y + Q[0][2]) / Q[0][0]; - dx = fabs(w.x-s.x); - cand = quadform(Q, w); - if (dx <= .5 && cand < min) { - min = cand; - xmin = w.x; - ymin = w.y; - } - } - fixx: - if (Q[1][1] == 0.0) { - goto corners; - } - for (z=0; z<2; z++) { /* value of the x-coordinate */ - w.x = s.x-0.5+z; - w.y = - (Q[1][0] * w.x + Q[1][2]) / Q[1][1]; - dy = fabs(w.y-s.y); - cand = quadform(Q, w); - if (dy <= .5 && cand < min) { - min = cand; - xmin = w.x; - ymin = w.y; - } - } - corners: - /* check four corners */ - for (l=0; l<2; l++) { - for (k=0; k<2; k++) { - w.x = s.x-0.5+l; - w.y = s.y-0.5+k; - cand = quadform(Q, w); - if (cand < min) { - min = cand; - xmin = w.x; - ymin = w.y; - } - } +#ifdef HAVE_GCC_LOOP_BUG + /* work around gcc bug #12243 */ + free(NULL); +#endif + + det = Q[0][0] * Q[1][1] - Q[0][1] * Q[1][0]; + if (det != 0.0) { + w.x = (-Q[0][2] * Q[1][1] + Q[1][2] * Q[0][1]) / det; + w.y = (Q[0][2] * Q[1][0] - Q[1][2] * Q[0][0]) / det; + break; + } + + /* matrix is singular - lines are parallel. Add another, + orthogonal axis, through the center of the unit square */ + if (Q[0][0] > Q[1][1]) { + v[0] = -Q[0][1]; + v[1] = Q[0][0]; + } else if (Q[1][1]) { + v[0] = -Q[1][1]; + v[1] = Q[1][0]; + } else { + v[0] = 1; + v[1] = 0; + } + d = sq(v[0]) + sq(v[1]); + v[2] = -v[1] * s.y - v[0] * s.x; + for (l = 0; l < 3; l++) { + for (k = 0; k < 3; k++) { + Q[l][k] += v[l] * v[k] / d; + } + } + } + dx = fabs(w.x - s.x); + dy = fabs(w.y - s.y); + if (dx <= .5 && dy <= .5) { + pp->curve.vertex[i].x = w.x + x0; + pp->curve.vertex[i].y = w.y + y0; + continue; + } + + /* the minimum was not in the unit square; now minimize quadratic + on boundary of square */ + min = quadform(Q, s); + xmin = s.x; + ymin = s.y; + + if (Q[0][0] == 0.0) { + goto fixx; + } + for (z = 0; z < 2; z++) { /* value of the y-coordinate */ + w.y = s.y - 0.5 + z; + w.x = -(Q[0][1] * w.y + Q[0][2]) / Q[0][0]; + dx = fabs(w.x - s.x); + cand = quadform(Q, w); + if (dx <= .5 && cand < min) { + min = cand; + xmin = w.x; + ymin = w.y; + } + } + fixx: + if (Q[1][1] == 0.0) { + goto corners; + } + for (z = 0; z < 2; z++) { /* value of the x-coordinate */ + w.x = s.x - 0.5 + z; + w.y = -(Q[1][0] * w.x + Q[1][2]) / Q[1][1]; + dy = fabs(w.y - s.y); + cand = quadform(Q, w); + if (dy <= .5 && cand < min) { + min = cand; + xmin = w.x; + ymin = w.y; + } + } + corners: + /* check four corners */ + for (l = 0; l < 2; l++) { + for (k = 0; k < 2; k++) { + w.x = s.x - 0.5 + l; + w.y = s.y - 0.5 + k; + cand = quadform(Q, w); + if (cand < min) { + min = cand; + xmin = w.x; + ymin = w.y; + } + } + } + + pp->curve.vertex[i].x = xmin + x0; + pp->curve.vertex[i].y = ymin + y0; + continue; } - pp->curve.vertex[i].x = xmin + x0; - pp->curve.vertex[i].y = ymin + y0; - continue; - } - - free(ctr); - free(dir); - free(q); - return 0; - - malloc_error: - free(ctr); - free(dir); - free(q); - return 1; + free(ctr); + free(dir); + free(q); + return 0; + +malloc_error: + free(ctr); + free(dir); + free(q); + return 1; } /* ---------------------------------------------------------------------- */ /* Stage 4: smoothing and corner analysis (Sec. 2.3.3) */ /* reverse orientation of a path */ -static void reverse(privcurve_t *curve) { - int m = curve->n; - int i, j; - dpoint_t tmp; - - for (i=0, j=m-1; ivertex[i]; - curve->vertex[i] = curve->vertex[j]; - curve->vertex[j] = tmp; - } +static void reverse(privcurve_t *curve) +{ + int m = curve->n; + int i, j; + dpoint_t tmp; + + for (i = 0, j = m - 1; i < j; i++, j--) { + tmp = curve->vertex[i]; + curve->vertex[i] = curve->vertex[j]; + curve->vertex[j] = tmp; + } } /* Always succeeds */ -static void smooth(privcurve_t *curve, double alphamax) { - int m = curve->n; - - int i, j, k; - double dd, denom, alpha; - dpoint_t p2, p3, p4; - - /* examine each vertex and find its best fit */ - for (i=0; ivertex[k], curve->vertex[j]); - - denom = ddenom(curve->vertex[i], curve->vertex[k]); - if (denom != 0.0) { - dd = dpara(curve->vertex[i], curve->vertex[j], curve->vertex[k]) / denom; - dd = fabs(dd); - alpha = dd>1 ? (1 - 1.0/dd) : 0; - alpha = alpha / 0.75; - } else { - alpha = 4/3.0; - } - curve->alpha0[j] = alpha; /* remember "original" value of alpha */ - - if (alpha > alphamax) { /* pointed corner */ - curve->tag[j] = POTRACE_CORNER; - curve->c[j][1] = curve->vertex[j]; - curve->c[j][2] = p4; - } else { - if (alpha < 0.55) { - alpha = 0.55; - } else if (alpha > 1) { - alpha = 1; - } - p2 = interval(.5+.5*alpha, curve->vertex[i], curve->vertex[j]); - p3 = interval(.5+.5*alpha, curve->vertex[k], curve->vertex[j]); - curve->tag[j] = POTRACE_CURVETO; - curve->c[j][0] = p2; - curve->c[j][1] = p3; - curve->c[j][2] = p4; +static void smooth(privcurve_t *curve, double alphamax) +{ + int m = curve->n; + + int i, j, k; + double dd, denom, alpha; + dpoint_t p2, p3, p4; + + /* examine each vertex and find its best fit */ + for (i = 0; i < m; i++) { + j = mod(i + 1, m); + k = mod(i + 2, m); + p4 = interval(1 / 2.0, curve->vertex[k], curve->vertex[j]); + + denom = ddenom(curve->vertex[i], curve->vertex[k]); + if (denom != 0.0) { + dd = dpara(curve->vertex[i], curve->vertex[j], curve->vertex[k]) / denom; + dd = fabs(dd); + alpha = dd > 1 ? (1 - 1.0 / dd) : 0; + alpha = alpha / 0.75; + } else { + alpha = 4 / 3.0; + } + curve->alpha0[j] = alpha; /* remember "original" value of alpha */ + + if (alpha > alphamax) { /* pointed corner */ + curve->tag[j] = POTRACE_CORNER; + curve->c[j][1] = curve->vertex[j]; + curve->c[j][2] = p4; + } else { + if (alpha < 0.55) { + alpha = 0.55; + } else if (alpha > 1) { + alpha = 1; + } + p2 = interval(.5 + .5 * alpha, curve->vertex[i], curve->vertex[j]); + p3 = interval(.5 + .5 * alpha, curve->vertex[k], curve->vertex[j]); + curve->tag[j] = POTRACE_CURVETO; + curve->c[j][0] = p2; + curve->c[j][1] = p3; + curve->c[j][2] = p4; + } + curve->alpha[j] = alpha; /* store the "cropped" value of alpha */ + curve->beta[j] = 0.5; } - curve->alpha[j] = alpha; /* store the "cropped" value of alpha */ - curve->beta[j] = 0.5; - } - curve->alphacurve = 1; + curve->alphacurve = 1; - return; + return; } /* ---------------------------------------------------------------------- */ @@ -905,341 +918,349 @@ static void smooth(privcurve_t *curve, double alphamax) { /* a private type for the result of opti_penalty */ struct opti_s { - double pen; /* penalty */ - dpoint_t c[2]; /* curve parameters */ - double t, s; /* curve parameters */ - double alpha; /* curve parameter */ + double pen; /* penalty */ + dpoint_t c[2]; /* curve parameters */ + double t, s; /* curve parameters */ + double alpha; /* curve parameter */ }; typedef struct opti_s opti_t; /* calculate best fit from i+.5 to j+.5. Assume icurve.n; - int k, k1, k2, conv, i1; - double area, alpha, d, d1, d2; - dpoint_t p0, p1, p2, p3, pt; - double A, R, A1, A2, A3, A4; - double s, t; - - /* check convexity, corner-freeness, and maximum bend < 179 degrees */ +static int opti_penalty(privpath_t *pp, int i, int j, opti_t *res, double opttolerance, int *convc, double *areac) +{ + int m = pp->curve.n; + int k, k1, k2, conv, i1; + double area, alpha, d, d1, d2; + dpoint_t p0, p1, p2, p3, pt; + double A, R, A1, A2, A3, A4; + double s, t; - if (i==j) { /* sanity - a full loop can never be an opticurve */ - return 1; - } + /* check convexity, corner-freeness, and maximum bend < 179 degrees */ - k = i; - i1 = mod(i+1, m); - k1 = mod(k+1, m); - conv = convc[k1]; - if (conv == 0) { - return 1; - } - d = ddist(pp->curve.vertex[i], pp->curve.vertex[i1]); - for (k=k1; k!=j; k=k1) { - k1 = mod(k+1, m); - k2 = mod(k+2, m); - if (convc[k1] != conv) { - return 1; + if (i == j) { /* sanity - a full loop can never be an opticurve */ + return 1; } - if (sign(cprod(pp->curve.vertex[i], pp->curve.vertex[i1], pp->curve.vertex[k1], pp->curve.vertex[k2])) != conv) { - return 1; + + k = i; + i1 = mod(i + 1, m); + k1 = mod(k + 1, m); + conv = convc[k1]; + if (conv == 0) { + return 1; } - if (iprod1(pp->curve.vertex[i], pp->curve.vertex[i1], pp->curve.vertex[k1], pp->curve.vertex[k2]) < d * ddist(pp->curve.vertex[k1], pp->curve.vertex[k2]) * COS179) { - return 1; + d = ddist(pp->curve.vertex[i], pp->curve.vertex[i1]); + for (k = k1; k != j; k = k1) { + k1 = mod(k + 1, m); + k2 = mod(k + 2, m); + if (convc[k1] != conv) { + return 1; + } + if (sign(cprod(pp->curve.vertex[i], pp->curve.vertex[i1], pp->curve.vertex[k1], pp->curve.vertex[k2])) != + conv) { + return 1; + } + if (iprod1(pp->curve.vertex[i], pp->curve.vertex[i1], pp->curve.vertex[k1], pp->curve.vertex[k2]) < + d * ddist(pp->curve.vertex[k1], pp->curve.vertex[k2]) * COS179) { + return 1; + } } - } - - /* the curve we're working in: */ - p0 = pp->curve.c[mod(i,m)][2]; - p1 = pp->curve.vertex[mod(i+1,m)]; - p2 = pp->curve.vertex[mod(j,m)]; - p3 = pp->curve.c[mod(j,m)][2]; - - /* determine its area */ - area = areac[j] - areac[i]; - area -= dpara(pp->curve.vertex[0], pp->curve.c[i][2], pp->curve.c[j][2])/2; - if (i>=j) { - area += areac[m]; - } - - /* find intersection o of p0p1 and p2p3. Let t,s such that o = - interval(t,p0,p1) = interval(s,p3,p2). Let A be the area of the - triangle (p0,o,p3). */ - - A1 = dpara(p0, p1, p2); - A2 = dpara(p0, p1, p3); - A3 = dpara(p0, p2, p3); - /* A4 = dpara(p1, p2, p3); */ - A4 = A1+A3-A2; - - if (A2 == A1) { /* this should never happen */ - return 1; - } - t = A3/(A3-A4); - s = A2/(A2-A1); - A = A2 * t / 2.0; - - if (A == 0.0) { /* this should never happen */ - return 1; - } + /* the curve we're working in: */ + p0 = pp->curve.c[mod(i, m)][2]; + p1 = pp->curve.vertex[mod(i + 1, m)]; + p2 = pp->curve.vertex[mod(j, m)]; + p3 = pp->curve.c[mod(j, m)][2]; + + /* determine its area */ + area = areac[j] - areac[i]; + area -= dpara(pp->curve.vertex[0], pp->curve.c[i][2], pp->curve.c[j][2]) / 2; + if (i >= j) { + area += areac[m]; + } - R = area / A; /* relative area */ - alpha = 2 - sqrt(4 - R / 0.3); /* overall alpha for p0-o-p3 curve */ + /* find intersection o of p0p1 and p2p3. Let t,s such that o = + interval(t,p0,p1) = interval(s,p3,p2). Let A be the area of the + triangle (p0,o,p3). */ - res->c[0] = interval(t * alpha, p0, p1); - res->c[1] = interval(s * alpha, p3, p2); - res->alpha = alpha; - res->t = t; - res->s = s; + A1 = dpara(p0, p1, p2); + A2 = dpara(p0, p1, p3); + A3 = dpara(p0, p2, p3); + /* A4 = dpara(p1, p2, p3); */ + A4 = A1 + A3 - A2; - p1 = res->c[0]; - p2 = res->c[1]; /* the proposed curve is now (p0,p1,p2,p3) */ + if (A2 == A1) { /* this should never happen */ + return 1; + } - res->pen = 0; + t = A3 / (A3 - A4); + s = A2 / (A2 - A1); + A = A2 * t / 2.0; - /* calculate penalty */ - /* check tangency with edges */ - for (k=mod(i+1,m); k!=j; k=k1) { - k1 = mod(k+1,m); - t = tangent(p0, p1, p2, p3, pp->curve.vertex[k], pp->curve.vertex[k1]); - if (t<-.5) { - return 1; - } - pt = bezier(t, p0, p1, p2, p3); - d = ddist(pp->curve.vertex[k], pp->curve.vertex[k1]); - if (d == 0.0) { /* this should never happen */ - return 1; + if (A == 0.0) { /* this should never happen */ + return 1; } - d1 = dpara(pp->curve.vertex[k], pp->curve.vertex[k1], pt) / d; - if (fabs(d1) > opttolerance) { - return 1; - } - if (iprod(pp->curve.vertex[k], pp->curve.vertex[k1], pt) < 0 || iprod(pp->curve.vertex[k1], pp->curve.vertex[k], pt) < 0) { - return 1; - } - res->pen += sq(d1); - } - - /* check corners */ - for (k=i; k!=j; k=k1) { - k1 = mod(k+1,m); - t = tangent(p0, p1, p2, p3, pp->curve.c[k][2], pp->curve.c[k1][2]); - if (t<-.5) { - return 1; - } - pt = bezier(t, p0, p1, p2, p3); - d = ddist(pp->curve.c[k][2], pp->curve.c[k1][2]); - if (d == 0.0) { /* this should never happen */ - return 1; - } - d1 = dpara(pp->curve.c[k][2], pp->curve.c[k1][2], pt) / d; - d2 = dpara(pp->curve.c[k][2], pp->curve.c[k1][2], pp->curve.vertex[k1]) / d; - d2 *= 0.75 * pp->curve.alpha[k1]; - if (d2 < 0) { - d1 = -d1; - d2 = -d2; - } - if (d1 < d2 - opttolerance) { - return 1; + + R = area / A; /* relative area */ + alpha = 2 - sqrt(4 - R / 0.3); /* overall alpha for p0-o-p3 curve */ + + res->c[0] = interval(t * alpha, p0, p1); + res->c[1] = interval(s * alpha, p3, p2); + res->alpha = alpha; + res->t = t; + res->s = s; + + p1 = res->c[0]; + p2 = res->c[1]; /* the proposed curve is now (p0,p1,p2,p3) */ + + res->pen = 0; + + /* calculate penalty */ + /* check tangency with edges */ + for (k = mod(i + 1, m); k != j; k = k1) { + k1 = mod(k + 1, m); + t = tangent(p0, p1, p2, p3, pp->curve.vertex[k], pp->curve.vertex[k1]); + if (t < -.5) { + return 1; + } + pt = bezier(t, p0, p1, p2, p3); + d = ddist(pp->curve.vertex[k], pp->curve.vertex[k1]); + if (d == 0.0) { /* this should never happen */ + return 1; + } + d1 = dpara(pp->curve.vertex[k], pp->curve.vertex[k1], pt) / d; + if (fabs(d1) > opttolerance) { + return 1; + } + if (iprod(pp->curve.vertex[k], pp->curve.vertex[k1], pt) < 0 || + iprod(pp->curve.vertex[k1], pp->curve.vertex[k], pt) < 0) { + return 1; + } + res->pen += sq(d1); } - if (d1 < d2) { - res->pen += sq(d1 - d2); + + /* check corners */ + for (k = i; k != j; k = k1) { + k1 = mod(k + 1, m); + t = tangent(p0, p1, p2, p3, pp->curve.c[k][2], pp->curve.c[k1][2]); + if (t < -.5) { + return 1; + } + pt = bezier(t, p0, p1, p2, p3); + d = ddist(pp->curve.c[k][2], pp->curve.c[k1][2]); + if (d == 0.0) { /* this should never happen */ + return 1; + } + d1 = dpara(pp->curve.c[k][2], pp->curve.c[k1][2], pt) / d; + d2 = dpara(pp->curve.c[k][2], pp->curve.c[k1][2], pp->curve.vertex[k1]) / d; + d2 *= 0.75 * pp->curve.alpha[k1]; + if (d2 < 0) { + d1 = -d1; + d2 = -d2; + } + if (d1 < d2 - opttolerance) { + return 1; + } + if (d1 < d2) { + res->pen += sq(d1 - d2); + } } - } - return 0; + return 0; } /* optimize the path p, replacing sequences of Bezier segments by a single segment when possible. Return 0 on success, 1 with errno set on failure. */ -static int opticurve(privpath_t *pp, double opttolerance) { - int m = pp->curve.n; - int *pt = NULL; /* pt[m+1] */ - double *pen = NULL; /* pen[m+1] */ - int *len = NULL; /* len[m+1] */ - opti_t *opt = NULL; /* opt[m+1] */ - int om; - int i,j,r; - opti_t o; - dpoint_t p0; - int i1; - double area; - double alpha; - double *s = NULL; - double *t = NULL; - - int *convc = NULL; /* conv[m]: pre-computed convexities */ - double *areac = NULL; /* cumarea[m+1]: cache for fast area computation */ - - SAFE_MALLOC(pt, m+1, int); - SAFE_MALLOC(pen, m+1, double); - SAFE_MALLOC(len, m+1, int); - SAFE_MALLOC(opt, m+1, opti_t); - SAFE_MALLOC(convc, m, int); - SAFE_MALLOC(areac, m+1, double); - - /* pre-calculate convexity: +1 = right turn, -1 = left turn, 0 = corner */ - for (i=0; icurve.tag[i] == POTRACE_CURVETO) { - convc[i] = sign(dpara(pp->curve.vertex[mod(i-1,m)], pp->curve.vertex[i], pp->curve.vertex[mod(i+1,m)])); - } else { - convc[i] = 0; +static int opticurve(privpath_t *pp, double opttolerance) +{ + int m = pp->curve.n; + int *pt = NULL; /* pt[m+1] */ + double *pen = NULL; /* pen[m+1] */ + int *len = NULL; /* len[m+1] */ + opti_t *opt = NULL; /* opt[m+1] */ + int om; + int i, j, r; + opti_t o; + dpoint_t p0; + int i1; + double area; + double alpha; + double *s = NULL; + double *t = NULL; + + int *convc = NULL; /* conv[m]: pre-computed convexities */ + double *areac = NULL; /* cumarea[m+1]: cache for fast area computation */ + + SAFE_MALLOC(pt, m + 1, int); + SAFE_MALLOC(pen, m + 1, double); + SAFE_MALLOC(len, m + 1, int); + SAFE_MALLOC(opt, m + 1, opti_t); + SAFE_MALLOC(convc, m, int); + SAFE_MALLOC(areac, m + 1, double); + + /* pre-calculate convexity: +1 = right turn, -1 = left turn, 0 = corner */ + for (i = 0; i < m; i++) { + if (pp->curve.tag[i] == POTRACE_CURVETO) { + convc[i] = + sign(dpara(pp->curve.vertex[mod(i - 1, m)], pp->curve.vertex[i], pp->curve.vertex[mod(i + 1, m)])); + } else { + convc[i] = 0; + } } - } - - /* pre-calculate areas */ - area = 0.0; - areac[0] = 0.0; - p0 = pp->curve.vertex[0]; - for (i=0; icurve.tag[i1] == POTRACE_CURVETO) { - alpha = pp->curve.alpha[i1]; - area += 0.3*alpha*(4-alpha)*dpara(pp->curve.c[i][2], pp->curve.vertex[i1], pp->curve.c[i1][2])/2; - area += dpara(p0, pp->curve.c[i][2], pp->curve.c[i1][2])/2; + + /* pre-calculate areas */ + area = 0.0; + areac[0] = 0.0; + p0 = pp->curve.vertex[0]; + for (i = 0; i < m; i++) { + i1 = mod(i + 1, m); + if (pp->curve.tag[i1] == POTRACE_CURVETO) { + alpha = pp->curve.alpha[i1]; + area += 0.3 * alpha * (4 - alpha) * dpara(pp->curve.c[i][2], pp->curve.vertex[i1], pp->curve.c[i1][2]) / 2; + area += dpara(p0, pp->curve.c[i][2], pp->curve.c[i1][2]) / 2; + } + areac[i + 1] = area; } - areac[i+1] = area; - } - - pt[0] = -1; - pen[0] = 0; - len[0] = 0; - - /* Fixme: we always start from a fixed point -- should find the best - curve cyclically */ - - for (j=1; j<=m; j++) { - /* calculate best path from 0 to j */ - pt[j] = j-1; - pen[j] = pen[j-1]; - len[j] = len[j-1]+1; - - for (i=j-2; i>=0; i--) { - r = opti_penalty(pp, i, mod(j,m), &o, opttolerance, convc, areac); - if (r) { - break; - } - if (len[j] > len[i]+1 || (len[j] == len[i]+1 && pen[j] > pen[i] + o.pen)) { - pt[j] = i; - pen[j] = pen[i] + o.pen; - len[j] = len[i] + 1; - opt[j] = o; - } + + pt[0] = -1; + pen[0] = 0; + len[0] = 0; + + /* Fixme: we always start from a fixed point -- should find the best + curve cyclically */ + + for (j = 1; j <= m; j++) { + /* calculate best path from 0 to j */ + pt[j] = j - 1; + pen[j] = pen[j - 1]; + len[j] = len[j - 1] + 1; + + for (i = j - 2; i >= 0; i--) { + r = opti_penalty(pp, i, mod(j, m), &o, opttolerance, convc, areac); + if (r) { + break; + } + if (len[j] > len[i] + 1 || (len[j] == len[i] + 1 && pen[j] > pen[i] + o.pen)) { + pt[j] = i; + pen[j] = pen[i] + o.pen; + len[j] = len[i] + 1; + opt[j] = o; + } + } } - } - om = len[m]; - r = privcurve_init(&pp->ocurve, om); - if (r) { - goto malloc_error; - } - SAFE_MALLOC(s, om, double); - SAFE_MALLOC(t, om, double); - - j = m; - for (i=om-1; i>=0; i--) { - if (pt[j]==j-1) { - pp->ocurve.tag[i] = pp->curve.tag[mod(j,m)]; - pp->ocurve.c[i][0] = pp->curve.c[mod(j,m)][0]; - pp->ocurve.c[i][1] = pp->curve.c[mod(j,m)][1]; - pp->ocurve.c[i][2] = pp->curve.c[mod(j,m)][2]; - pp->ocurve.vertex[i] = pp->curve.vertex[mod(j,m)]; - pp->ocurve.alpha[i] = pp->curve.alpha[mod(j,m)]; - pp->ocurve.alpha0[i] = pp->curve.alpha0[mod(j,m)]; - pp->ocurve.beta[i] = pp->curve.beta[mod(j,m)]; - s[i] = t[i] = 1.0; - } else { - pp->ocurve.tag[i] = POTRACE_CURVETO; - pp->ocurve.c[i][0] = opt[j].c[0]; - pp->ocurve.c[i][1] = opt[j].c[1]; - pp->ocurve.c[i][2] = pp->curve.c[mod(j,m)][2]; - pp->ocurve.vertex[i] = interval(opt[j].s, pp->curve.c[mod(j,m)][2], pp->curve.vertex[mod(j,m)]); - pp->ocurve.alpha[i] = opt[j].alpha; - pp->ocurve.alpha0[i] = opt[j].alpha; - s[i] = opt[j].s; - t[i] = opt[j].t; + om = len[m]; + r = privcurve_init(&pp->ocurve, om); + if (r) { + goto malloc_error; } - j = pt[j]; - } - - /* calculate beta parameters */ - for (i=0; iocurve.beta[i] = s[i] / (s[i] + t[i1]); - } - pp->ocurve.alphacurve = 1; - - free(pt); - free(pen); - free(len); - free(opt); - free(s); - free(t); - free(convc); - free(areac); - return 0; - - malloc_error: - free(pt); - free(pen); - free(len); - free(opt); - free(s); - free(t); - free(convc); - free(areac); - return 1; + SAFE_MALLOC(s, om, double); + SAFE_MALLOC(t, om, double); + + j = m; + for (i = om - 1; i >= 0; i--) { + if (pt[j] == j - 1) { + pp->ocurve.tag[i] = pp->curve.tag[mod(j, m)]; + pp->ocurve.c[i][0] = pp->curve.c[mod(j, m)][0]; + pp->ocurve.c[i][1] = pp->curve.c[mod(j, m)][1]; + pp->ocurve.c[i][2] = pp->curve.c[mod(j, m)][2]; + pp->ocurve.vertex[i] = pp->curve.vertex[mod(j, m)]; + pp->ocurve.alpha[i] = pp->curve.alpha[mod(j, m)]; + pp->ocurve.alpha0[i] = pp->curve.alpha0[mod(j, m)]; + pp->ocurve.beta[i] = pp->curve.beta[mod(j, m)]; + s[i] = t[i] = 1.0; + } else { + pp->ocurve.tag[i] = POTRACE_CURVETO; + pp->ocurve.c[i][0] = opt[j].c[0]; + pp->ocurve.c[i][1] = opt[j].c[1]; + pp->ocurve.c[i][2] = pp->curve.c[mod(j, m)][2]; + pp->ocurve.vertex[i] = interval(opt[j].s, pp->curve.c[mod(j, m)][2], pp->curve.vertex[mod(j, m)]); + pp->ocurve.alpha[i] = opt[j].alpha; + pp->ocurve.alpha0[i] = opt[j].alpha; + s[i] = opt[j].s; + t[i] = opt[j].t; + } + j = pt[j]; + } + + /* calculate beta parameters */ + for (i = 0; i < om; i++) { + i1 = mod(i + 1, om); + pp->ocurve.beta[i] = s[i] / (s[i] + t[i1]); + } + pp->ocurve.alphacurve = 1; + + free(pt); + free(pen); + free(len); + free(opt); + free(s); + free(t); + free(convc); + free(areac); + return 0; + +malloc_error: + free(pt); + free(pen); + free(len); + free(opt); + free(s); + free(t); + free(convc); + free(areac); + return 1; } /* ---------------------------------------------------------------------- */ -#define TRY(x) if (x) goto try_error +#define TRY(x) \ + if (x) \ + goto try_error /* return 0 on success, 1 on error with errno set. */ -int process_path(path_t *plist, const potrace_param_t *param, progress_t *progress) { - path_t *p; - double nn = 0, cn = 0; - - if (progress->callback) { - /* precompute task size for progress estimates */ - nn = 0; - list_forall (p, plist) { - nn += p->priv->len; - } - cn = 0; - } - - /* call downstream function with each path */ - list_forall (p, plist) { - TRY(calc_sums(p->priv)); - TRY(calc_lon(p->priv)); - TRY(bestpolygon(p->priv)); - TRY(adjust_vertices(p->priv)); - if (p->sign == '-') { /* reverse orientation of negative paths */ - reverse(&p->priv->curve); - } - smooth(&p->priv->curve, param->alphamax); - if (param->opticurve) { - TRY(opticurve(p->priv, param->opttolerance)); - p->priv->fcurve = &p->priv->ocurve; - } else { - p->priv->fcurve = &p->priv->curve; - } - privcurve_to_curve(p->priv->fcurve, &p->curve); +int process_path(path_t *plist, const potrace_param_t *param, progress_t *progress) +{ + path_t *p; + double nn = 0, cn = 0; if (progress->callback) { - cn += p->priv->len; - progress_update(cn/nn, progress); + /* precompute task size for progress estimates */ + nn = 0; + list_forall(p, plist) { nn += p->priv->len; } + cn = 0; + } + + /* call downstream function with each path */ + list_forall(p, plist) + { + TRY(calc_sums(p->priv)); + TRY(calc_lon(p->priv)); + TRY(bestpolygon(p->priv)); + TRY(adjust_vertices(p->priv)); + if (p->sign == '-') { /* reverse orientation of negative paths */ + reverse(&p->priv->curve); + } + smooth(&p->priv->curve, param->alphamax); + if (param->opticurve) { + TRY(opticurve(p->priv, param->opttolerance)); + p->priv->fcurve = &p->priv->ocurve; + } else { + p->priv->fcurve = &p->priv->curve; + } + privcurve_to_curve(p->priv->fcurve, &p->curve); + + if (progress->callback) { + cn += p->priv->len; + progress_update(cn / nn, progress); + } } - } - progress_update(1.0, progress); + progress_update(1.0, progress); - return 0; + return 0; - try_error: - return 1; +try_error: + return 1; } -- cgit v1.2.3 From a2dca92e95dabcdaa66310f489dba5178e3e1914 Mon Sep 17 00:00:00 2001 From: "Johan B. C. Engelen" Date: Sun, 12 Oct 2014 13:16:16 +0200 Subject: undo whitespace edit. it's an external dependency, and let's assume that the potential bug is not a bug (line 611, possibly uninitialized used of prev[i]) (bzr r13594) --- src/trace/potrace/trace.cpp | 2059 +++++++++++++++++++++---------------------- 1 file changed, 1019 insertions(+), 1040 deletions(-) (limited to 'src') diff --git a/src/trace/potrace/trace.cpp b/src/trace/potrace/trace.cpp index 2bdefb90a..f1e88a908 100644 --- a/src/trace/potrace/trace.cpp +++ b/src/trace/potrace/trace.cpp @@ -16,129 +16,124 @@ #include "trace.h" #include "progress.h" -#define INFTY 10000000 // it suffices that this is longer than any - // path; it need not be really infinite -#define COS179 -0.999847695156 // the cosine of 179 degrees +#define INFTY 10000000 /* it suffices that this is longer than any + path; it need not be really infinite */ +#define COS179 -0.999847695156 /* the cosine of 179 degrees */ /* ---------------------------------------------------------------------- */ #define SAFE_MALLOC(var, n, typ) \ - if ((var = (typ *)malloc((n)*sizeof(typ))) == NULL) goto malloc_error + if ((var = (typ *)malloc((n)*sizeof(typ))) == NULL) goto malloc_error /* ---------------------------------------------------------------------- */ /* auxiliary functions */ /* return a direction that is 90 degrees counterclockwise from p2-p0, but then restricted to one of the major wind directions (n, nw, w, etc) */ -static inline point_t dorth_infty(dpoint_t p0, dpoint_t p2) -{ - point_t r; - - r.y = sign(p2.x - p0.x); - r.x = -sign(p2.y - p0.y); +static inline point_t dorth_infty(dpoint_t p0, dpoint_t p2) { + point_t r; + + r.y = sign(p2.x-p0.x); + r.x = -sign(p2.y-p0.y); - return r; + return r; } /* return (p1-p0)x(p2-p0), the area of the parallelogram */ -static inline double dpara(dpoint_t p0, dpoint_t p1, dpoint_t p2) -{ - double x1, y1, x2, y2; +static inline double dpara(dpoint_t p0, dpoint_t p1, dpoint_t p2) { + double x1, y1, x2, y2; - x1 = p1.x - p0.x; - y1 = p1.y - p0.y; - x2 = p2.x - p0.x; - y2 = p2.y - p0.y; + x1 = p1.x-p0.x; + y1 = p1.y-p0.y; + x2 = p2.x-p0.x; + y2 = p2.y-p0.y; - return x1 * y2 - x2 * y1; + return x1*y2 - x2*y1; } /* ddenom/dpara have the property that the square of radius 1 centered at p1 intersects the line p0p2 iff |dpara(p0,p1,p2)| <= ddenom(p0,p2) */ -static inline double ddenom(dpoint_t p0, dpoint_t p2) -{ - point_t r = dorth_infty(p0, p2); +static inline double ddenom(dpoint_t p0, dpoint_t p2) { + point_t r = dorth_infty(p0, p2); - return r.y * (p2.x - p0.x) - r.x * (p2.y - p0.y); + return r.y*(p2.x-p0.x) - r.x*(p2.y-p0.y); } /* return 1 if a <= b < c < a, in a cyclic sense (mod n) */ -static inline int cyclic(int a, int b, int c) -{ - if (a <= c) { - return (a <= b && b < c); - } else { - return (a <= b || b < c); - } +static inline int cyclic(int a, int b, int c) { + if (a<=c) { + return (a<=b && blen; - sums_t *sums = pp->sums; - - double x, y, x2, xy, y2; - double k; - double a, b, c, lambda2, l; - int r = 0; /* rotations from i to j */ - - while (j >= n) { - j -= n; - r += 1; - } - while (i >= n) { - i -= n; - r -= 1; - } - while (j < 0) { - j += n; - r -= 1; - } - while (i < 0) { - i += n; - r += 1; - } - - x = sums[j + 1].x - sums[i].x + r * sums[n].x; - y = sums[j + 1].y - sums[i].y + r * sums[n].y; - x2 = sums[j + 1].x2 - sums[i].x2 + r * sums[n].x2; - xy = sums[j + 1].xy - sums[i].xy + r * sums[n].xy; - y2 = sums[j + 1].y2 - sums[i].y2 + r * sums[n].y2; - k = j + 1 - i + r * n; - - ctr->x = x / k; - ctr->y = y / k; - - a = (x2 - (double)x * x / k) / k; - b = (xy - (double)x * y / k) / k; - c = (y2 - (double)y * y / k) / k; - - lambda2 = (a + c + sqrt((a - c) * (a - c) + 4 * b * b)) / 2; /* larger e.value */ - - /* now find e.vector for lambda2 */ - a -= lambda2; - c -= lambda2; - - if (fabs(a) >= fabs(c)) { - l = sqrt(a * a + b * b); - if (l != 0) { - dir->x = -b / l; - dir->y = a / l; - } - } else { - l = sqrt(c * c + b * b); - if (l != 0) { - dir->x = -c / l; - dir->y = b / l; - } +static void pointslope(privpath_t *pp, int i, int j, dpoint_t *ctr, dpoint_t *dir) { + /* assume ilen; + sums_t *sums = pp->sums; + + double x, y, x2, xy, y2; + double k; + double a, b, c, lambda2, l; + int r=0; /* rotations from i to j */ + + while (j>=n) { + j-=n; + r+=1; + } + while (i>=n) { + i-=n; + r-=1; + } + while (j<0) { + j+=n; + r-=1; + } + while (i<0) { + i+=n; + r+=1; + } + + x = sums[j+1].x-sums[i].x+r*sums[n].x; + y = sums[j+1].y-sums[i].y+r*sums[n].y; + x2 = sums[j+1].x2-sums[i].x2+r*sums[n].x2; + xy = sums[j+1].xy-sums[i].xy+r*sums[n].xy; + y2 = sums[j+1].y2-sums[i].y2+r*sums[n].y2; + k = j+1-i+r*n; + + ctr->x = x/k; + ctr->y = y/k; + + a = (x2-(double)x*x/k)/k; + b = (xy-(double)x*y/k)/k; + c = (y2-(double)y*y/k)/k; + + lambda2 = (a+c+sqrt((a-c)*(a-c)+4*b*b))/2; /* larger e.value */ + + /* now find e.vector for lambda2 */ + a -= lambda2; + c -= lambda2; + + if (fabs(a) >= fabs(c)) { + l = sqrt(a*a+b*b); + if (l!=0) { + dir->x = -b/l; + dir->y = a/l; } - if (l == 0) { - dir->x = dir->y = 0; /* sometimes this can happen when k=4: - the two eigenvalues coincide */ + } else { + l = sqrt(c*c+b*b); + if (l!=0) { + dir->x = -c/l; + dir->y = b/l; } + } + if (l==0) { + dir->x = dir->y = 0; /* sometimes this can happen when k=4: + the two eigenvalues coincide */ + } } /* the type of (affine) quadratic forms, represented as symmetric 3x3 @@ -147,158 +142,155 @@ static void pointslope(privpath_t *pp, int i, int j, dpoint_t *ctr, dpoint_t *di typedef double quadform_t[3][3]; /* Apply quadratic form Q to vector w = (w.x,w.y) */ -static inline double quadform(quadform_t Q, dpoint_t w) -{ - double v[3]; - int i, j; - double sum; - - v[0] = w.x; - v[1] = w.y; - v[2] = 1; - sum = 0.0; - - for (i = 0; i < 3; i++) { - for (j = 0; j < 3; j++) { - sum += v[i] * Q[i][j] * v[j]; - } +static inline double quadform(quadform_t Q, dpoint_t w) { + double v[3]; + int i, j; + double sum; + + v[0] = w.x; + v[1] = w.y; + v[2] = 1; + sum = 0.0; + + for (i=0; i<3; i++) { + for (j=0; j<3; j++) { + sum += v[i] * Q[i][j] * v[j]; } - return sum; + } + return sum; } /* calculate p1 x p2 */ -static inline int xprod(point_t p1, point_t p2) { return p1.x * p2.y - p1.y * p2.x; } +static inline int xprod(point_t p1, point_t p2) { + return p1.x*p2.y - p1.y*p2.x; +} /* calculate (p1-p0)x(p3-p2) */ -static inline double cprod(dpoint_t p0, dpoint_t p1, dpoint_t p2, dpoint_t p3) -{ - double x1, y1, x2, y2; +static inline double cprod(dpoint_t p0, dpoint_t p1, dpoint_t p2, dpoint_t p3) { + double x1, y1, x2, y2; - x1 = p1.x - p0.x; - y1 = p1.y - p0.y; - x2 = p3.x - p2.x; - y2 = p3.y - p2.y; + x1 = p1.x - p0.x; + y1 = p1.y - p0.y; + x2 = p3.x - p2.x; + y2 = p3.y - p2.y; - return x1 * y2 - x2 * y1; + return x1*y2 - x2*y1; } /* calculate (p1-p0)*(p2-p0) */ -static inline double iprod(dpoint_t p0, dpoint_t p1, dpoint_t p2) -{ - double x1, y1, x2, y2; +static inline double iprod(dpoint_t p0, dpoint_t p1, dpoint_t p2) { + double x1, y1, x2, y2; - x1 = p1.x - p0.x; - y1 = p1.y - p0.y; - x2 = p2.x - p0.x; - y2 = p2.y - p0.y; + x1 = p1.x - p0.x; + y1 = p1.y - p0.y; + x2 = p2.x - p0.x; + y2 = p2.y - p0.y; - return x1 * x2 + y1 * y2; + return x1*x2 + y1*y2; } /* calculate (p1-p0)*(p3-p2) */ -static inline double iprod1(dpoint_t p0, dpoint_t p1, dpoint_t p2, dpoint_t p3) -{ - double x1, y1, x2, y2; +static inline double iprod1(dpoint_t p0, dpoint_t p1, dpoint_t p2, dpoint_t p3) { + double x1, y1, x2, y2; - x1 = p1.x - p0.x; - y1 = p1.y - p0.y; - x2 = p3.x - p2.x; - y2 = p3.y - p2.y; + x1 = p1.x - p0.x; + y1 = p1.y - p0.y; + x2 = p3.x - p2.x; + y2 = p3.y - p2.y; - return x1 * x2 + y1 * y2; + return x1*x2 + y1*y2; } /* calculate distance between two points */ -static inline double ddist(dpoint_t p, dpoint_t q) { return sqrt(sq(p.x - q.x) + sq(p.y - q.y)); } +static inline double ddist(dpoint_t p, dpoint_t q) { + return sqrt(sq(p.x-q.x)+sq(p.y-q.y)); +} /* calculate point of a bezier curve */ -static inline dpoint_t bezier(double t, dpoint_t p0, dpoint_t p1, dpoint_t p2, dpoint_t p3) -{ - double s = 1 - t; - dpoint_t res; +static inline dpoint_t bezier(double t, dpoint_t p0, dpoint_t p1, dpoint_t p2, dpoint_t p3) { + double s = 1-t; + dpoint_t res; - /* Note: a good optimizing compiler (such as gcc-3) reduces the - following to 16 multiplications, using common subexpression - elimination. */ + /* Note: a good optimizing compiler (such as gcc-3) reduces the + following to 16 multiplications, using common subexpression + elimination. */ - res.x = s * s * s * p0.x + 3 * (s * s * t) * p1.x + 3 * (t * t * s) * p2.x + t * t * t * p3.x; - res.y = s * s * s * p0.y + 3 * (s * s * t) * p1.y + 3 * (t * t * s) * p2.y + t * t * t * p3.y; + res.x = s*s*s*p0.x + 3*(s*s*t)*p1.x + 3*(t*t*s)*p2.x + t*t*t*p3.x; + res.y = s*s*s*p0.y + 3*(s*s*t)*p1.y + 3*(t*t*s)*p2.y + t*t*t*p3.y; - return res; + return res; } /* calculate the point t in [0..1] on the (convex) bezier curve (p0,p1,p2,p3) which is tangent to q1-q0. Return -1.0 if there is no solution in [0..1]. */ -static double tangent(dpoint_t p0, dpoint_t p1, dpoint_t p2, dpoint_t p3, dpoint_t q0, dpoint_t q1) -{ - double A, B, C; /* (1-t)^2 A + 2(1-t)t B + t^2 C = 0 */ - double a, b, c; /* a t^2 + b t + c = 0 */ - double d, s, r1, r2; - - A = cprod(p0, p1, q0, q1); - B = cprod(p1, p2, q0, q1); - C = cprod(p2, p3, q0, q1); - - a = A - 2 * B + C; - b = -2 * A + 2 * B; - c = A; - - d = b * b - 4 * a * c; - - if (a == 0 || d < 0) { - return -1.0; - } - - s = sqrt(d); - - r1 = (-b + s) / (2 * a); - r2 = (-b - s) / (2 * a); - - if (r1 >= 0 && r1 <= 1) { - return r1; - } else if (r2 >= 0 && r2 <= 1) { - return r2; - } else { - return -1.0; - } +static double tangent(dpoint_t p0, dpoint_t p1, dpoint_t p2, dpoint_t p3, dpoint_t q0, dpoint_t q1) { + double A, B, C; /* (1-t)^2 A + 2(1-t)t B + t^2 C = 0 */ + double a, b, c; /* a t^2 + b t + c = 0 */ + double d, s, r1, r2; + + A = cprod(p0, p1, q0, q1); + B = cprod(p1, p2, q0, q1); + C = cprod(p2, p3, q0, q1); + + a = A - 2*B + C; + b = -2*A + 2*B; + c = A; + + d = b*b - 4*a*c; + + if (a==0 || d<0) { + return -1.0; + } + + s = sqrt(d); + + r1 = (-b + s) / (2 * a); + r2 = (-b - s) / (2 * a); + + if (r1 >= 0 && r1 <= 1) { + return r1; + } else if (r2 >= 0 && r2 <= 1) { + return r2; + } else { + return -1.0; + } } /* ---------------------------------------------------------------------- */ /* Preparation: fill in the sum* fields of a path (used for later rapid summing). Return 0 on success, 1 with errno set on failure. */ -static int calc_sums(privpath_t *pp) -{ - int i, x, y; - int n = pp->len; - - SAFE_MALLOC(pp->sums, pp->len + 1, sums_t); - - /* origin */ - pp->x0 = pp->pt[0].x; - pp->y0 = pp->pt[0].y; - - /* preparatory computation for later fast summing */ - pp->sums[0].x2 = pp->sums[0].xy = pp->sums[0].y2 = pp->sums[0].x = pp->sums[0].y = 0; - for (i = 0; i < n; i++) { - x = pp->pt[i].x - pp->x0; - y = pp->pt[i].y - pp->y0; - pp->sums[i + 1].x = pp->sums[i].x + x; - pp->sums[i + 1].y = pp->sums[i].y + y; - pp->sums[i + 1].x2 = pp->sums[i].x2 + x * x; - pp->sums[i + 1].xy = pp->sums[i].xy + x * y; - pp->sums[i + 1].y2 = pp->sums[i].y2 + y * y; - } - return 0; - -malloc_error: - return 1; +static int calc_sums(privpath_t *pp) { + int i, x, y; + int n = pp->len; + + SAFE_MALLOC(pp->sums, pp->len+1, sums_t); + + /* origin */ + pp->x0 = pp->pt[0].x; + pp->y0 = pp->pt[0].y; + + /* preparatory computation for later fast summing */ + pp->sums[0].x2 = pp->sums[0].xy = pp->sums[0].y2 = pp->sums[0].x = pp->sums[0].y = 0; + for (i=0; ipt[i].x - pp->x0; + y = pp->pt[i].y - pp->y0; + pp->sums[i+1].x = pp->sums[i].x + x; + pp->sums[i+1].y = pp->sums[i].y + y; + pp->sums[i+1].x2 = pp->sums[i].x2 + x*x; + pp->sums[i+1].xy = pp->sums[i].xy + x*y; + pp->sums[i+1].y2 = pp->sums[i].y2 + y*y; + } + return 0; + + malloc_error: + return 1; } /* ---------------------------------------------------------------------- */ /* Stage 1: determine the straight subpaths (Sec. 2.2.1). Fill in the - "lon" component of a path object (based on pt/len). For each i, + "lon" component of a path object (based on pt/len). For each i, lon[i] is the furthest index such that a straight line can be drawn from i to lon[i]. Return 1 on error with errno set, else 0. */ @@ -326,208 +318,206 @@ malloc_error: substantial. */ /* returns 0 on success, 1 on error with errno set */ -static int calc_lon(privpath_t *pp) -{ - point_t *pt = pp->pt; - int n = pp->len; - int i, j, k, k1; - int ct[4], dir; - point_t constraint[2]; - point_t cur; - point_t off; - int *pivk = NULL; /* pivk[n] */ - int *nc = NULL; /* nc[n]: next corner */ - point_t dk; /* direction of k-k1 */ - int a, b, c, d; - - SAFE_MALLOC(pivk, n, int); - SAFE_MALLOC(nc, n, int); - - /* initialize the nc data structure. Point from each point to the - furthest future point to which it is connected by a vertical or - horizontal segment. We take advantage of the fact that there is - always a direction change at 0 (due to the path decomposition - algorithm). But even if this were not so, there is no harm, as - in practice, correctness does not depend on the word "furthest" - above. */ - k = 0; - for (i = n - 1; i >= 0; i--) { - if (pt[i].x != pt[k].x && pt[i].y != pt[k].y) { - k = i + 1; /* necessarily ipt; + int n = pp->len; + int i, j, k, k1; + int ct[4], dir; + point_t constraint[2]; + point_t cur; + point_t off; + int *pivk = NULL; /* pivk[n] */ + int *nc = NULL; /* nc[n]: next corner */ + point_t dk; /* direction of k-k1 */ + int a, b, c, d; + + SAFE_MALLOC(pivk, n, int); + SAFE_MALLOC(nc, n, int); + + /* initialize the nc data structure. Point from each point to the + furthest future point to which it is connected by a vertical or + horizontal segment. We take advantage of the fact that there is + always a direction change at 0 (due to the path decomposition + algorithm). But even if this were not so, there is no harm, as + in practice, correctness does not depend on the word "furthest" + above. */ + k = 0; + for (i=n-1; i>=0; i--) { + if (pt[i].x != pt[k].x && pt[i].y != pt[k].y) { + k = i+1; /* necessarily ilon, n, int); - - /* determine pivot points: for each i, let pivk[i] be the furthest k - such that all j with i= 0; i--) { - ct[0] = ct[1] = ct[2] = ct[3] = 0; - - /* keep track of "directions" that have occurred */ - dir = (3 + 3 * (pt[mod(i + 1, n)].x - pt[i].x) + (pt[mod(i + 1, n)].y - pt[i].y)) / 2; - ct[dir]++; - - constraint[0].x = 0; - constraint[0].y = 0; - constraint[1].x = 0; - constraint[1].y = 0; - - /* find the next k such that no straight line from i to k */ - k = nc[i]; - k1 = i; - while (1) { - - dir = (3 + 3 * sign(pt[k].x - pt[k1].x) + sign(pt[k].y - pt[k1].y)) / 2; - ct[dir]++; - - /* if all four "directions" have occurred, cut this path */ - if (ct[0] && ct[1] && ct[2] && ct[3]) { - pivk[i] = k1; - goto foundk; - } - - cur.x = pt[k].x - pt[i].x; - cur.y = pt[k].y - pt[i].y; - - /* see if current constraint is violated */ - if (xprod(constraint[0], cur) < 0 || xprod(constraint[1], cur) > 0) { - goto constraint_viol; - } - - /* else, update constraint */ - if (abs(cur.x) <= 1 && abs(cur.y) <= 1) { - /* no constraint */ - } else { - off.x = cur.x + ((cur.y >= 0 && (cur.y > 0 || cur.x < 0)) ? 1 : -1); - off.y = cur.y + ((cur.x <= 0 && (cur.x < 0 || cur.y < 0)) ? 1 : -1); - if (xprod(constraint[0], off) >= 0) { - constraint[0] = off; - } - off.x = cur.x + ((cur.y <= 0 && (cur.y < 0 || cur.x < 0)) ? 1 : -1); - off.y = cur.y + ((cur.x >= 0 && (cur.x > 0 || cur.y < 0)) ? 1 : -1); - if (xprod(constraint[1], off) <= 0) { - constraint[1] = off; - } - } - k1 = k; - k = nc[k1]; - if (!cyclic(k, i, k1)) { - break; - } - } - constraint_viol: - /* k1 was the last "corner" satisfying the current constraint, and - k is the first one violating it. We now need to find the last - point along k1..k which satisfied the constraint. */ - dk.x = sign(pt[k].x - pt[k1].x); - dk.y = sign(pt[k].y - pt[k1].y); - cur.x = pt[k1].x - pt[i].x; - cur.y = pt[k1].y - pt[i].y; - /* find largest integer j such that xprod(constraint[0], cur+j*dk) - >= 0 and xprod(constraint[1], cur+j*dk) <= 0. Use bilinearity - of xprod. */ - a = xprod(constraint[0], cur); - b = xprod(constraint[0], dk); - c = xprod(constraint[1], cur); - d = xprod(constraint[1], dk); - /* find largest integer j such that a+j*b>=0 and c+j*d<=0. This - can be solved with integer arithmetic. */ - j = INFTY; - if (b < 0) { - j = floordiv(a, -b); - } - if (d > 0) { - j = min(j, floordiv(-c, d)); - } - pivk[i] = mod(k1 + j, n); - foundk: - ; - } /* for i */ - - /* clean up: for each i, let lon[i] be the largest k such that for - all i' with i<=i'lon[n - 1] = j; - for (i = n - 2; i >= 0; i--) { - if (cyclic(i + 1, pivk[i], j)) { - j = pivk[i]; - } - pp->lon[i] = j; + nc[i] = k; + } + + SAFE_MALLOC(pp->lon, n, int); + + /* determine pivot points: for each i, let pivk[i] be the furthest k + such that all j with i=0; i--) { + ct[0] = ct[1] = ct[2] = ct[3] = 0; + + /* keep track of "directions" that have occurred */ + dir = (3+3*(pt[mod(i+1,n)].x-pt[i].x)+(pt[mod(i+1,n)].y-pt[i].y))/2; + ct[dir]++; + + constraint[0].x = 0; + constraint[0].y = 0; + constraint[1].x = 0; + constraint[1].y = 0; + + /* find the next k such that no straight line from i to k */ + k = nc[i]; + k1 = i; + while (1) { + + dir = (3+3*sign(pt[k].x-pt[k1].x)+sign(pt[k].y-pt[k1].y))/2; + ct[dir]++; + + /* if all four "directions" have occurred, cut this path */ + if (ct[0] && ct[1] && ct[2] && ct[3]) { + pivk[i] = k1; + goto foundk; + } + + cur.x = pt[k].x - pt[i].x; + cur.y = pt[k].y - pt[i].y; + + /* see if current constraint is violated */ + if (xprod(constraint[0], cur) < 0 || xprod(constraint[1], cur) > 0) { + goto constraint_viol; + } + + /* else, update constraint */ + if (abs(cur.x) <= 1 && abs(cur.y) <= 1) { + /* no constraint */ + } else { + off.x = cur.x + ((cur.y>=0 && (cur.y>0 || cur.x<0)) ? 1 : -1); + off.y = cur.y + ((cur.x<=0 && (cur.x<0 || cur.y<0)) ? 1 : -1); + if (xprod(constraint[0], off) >= 0) { + constraint[0] = off; + } + off.x = cur.x + ((cur.y<=0 && (cur.y<0 || cur.x<0)) ? 1 : -1); + off.y = cur.y + ((cur.x>=0 && (cur.x>0 || cur.y<0)) ? 1 : -1); + if (xprod(constraint[1], off) <= 0) { + constraint[1] = off; + } + } + k1 = k; + k = nc[k1]; + if (!cyclic(k,i,k1)) { + break; + } } - - for (i = n - 1; cyclic(mod(i + 1, n), j, pp->lon[i]); i--) { - pp->lon[i] = j; + constraint_viol: + /* k1 was the last "corner" satisfying the current constraint, and + k is the first one violating it. We now need to find the last + point along k1..k which satisfied the constraint. */ + dk.x = sign(pt[k].x-pt[k1].x); + dk.y = sign(pt[k].y-pt[k1].y); + cur.x = pt[k1].x - pt[i].x; + cur.y = pt[k1].y - pt[i].y; + /* find largest integer j such that xprod(constraint[0], cur+j*dk) + >= 0 and xprod(constraint[1], cur+j*dk) <= 0. Use bilinearity + of xprod. */ + a = xprod(constraint[0], cur); + b = xprod(constraint[0], dk); + c = xprod(constraint[1], cur); + d = xprod(constraint[1], dk); + /* find largest integer j such that a+j*b>=0 and c+j*d<=0. This + can be solved with integer arithmetic. */ + j = INFTY; + if (b<0) { + j = floordiv(a,-b); + } + if (d>0) { + j = min(j, floordiv(-c,d)); } + pivk[i] = mod(k1+j,n); + foundk: + ; + } /* for i */ + + /* clean up: for each i, let lon[i] be the largest k such that for + all i' with i<=i'lon[n-1]=j; + for (i=n-2; i>=0; i--) { + if (cyclic(i+1,pivk[i],j)) { + j=pivk[i]; + } + pp->lon[i]=j; + } - free(pivk); - free(nc); - return 0; + for (i=n-1; cyclic(mod(i+1,n),j,pp->lon[i]); i--) { + pp->lon[i] = j; + } -malloc_error: - free(pivk); - free(nc); - return 1; + free(pivk); + free(nc); + return 0; + + malloc_error: + free(pivk); + free(nc); + return 1; } /* ---------------------------------------------------------------------- */ -/* Stage 2: calculate the optimal polygon (Sec. 2.2.2-2.2.4). */ +/* Stage 2: calculate the optimal polygon (Sec. 2.2.2-2.2.4). */ /* Auxiliary function: calculate the penalty of an edge from i to j in the given path. This needs the "lon" and "sum*" data. */ -static double penalty3(privpath_t *pp, int i, int j) -{ - int n = pp->len; - point_t *pt = pp->pt; - sums_t *sums = pp->sums; - - /* assume 0<=i= n) { - j -= n; - r = 1; - } - - /* critical inner loop: the "if" gives a 4.6 percent speedup */ - if (r == 0) { - x = sums[j + 1].x - sums[i].x; - y = sums[j + 1].y - sums[i].y; - x2 = sums[j + 1].x2 - sums[i].x2; - xy = sums[j + 1].xy - sums[i].xy; - y2 = sums[j + 1].y2 - sums[i].y2; - k = j + 1 - i; - } else { - x = sums[j + 1].x - sums[i].x + sums[n].x; - y = sums[j + 1].y - sums[i].y + sums[n].y; - x2 = sums[j + 1].x2 - sums[i].x2 + sums[n].x2; - xy = sums[j + 1].xy - sums[i].xy + sums[n].xy; - y2 = sums[j + 1].y2 - sums[i].y2 + sums[n].y2; - k = j + 1 - i + n; - } - - px = (pt[i].x + pt[j].x) / 2.0 - pt[0].x; - py = (pt[i].y + pt[j].y) / 2.0 - pt[0].y; - ey = (pt[j].x - pt[i].x); - ex = -(pt[j].y - pt[i].y); - - a = ((x2 - 2 * x * px) / k + px * px); - b = ((xy - x * py - y * px) / k + px * py); - c = ((y2 - 2 * y * py) / k + py * py); - - s = ex * ex * a + 2 * ex * ey * b + ey * ey * c; - - return sqrt(s); +static double penalty3(privpath_t *pp, int i, int j) { + int n = pp->len; + point_t *pt = pp->pt; + sums_t *sums = pp->sums; + + /* assume 0<=i=n) { + j -= n; + r = 1; + } + + /* critical inner loop: the "if" gives a 4.6 percent speedup */ + if (r == 0) { + x = sums[j+1].x - sums[i].x; + y = sums[j+1].y - sums[i].y; + x2 = sums[j+1].x2 - sums[i].x2; + xy = sums[j+1].xy - sums[i].xy; + y2 = sums[j+1].y2 - sums[i].y2; + k = j+1 - i; + } else { + x = sums[j+1].x - sums[i].x + sums[n].x; + y = sums[j+1].y - sums[i].y + sums[n].y; + x2 = sums[j+1].x2 - sums[i].x2 + sums[n].x2; + xy = sums[j+1].xy - sums[i].xy + sums[n].xy; + y2 = sums[j+1].y2 - sums[i].y2 + sums[n].y2; + k = j+1 - i + n; + } + + px = (pt[i].x + pt[j].x) / 2.0 - pt[0].x; + py = (pt[i].y + pt[j].y) / 2.0 - pt[0].y; + ey = (pt[j].x - pt[i].x); + ex = -(pt[j].y - pt[i].y); + + a = ((x2 - 2*x*px) / k + px*px); + b = ((xy - x*py - y*px) / k + px*py); + c = ((y2 - 2*y*py) / k + py*py); + + s = ex*ex*a + 2*ex*ey*b + ey*ey*c; + + return sqrt(s); } /* find the optimal polygon. Fill in the m and po components. Return 1 @@ -535,109 +525,109 @@ static double penalty3(privpath_t *pp, int i, int j) is in the polygon. Fixme: implement cyclic version. */ static int bestpolygon(privpath_t *pp) { - int i, j, m, k; - int n = pp->len; - double *pen = NULL; /* pen[n+1]: penalty vector */ - int *prev = NULL; /* prev[n+1]: best path pointer vector */ - int *clip0 = NULL; /* clip0[n]: longest segment pointer, non-cyclic */ - int *clip1 = NULL; /* clip1[n+1]: backwards segment pointer, non-cyclic */ - int *seg0 = NULL; /* seg0[m+1]: forward segment bounds, m<=n */ - int *seg1 = NULL; /* seg1[m+1]: backward segment bounds, m<=n */ - double thispen; - double best; - int c; - - SAFE_MALLOC(pen, n + 1, double); - SAFE_MALLOC(prev, n + 1, int); - SAFE_MALLOC(clip0, n, int); - SAFE_MALLOC(clip1, n + 1, int); - SAFE_MALLOC(seg0, n + 1, int); - SAFE_MALLOC(seg1, n + 1, int); - - /* calculate clipped paths */ - for (i = 0; i < n; i++) { - c = mod(pp->lon[mod(i - 1, n)] - 1, n); - if (c == i) { - c = mod(i + 1, n); - } - if (c < i) { - clip0[i] = n; - } else { - clip0[i] = c; - } - } - - /* calculate backwards path clipping, non-cyclic. j <= clip0[i] iff - clip1[j] <= i, for i,j=0..n. */ - j = 1; - for (i = 0; i < n; i++) { - while (j <= clip0[i]) { - clip1[j] = i; - j++; - } - } - - /* calculate seg0[j] = longest path from 0 with j segments */ - i = 0; - for (j = 0; i < n; j++) { - seg0[j] = i; - i = clip0[i]; + int i, j, m, k; + int n = pp->len; + double *pen = NULL; /* pen[n+1]: penalty vector */ + int *prev = NULL; /* prev[n+1]: best path pointer vector */ + int *clip0 = NULL; /* clip0[n]: longest segment pointer, non-cyclic */ + int *clip1 = NULL; /* clip1[n+1]: backwards segment pointer, non-cyclic */ + int *seg0 = NULL; /* seg0[m+1]: forward segment bounds, m<=n */ + int *seg1 = NULL; /* seg1[m+1]: backward segment bounds, m<=n */ + double thispen; + double best; + int c; + + SAFE_MALLOC(pen, n+1, double); + SAFE_MALLOC(prev, n+1, int); + SAFE_MALLOC(clip0, n, int); + SAFE_MALLOC(clip1, n+1, int); + SAFE_MALLOC(seg0, n+1, int); + SAFE_MALLOC(seg1, n+1, int); + + /* calculate clipped paths */ + for (i=0; ilon[mod(i-1,n)]-1,n); + if (c == i) { + c = mod(i+1,n); } - seg0[j] = n; - m = j; - - /* calculate seg1[j] = longest path to n with m-j segments */ - i = n; - for (j = m; j > 0; j--) { - seg1[j] = i; - i = clip1[i]; + if (c < i) { + clip0[i] = n; + } else { + clip0[i] = c; } - seg1[0] = 0; - - /* now find the shortest path with m segments, based on penalty3 */ - /* note: the outer 2 loops jointly have at most n iterations, thus - the worst-case behavior here is quadratic. In practice, it is - close to linear since the inner loop tends to be short. */ - pen[0] = 0; - for (j = 1; j <= m; j++) { - for (i = seg1[j]; i <= seg0[j]; i++) { - best = -1; - for (k = seg0[j - 1]; k >= clip1[i]; k--) { - thispen = penalty3(pp, k, i) + pen[k]; - if (best < 0 || thispen < best) { - prev[i] = k; - best = thispen; - } - } - pen[i] = best; - } + } + + /* calculate backwards path clipping, non-cyclic. j <= clip0[i] iff + clip1[j] <= i, for i,j=0..n. */ + j = 1; + for (i=0; im = m; - SAFE_MALLOC(pp->po, m, int); - - /* read off shortest path */ - for (i = n, j = m - 1; i > 0; j--) { - i = prev[i]; - pp->po[j] = i; + } + + /* calculate seg0[j] = longest path from 0 with j segments */ + i = 0; + for (j=0; i0; j--) { + seg1[j] = i; + i = clip1[i]; + } + seg1[0] = 0; + + /* now find the shortest path with m segments, based on penalty3 */ + /* note: the outer 2 loops jointly have at most n iterations, thus + the worst-case behavior here is quadratic. In practice, it is + close to linear since the inner loop tends to be short. */ + pen[0]=0; + for (j=1; j<=m; j++) { + for (i=seg1[j]; i<=seg0[j]; i++) { + best = -1; + for (k=seg0[j-1]; k>=clip1[i]; k--) { + thispen = penalty3(pp, k, i) + pen[k]; + if (best < 0 || thispen < best) { + prev[i] = k; + best = thispen; + } + } + pen[i] = best; } - - free(pen); - free(prev); - free(clip0); - free(clip1); - free(seg0); - free(seg1); - return 0; - -malloc_error: - free(pen); - free(prev); - free(clip0); - free(clip1); - free(seg0); - free(seg1); - return 1; + } + + pp->m = m; + SAFE_MALLOC(pp->po, m, int); + + /* read off shortest path */ + for (i=n, j=m-1; i>0; j--) { + i = prev[i]; + pp->po[j] = i; + } + + free(pen); + free(prev); + free(clip0); + free(clip1); + free(seg0); + free(seg1); + return 0; + + malloc_error: + free(pen); + free(prev); + free(clip0); + free(clip1); + free(seg0); + free(seg1); + return 1; } /* ---------------------------------------------------------------------- */ @@ -648,269 +638,266 @@ malloc_error: if it lies outside. Return 1 with errno set on error; 0 on success. */ -static int adjust_vertices(privpath_t *pp) -{ - int m = pp->m; - int *po = pp->po; - int n = pp->len; - point_t *pt = pp->pt; - int x0 = pp->x0; - int y0 = pp->y0; - - dpoint_t *ctr = NULL; /* ctr[m] */ - dpoint_t *dir = NULL; /* dir[m] */ - quadform_t *q = NULL; /* q[m] */ - double v[3]; - double d; - int i, j, k, l; - dpoint_t s; - int r; - - SAFE_MALLOC(ctr, m, dpoint_t); - SAFE_MALLOC(dir, m, dpoint_t); - SAFE_MALLOC(q, m, quadform_t); - - r = privcurve_init(&pp->curve, m); - if (r) { - goto malloc_error; - } - - /* calculate "optimal" point-slope representation for each line - segment */ - for (i = 0; i < m; i++) { - j = po[mod(i + 1, m)]; - j = mod(j - po[i], n) + po[i]; - pointslope(pp, po[i], j, &ctr[i], &dir[i]); +static int adjust_vertices(privpath_t *pp) { + int m = pp->m; + int *po = pp->po; + int n = pp->len; + point_t *pt = pp->pt; + int x0 = pp->x0; + int y0 = pp->y0; + + dpoint_t *ctr = NULL; /* ctr[m] */ + dpoint_t *dir = NULL; /* dir[m] */ + quadform_t *q = NULL; /* q[m] */ + double v[3]; + double d; + int i, j, k, l; + dpoint_t s; + int r; + + SAFE_MALLOC(ctr, m, dpoint_t); + SAFE_MALLOC(dir, m, dpoint_t); + SAFE_MALLOC(q, m, quadform_t); + + r = privcurve_init(&pp->curve, m); + if (r) { + goto malloc_error; + } + + /* calculate "optimal" point-slope representation for each line + segment */ + for (i=0; i Q[1][1]) { - v[0] = -Q[0][1]; - v[1] = Q[0][0]; - } else if (Q[1][1]) { - v[0] = -Q[1][1]; - v[1] = Q[1][0]; - } else { - v[0] = 1; - v[1] = 0; - } - d = sq(v[0]) + sq(v[1]); - v[2] = -v[1] * s.y - v[0] * s.x; - for (l = 0; l < 3; l++) { - for (k = 0; k < 3; k++) { - Q[l][k] += v[l] * v[k] / d; - } - } - } - dx = fabs(w.x - s.x); - dy = fabs(w.y - s.y); - if (dx <= .5 && dy <= .5) { - pp->curve.vertex[i].x = w.x + x0; - pp->curve.vertex[i].y = w.y + y0; - continue; - } - - /* the minimum was not in the unit square; now minimize quadratic - on boundary of square */ - min = quadform(Q, s); - xmin = s.x; - ymin = s.y; - - if (Q[0][0] == 0.0) { - goto fixx; - } - for (z = 0; z < 2; z++) { /* value of the y-coordinate */ - w.y = s.y - 0.5 + z; - w.x = -(Q[0][1] * w.y + Q[0][2]) / Q[0][0]; - dx = fabs(w.x - s.x); - cand = quadform(Q, w); - if (dx <= .5 && cand < min) { - min = cand; - xmin = w.x; - ymin = w.y; - } - } - fixx: - if (Q[1][1] == 0.0) { - goto corners; - } - for (z = 0; z < 2; z++) { /* value of the x-coordinate */ - w.x = s.x - 0.5 + z; - w.y = -(Q[1][0] * w.x + Q[1][2]) / Q[1][1]; - dy = fabs(w.y - s.y); - cand = quadform(Q, w); - if (dy <= .5 && cand < min) { - min = cand; - xmin = w.x; - ymin = w.y; - } - } - corners: - /* check four corners */ - for (l = 0; l < 2; l++) { - for (k = 0; k < 2; k++) { - w.x = s.x - 0.5 + l; - w.y = s.y - 0.5 + k; - cand = quadform(Q, w); - if (cand < min) { - min = cand; - xmin = w.x; - ymin = w.y; - } - } - } - - pp->curve.vertex[i].x = xmin + x0; - pp->curve.vertex[i].y = ymin + y0; - continue; + + det = Q[0][0]*Q[1][1] - Q[0][1]*Q[1][0]; + if (det != 0.0) { + w.x = (-Q[0][2]*Q[1][1] + Q[1][2]*Q[0][1]) / det; + w.y = ( Q[0][2]*Q[1][0] - Q[1][2]*Q[0][0]) / det; + break; + } + + /* matrix is singular - lines are parallel. Add another, + orthogonal axis, through the center of the unit square */ + if (Q[0][0]>Q[1][1]) { + v[0] = -Q[0][1]; + v[1] = Q[0][0]; + } else if (Q[1][1]) { + v[0] = -Q[1][1]; + v[1] = Q[1][0]; + } else { + v[0] = 1; + v[1] = 0; + } + d = sq(v[0]) + sq(v[1]); + v[2] = - v[1] * s.y - v[0] * s.x; + for (l=0; l<3; l++) { + for (k=0; k<3; k++) { + Q[l][k] += v[l] * v[k] / d; + } + } + } + dx = fabs(w.x-s.x); + dy = fabs(w.y-s.y); + if (dx <= .5 && dy <= .5) { + pp->curve.vertex[i].x = w.x+x0; + pp->curve.vertex[i].y = w.y+y0; + continue; } - free(ctr); - free(dir); - free(q); - return 0; + /* the minimum was not in the unit square; now minimize quadratic + on boundary of square */ + min = quadform(Q, s); + xmin = s.x; + ymin = s.y; -malloc_error: - free(ctr); - free(dir); - free(q); - return 1; + if (Q[0][0] == 0.0) { + goto fixx; + } + for (z=0; z<2; z++) { /* value of the y-coordinate */ + w.y = s.y-0.5+z; + w.x = - (Q[0][1] * w.y + Q[0][2]) / Q[0][0]; + dx = fabs(w.x-s.x); + cand = quadform(Q, w); + if (dx <= .5 && cand < min) { + min = cand; + xmin = w.x; + ymin = w.y; + } + } + fixx: + if (Q[1][1] == 0.0) { + goto corners; + } + for (z=0; z<2; z++) { /* value of the x-coordinate */ + w.x = s.x-0.5+z; + w.y = - (Q[1][0] * w.x + Q[1][2]) / Q[1][1]; + dy = fabs(w.y-s.y); + cand = quadform(Q, w); + if (dy <= .5 && cand < min) { + min = cand; + xmin = w.x; + ymin = w.y; + } + } + corners: + /* check four corners */ + for (l=0; l<2; l++) { + for (k=0; k<2; k++) { + w.x = s.x-0.5+l; + w.y = s.y-0.5+k; + cand = quadform(Q, w); + if (cand < min) { + min = cand; + xmin = w.x; + ymin = w.y; + } + } + } + + pp->curve.vertex[i].x = xmin + x0; + pp->curve.vertex[i].y = ymin + y0; + continue; + } + + free(ctr); + free(dir); + free(q); + return 0; + + malloc_error: + free(ctr); + free(dir); + free(q); + return 1; } /* ---------------------------------------------------------------------- */ /* Stage 4: smoothing and corner analysis (Sec. 2.3.3) */ /* reverse orientation of a path */ -static void reverse(privcurve_t *curve) -{ - int m = curve->n; - int i, j; - dpoint_t tmp; - - for (i = 0, j = m - 1; i < j; i++, j--) { - tmp = curve->vertex[i]; - curve->vertex[i] = curve->vertex[j]; - curve->vertex[j] = tmp; - } +static void reverse(privcurve_t *curve) { + int m = curve->n; + int i, j; + dpoint_t tmp; + + for (i=0, j=m-1; ivertex[i]; + curve->vertex[i] = curve->vertex[j]; + curve->vertex[j] = tmp; + } } /* Always succeeds */ -static void smooth(privcurve_t *curve, double alphamax) -{ - int m = curve->n; - - int i, j, k; - double dd, denom, alpha; - dpoint_t p2, p3, p4; - - /* examine each vertex and find its best fit */ - for (i = 0; i < m; i++) { - j = mod(i + 1, m); - k = mod(i + 2, m); - p4 = interval(1 / 2.0, curve->vertex[k], curve->vertex[j]); - - denom = ddenom(curve->vertex[i], curve->vertex[k]); - if (denom != 0.0) { - dd = dpara(curve->vertex[i], curve->vertex[j], curve->vertex[k]) / denom; - dd = fabs(dd); - alpha = dd > 1 ? (1 - 1.0 / dd) : 0; - alpha = alpha / 0.75; - } else { - alpha = 4 / 3.0; - } - curve->alpha0[j] = alpha; /* remember "original" value of alpha */ - - if (alpha > alphamax) { /* pointed corner */ - curve->tag[j] = POTRACE_CORNER; - curve->c[j][1] = curve->vertex[j]; - curve->c[j][2] = p4; - } else { - if (alpha < 0.55) { - alpha = 0.55; - } else if (alpha > 1) { - alpha = 1; - } - p2 = interval(.5 + .5 * alpha, curve->vertex[i], curve->vertex[j]); - p3 = interval(.5 + .5 * alpha, curve->vertex[k], curve->vertex[j]); - curve->tag[j] = POTRACE_CURVETO; - curve->c[j][0] = p2; - curve->c[j][1] = p3; - curve->c[j][2] = p4; - } - curve->alpha[j] = alpha; /* store the "cropped" value of alpha */ - curve->beta[j] = 0.5; +static void smooth(privcurve_t *curve, double alphamax) { + int m = curve->n; + + int i, j, k; + double dd, denom, alpha; + dpoint_t p2, p3, p4; + + /* examine each vertex and find its best fit */ + for (i=0; ivertex[k], curve->vertex[j]); + + denom = ddenom(curve->vertex[i], curve->vertex[k]); + if (denom != 0.0) { + dd = dpara(curve->vertex[i], curve->vertex[j], curve->vertex[k]) / denom; + dd = fabs(dd); + alpha = dd>1 ? (1 - 1.0/dd) : 0; + alpha = alpha / 0.75; + } else { + alpha = 4/3.0; + } + curve->alpha0[j] = alpha; /* remember "original" value of alpha */ + + if (alpha > alphamax) { /* pointed corner */ + curve->tag[j] = POTRACE_CORNER; + curve->c[j][1] = curve->vertex[j]; + curve->c[j][2] = p4; + } else { + if (alpha < 0.55) { + alpha = 0.55; + } else if (alpha > 1) { + alpha = 1; + } + p2 = interval(.5+.5*alpha, curve->vertex[i], curve->vertex[j]); + p3 = interval(.5+.5*alpha, curve->vertex[k], curve->vertex[j]); + curve->tag[j] = POTRACE_CURVETO; + curve->c[j][0] = p2; + curve->c[j][1] = p3; + curve->c[j][2] = p4; } - curve->alphacurve = 1; + curve->alpha[j] = alpha; /* store the "cropped" value of alpha */ + curve->beta[j] = 0.5; + } + curve->alphacurve = 1; - return; + return; } /* ---------------------------------------------------------------------- */ @@ -918,349 +905,341 @@ static void smooth(privcurve_t *curve, double alphamax) /* a private type for the result of opti_penalty */ struct opti_s { - double pen; /* penalty */ - dpoint_t c[2]; /* curve parameters */ - double t, s; /* curve parameters */ - double alpha; /* curve parameter */ + double pen; /* penalty */ + dpoint_t c[2]; /* curve parameters */ + double t, s; /* curve parameters */ + double alpha; /* curve parameter */ }; typedef struct opti_s opti_t; /* calculate best fit from i+.5 to j+.5. Assume icurve.n; - int k, k1, k2, conv, i1; - double area, alpha, d, d1, d2; - dpoint_t p0, p1, p2, p3, pt; - double A, R, A1, A2, A3, A4; - double s, t; +static int opti_penalty(privpath_t *pp, int i, int j, opti_t *res, double opttolerance, int *convc, double *areac) { + int m = pp->curve.n; + int k, k1, k2, conv, i1; + double area, alpha, d, d1, d2; + dpoint_t p0, p1, p2, p3, pt; + double A, R, A1, A2, A3, A4; + double s, t; - /* check convexity, corner-freeness, and maximum bend < 179 degrees */ + /* check convexity, corner-freeness, and maximum bend < 179 degrees */ - if (i == j) { /* sanity - a full loop can never be an opticurve */ - return 1; - } + if (i==j) { /* sanity - a full loop can never be an opticurve */ + return 1; + } - k = i; - i1 = mod(i + 1, m); - k1 = mod(k + 1, m); - conv = convc[k1]; - if (conv == 0) { - return 1; + k = i; + i1 = mod(i+1, m); + k1 = mod(k+1, m); + conv = convc[k1]; + if (conv == 0) { + return 1; + } + d = ddist(pp->curve.vertex[i], pp->curve.vertex[i1]); + for (k=k1; k!=j; k=k1) { + k1 = mod(k+1, m); + k2 = mod(k+2, m); + if (convc[k1] != conv) { + return 1; } - d = ddist(pp->curve.vertex[i], pp->curve.vertex[i1]); - for (k = k1; k != j; k = k1) { - k1 = mod(k + 1, m); - k2 = mod(k + 2, m); - if (convc[k1] != conv) { - return 1; - } - if (sign(cprod(pp->curve.vertex[i], pp->curve.vertex[i1], pp->curve.vertex[k1], pp->curve.vertex[k2])) != - conv) { - return 1; - } - if (iprod1(pp->curve.vertex[i], pp->curve.vertex[i1], pp->curve.vertex[k1], pp->curve.vertex[k2]) < - d * ddist(pp->curve.vertex[k1], pp->curve.vertex[k2]) * COS179) { - return 1; - } + if (sign(cprod(pp->curve.vertex[i], pp->curve.vertex[i1], pp->curve.vertex[k1], pp->curve.vertex[k2])) != conv) { + return 1; } - - /* the curve we're working in: */ - p0 = pp->curve.c[mod(i, m)][2]; - p1 = pp->curve.vertex[mod(i + 1, m)]; - p2 = pp->curve.vertex[mod(j, m)]; - p3 = pp->curve.c[mod(j, m)][2]; - - /* determine its area */ - area = areac[j] - areac[i]; - area -= dpara(pp->curve.vertex[0], pp->curve.c[i][2], pp->curve.c[j][2]) / 2; - if (i >= j) { - area += areac[m]; + if (iprod1(pp->curve.vertex[i], pp->curve.vertex[i1], pp->curve.vertex[k1], pp->curve.vertex[k2]) < d * ddist(pp->curve.vertex[k1], pp->curve.vertex[k2]) * COS179) { + return 1; } + } + + /* the curve we're working in: */ + p0 = pp->curve.c[mod(i,m)][2]; + p1 = pp->curve.vertex[mod(i+1,m)]; + p2 = pp->curve.vertex[mod(j,m)]; + p3 = pp->curve.c[mod(j,m)][2]; + + /* determine its area */ + area = areac[j] - areac[i]; + area -= dpara(pp->curve.vertex[0], pp->curve.c[i][2], pp->curve.c[j][2])/2; + if (i>=j) { + area += areac[m]; + } + + /* find intersection o of p0p1 and p2p3. Let t,s such that o = + interval(t,p0,p1) = interval(s,p3,p2). Let A be the area of the + triangle (p0,o,p3). */ + + A1 = dpara(p0, p1, p2); + A2 = dpara(p0, p1, p3); + A3 = dpara(p0, p2, p3); + /* A4 = dpara(p1, p2, p3); */ + A4 = A1+A3-A2; + + if (A2 == A1) { /* this should never happen */ + return 1; + } - /* find intersection o of p0p1 and p2p3. Let t,s such that o = - interval(t,p0,p1) = interval(s,p3,p2). Let A be the area of the - triangle (p0,o,p3). */ + t = A3/(A3-A4); + s = A2/(A2-A1); + A = A2 * t / 2.0; + + if (A == 0.0) { /* this should never happen */ + return 1; + } - A1 = dpara(p0, p1, p2); - A2 = dpara(p0, p1, p3); - A3 = dpara(p0, p2, p3); - /* A4 = dpara(p1, p2, p3); */ - A4 = A1 + A3 - A2; + R = area / A; /* relative area */ + alpha = 2 - sqrt(4 - R / 0.3); /* overall alpha for p0-o-p3 curve */ - if (A2 == A1) { /* this should never happen */ - return 1; - } + res->c[0] = interval(t * alpha, p0, p1); + res->c[1] = interval(s * alpha, p3, p2); + res->alpha = alpha; + res->t = t; + res->s = s; - t = A3 / (A3 - A4); - s = A2 / (A2 - A1); - A = A2 * t / 2.0; + p1 = res->c[0]; + p2 = res->c[1]; /* the proposed curve is now (p0,p1,p2,p3) */ - if (A == 0.0) { /* this should never happen */ - return 1; - } + res->pen = 0; - R = area / A; /* relative area */ - alpha = 2 - sqrt(4 - R / 0.3); /* overall alpha for p0-o-p3 curve */ - - res->c[0] = interval(t * alpha, p0, p1); - res->c[1] = interval(s * alpha, p3, p2); - res->alpha = alpha; - res->t = t; - res->s = s; - - p1 = res->c[0]; - p2 = res->c[1]; /* the proposed curve is now (p0,p1,p2,p3) */ - - res->pen = 0; - - /* calculate penalty */ - /* check tangency with edges */ - for (k = mod(i + 1, m); k != j; k = k1) { - k1 = mod(k + 1, m); - t = tangent(p0, p1, p2, p3, pp->curve.vertex[k], pp->curve.vertex[k1]); - if (t < -.5) { - return 1; - } - pt = bezier(t, p0, p1, p2, p3); - d = ddist(pp->curve.vertex[k], pp->curve.vertex[k1]); - if (d == 0.0) { /* this should never happen */ - return 1; - } - d1 = dpara(pp->curve.vertex[k], pp->curve.vertex[k1], pt) / d; - if (fabs(d1) > opttolerance) { - return 1; - } - if (iprod(pp->curve.vertex[k], pp->curve.vertex[k1], pt) < 0 || - iprod(pp->curve.vertex[k1], pp->curve.vertex[k], pt) < 0) { - return 1; - } - res->pen += sq(d1); + /* calculate penalty */ + /* check tangency with edges */ + for (k=mod(i+1,m); k!=j; k=k1) { + k1 = mod(k+1,m); + t = tangent(p0, p1, p2, p3, pp->curve.vertex[k], pp->curve.vertex[k1]); + if (t<-.5) { + return 1; } - - /* check corners */ - for (k = i; k != j; k = k1) { - k1 = mod(k + 1, m); - t = tangent(p0, p1, p2, p3, pp->curve.c[k][2], pp->curve.c[k1][2]); - if (t < -.5) { - return 1; - } - pt = bezier(t, p0, p1, p2, p3); - d = ddist(pp->curve.c[k][2], pp->curve.c[k1][2]); - if (d == 0.0) { /* this should never happen */ - return 1; - } - d1 = dpara(pp->curve.c[k][2], pp->curve.c[k1][2], pt) / d; - d2 = dpara(pp->curve.c[k][2], pp->curve.c[k1][2], pp->curve.vertex[k1]) / d; - d2 *= 0.75 * pp->curve.alpha[k1]; - if (d2 < 0) { - d1 = -d1; - d2 = -d2; - } - if (d1 < d2 - opttolerance) { - return 1; - } - if (d1 < d2) { - res->pen += sq(d1 - d2); - } + pt = bezier(t, p0, p1, p2, p3); + d = ddist(pp->curve.vertex[k], pp->curve.vertex[k1]); + if (d == 0.0) { /* this should never happen */ + return 1; } + d1 = dpara(pp->curve.vertex[k], pp->curve.vertex[k1], pt) / d; + if (fabs(d1) > opttolerance) { + return 1; + } + if (iprod(pp->curve.vertex[k], pp->curve.vertex[k1], pt) < 0 || iprod(pp->curve.vertex[k1], pp->curve.vertex[k], pt) < 0) { + return 1; + } + res->pen += sq(d1); + } + + /* check corners */ + for (k=i; k!=j; k=k1) { + k1 = mod(k+1,m); + t = tangent(p0, p1, p2, p3, pp->curve.c[k][2], pp->curve.c[k1][2]); + if (t<-.5) { + return 1; + } + pt = bezier(t, p0, p1, p2, p3); + d = ddist(pp->curve.c[k][2], pp->curve.c[k1][2]); + if (d == 0.0) { /* this should never happen */ + return 1; + } + d1 = dpara(pp->curve.c[k][2], pp->curve.c[k1][2], pt) / d; + d2 = dpara(pp->curve.c[k][2], pp->curve.c[k1][2], pp->curve.vertex[k1]) / d; + d2 *= 0.75 * pp->curve.alpha[k1]; + if (d2 < 0) { + d1 = -d1; + d2 = -d2; + } + if (d1 < d2 - opttolerance) { + return 1; + } + if (d1 < d2) { + res->pen += sq(d1 - d2); + } + } - return 0; + return 0; } /* optimize the path p, replacing sequences of Bezier segments by a single segment when possible. Return 0 on success, 1 with errno set on failure. */ -static int opticurve(privpath_t *pp, double opttolerance) -{ - int m = pp->curve.n; - int *pt = NULL; /* pt[m+1] */ - double *pen = NULL; /* pen[m+1] */ - int *len = NULL; /* len[m+1] */ - opti_t *opt = NULL; /* opt[m+1] */ - int om; - int i, j, r; - opti_t o; - dpoint_t p0; - int i1; - double area; - double alpha; - double *s = NULL; - double *t = NULL; - - int *convc = NULL; /* conv[m]: pre-computed convexities */ - double *areac = NULL; /* cumarea[m+1]: cache for fast area computation */ - - SAFE_MALLOC(pt, m + 1, int); - SAFE_MALLOC(pen, m + 1, double); - SAFE_MALLOC(len, m + 1, int); - SAFE_MALLOC(opt, m + 1, opti_t); - SAFE_MALLOC(convc, m, int); - SAFE_MALLOC(areac, m + 1, double); - - /* pre-calculate convexity: +1 = right turn, -1 = left turn, 0 = corner */ - for (i = 0; i < m; i++) { - if (pp->curve.tag[i] == POTRACE_CURVETO) { - convc[i] = - sign(dpara(pp->curve.vertex[mod(i - 1, m)], pp->curve.vertex[i], pp->curve.vertex[mod(i + 1, m)])); - } else { - convc[i] = 0; - } - } - - /* pre-calculate areas */ - area = 0.0; - areac[0] = 0.0; - p0 = pp->curve.vertex[0]; - for (i = 0; i < m; i++) { - i1 = mod(i + 1, m); - if (pp->curve.tag[i1] == POTRACE_CURVETO) { - alpha = pp->curve.alpha[i1]; - area += 0.3 * alpha * (4 - alpha) * dpara(pp->curve.c[i][2], pp->curve.vertex[i1], pp->curve.c[i1][2]) / 2; - area += dpara(p0, pp->curve.c[i][2], pp->curve.c[i1][2]) / 2; - } - areac[i + 1] = area; - } - - pt[0] = -1; - pen[0] = 0; - len[0] = 0; - - /* Fixme: we always start from a fixed point -- should find the best - curve cyclically */ - - for (j = 1; j <= m; j++) { - /* calculate best path from 0 to j */ - pt[j] = j - 1; - pen[j] = pen[j - 1]; - len[j] = len[j - 1] + 1; - - for (i = j - 2; i >= 0; i--) { - r = opti_penalty(pp, i, mod(j, m), &o, opttolerance, convc, areac); - if (r) { - break; - } - if (len[j] > len[i] + 1 || (len[j] == len[i] + 1 && pen[j] > pen[i] + o.pen)) { - pt[j] = i; - pen[j] = pen[i] + o.pen; - len[j] = len[i] + 1; - opt[j] = o; - } - } +static int opticurve(privpath_t *pp, double opttolerance) { + int m = pp->curve.n; + int *pt = NULL; /* pt[m+1] */ + double *pen = NULL; /* pen[m+1] */ + int *len = NULL; /* len[m+1] */ + opti_t *opt = NULL; /* opt[m+1] */ + int om; + int i,j,r; + opti_t o; + dpoint_t p0; + int i1; + double area; + double alpha; + double *s = NULL; + double *t = NULL; + + int *convc = NULL; /* conv[m]: pre-computed convexities */ + double *areac = NULL; /* cumarea[m+1]: cache for fast area computation */ + + SAFE_MALLOC(pt, m+1, int); + SAFE_MALLOC(pen, m+1, double); + SAFE_MALLOC(len, m+1, int); + SAFE_MALLOC(opt, m+1, opti_t); + SAFE_MALLOC(convc, m, int); + SAFE_MALLOC(areac, m+1, double); + + /* pre-calculate convexity: +1 = right turn, -1 = left turn, 0 = corner */ + for (i=0; icurve.tag[i] == POTRACE_CURVETO) { + convc[i] = sign(dpara(pp->curve.vertex[mod(i-1,m)], pp->curve.vertex[i], pp->curve.vertex[mod(i+1,m)])); + } else { + convc[i] = 0; } - om = len[m]; - r = privcurve_init(&pp->ocurve, om); - if (r) { - goto malloc_error; + } + + /* pre-calculate areas */ + area = 0.0; + areac[0] = 0.0; + p0 = pp->curve.vertex[0]; + for (i=0; icurve.tag[i1] == POTRACE_CURVETO) { + alpha = pp->curve.alpha[i1]; + area += 0.3*alpha*(4-alpha)*dpara(pp->curve.c[i][2], pp->curve.vertex[i1], pp->curve.c[i1][2])/2; + area += dpara(p0, pp->curve.c[i][2], pp->curve.c[i1][2])/2; } - SAFE_MALLOC(s, om, double); - SAFE_MALLOC(t, om, double); - - j = m; - for (i = om - 1; i >= 0; i--) { - if (pt[j] == j - 1) { - pp->ocurve.tag[i] = pp->curve.tag[mod(j, m)]; - pp->ocurve.c[i][0] = pp->curve.c[mod(j, m)][0]; - pp->ocurve.c[i][1] = pp->curve.c[mod(j, m)][1]; - pp->ocurve.c[i][2] = pp->curve.c[mod(j, m)][2]; - pp->ocurve.vertex[i] = pp->curve.vertex[mod(j, m)]; - pp->ocurve.alpha[i] = pp->curve.alpha[mod(j, m)]; - pp->ocurve.alpha0[i] = pp->curve.alpha0[mod(j, m)]; - pp->ocurve.beta[i] = pp->curve.beta[mod(j, m)]; - s[i] = t[i] = 1.0; - } else { - pp->ocurve.tag[i] = POTRACE_CURVETO; - pp->ocurve.c[i][0] = opt[j].c[0]; - pp->ocurve.c[i][1] = opt[j].c[1]; - pp->ocurve.c[i][2] = pp->curve.c[mod(j, m)][2]; - pp->ocurve.vertex[i] = interval(opt[j].s, pp->curve.c[mod(j, m)][2], pp->curve.vertex[mod(j, m)]); - pp->ocurve.alpha[i] = opt[j].alpha; - pp->ocurve.alpha0[i] = opt[j].alpha; - s[i] = opt[j].s; - t[i] = opt[j].t; - } - j = pt[j]; + areac[i+1] = area; + } + + pt[0] = -1; + pen[0] = 0; + len[0] = 0; + + /* Fixme: we always start from a fixed point -- should find the best + curve cyclically */ + + for (j=1; j<=m; j++) { + /* calculate best path from 0 to j */ + pt[j] = j-1; + pen[j] = pen[j-1]; + len[j] = len[j-1]+1; + + for (i=j-2; i>=0; i--) { + r = opti_penalty(pp, i, mod(j,m), &o, opttolerance, convc, areac); + if (r) { + break; + } + if (len[j] > len[i]+1 || (len[j] == len[i]+1 && pen[j] > pen[i] + o.pen)) { + pt[j] = i; + pen[j] = pen[i] + o.pen; + len[j] = len[i] + 1; + opt[j] = o; + } } - - /* calculate beta parameters */ - for (i = 0; i < om; i++) { - i1 = mod(i + 1, om); - pp->ocurve.beta[i] = s[i] / (s[i] + t[i1]); + } + om = len[m]; + r = privcurve_init(&pp->ocurve, om); + if (r) { + goto malloc_error; + } + SAFE_MALLOC(s, om, double); + SAFE_MALLOC(t, om, double); + + j = m; + for (i=om-1; i>=0; i--) { + if (pt[j]==j-1) { + pp->ocurve.tag[i] = pp->curve.tag[mod(j,m)]; + pp->ocurve.c[i][0] = pp->curve.c[mod(j,m)][0]; + pp->ocurve.c[i][1] = pp->curve.c[mod(j,m)][1]; + pp->ocurve.c[i][2] = pp->curve.c[mod(j,m)][2]; + pp->ocurve.vertex[i] = pp->curve.vertex[mod(j,m)]; + pp->ocurve.alpha[i] = pp->curve.alpha[mod(j,m)]; + pp->ocurve.alpha0[i] = pp->curve.alpha0[mod(j,m)]; + pp->ocurve.beta[i] = pp->curve.beta[mod(j,m)]; + s[i] = t[i] = 1.0; + } else { + pp->ocurve.tag[i] = POTRACE_CURVETO; + pp->ocurve.c[i][0] = opt[j].c[0]; + pp->ocurve.c[i][1] = opt[j].c[1]; + pp->ocurve.c[i][2] = pp->curve.c[mod(j,m)][2]; + pp->ocurve.vertex[i] = interval(opt[j].s, pp->curve.c[mod(j,m)][2], pp->curve.vertex[mod(j,m)]); + pp->ocurve.alpha[i] = opt[j].alpha; + pp->ocurve.alpha0[i] = opt[j].alpha; + s[i] = opt[j].s; + t[i] = opt[j].t; } - pp->ocurve.alphacurve = 1; - - free(pt); - free(pen); - free(len); - free(opt); - free(s); - free(t); - free(convc); - free(areac); - return 0; - -malloc_error: - free(pt); - free(pen); - free(len); - free(opt); - free(s); - free(t); - free(convc); - free(areac); - return 1; + j = pt[j]; + } + + /* calculate beta parameters */ + for (i=0; iocurve.beta[i] = s[i] / (s[i] + t[i1]); + } + pp->ocurve.alphacurve = 1; + + free(pt); + free(pen); + free(len); + free(opt); + free(s); + free(t); + free(convc); + free(areac); + return 0; + + malloc_error: + free(pt); + free(pen); + free(len); + free(opt); + free(s); + free(t); + free(convc); + free(areac); + return 1; } /* ---------------------------------------------------------------------- */ -#define TRY(x) \ - if (x) \ - goto try_error +#define TRY(x) if (x) goto try_error /* return 0 on success, 1 on error with errno set. */ -int process_path(path_t *plist, const potrace_param_t *param, progress_t *progress) -{ - path_t *p; - double nn = 0, cn = 0; - - if (progress->callback) { - /* precompute task size for progress estimates */ - nn = 0; - list_forall(p, plist) { nn += p->priv->len; } - cn = 0; +int process_path(path_t *plist, const potrace_param_t *param, progress_t *progress) { + path_t *p; + double nn = 0, cn = 0; + + if (progress->callback) { + /* precompute task size for progress estimates */ + nn = 0; + list_forall (p, plist) { + nn += p->priv->len; + } + cn = 0; + } + + /* call downstream function with each path */ + list_forall (p, plist) { + TRY(calc_sums(p->priv)); + TRY(calc_lon(p->priv)); + TRY(bestpolygon(p->priv)); + TRY(adjust_vertices(p->priv)); + if (p->sign == '-') { /* reverse orientation of negative paths */ + reverse(&p->priv->curve); } + smooth(&p->priv->curve, param->alphamax); + if (param->opticurve) { + TRY(opticurve(p->priv, param->opttolerance)); + p->priv->fcurve = &p->priv->ocurve; + } else { + p->priv->fcurve = &p->priv->curve; + } + privcurve_to_curve(p->priv->fcurve, &p->curve); - /* call downstream function with each path */ - list_forall(p, plist) - { - TRY(calc_sums(p->priv)); - TRY(calc_lon(p->priv)); - TRY(bestpolygon(p->priv)); - TRY(adjust_vertices(p->priv)); - if (p->sign == '-') { /* reverse orientation of negative paths */ - reverse(&p->priv->curve); - } - smooth(&p->priv->curve, param->alphamax); - if (param->opticurve) { - TRY(opticurve(p->priv, param->opttolerance)); - p->priv->fcurve = &p->priv->ocurve; - } else { - p->priv->fcurve = &p->priv->curve; - } - privcurve_to_curve(p->priv->fcurve, &p->curve); - - if (progress->callback) { - cn += p->priv->len; - progress_update(cn / nn, progress); - } + if (progress->callback) { + cn += p->priv->len; + progress_update(cn/nn, progress); } + } - progress_update(1.0, progress); + progress_update(1.0, progress); - return 0; + return 0; -try_error: - return 1; + try_error: + return 1; } -- cgit v1.2.3 From 61da6fdd1620ae1d4d47e7bcc8c0a5671fcfe122 Mon Sep 17 00:00:00 2001 From: Tavmjong Bah Date: Sun, 12 Oct 2014 14:11:21 +0200 Subject: Fix "Access to field 'image_out' results in a dereference of a null pointer" (bzr r13595) --- src/filters/blend.cpp | 26 ++++++++++++++++---------- 1 file changed, 16 insertions(+), 10 deletions(-) (limited to 'src') diff --git a/src/filters/blend.cpp b/src/filters/blend.cpp index fd5a9871e..ca1d5bf96 100644 --- a/src/filters/blend.cpp +++ b/src/filters/blend.cpp @@ -183,6 +183,8 @@ void SPFeBlend::update(SPCtx *ctx, guint flags) { /* Unlike normal in, in2 is required attribute. Make sure, we can call * it by some name. */ + /* This may not be true.... see issue at + * http://www.w3.org/TR/filter-effects/#feBlendElement (but it doesn't hurt). */ if (this->in2 == Inkscape::Filters::NR_FILTER_SLOT_NOT_SET || this->in2 == Inkscape::Filters::NR_FILTER_UNNAMED_SLOT) { @@ -206,26 +208,30 @@ Inkscape::XML::Node* SPFeBlend::write(Inkscape::XML::Document *doc, Inkscape::XM repr = doc->createElement("svg:feBlend"); } - gchar const *out_name = sp_filter_name_for_image(parent, this->in2); + gchar const *in2_name = sp_filter_name_for_image(parent, this->in2); - if (out_name) { - repr->setAttribute("in2", out_name); - } else { + if( !in2_name ) { + + // This code is very similar to sp_filter_primtive_name_previous_out() SPObject *i = parent->children; + // Find previous filter primitive while (i && i->next != this) { i = i->next; } - SPFilterPrimitive *i_prim = SP_FILTER_PRIMITIVE(i); - out_name = sp_filter_name_for_image(parent, i_prim->image_out); - repr->setAttribute("in2", out_name); - - if (!out_name) { - g_warning("Unable to set in2 for feBlend"); + if( i ) { + SPFilterPrimitive *i_prim = SP_FILTER_PRIMITIVE(i); + in2_name = sp_filter_name_for_image(parent, i_prim->image_out); } } + if (in2_name) { + repr->setAttribute("in2", in2_name); + } else { + g_warning("Unable to set in2 for feBlend"); + } + char const *mode; switch(this->blend_mode) { case Inkscape::Filters::BLEND_NORMAL: -- cgit v1.2.3 From 17cb730370ab87bdb323c1e6972d09fe75b59f5a Mon Sep 17 00:00:00 2001 From: Tavmjong Bah Date: Sun, 12 Oct 2014 14:16:05 +0200 Subject: Fix "Access to field 'image_out' results in a dereference of a null pointer" (bzr r13596) --- src/filters/composite.cpp | 26 ++++++++++++++++---------- 1 file changed, 16 insertions(+), 10 deletions(-) (limited to 'src') diff --git a/src/filters/composite.cpp b/src/filters/composite.cpp index 257292f12..e600b6d20 100644 --- a/src/filters/composite.cpp +++ b/src/filters/composite.cpp @@ -205,6 +205,8 @@ void SPFeComposite::update(SPCtx *ctx, guint flags) { /* Unlike normal in, in2 is required attribute. Make sure, we can call * it by some name. */ + /* This may not be true.... see issue at + * http://www.w3.org/TR/filter-effects/#feBlendElement (but it doesn't hurt). */ if (this->in2 == Inkscape::Filters::NR_FILTER_SLOT_NOT_SET || this->in2 == Inkscape::Filters::NR_FILTER_UNNAMED_SLOT) { @@ -228,26 +230,30 @@ Inkscape::XML::Node* SPFeComposite::write(Inkscape::XML::Document *doc, Inkscape repr = doc->createElement("svg:feComposite"); } - gchar const *out_name = sp_filter_name_for_image(parent, this->in2); + gchar const *in2_name = sp_filter_name_for_image(parent, this->in2); - if (out_name) { - repr->setAttribute("in2", out_name); - } else { + if( !in2_name ) { + + // This code is very similar to sp_filter_primitive_name_previous_out() SPObject *i = parent->children; + // Find previous filter primitive while (i && i->next != this) { i = i->next; } - SPFilterPrimitive *i_prim = SP_FILTER_PRIMITIVE(i); - out_name = sp_filter_name_for_image(parent, i_prim->image_out); - repr->setAttribute("in2", out_name); - - if (!out_name) { - g_warning("Unable to set in2 for feComposite"); + if( i ) { + SPFilterPrimitive *i_prim = SP_FILTER_PRIMITIVE(i); + in2_name = sp_filter_name_for_image(parent, i_prim->image_out); } } + if (in2_name) { + repr->setAttribute("in2", in2_name); + } else { + g_warning("Unable to set in2 for feComposite"); + } + char const *comp_op; switch (this->composite_operator) { -- cgit v1.2.3 From be65f03b5b02746d3d21b5e905663d3264226321 Mon Sep 17 00:00:00 2001 From: Tavmjong Bah Date: Sun, 12 Oct 2014 14:19:49 +0200 Subject: Fix "Access to field 'image_out' results in a dereference of a null pointer" (bzr r13597) --- src/filters/displacementmap.cpp | 25 +++++++++++++++---------- 1 file changed, 15 insertions(+), 10 deletions(-) (limited to 'src') diff --git a/src/filters/displacementmap.cpp b/src/filters/displacementmap.cpp index 473e4913f..7dedfc031 100644 --- a/src/filters/displacementmap.cpp +++ b/src/filters/displacementmap.cpp @@ -202,25 +202,30 @@ Inkscape::XML::Node* SPFeDisplacementMap::write(Inkscape::XML::Document *doc, In repr = doc->createElement("svg:feDisplacementMap"); } - gchar const *out_name = sp_filter_name_for_image(parent, this->in2); - if (out_name) { - repr->setAttribute("in2", out_name); - } else { + gchar const *in2_name = sp_filter_name_for_image(parent, this->in2); + + if( !in2_name ) { + + // This code is very similar to sp_filter_primtive_name_previous_out() SPObject *i = parent->children; + // Find previous filter primitive while (i && i->next != this) { i = i->next; } - SPFilterPrimitive *i_prim = SP_FILTER_PRIMITIVE(i); - out_name = sp_filter_name_for_image(parent, i_prim->image_out); - repr->setAttribute("in2", out_name); - - if (!out_name) { - g_warning("Unable to set in2 for feDisplacementMap"); + if( i ) { + SPFilterPrimitive *i_prim = SP_FILTER_PRIMITIVE(i); + in2_name = sp_filter_name_for_image(parent, i_prim->image_out); } } + if (in2_name) { + repr->setAttribute("in2", in2_name); + } else { + g_warning("Unable to set in2 for feDisplacementMap"); + } + sp_repr_set_svg_double(repr, "scale", this->scale); repr->setAttribute("xChannelSelector", get_channelselector_name(this->xChannelSelector)); -- cgit v1.2.3 From 7fdadb51f279180e7613a7ef703cda1818fe007d Mon Sep 17 00:00:00 2001 From: Tavmjong Bah Date: Sun, 12 Oct 2014 14:27:21 +0200 Subject: Fix "Value stored to 'child' is never read" (bzr r13598) --- src/xml/repr-io.cpp | 1 - 1 file changed, 1 deletion(-) (limited to 'src') diff --git a/src/xml/repr-io.cpp b/src/xml/repr-io.cpp index 0319bb5e3..a4146f215 100644 --- a/src/xml/repr-io.cpp +++ b/src/xml/repr-io.cpp @@ -623,7 +623,6 @@ static Node *sp_repr_svg_read_node (Document *xml_doc, xmlNodePtr node, const gc repr->setContent(reinterpret_cast(node->content)); } - child = node->xmlChildrenNode; for (child = node->xmlChildrenNode; child != NULL; child = child->next) { Node *crepr = sp_repr_svg_read_node (xml_doc, child, default_ns, prefix_map); if (crepr) { -- cgit v1.2.3 From e85f0be00539e62796c7160d2a18f20f5d264b7a Mon Sep 17 00:00:00 2001 From: "Johan B. C. Engelen" Date: Sun, 12 Oct 2014 15:50:17 +0200 Subject: Extension enumeration/dropdownbox parameter: fix potential NULL-deref crash on xml == nullptr. (bzr r13599) --- src/extension/param/enum.cpp | 55 ++++++++++++++++++++++---------------------- 1 file changed, 28 insertions(+), 27 deletions(-) (limited to 'src') diff --git a/src/extension/param/enum.cpp b/src/extension/param/enum.cpp index bb50c06e1..74b2a75ad 100644 --- a/src/extension/param/enum.cpp +++ b/src/extension/param/enum.cpp @@ -41,8 +41,8 @@ namespace Extension { class enumentry { public: enumentry (Glib::ustring &val, Glib::ustring &text) : - value(val), - guitext(text) + value(val), + guitext(text) {} Glib::ustring value; @@ -50,16 +50,19 @@ public: }; -ParamComboBox::ParamComboBox (const gchar * name, const gchar * guitext, const gchar * desc, const Parameter::_scope_t scope, bool gui_hidden, const gchar * gui_tip, Inkscape::Extension::Extension * ext, Inkscape::XML::Node * xml) : - Parameter(name, guitext, desc, scope, gui_hidden, gui_tip, ext), _indent(0) +ParamComboBox::ParamComboBox(const gchar *name, const gchar *guitext, const gchar *desc, + const Parameter::_scope_t scope, bool gui_hidden, const gchar *gui_tip, + Inkscape::Extension::Extension *ext, Inkscape::XML::Node *xml) + : Parameter(name, guitext, desc, scope, gui_hidden, gui_tip, ext) + , _value(NULL) + , _indent(0) + , choices(NULL) { - choices = NULL; - _value = NULL; + const char *xmlval = NULL; // the value stored in XML - // Read XML tree to add enumeration items: - // printf("Extension Constructor: "); if (xml != NULL) { - for (Inkscape::XML::Node *node = xml->firstChild(); node; node = node->next()) { + // Read XML tree to add enumeration items: + for (Inkscape::XML::Node *node = xml->firstChild(); node; node = node->next()) { char const * chname = node->name(); if (!strcmp(chname, INKSCAPE_EXTENSION_NS "item") || !strcmp(chname, INKSCAPE_EXTENSION_NS "_item")) { Glib::ustring newguitext, newvalue; @@ -69,8 +72,8 @@ ParamComboBox::ParamComboBox (const gchar * name, const gchar * guitext, const g } if (contents != NULL) { // don't translate when 'item' but do translate when '_item' - // NOTE: internal extensions use build_from_mem and don't need _item but - // still need to include if are to be localized + // NOTE: internal extensions use build_from_mem and don't need _item but + // still need to include if are to be localized if (!strcmp(chname, INKSCAPE_EXTENSION_NS "_item")) { if (node->attribute("msgctxt") != NULL) { newguitext = g_dpgettext2(NULL, node->attribute("msgctxt"), contents); @@ -95,30 +98,28 @@ ParamComboBox::ParamComboBox (const gchar * name, const gchar * guitext, const g } } } - } - - // Initialize _value with the default value from xml - // for simplicity : default to the contents of the first xml-child - const char * defaultval = NULL; - if (xml->firstChild() && xml->firstChild()->firstChild()) { - defaultval = xml->firstChild()->attribute("value"); - } + + // Initialize _value with the default value from xml + // for simplicity : default to the contents of the first xml-child + if (xml->firstChild() && xml->firstChild()->firstChild()) { + xmlval = xml->firstChild()->attribute("value"); + } - const char * indent = xml->attribute("indent"); - if (indent != NULL) { - _indent = atoi(indent) * 12; + const char *indent = xml->attribute("indent"); + if (indent != NULL) { + _indent = atoi(indent) * 12; + } } gchar * pref_name = this->pref_name(); Inkscape::Preferences *prefs = Inkscape::Preferences::get(); - Glib::ustring paramval = prefs->getString(extension_pref_root + pref_name); + Glib::ustring paramval = prefs ? prefs->getString(extension_pref_root + pref_name) : ""; g_free(pref_name); if (!paramval.empty()) { - defaultval = paramval.data(); - } - if (defaultval != NULL) { - _value = g_strdup(defaultval); + _value = g_strdup(paramval.data()); + } else if (xmlval) { + _value = g_strdup(xmlval); } } -- cgit v1.2.3 From b021cb56b6da1e6e08da95209e75bb3d102d5267 Mon Sep 17 00:00:00 2001 From: "Johan B. C. Engelen" Date: Sun, 12 Oct 2014 17:55:41 +0200 Subject: noop: remove duplicate assignment in sp-lpe-item.cpp (bzr r13600) --- src/sp-lpe-item.cpp | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) (limited to 'src') diff --git a/src/sp-lpe-item.cpp b/src/sp-lpe-item.cpp index 321d2fc42..17cacbae2 100644 --- a/src/sp-lpe-item.cpp +++ b/src/sp-lpe-item.cpp @@ -655,9 +655,8 @@ sp_lpe_item_apply_to_mask(SPItem * item) { SPMask *mask = item->mask_ref->getObject(); if(mask) { - SPObject * mask_data = mask->firstChild(); + SPObject *mask_data = mask->firstChild(); SPCurve * mask_curve = NULL; - mask_data = mask->firstChild(); if (SP_IS_PATH(mask_data)) { mask_curve = SP_PATH(mask_data)->get_original_curve(); } else if(SP_IS_SHAPE(mask_data)) { -- cgit v1.2.3 From acb22a496132ae58d7a500b4f92e16786d8a27c9 Mon Sep 17 00:00:00 2001 From: "Johan B. C. Engelen" Date: Sun, 12 Oct 2014 17:59:21 +0200 Subject: object-snapper.cpp : make the logic easier to read (use else-clause as default initialization value) (bzr r13601) --- src/object-snapper.cpp | 5 +---- 1 file changed, 1 insertion(+), 4 deletions(-) (limited to 'src') diff --git a/src/object-snapper.cpp b/src/object-snapper.cpp index 9ff6df3bf..2b29814b0 100644 --- a/src/object-snapper.cpp +++ b/src/object-snapper.cpp @@ -289,7 +289,7 @@ void Inkscape::ObjectSnapper::_snapNodes(IntermSnapResults &isr, for (std::vector::const_iterator k = _points_to_snap_to->begin(); k != _points_to_snap_to->end(); ++k) { if (_allowSourceToSnapToTarget(p.getSourceType(), (*k).getTargetType(), strict_snapping)) { Geom::Point target_pt = (*k).getPoint(); - Geom::Coord dist = Geom::infinity(); + Geom::Coord dist = Geom::L2(target_pt - p.getPoint()); // Default: free (unconstrained) snapping if (!c.isUndefined()) { // We're snapping to nodes along a constraint only, so find out if this node // is at the constraint, while allowing for a small margin @@ -299,9 +299,6 @@ void Inkscape::ObjectSnapper::_snapNodes(IntermSnapResults &isr, continue; } dist = Geom::L2(target_pt - p_proj_on_constraint); - } else { - // Free (unconstrained) snapping - dist = Geom::L2(target_pt - p.getPoint()); } if (dist < getSnapperTolerance() && dist < s.getSnapDistance()) { -- cgit v1.2.3 From 1d31429ac10d5fa47a611cfc2e89f4c313c1ad16 Mon Sep 17 00:00:00 2001 From: "Johan B. C. Engelen" Date: Sun, 12 Oct 2014 18:19:14 +0200 Subject: path combine: add an assert to indicate that item list cannot be empty (bzr r13602) --- src/path-chemistry.cpp | 1 + 1 file changed, 1 insertion(+) (limited to 'src') diff --git a/src/path-chemistry.cpp b/src/path-chemistry.cpp index d355b49fe..787193060 100644 --- a/src/path-chemistry.cpp +++ b/src/path-chemistry.cpp @@ -78,6 +78,7 @@ sp_selected_path_combine(SPDesktop *desktop) items = g_slist_sort(items, (GCompareFunc) sp_item_repr_compare_position); items = g_slist_reverse(items); + assert(items); // cannot be NULL because of list length check at top of function // remember the position, id, transform and style of the topmost path, they will be assigned to the combined one gint position = 0; -- cgit v1.2.3 From 2f2f164bbe75b4a4bb6d553242cecd587ed03456 Mon Sep 17 00:00:00 2001 From: "Johan B. C. Engelen" Date: Sun, 12 Oct 2014 19:25:47 +0200 Subject: DrawingItem: add NULL check (bug found by clang static analyzer) (bzr r13603) --- src/display/drawing-item.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) (limited to 'src') diff --git a/src/display/drawing-item.cpp b/src/display/drawing-item.cpp index 80eb69546..2ea3641dc 100644 --- a/src/display/drawing-item.cpp +++ b/src/display/drawing-item.cpp @@ -967,7 +967,7 @@ DrawingItem::_setStyleCommon(SPStyle *&_style, SPStyle *style) if (_style) sp_style_unref(_style); _style = style; - if (style->filter.set && style->getFilter()) { + if (style && style->filter.set && style->getFilter()) { if (!_filter) { int primitives = sp_filter_primitive_count(SP_FILTER(style->getFilter())); _filter = new Inkscape::Filters::Filter(primitives); -- cgit v1.2.3 From 932ab442403e4007b7dfb8e4c587e5d979379df1 Mon Sep 17 00:00:00 2001 From: Tavmjong Bah Date: Mon, 13 Oct 2014 10:05:35 +0200 Subject: Fix: "Forming reference to null pointer" (bzr r13606) --- src/style-internal.cpp | 11 ++++++++++- 1 file changed, 10 insertions(+), 1 deletion(-) (limited to 'src') diff --git a/src/style-internal.cpp b/src/style-internal.cpp index a55a11403..f8a4a43f7 100644 --- a/src/style-internal.cpp +++ b/src/style-internal.cpp @@ -981,7 +981,16 @@ SPIPaint::read( gchar const *str ) { if (streq(str, "currentColor")) { set = true; currentcolor = true; - setColor( style->color.value.color ); + if (style) { + setColor( style->color.value.color ); + } else { + // Normally an SPIPaint is part of an SPStyle and the value of 'color' is + // available. SPIPaint can be used 'stand-alone' (e.g. to parse color values) in + // which case a value of 'currentColor' is meaningless, thus we shouldn't reach + // here. + std::cerr << "SPIPaint::read(): value is 'currentColor' but 'color' not available." << std::endl; + setColor( 0 ); + } } else if (streq(str, "none")) { set = true; noneSet = true; -- cgit v1.2.3 From fd587105dcf2e883a6712013b7ea7ce6e3b3b11d Mon Sep 17 00:00:00 2001 From: Tavmjong Bah Date: Mon, 13 Oct 2014 10:21:07 +0200 Subject: Fix "Called C++ object pointer is null" (bzr r13607) --- src/style-internal.cpp | 15 ++++++++++++--- 1 file changed, 12 insertions(+), 3 deletions(-) (limited to 'src') diff --git a/src/style-internal.cpp b/src/style-internal.cpp index f8a4a43f7..a81a2799d 100644 --- a/src/style-internal.cpp +++ b/src/style-internal.cpp @@ -1405,11 +1405,20 @@ SPIFilter::read( gchar const *str ) { set = true; // Create href if not already done. - if (!href && style->object) { - href = new SPFilterReference(style->object); - href->changedSignal().connect(sigc::bind(sigc::ptr_fun(sp_style_filter_ref_changed), style)); + if (!href) { + if (style->object) { + href = new SPFilterReference(style->object); + } + // Do we have href now? + if ( href ) { + href->changedSignal().connect(sigc::bind(sigc::ptr_fun(sp_style_filter_ref_changed), style)); + } else { + std::cerr << "SPIFilter::read(): Could not allocate 'href'" << std::endl; + return; + } } + // We have href try { href->attach(Inkscape::URI(uri)); } catch (Inkscape::BadURIException &e) { -- cgit v1.2.3 From 25f2a039c44b34fce5f0f4eb84259c4c1755aa62 Mon Sep 17 00:00:00 2001 From: Tavmjong Bah Date: Mon, 13 Oct 2014 10:50:39 +0200 Subject: Fix multiple: "Called C++ object pointer is null" (bzr r13608) --- src/widgets/stroke-style.cpp | 6 ++++++ 1 file changed, 6 insertions(+) (limited to 'src') diff --git a/src/widgets/stroke-style.cpp b/src/widgets/stroke-style.cpp index 0e0a4fd72..9c8d45908 100644 --- a/src/widgets/stroke-style.cpp +++ b/src/widgets/stroke-style.cpp @@ -763,6 +763,9 @@ StrokeStyle::setJoinType (unsigned const jointype) tb = joinBevel; break; default: + // Should not happen + std::cerr << "StrokeStyle::setJoinType(): Invalid value" << std::endl; + tb = joinMiter; break; } setJoinButtons(tb); @@ -786,6 +789,9 @@ StrokeStyle::setCapType (unsigned const captype) tb = capSquare; break; default: + // Should not happen + std::cerr << "StrokeStyle::setCapType(): Invalid value" << std::endl; + tb = capButt; break; } setCapButtons(tb); -- cgit v1.2.3 From a2c47e60a02938013ab913bfee2ef9c9fa5ea870 Mon Sep 17 00:00:00 2001 From: Tavmjong Bah Date: Mon, 13 Oct 2014 11:11:01 +0200 Subject: Fixed: "Called C++ object pointer is null" (bzr r13609) --- src/ui/dialog/clonetiler.cpp | 7 ++++++- 1 file changed, 6 insertions(+), 1 deletion(-) (limited to 'src') diff --git a/src/ui/dialog/clonetiler.cpp b/src/ui/dialog/clonetiler.cpp index e5f18216c..9b67132c0 100644 --- a/src/ui/dialog/clonetiler.cpp +++ b/src/ui/dialog/clonetiler.cpp @@ -2235,7 +2235,11 @@ void CloneTiler::clonetiler_apply(GtkWidget */*widget*/, GtkWidget *dlg) gdk_window_process_all_updates(); SPObject *obj = selection->singleItem(); - SPItem *item = SP_IS_ITEM(obj) ? SP_ITEM(obj) : 0; + if (!obj) { + // Should never happen (empty selection checked above). + std::cerr << "CloneTiler::clonetile_apply(): No object in single item selection!!!" << std::endl; + return; + } Inkscape::XML::Node *obj_repr = obj->getRepr(); const char *id_href = g_strdup_printf("#%s", obj_repr->attribute("id")); SPObject *parent = obj->parent; @@ -2326,6 +2330,7 @@ void CloneTiler::clonetiler_apply(GtkWidget */*widget*/, GtkWidget *dlg) bool invert_picked = prefs->getBool(prefs_path + "invert_picked"); double gamma_picked = prefs->getDoubleLimited(prefs_path + "gamma_picked", 0, -10, 10); + SPItem *item = SP_IS_ITEM(obj) ? SP_ITEM(obj) : 0; if (dotrace) { clonetiler_trace_setup (sp_desktop_document(desktop), 1.0, item); } -- cgit v1.2.3 From 9c35de31702bb501f596555d347e9db76b248f14 Mon Sep 17 00:00:00 2001 From: Bryce Harrington Date: Mon, 13 Oct 2014 19:33:59 -0700 Subject: 2geom: Re-fix IS_NAN and IS_FINITE for solaris This patch should go upstream to 2geom... (bzr r13610) --- src/2geom/math-utils.h | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) (limited to 'src') diff --git a/src/2geom/math-utils.h b/src/2geom/math-utils.h index 83c7b4f5e..0d5d3667b 100644 --- a/src/2geom/math-utils.h +++ b/src/2geom/math-utils.h @@ -111,7 +111,7 @@ inline void sincos(double angle, double &sin_, double &cos_) { # define IS_NAN(_a) (_isnan(_a)) /* Win32 definition */ #elif defined(isnan) || defined(__FreeBSD__) || defined(__osf__) # define IS_NAN(_a) (isnan(_a)) /* GNU definition */ -#elif defined (SOLARIS_2_8) && __GNUC__ == 3 && __GNUC_MINOR__ == 2 +#elif defined (SOLARIS) # define IS_NAN(_a) (isnan(_a)) /* GNU definition */ #else # define IS_NAN(_a) (boost::math::isnan(_a)) @@ -129,7 +129,7 @@ inline void sincos(double angle, double &sin_, double &cos_) { # define IS_FINITE(_a) (isfinite(_a)) #elif defined(__osf__) # define IS_FINITE(_a) (finite(_a) && !IS_NAN(_a)) -#elif defined (SOLARIS_2_8) && __GNUC__ == 3 && __GNUC_MINOR__ == 2 +#elif defined (SOLARIS) #include #define IS_FINITE(_a) (finite(_a) && !IS_NAN(_a)) #else -- cgit v1.2.3