summaryrefslogtreecommitdiffstats
path: root/src/path
diff options
context:
space:
mode:
Diffstat (limited to 'src/path')
-rw-r--r--src/path/CMakeLists.txt27
-rw-r--r--src/path/README11
-rw-r--r--src/path/path-boolop.cpp936
-rw-r--r--src/path/path-boolop.h34
-rw-r--r--src/path/path-object-set.cpp166
-rw-r--r--src/path/path-offset.cpp468
-rw-r--r--src/path/path-offset.h41
-rw-r--r--src/path/path-outline.cpp792
-rw-r--r--src/path/path-outline.h60
-rw-r--r--src/path/path-simplify.cpp126
-rw-r--r--src/path/path-simplify.h28
-rw-r--r--src/path/path-util.cpp184
-rw-r--r--src/path/path-util.h162
-rw-r--r--src/path/splinefit/bezier-fit.cpp90
-rw-r--r--src/path/splinefit/bezier-fit.h20
-rw-r--r--src/path/splinefit/splinefit.c1427
-rw-r--r--src/path/splinefit/splinefit.h78
-rw-r--r--src/path/splinefit/splinefont.c1174
-rw-r--r--src/path/splinefit/splinefont.h191
-rw-r--r--src/path/splinefit/splinerefigure.c117
-rw-r--r--src/path/splinefit/splinerefigure.h9
21 files changed, 6141 insertions, 0 deletions
diff --git a/src/path/CMakeLists.txt b/src/path/CMakeLists.txt
new file mode 100644
index 0000000..cbdac41
--- /dev/null
+++ b/src/path/CMakeLists.txt
@@ -0,0 +1,27 @@
+# SPDX-License-Identifier: GPL-2.0-or-later
+
+
+set(path_SRC
+ path-boolop.cpp
+ path-object-set.cpp
+ path-offset.cpp
+ path-outline.cpp
+ path-simplify.cpp
+ path-util.cpp
+ splinefit/bezier-fit.cpp
+ splinefit/splinefit.c
+ splinefit/splinefont.c
+ splinefit/splinerefigure.c
+
+ # -------
+ # Headers
+ path-boolop.h
+ path-offset.h
+ path-outline.h
+ path-simplify.h
+ path-util.h
+ splinefit/splinefit.h
+ splinefit/splinefont.h
+)
+
+add_inkscape_source("${path_SRC}")
diff --git a/src/path/README b/src/path/README
new file mode 100644
index 0000000..309e2d6
--- /dev/null
+++ b/src/path/README
@@ -0,0 +1,11 @@
+
+
+This directory contains code related to path manipulations.
+
+To do:
+
+* Move all path manipulation related code here (livarot, etc.).
+* Merge code into lib2geom if possible.
+* Get rid of livarot
+* Remove GUI dependencies.
+* Find a better place to put directory.
diff --git a/src/path/path-boolop.cpp b/src/path/path-boolop.cpp
new file mode 100644
index 0000000..b4c4922
--- /dev/null
+++ b/src/path/path-boolop.cpp
@@ -0,0 +1,936 @@
+// SPDX-License-Identifier: GPL-2.0-or-later
+/** @file
+ * Boolean operations.
+ *//*
+ * Authors:
+ * see git history
+ * Created by fred on Fri Dec 05 2003.
+ * tweaked endlessly by bulia byak <buliabyak@users.sf.net>
+ *
+ * Copyright (C) 2018 Authors
+ * Released under GNU GPL v2+, read the file 'COPYING' for more information.
+ */
+
+#include <vector>
+
+#include <glibmm/i18n.h>
+#include <2geom/intersection-graph.h>
+#include <2geom/svg-path-parser.h> // to get from SVG on boolean to Geom::Path
+#include <2geom/utils.h>
+
+#include "path-boolop.h"
+#include "path-util.h"
+
+#include "message-stack.h"
+#include "path-chemistry.h" // copy_object_properties()
+
+#include "helper/geom.h" // pathv_to_linear_and_cubic_beziers()
+
+#include "livarot/Path.h"
+#include "livarot/Shape.h"
+
+#include "object/sp-flowtext.h"
+#include "object/sp-shape.h"
+#include "object/sp-text.h"
+
+#include "display/curve.h"
+
+#include "ui/icon-names.h"
+
+#include "svg/svg.h"
+
+#include "xml/repr-sorting.h"
+
+using Inkscape::DocumentUndo;
+
+// fonctions utilitaires
+bool
+Ancetre(Inkscape::XML::Node *a, Inkscape::XML::Node *who)
+{
+ if (who == nullptr || a == nullptr)
+ return false;
+ if (who == a)
+ return true;
+ return Ancetre(a->parent(), who);
+}
+
+
+bool Inkscape::ObjectSet::pathUnion(const bool skip_undo, bool silent) {
+ BoolOpErrors result = pathBoolOp(bool_op_union, skip_undo, false, INKSCAPE_ICON("path-union"),
+ _("Union"), silent);
+ return DONE == result;
+}
+
+bool
+Inkscape::ObjectSet::pathIntersect(const bool skip_undo, bool silent)
+{
+ BoolOpErrors result = pathBoolOp(bool_op_inters, skip_undo, false, INKSCAPE_ICON("path-intersection"),
+ _("Intersection"), silent);
+ return DONE == result;
+}
+
+bool
+Inkscape::ObjectSet::pathDiff(const bool skip_undo, bool silent)
+{
+ BoolOpErrors result = pathBoolOp(bool_op_diff, skip_undo, false, INKSCAPE_ICON("path-difference"),
+ _("Difference"), silent);
+ return DONE == result;
+}
+
+bool
+Inkscape::ObjectSet::pathSymDiff(const bool skip_undo, bool silent)
+{
+ BoolOpErrors result = pathBoolOp(bool_op_symdiff, skip_undo, false, INKSCAPE_ICON("path-exclusion"),
+ _("Exclusion"), silent);
+ return DONE == result;
+}
+
+bool
+Inkscape::ObjectSet::pathCut(const bool skip_undo, bool silent)
+{
+ BoolOpErrors result = pathBoolOp(bool_op_cut, skip_undo, false, INKSCAPE_ICON("path-division"),
+ _("Division"), silent);
+ return DONE == result;
+}
+
+bool
+Inkscape::ObjectSet::pathSlice(const bool skip_undo, bool silent)
+{
+ BoolOpErrors result = pathBoolOp(bool_op_slice, skip_undo, false, INKSCAPE_ICON("path-cut"),
+ _("Cut path"), silent);
+ return DONE == result;
+}
+
+// helper for printing error messages, regardless of whether we have a GUI or not
+// If desktop == NULL, errors will be shown on stderr
+static void
+boolop_display_error_message(SPDesktop *desktop, Glib::ustring const &msg)
+{
+ if (desktop) {
+ desktop->messageStack()->flash(Inkscape::ERROR_MESSAGE, msg);
+ } else {
+ g_printerr("%s\n", msg.c_str());
+ }
+}
+
+/**
+ * Calculate the threshold for the given PathVector based
+ * on it's bounding box.
+ *
+ * @param path - The PathVector to calculate the threshold for.
+ * @param threshold - The starting threshold, usually 0.1
+ */
+double get_threshold(Geom::PathVector const &path, double threshold)
+{
+ auto maybe_box = path.boundsFast();
+ if (!maybe_box)
+ return 0.0;
+ Geom::Rect box = *maybe_box;
+ double diagonal = Geom::distance(
+ Geom::Point(box[Geom::X].min(), box[Geom::Y].min()),
+ Geom::Point(box[Geom::X].max(), box[Geom::Y].max())
+ );
+ return threshold * (diagonal / 100);
+}
+
+void
+sp_flatten(Geom::PathVector &pathvector, FillRule fillkind)
+{
+ Path *orig = new Path;
+ orig->LoadPathVector(pathvector);
+ Shape *theShape = new Shape;
+ Shape *theRes = new Shape;
+ orig->ConvertWithBackData(get_threshold(pathvector, 0.1));
+ orig->Fill(theShape, 0);
+ theRes->ConvertToShape(theShape, fillkind);
+ Path *originaux[1];
+ originaux[0] = orig;
+ Path *res = new Path;
+ theRes->ConvertToForme(res, 1, originaux, true);
+
+ delete theShape;
+ delete theRes;
+ pathvector = res->MakePathVector();
+ delete res;
+ delete orig;
+}
+
+// boolean operations PathVectors A,B -> PathVector result.
+// This is derived from sp_selected_path_boolop
+// take the source paths from the file, do the operation, delete the originals and add the results
+// fra,fra are fill_rules for PathVectors a,b
+Geom::PathVector
+sp_pathvector_boolop(Geom::PathVector const &pathva, Geom::PathVector const &pathvb, bool_op bop,
+ fill_typ fra, fill_typ frb, bool livarotonly, bool flattenbefore)
+{
+ int error = 0;
+ return sp_pathvector_boolop(pathva, pathvb, bop, fra, frb, livarotonly, flattenbefore, error);
+}
+
+Geom::PathVector
+sp_pathvector_boolop(Geom::PathVector const &pathva, Geom::PathVector const &pathvb, bool_op bop,
+ fill_typ fra, fill_typ frb, bool livarotonly, bool flattenbefore, int &error)
+{
+ if (!livarotonly) {
+ try {
+ Geom::PathVector a = pathv_to_linear_and_cubic_beziers(pathva);
+ Geom::PathVector b = pathv_to_linear_and_cubic_beziers(pathvb);
+ if (flattenbefore) {
+ sp_flatten(a, fra);
+ sp_flatten(b, frb);
+ }
+ Geom::PathVector out;
+ // dont change tolerande give errors on boolops
+ auto pig = Geom::PathIntersectionGraph(a, b, Geom::EPSILON);
+ if (bop == bool_op_inters) {
+ out = pig.getIntersection();
+ } else if (bop == bool_op_union) {
+ out = pig.getUnion();
+
+ } else if (bop == bool_op_symdiff) {
+ out = pig.getXOR();
+ } else if (bop == bool_op_diff) {
+ out = pig.getBminusA(); //livarot order...
+ } else if (bop == bool_op_cut) {
+ out = pig.getBminusA();
+ auto tmp = pig.getIntersection();
+ out.insert(out.end(), tmp.begin(), tmp.end());
+ } else if (bop == bool_op_slice) {
+ // go to livarot
+ livarotonly = true;
+ }
+ if (!livarotonly) {
+ return out;
+ }
+ } catch (...) {
+ g_debug("Path Intersection Graph failed boolops, fallback to livarot");
+ }
+ }
+ error = 1;
+ // extract the livarot Paths from the source objects
+ // also get the winding rule specified in the style
+ int nbOriginaux = 2;
+ std::vector<Path *> originaux(nbOriginaux);
+ std::vector<FillRule> origWind(nbOriginaux);
+ origWind[0]=fra;
+ origWind[1]=frb;
+ Geom::PathVector patht;
+ // Livarot's outline of arcs is broken. So convert the path to linear and cubics only, for which the outline is created correctly.
+ originaux[0] = Path_for_pathvector(pathv_to_linear_and_cubic_beziers( pathva));
+ originaux[1] = Path_for_pathvector(pathv_to_linear_and_cubic_beziers( pathvb));
+
+ // some temporary instances, first
+ Shape *theShapeA = new Shape;
+ Shape *theShapeB = new Shape;
+ Shape *theShape = new Shape;
+ Path *res = new Path;
+ res->SetBackData(false);
+
+ Path::cut_position *toCut=nullptr;
+ int nbToCut = 0;
+
+ if ( bop == bool_op_inters || bop == bool_op_union || bop == bool_op_diff || bop == bool_op_symdiff ) {
+ // true boolean op
+ // get the polygons of each path, with the winding rule specified, and apply the operation iteratively
+ originaux[0]->ConvertWithBackData(get_threshold(pathva, 0.1));
+
+ originaux[0]->Fill(theShape, 0);
+
+ theShapeA->ConvertToShape(theShape, origWind[0]);
+
+ originaux[1]->ConvertWithBackData(get_threshold(pathvb, 0.1));
+
+ originaux[1]->Fill(theShape, 1);
+
+ theShapeB->ConvertToShape(theShape, origWind[1]);
+
+ theShape->Booleen(theShapeB, theShapeA, bop);
+
+ } else if ( bop == bool_op_cut ) {
+ // cuts= sort of a bastard boolean operation, thus not the axact same modus operandi
+ // technically, the cut path is not necessarily a polygon (thus has no winding rule)
+ // it is just uncrossed, and cleaned from duplicate edges and points
+ // then it's fed to Booleen() which will uncross it against the other path
+ // then comes the trick: each edge of the cut path is duplicated (one in each direction),
+ // thus making a polygon. the weight of the edges of the cut are all 0, but
+ // the Booleen need to invert the ones inside the source polygon (for the subsequent
+ // ConvertToForme)
+
+ // the cut path needs to have the highest pathID in the back data
+ // that's how the Booleen() function knows it's an edge of the cut
+ {
+ Path* swap=originaux[0];originaux[0]=originaux[1];originaux[1]=swap;
+ int swai=origWind[0];origWind[0]=origWind[1];origWind[1]=(fill_typ)swai;
+ }
+ originaux[0]->ConvertWithBackData(get_threshold(pathvb, 0.1));
+
+ originaux[0]->Fill(theShape, 0);
+
+ theShapeA->ConvertToShape(theShape, origWind[0]);
+
+ originaux[1]->ConvertWithBackData(get_threshold(pathva, 0.1));
+
+ if ((originaux[1]->pts.size() == 2) && originaux[1]->pts[0].isMoveTo && !originaux[1]->pts[1].isMoveTo)
+ originaux[1]->Fill(theShape, 1,false,true,false); // see LP Bug 177956
+ else
+ originaux[1]->Fill(theShape, 1,false,false,false); //do not closeIfNeeded
+
+ theShapeB->ConvertToShape(theShape, fill_justDont); // fill_justDont doesn't computes winding numbers
+
+ // les elements arrivent en ordre inverse dans la liste
+ theShape->Booleen(theShapeB, theShapeA, bool_op_cut, 1);
+
+ } else if ( bop == bool_op_slice ) {
+ // slice is not really a boolean operation
+ // you just put the 2 shapes in a single polygon, uncross it
+ // the points where the degree is > 2 are intersections
+ // just check it's an intersection on the path you want to cut, and keep it
+ // the intersections you have found are then fed to ConvertPositionsToMoveTo() which will
+ // make new subpath at each one of these positions
+ // inversion pour l'opération
+ {
+ Path* swap=originaux[0];originaux[0]=originaux[1];originaux[1]=swap;
+ int swai=origWind[0];origWind[0]=origWind[1];origWind[1]=(fill_typ)swai;
+ }
+ originaux[0]->ConvertWithBackData(get_threshold(pathvb, 0.1));
+
+ originaux[0]->Fill(theShapeA, 0,false,false,false); // don't closeIfNeeded
+
+ originaux[1]->ConvertWithBackData(get_threshold(pathva, 0.1));
+
+ originaux[1]->Fill(theShapeA, 1,true,false,false);// don't closeIfNeeded and just dump in the shape, don't reset it
+
+ theShape->ConvertToShape(theShapeA, fill_justDont);
+
+ if ( theShape->hasBackData() ) {
+ // should always be the case, but ya never know
+ {
+ for (int i = 0; i < theShape->numberOfPoints(); i++) {
+ if ( theShape->getPoint(i).totalDegree() > 2 ) {
+ // possibly an intersection
+ // we need to check that at least one edge from the source path is incident to it
+ // before we declare it's an intersection
+ int cb = theShape->getPoint(i).incidentEdge[FIRST];
+ int nbOrig=0;
+ int nbOther=0;
+ int piece=-1;
+ float t=0.0;
+ while ( cb >= 0 && cb < theShape->numberOfEdges() ) {
+ if ( theShape->ebData[cb].pathID == 0 ) {
+ // the source has an edge incident to the point, get its position on the path
+ piece=theShape->ebData[cb].pieceID;
+ if ( theShape->getEdge(cb).st == i ) {
+ t=theShape->ebData[cb].tSt;
+ } else {
+ t=theShape->ebData[cb].tEn;
+ }
+ nbOrig++;
+ }
+ if ( theShape->ebData[cb].pathID == 1 ) nbOther++; // the cut is incident to this point
+ cb=theShape->NextAt(i, cb);
+ }
+ if ( nbOrig > 0 && nbOther > 0 ) {
+ // point incident to both path and cut: an intersection
+ // note that you only keep one position on the source; you could have degenerate
+ // cases where the source crosses itself at this point, and you wouyld miss an intersection
+ toCut=(Path::cut_position*)realloc(toCut, (nbToCut+1)*sizeof(Path::cut_position));
+ toCut[nbToCut].piece=piece;
+ toCut[nbToCut].t=t;
+ nbToCut++;
+ }
+ }
+ }
+ }
+ {
+ // i think it's useless now
+ int i = theShape->numberOfEdges() - 1;
+ for (;i>=0;i--) {
+ if ( theShape->ebData[i].pathID == 1 ) {
+ theShape->SubEdge(i);
+ }
+ }
+ }
+
+ }
+ }
+
+ // pour compenser le swap juste avant
+ if ( bop == bool_op_slice ) {
+ res->Copy(originaux[0]);
+ res->ConvertPositionsToMoveTo(nbToCut, toCut); // cut where you found intersections
+ free(toCut);
+ } else {
+ theShape->ConvertToForme(res, nbOriginaux, &originaux[0]);
+ }
+
+ delete theShape;
+ delete theShapeA;
+ delete theShapeB;
+ delete originaux[0];
+ delete originaux[1];
+
+ auto outres = res->MakePathVector();
+
+ delete res;
+ return outres;
+}
+
+/**
+ * Workaround for buggy Path::Transform() which incorrectly transforms arc commands.
+ *
+ * TODO: Fix PathDescrArcTo::transform() and then remove this workaround.
+ */
+static void transformLivarotPath(Path *res, Geom::Affine const &affine)
+{
+ res->LoadPathVector(res->MakePathVector() * affine);
+}
+
+// boolean operations on the desktop
+// take the source paths from the file, do the operation, delete the originals and add the results
+BoolOpErrors Inkscape::ObjectSet::pathBoolOp(bool_op bop, const bool skip_undo, const bool checked,
+ const Glib::ustring icon_name, const Glib::ustring description,
+ bool const silent)
+{
+ if (nullptr != desktop() && !checked) {
+ SPDocument *doc = desktop()->getDocument();
+ BoolOpErrors returnCode = ObjectSet::pathBoolOp(bop, true, true,icon_name);
+
+ switch(returnCode) {
+ case ERR_TOO_LESS_PATHS_1:
+ if (!silent) {
+ boolop_display_error_message(desktop(),
+ _("Select <b>at least 1 path</b> to perform a boolean union."));
+ }
+ break;
+ case ERR_TOO_LESS_PATHS_2:
+ if (!silent) {
+ boolop_display_error_message(desktop(),
+ _("Select <b>at least 2 paths</b> to perform a boolean operation."));
+ }
+ break;
+ case ERR_NO_PATHS:
+ if (!silent) {
+ boolop_display_error_message(desktop(),
+ _("One of the objects is <b>not a path</b>, cannot perform boolean operation."));
+ }
+ break;
+ case ERR_Z_ORDER:
+ if (!silent) {
+ boolop_display_error_message(desktop(),
+ _("Unable to determine the <b>z-order</b> of the objects selected for difference, XOR, division, or path cut."));
+ }
+ break;
+ case DONE_NO_PATH:
+ if (!skip_undo) {
+ DocumentUndo::done(doc, description, "");
+ }
+ break;
+ case DONE:
+ if (!skip_undo) {
+ DocumentUndo::done(doc, description, icon_name);
+ }
+ break;
+ case DONE_NO_ACTION:
+ // Do nothing (?)
+ break;
+ }
+ return returnCode;
+ }
+
+ SPDocument *doc = document();
+ std::vector<SPItem*> il(items().begin(), items().end());
+
+ // allow union on a single object for the purpose of removing self overlapse (svn log, revision 13334)
+ if (il.size() < 2 && bop != bool_op_union) {
+ return ERR_TOO_LESS_PATHS_2;
+ }
+ else if (il.size() < 1) {
+ return ERR_TOO_LESS_PATHS_1;
+ }
+
+ g_assert(!il.empty());
+
+ // reverseOrderForOp marks whether the order of the list is the top->down order
+ // it's only used when there are 2 objects, and for operations who need to know the
+ // topmost object (differences, cuts)
+ bool reverseOrderForOp = false;
+
+ if (bop == bool_op_diff || bop == bool_op_cut || bop == bool_op_slice) {
+ // check in the tree to find which element of the selection list is topmost (for 2-operand commands only)
+ Inkscape::XML::Node *a = il.front()->getRepr();
+ Inkscape::XML::Node *b = il.back()->getRepr();
+
+ if (a == nullptr || b == nullptr) {
+ return ERR_Z_ORDER;
+ }
+
+ if (Ancetre(a, b)) {
+ // a is the parent of b, already in the proper order
+ } else if (Ancetre(b, a)) {
+ // reverse order
+ reverseOrderForOp = true;
+ } else {
+
+ // objects are not in parent/child relationship;
+ // find their lowest common ancestor
+ Inkscape::XML::Node *parent = LCA(a, b);
+ if (parent == nullptr) {
+ return ERR_Z_ORDER;
+ }
+
+ // find the children of the LCA that lead from it to the a and b
+ Inkscape::XML::Node *as = AncetreFils(a, parent);
+ Inkscape::XML::Node *bs = AncetreFils(b, parent);
+
+ // find out which comes first
+ for (Inkscape::XML::Node *child = parent->firstChild(); child; child = child->next()) {
+ if (child == as) {
+ /* a first, so reverse. */
+ reverseOrderForOp = true;
+ break;
+ }
+ if (child == bs)
+ break;
+ }
+ }
+ }
+
+ g_assert(!il.empty());
+
+ // first check if all the input objects have shapes
+ // otherwise bail out
+ for (auto item : il)
+ {
+ if (!is<SPShape>(item) && !is<SPText>(item) && !is<SPFlowtext>(item))
+ {
+ return ERR_NO_PATHS;
+ }
+ }
+
+ // extract the livarot Paths from the source objects
+ // also get the winding rule specified in the style
+ int nbOriginaux = il.size();
+ std::vector<Path *> originaux(nbOriginaux);
+ std::vector<FillRule> origWind(nbOriginaux);
+ std::vector<double> origThresh(nbOriginaux);
+ int curOrig;
+ {
+ curOrig = 0;
+ for (auto item : il)
+ {
+ // apply live path effects prior to performing boolean operation
+ char const *id = item->getAttribute("id");
+ auto lpeitem = cast<SPLPEItem>(item);
+ if (lpeitem) {
+ SPDocument * document = item->document;
+ lpeitem->removeAllPathEffects(true);
+ SPObject *elemref = document->getObjectById(id);
+ if (elemref && elemref != item) {
+ // If the LPE item is a shape, it is converted to a path
+ // so we need to reupdate the item
+ item = cast<SPItem>(elemref);
+ }
+ }
+ SPCSSAttr *css = sp_repr_css_attr(reinterpret_cast<SPObject *>(il[0])->getRepr(), "style");
+ gchar const *val = sp_repr_css_property(css, "fill-rule", nullptr);
+ if (val && strcmp(val, "nonzero") == 0) {
+ origWind[curOrig]= fill_nonZero;
+ } else if (val && strcmp(val, "evenodd") == 0) {
+ origWind[curOrig]= fill_oddEven;
+ } else {
+ origWind[curOrig]= fill_nonZero;
+ }
+ sp_repr_css_attr_unref(css);
+
+ if (auto curve = curve_for_item(item)) {
+ auto pathv = curve->get_pathvector() * item->i2doc_affine();
+ originaux[curOrig] = Path_for_pathvector(pathv);
+ origThresh[curOrig] = get_threshold(pathv, 0.1);
+ } else {
+ originaux[curOrig] = nullptr;
+ }
+
+ if (originaux[curOrig] == nullptr || originaux[curOrig]->descr_cmd.size() <= 1)
+ {
+ for (int i = curOrig; i >= 0; i--) delete originaux[i];
+ return DONE_NO_ACTION;
+ }
+ curOrig++;
+ }
+ }
+ // reverse if needed
+ // note that the selection list keeps its order
+ if ( reverseOrderForOp ) {
+ std::swap(originaux[0], originaux[1]);
+ std::swap(origWind[0], origWind[1]);
+ std::swap(origThresh[0], origThresh[1]);
+ }
+
+ // and work
+ // some temporary instances, first
+ Shape *theShapeA = new Shape;
+ Shape *theShapeB = new Shape;
+ Shape *theShape = new Shape;
+ Path *res = new Path;
+ res->SetBackData(false);
+ Path::cut_position *toCut=nullptr;
+ int nbToCut=0;
+
+ if ( bop == bool_op_inters || bop == bool_op_union || bop == bool_op_diff || bop == bool_op_symdiff ) {
+ // true boolean op
+ // get the polygons of each path, with the winding rule specified, and apply the operation iteratively
+ originaux[0]->ConvertWithBackData(origThresh[0]);
+
+ originaux[0]->Fill(theShape, 0);
+
+ theShapeA->ConvertToShape(theShape, origWind[0]);
+
+ curOrig = 1;
+ for (auto item : il){
+ if(item==il[0])continue;
+ originaux[curOrig]->ConvertWithBackData(origThresh[curOrig]);
+
+ originaux[curOrig]->Fill(theShape, curOrig);
+
+ theShapeB->ConvertToShape(theShape, origWind[curOrig]);
+
+ /* Due to quantization of the input shape coordinates, we may end up with A or B being empty.
+ * If this is a union or symdiff operation, we just use the non-empty shape as the result:
+ * A=0 => (0 or B) == B
+ * B=0 => (A or 0) == A
+ * A=0 => (0 xor B) == B
+ * B=0 => (A xor 0) == A
+ * If this is an intersection operation, we just use the empty shape as the result:
+ * A=0 => (0 and B) == 0 == A
+ * B=0 => (A and 0) == 0 == B
+ * If this a difference operation, and the upper shape (A) is empty, we keep B.
+ * If the lower shape (B) is empty, we still keep B, as it's empty:
+ * A=0 => (B - 0) == B
+ * B=0 => (0 - A) == 0 == B
+ *
+ * In any case, the output from this operation is stored in shape A, so we may apply
+ * the above rules simply by judicious use of swapping A and B where necessary.
+ */
+ bool zeroA = theShapeA->numberOfEdges() == 0;
+ bool zeroB = theShapeB->numberOfEdges() == 0;
+ if (zeroA || zeroB) {
+ // We might need to do a swap. Apply the above rules depending on operation type.
+ bool resultIsB = ((bop == bool_op_union || bop == bool_op_symdiff) && zeroA)
+ || ((bop == bool_op_inters) && zeroB)
+ || (bop == bool_op_diff);
+ if (resultIsB) {
+ // Swap A and B to use B as the result
+ Shape *swap = theShapeB;
+ theShapeB = theShapeA;
+ theShapeA = swap;
+ }
+ } else {
+ // Just do the Boolean operation as usual
+ // les elements arrivent en ordre inverse dans la liste
+ theShape->Booleen(theShapeB, theShapeA, bop);
+ Shape *swap = theShape;
+ theShape = theShapeA;
+ theShapeA = swap;
+ }
+ curOrig++;
+ }
+
+ {
+ Shape *swap = theShape;
+ theShape = theShapeA;
+ theShapeA = swap;
+ }
+
+ } else if ( bop == bool_op_cut ) {
+ // cuts= sort of a bastard boolean operation, thus not the axact same modus operandi
+ // technically, the cut path is not necessarily a polygon (thus has no winding rule)
+ // it is just uncrossed, and cleaned from duplicate edges and points
+ // then it's fed to Booleen() which will uncross it against the other path
+ // then comes the trick: each edge of the cut path is duplicated (one in each direction),
+ // thus making a polygon. the weight of the edges of the cut are all 0, but
+ // the Booleen need to invert the ones inside the source polygon (for the subsequent
+ // ConvertToForme)
+
+ // the cut path needs to have the highest pathID in the back data
+ // that's how the Booleen() function knows it's an edge of the cut
+ {
+ Path* swap=originaux[0];originaux[0]=originaux[1];originaux[1]=swap;
+ int swai=origWind[0];origWind[0]=origWind[1];origWind[1]=(fill_typ)swai;
+ std::swap(origThresh[0], origThresh[1]);
+ }
+ originaux[0]->ConvertWithBackData(origThresh[0]);
+
+ originaux[0]->Fill(theShape, 0);
+
+ theShapeA->ConvertToShape(theShape, origWind[0]);
+
+ originaux[1]->ConvertWithBackData(origThresh[1]);
+
+ if ((originaux[1]->pts.size() == 2) && originaux[1]->pts[0].isMoveTo && !originaux[1]->pts[1].isMoveTo)
+ originaux[1]->Fill(theShape, 1,false,true,false); // see LP Bug 177956
+ else
+ originaux[1]->Fill(theShape, 1,false,false,false); //do not closeIfNeeded
+
+ theShapeB->ConvertToShape(theShape, fill_justDont); // fill_justDont doesn't computes winding numbers
+
+ // les elements arrivent en ordre inverse dans la liste
+ theShape->Booleen(theShapeB, theShapeA, bool_op_cut, 1);
+
+ } else if ( bop == bool_op_slice ) {
+ // slice is not really a boolean operation
+ // you just put the 2 shapes in a single polygon, uncross it
+ // the points where the degree is > 2 are intersections
+ // just check it's an intersection on the path you want to cut, and keep it
+ // the intersections you have found are then fed to ConvertPositionsToMoveTo() which will
+ // make new subpath at each one of these positions
+ // inversion pour l'opération
+ {
+ Path* swap=originaux[0];originaux[0]=originaux[1];originaux[1]=swap;
+ int swai=origWind[0];origWind[0]=origWind[1];origWind[1]=(fill_typ)swai;
+ std::swap(origThresh[0], origThresh[1]);
+ }
+ originaux[0]->ConvertWithBackData(origThresh[0]);
+
+ originaux[0]->Fill(theShapeA, 0,false,false,false); // don't closeIfNeeded
+
+ originaux[1]->ConvertWithBackData(origThresh[1]);
+
+ originaux[1]->Fill(theShapeA, 1,true,false,false);// don't closeIfNeeded and just dump in the shape, don't reset it
+
+ theShape->ConvertToShape(theShapeA, fill_justDont);
+
+ if ( theShape->hasBackData() ) {
+ // should always be the case, but ya never know
+ {
+ for (int i = 0; i < theShape->numberOfPoints(); i++) {
+ if ( theShape->getPoint(i).totalDegree() > 2 ) {
+ // possibly an intersection
+ // we need to check that at least one edge from the source path is incident to it
+ // before we declare it's an intersection
+ int cb = theShape->getPoint(i).incidentEdge[FIRST];
+ int nbOrig=0;
+ int nbOther=0;
+ int piece=-1;
+ float t=0.0;
+ while ( cb >= 0 && cb < theShape->numberOfEdges() ) {
+ if ( theShape->ebData[cb].pathID == 0 ) {
+ // the source has an edge incident to the point, get its position on the path
+ piece=theShape->ebData[cb].pieceID;
+ if ( theShape->getEdge(cb).st == i ) {
+ t=theShape->ebData[cb].tSt;
+ } else {
+ t=theShape->ebData[cb].tEn;
+ }
+ nbOrig++;
+ }
+ if ( theShape->ebData[cb].pathID == 1 ) nbOther++; // the cut is incident to this point
+ cb=theShape->NextAt(i, cb);
+ }
+ if ( nbOrig > 0 && nbOther > 0 ) {
+ // point incident to both path and cut: an intersection
+ // note that you only keep one position on the source; you could have degenerate
+ // cases where the source crosses itself at this point, and you wouyld miss an intersection
+ toCut=(Path::cut_position*)realloc(toCut, (nbToCut+1)*sizeof(Path::cut_position));
+ toCut[nbToCut].piece=piece;
+ toCut[nbToCut].t=t;
+ nbToCut++;
+ }
+ }
+ }
+ }
+ {
+ // i think it's useless now
+ int i = theShape->numberOfEdges() - 1;
+ for (;i>=0;i--) {
+ if ( theShape->ebData[i].pathID == 1 ) {
+ theShape->SubEdge(i);
+ }
+ }
+ }
+
+ }
+ }
+
+ int* nesting=nullptr;
+ int* conts=nullptr;
+ int nbNest=0;
+ // pour compenser le swap juste avant
+ if ( bop == bool_op_slice ) {
+// theShape->ConvertToForme(res, nbOriginaux, originaux, true);
+// res->ConvertForcedToMoveTo();
+ res->Copy(originaux[0]);
+ res->ConvertPositionsToMoveTo(nbToCut, toCut); // cut where you found intersections
+ free(toCut);
+ } else if ( bop == bool_op_cut ) {
+ // il faut appeler pour desallouer PointData (pas vital, mais bon)
+ // the Booleen() function did not deallocate the point_data array in theShape, because this
+ // function needs it.
+ // this function uses the point_data to get the winding number of each path (ie: is a hole or not)
+ // for later reconstruction in objects, you also need to extract which path is parent of holes (nesting info)
+ theShape->ConvertToFormeNested(res, nbOriginaux, &originaux[0], 1, nbNest, nesting, conts);
+ } else {
+ theShape->ConvertToForme(res, nbOriginaux, &originaux[0]);
+ }
+
+ delete theShape;
+ delete theShapeA;
+ delete theShapeB;
+ for (int i = 0; i < nbOriginaux; i++) delete originaux[i];
+
+ if (res->descr_cmd.size() <= 1)
+ {
+ // only one command, presumably a moveto: it isn't a path
+ for (auto l : il){
+ l->deleteObject();
+ }
+ clear();
+
+ delete res;
+ return DONE_NO_PATH;
+ }
+
+ // get the source path object
+ SPObject *source;
+ if ( bop == bool_op_diff || bop == bool_op_cut || bop == bool_op_slice ) {
+ if (reverseOrderForOp) {
+ source = il[0];
+ } else {
+ source = il.back();
+ }
+ } else {
+ // find out the bottom object
+ std::vector<Inkscape::XML::Node*> sorted(xmlNodes().begin(), xmlNodes().end());
+
+ sort(sorted.begin(),sorted.end(),sp_repr_compare_position_bool);
+
+ source = doc->getObjectByRepr(sorted.front());
+ }
+
+ // adjust style properties that depend on a possible transform in the source object in order
+ // to get a correct style attribute for the new path
+ auto item_source = cast<SPItem>(source);
+ Geom::Affine i2doc(item_source->i2doc_affine());
+
+ Inkscape::XML::Node *repr_source = source->getRepr();
+
+ // remember important aspects of the source path, to be restored
+ gint pos = repr_source->position();
+ Inkscape::XML::Node *parent = repr_source->parent();
+ // remove source paths
+ clear();
+ for (auto l : il){
+ if (l != item_source) {
+ // delete the object for real, so that its clones can take appropriate action
+ l->deleteObject();
+ }
+ }
+
+ auto const source2doc_inverse = i2doc.inverse();
+ char const *const old_transform_attibute = repr_source->attribute("transform");
+
+ // now that we have the result, add it on the canvas
+ if ( bop == bool_op_cut || bop == bool_op_slice ) {
+ int nbRP=0;
+ Path** resPath;
+ if ( bop == bool_op_slice ) {
+ // there are moveto's at each intersection, but it's still one unique path
+ // so break it down and add each subpath independently
+ // we could call break_apart to do this, but while we have the description...
+ resPath=res->SubPaths(nbRP, false);
+ } else {
+ // cut operation is a bit wicked: you need to keep holes
+ // that's why you needed the nesting
+ // ConvertToFormeNested() dumped all the subpath in a single Path "res", so we need
+ // to get the path for each part of the polygon. that's why you need the nesting info:
+ // to know in which subpath to add a subpath
+ resPath=res->SubPathsWithNesting(nbRP, true, nbNest, nesting, conts);
+
+ // cleaning
+ if ( conts ) free(conts);
+ if ( nesting ) free(nesting);
+ }
+
+ // add all the pieces resulting from cut or slice
+ std::vector <Inkscape::XML::Node*> selection;
+ for (int i=0;i<nbRP;i++) {
+ transformLivarotPath(resPath[i], source2doc_inverse);
+ gchar *d = resPath[i]->svg_dump_path();
+
+ Inkscape::XML::Document *xml_doc = doc->getReprDoc();
+ Inkscape::XML::Node *repr = xml_doc->createElement("svg:path");
+
+ Inkscape::copy_object_properties(repr, repr_source);
+
+ // Delete source on last iteration (after we don't need repr_source anymore). As a consequence, the last
+ // item will inherit the original's id.
+ if (i + 1 == nbRP) {
+ item_source->deleteObject(false);
+ }
+
+ repr->setAttribute("d", d);
+ g_free(d);
+
+ // for slice, remove fill
+ if (bop == bool_op_slice) {
+ SPCSSAttr *css;
+
+ css = sp_repr_css_attr_new();
+ sp_repr_css_set_property(css, "fill", "none");
+
+ sp_repr_css_change(repr, css, "style");
+
+ sp_repr_css_attr_unref(css);
+ }
+
+ repr->setAttributeOrRemoveIfEmpty("transform", old_transform_attibute);
+
+ // add the new repr to the parent
+ // move to the saved position
+ parent->addChildAtPos(repr, pos);
+
+ selection.push_back(repr);
+ Inkscape::GC::release(repr);
+
+ delete resPath[i];
+ }
+ setReprList(selection);
+ if ( resPath ) free(resPath);
+
+ } else {
+ transformLivarotPath(res, source2doc_inverse);
+ gchar *d = res->svg_dump_path();
+
+ Inkscape::XML::Document *xml_doc = doc->getReprDoc();
+ Inkscape::XML::Node *repr = xml_doc->createElement("svg:path");
+
+ Inkscape::copy_object_properties(repr, repr_source);
+
+ // delete it so that its clones don't get alerted; this object will be restored shortly, with the same id
+ item_source->deleteObject(false);
+
+ repr->setAttribute("d", d);
+ g_free(d);
+
+ repr->setAttributeOrRemoveIfEmpty("transform", old_transform_attibute);
+
+ parent->addChildAtPos(repr, pos);
+
+ set(repr);
+ Inkscape::GC::release(repr);
+ }
+
+ delete res;
+
+ return DONE;
+}
+
+/*
+ Local Variables:
+ mode:c++
+ c-file-style:"stroustrup"
+ c-file-offsets:((innamespace . 0)(inline-open . 0)(case-label . +))
+ indent-tabs-mode:nil
+ fill-column:99
+ End:
+*/
+// vim: filetype=cpp:expandtab:shiftwidth=4:tabstop=8:softtabstop=4:fileencoding=utf-8:textwidth=99 :
diff --git a/src/path/path-boolop.h b/src/path/path-boolop.h
new file mode 100644
index 0000000..f587590
--- /dev/null
+++ b/src/path/path-boolop.h
@@ -0,0 +1,34 @@
+// SPDX-License-Identifier: GPL-2.0-or-later
+/** @file
+ * Boolean operations.
+ *//*
+ * Authors: see git history
+ *
+ * Copyright (C) 2018 Authors
+ * Released under GNU GPL v2+, read the file 'COPYING' for more information.
+ */
+#ifndef PATH_BOOLOP_H
+#define PATH_BOOLOP_H
+
+#include <2geom/path.h>
+#include "livarot/Path.h" // FillRule
+#include "object/object-set.h" // bool_op
+
+void sp_flatten(Geom::PathVector &pathvector, FillRule fillkind);
+Geom::PathVector sp_pathvector_boolop(Geom::PathVector const &pathva, Geom::PathVector const &pathvb, bool_op bop,
+ FillRule fra, FillRule frb, bool livarotonly, bool flattenbefore, int &error);
+Geom::PathVector sp_pathvector_boolop(Geom::PathVector const &pathva, Geom::PathVector const &pathvb, bool_op bop,
+ FillRule fra, FillRule frb, bool livarotonly = false, bool flattenbefore = true);
+
+#endif // PATH_BOOLOP_H
+
+/*
+ Local Variables:
+ mode:c++
+ c-file-style:"stroustrup"
+ c-file-offsets:((innamespace . 0)(inline-open . 0)(case-label . +))
+ indent-tabs-mode:nil
+ fill-column:99
+ End:
+*/
+// vim: filetype=cpp:expandtab:shiftwidth=4:tabstop=8:softtabstop=4:fileencoding=utf-8:textwidth=99 :
diff --git a/src/path/path-object-set.cpp b/src/path/path-object-set.cpp
new file mode 100644
index 0000000..f5d0770
--- /dev/null
+++ b/src/path/path-object-set.cpp
@@ -0,0 +1,166 @@
+// SPDX-License-Identifier: GPL-2.0-or-later
+
+/** @file
+ *
+ * Path related functions for ObjectSet.
+ *
+ * Copyright (C) 2020 Tavmjong Bah
+ *
+ * Released under GNU GPL v2+, read the file 'COPYING' for more information.
+ *
+ * TODO: Move related code from path-chemistry.cpp
+ *
+ */
+
+#include <glibmm/i18n.h>
+
+#include "document-undo.h"
+#include "message-stack.h"
+
+#include "attribute-rel-util.h"
+
+#include "object/object-set.h"
+#include "path/path-outline.h"
+#include "path/path-simplify.h"
+#include "ui/icon-names.h"
+
+using Inkscape::ObjectSet;
+
+bool
+ObjectSet::strokesToPaths(bool legacy, bool skip_undo)
+{
+ if (desktop() && isEmpty()) {
+ desktop()->messageStack()->flash(Inkscape::WARNING_MESSAGE, _("Select <b>stroked path(s)</b> to convert stroke to path."));
+ return false;
+ }
+
+ bool did = false;
+
+ Inkscape::Preferences *prefs = Inkscape::Preferences::get();
+ if (prefs->getBool("/options/pathoperationsunlink/value", true)) {
+ did = unlinkRecursive(true);
+ }
+
+ // Need to turn on stroke scaling to ensure stroke is scaled when transformed!
+ bool scale_stroke = prefs->getBool("/options/transform/stroke", true);
+ prefs->setBool("/options/transform/stroke", true);
+
+
+ std::vector<SPItem *> my_items(items().begin(), items().end());
+
+ for (auto item : my_items) {
+ // Do not remove the object from the selection here
+ // as we want to keep it selected if the whole operation fails
+ Inkscape::XML::Node *new_node = item_to_paths(item, legacy);
+ if (new_node) {
+ SPObject* new_item = document()->getObjectByRepr(new_node);
+
+ // Markers don't inherit properties from outside the
+ // marker. When converted to paths objects they need to be
+ // protected from inheritance. This is why (probably) the stroke
+ // to path code uses SP_STYLE_FLAG_ALWAYS when defining the
+ // style of the fill and stroke during the conversion. This
+ // means the style contains every possible property. Once we've
+ // finished the stroke to path conversion, we can eliminate
+ // unneeded properties from the style element.
+ sp_attribute_clean_recursive(new_node, SP_ATTRCLEAN_STYLE_REMOVE | SP_ATTRCLEAN_DEFAULT_REMOVE);
+
+ add(new_item); // Add to selection.
+ did = true;
+ }
+ }
+
+ // Reset
+ prefs->setBool("/options/transform/stroke", scale_stroke);
+
+ if (desktop() && !did) {
+ desktop()->messageStack()->flash(Inkscape::ERROR_MESSAGE, _("<b>No stroked paths</b> in the selection."));
+ }
+
+ if (did && !skip_undo) {
+ Inkscape::DocumentUndo::done(document(), _("Convert stroke to path"), "");
+ } else if (!did && !skip_undo) {
+ Inkscape::DocumentUndo::cancel(document());
+ }
+
+ return did;
+}
+
+bool
+ObjectSet::simplifyPaths(bool skip_undo)
+{
+ if (desktop() && isEmpty()) {
+ desktop()->messageStack()->flash(Inkscape::WARNING_MESSAGE, _("Select <b>path(s)</b> to simplify."));
+ return false;
+ }
+
+ Inkscape::Preferences *prefs = Inkscape::Preferences::get();
+ double threshold = prefs->getDouble("/options/simplifythreshold/value", 0.003);
+ bool justCoalesce = prefs->getBool( "/options/simplifyjustcoalesce/value", false);
+
+ // Keep track of accelerated simplify
+ static gint64 previous_time = 0;
+ static gdouble multiply = 1.0;
+
+ // Get the current time
+ gint64 current_time = g_get_monotonic_time();
+
+ // Was the previous call to this function recent? (<0.5 sec)
+ if (previous_time > 0 && current_time - previous_time < 500000) {
+
+ // add to the threshold 1/2 of its original value
+ multiply += 0.5;
+ threshold *= multiply;
+
+ } else {
+ // reset to the default
+ multiply = 1;
+ }
+
+ // Remember time for next call
+ previous_time = current_time;
+
+ // set "busy" cursor
+ if (desktop()) {
+ desktop()->setWaitingCursor();
+ }
+
+ Geom::OptRect selectionBbox = visualBounds();
+ if (!selectionBbox) {
+ std::cerr << "ObjectSet::: selection has no visual bounding box!" << std::endl;
+ return false;
+ }
+ double size = L2(selectionBbox->dimensions());
+
+ int pathsSimplified = 0;
+ std::vector<SPItem *> my_items(items().begin(), items().end());
+ for (auto item : my_items) {
+ pathsSimplified += path_simplify(item, threshold, justCoalesce, size);
+ }
+
+ if (pathsSimplified > 0 && !skip_undo) {
+ DocumentUndo::done(document(), _("Simplify"), INKSCAPE_ICON("path-simplify"));
+ }
+
+ if (desktop()) {
+ desktop()->clearWaitingCursor();
+ if (pathsSimplified > 0) {
+ desktop()->messageStack()->flashF(Inkscape::NORMAL_MESSAGE, _("<b>%d</b> paths simplified."), pathsSimplified);
+ } else {
+ desktop()->messageStack()->flash(Inkscape::ERROR_MESSAGE, _("<b>No paths</b> to simplify in the selection."));
+ }
+ }
+
+ return (pathsSimplified > 0);
+}
+
+/*
+ Local Variables:
+ mode:c++
+ c-file-style:"stroustrup"
+ c-file-offsets:((innamespace . 0)(inline-open . 0)(case-label . +))
+ indent-tabs-mode:nil
+ fill-column:99
+ End:
+*/
+// vim: filetype=cpp:expandtab:shiftwidth=4:tabstop=8:softtabstop=4:fileencoding=utf-8:textwidth=99 :
diff --git a/src/path/path-offset.cpp b/src/path/path-offset.cpp
new file mode 100644
index 0000000..f34c61a
--- /dev/null
+++ b/src/path/path-offset.cpp
@@ -0,0 +1,468 @@
+// SPDX-License-Identifier: GPL-2.0-or-later
+/** @file
+ * Path offsets.
+ *//*
+ * Authors:
+ * see git history
+ * Created by fred on Fri Dec 05 2003.
+ * tweaked endlessly by bulia byak <buliabyak@users.sf.net>
+ *
+ * Copyright (C) 2018 Authors
+ * Released under GNU GPL v2+, read the file 'COPYING' for more information.
+ */
+
+/*
+ * contains lots of stitched pieces of path-chemistry.c
+ */
+
+#include <vector>
+
+#include <glibmm/i18n.h>
+
+#include "path-offset.h"
+#include "path-util.h"
+
+#include "message-stack.h"
+#include "path-chemistry.h" // copy_object_properties()
+#include "selection.h"
+
+#include "display/curve.h"
+
+#include "livarot/Path.h"
+#include "livarot/Shape.h"
+
+#include "object/sp-flowtext.h"
+#include "object/sp-path.h"
+#include "object/sp-text.h"
+
+#include "ui/icon-names.h"
+
+#include "style.h"
+
+#define MIN_OFFSET 0.01
+
+using Inkscape::DocumentUndo;
+
+void sp_selected_path_do_offset(SPDesktop *desktop, bool expand, double prefOffset);
+void sp_selected_path_create_offset_object(SPDesktop *desktop, int expand, bool updating);
+
+void
+sp_selected_path_offset(SPDesktop *desktop)
+{
+ Inkscape::Preferences *prefs = Inkscape::Preferences::get();
+ double prefOffset = prefs->getDouble("/options/defaultoffsetwidth/value", 1.0, "px");
+
+ sp_selected_path_do_offset(desktop, true, prefOffset);
+}
+void
+sp_selected_path_inset(SPDesktop *desktop)
+{
+ Inkscape::Preferences *prefs = Inkscape::Preferences::get();
+ double prefOffset = prefs->getDouble("/options/defaultoffsetwidth/value", 1.0, "px");
+
+ sp_selected_path_do_offset(desktop, false, prefOffset);
+}
+
+void
+sp_selected_path_offset_screen(SPDesktop *desktop, double pixels)
+{
+ sp_selected_path_do_offset(desktop, true, pixels / desktop->current_zoom());
+}
+
+void
+sp_selected_path_inset_screen(SPDesktop *desktop, double pixels)
+{
+ sp_selected_path_do_offset(desktop, false, pixels / desktop->current_zoom());
+}
+
+
+void sp_selected_path_create_offset_object_zero(SPDesktop *desktop)
+{
+ sp_selected_path_create_offset_object(desktop, 0, false);
+}
+
+void sp_selected_path_create_offset(SPDesktop *desktop)
+{
+ sp_selected_path_create_offset_object(desktop, 1, false);
+}
+void sp_selected_path_create_inset(SPDesktop *desktop)
+{
+ sp_selected_path_create_offset_object(desktop, -1, false);
+}
+
+void sp_selected_path_create_updating_offset_object_zero(SPDesktop *desktop)
+{
+ sp_selected_path_create_offset_object(desktop, 0, true);
+}
+
+void sp_selected_path_create_updating_offset(SPDesktop *desktop)
+{
+ sp_selected_path_create_offset_object(desktop, 1, true);
+}
+
+void sp_selected_path_create_updating_inset(SPDesktop *desktop)
+{
+ sp_selected_path_create_offset_object(desktop, -1, true);
+}
+
+void sp_selected_path_create_offset_object(SPDesktop *desktop, int expand, bool updating)
+{
+ Inkscape::Selection *selection = desktop->getSelection();
+ SPItem *item = selection->singleItem();
+
+ if (auto shape = cast<SPShape>(item)) {
+ if (!shape->curve())
+ return;
+ } else if (is<SPText>(item)) {
+ } else {
+ desktop->messageStack()->flash(Inkscape::ERROR_MESSAGE, _("Selected object is <b>not a path</b>, cannot inset/outset."));
+ return;
+ }
+
+ Geom::Affine const transform(item->transform);
+ auto scaling_factor = item->i2doc_affine().descrim();
+
+ item->doWriteTransform(Geom::identity());
+
+ // remember the position of the item
+ gint pos = item->getRepr()->position();
+
+ // remember parent
+ Inkscape::XML::Node *parent = item->getRepr()->parent();
+
+ float o_width = 0;
+ {
+ Inkscape::Preferences *prefs = Inkscape::Preferences::get();
+ o_width = prefs->getDouble("/options/defaultoffsetwidth/value", 1.0, "px");
+ o_width /= scaling_factor;
+ if (scaling_factor == 0 || o_width < MIN_OFFSET) {
+ o_width = MIN_OFFSET;
+ }
+ }
+
+ Path *orig = Path_for_item(item, true, false);
+ if (orig == nullptr)
+ {
+ return;
+ }
+
+ Path *res = new Path;
+ res->SetBackData(false);
+
+ {
+ Shape *theShape = new Shape;
+ Shape *theRes = new Shape;
+
+ orig->ConvertWithBackData(1.0);
+ orig->Fill(theShape, 0);
+
+ SPCSSAttr *css = sp_repr_css_attr(item->getRepr(), "style");
+ gchar const *val = sp_repr_css_property(css, "fill-rule", nullptr);
+ if (val && strcmp(val, "nonzero") == 0)
+ {
+ theRes->ConvertToShape(theShape, fill_nonZero);
+ }
+ else if (val && strcmp(val, "evenodd") == 0)
+ {
+ theRes->ConvertToShape(theShape, fill_oddEven);
+ }
+ else
+ {
+ theRes->ConvertToShape(theShape, fill_nonZero);
+ }
+
+ Path *originaux[1];
+ originaux[0] = orig;
+ theRes->ConvertToForme(res, 1, originaux);
+
+ delete theShape;
+ delete theRes;
+ }
+
+ if (res->descr_cmd.size() <= 1)
+ {
+ // pas vraiment de points sur le resultat
+ // donc il ne reste rien
+ DocumentUndo::done(desktop->getDocument(),
+ (updating ? _("Create linked offset")
+ : _("Create dynamic offset")),
+ (updating ? INKSCAPE_ICON("path-offset-linked")
+ : INKSCAPE_ICON("path-offset-dynamic")));
+ selection->clear();
+
+ delete res;
+ delete orig;
+ return;
+ }
+
+ {
+ Inkscape::XML::Document *xml_doc = desktop->doc()->getReprDoc();
+ Inkscape::XML::Node *repr = xml_doc->createElement("svg:path");
+
+ if (!updating) {
+ Inkscape::copy_object_properties(repr, item->getRepr());
+ } else {
+ gchar const *style = item->getRepr()->attribute("style");
+ repr->setAttribute("style", style);
+ }
+
+ repr->setAttribute("sodipodi:type", "inkscape:offset");
+ repr->setAttributeSvgDouble("inkscape:radius", ( expand > 0
+ ? o_width
+ : expand < 0
+ ? -o_width
+ : 0 ));
+
+ gchar *str = res->svg_dump_path();
+ repr->setAttribute("inkscape:original", str);
+ g_free(str);
+ str = nullptr;
+
+ if ( updating ) {
+
+ //XML Tree being used directly here while it shouldn't be
+ item->doWriteTransform(transform);
+ char const *id = item->getRepr()->attribute("id");
+ char const *uri = g_strdup_printf("#%s", id);
+ repr->setAttribute("xlink:href", uri);
+ g_free((void *) uri);
+ } else {
+ repr->removeAttribute("inkscape:href");
+ // delete original
+ item->deleteObject(false);
+ }
+
+ // add the new repr to the parent
+ // move to the saved position
+ parent->addChildAtPos(repr, pos);
+
+ SPItem *nitem = reinterpret_cast<SPItem *>(desktop->getDocument()->getObjectByRepr(repr));
+
+ if ( !updating ) {
+ // apply the transform to the offset
+ nitem->doWriteTransform(transform);
+ }
+
+ // The object just created from a temporary repr is only a seed.
+ // We need to invoke its write which will update its real repr (in particular adding d=)
+ nitem->updateRepr();
+
+ Inkscape::GC::release(repr);
+
+ selection->set(nitem);
+ }
+
+ DocumentUndo::done(desktop->getDocument(),
+ (updating ? _("Create linked offset")
+ : _("Create dynamic offset")),
+ (updating ? INKSCAPE_ICON("path-offset-linked")
+ : INKSCAPE_ICON("path-offset-dynamic")));
+
+ delete res;
+ delete orig;
+}
+
+/**
+ * Apply offset to selected paths
+ * @param desktop Targeted desktop
+ * @param expand True if offset expands, False if it shrinks paths
+ * @param prefOffset Size of offset in pixels
+ */
+void
+sp_selected_path_do_offset(SPDesktop *desktop, bool expand, double prefOffset)
+{
+ Inkscape::Selection *selection = desktop->getSelection();
+
+ if (selection->isEmpty()) {
+ desktop->messageStack()->flash(Inkscape::WARNING_MESSAGE, _("Select <b>path(s)</b> to inset/outset."));
+ return;
+ }
+
+ bool did = false;
+ std::vector<SPItem*> il(selection->items().begin(), selection->items().end());
+ for (auto item : il){
+ if (auto shape = cast<SPShape>(item)) {
+ if (!shape->curve())
+ continue;
+ } else if (is<SPText>(item) || is<SPFlowtext>(item)) {
+ } else {
+ continue;
+ }
+
+ Geom::Affine const transform(item->transform);
+ auto scaling_factor = item->i2doc_affine().descrim();
+
+ item->doWriteTransform(Geom::identity());
+
+ float o_width = 0;
+ float o_miter = 0;
+ JoinType o_join = join_straight;
+ //ButtType o_butt = butt_straight;
+
+ {
+ SPStyle *i_style = item->style;
+ int jointype = i_style->stroke_linejoin.value;
+
+ switch (jointype) {
+ case SP_STROKE_LINEJOIN_MITER:
+ o_join = join_pointy;
+ break;
+ case SP_STROKE_LINEJOIN_ROUND:
+ o_join = join_round;
+ break;
+ default:
+ o_join = join_straight;
+ break;
+ }
+
+ // scale to account for transforms and document units
+ o_width = prefOffset / scaling_factor;
+
+ if (scaling_factor == 0 || o_width < MIN_OFFSET) {
+ o_width = MIN_OFFSET;
+ }
+ o_miter = i_style->stroke_miterlimit.value * o_width;
+ }
+
+ Path *orig = Path_for_item(item, false);
+ if (orig == nullptr) {
+ continue;
+ }
+
+ Path *res = new Path;
+ res->SetBackData(false);
+
+ {
+ Shape *theShape = new Shape;
+ Shape *theRes = new Shape;
+
+ orig->ConvertWithBackData(0.03);
+ orig->Fill(theShape, 0);
+
+ SPCSSAttr *css = sp_repr_css_attr(item->getRepr(), "style");
+ gchar const *val = sp_repr_css_property(css, "fill-rule", nullptr);
+ if (val && strcmp(val, "nonzero") == 0)
+ {
+ theRes->ConvertToShape(theShape, fill_nonZero);
+ }
+ else if (val && strcmp(val, "evenodd") == 0)
+ {
+ theRes->ConvertToShape(theShape, fill_oddEven);
+ }
+ else
+ {
+ theRes->ConvertToShape(theShape, fill_nonZero);
+ }
+
+ // et maintenant: offset
+ // methode inexacte
+/* Path *originaux[1];
+ originaux[0] = orig;
+ theRes->ConvertToForme(res, 1, originaux);
+
+ if (expand) {
+ res->OutsideOutline(orig, 0.5 * o_width, o_join, o_butt, o_miter);
+ } else {
+ res->OutsideOutline(orig, -0.5 * o_width, o_join, o_butt, o_miter);
+ }
+
+ orig->ConvertWithBackData(1.0);
+ orig->Fill(theShape, 0);
+ theRes->ConvertToShape(theShape, fill_positive);
+ originaux[0] = orig;
+ theRes->ConvertToForme(res, 1, originaux);
+
+ if (o_width >= 0.5) {
+ // res->Coalesce(1.0);
+ res->ConvertEvenLines(1.0);
+ res->Simplify(1.0);
+ } else {
+ // res->Coalesce(o_width);
+ res->ConvertEvenLines(1.0*o_width);
+ res->Simplify(1.0 * o_width);
+ } */
+ // methode par makeoffset
+
+ if (expand)
+ {
+ theShape->MakeOffset(theRes, o_width, o_join, o_miter);
+ }
+ else
+ {
+ theShape->MakeOffset(theRes, -o_width, o_join, o_miter);
+ }
+ theRes->ConvertToShape(theShape, fill_positive);
+
+ res->Reset();
+ theRes->ConvertToForme(res);
+
+ res->ConvertEvenLines(0.1);
+ res->Simplify(0.1);
+
+ delete theShape;
+ delete theRes;
+ }
+
+ did = true;
+
+ // remember the position of the item
+ gint pos = item->getRepr()->position();
+ // remember parent
+ Inkscape::XML::Node *parent = item->getRepr()->parent();
+
+ selection->remove(item);
+
+ Inkscape::XML::Node *repr = nullptr;
+
+ if (res->descr_cmd.size() > 1) { // if there's 0 or 1 node left, drop this path altogether
+ Inkscape::XML::Document *xml_doc = desktop->doc()->getReprDoc();
+ repr = xml_doc->createElement("svg:path");
+
+ Inkscape::copy_object_properties(repr, item->getRepr());
+ }
+
+ item->deleteObject(false);
+
+ if (repr) {
+ gchar *str = res->svg_dump_path();
+ repr->setAttribute("d", str);
+ g_free(str);
+
+ // add the new repr to the parent
+ // move to the saved position
+ parent->addChildAtPos(repr, pos);
+
+ SPItem *newitem = (SPItem *) desktop->getDocument()->getObjectByRepr(repr);
+
+ // reapply the transform
+ newitem->doWriteTransform(transform);
+
+ selection->add(repr);
+
+ Inkscape::GC::release(repr);
+ }
+
+ delete orig;
+ delete res;
+ }
+
+ if (did) {
+ DocumentUndo::done(desktop->getDocument(),
+ (expand ? _("Outset path") : _("Inset path")),
+ (expand ? INKSCAPE_ICON("path-outset") : INKSCAPE_ICON("path-inset")));
+ } else {
+ desktop->messageStack()->flash(Inkscape::ERROR_MESSAGE, _("<b>No paths</b> to inset/outset in the selection."));
+ return;
+ }
+}
+
+/*
+ Local Variables:
+ mode:c++
+ c-file-style:"stroustrup"
+ c-file-offsets:((innamespace . 0)(inline-open . 0)(case-label . +))
+ indent-tabs-mode:nil
+ fill-column:99
+ End:
+*/
+// vim: filetype=cpp:expandtab:shiftwidth=4:tabstop=8:softtabstop=4:fileencoding=utf-8:textwidth=99 :
diff --git a/src/path/path-offset.h b/src/path/path-offset.h
new file mode 100644
index 0000000..d4996bc
--- /dev/null
+++ b/src/path/path-offset.h
@@ -0,0 +1,41 @@
+// SPDX-License-Identifier: GPL-2.0-or-later
+/** @file
+ * Path offsets.
+ *//*
+ * Authors: see git history
+ *
+ * Copyright (C) 2018 Authors
+ * Released under GNU GPL v2+, read the file 'COPYING' for more information.
+ */
+#ifndef PATH_OFFSET_H
+#define PATH_OFFSET_H
+
+class SPDesktop;
+
+// offset/inset of a curve
+// takes the fill-rule in consideration
+// offset amount is the stroke-width of the curve
+void sp_selected_path_offset (SPDesktop *desktop);
+void sp_selected_path_offset_screen (SPDesktop *desktop, double pixels);
+void sp_selected_path_inset (SPDesktop *desktop);
+void sp_selected_path_inset_screen (SPDesktop *desktop, double pixels);
+void sp_selected_path_create_offset (SPDesktop *desktop);
+void sp_selected_path_create_inset (SPDesktop *desktop);
+void sp_selected_path_create_updating_offset (SPDesktop *desktop);
+void sp_selected_path_create_updating_inset (SPDesktop *desktop);
+
+void sp_selected_path_create_offset_object_zero (SPDesktop *desktop);
+void sp_selected_path_create_updating_offset_object_zero (SPDesktop *desktop);
+
+#endif // PATH_OFFSET_H
+
+/*
+ Local Variables:
+ mode:c++
+ c-file-style:"stroustrup"
+ c-file-offsets:((innamespace . 0)(inline-open . 0)(case-label . +))
+ indent-tabs-mode:nil
+ fill-column:99
+ End:
+*/
+// vim: filetype=cpp:expandtab:shiftwidth=4:tabstop=8:softtabstop=4:fileencoding=utf-8:textwidth=99 :
diff --git a/src/path/path-outline.cpp b/src/path/path-outline.cpp
new file mode 100644
index 0000000..47c3dee
--- /dev/null
+++ b/src/path/path-outline.cpp
@@ -0,0 +1,792 @@
+// SPDX-License-Identifier: GPL-2.0-or-later
+
+/** @file
+ *
+ * Two related object to path operations:
+ *
+ * 1. Find a path that includes fill, stroke, and markers. Useful for finding a visual bounding box.
+ * 2. Take a set of objects and find an identical visual representation using only paths.
+ *
+ * Copyright (C) 2020 Tavmjong Bah
+ * Copyright (C) 2018 Authors
+ *
+ * Released under GNU GPL v2+, read the file 'COPYING' for more information.
+ *
+ * Code moved from splivarot.cpp
+ *
+ */
+
+#include "path-outline.h"
+
+#include <vector>
+
+#include "path-chemistry.h" // Should be moved to path directory
+#include "message-stack.h" // Should be removed.
+#include "selection.h"
+#include "style.h"
+
+#include "display/curve.h" // Should be moved to path directory
+
+#include "helper/geom.h" // pathv_to_linear_and_cubic()
+
+#include "livarot/LivarotDefs.h"
+#include "livarot/Path.h"
+#include "livarot/Shape.h"
+
+#include "object/object-set.h"
+#include "object/box3d.h"
+#include "object/sp-item.h"
+#include "object/sp-marker.h"
+#include "object/sp-shape.h"
+#include "object/sp-text.h"
+#include "object/sp-flowtext.h"
+
+#include "svg/svg.h"
+
+/**
+ * Given an item, find a path representing the fill and a path representing the stroke.
+ * Returns true if fill path found. Item may not have a stroke in which case stroke path is empty.
+ * bbox_only==true skips cleaning up the stroke path.
+ * Encapsulates use of livarot.
+ */
+bool
+item_find_paths(const SPItem *item, Geom::PathVector& fill, Geom::PathVector& stroke, bool bbox_only = false)
+{
+ auto shape = cast<SPShape>(item);
+ auto text = cast<SPText>(item);
+
+ if (!shape && !text) {
+ return false;
+ }
+
+ std::optional<SPCurve> curve;
+ if (shape) {
+ curve = SPCurve::ptr_to_opt(shape->curve());
+ } else if (text) {
+ curve = text->getNormalizedBpath();
+ } else {
+ std::cerr << "item_find_paths: item not shape or text!" << std::endl;
+ return false;
+ }
+
+ if (!curve) {
+ std::cerr << "item_find_paths: no curve!" << std::endl;
+ return false;
+ }
+
+ if (curve->get_pathvector().empty()) {
+ std::cerr << "item_find_paths: curve empty!" << std::endl;
+ return false;
+ }
+
+ fill = curve->get_pathvector();
+
+ if (!item->style) {
+ // Should never happen
+ std::cerr << "item_find_paths: item with no style!" << std::endl;
+ return false;
+ }
+
+ if (item->style->stroke.isNone()) {
+ // No stroke, no chocolate!
+ return true;
+ }
+
+ // Now that we have a valid curve with stroke, do offset. We use Livarot for this as
+ // lib2geom does not yet handle offsets correctly.
+
+ // Livarot's outline of arcs is broken. So convert the path to linear and cubics only, for
+ // which the outline is created correctly.
+ Geom::PathVector pathv = pathv_to_linear_and_cubic_beziers( fill );
+
+ SPStyle *style = item->style;
+
+ double stroke_width = style->stroke_width.computed;
+ if (stroke_width < Geom::EPSILON) {
+ // https://bugs.launchpad.net/inkscape/+bug/1244861
+ stroke_width = Geom::EPSILON;
+ }
+ double miter = style->stroke_miterlimit.value * stroke_width;
+
+ JoinType join;
+ switch (style->stroke_linejoin.computed) {
+ case SP_STROKE_LINEJOIN_MITER:
+ join = join_pointy;
+ break;
+ case SP_STROKE_LINEJOIN_ROUND:
+ join = join_round;
+ break;
+ default:
+ join = join_straight;
+ break;
+ }
+
+ ButtType butt;
+ switch (style->stroke_linecap.computed) {
+ case SP_STROKE_LINECAP_SQUARE:
+ butt = butt_square;
+ break;
+ case SP_STROKE_LINECAP_ROUND:
+ butt = butt_round;
+ break;
+ default:
+ butt = butt_straight;
+ break;
+ }
+
+ Path *origin = new Path; // Fill
+ Path *offset = new Path;
+
+ Geom::Affine const transform(item->transform);
+ double const scale = transform.descrim();
+
+ origin->LoadPathVector(pathv);
+ offset->SetBackData(false);
+
+ if (!style->stroke_dasharray.values.empty() && style->stroke_dasharray.is_valid()) {
+ // We have dashes!
+ origin->ConvertWithBackData(0.005); // Approximate by polyline
+ origin->DashPolylineFromStyle(style, scale, 0);
+ auto bounds = Geom::bounds_fast(pathv);
+ if (bounds) {
+ double size = Geom::L2(bounds->dimensions());
+ origin->Simplify(size * 0.000005); // Polylines to Beziers
+ }
+ }
+
+ // Finally do offset!
+ origin->Outline(offset, 0.5 * stroke_width, join, butt, 0.5 * miter);
+
+ if (bbox_only) {
+ stroke = offset->MakePathVector();
+ } else {
+ // Clean-up shape
+
+ offset->ConvertWithBackData(1.0); // Approximate by polyline
+
+ Shape *theShape = new Shape;
+ offset->Fill(theShape, 0); // Convert polyline to shape, step 1.
+
+ Shape *theOffset = new Shape;
+ theOffset->ConvertToShape(theShape, fill_positive); // Create an intersection free polygon (theOffset), step2.
+ theOffset->ConvertToForme(origin, 1, &offset); // Turn shape into contour (stored in origin).
+
+ stroke = origin->MakePathVector(); // Note origin was replaced above by stroke!
+ }
+
+ delete origin;
+ delete offset;
+
+ // std::cout << " fill: " << sp_svg_write_path(fill) << " count: " << fill.curveCount() << std::endl;
+ // std::cout << " stroke: " << sp_svg_write_path(stroke) << " count: " << stroke.curveCount() << std::endl;
+ return true;
+}
+
+
+// ======================== Item to Outline ===================== //
+
+static
+void item_to_outline_add_marker_child( SPItem const *item, Geom::Affine marker_transform, Geom::PathVector* pathv_in )
+{
+ Geom::Affine tr(marker_transform);
+ tr = item->transform * tr;
+
+ // note: a marker child item can be an item group!
+ if (is<SPGroup>(item)) {
+ // recurse through all childs:
+ for (auto& o: item->children) {
+ if (auto childitem = cast<SPItem>(&o)) {
+ item_to_outline_add_marker_child(childitem, tr, pathv_in);
+ }
+ }
+ } else {
+ Geom::PathVector* marker_pathv = item_to_outline(item);
+
+ if (marker_pathv) {
+ for (const auto & j : *marker_pathv) {
+ pathv_in->push_back(j * tr);
+ }
+ delete marker_pathv;
+ }
+ }
+}
+
+static
+void item_to_outline_add_marker( SPObject const *marker_object, Geom::Affine marker_transform,
+ Geom::Scale stroke_scale, Geom::PathVector* pathv_in )
+{
+ SPMarker const * marker = cast<SPMarker>(marker_object);
+
+ Geom::Affine tr(marker_transform);
+ if (marker->markerUnits == SP_MARKER_UNITS_STROKEWIDTH) {
+ tr = stroke_scale * tr;
+ }
+ // total marker transform
+ tr = marker->c2p * tr;
+
+ SPItem const * marker_item = sp_item_first_item_child(marker_object); // why only consider the first item? can a marker only consist of a single item (that may be a group)?
+ if (marker_item) {
+ item_to_outline_add_marker_child(marker_item, tr, pathv_in);
+ }
+}
+
+
+/**
+ * Returns a pathvector that is the outline of the stroked item, with markers.
+ * item must be an SPShape or an SPText.
+ * The only current use of this function has exclude_markers true! (SPShape::either_bbox).
+ * TODO: See if SPShape::either_bbox's union with markers is the same as one would get
+ * with bbox_only false.
+ */
+Geom::PathVector* item_to_outline(SPItem const *item, bool exclude_markers)
+{
+ Geom::PathVector fill; // Used for locating markers.
+ Geom::PathVector stroke; // Used for creating outline (and finding bbox).
+ item_find_paths(item, fill, stroke, true); // Skip cleaning up stroke shape.
+
+ Geom::PathVector *ret_pathv = nullptr;
+
+ if (fill.curveCount() == 0) {
+ std::cerr << "item_to_outline: fill path has no segments!" << std::endl;
+ return ret_pathv;
+ }
+
+ if (stroke.size() > 0) {
+ ret_pathv = new Geom::PathVector(stroke);
+ } else {
+ // No stroke, use fill path.
+ ret_pathv = new Geom::PathVector(fill);
+ }
+
+ if (exclude_markers) {
+ return ret_pathv;
+ }
+
+ auto shape = cast<SPShape>(item);
+ if (shape && shape->hasMarkers()) {
+
+ SPStyle *style = shape->style;
+ Geom::Scale scale(style->stroke_width.computed);
+
+ // START marker
+ for (int i = 0; i < 2; i++) { // SP_MARKER_LOC and SP_MARKER_LOC_START
+ if ( SPObject *marker_obj = shape->_marker[i] ) {
+ Geom::Affine const m (sp_shape_marker_get_transform_at_start(fill.front().front()));
+ item_to_outline_add_marker( marker_obj, m, scale, ret_pathv );
+ }
+ }
+
+ // MID marker
+ for (int i = 0; i < 3; i += 2) { // SP_MARKER_LOC and SP_MARKER_LOC_MID
+ SPObject *midmarker_obj = shape->_marker[i];
+ if (!midmarker_obj) continue;
+ for(Geom::PathVector::const_iterator path_it = fill.begin(); path_it != fill.end(); ++path_it) {
+
+ // START position
+ if ( path_it != fill.begin() &&
+ ! ((path_it == (fill.end()-1)) && (path_it->size_default() == 0)) ) // if this is the last path and it is a moveto-only, there is no mid marker there
+ {
+ Geom::Affine const m (sp_shape_marker_get_transform_at_start(path_it->front()));
+ item_to_outline_add_marker( midmarker_obj, m, scale, ret_pathv);
+ }
+
+ // MID position
+ if (path_it->size_default() > 1) {
+ Geom::Path::const_iterator curve_it1 = path_it->begin(); // incoming curve
+ Geom::Path::const_iterator curve_it2 = ++(path_it->begin()); // outgoing curve
+ while (curve_it2 != path_it->end_default())
+ {
+ /* Put marker between curve_it1 and curve_it2.
+ * Loop to end_default (so including closing segment), because when a path is closed,
+ * there should be a midpoint marker between last segment and closing straight line segment
+ */
+ Geom::Affine const m (sp_shape_marker_get_transform(*curve_it1, *curve_it2));
+ item_to_outline_add_marker( midmarker_obj, m, scale, ret_pathv);
+
+ ++curve_it1;
+ ++curve_it2;
+ }
+ }
+
+ // END position
+ if ( path_it != (fill.end()-1) && !path_it->empty()) {
+ Geom::Curve const &lastcurve = path_it->back_default();
+ Geom::Affine const m = sp_shape_marker_get_transform_at_end(lastcurve);
+ item_to_outline_add_marker( midmarker_obj, m, scale, ret_pathv );
+ }
+ }
+ }
+
+ // END marker
+ for (int i = 0; i < 4; i += 3) { // SP_MARKER_LOC and SP_MARKER_LOC_END
+ if ( SPObject *marker_obj = shape->_marker[i] ) {
+ /* Get reference to last curve in the path.
+ * For moveto-only path, this returns the "closing line segment". */
+ Geom::Path const &path_last = fill.back();
+ unsigned int index = path_last.size_default();
+ if (index > 0) {
+ index--;
+ }
+ Geom::Curve const &lastcurve = path_last[index];
+
+ Geom::Affine const m = sp_shape_marker_get_transform_at_end(lastcurve);
+ item_to_outline_add_marker( marker_obj, m, scale, ret_pathv );
+ }
+ }
+ }
+
+ return ret_pathv;
+}
+
+
+
+// ========================= Stroke to Path ====================== //
+
+static
+void item_to_paths_add_marker( SPItem *context,
+ SPObject *marker_object,
+ Geom::Affine marker_transform,
+ double linewidth,
+ bool start_marker,
+ Inkscape::XML::Node *g_repr,
+ bool legacy)
+{
+ auto doc = context->document;
+ auto marker = cast<SPMarker>(marker_object);
+ SPItem* marker_item = sp_item_first_item_child(marker_object);
+ if (!marker_item) {
+ return;
+ }
+
+
+ // total marker transform
+ auto tr = marker->get_marker_transform(marker_transform, linewidth, start_marker);
+ tr = marker_item->transform * marker->c2p * tr;
+
+ if (marker_item->getRepr()) {
+ Inkscape::XML::Node *m_repr = marker_item->getRepr()->duplicate(doc->getReprDoc());
+ g_repr->addChildAtPos(m_repr, 0);
+ SPItem *marker_item = (SPItem *) doc->getObjectByRepr(m_repr);
+ marker_item->doWriteTransform(tr);
+ if (!legacy) {
+ item_to_paths(marker_item, legacy, context);
+ }
+ }
+}
+
+
+/*
+ * Find an outline that represents an item.
+ * If legacy, text will not be handled as it is not a shape.
+ * If a new item is created it is returned.
+ * If the input item is a group and that group contains a changed item, the group node is returned
+ * (marking a change).
+ *
+ * The return value is used externally to update a selection. It is nullptr if no change is made.
+ */
+Inkscape::XML::Node*
+item_to_paths(SPItem *item, bool legacy, SPItem *context)
+{
+ char const *id = item->getAttribute("id");
+ SPDocument *doc = item->document;
+ bool flatten = false;
+ // flatten all paths effects
+ auto lpeitem = cast<SPLPEItem>(item);
+ if (lpeitem && lpeitem->hasPathEffect()) {
+ lpeitem->removeAllPathEffects(true);
+ SPObject *elemref = doc->getObjectById(id);
+ if (elemref && elemref != item) {
+ // If the LPE item is a shape, it is converted to a path
+ // so we need to reupdate the item
+ item = cast<SPItem>(elemref);
+ }
+ auto flat_item = cast<SPLPEItem>(elemref);
+ if (!flat_item || !flat_item->hasPathEffect()) {
+ flatten = true;
+ }
+ }
+ // convert text/3dbox to path
+ if (is<SPText>(item) || is<SPFlowtext>(item) || is<SPBox3D>(item)) {
+ if (legacy) {
+ return nullptr;
+ }
+
+ Inkscape::ObjectSet original_objects {doc}; // doc or desktop shouldn't be necessary
+ original_objects.add(item);
+ original_objects.toCurves(true);
+ SPItem * new_item = original_objects.singleItem();
+ if (new_item && new_item != item) {
+ flatten = true;
+ item = new_item;
+ } else {
+ g_warning("item_to_paths: flattening text or 3D box failed.");
+ return nullptr;
+ }
+ }
+ // if group, recurse
+ auto group = cast<SPGroup>(item);
+ if (group) {
+ if (legacy) {
+ return nullptr;
+ }
+ std::vector<SPItem*> const item_list = group->item_list();
+ bool did = false;
+ for (auto subitem : item_list) {
+ if (item_to_paths(subitem, legacy)) {
+ did = true;
+ }
+ }
+ if (did || flatten) {
+ // This indicates that at least one thing was changed inside the group.
+ return group->getRepr();
+ } else {
+ return nullptr;
+ }
+ }
+
+ auto shape = cast<SPShape>(item);
+ if (!shape) {
+ return nullptr;
+ }
+
+ Geom::PathVector fill_path;
+ Geom::PathVector stroke_path;
+ bool status = item_find_paths(item, fill_path, stroke_path);
+
+ if (!status) {
+ // Was not a well structured shape (or text).
+ return nullptr;
+ }
+
+ // The styles ------------------------
+
+ // Copying stroke style to fill will fail for properties not defined by style attribute
+ // (i.e., properties defined in style sheet or by attributes).
+ SPStyle *style = item->style;
+ SPCSSAttr *ncss = sp_css_attr_from_style(style, SP_STYLE_FLAG_ALWAYS);
+ SPCSSAttr *ncsf = sp_css_attr_from_style(style, SP_STYLE_FLAG_ALWAYS);
+
+ if (context) {
+ SPCSSAttr *ctxt_style = sp_css_attr_from_style(context->style, SP_STYLE_FLAG_ALWAYS);
+
+ // TODO: browsers have different behaviours with context on markers
+ // we need to revisit in the future for best matching
+ // also dont know if opacity is or should be included in context
+ gchar const *s_val = sp_repr_css_property(ctxt_style, "stroke", nullptr);
+ gchar const *f_val = sp_repr_css_property(ctxt_style, "fill", nullptr);
+ if (style->fill.paintOrigin == SP_CSS_PAINT_ORIGIN_CONTEXT_STROKE ||
+ style->fill.paintOrigin == SP_CSS_PAINT_ORIGIN_CONTEXT_FILL)
+ {
+ gchar const *fill_value = (style->fill.paintOrigin == SP_CSS_PAINT_ORIGIN_CONTEXT_STROKE) ? s_val : f_val;
+ sp_repr_css_set_property(ncss, "fill", fill_value);
+ sp_repr_css_set_property(ncsf, "fill", fill_value);
+ }
+ if (style->stroke.paintOrigin == SP_CSS_PAINT_ORIGIN_CONTEXT_STROKE ||
+ style->stroke.paintOrigin == SP_CSS_PAINT_ORIGIN_CONTEXT_FILL)
+ {
+ gchar const *stroke_value = (style->stroke.paintOrigin == SP_CSS_PAINT_ORIGIN_CONTEXT_FILL) ? f_val : s_val;
+ sp_repr_css_set_property(ncss, "stroke", stroke_value);
+ sp_repr_css_set_property(ncsf, "stroke", stroke_value);
+ }
+ }
+ // Stroke
+
+ gchar const *s_val = sp_repr_css_property(ncss, "stroke", nullptr);
+ gchar const *s_opac = sp_repr_css_property(ncss, "stroke-opacity", nullptr);
+ gchar const *f_val = sp_repr_css_property(ncss, "fill", nullptr);
+ gchar const *opacity = sp_repr_css_property(ncss, "opacity", nullptr); // Also for markers
+ gchar const *filter = sp_repr_css_property(ncss, "filter", nullptr); // Also for markers
+
+ sp_repr_css_set_property(ncss, "stroke", "none");
+ sp_repr_css_set_property(ncss, "stroke-width", nullptr);
+ sp_repr_css_set_property(ncss, "stroke-opacity", "1.0");
+ sp_repr_css_set_property(ncss, "filter", nullptr);
+ sp_repr_css_set_property(ncss, "opacity", nullptr);
+ sp_repr_css_unset_property(ncss, "marker-start");
+ sp_repr_css_unset_property(ncss, "marker-mid");
+ sp_repr_css_unset_property(ncss, "marker-end");
+
+ // we change the stroke to fill on ncss to create the filled stroke
+ sp_repr_css_set_property(ncss, "fill", s_val);
+ if ( s_opac ) {
+ sp_repr_css_set_property(ncss, "fill-opacity", s_opac);
+ } else {
+ sp_repr_css_set_property(ncss, "fill-opacity", "1.0");
+ }
+
+ sp_repr_css_set_property(ncsf, "stroke", "none");
+ sp_repr_css_set_property(ncsf, "stroke-width", nullptr);
+ sp_repr_css_set_property(ncsf, "stroke-opacity", "1.0");
+ sp_repr_css_set_property(ncsf, "filter", nullptr);
+ sp_repr_css_set_property(ncsf, "opacity", nullptr);
+ sp_repr_css_unset_property(ncsf, "marker-start");
+ sp_repr_css_unset_property(ncsf, "marker-mid");
+ sp_repr_css_unset_property(ncsf, "marker-end");
+
+ // The object tree -------------------
+
+ // Remember the position of the item
+ gint pos = item->getRepr()->position();
+
+ // Remember parent
+ Inkscape::XML::Node *parent = item->getRepr()->parent();
+
+ Inkscape::XML::Document *xml_doc = doc->getReprDoc();
+
+ // Create a group to put everything in.
+ Inkscape::XML::Node *g_repr = xml_doc->createElement("svg:g");
+
+ Inkscape::copy_object_properties(g_repr, item->getRepr());
+ // drop copied style, children will be re-styled (stroke becomes fill)
+ g_repr->removeAttribute("style");
+
+ // Add the group to the parent, move to the saved position
+ parent->addChildAtPos(g_repr, pos);
+
+ // The stroke ------------------------
+ Inkscape::XML::Node *stroke = nullptr;
+ if (s_val && g_strcmp0(s_val,"none") != 0 && stroke_path.size() > 0) {
+ stroke = xml_doc->createElement("svg:path");
+ sp_repr_css_change(stroke, ncss, "style");
+
+ stroke->setAttribute("d", sp_svg_write_path(stroke_path));
+ }
+ sp_repr_css_attr_unref(ncss);
+
+ // The fill --------------------------
+ Inkscape::XML::Node *fill = nullptr;
+ if (f_val && g_strcmp0(f_val,"none") != 0 && !legacy) {
+ fill = xml_doc->createElement("svg:path");
+ sp_repr_css_change(fill, ncsf, "style");
+
+ fill->setAttribute("d", sp_svg_write_path(fill_path));
+ }
+ sp_repr_css_attr_unref(ncsf);
+
+ // The markers -----------------------
+ Inkscape::XML::Node *markers = nullptr;
+
+ if (shape->hasMarkers()) {
+ auto linewidth = style->stroke_width.computed;
+ if (!legacy) {
+ markers = xml_doc->createElement("svg:g");
+ g_repr->addChildAtPos(markers, pos);
+ } else {
+ markers = g_repr;
+ }
+
+ // START marker
+ for (int i = 0; i < 2; i++) { // SP_MARKER_LOC and SP_MARKER_LOC_START
+ if ( SPObject *marker_obj = shape->_marker[i] ) {
+ Geom::Affine const m (sp_shape_marker_get_transform_at_start(fill_path.front().front()));
+ item_to_paths_add_marker(item, marker_obj, m, linewidth, true, markers, legacy);
+ }
+ }
+
+ // MID marker
+ for (int i = 0; i < 3; i += 2) { // SP_MARKER_LOC and SP_MARKER_LOC_MID
+ SPObject *midmarker_obj = shape->_marker[i];
+ if (!midmarker_obj) continue; // TODO use auto below
+ for(Geom::PathVector::const_iterator path_it = fill_path.begin(); path_it != fill_path.end(); ++path_it) {
+
+ // START position
+ if ( path_it != fill_path.begin() &&
+ ! ((path_it == (fill_path.end()-1)) && (path_it->size_default() == 0)) ) // if this is the last path and it is a moveto-only, there is no mid marker there
+ {
+ Geom::Affine const m (sp_shape_marker_get_transform_at_start(path_it->front()));
+ item_to_paths_add_marker(item, midmarker_obj, m, linewidth, false, markers, legacy);
+ }
+
+ // MID position
+ if (path_it->size_default() > 1) {
+ Geom::Path::const_iterator curve_it1 = path_it->begin(); // incoming curve
+ Geom::Path::const_iterator curve_it2 = ++(path_it->begin()); // outgoing curve
+ while (curve_it2 != path_it->end_default()) {
+ /* Put marker between curve_it1 and curve_it2.
+ * Loop to end_default (so including closing segment), because when a path is closed,
+ * there should be a midpoint marker between last segment and closing straight line segment
+ */
+ Geom::Affine const m (sp_shape_marker_get_transform(*curve_it1, *curve_it2));
+ item_to_paths_add_marker(item, midmarker_obj, m, linewidth, false, markers, legacy);
+
+ ++curve_it1;
+ ++curve_it2;
+ }
+ }
+
+ // END position
+ if ( path_it != (fill_path.end()-1) && !path_it->empty()) {
+ Geom::Curve const &lastcurve = path_it->back_default();
+ Geom::Affine const m = sp_shape_marker_get_transform_at_end(lastcurve);
+ item_to_paths_add_marker(item, midmarker_obj, m, linewidth, false, markers, legacy);
+ }
+ }
+ }
+
+ // END marker
+ for (int i = 0; i < 4; i += 3) { // SP_MARKER_LOC and SP_MARKER_LOC_END
+ if ( SPObject *marker_obj = shape->_marker[i] ) {
+ /* Get reference to last curve in the path.
+ * For moveto-only path, this returns the "closing line segment". */
+ Geom::Path const &path_last = fill_path.back();
+ unsigned int index = path_last.size_default();
+ if (index > 0) {
+ index--;
+ }
+ Geom::Curve const &lastcurve = path_last[index];
+
+ Geom::Affine const m = sp_shape_marker_get_transform_at_end(lastcurve);
+ item_to_paths_add_marker(item, marker_obj, m, linewidth, false, markers, legacy);
+ }
+ }
+ }
+
+ gchar const *paint_order = sp_repr_css_property(ncss, "paint-order", nullptr);
+ SPIPaintOrder temp;
+ temp.read( paint_order );
+ bool unique = false;
+ if ((!fill && !markers) || (!fill && !stroke) || (!markers && !stroke)) {
+ unique = true;
+ }
+ if (temp.layer[0] != SP_CSS_PAINT_ORDER_NORMAL && !legacy && !unique) {
+
+ if (temp.layer[0] == SP_CSS_PAINT_ORDER_FILL) {
+ if (temp.layer[1] == SP_CSS_PAINT_ORDER_STROKE) {
+ if ( fill ) {
+ g_repr->appendChild(fill);
+ }
+ if ( stroke ) {
+ g_repr->appendChild(stroke);
+ }
+ if ( markers ) {
+ markers->setPosition(2);
+ }
+ } else {
+ if ( fill ) {
+ g_repr->appendChild(fill);
+ }
+ if ( markers ) {
+ markers->setPosition(1);
+ }
+ if ( stroke ) {
+ g_repr->appendChild(stroke);
+ }
+ }
+ } else if (temp.layer[0] == SP_CSS_PAINT_ORDER_STROKE) {
+ if (temp.layer[1] == SP_CSS_PAINT_ORDER_FILL) {
+ if ( stroke ) {
+ g_repr->appendChild(stroke);
+ }
+ if ( fill ) {
+ g_repr->appendChild(fill);
+ }
+ if ( markers ) {
+ markers->setPosition(2);
+ }
+ } else {
+ if ( stroke ) {
+ g_repr->appendChild(stroke);
+ }
+ if ( markers ) {
+ markers->setPosition(1);
+ }
+ if ( fill ) {
+ g_repr->appendChild(fill);
+ }
+ }
+ } else {
+ if (temp.layer[1] == SP_CSS_PAINT_ORDER_STROKE) {
+ if ( markers ) {
+ markers->setPosition(0);
+ }
+ if ( stroke ) {
+ g_repr->appendChild(stroke);
+ }
+ if ( fill ) {
+ g_repr->appendChild(fill);
+ }
+ } else {
+ if ( markers ) {
+ markers->setPosition(0);
+ }
+ if ( fill ) {
+ g_repr->appendChild(fill);
+ }
+ if ( stroke ) {
+ g_repr->appendChild(stroke);
+ }
+ }
+ }
+
+ } else if (!unique) {
+ if ( fill ) {
+ g_repr->appendChild(fill);
+ }
+ if ( stroke ) {
+ g_repr->appendChild(stroke);
+ }
+ if ( markers ) {
+ markers->setPosition(2);
+ }
+ }
+
+ bool did = false;
+ // only consider it a change if more than a fill is created.
+ if (stroke || markers) {
+ did = true;
+ }
+
+ Inkscape::XML::Node *out = nullptr;
+
+ if (!fill && !markers && did) {
+ out = stroke;
+ } else if (!fill && !stroke && did) {
+ out = markers;
+ } else if(did) {
+ out = g_repr;
+ } else {
+ parent->removeChild(g_repr);
+ Inkscape::GC::release(g_repr);
+ if (fill) {
+ // Copy the style, to preserve context-fill cascade
+ if (context) {
+ item->setAttribute("style", fill->attribute("style"));
+ }
+ Inkscape::GC::release(fill);
+ }
+ return (flatten ? item->getRepr() : nullptr);
+ }
+
+ SPCSSAttr *r_style = sp_repr_css_attr_new();
+ sp_repr_css_set_property(r_style, "opacity", opacity);
+ sp_repr_css_set_property(r_style, "filter", filter);
+ sp_repr_css_change(out, r_style, "style");
+
+ sp_repr_css_attr_unref(r_style);
+ if (unique && out != markers) { // markers are already a child of g_repr
+ g_assert(out != g_repr);
+ parent->addChild(out, g_repr);
+ parent->removeChild(g_repr);
+ Inkscape::GC::release(g_repr);
+ }
+ out->setAttribute("transform", item->getRepr()->attribute("transform"));
+
+ // We're replacing item, delete it.
+ item->deleteObject(false);
+
+ out->setAttribute("id",id);
+ Inkscape::GC::release(out);
+
+ return out;
+}
+
+/*
+ Local Variables:
+ mode:c++
+ c-file-style:"stroustrup"
+ c-file-offsets:((innamespace . 0)(inline-open . 0)(case-label . +))
+ indent-tabs-mode:nil
+ fill-column:99
+ End:
+*/
+// vim: filetype=cpp:expandtab:shiftwidth=4:tabstop=8:softtabstop=4:fileencoding=utf-8:textwidth=99 :
diff --git a/src/path/path-outline.h b/src/path/path-outline.h
new file mode 100644
index 0000000..10df017
--- /dev/null
+++ b/src/path/path-outline.h
@@ -0,0 +1,60 @@
+// SPDX-License-Identifier: GPL-2.0-or-later
+
+/** @file
+ *
+ * Two related object to path operations:
+ *
+ * 1. Find a path that includes fill, stroke, and markers. Useful for finding a visual bounding box.
+ * 2. Take a set of objects and find an identical visual representation using only paths.
+ *
+ * Copyright (C) 2020 Tavmjong Bah
+ * Copyright (C) 2018 Authors
+ *
+ * Released under GNU GPL v2+, read the file 'COPYING' for more information.
+ */
+
+#ifndef SEEN_PATH_OUTLINE_H
+#define SEEN_PATH_OUTLINE_H
+
+class SPDesktop;
+class SPItem;
+
+namespace Geom {
+ class PathVector;
+}
+
+namespace Inkscape {
+namespace XML {
+ class Node;
+}
+}
+
+/**
+ * Find an outline that represents an item.
+ */
+Geom::PathVector* item_to_outline (SPItem const *item, bool exclude_markers = false);
+
+/**
+ * Replace item by path objects (a.k.a. stroke to path).
+ */
+Inkscape::XML::Node* item_to_paths(SPItem *item, bool legacy = false, SPItem *context = nullptr);
+
+/**
+ * Replace selected items by path objects (a.k.a. stroke to >path).
+ * TODO: remove desktop dependency.
+ */
+void selection_to_paths (SPDesktop *desktop, bool legacy = false);
+
+
+#endif // SEEN_PATH_OUTLINE_H
+
+/*
+ Local Variables:
+ mode:c++
+ c-file-style:"stroustrup"
+ c-file-offsets:((innamespace . 0)(inline-open . 0)(case-label . +))
+ indent-tabs-mode:nil
+ fill-column:99
+ End:
+*/
+// vim: filetype=cpp:expandtab:shiftwidth=4:tabstop=8:softtabstop=4:fileencoding=utf-8:textwidth=99 :
diff --git a/src/path/path-simplify.cpp b/src/path/path-simplify.cpp
new file mode 100644
index 0000000..2625866
--- /dev/null
+++ b/src/path/path-simplify.cpp
@@ -0,0 +1,126 @@
+// SPDX-License-Identifier: GPL-2.0-or-later
+/** @file
+ * Simplify paths (reduce node count).
+ *//*
+ * Authors:
+ * see git history
+ * Created by fred on Fri Dec 05 2003.
+ * tweaked endlessly by bulia byak <buliabyak@users.sf.net>
+ *
+ * Copyright (C) 2018 Authors
+ * Released under GNU GPL v2+, read the file 'COPYING' for more information.
+ */
+
+#ifdef HAVE_CONFIG_H
+#endif
+
+#include <vector>
+
+#include "path-simplify.h"
+#include "path-util.h"
+
+#include "document-undo.h"
+#include "preferences.h"
+
+#include "livarot/Path.h"
+
+#include "object/sp-item-group.h"
+#include "object/sp-path.h"
+
+using Inkscape::DocumentUndo;
+
+// Return number of paths simplified (can be greater than one if group).
+int
+path_simplify(SPItem *item, float threshold, bool justCoalesce, double size)
+{
+ //If this is a group, do the children instead
+ auto group = cast<SPGroup>(item);
+ if (group) {
+ int pathsSimplified = 0;
+ std::vector<SPItem*> items = group->item_list();
+ for (auto item : items) {
+ pathsSimplified += path_simplify(item, threshold, justCoalesce, size);
+ }
+ return pathsSimplified;
+ }
+
+ auto path = cast<SPPath>(item);
+ if (!path) {
+ return 0;
+ }
+
+ // There is actually no option in the preferences dialog for this!
+ Inkscape::Preferences *prefs = Inkscape::Preferences::get();
+ bool simplifyIndividualPaths = prefs->getBool("/options/simplifyindividualpaths/value");
+ if (simplifyIndividualPaths) {
+ Geom::OptRect itemBbox = item->documentVisualBounds();
+ if (itemBbox) {
+ size = L2(itemBbox->dimensions());
+ } else {
+ size = 0;
+ }
+ }
+
+ // Correct virtual size by full transform (bug #166937).
+ size /= item->i2doc_affine().descrim();
+
+ // Save the transform, to re-apply it after simplification.
+ Geom::Affine const transform(item->transform);
+
+ /*
+ reset the transform, effectively transforming the item by transform.inverse();
+ this is necessary so that the item is transformed twice back and forth,
+ allowing all compensations to cancel out regardless of the preferences
+ */
+ item->doWriteTransform(Geom::identity());
+
+ // SPLivarot: Start -----------------
+
+ // Get path to simplify (note that the path *before* LPE calculation is needed)
+ Path *orig = Path_for_item_before_LPE(item, false);
+ if (orig == nullptr) {
+ return 0;
+ }
+
+ if ( justCoalesce ) {
+ orig->Coalesce(threshold * size);
+ } else {
+ orig->ConvertEvenLines(threshold * size);
+ orig->Simplify(threshold * size);
+ }
+
+ // Path
+ gchar *str = orig->svg_dump_path();
+
+ // SPLivarot: End -------------------
+
+ char const *patheffect = item->getRepr()->attribute("inkscape:path-effect");
+ if (patheffect) {
+ item->setAttribute("inkscape:original-d", str);
+ } else {
+ item->setAttribute("d", str);
+ }
+ g_free(str);
+
+ // reapply the transform
+ item->doWriteTransform(transform);
+
+ // remove irrelevant old nodetypes attibute
+ item->removeAttribute("sodipodi:nodetypes");
+
+ // clean up
+ if (orig) delete orig;
+
+ return 1;
+}
+
+/*
+ Local Variables:
+ mode:c++
+ c-file-style:"stroustrup"
+ c-file-offsets:((innamespace . 0)(inline-open . 0)(case-label . +))
+ indent-tabs-mode:nil
+ fill-column:99
+ End:
+*/
+// vim: filetype=cpp:expandtab:shiftwidth=4:tabstop=8:softtabstop=4:fileencoding=utf-8:textwidth=99 :
diff --git a/src/path/path-simplify.h b/src/path/path-simplify.h
new file mode 100644
index 0000000..0beda44
--- /dev/null
+++ b/src/path/path-simplify.h
@@ -0,0 +1,28 @@
+// SPDX-License-Identifier: GPL-2.0-or-later
+/** @file
+ * Simplify paths (reduce node count).
+ *//*
+ * Authors: see git history
+ *
+ * Copyright (C) 2018 Authors
+ * Released under GNU GPL v2+, read the file 'COPYING' for more information.
+ */
+#ifndef PATH_SIMPLIFY_H
+#define PATH_SIMPLIFY_H
+
+class SPItem;
+
+int path_simplify(SPItem *item, float threshold, bool justCoalesce, double size);
+
+#endif // PATH_SIMPLIFY_H
+
+/*
+ Local Variables:
+ mode:c++
+ c-file-style:"stroustrup"
+ c-file-offsets:((innamespace . 0)(inline-open . 0)(case-label . +))
+ indent-tabs-mode:nil
+ fill-column:99
+ End:
+*/
+// vim: filetype=cpp:expandtab:shiftwidth=4:tabstop=8:softtabstop=4:fileencoding=utf-8:textwidth=99 :
diff --git a/src/path/path-util.cpp b/src/path/path-util.cpp
new file mode 100644
index 0000000..b1df9d3
--- /dev/null
+++ b/src/path/path-util.cpp
@@ -0,0 +1,184 @@
+// SPDX-License-Identifier: GPL-2.0-or-later
+/** @file
+ * Path utilities.
+ *//*
+ * Authors:
+ * see git history
+ * Created by fred on Fri Dec 05 2003.
+ * tweaked endlessly by bulia byak <buliabyak@users.sf.net>
+ *
+ * Copyright (C) 2018 Authors
+ * Released under GNU GPL v2+, read the file 'COPYING' for more information.
+ */
+
+#ifdef HAVE_CONFIG_H
+#endif
+
+#include <vector>
+
+#include "path-util.h"
+
+#include "text-editing.h"
+
+#include "livarot/Path.h"
+#include "livarot/Shape.h"
+
+#include "object/sp-flowtext.h"
+#include "object/sp-image.h"
+#include "object/sp-marker.h"
+#include "object/sp-path.h"
+#include "object/sp-text.h"
+
+#include "display/curve.h"
+
+Path *
+Path_for_pathvector(Geom::PathVector const &epathv)
+{
+ /*std::cout << "converting to Livarot path" << std::endl;
+
+ Geom::SVGPathWriter wr;
+ wr.feed(epathv);
+ std::cout << wr.str() << std::endl;*/
+
+ Path *dest = new Path;
+ dest->LoadPathVector(epathv);
+ return dest;
+}
+
+Path *
+Path_for_item(SPItem *item, bool doTransformation, bool transformFull)
+{
+ auto curve = curve_for_item(item);
+
+ if (!curve) {
+ return nullptr;
+ }
+
+ Geom::PathVector *pathv =
+ pathvector_for_curve(item, &*curve, doTransformation, transformFull, Geom::identity(), Geom::identity());
+
+ /*std::cout << "converting to Livarot path" << std::endl;
+
+ Geom::SVGPathWriter wr;
+ if (pathv) {
+ wr.feed(*pathv);
+ }
+ std::cout << wr.str() << std::endl;*/
+
+ Path *dest = new Path;
+ dest->LoadPathVector(*pathv);
+ delete pathv;
+
+ /*gchar *str = dest->svg_dump_path();
+ std::cout << "After conversion:\n" << str << std::endl;
+ g_free(str);*/
+
+ return dest;
+}
+
+Path *
+Path_for_item_before_LPE(SPItem *item, bool doTransformation, bool transformFull)
+{
+ auto curve = curve_for_item_before_LPE(item);
+
+ if (!curve) {
+ return nullptr;
+ }
+
+ Geom::PathVector *pathv =
+ pathvector_for_curve(item, &*curve, doTransformation, transformFull, Geom::identity(), Geom::identity());
+
+ Path *dest = new Path;
+ dest->LoadPathVector(*pathv);
+ delete pathv;
+
+ return dest;
+}
+
+Geom::PathVector*
+pathvector_for_curve(SPItem *item, SPCurve *curve, bool doTransformation, bool transformFull, Geom::Affine extraPreAffine, Geom::Affine extraPostAffine)
+{
+ if (curve == nullptr)
+ return nullptr;
+
+ Geom::PathVector *dest = new Geom::PathVector;
+ *dest = curve->get_pathvector(); // Make a copy; must be freed by the caller!
+
+ if (doTransformation) {
+ if (transformFull) {
+ *dest *= extraPreAffine * item->i2doc_affine() * extraPostAffine;
+ } else {
+ *dest *= extraPreAffine * (Geom::Affine)item->transform * extraPostAffine;
+ }
+ } else {
+ *dest *= extraPreAffine * extraPostAffine;
+ }
+
+ return dest;
+}
+
+std::optional<SPCurve> curve_for_item(SPItem *item)
+{
+ if (!item) {
+ return {};
+ }
+
+ if (auto path = cast<SPPath>(item)) {
+ return SPCurve::ptr_to_opt(path->curveForEdit());
+ } else if (auto shape = cast<SPShape>(item)) {
+ return SPCurve::ptr_to_opt(shape->curve());
+ } else if (is<SPText>(item) || is<SPFlowtext>(item)) {
+ return te_get_layout(item)->convertToCurves();
+ } else if (auto image = cast<SPImage>(item)) {
+ return SPCurve::ptr_to_opt(image->get_curve());
+ }
+
+ return {};
+}
+
+std::optional<SPCurve> curve_for_item_before_LPE(SPItem *item)
+{
+ if (!item) {
+ return {};
+ }
+
+ if (auto shape = cast<SPShape>(item)) {
+ return SPCurve::ptr_to_opt(shape->curveForEdit());
+ } else if (is<SPText>(item) || is<SPFlowtext>(item)) {
+ return te_get_layout(item)->convertToCurves();
+ } else if (auto image = cast<SPImage>(item)) {
+ return SPCurve::ptr_to_opt(image->get_curve());
+ }
+
+ return {};
+}
+
+std::optional<Path::cut_position> get_nearest_position_on_Path(Path *path, Geom::Point p, unsigned seg)
+{
+ std::optional<Path::cut_position> result;
+ if (!path) {
+ return result; // returns empty std::optional
+ }
+ //get nearest position on path
+ result = path->PointToCurvilignPosition(p, seg);
+ return result;
+}
+
+Geom::Point get_point_on_Path(Path *path, int piece, double t)
+{
+ Geom::Point p;
+ path->PointAt(piece, t, p);
+ return p;
+}
+
+
+/*
+ Local Variables:
+ mode:c++
+ c-file-style:"stroustrup"
+ c-file-offsets:((innamespace . 0)(inline-open . 0)(case-label . +))
+ indent-tabs-mode:nil
+ fill-column:99
+ End:
+*/
+// vim: filetype=cpp:expandtab:shiftwidth=4:tabstop=8:softtabstop=4:fileencoding=utf-8:textwidth=99 :
diff --git a/src/path/path-util.h b/src/path/path-util.h
new file mode 100644
index 0000000..77e3fcb
--- /dev/null
+++ b/src/path/path-util.h
@@ -0,0 +1,162 @@
+// SPDX-License-Identifier: GPL-2.0-or-later
+/** @file
+ * Path utilities.
+ *//*
+ * Authors: see git history
+ *
+ * Copyright (C) 2018 Authors
+ * Released under GNU GPL v2+, read the file 'COPYING' for more information.
+ */
+#ifndef PATH_UTIL_H
+#define PATH_UTIL_H
+
+#include <2geom/forward.h>
+#include <2geom/path.h>
+
+#include "livarot/Path.h"
+#include "display/curve.h"
+
+#include <optional>
+
+class SPItem;
+
+/**
+ * Creates a Livarot Path object from the Geom::PathVector.
+ *
+ * @param pathv A reference to the pathvector to convert.
+ *
+ * @return A pointer to the livarot Path object. The caller would need
+ * to free the object after using it.
+ */
+Path *Path_for_pathvector(Geom::PathVector const &pathv);
+
+/**
+ * Creates a Livarot Path object from an SPItem. The Geom::PathVector extracted from
+ * the item would be before applying LPEs for SPPath and after applying LPEs for other shapes.
+ *
+ * @param item A pointer to the SPItem object.
+ * @param doTransformation See the same parameter in pathvector_for_curve().
+ * @param transformFull See the same parameter in pathvector_for_curve().
+ *
+ * @return A pointer to the livarot Path object. The caller would need
+ * to free the object after using it.
+ */
+Path *Path_for_item(SPItem *item, bool doTransformation, bool transformFull = true);
+
+/**
+ * Creates a Livarot Path object from the SPItem. This function ensures that the Geom::PathVector
+ * extracted is the one before applying the LPE stack.
+ *
+ * @param item A pointer to the SPItem object.
+ * @param doTransformation See the same parameter in pathvector_for_curve().
+ * @param transformFull See the same parameter in pathvector_for_curve().
+ *
+ * @return A pointer to the livarot Path object. The caller would need
+ * to free the object after using it.
+ */
+Path *Path_for_item_before_LPE(SPItem *item, bool doTransformation, bool transformFull = true);
+
+/**
+ * Gets a Geom::PathVector from the SPCurve object.
+ *
+ * TODO: see if calling this method can be optimized. All the pathvector copying might be slow.
+ *
+ * @param item A pointer to the original SPItem. This is required to get the transformation
+ * information.
+ * @param curve A pointer to the SPCurve. If this pointer is null, the function returns nullptr too.
+ * @param doTransformation If set to true, the transformation in the SPItem is applied.
+ * @param transformFull If doTransformation and transformFull are both set to true, the
+ * i2doc_affine transformation is applied, which includes the all the transformations of the
+ * ancestors and even the document viewport. If doTransformation is true but transformFull is false
+ * only the item's own transformation gets applied.
+ * @param extraPreAffine This is a Geom::Affine transformation that gets applied before any
+ * transformations from SPItem.
+ * @param extraPostAffine This is a Geom::Affine transformation that gets applied after any
+ * transformations from SPItem.
+ *
+ * @return A pointer to the Geom::PathVector. Must be freed by the caller.
+ */
+Geom::PathVector* pathvector_for_curve(SPItem *item, SPCurve *curve, bool doTransformation, bool transformFull, Geom::Affine extraPreAffine, Geom::Affine extraPostAffine);
+
+/**
+ * Gets an SPCurve from the SPItem.
+ *
+ * An SPCurve is basically a wrapper around a Geom::PathVector. This function extracts
+ * an SPCurve object from the SPItem. The behavior is a little ill-defined though. It
+ * returns the path before applying LPE stack if the item is an SPPath and the path
+ * after applying the LPE stack for all other types.
+ *
+ * For SPText, SPFlowtext and SPImage LPEs don't matter anyways (AFAIK). So, curve_for_item
+ * and curve_for_item_before_LPE behave identically. For SPShape objects, there is a difference.
+ * curve_for_item returns path before LPE only if it's SPPath and path after LPE otherwise.
+ * curve_for_item_before_LPE returns path before LPE all the time. See the inheritance diagram
+ * for these classes if this doesn't make sense.
+ *
+ * @param item Pointer to the SPItem object.
+ *
+ * @returns The extracted SPCurve
+ */
+std::optional<SPCurve> curve_for_item(SPItem *item);
+
+/**
+ * Gets an SPCurve from the SPItem before any LPE.
+ *
+ * An SPCurve is basically a wrapper around a Geom::PathVector. This function extracts
+ * an SPCurve object from the SPItem but it ensures that if LPEs are supported on that
+ * SPItem the path before the LPE is returned.
+ *
+ * For SPText, SPFlowtext and SPImage LPEs don't matter anyways (AFAIK). So, curve_for_item
+ * and curve_for_item_before_LPE behave identically. For SPShape objects, there is a difference.
+ * curve_for_item returns path before LPE only if it's SPPath and path after LPE otherwise.
+ * curve_for_item_before_LPE returns path before LPE all the time. See the inheritance diagram
+ * for these classes if this doesn't make sense.
+ *
+ * @param item Pointer to the SPItem object.
+ *
+ * @returns The extracted SPCurve
+ */
+std::optional<SPCurve> curve_for_item_before_LPE(SPItem *item);
+
+/**
+ * Get the nearest position given a Livarot Path and a point.
+ *
+ * I've not properly read the internals of the functions that get called for this function.
+ * However, if the name is correct, it seems the function returns the point in time and the
+ * piece where the path gets closest to the given point.
+ *
+ * TODO: Confirm this by checking the internals of the function that gets called.
+ *
+ * @param path A pointer to the Livarot Path object.
+ * @param p The point to which you wanna find the nearest point on the Path.
+ * @param seg If it is set to 0, all segments are considered. If set to any other number, only
+ * that segment will be considered.
+ *
+ * @return The point at time t.
+ */
+std::optional<Path::cut_position> get_nearest_position_on_Path(Path *path, Geom::Point p, unsigned seg = 0);
+
+/**
+ * Gets the point at a particular time in a particular piece in a path description.
+ *
+ * @param path Pointer to the Livarot Path object.
+ * @param piece The index of the path description. If this index is invalid, a point
+ * at the origin will be returned.
+ * @param t The time on the path description. Should be between 0 and 1.
+ *
+ * @return The cut position. This structure has both the piece index and the time for the nearest
+ * point. This cut position object is wrapped in boost::optional.
+ */
+Geom::Point get_point_on_Path(Path *path, int piece, double t);
+
+#endif // PATH_UTIL_H
+
+/*
+ Local Variables:
+ mode:c++
+ c-file-style:"stroustrup"
+ c-file-offsets:((innamespace . 0)(inline-open . 0)(case-label . +))
+ indent-tabs-mode:nil
+ fill-column:99
+ End:
+*/
+// vim: filetype=cpp:expandtab:shiftwidth=4:tabstop=8:softtabstop=4:fileencoding=utf-8:textwidth=99 :
diff --git a/src/path/splinefit/bezier-fit.cpp b/src/path/splinefit/bezier-fit.cpp
new file mode 100644
index 0000000..111f538
--- /dev/null
+++ b/src/path/splinefit/bezier-fit.cpp
@@ -0,0 +1,90 @@
+// SPDX-License-Identifier: GPL-2.0-or-later
+
+#include <iostream>
+#include <vector>
+#include "bezier-fit.h"
+#include <2geom/bezier-utils.h>
+#include <2geom/point.h>
+
+extern "C" {
+ #include "splinefit.h"
+ #include "splinefont.h"
+}
+
+int bezier_fit(Geom::Point bezier[4], const std::vector<InputPoint>& data) {
+
+ if (data.size() <= 2) return 0;
+
+ int order2 = false; // not 2nd order, so cubic
+ // "Fitting cubic BĂ©zier curves"
+ // https://raphlinus.github.io/curves/2021/03/11/bezier-fitting.html
+ mergetype mt = mt_levien;
+ auto len = data.size();
+
+ std::vector<FitPoint> fit;
+ for (int i = 0; i < len; ++i) {
+ fit.push_back({});
+ auto& fp = fit.back();
+ fp.p.x = data[i].x();
+ fp.p.y = data[i].y();
+ fp.t = data[i].t;
+ fp.ut.x = fp.ut.y = 0;
+ }
+
+ // transform data into spline set format
+
+ auto input = (SplineSet*)chunkalloc(sizeof(SplineSet));
+
+ for (int i = 0; i < len; ++i) {
+ auto& d = data[i];
+ auto sp = SplinePointCreate(d.x(), d.y());
+ if (d.have_slope) {
+ sp->nextcp.x = d.front.x();
+ sp->nextcp.y = d.front.y();
+ sp->nonextcp = false;
+ sp->prevcp.x = d.back.x();
+ sp->prevcp.y = d.back.y();
+ sp->noprevcp = false;
+ }
+
+ if (i == 0) {
+ input->first = input->last = sp;
+ }
+ else {
+ SplineMake(input->last, sp, order2);
+ input->last = sp;
+ }
+ }
+
+ Spline* spline = ApproximateSplineFromPointsSlopes(input->first, input->last, fit.data(), fit.size(), order2, mt);
+ bool ok = spline != nullptr;
+
+ if (!spline) {
+ std::vector<Geom::Point> inp;
+ inp.reserve(data.size());
+ for (auto& pt : data) {
+ inp.push_back(pt);
+ }
+ ok = bezier_fit_cubic(bezier, inp.data(), inp.size(), 0.5) > 0;
+ }
+
+ if (spline) {
+ bezier[0].x() = spline->from->me.x;
+ bezier[0].y() = spline->from->me.y;
+
+ bezier[1].x() = spline->from->nextcp.x;
+ bezier[1].y() = spline->from->nextcp.y;
+
+ bezier[2].x() = spline->to->prevcp.x;
+ bezier[2].y() = spline->to->prevcp.y;
+
+ bezier[3].x() = spline->to->me.x;
+ bezier[3].y() = spline->to->me.y;
+ }
+
+ SplinePointListFree(input);
+ //TODO: verify that all C structs are freed up
+ // SplineFree(spline);
+
+ return ok;
+}
diff --git a/src/path/splinefit/bezier-fit.h b/src/path/splinefit/bezier-fit.h
new file mode 100644
index 0000000..7eb0440
--- /dev/null
+++ b/src/path/splinefit/bezier-fit.h
@@ -0,0 +1,20 @@
+// SPDX-License-Identifier: GPL-2.0-or-later
+
+#include <vector>
+#include "2geom/point.h"
+
+struct InputPoint : Geom::Point {
+ InputPoint() {}
+ InputPoint(const Geom::Point& pt) : Point(pt) {}
+ InputPoint(const Geom::Point& pt, double t) : Point(pt), t(t) {}
+ InputPoint(const Geom::Point& pt, const Geom::Point& front, const Geom::Point& back, double t)
+ : Point(pt), front(front), back(back), t(t), have_slope(true) {}
+
+ Geom::Point front;
+ Geom::Point back;
+ double t = 0;
+ bool have_slope = false;
+};
+
+// Fit cubic Bezier to input points; use slope of the first and last points to find a fit
+int bezier_fit(Geom::Point bezier[4], const std::vector<InputPoint>& data);
diff --git a/src/path/splinefit/splinefit.c b/src/path/splinefit/splinefit.c
new file mode 100644
index 0000000..9e3039a
--- /dev/null
+++ b/src/path/splinefit/splinefit.c
@@ -0,0 +1,1427 @@
+// SPDX-License-Identifier: GPL-2.0-or-later
+/* -*- coding: utf-8 -*- */
+/* Copyright (C) 2000-2012 by George Williams, 2021 by Linus Romer */
+/*
+ * Redistribution and use in source and binary forms, with or without
+ * modification, are permitted provided that the following conditions are met:
+
+ * Redistributions of source code must retain the above copyright notice, this
+ * list of conditions and the following disclaimer.
+
+ * Redistributions in binary form must reproduce the above copyright notice,
+ * this list of conditions and the following disclaimer in the documentation
+ * and/or other materials provided with the distribution.
+
+ * The name of the author may not be used to endorse or promote products
+ * derived from this software without specific prior written permission.
+
+ * THIS SOFTWARE IS PROVIDED BY THE AUTHOR ``AS IS'' AND ANY EXPRESS OR IMPLIED
+ * WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF
+ * MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO
+ * EVENT SHALL THE AUTHOR BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
+ * SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
+ * PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS;
+ * OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY,
+ * WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR
+ * OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF
+ * ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+ */
+
+/* Incorporated into Inkscape sources; splinefit.c as of 2023-02-15.
+
+ * Note: The only changes to this file in Inkscape codebase are modification of includes
+ * and adding this note. Formatting is intact to facilitate future merges.
+ */
+
+#include <memory.h>
+#include <stdbool.h>
+#include <stdint.h>
+#include "splinefit.h"
+
+#include <math.h>
+
+static Spline *IsLinearApprox(SplinePoint *from, SplinePoint *to,
+ FitPoint *mid, int cnt, int order2) {
+ bigreal vx, vy, slope;
+ int i;
+
+ vx = to->me.x-from->me.x; vy = to->me.y-from->me.y;
+ if ( vx==0 && vy==0 ) {
+ for ( i=0; i<cnt; ++i )
+ if ( mid[i].p.x != from->me.x || mid[i].p.y != from->me.y )
+return( NULL );
+ } else if ( fabs(vx)>fabs(vy) ) {
+ slope = vy/vx;
+ for ( i=0; i<cnt; ++i )
+ if ( !RealWithin(mid[i].p.y,from->me.y+slope*(mid[i].p.x-from->me.x),.7) )
+return( NULL );
+ } else {
+ slope = vx/vy;
+ for ( i=0; i<cnt; ++i )
+ if ( !RealWithin(mid[i].p.x,from->me.x+slope*(mid[i].p.y-from->me.y),.7) )
+return( NULL );
+ }
+ from->nextcp = from->me;
+ to->prevcp = to->me;
+return( SplineMake(from,to,order2) );
+}
+
+/* This routine should almost never be called now. It uses a flawed algorithm */
+/* which won't produce the best results. It gets called only when the better */
+/* approach doesn't work (singular matrices, etc.) */
+/* Old comment, back when I was confused... */
+/* Least squares tells us that:
+ | S(xi*ti^3) | | S(ti^6) S(ti^5) S(ti^4) S(ti^3) | | a |
+ | S(xi*ti^2) | = | S(ti^5) S(ti^4) S(ti^3) S(ti^2) | * | b |
+ | S(xi*ti) | | S(ti^4) S(ti^3) S(ti^2) S(ti) | | c |
+ | S(xi) | | S(ti^3) S(ti^2) S(ti) n | | d |
+ and the definition of a spline tells us:
+ | x1 | = | 1 1 1 1 | * (a b c d)
+ | x0 | = | 0 0 0 1 | * (a b c d)
+So we're a bit over specified. Let's use the last two lines of least squares
+and the 2 from the spline defn. So d==x0. Now we've got three unknowns
+and only three equations...
+
+For order2 splines we've got
+ | S(xi*ti^2) | | S(ti^4) S(ti^3) S(ti^2) | | b |
+ | S(xi*ti) | = | S(ti^3) S(ti^2) S(ti) | * | c |
+ | S(xi) | | S(ti^2) S(ti) n | | d |
+ and the definition of a spline tells us:
+ | x1 | = | 1 1 1 | * (b c d)
+ | x0 | = | 0 0 1 | * (b c d)
+=>
+ d = x0
+ b+c = x1-x0
+ S(ti^2)*b + S(ti)*c = S(xi)-n*x0
+ S(ti^2)*b + S(ti)*(x1-x0-b) = S(xi)-n*x0
+ [ S(ti^2)-S(ti) ]*b = S(xi)-S(ti)*(x1-x0) - n*x0
+*/
+static int _ApproximateSplineFromPoints(SplinePoint *from, SplinePoint *to,
+ FitPoint *mid, int cnt, BasePoint *nextcp, BasePoint *prevcp,
+ int order2) {
+ bigreal tt, ttn;
+ int i, j, ret;
+ bigreal vx[3], vy[3], m[3][3];
+ bigreal ts[7], xts[4], yts[4];
+ BasePoint nres, pres;
+ int nrescnt=0, prescnt=0;
+ bigreal nmin, nmax, pmin, pmax, test, ptest;
+ bigreal bx, by, cx, cy;
+
+ memset(&nres,0,sizeof(nres)); memset(&pres,0,sizeof(pres));
+
+ /* Add the initial and end points */
+ ts[0] = 2; for ( i=1; i<7; ++i ) ts[i] = 1;
+ xts[0] = from->me.x+to->me.x; yts[0] = from->me.y+to->me.y;
+ xts[3] = xts[2] = xts[1] = to->me.x; yts[3] = yts[2] = yts[1] = to->me.y;
+ nmin = pmin = 0; nmax = pmax = (to->me.x-from->me.x)*(to->me.x-from->me.x)+(to->me.y-from->me.y)*(to->me.y-from->me.y);
+ for ( i=0; i<cnt; ++i ) {
+ xts[0] += mid[i].p.x;
+ yts[0] += mid[i].p.y;
+ ++ts[0];
+ tt = mid[i].t;
+ xts[1] += tt*mid[i].p.x;
+ yts[1] += tt*mid[i].p.y;
+ ts[1] += tt;
+ ts[2] += (ttn=tt*tt);
+ xts[2] += ttn*mid[i].p.x;
+ yts[2] += ttn*mid[i].p.y;
+ ts[3] += (ttn*=tt);
+ xts[3] += ttn*mid[i].p.x;
+ yts[3] += ttn*mid[i].p.y;
+ ts[4] += (ttn*=tt);
+ ts[5] += (ttn*=tt);
+ ts[6] += (ttn*=tt);
+
+ test = (mid[i].p.x-from->me.x)*(to->me.x-from->me.x) + (mid[i].p.y-from->me.y)*(to->me.y-from->me.y);
+ if ( test<nmin ) nmin=test;
+ if ( test>nmax ) nmax=test;
+ test = (mid[i].p.x-to->me.x)*(from->me.x-to->me.x) + (mid[i].p.y-to->me.y)*(from->me.y-to->me.y);
+ if ( test<pmin ) pmin=test;
+ if ( test>pmax ) pmax=test;
+ }
+ pmin *= 1.2; pmax *= 1.2; nmin *= 1.2; nmax *= 1.2;
+
+ for ( j=0; j<3; ++j ) {
+ if ( order2 ) {
+ if ( RealNear(ts[j+2],ts[j+1]) )
+ continue;
+ /* This produces really bad results!!!! But I don't see what I can do to improve it */
+ bx = (xts[j]-ts[j+1]*(to->me.x-from->me.x) - ts[j]*from->me.x) / (ts[j+2]-ts[j+1]);
+ by = (yts[j]-ts[j+1]*(to->me.y-from->me.y) - ts[j]*from->me.y) / (ts[j+2]-ts[j+1]);
+ cx = to->me.x-from->me.x-bx;
+ cy = to->me.y-from->me.y-by;
+
+ nextcp->x = from->me.x + cx/2;
+ nextcp->y = from->me.y + cy/2;
+ *prevcp = *nextcp;
+ } else {
+ vx[0] = xts[j+1]-ts[j+1]*from->me.x;
+ vx[1] = xts[j]-ts[j]*from->me.x;
+ vx[2] = to->me.x-from->me.x; /* always use the defn of spline */
+
+ vy[0] = yts[j+1]-ts[j+1]*from->me.y;
+ vy[1] = yts[j]-ts[j]*from->me.y;
+ vy[2] = to->me.y-from->me.y;
+
+ m[0][0] = ts[j+4]; m[0][1] = ts[j+3]; m[0][2] = ts[j+2];
+ m[1][0] = ts[j+3]; m[1][1] = ts[j+2]; m[1][2] = ts[j+1];
+ m[2][0] = 1; m[2][1] = 1; m[2][2] = 1;
+
+ /* Remove a terms from rows 0 and 1 */
+ vx[0] -= ts[j+4]*vx[2];
+ vy[0] -= ts[j+4]*vy[2];
+ m[0][0] = 0; m[0][1] -= ts[j+4]; m[0][2] -= ts[j+4];
+ vx[1] -= ts[j+3]*vx[2];
+ vy[1] -= ts[j+3]*vy[2];
+ m[1][0] = 0; m[1][1] -= ts[j+3]; m[1][2] -= ts[j+3];
+
+ if ( fabs(m[1][1])<fabs(m[0][1]) ) {
+ bigreal temp;
+ temp = vx[1]; vx[1] = vx[0]; vx[0] = temp;
+ temp = vy[1]; vy[1] = vy[0]; vy[0] = temp;
+ temp = m[1][1]; m[1][1] = m[0][1]; m[0][1] = temp;
+ temp = m[1][2]; m[1][2] = m[0][2]; m[0][2] = temp;
+ }
+ /* remove b terms from rows 0 and 2 (first normalize row 1 so m[1][1] is 1*/
+ vx[1] /= m[1][1];
+ vy[1] /= m[1][1];
+ m[1][2] /= m[1][1];
+ m[1][1] = 1;
+ vx[0] -= m[0][1]*vx[1];
+ vy[0] -= m[0][1]*vy[1];
+ m[0][2] -= m[0][1]*m[1][2]; m[0][1] = 0;
+ vx[2] -= m[2][1]*vx[1];
+ vy[2] -= m[2][1]*vy[1];
+ m[2][2] -= m[2][1]*m[1][2]; m[2][1] = 0;
+
+ vx[0] /= m[0][2]; /* This is cx */
+ vy[0] /= m[0][2]; /* This is cy */
+ /*m[0][2] = 1;*/
+
+ vx[1] -= m[1][2]*vx[0]; /* This is bx */
+ vy[1] -= m[1][2]*vy[0]; /* This is by */
+ /* m[1][2] = 0; */
+ vx[2] -= m[2][2]*vx[0]; /* This is ax */
+ vy[2] -= m[2][2]*vy[0]; /* This is ay */
+ /* m[2][2] = 0; */
+
+ nextcp->x = from->me.x + vx[0]/3;
+ nextcp->y = from->me.y + vy[0]/3;
+ prevcp->x = nextcp->x + (vx[0]+vx[1])/3;
+ prevcp->y = nextcp->y + (vy[0]+vy[1])/3;
+ }
+
+ test = (nextcp->x-from->me.x)*(to->me.x-from->me.x) +
+ (nextcp->y-from->me.y)*(to->me.y-from->me.y);
+ ptest = (prevcp->x-to->me.x)*(from->me.x-to->me.x) +
+ (prevcp->y-to->me.y)*(from->me.y-to->me.y);
+ if ( order2 &&
+ (test<nmin || test>nmax || ptest<pmin || ptest>pmax))
+ continue;
+ if ( test>=nmin && test<=nmax ) {
+ nres.x += nextcp->x; nres.y += nextcp->y;
+ ++nrescnt;
+ }
+ if ( test>=pmin && test<=pmax ) {
+ pres.x += prevcp->x; pres.y += prevcp->y;
+ ++prescnt;
+ }
+ if ( nrescnt==1 && prescnt==1 )
+ break;
+ }
+
+ ret = 0;
+ if ( nrescnt>0 ) {
+ ret |= 1;
+ nextcp->x = nres.x/nrescnt;
+ nextcp->y = nres.y/nrescnt;
+ } else
+ *nextcp = from->nextcp;
+ if ( prescnt>0 ) {
+ ret |= 2;
+ prevcp->x = pres.x/prescnt;
+ prevcp->y = pres.y/prescnt;
+ } else
+ *prevcp = to->prevcp;
+ if ( order2 && ret!=3 ) {
+ nextcp->x = (nextcp->x + prevcp->x)/2;
+ nextcp->y = (nextcp->y + prevcp->y)/2;
+ }
+ if ( order2 )
+ *prevcp = *nextcp;
+return( ret );
+}
+
+static void TestForLinear(SplinePoint *from,SplinePoint *to) {
+ BasePoint off, cpoff, cpoff2;
+ bigreal len, co, co2;
+
+ /* Did we make a line? */
+ off.x = to->me.x-from->me.x; off.y = to->me.y-from->me.y;
+ len = sqrt(off.x*off.x + off.y*off.y);
+ if ( len!=0 ) {
+ off.x /= len; off.y /= len;
+ cpoff.x = from->nextcp.x-from->me.x; cpoff.y = from->nextcp.y-from->me.y;
+ len = sqrt(cpoff.x*cpoff.x + cpoff.y*cpoff.y);
+ if ( len!=0 ) {
+ cpoff.x /= len; cpoff.y /= len;
+ }
+ cpoff2.x = to->prevcp.x-from->me.x; cpoff2.y = to->prevcp.y-from->me.y;
+ len = sqrt(cpoff2.x*cpoff2.x + cpoff2.y*cpoff2.y);
+ if ( len!=0 ) {
+ cpoff2.x /= len; cpoff2.y /= len;
+ }
+ co = cpoff.x*off.y - cpoff.y*off.x; co2 = cpoff2.x*off.y - cpoff2.y*off.x;
+ if ( co<.05 && co>-.05 && co2<.05 && co2>-.05 ) {
+ from->nextcp = from->me;
+ to->prevcp = to->me;
+ } else {
+ Spline temp;
+ memset(&temp,0,sizeof(temp));
+ temp.from = from; temp.to = to;
+ SplineRefigure(&temp);
+ if ( SplineIsLinear(&temp)) {
+ from->nextcp = from->me;
+ to->prevcp = to->me;
+ }
+ }
+ }
+}
+
+/* Find a spline which best approximates the list of intermediate points we */
+/* are given. No attempt is made to use fixed slope angles */
+/* given a set of points (x,y,t) */
+/* find the bezier spline which best fits those points */
+
+/* OK, we know the end points, so all we really need are the control points */
+ /* For cubics.... */
+/* Pf = point from */
+/* CPf = control point, from nextcp */
+/* CPt = control point, to prevcp */
+/* Pt = point to */
+/* S(t) = Pf + 3*(CPf-Pf)*t + 3*(CPt-2*CPf+Pf)*t^2 + (Pt-3*CPt+3*CPf-Pf)*t^3 */
+/* S(t) = Pf - 3*Pf*t + 3*Pf*t^2 - Pf*t^3 + Pt*t^3 + */
+/* 3*(t-2*t^2+t^3)*CPf + */
+/* 3*(t^2-t^3)*CPt */
+/* We want to minimize ÎŁ [S(ti)-Pi]^2 */
+/* There are four variables CPf.x, CPf.y, CPt.x, CPt.y */
+/* When we take the derivative of the error term above with each of these */
+/* variables, we find that the two coordinates are separate. So I shall only */
+/* work through the equations once, leaving off the coordinate */
+/* d error/dCPf = ÎŁ 2*3*(t-2*t^2+t^3) * [S(ti)-Pi] = 0 */
+/* d error/dCPt = ÎŁ 2*3*(t^2-t^3) * [S(ti)-Pi] = 0 */
+ /* For quadratics.... */
+/* CP = control point, there's only one */
+/* S(t) = Pf + 2*(CP-Pf)*t + (Pt-2*CP+Pf)*t^2 */
+/* S(t) = Pf - 2*Pf*t + Pf*t^2 + Pt*t^2 + */
+/* 2*(t-2*t^2)*CP */
+/* We want to minimize ÎŁ [S(ti)-Pi]^2 */
+/* There are two variables CP.x, CP.y */
+/* d error/dCP = ÎŁ 2*2*(t-2*t^2) * [S(ti)-Pi] = 0 */
+/* ÎŁ (t-2*t^2) * [Pf - 2*Pf*t + Pf*t^2 + Pt*t^2 - Pi + */
+/* 2*(t-2*t^2)*CP] = 0 */
+/* CP * (ÎŁ 2*(t-2*t^2)*(t-2*t^2)) = ÎŁ (t-2*t^2) * [Pf - 2*Pf*t + Pf*t^2 + Pt*t^2 - Pi] */
+
+/* ÎŁ (t-2*t^2) * [Pf - 2*Pf*t + Pf*t^2 + Pt*t^2 - Pi] */
+/* CP = ----------------------------------------------------- */
+/* ÎŁ 2*(t-2*t^2)*(t-2*t^2) */
+Spline *ApproximateSplineFromPoints(SplinePoint *from, SplinePoint *to,
+ FitPoint *mid, int cnt, int order2) {
+ int ret;
+ Spline *spline;
+ BasePoint nextcp, prevcp;
+ int i;
+
+ if ( order2 ) {
+ bigreal xconst, yconst, term /* Same for x and y */;
+ xconst = yconst = term = 0;
+ for ( i=0; i<cnt; ++i ) {
+ bigreal t = mid[i].t, t2 = t*t;
+ bigreal tfactor = (t-2*t2);
+ term += 2*tfactor*tfactor;
+ xconst += tfactor*(from->me.x*(1-2*t+t2) + to->me.x*t2 - mid[i].p.x);
+ yconst += tfactor*(from->me.y*(1-2*t+t2) + to->me.y*t2 - mid[i].p.y);
+ }
+ if ( term!=0 ) {
+ BasePoint cp;
+ cp.x = xconst/term; cp.y = yconst/term;
+ from->nextcp = to->prevcp = cp;
+return( SplineMake2(from,to));
+ }
+ } else {
+ bigreal xconst[2], yconst[2], f_term[2], t_term[2] /* Same for x and y */;
+ bigreal tfactor[2], determinant;
+ xconst[0] = xconst[1] = yconst[0] = yconst[1] =
+ f_term[0] = f_term[1] = t_term[0] = t_term[1] = 0;
+ for ( i=0; i<cnt; ++i ) {
+ bigreal t = mid[i].t, t2 = t*t, t3=t*t2;
+ bigreal xc = (from->me.x*(1-3*t+3*t2-t3) + to->me.x*t3 - mid[i].p.x);
+ bigreal yc = (from->me.y*(1-3*t+3*t2-t3) + to->me.y*t3 - mid[i].p.y);
+ tfactor[0] = (t-2*t2+t3); tfactor[1]=(t2-t3);
+ xconst[0] += tfactor[0]*xc;
+ xconst[1] += tfactor[1]*xc;
+ yconst[0] += tfactor[0]*yc;
+ yconst[1] += tfactor[1]*yc;
+ f_term[0] += 3*tfactor[0]*tfactor[0];
+ f_term[1] += 3*tfactor[0]*tfactor[1];
+ /*t_term[0] += 3*tfactor[1]*tfactor[0];*/
+ t_term[1] += 3*tfactor[1]*tfactor[1];
+ }
+ t_term[0] = f_term[1];
+ determinant = f_term[1]*t_term[0] - f_term[0]*t_term[1];
+ if ( determinant!=0 ) {
+ to->prevcp.x = -(xconst[0]*f_term[1]-xconst[1]*f_term[0])/determinant;
+ to->prevcp.y = -(yconst[0]*f_term[1]-yconst[1]*f_term[0])/determinant;
+ if ( f_term[0]!=0 ) {
+ from->nextcp.x = (-xconst[0]-t_term[0]*to->prevcp.x)/f_term[0];
+ from->nextcp.y = (-yconst[0]-t_term[0]*to->prevcp.y)/f_term[0];
+ } else {
+ from->nextcp.x = (-xconst[1]-t_term[1]*to->prevcp.x)/f_term[1];
+ from->nextcp.y = (-yconst[1]-t_term[1]*to->prevcp.y)/f_term[1];
+ }
+return( SplineMake3(from,to));
+ }
+ }
+
+ if ( (spline = IsLinearApprox(from,to,mid,cnt,order2))!=NULL )
+return( spline );
+
+ ret = _ApproximateSplineFromPoints(from,to,mid,cnt,&nextcp,&prevcp,order2);
+
+ if ( ret&1 ) {
+ from->nextcp = nextcp;
+ } else {
+ from->nextcp = from->me;
+ }
+ if ( ret&2 ) {
+ to->prevcp = prevcp;
+ } else {
+ to->prevcp = to->me;
+ }
+ TestForLinear(from,to);
+ spline = SplineMake(from,to,order2);
+return( spline );
+}
+
+static bigreal ClosestSplineSolve(Spline1D *sp,bigreal sought,bigreal close_to_t) {
+ /* We want to find t so that spline(t) = sought */
+ /* find the value which is closest to close_to_t */
+ /* on error return closetot */
+ extended ts[3];
+ int i;
+ bigreal t, best, test;
+
+ _CubicSolve(sp,sought,ts);
+ best = 9e20; t= close_to_t;
+ for ( i=0; i<3; ++i ) if ( ts[i]>-.0001 && ts[i]<1.0001 ) {
+ if ( (test=ts[i]-close_to_t)<0 ) test = -test;
+ if ( test<best ) {
+ best = test;
+ t = ts[i];
+ }
+ }
+
+return( t );
+}
+
+struct dotbounds {
+ BasePoint unit;
+ BasePoint base;
+ bigreal len;
+ bigreal min,max; /* If min<0 || max>len the spline extends beyond its endpoints */
+};
+
+static bigreal SigmaDeltas(Spline *spline, FitPoint *mid, int cnt, DBounds *b, struct dotbounds *db) {
+ int i;
+ bigreal xdiff, ydiff, sum, temp, t;
+ SplinePoint *to = spline->to, *from = spline->from;
+ extended ts[2], x,y;
+ struct dotbounds db2;
+ bigreal dot;
+ int near_vert, near_horiz;
+
+ if ( (xdiff = to->me.x-from->me.x)<0 ) xdiff = -xdiff;
+ if ( (ydiff = to->me.y-from->me.y)<0 ) ydiff = -ydiff;
+ near_vert = ydiff>2*xdiff;
+ near_horiz = xdiff>2*ydiff;
+
+ sum = 0;
+ for ( i=0; i<cnt; ++i ) {
+ if ( near_vert ) {
+ t = ClosestSplineSolve(&spline->splines[1],mid[i].p.y,mid[i].t);
+ } else if ( near_horiz ) {
+ t = ClosestSplineSolve(&spline->splines[0],mid[i].p.x,mid[i].t);
+ } else {
+ t = (ClosestSplineSolve(&spline->splines[1],mid[i].p.y,mid[i].t) +
+ ClosestSplineSolve(&spline->splines[0],mid[i].p.x,mid[i].t))/2;
+ }
+ temp = mid[i].p.x - ( ((spline->splines[0].a*t+spline->splines[0].b)*t+spline->splines[0].c)*t + spline->splines[0].d );
+ sum += temp*temp;
+ temp = mid[i].p.y - ( ((spline->splines[1].a*t+spline->splines[1].b)*t+spline->splines[1].c)*t + spline->splines[1].d );
+ sum += temp*temp;
+ }
+
+ /* Ok, we've got distances from a set of points on the old spline to the */
+ /* new one. Let's do the reverse: find the extrema of the new and see how*/
+ /* close they get to the bounding box of the old */
+ /* And get really unhappy if the spline extends beyond the end-points */
+ db2.min = 0; db2.max = db->len;
+ SplineFindExtrema(&spline->splines[0],&ts[0],&ts[1]);
+ for ( i=0; i<2; ++i ) {
+ if ( ts[i]!=-1 ) {
+ x = ((spline->splines[0].a*ts[i]+spline->splines[0].b)*ts[i]+spline->splines[0].c)*ts[i] + spline->splines[0].d;
+ y = ((spline->splines[1].a*ts[i]+spline->splines[1].b)*ts[i]+spline->splines[1].c)*ts[i] + spline->splines[1].d;
+ if ( x<b->minx )
+ sum += (x-b->minx)*(x-b->minx);
+ else if ( x>b->maxx )
+ sum += (x-b->maxx)*(x-b->maxx);
+ dot = (x-db->base.x)*db->unit.x + (y-db->base.y)*db->unit.y;
+ if ( dot<db2.min ) db2.min = dot;
+ if ( dot>db2.max ) db2.max = dot;
+ }
+ }
+ SplineFindExtrema(&spline->splines[1],&ts[0],&ts[1]);
+ for ( i=0; i<2; ++i ) {
+ if ( ts[i]!=-1 ) {
+ x = ((spline->splines[0].a*ts[i]+spline->splines[0].b)*ts[i]+spline->splines[0].c)*ts[i] + spline->splines[0].d;
+ y = ((spline->splines[1].a*ts[i]+spline->splines[1].b)*ts[i]+spline->splines[1].c)*ts[i] + spline->splines[1].d;
+ if ( y<b->miny )
+ sum += (y-b->miny)*(y-b->miny);
+ else if ( y>b->maxy )
+ sum += (y-b->maxy)*(y-b->maxy);
+ dot = (x-db->base.x)*db->unit.x + (y-db->base.y)*db->unit.y;
+ if ( dot<db2.min ) db2.min = dot;
+ if ( dot>db2.max ) db2.max = dot;
+ }
+ }
+
+ /* Big penalty for going beyond the range of the desired spline */
+ if ( db->min==0 && db2.min<0 )
+ sum += 10000 + db2.min*db2.min;
+ else if ( db2.min<db->min )
+ sum += 100 + (db2.min-db->min)*(db2.min-db->min);
+ if ( db->max==db->len && db2.max>db->len )
+ sum += 10000 + (db2.max-db->max)*(db2.max-db->max);
+ else if ( db2.max>db->max )
+ sum += 100 + (db2.max-db->max)*(db2.max-db->max);
+
+return( sum );
+}
+
+static void ApproxBounds(DBounds *b, FitPoint *mid, int cnt, struct dotbounds *db) {
+ int i;
+ bigreal dot;
+
+ b->minx = b->maxx = mid[0].p.x;
+ b->miny = b->maxy = mid[0].p.y;
+ db->min = 0; db->max = db->len;
+ for ( i=1; i<cnt; ++i ) {
+ if ( mid[i].p.x>b->maxx ) b->maxx = mid[i].p.x;
+ if ( mid[i].p.x<b->minx ) b->minx = mid[i].p.x;
+ if ( mid[i].p.y>b->maxy ) b->maxy = mid[i].p.y;
+ if ( mid[i].p.y<b->miny ) b->miny = mid[i].p.y;
+ dot = (mid[i].p.x-db->base.x)*db->unit.x + (mid[i].p.y-db->base.y)*db->unit.y;
+ if ( dot<db->min ) db->min = dot;
+ if ( dot>db->max ) db->max = dot;
+ }
+}
+
+/* pf == point from (start point) */
+/* Δf == slope from (cp(from) - from) */
+/* pt == point to (end point, t==1) */
+/* Δt == slope to (cp(to) - to) */
+
+/* A spline from pf to pt with slope vectors rf*Δf, rt*Δt is: */
+/* p(t) = pf + [ 3*rf*Δf ]*t + 3*[pt-pf+rt*Δt-2*rf*Δf] *t^2 + */
+/* [2*pf-2*pt+3*rf*Δf-3*rt*Δt]*t^3 */
+
+/* So I want */
+/* d ÎŁ (p(t(i))-p(i))^2/ d rf == 0 */
+/* d ÎŁ (p(t(i))-p(i))^2/ d rt == 0 */
+/* now... */
+/* d ÎŁ (p(t(i))-p(i))^2/ d rf == 0 */
+/* => Σ 3*t*Δf*(1-2*t+t^2)*
+ * [pf-pi+ 3*(pt-pf)*t^2 + 2*(pf-pt)*t^3] +
+ * 3*[t - 2*t^2 + t^3]*Δf*rf +
+ * 3*[t^2-t^3]*Δt*rt */
+/* and... */
+/* d ÎŁ (p(t(i))-p(i))^2/ d rt == 0 */
+/* => Σ 3*t^2*Δt*(1-t)*
+ * [pf-pi+ 3*(pt-pf)*t^2 + 2*(pf-pt)*t^3] +
+ * 3*[t - 2*t^2 + t^3]*Δf*rf +
+ * 3*[t^2-t^3]*Δt*rt */
+
+/* Now for a long time I looked at that and saw four equations and two unknowns*/
+/* That was I was trying to solve for x and y separately, and that doesn't work. */
+/* There are really just two equations and each sums over both x and y components */
+
+/* Old comment: */
+/* I used to do a least squares approach adding two more to the above set of equations */
+/* which held the slopes constant. But that didn't work very well. So instead*/
+/* Then I tried doing the approximation, and then forcing the control points */
+/* to be in line (with the original slopes), getting a better approximation */
+/* to "t" for each data point and then calculating an error array, approximating*/
+/* it, and using that to fix up the final result */
+/* Then I tried checking various possible cp lengths in the desired directions*/
+/* finding the best one or two, and doing a 2D binary search using that as a */
+/* starting point. */
+/* And sometimes a least squares approach will give us the right answer, so */
+/* try that too. */
+/* This still isn't as good as I'd like it... But I haven't been able to */
+/* improve it further yet */
+/* The mergetype mt is either of: */
+/* mt_matrix; original, fast, all-purpose (relies on matrix calculations) */
+/* mt_levien; by Raph Levien (implemented by Linus Romer), fast, accurate, use only if mid is on spline */
+/* mt_bruteforce; slow, all-purpose, normally more accurate than mt_matrix.*/
+/* The mt_levien algorithm is explained here: */
+/* raphlinus.github.io/curves/2021/03/11/bezier-fitting.html */
+/* The notation used here is a bit different: Instead of theta1, theta2, */
+/* delta1, delta2, momentx, area we use alpha,beta,a,b,m,f: */
+/* Here is to complete math that we are using: */
+/* Signed area of the cubic bezier spline a .. controls b and c .. d to the x-axis */
+/* f = ((xb-xa)*(10*ya+6*yb+3*yc+yd)+(xc-xb)*(4*ya+6*yb+6*yc+4*yd)+(xd-xc)*(ya+3*yb+6*yc+10*yd))/20; */
+/* simplified for the normed case */
+/* f = 3/20*(2*a*sin(alpha)+2*b*sin(beta)-a*b*sin(alpha+beta)); */
+/* solved for b */
+/* b = (20*f-6*a*sin(alpha))/(6*sin(beta)-3*a*sin(alpha+beta)). */
+/* Signed area of the cubic bezier spline a .. controls b and c .. d to the x-axis */
+/* from point a up to the bezier point at time t */
+/* f(t) = ((((1-t)*xa+xb*t)-xa)*(10*ya+6*((1-t)*ya+yb*t)+3*((1-t)^2*ya+2*(1-t)*t*yb+t^2*yc) */
+/* +((1-t)^3*ya+3*(1-t)^2*t*yb+3*(1-t)*t^2*yc+t^3*yd))+(((1-t)^2*xa+2*(1-t)*t*xb+t^2*xc) */
+/* -((1-t)*xa+xb*t))*(4*ya+6*((1-t)*ya+yb*t)+6*((1-t)^2*ya+2*(1-t)*t*yb+t^2*yc) */
+/* +4*((1-t)^3*ya+3*(1-t)^2*t*yb+3*(1-t)*t^2*yc+t^3*yd))+(((1-t)^3*xa+3*(1-t)^2*t*xb */
+/* +3*(1-t)*t^2*xc+t^3*xd)-((1-t)^2*xa+2*(1-t)*t*xb+t^2*xc))*(ya+3*((1-t)*ya+yb*t) */
+/* +6*((1-t)^2*ya+2*(1-t)*t*yb+t^2*yc)+10*((1-t)^3*ya+3*(1-t)^2*t*yb+3*(1-t)*t^2*yc+t^3*yd)))/20; */
+/* simplified for the normed case: */
+/* f(t) = -(3*(30*a*b*sin(beta-alpha)*t^6+15*b^2*sin(2*beta)*t^6-20*b*sin(beta)*t^6 */
+/* -15*a^2*sin(2*alpha)*t^6+20*a*sin(alpha)*t^6+6*a*b*sin(beta+alpha)*t^5 */
+/* -90*a*b*sin(beta-alpha)*t^5-30*b^2*sin(2*beta)*t^5+48*b*sin(beta)*t^5 */
+/* +60*a^2*sin(2*alpha)*t^5-72*a*sin(alpha)*t^5-15*a*b*sin(beta+alpha)*t^4 */
+/* +90*a*b*sin(beta-alpha)*t^4+15*b^2*sin(2*beta)*t^4-30*b*sin(beta)*t^4 */
+/* -90*a^2*sin(2*alpha)*t^4+90*a*sin(alpha)*t^4+10*a*b*sin(beta+alpha)*t^3 */
+/* -30*a*b*sin(beta-alpha)*t^3+60*a^2*sin(2*alpha)*t^3-40*a*sin(alpha)*t^3 */
+/* -15*a^2*sin(2*alpha)*t^2))/20. */
+/* First moment about y-axis = \int x dA = \int x dA/dt dt for a cubic bezier */
+/* path a .. controls b and c .. d */
+/* m = (280*xd^2*yd-105*xc*xd*yd-30*xb*xd*yd-5*xa*xd*yd-45*xc^2*yd-45*xb*xc*yd */
+/* -12*xa*xc*yd-18*xb^2*yd-15*xa*xb*yd-5*xa^2*yd+105*xd^2*yc+45*xc*xd*yc */
+/* -3*xa*xd*yc-27*xb*xc*yc-18*xa*xc*yc-27*xb^2*yc-45*xa*xb*yc-30*xa^2*yc */
+/* +30*xd^2*yb+45*xc*xd*yb+18*xb*xd*yb+3*xa*xd*yb+27*xc^2*yb+27*xb*xc*yb */
+/* -45*xa*xb*yb-105*xa^2*yb+5*xd^2*ya+15*xc*xd*ya+12*xb*xd*ya+5*xa*xd*ya */
+/* +18*xc^2*ya+45*xb*xc*ya+30*xa*xc*ya+45*xb^2*ya+105*xa*xb*ya-280*xa^2*ya)/840; */
+/* simplified for the normed case */
+/* m = (9*a*cos(alpha)*b^2*cos(beta)*sin(beta)-15*b^2*cos(beta)*sin(beta) */
+/* -9*a^2*cos(alpha)^2*b*sin(beta)-9*a*cos(alpha)*b*sin(beta)+50*b*sin(beta) */
+/* +9*a*sin(alpha)*b^2*cos(beta)^2-9*a^2*cos(alpha)*sin(alpha)*b*cos(beta) */
+/* -33*a*sin(alpha)*b*cos(beta)+15*a^2*cos(alpha)*sin(alpha)+34*a*sin(alpha))/280; */
+/* normed case combined with the formula for b depending on the area (see above): */
+/* m = (34*a*sin(alpha)+50*(20*f-6*a*sin(alpha))/(6*sin(beta)-3*a*sin(beta+alpha))*sin(beta) */
+/* +15*a^2*sin(alpha)*cos(alpha)-15*(20*f-6*a*sin(alpha))/(6*sin(beta) */
+/* -3*a*sin(beta+alpha))^2*sin(beta)*cos(beta)-a*(20*f-6*a*sin(alpha))/(6*sin(beta) */
+/* -3*a*sin(beta+alpha))*(33*sin(alpha)*cos(beta)+9*cos(alpha)*sin(beta)) */
+/* -9*a^2*(20*f-6*a*sin(alpha))/(6*sin(beta)-3*a*sin(beta+alpha))*sin(alpha+beta)*cos(alpha) */
+/* +9*a*(20*f-6*a*sin(alpha))/(6*sin(beta)-3*a*sin(beta+alpha))^2*sin(alpha+beta)*cos(beta))/280; */
+/* and reduced to a quartic equation with sa = sin(alpha), sb = sin(beta), ca = cos(alpha), cb = cos(beta) */
+/* 0 = -9*ca*(((2*sb*cb*ca+sa*(2*cb*cb-1))*ca-2*sb*cb)*ca-cb*cb*sa) * a^4 */
+/* + 12*((((cb*(30*f*cb-sb)-15*f)*ca+2*sa-cb*sa*(cb+30*f*sb))*ca+cb*(sb-15*f*cb))*ca-sa*cb*cb) * a^3 */
+/* + 12*((((70*m+15*f)*sb^2+cb*(9*sb-70*cb*m-5*cb*f))*ca-5*sa*sb*(3*sb-4*cb*(7*m+f)))*ca-cb*(9*sb-70*cb*m-5*cb*f)) * a^2 */
+/* + 16*(((12*sa-5*ca*(42*m-17*f))*sb-70*cb*(3*m-f)*sa-75*ca*cb*f*f)*sb-75*cb^2*f^2*sa) * a */
+/* + 80*sb*(42*sb*m-25*f*(sb-cb*f)); */
+/* this quartic equation reduces to a quadratic for the special case beta = pi - alpha or beta = -alpha */
+/* 0 = -9*ca*sa^2 * a^3 */
+/* + 6*sa*(4*sa+5*ca*f) * a^2 */
+/* + 10*((42*m-25*f)*sa-25*ca*f^2). */
+/* The derivative of the first moment (not the quartic) = 0 results in a quartic as well: */
+/* 0 = -9*ca*sa*sab^3 * a^4 */
+/* -3*sab^2*(9*ca*sa*sb-(17*sa+30*ca*f)*sab+15*cb*sa^2) * a^3 */
+/* +18*sab*sb*(21*ca*sa*sb-(17*sa+30*ca*f)*sab+15*cb*sa^2) * a^2 */
+/* -4*(144*ca*sa*sb^3+((-78*sa-135*ca*f)*sab+108*cb*sa^2)*sb^2+(-125*f*sab^2-45*cb*f*sa*sab)*sb+150*cb*f^2*sab^2) * a */
+/* +8*sb*((24*sa+45*ca*f)*sb^2+(15*cb*f*sa-125*f*sab)*sb+100*cb*f^2*sab) */
+/* this quartic equation reduces to a linear for the special case beta = pi - alpha or beta = -alpha */
+/* 0 = -3*ca*sa * a */
+/* +4*sa+5*ca*f */
+#define TRY_CNT 2
+#define DECIMATION 5
+Spline *ApproximateSplineFromPointsSlopes(SplinePoint *from, SplinePoint *to,
+ FitPoint *mid, int cnt, int order2, enum mergetype mt) {
+ BasePoint tounit, fromunit, ftunit;
+ bigreal flen,tlen,ftlen,dot;
+ Spline *spline, temp;
+ BasePoint nextcp;
+ int bettern, betterp;
+ bigreal offn, offp, incrn, incrp, trylen;
+ int nocnt = 0, totcnt;
+ bigreal curdiff, bestdiff[TRY_CNT];
+ int i,j,besti[TRY_CNT],bestj[TRY_CNT],k,l;
+ bigreal fdiff, tdiff, fmax, tmax, fdotft, tdotft;
+ DBounds b;
+ struct dotbounds db;
+ bigreal offn_, offp_, finaldiff;
+ bigreal pt_pf_x, pt_pf_y, determinant;
+ bigreal consts[2], rt_terms[2], rf_terms[2];
+
+ /* If all the selected points are at the same spot, and one of the */
+ /* end-points is also at that spot, then just copy the control point */
+ /* But our caller seems to have done that for us */
+
+ /* If the two end-points are corner points then allow the slope to vary */
+ /* Or if one end-point is a tangent but the point defining the tangent's */
+ /* slope is being removed then allow the slope to vary */
+ /* Except if the slope is horizontal or vertical then keep it fixed */
+ if ( ( !from->nonextcp && ( from->nextcp.x==from->me.x || from->nextcp.y==from->me.y)) ||
+ (!to->noprevcp && ( to->prevcp.x==to->me.x || to->prevcp.y==to->me.y)) )
+ /* Preserve the slope */;
+ else if ( ((from->pointtype == pt_corner && from->nonextcp) ||
+ (from->pointtype == pt_tangent &&
+ ((from->nonextcp && from->noprevcp) || !from->noprevcp))) &&
+ ((to->pointtype == pt_corner && to->noprevcp) ||
+ (to->pointtype == pt_tangent &&
+ ((to->nonextcp && to->noprevcp) || !to->nonextcp))) ) {
+ from->pointtype = to->pointtype = pt_corner;
+return( ApproximateSplineFromPoints(from,to,mid,cnt,order2) );
+ }
+
+ /* If we are going to honour the slopes of a quadratic spline, there is */
+ /* only one possibility */
+ if ( order2 ) {
+ if ( from->nonextcp )
+ from->nextcp = from->next->to->me;
+ if ( to->noprevcp )
+ to->prevcp = to->prev->from->me;
+ fromunit.x = from->nextcp.x-from->me.x; fromunit.y = from->nextcp.y-from->me.y;
+ tounit.x = to->prevcp.x-to->me.x; tounit.y = to->prevcp.y-to->me.y;
+
+ if ( !IntersectLines(&nextcp,&from->nextcp,&from->me,&to->prevcp,&to->me) ||
+ (nextcp.x-from->me.x)*fromunit.x + (nextcp.y-from->me.y)*fromunit.y < 0 ||
+ (nextcp.x-to->me.x)*tounit.x + (nextcp.y-to->me.y)*tounit.y < 0 ) {
+ /* If the slopes don't intersect then use a line */
+ /* (or if the intersection is patently absurd) */
+ from->nextcp = from->me;
+ to->prevcp = to->me;
+ TestForLinear(from,to);
+ } else {
+ from->nextcp = to->prevcp = nextcp;
+ }
+return( SplineMake2(from,to));
+ }
+ /* From here down we are only working with cubic splines */
+
+ if ( cnt==0 ) {
+ /* Just use what we've got, no data to improve it */
+ /* But we do sometimes get some cps which are too big */
+ bigreal len = sqrt((to->me.x-from->me.x)*(to->me.x-from->me.x) + (to->me.y-from->me.y)*(to->me.y-from->me.y));
+ if ( len==0 ) {
+ from->nextcp = from->me;
+ to->prevcp = to->me;
+ } else {
+ BasePoint noff, poff;
+ bigreal nlen, plen;
+ noff.x = from->nextcp.x-from->me.x; noff.y = from->nextcp.y-from->me.y;
+ poff.x = to->me.x-to->prevcp.x; poff.y = to->me.y-to->prevcp.y;
+ nlen = sqrt(noff.x*noff.x + noff.y+noff.y);
+ plen = sqrt(poff.x*poff.x + poff.y+poff.y);
+ if ( nlen>len/3 ) {
+ noff.x *= len/3/nlen; noff.y *= len/3/nlen;
+ from->nextcp.x = from->me.x + noff.x;
+ from->nextcp.y = from->me.y + noff.y;
+ }
+ if ( plen>len/3 ) {
+ poff.x *= len/3/plen; poff.y *= len/3/plen;
+ to->prevcp.x = to->me.x + poff.x;
+ to->prevcp.y = to->me.y + poff.y;
+ }
+ }
+return( SplineMake3(from,to));
+ }
+
+ if ( to->prev!=NULL && (( to->noprevcp && to->nonextcp ) || to->prev->knownlinear )) {
+ tounit.x = to->prev->from->me.x-to->me.x; tounit.y = to->prev->from->me.y-to->me.y;
+// g_message("1 tu: %f %f, %d %d %d", tounit.x, tounit.y, to->noprevcp, to->nonextcp, to->prev->knownlinear);
+ } else if ( !to->noprevcp || to->pointtype == pt_corner ) {
+ tounit.x = to->prevcp.x-to->me.x; tounit.y = to->prevcp.y-to->me.y;
+// g_message("2 tu: %f %f", tounit.x, tounit.y);
+ } else {
+ tounit.x = to->me.x-to->nextcp.x; tounit.y = to->me.y-to->nextcp.y;
+// g_message("3 tu: %f %f", tounit.x, tounit.y);
+ }
+ tlen = sqrt(tounit.x*tounit.x + tounit.y*tounit.y);
+ if ( from->next!=NULL && (( from->noprevcp && from->nonextcp ) || from->next->knownlinear) ) {
+ fromunit.x = from->next->to->me.x-from->me.x; fromunit.y = from->next->to->me.y-from->me.y;
+// g_message("1 fu");
+ } else if ( !from->nonextcp || from->pointtype == pt_corner ) {
+ fromunit.x = from->nextcp.x-from->me.x; fromunit.y = from->nextcp.y-from->me.y;
+// g_message("2 fu");
+ } else {
+ fromunit.x = from->me.x-from->prevcp.x; fromunit.y = from->me.y-from->prevcp.y;
+// g_message("3 fu");
+ }
+ flen = sqrt(fromunit.x*fromunit.x + fromunit.y*fromunit.y);
+// g_message("tu: %f %f, fu: %f, %f", tounit.x, tounit.y, fromunit.x, fromunit.y);
+ if ( (tlen==0 || flen==0) && (from->next==NULL || to->prev==NULL) ) {
+ memset(&temp,0,sizeof(temp));
+ temp.from = from; temp.to = to;
+ SplineRefigure(&temp);
+ from->next = to->prev = NULL;
+ }
+ if ( tlen==0 ) {
+ if ( to->prev!=NULL ) {
+ temp = *to->prev;
+ }
+ if ( (to->pointtype==pt_curve || to->pointtype==pt_hvcurve) &&
+ to->next && !to->nonextcp ) {
+ tounit.x = to->me.x-to->nextcp.x; tounit.y = to->me.y-to->nextcp.y;
+ } else {
+ tounit.x = -( (3*temp.splines[0].a*.9999+2*temp.splines[0].b)*.9999+temp.splines[0].c );
+ tounit.y = -( (3*temp.splines[1].a*.9999+2*temp.splines[1].b)*.9999+temp.splines[1].c );
+ }
+ tlen = sqrt(tounit.x*tounit.x + tounit.y*tounit.y);
+ }
+ tounit.x /= tlen; tounit.y /= tlen;
+
+ if ( flen==0 ) {
+ if ( from->next!=NULL ) {
+ temp = *from->next;
+ }
+ if ( (from->pointtype==pt_curve || from->pointtype==pt_hvcurve) &&
+ from->prev && !from->noprevcp ) {
+ fromunit.x = from->me.x-from->prevcp.x; fromunit.y = from->me.y-from->prevcp.y;
+ } else {
+ fromunit.x = ( (3*temp.splines[0].a*.0001+2*temp.splines[0].b)*.0001+temp.splines[0].c );
+ fromunit.y = ( (3*temp.splines[1].a*.0001+2*temp.splines[1].b)*.0001+temp.splines[1].c );
+ }
+ flen = sqrt(fromunit.x*fromunit.x + fromunit.y*fromunit.y);
+ }
+ fromunit.x /= flen; fromunit.y /= flen;
+
+ ftunit.x = (to->me.x-from->me.x); ftunit.y = (to->me.y-from->me.y);
+ ftlen = sqrt(ftunit.x*ftunit.x + ftunit.y*ftunit.y);
+ if ( ftlen!=0 ) {
+ ftunit.x /= ftlen; ftunit.y /= ftlen;
+ }
+
+ if ( (dot=fromunit.x*tounit.y - fromunit.y*tounit.x)<.0001 && dot>-.0001 &&
+ (dot=ftunit.x*tounit.y - ftunit.y*tounit.x)<.0001 && dot>-.0001 ) {
+ /* It's a line. Slopes are parallel, and parallel to vector between (from,to) */
+ from->nextcp = from->me; to->prevcp = to->me;
+return( SplineMake3(from,to));
+ }
+ /* This is the generic case, where a generic part is approximated by a cubic */
+ /* bezier spline. */
+ if ( ( ftlen == 0 ) && ( mt != mt_matrix ) )
+ mt = mt_matrix;
+ if ( mt == mt_levien ) {
+ bigreal f,m,xa,ya,xb,yb,xc,yc,xd,yd,sasa,sab;
+ int numberOfSolutions;
+ SplinePoint *frompoint,*topoint;
+ f = 0; /* area */
+ m = 0; /* first area moment about y (along x) */
+ frompoint = from;
+ if ( from->next==NULL )
+ topoint=to;
+ else
+ topoint=from->next->to;
+ for ( ; ; frompoint = topoint->next->from, topoint = topoint->next->to ) {
+ /* normalizing transformation (chord to x axis and length 1) */
+ xa = ((frompoint->me.x-from->me.x)*ftunit.x+(frompoint->me.y-from->me.y)*ftunit.y)/ftlen;
+ ya = (-(frompoint->me.x-from->me.x)*ftunit.y+(frompoint->me.y-from->me.y)*ftunit.x)/ftlen;
+ xb = ((frompoint->nextcp.x-from->me.x)*ftunit.x+(frompoint->nextcp.y-from->me.y)*ftunit.y)/ftlen;
+ yb = (-(frompoint->nextcp.x-from->me.x)*ftunit.y+(frompoint->nextcp.y-from->me.y)*ftunit.x)/ftlen;
+ xc = ((topoint->prevcp.x-from->me.x)*ftunit.x+(topoint->prevcp.y-from->me.y)*ftunit.y)/ftlen;
+ yc = (-(topoint->prevcp.x-from->me.x)*ftunit.y+(topoint->prevcp.y-from->me.y)*ftunit.x)/ftlen;
+ xd = ((topoint->me.x-from->me.x)*ftunit.x+(topoint->me.y-from->me.y)*ftunit.y)/ftlen;
+ yd = (-(topoint->me.x-from->me.x)*ftunit.y+(topoint->me.y-from->me.y)*ftunit.x)/ftlen;
+ f += ((xb-xa)*(10*ya+6*yb+3*yc+yd)+(xc-xb)*(4*ya+6*yb+6*yc+4*yd)+(xd-xc)*(ya+3*yb+6*yc+10*yd))/20;
+ m += (280*xd*xd*yd-105*xc*xd*yd-30*xb*xd*yd-5*xa*xd*yd-45*xc*xc*yd-45*xb*xc*yd-12*xa*xc*yd-18*xb*xb*yd
+ -15*xa*xb*yd-5*xa*xa*yd+105*xd*xd*yc+45*xc*xd*yc-3*xa*xd*yc-27*xb*xc*yc-18*xa*xc*yc-27*xb*xb*yc
+ -45*xa*xb*yc-30*xa*xa*yc+30*xd*xd*yb+45*xc*xd*yb+18*xb*xd*yb+3*xa*xd*yb+27*xc*xc*yb+27*xb*xc*yb
+ -45*xa*xb*yb-105*xa*xa*yb+5*xd*xd*ya+15*xc*xd*ya+12*xb*xd*ya+5*xa*xd*ya+18*xc*xc*ya+45*xb*xc*ya
+ +30*xa*xc*ya+45*xb*xb*ya+105*xa*xb*ya-280*xa*xa*ya)/840;
+ if ( topoint==to )
+ break;
+ }
+ BasePoint aunit = (BasePoint) { BPDot(ftunit, fromunit), BPCross(ftunit, fromunit) }; /* normed direction at "from" */
+ BasePoint bunit = (BasePoint) { BPDot(BPRev(ftunit), tounit),BPCross(ftunit, tounit) }; /* normed direction at "to" */
+ if ( aunit.y < 0 ) { /* normalize aunit.y to >= 0: */
+ aunit.y = -aunit.y;
+ bunit.y = -bunit.y;
+ m = -m;
+ f = -f;
+ }
+ /* calculate the Tunni point (where the tangents at "from" and "to" intersect) */
+ bigreal aMax = 100; /* maximal value that the handle a can reach up to the Tunni point, 100 is really long */
+ bigreal bMax = 100; /* maximal value that the handle b can reach up to the Tunni point, 100 is really long */
+ sab = aunit.y*bunit.x+aunit.x*bunit.y;
+ if (sab != 0) { /* if handles not parallel */
+ aMax = bunit.y/sab;
+ bMax = aunit.y/sab;
+ if ( aMax < 0 ) {
+ aMax = 100;
+ }
+ if ( bMax < 0 ) {
+ bMax = 100;
+ }
+ }
+ /* start approximation by solving the quartic equation */
+ sasa = aunit.y*aunit.y; /* reducing the multiplications */
+ Quartic quad;
+ if ( (aunit.x == -bunit.x && aunit.y == bunit.y) || (aunit.x == bunit.x && aunit.y == -bunit.y) ) { /* handles are parallel */
+ quad.a = 0;
+ quad.b = 0;
+ quad.c = -9*aunit.x*sasa;
+ quad.d = 6*aunit.y*(4*aunit.y+5*aunit.x*f);
+ quad.e = 10*((42*m-25*f)*aunit.y-25*aunit.x*f*f);
+ } else { /* generic situation */
+ quad.a = -9*aunit.x*(((2*bunit.y*bunit.x*aunit.x+aunit.y
+ *(2*bunit.x*bunit.x-1))*aunit.x-2*bunit.y*bunit.x)
+ *aunit.x-bunit.x*bunit.x*aunit.y);
+ quad.b = 12*((((bunit.x*(30*f*bunit.x-bunit.y)-15*f)
+ *aunit.x+2*aunit.y-bunit.x*aunit.y*(bunit.x+30*f*bunit.y))
+ *aunit.x+bunit.x*(bunit.y-15*f*bunit.x))
+ *aunit.x-aunit.y*bunit.x*bunit.x);
+ quad.c = 12*((((70*m+15*f)*bunit.y*bunit.y+bunit.x
+ *(9*bunit.y-70*bunit.x*m-5*bunit.x*f))
+ *aunit.x-5*aunit.y*bunit.y*(3*bunit.y-4*bunit.x
+ *(7*m+f)))*aunit.x-bunit.x*(9*bunit.y-70*bunit.x*m-5*bunit.x*f));
+ quad.d = 16*(((12*aunit.y-5*aunit.x*(42*m-17*f))*bunit.y
+ -70*bunit.x*(3*m-f)*aunit.y-75*aunit.x*bunit.x*f*f)
+ *bunit.y-75*bunit.x*bunit.x*f*f*aunit.y);
+ quad.e = 80*bunit.y*(42*bunit.y*m-25*f*(bunit.y-bunit.x*f));
+ }
+ extended solutions[4] = {-999999,-999999,-999999,-999999};
+ _QuarticSolve(&quad,solutions);
+ extended abSolutions[10][2]; /* there are at most 4+4+1+1=10 solutions of pairs of a and b (quartic=0,derivative=0,b=0.01,a=0.01) */
+ numberOfSolutions = 0;
+ extended a,b;
+ for( int i = 0; i < 4; i++ ){
+ a = solutions[i];
+ if ( a >= 0 && a < aMax ) {
+ b = (20*f-6*a*aunit.y)/(3*(2*bunit.y-a*sab));
+ if ( b >= 0 && b < bMax ) {
+ abSolutions[numberOfSolutions][0] = a;
+ abSolutions[numberOfSolutions++][1] = b;
+ }
+ }
+ }
+ /* and now again for the derivative (of m not of the upper quartic): */
+ if ( (aunit.x == -bunit.x && aunit.y == bunit.y) || (aunit.x == bunit.x && aunit.y == -bunit.y) ) { /* handles are parallel */
+ quad.a = 0;
+ quad.b = 0;
+ quad.c = 0;
+ quad.d = -3*aunit.x*aunit.y;
+ quad.e = 4*aunit.y+5*aunit.x*f;
+ } else { /* generic situation */
+ bigreal sbsb = bunit.y*bunit.y;
+ bigreal sabsab = sab*sab;
+ quad.a = -9*aunit.x*aunit.y*sabsab*sab;
+ quad.b = -3*sabsab*(9*aunit.x*aunit.y*bunit.y-(17*aunit.y
+ +30*aunit.x*f)*sab+15*bunit.x*sasa);
+ quad.c = 18*sab*bunit.y*(21*aunit.x*aunit.y*bunit.y-(17*aunit.y
+ +30*aunit.x*f)*sab+15*bunit.x*sasa);
+ quad.d = -4*(144*aunit.x*aunit.y*sbsb*bunit.y+((-78*aunit.y-
+ 135*aunit.x*f)*sab+108*bunit.x*sasa)*sbsb+(-125*f*sabsab
+ -45*bunit.x*f*aunit.y*sab)*bunit.y+150*bunit.x*f*f*sabsab);
+ quad.e = 8*bunit.y*((24*aunit.y+45*aunit.x*f)*sbsb
+ +(15*bunit.x*f*aunit.y-125*f*sab)*bunit.y+100*bunit.x*f*f*sab);
+ }
+ for( int i = 0; i < 4; i++ ) /* overwriting (reusing) */
+ solutions[i] = -999999;
+ _QuarticSolve(&quad,solutions);
+ for( int i = 0; i < 4; i++ ){
+ a = solutions[i];
+ if ( a >= 0 && a < aMax ) {
+ b = (20*f-6*a*aunit.y)/(3*(2*bunit.y-a*sab));
+ if ( b >= 0 && b < bMax ) {
+ abSolutions[numberOfSolutions][0] = a;
+ abSolutions[numberOfSolutions++][1] = b;
+ }
+ }
+ }
+ /* Add the solution of b = 0.01 (approximately 0 but above because of direction). */
+#if 0 /* those solutions lead to subpar approximations */
+ /* This solution is not part of the original algorithm by Raph Levien. */
+ a = (2000*f-6*bunit.y)/(600*aunit.y-3*sab);
+ if ( a >= 0 && a < aMax ) {
+ abSolutions[numberOfSolutions][0] = a;
+ abSolutions[numberOfSolutions++][1] = 0.01;
+ }
+ /* Add the solution of a = 0.01 (approximately 0 but above because of direction). */
+ /* This solution is not part of the original algorithm by Raph Levien. */
+ b = (2000*f-6*aunit.y)/(600*bunit.y-3*sab);
+ if ( b >= 0 && b < bMax ) {
+ abSolutions[numberOfSolutions][0] = 0.01;
+ abSolutions[numberOfSolutions++][1] = b;
+ }
+ if ( numberOfSolutions == 0) { /* add solutions that extend up to the Tunni point */
+ /* try solution with a = aMax and b area-equal*/
+ b = (20*f-6*aMax*aunit.y)/(3*(2*bunit.y-aMax*sab));
+ if ( b >= 0 && b < bMax ) {
+ abSolutions[numberOfSolutions][0] = aMax;
+ abSolutions[numberOfSolutions++][1] = b;
+ }
+ /* try solution with b = bMax and a area-equal*/
+ a = (20*f-6*bMax*bunit.y)/(3*(2*aunit.y-bMax*sab));
+ if ( a >= 0 && a < aMax ) {
+ abSolutions[numberOfSolutions][0] = a;
+ abSolutions[numberOfSolutions++][1] = bMax;
+ }
+ }
+#endif
+ if ( numberOfSolutions == 0) {
+ /* no solutions found, quit */
+return NULL;
+
+ /* solution with a = aMax and b = bMax*/
+ abSolutions[numberOfSolutions][0] = aMax;
+ abSolutions[numberOfSolutions++][1] = bMax;
+ }
+ if ( numberOfSolutions == 1) {
+ from->nextcp.x = from->me.x+ftlen*fromunit.x*abSolutions[0][0];
+ from->nextcp.y = from->me.y+ftlen*fromunit.y*abSolutions[0][0];
+ to->prevcp.x = to->me.x+ftlen*tounit.x*abSolutions[0][1];
+ to->prevcp.y = to->me.y+ftlen*tounit.y*abSolutions[0][1];
+ } else { /* compare L2 errors to choose the best solution */
+ bigreal bestError = 1e30;
+ bigreal t,error,errorsum,dist;
+ BasePoint prevcp,nextcp,coeff1,coeff2,coeff3;
+ int last_best_j;
+ for (int k=0; k<numberOfSolutions; k++) {
+ nextcp.x = from->me.x+ftlen*fromunit.x*abSolutions[k][0];
+ nextcp.y = from->me.y+ftlen*fromunit.y*abSolutions[k][0];
+ prevcp.x = to->me.x+ftlen*tounit.x*abSolutions[k][1];
+ prevcp.y = to->me.y+ftlen*tounit.y*abSolutions[k][1];
+ /* Calculate the error of the cubic bezier path from,nextcp,prevcp,to: */
+ /* In order to do that we calculate 99 points on the bezier path. */
+ coeff3.x = -from->me.x+3*nextcp.x-3*prevcp.x+to->me.x;
+ coeff3.y = -from->me.y+3*nextcp.y-3*prevcp.y+to->me.y;
+ coeff2.x = 3*from->me.x-6*nextcp.x+3*prevcp.x;
+ coeff2.y = 3*from->me.y-6*nextcp.y+3*prevcp.y;
+ coeff1.x = -3*from->me.x+3*nextcp.x;
+ coeff1.y = -3*from->me.y+3*nextcp.y;
+ BasePoint approx[99];
+ for (int i=0; i<99; i++) {
+ t = (i+1)/100.0;
+ approx[i].x = from->me.x+t*(coeff1.x+t*(coeff2.x+t*coeff3.x));
+ approx[i].y = from->me.y+t*(coeff1.y+t*(coeff2.y+t*coeff3.y));
+ }
+ /* Now we calculate the error by determining the minimal quadratic distance to the mid points. */
+ errorsum = 0.0;
+ last_best_j = 0;
+ for (int i=0; i<cnt; i++) { /* Going through the mid points */
+ error = 1e30;
+ /* For this mid point, find the distance to the closest one of the */
+ /* 99 points on the approximate cubic bezier. */
+ /* To not favour approximations which trace the original multiple times */
+ /* by going back and forth, only consider monotonic mappings. */
+ /* I.e., start from the point that was closest to the previous mid point: */
+ for (int j=last_best_j; j<99; j++) {
+ dist = (mid[i].p.x-approx[j].x)*(mid[i].p.x-approx[j].x)
+ +(mid[i].p.y-approx[j].y)*(mid[i].p.y-approx[j].y);
+ if (dist < error) {
+ error = dist;
+ last_best_j = j;
+ }
+ }
+ errorsum += error;
+ if (errorsum > bestError)
+ break;
+ }
+ if (errorsum < bestError) {
+ bestError = errorsum;
+ from->nextcp = nextcp;
+ to->prevcp = prevcp;
+ }
+ }
+ }
+ return( SplineMake3(from,to));
+ } else if ( mt == mt_bruteforce ) {
+ bigreal best_error = 1e30;
+ bigreal t,error,errorsum,dist;
+ BasePoint prevcp,coeff1,coeff2,coeff3;
+ bigreal best_fromhandle = 0.0;
+ bigreal best_tohandle = 0.0;
+ BasePoint approx[99]; /* The 99 points on the approximate cubic bezier */
+ /* We make 2 runs: The first run to narrow the variation range, the second run to finetune */
+ /* The optimal length of the two handles are determined by brute force. */
+ for (int run=0; run<2; ++run) {
+ for (int fromhandle=((run==0)?1:-29); fromhandle<=((run==0)?60:29); ++fromhandle) {
+ for (int tohandle=((run==0)?1:-29); tohandle<=((run==0)?60:29); ++tohandle) {
+ nextcp.x = from->me.x+ftlen*fromunit.x*( (run==0)?fromhandle:best_fromhandle+fromhandle/30.0 )/60.0;
+ nextcp.y = from->me.y+ftlen*fromunit.y*( (run==0)?fromhandle:best_fromhandle+fromhandle/30.0 )/60.0;
+ prevcp.x = to->me.x+ftlen*tounit.x*( (run==0)?tohandle:best_tohandle+tohandle/30.0 )/60.0;
+ prevcp.y = to->me.y+ftlen*tounit.y*( (run==0)?tohandle:best_tohandle+tohandle/30.0 )/60.0;
+ /* Calculate the error of the cubic bezier path from,nextcp,prevcp,to: */
+ /* In order to do that we calculate 99 points on the bezier path. */
+ coeff3.x = -from->me.x+3*nextcp.x-3*prevcp.x+to->me.x;
+ coeff3.y = -from->me.y+3*nextcp.y-3*prevcp.y+to->me.y;
+ coeff2.x = 3*from->me.x-6*nextcp.x+3*prevcp.x;
+ coeff2.y = 3*from->me.y-6*nextcp.y+3*prevcp.y;
+ coeff1.x = -3*from->me.x+3*nextcp.x;
+ coeff1.y = -3*from->me.y+3*nextcp.y;
+ for (int i=0; i<99; ++i) {
+ t = (i+1)/100.0;
+ approx[i].x = from->me.x+t*(coeff1.x+t*(coeff2.x+t*coeff3.x));
+ approx[i].y = from->me.y+t*(coeff1.y+t*(coeff2.y+t*coeff3.y));
+ }
+ /* Now we calculate the error by determining the minimal quadratic distance to the mid points. */
+ errorsum = 0.0;
+ for (int i=0; i<cnt; ++i) { /* Going through the mid points */
+ error = (mid[i].p.x-approx[0].x)*(mid[i].p.x-approx[0].x)
+ +(mid[i].p.y-approx[0].y)*(mid[i].p.y-approx[0].y);
+ /* Above we have just initialized the error and */
+ /* now we are going through the remaining 98 of */
+ /* 99 points on the approximate cubic bezier: */
+ for (int j=1; j<99; ++j) {
+ dist = (mid[i].p.x-approx[j].x)*(mid[i].p.x-approx[j].x)
+ +(mid[i].p.y-approx[j].y)*(mid[i].p.y-approx[j].y);
+ if (dist < error)
+ error = dist;
+ }
+ errorsum += error;
+ if (errorsum > best_error)
+ break;
+ }
+ if (errorsum < best_error) {
+ best_error = errorsum;
+ if (run == 0) {
+ best_fromhandle = fromhandle;
+ best_tohandle = tohandle;
+ }
+ from->nextcp = nextcp;
+ to->prevcp = prevcp;
+ }
+ }
+ }
+ }
+ return( SplineMake3(from,to));
+ }
+ else { /* mergetype mt_matrix (original algorithm) */
+ pt_pf_x = to->me.x - from->me.x;
+ pt_pf_y = to->me.y - from->me.y;
+ consts[0] = consts[1] = rt_terms[0] = rt_terms[1] = rf_terms[0] = rf_terms[1] = 0;
+ for ( i=0; i<cnt; ++i ) {
+ bigreal t = mid[i].t, t2 = t*t, t3=t2*t;
+ bigreal factor_from = t-2*t2+t3;
+ bigreal factor_to = t2-t3;
+ bigreal const_x = from->me.x-mid[i].p.x + 3*pt_pf_x*t2 - 2*pt_pf_x*t3;
+ bigreal const_y = from->me.y-mid[i].p.y + 3*pt_pf_y*t2 - 2*pt_pf_y*t3;
+ bigreal temp1 = 3*(t-2*t2+t3);
+ bigreal rf_term_x = temp1*fromunit.x;
+ bigreal rf_term_y = temp1*fromunit.y;
+ bigreal temp2 = 3*(t2-t3);
+ bigreal rt_term_x = -temp2*tounit.x;
+ bigreal rt_term_y = -temp2*tounit.y;
+
+ consts[0] += factor_from*( fromunit.x*const_x + fromunit.y*const_y );
+ consts[1] += factor_to *( -tounit.x*const_x + -tounit.y*const_y);
+ rf_terms[0] += factor_from*( fromunit.x*rf_term_x + fromunit.y*rf_term_y);
+ rf_terms[1] += factor_to*( -tounit.x*rf_term_x + -tounit.y*rf_term_y);
+ rt_terms[0] += factor_from*( fromunit.x*rt_term_x + fromunit.y*rt_term_y);
+ rt_terms[1] += factor_to*( -tounit.x*rt_term_x + -tounit.y*rt_term_y);
+ }
+
+ /* I've only seen singular matrices (determinant==0) when cnt==1 */
+ /* but even with cnt==1 the determinant is usually non-0 (16 times out of 17)*/
+ determinant = (rt_terms[0]*rf_terms[1]-rt_terms[1]*rf_terms[0]);
+ if ( determinant!=0 ) {
+ bigreal rt, rf;
+ rt = (consts[1]*rf_terms[0]-consts[0]*rf_terms[1])/determinant;
+ if ( rf_terms[0]!=0 )
+ rf = -(consts[0]+rt*rt_terms[0])/rf_terms[0];
+ else /* if ( rf_terms[1]!=0 ) This can't happen, otherwise the determinant would be 0 */
+ rf = -(consts[1]+rt*rt_terms[1])/rf_terms[1];
+ /* If we get bad values (ones that point diametrically opposed to what*/
+ /* we need), then fix that factor at 0, and see what we get for the */
+ /* other */
+ if ( rf>=0 && rt>0 && rf_terms[0]!=0 &&
+ (rf = -consts[0]/rf_terms[0])>0 ) {
+ rt = 0;
+ } else if ( rf<0 && rt<=0 && rt_terms[1]!=0 &&
+ (rt = -consts[1]/rt_terms[1])<0 ) {
+ rf = 0;
+ }
+ if ( rt<=0 && rf>=0 ) {
+ from->nextcp.x = from->me.x + rf*fromunit.x;
+ from->nextcp.y = from->me.y + rf*fromunit.y;
+ to->prevcp.x = to->me.x - rt*tounit.x;
+ to->prevcp.y = to->me.y - rt*tounit.y;
+return( SplineMake3(from,to));
+ }
+ }
+
+ trylen = (to->me.x-from->me.x)*fromunit.x + (to->me.y-from->me.y)*fromunit.y;
+ if ( trylen>flen ) flen = trylen;
+
+ trylen = (from->me.x-to->me.x)*tounit.x + (from->me.y-to->me.y)*tounit.y;
+ if ( trylen>tlen ) tlen = trylen;
+
+ for ( i=0; i<cnt; ++i ) {
+ trylen = (mid[i].p.x-from->me.x)*fromunit.x + (mid[i].p.y-from->me.y)*fromunit.y;
+ if ( trylen>flen ) flen = trylen;
+ trylen = (mid[i].p.x-to->me.x)*tounit.x + (mid[i].p.y-to->me.y)*tounit.y;
+ if ( trylen>tlen ) tlen = trylen;
+ }
+
+ fdotft = fromunit.x*ftunit.x + fromunit.y*ftunit.y;
+ fmax = fdotft>0 ? ftlen/fdotft : 1e10;
+ tdotft = -tounit.x*ftunit.x - tounit.y*ftunit.y;
+ tmax = tdotft>0 ? ftlen/tdotft : 1e10;
+ /* At fmax, tmax the control points will stretch beyond the other endpoint*/
+ /* when projected along the line between the two endpoints */
+
+ db.base = from->me;
+ db.unit = ftunit;
+ db.len = ftlen;
+ ApproxBounds(&b,mid,cnt,&db);
+
+ for ( k=0; k<TRY_CNT; ++k ) {
+ bestdiff[k] = 1e20;
+ besti[k] = -1; bestj[k] = -1;
+ }
+ fdiff = flen/DECIMATION;
+ tdiff = tlen/DECIMATION;
+ from->nextcp = from->me;
+ memset(&temp,0,sizeof(Spline));
+ temp.from = from; temp.to = to;
+ for ( i=1; i<DECIMATION; ++i ) {
+ from->nextcp.x += fdiff*fromunit.x; from->nextcp.y += fdiff*fromunit.y;
+ to->prevcp = to->me;
+ for ( j=1; j<DECIMATION; ++j ) {
+ to->prevcp.x += tdiff*tounit.x; to->prevcp.y += tdiff*tounit.y;
+ SplineRefigure(&temp);
+ curdiff = SigmaDeltas(&temp,mid,cnt,&b,&db);
+ for ( k=0; k<TRY_CNT; ++k ) {
+ if ( curdiff<bestdiff[k] ) {
+ for ( l=TRY_CNT-1; l>k; --l ) {
+ bestdiff[l] = bestdiff[l-1];
+ besti[l] = besti[l-1];
+ bestj[l] = bestj[l-1];
+ }
+ bestdiff[k] = curdiff;
+ besti[k] = i; bestj[k]=j;
+ break;
+ }
+ }
+ }
+ }
+
+ finaldiff = 1e20;
+ offn_ = offp_ = -1;
+ spline = SplineMake(from,to,false);
+ for ( k=-1; k<TRY_CNT; ++k ) {
+ if ( k<0 ) {
+ BasePoint nextcp, prevcp;
+ bigreal temp1, temp2;
+ int ret = _ApproximateSplineFromPoints(from,to,mid,cnt,&nextcp,&prevcp,false);
+ /* sometimes least squares gives us the right answer */
+ if ( !(ret&1) || !(ret&2))
+ continue;
+ temp1 = (prevcp.x-to->me.x)*tounit.x + (prevcp.y-to->me.y)*tounit.y;
+ temp2 = (nextcp.x-from->me.x)*fromunit.x + (nextcp.y-from->me.y)*fromunit.y;
+ if ( temp1<=0 || temp2<=0 ) /* A nice solution... but the control points are diametrically opposed to what they should be */
+ continue;
+ tlen = temp1; flen = temp2;
+ } else {
+ if ( bestj[k]<0 || besti[k]<0 )
+ continue;
+ tlen = bestj[k]*tdiff; flen = besti[k]*fdiff;
+ }
+ to->prevcp.x = to->me.x + tlen*tounit.x; to->prevcp.y = to->me.y + tlen*tounit.y;
+ from->nextcp.x = from->me.x + flen*fromunit.x; from->nextcp.y = from->me.y + flen*fromunit.y;
+ SplineRefigure(spline);
+
+ bettern = betterp = false;
+ incrn = tdiff/2.0; incrp = fdiff/2.0;
+ offn = flen; offp = tlen;
+ nocnt = 0;
+ curdiff = SigmaDeltas(spline,mid,cnt,&b,&db);
+ totcnt = 0;
+ for (;;) {
+ bigreal fadiff, fsdiff;
+ bigreal tadiff, tsdiff;
+
+ from->nextcp.x = from->me.x + (offn+incrn)*fromunit.x; from->nextcp.y = from->me.y + (offn+incrn)*fromunit.y;
+ to->prevcp.x = to->me.x + offp*tounit.x; to->prevcp.y = to->me.y + offp*tounit.y;
+ SplineRefigure(spline);
+ fadiff = SigmaDeltas(spline,mid,cnt,&b,&db);
+ from->nextcp.x = from->me.x + (offn-incrn)*fromunit.x; from->nextcp.y = from->me.y + (offn-incrn)*fromunit.y;
+ SplineRefigure(spline);
+ fsdiff = SigmaDeltas(spline,mid,cnt,&b,&db);
+ from->nextcp.x = from->me.x + offn*fromunit.x; from->nextcp.y = from->me.y + offn*fromunit.y;
+ if ( offn-incrn<=0 )
+ fsdiff += 1e10;
+
+ to->prevcp.x = to->me.x + (offp+incrp)*tounit.x; to->prevcp.y = to->me.y + (offp+incrp)*tounit.y;
+ SplineRefigure(spline);
+ tadiff = SigmaDeltas(spline,mid,cnt,&b,&db);
+ to->prevcp.x = to->me.x + (offp-incrp)*tounit.x; to->prevcp.y = to->me.y + (offp-incrp)*tounit.y;
+ SplineRefigure(spline);
+ tsdiff = SigmaDeltas(spline,mid,cnt,&b,&db);
+ to->prevcp.x = to->me.x + offp*tounit.x; to->prevcp.y = to->me.y + offp*tounit.y;
+ if ( offp-incrp<=0 )
+ tsdiff += 1e10;
+
+ if ( offn>=incrn && fsdiff<curdiff &&
+ (fsdiff<fadiff && fsdiff<tsdiff && fsdiff<tadiff)) {
+ offn -= incrn;
+ if ( bettern>0 )
+ incrn /= 2;
+ bettern = -1;
+ nocnt = 0;
+ curdiff = fsdiff;
+ } else if ( offn+incrn<fmax && fadiff<curdiff &&
+ (fadiff<=fsdiff && fadiff<tsdiff && fadiff<tadiff)) {
+ offn += incrn;
+ if ( bettern<0 )
+ incrn /= 2;
+ bettern = 1;
+ nocnt = 0;
+ curdiff = fadiff;
+ } else if ( offp>=incrp && tsdiff<curdiff &&
+ (tsdiff<=fsdiff && tsdiff<=fadiff && tsdiff<tadiff)) {
+ offp -= incrp;
+ if ( betterp>0 )
+ incrp /= 2;
+ betterp = -1;
+ nocnt = 0;
+ curdiff = tsdiff;
+ } else if ( offp+incrp<tmax && tadiff<curdiff &&
+ (tadiff<=fsdiff && tadiff<=fadiff && tadiff<=tsdiff)) {
+ offp += incrp;
+ if ( betterp<0 )
+ incrp /= 2;
+ betterp = 1;
+ nocnt = 0;
+ curdiff = tadiff;
+ } else {
+ if ( ++nocnt > 6 )
+ break;
+ incrn /= 2;
+ incrp /= 2;
+ }
+ if ( curdiff<1 )
+ break;
+ if ( incrp<tdiff/1024 || incrn<fdiff/1024 )
+ break;
+ if ( ++totcnt>200 )
+ break;
+ if ( offn<0 || offp<0 ) {
+ IError("Approximation got inverse control points");
+ break;
+ }
+ }
+ if ( curdiff<finaldiff ) {
+ finaldiff = curdiff;
+ offn_ = offn;
+ offp_ = offp;
+ }
+ }
+
+ to->noprevcp = offp_==0;
+ from->nonextcp = offn_==0;
+ to->prevcp.x = to->me.x + offp_*tounit.x; to->prevcp.y = to->me.y + offp_*tounit.y;
+ from->nextcp.x = from->me.x + offn_*fromunit.x; from->nextcp.y = from->me.y + offn_*fromunit.y;
+ /* I used to check for a spline very close to linear (and if so, make it */
+ /* linear). But in when stroking a path with an elliptical pen we transform*/
+ /* the coordinate system and our normal definitions of "close to linear" */
+ /* don't apply */
+ /*TestForLinear(from,to);*/
+ SplineRefigure(spline);
+
+return( spline );
+ }
+}
+#undef TRY_CNT
+#undef DECIMATION
+
+SplinePoint *_ApproximateSplineSetFromGen(SplinePoint *from, SplinePoint *to,
+ bigreal start_t, bigreal end_t,
+ bigreal toler, int toler_is_sumsq,
+ GenPointsP genp, void *tok,
+ int order2, int depth) {
+ int cnt, i, maxerri=0, created = false;
+ bigreal errsum=0, maxerr=0, d, mid_t;
+ FitPoint *fp;
+ SplinePoint *mid, *r;
+
+ cnt = (*genp)(tok, start_t, end_t, &fp);
+ if ( cnt < 2 )
+ return NULL;
+
+ // Rescale zero to one
+ for ( i=1; i<(cnt-1); ++i )
+ fp[i].t = (fp[i].t-fp[0].t)/(fp[cnt-1].t-fp[0].t);
+ fp[0].t = 0.0;
+ fp[cnt-1].t = 1.0;
+
+ from->nextcp.x = from->me.x + fp[0].ut.x;
+ from->nextcp.y = from->me.y + fp[0].ut.y;
+ from->nonextcp = false;
+ if ( to!=NULL )
+ to->me = fp[cnt-1].p;
+ else {
+ to = SplinePointCreate(fp[cnt-1].p.x, fp[cnt-1].p.y);
+ created = true;
+ }
+ to->prevcp.x = to->me.x - fp[cnt-1].ut.x;
+ to->prevcp.y = to->me.y - fp[cnt-1].ut.y;
+ to->noprevcp = false;
+ ApproximateSplineFromPointsSlopes(from,to,fp+1,cnt-2,order2,mt_matrix);
+
+ for ( i=0; i<cnt; ++i ) {
+ d = SplineMinDistanceToPoint(from->next, &fp[i].p);
+ errsum += d*d;
+ if ( d>maxerr ) {
+ maxerr = d;
+ maxerri = i;
+ }
+ }
+ // printf(" Error sum %lf, max error %lf at depth %d\n", errsum, maxerr, depth);
+
+ if ( (toler_is_sumsq ? errsum : maxerr) > toler && depth < 6 ) {
+ mid_t = fp[maxerri].t * (end_t-start_t) + start_t;
+ free(fp);
+ SplineFree(from->next);
+ from->next = NULL;
+ to->prev = NULL;
+ mid = _ApproximateSplineSetFromGen(from, NULL, start_t, mid_t, toler,
+ toler_is_sumsq, genp, tok, order2,
+ depth+1);
+ if ( mid ) {
+ r = _ApproximateSplineSetFromGen(mid, to, mid_t, end_t, toler,
+ toler_is_sumsq, genp, tok,
+ order2, depth+1);
+ if ( r )
+ return r;
+ else {
+ if ( created )
+ SplinePointFree(to);
+ else
+ to->prev = NULL;
+ SplinePointFree(mid);
+ SplineFree(from->next);
+ from->next = NULL;
+ return NULL;
+ }
+ } else {
+ if ( created )
+ SplinePointFree(to);
+ return NULL;
+ }
+ } else if ( (toler_is_sumsq ? errsum : maxerr) > toler ) {
+ TRACE("%s %lf exceeds %lf at maximum depth %d\n",
+ toler_is_sumsq ? "Sum of squared errors" : "Maximum error length",
+ toler_is_sumsq ? errsum : maxerr, toler, depth);
+ }
+ free(fp);
+ return to;
+}
+
+SplinePoint *ApproximateSplineSetFromGen(SplinePoint *from, SplinePoint *to,
+ bigreal start_t, bigreal end_t,
+ bigreal toler, int toler_is_sumsq,
+ GenPointsP genp, void *tok,
+ int order2) {
+ return _ApproximateSplineSetFromGen(from, to, start_t, end_t, toler,
+ toler_is_sumsq, genp, tok, order2, 0);
+}
diff --git a/src/path/splinefit/splinefit.h b/src/path/splinefit/splinefit.h
new file mode 100644
index 0000000..22e21fc
--- /dev/null
+++ b/src/path/splinefit/splinefit.h
@@ -0,0 +1,78 @@
+// SPDX-License-Identifier: GPL-2.0-or-later
+/* Copyright (C) 2000-2012 by George Williams, 2019 by Skef Iterum, 2021 by Linus Romer */
+/*
+ * Redistribution and use in source and binary forms, with or without
+ * modification, are permitted provided that the following conditions are met:
+
+ * Redistributions of source code must retain the above copyright notice, this
+ * list of conditions and the following disclaimer.
+
+ * Redistributions in binary form must reproduce the above copyright notice,
+ * this list of conditions and the following disclaimer in the documentation
+ * and/or other materials provided with the distribution.
+
+ * The name of the author may not be used to endorse or promote products
+ * derived from this software without specific prior written permission.
+
+ * THIS SOFTWARE IS PROVIDED BY THE AUTHOR ``AS IS'' AND ANY EXPRESS OR IMPLIED
+ * WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF
+ * MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO
+ * EVENT SHALL THE AUTHOR BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
+ * SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
+ * PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS;
+ * OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY,
+ * WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR
+ * OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF
+ * ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+ */
+
+#ifndef FONTFORGE_SPLINEFIT_H
+#define FONTFORGE_SPLINEFIT_H
+
+// #include <fontforge-config.h>
+#include "splinefont.h"
+
+typedef struct fitpoint {
+ BasePoint p;
+ BasePoint ut;
+ bigreal t;
+} FitPoint;
+
+#define FITPOINT_EMPTY { {0.0, 0.0}, {0.0, 0.0}, 0.0 }
+
+enum mergetype { mt_matrix, mt_levien, mt_bruteforce };
+
+extern Spline *ApproximateSplineFromPoints(SplinePoint *from, SplinePoint *to,
+ FitPoint *mid, int cnt, int order2);
+extern Spline *ApproximateSplineFromPointsSlopes(SplinePoint *from, SplinePoint *to,
+ FitPoint *mid, int cnt, int order2, enum mergetype mt);
+
+/* ApproximateSplineSetFromGen() fits a one or more splines to data
+ * generated by calls to genp, within the tolerance toler. The data
+ * are treated as noiseless although some noise may be harmless. When
+ * an interval cannot be fit with a single spline it is divided at
+ * the point of largest error, where the position and slope are fixed
+ * according to the FitPoint entries.
+ *
+ * A GenPointsP function should allocate a free-able array of FitPoints,
+ * set *fpp to point to it, and return the number of points allocated.
+ * The points should correspond to values chosen between the interval
+ * t_start and t_end, but just what this interval corresponds to depends
+ * on vinfo and is therefore opaque to the tracing system. It is up to the
+ * implementation how many points to return on the interval and how they are
+ * "t-spaced".
+ *
+ * The function can also return 0 points, in which case *fpp should be set
+ * to NULL. This will cancel the tracing. This semantic is useful for
+ * communicating information found during tracing back to the caller
+ * of ApproximateSplineSetFromGen(), presumably via the vinfo structure,
+ * so that it can retrace with a different choice of interval.
+ */
+typedef int (*GenPointsP)(void *vinfo, bigreal t_start, bigreal t_end, FitPoint **fpp);
+
+extern SplinePoint *ApproximateSplineSetFromGen(SplinePoint *from, SplinePoint *to,
+ bigreal start_t, bigreal end_t,
+ bigreal toler, int toler_is_sumsq,
+ GenPointsP genp, void *tok, int order2);
+
+#endif /* FONTFORGE_SPLINEFIT_H */
diff --git a/src/path/splinefit/splinefont.c b/src/path/splinefit/splinefont.c
new file mode 100644
index 0000000..4681628
--- /dev/null
+++ b/src/path/splinefit/splinefont.c
@@ -0,0 +1,1174 @@
+// SPDX-License-Identifier: GPL-2.0-or-later
+
+#include <math.h>
+#include <stdbool.h>
+#include <stdint.h>
+#include "splinefont.h"
+#include "splinefit.h"
+
+#define FONTFORGE_CONFIG_USE_DOUBLE 1
+
+bigreal BPDot(BasePoint v1, BasePoint v2) {
+ return v1.x * v2.x + v1.y * v2.y;
+}
+
+bigreal BPCross(BasePoint v1, BasePoint v2) {
+ return v1.x * v2.y - v1.y * v2.x;
+}
+
+BasePoint BPRev(BasePoint v) {
+ return (BasePoint) { -v.x, -v.y };
+}
+
+int RealWithin(real a,real b,real fudge) {
+ return( b>=a-fudge && b<=a+fudge );
+}
+
+BOOL RealNear(real a,real b) {
+ real d = a-b;
+#ifdef FONTFORGE_CONFIG_USE_DOUBLE
+ // These tighter equals-zero tests are retained for code tuned when
+ // passing zero as a constant
+ if ( a==0 )
+ return b>-1e-8 && b<1e-8;
+ if ( b==0 )
+ return a>-1e-8 && a<1e-8;
+
+ return d>-1e-6 && d<1e-6;
+#else /* For floats */
+ return d>-1e-5 && d<1e-5
+#endif
+}
+
+int RealApprox(real a,real b) {
+
+ if ( a==0 ) {
+ if ( b<.0001 && b>-.0001 )
+return( true );
+ } else if ( b==0 ) {
+ if ( a<.0001 && a>-.0001 )
+return( true );
+ } else {
+ a /= b;
+ if ( a>=.95 && a<=1.05 )
+return( true );
+ }
+return( false );
+}
+
+void LineListFree(LineList *ll) {
+ LineList *next;
+
+ while ( ll!=NULL ) {
+ next = ll->next;
+ chunkfree(ll,sizeof(LineList));
+ ll = next;
+ }
+}
+
+void LinearApproxFree(LinearApprox *la) {
+ LinearApprox *next;
+
+ while ( la!=NULL ) {
+ next = la->next;
+ LineListFree(la->lines);
+ chunkfree(la,sizeof(LinearApprox));
+ la = next;
+ }
+}
+
+void SplineFree(Spline *spline) {
+ LinearApproxFree(spline->approx);
+ chunkfree(spline,sizeof(Spline));
+}
+
+SplinePoint *SplinePointCreate(real x, real y) {
+ SplinePoint *sp;
+ if ( (sp=chunkalloc(sizeof(SplinePoint)))!=NULL ) {
+ sp->me.x = x; sp->me.y = y;
+ sp->nextcp = sp->prevcp = sp->me;
+ sp->nonextcp = sp->noprevcp = true;
+ sp->nextcpdef = sp->prevcpdef = false;
+ sp->ttfindex = sp->nextcpindex = 0xfffe;
+ sp->name = NULL;
+ }
+ return( sp );
+}
+
+void SplinePointsFree(SplinePointList *spl) {
+ Spline *first, *spline, *next;
+ int nonext;
+
+ if ( spl==NULL )
+ return;
+ if ( spl->first!=NULL ) {
+ nonext = spl->first->next==NULL; // If there is no spline, we set a flag.
+ first = NULL;
+ // We start on the first spline if it exists.
+ for ( spline = spl->first->next; spline!=NULL && spline!=first; spline = next ) {
+ next = spline->to->next; // Cache the location of the next spline.
+ SplinePointFree(spline->to); // Free the destination point.
+ SplineFree(spline); // Free the spline.
+ if ( first==NULL ) first = spline; // We want to avoid repeating the circuit.
+ }
+ // If the path is open or has no splines, free the starting point.
+ if ( spl->last!=spl->first || nonext )
+ SplinePointFree(spl->first);
+ }
+}
+
+void SplinePointListFree(SplinePointList *spl) {
+
+ if ( spl==NULL ) return;
+ SplinePointsFree(spl);
+ // free(spl->spiros);
+ free(spl->contour_name);
+ chunkfree(spl,sizeof(SplinePointList));
+}
+
+void SplineRefigure2(Spline *spline) {
+ SplinePoint *from = spline->from, *to = spline->to;
+ Spline1D *xsp = &spline->splines[0], *ysp = &spline->splines[1];
+ Spline old;
+
+#ifdef DEBUG
+ if ( RealNear(from->me.x,to->me.x) && RealNear(from->me.y,to->me.y))
+ IError("Zero length spline created");
+#endif
+ if ( spline->acceptableextrema )
+ old = *spline;
+
+ if ( ( from->nextcp.x==from->me.x && from->nextcp.y==from->me.y && from->nextcpindex>=0xfffe )
+ || ( to->prevcp.x==to->me.x && to->prevcp.y==to->me.y && from->nextcpindex>=0xfffe ) ) {
+ from->nonextcp = to->noprevcp = true;
+ from->nextcp = from->me;
+ to->prevcp = to->me;
+ } else {
+ from->nonextcp = to->noprevcp = false;
+ if ( from->nextcp.x==from->me.x && from->nextcp.y==from->me.y )
+ to->prevcp = from->me;
+ else if ( to->prevcp.x==to->me.x && to->prevcp.y==to->me.y )
+ from->nextcp = to->me;
+ }
+
+ if ( from->nonextcp && to->noprevcp )
+ /* Ok */;
+ else if ( from->nextcp.x!=to->prevcp.x || from->nextcp.y!=to->prevcp.y ) {
+ if ( RealNear(from->nextcp.x,to->prevcp.x) &&
+ RealNear(from->nextcp.y,to->prevcp.y)) {
+ from->nextcp.x = to->prevcp.x = (from->nextcp.x+to->prevcp.x)/2;
+ from->nextcp.y = to->prevcp.y = (from->nextcp.y+to->prevcp.y)/2;
+ } else {
+ IError("Invalid 2nd order spline in SplineRefigure2" );
+#ifndef GWW_TEST
+ /* I don't want these to go away when I'm debugging. I want to */
+ /* know how I got them */
+ from->nextcp.x = to->prevcp.x = (from->nextcp.x+to->prevcp.x)/2;
+ from->nextcp.y = to->prevcp.y = (from->nextcp.y+to->prevcp.y)/2;
+#endif
+ }
+ }
+
+ xsp->d = from->me.x; ysp->d = from->me.y;
+ if ( from->nonextcp && to->noprevcp ) {
+ spline->islinear = true;
+ xsp->c = to->me.x-from->me.x;
+ ysp->c = to->me.y-from->me.y;
+ xsp->a = xsp->b = 0;
+ ysp->a = ysp->b = 0;
+ } else {
+ /* from p. 393 (Operator Details, curveto) PostScript Lang. Ref. Man. (Red book) */
+ xsp->c = 2*(from->nextcp.x-from->me.x);
+ ysp->c = 2*(from->nextcp.y-from->me.y);
+ xsp->b = to->me.x-from->me.x-xsp->c;
+ ysp->b = to->me.y-from->me.y-ysp->c;
+ xsp->a = 0;
+ ysp->a = 0;
+ if ( RealNear(xsp->c,0)) xsp->c=0;
+ if ( RealNear(ysp->c,0)) ysp->c=0;
+ if ( RealNear(xsp->b,0)) xsp->b=0;
+ if ( RealNear(ysp->b,0)) ysp->b=0;
+ spline->islinear = false;
+ if ( ysp->b==0 && xsp->b==0 )
+ spline->islinear = true; /* This seems extremely unlikely... */
+ if ( from->nextcpselected || to->prevcpselected ) {
+ // The convention for tracking selection of quadratic control
+ // points is to use nextcpselected except at the tail of the
+ // list, where it's prevcpselected on the first point.
+ from->nextcpselected = true;
+ to->prevcpselected = false;
+ }
+ }
+ if ( isnan(ysp->b) || isnan(xsp->b) )
+ IError("NaN value in spline creation");
+ LinearApproxFree(spline->approx);
+ spline->approx = NULL;
+ spline->knowncurved = false;
+ spline->knownlinear = spline->islinear;
+ SplineIsLinear(spline);
+ spline->isquadratic = !spline->knownlinear;
+ spline->order2 = true;
+
+ if ( spline->acceptableextrema ) {
+ /* I don't check "d", because changes to that reflect simple */
+ /* translations which will not affect the shape of the spline */
+ /* (I don't check "a" because it is always 0 in a quadratic spline) */
+ if ( !RealNear(old.splines[0].b,spline->splines[0].b) ||
+ !RealNear(old.splines[0].c,spline->splines[0].c) ||
+ !RealNear(old.splines[1].b,spline->splines[1].b) ||
+ !RealNear(old.splines[1].c,spline->splines[1].c) )
+ spline->acceptableextrema = false;
+ }
+}
+
+Spline *SplineMake(SplinePoint *from, SplinePoint *to, int order2) {
+ if (order2 > 0)
+return( SplineMake2(from,to));
+ else
+return( SplineMake3(from,to));
+}
+
+Spline *SplineMake2(SplinePoint *from, SplinePoint *to) {
+ Spline *spline = chunkalloc(sizeof(Spline));
+
+ spline->from = from; spline->to = to;
+ from->next = to->prev = spline;
+ spline->order2 = true;
+ SplineRefigure2(spline);
+return( spline );
+}
+
+Spline *SplineMake3(SplinePoint *from, SplinePoint *to) {
+ Spline *spline = chunkalloc(sizeof(Spline));
+
+ spline->from = from; spline->to = to;
+ from->next = to->prev = spline;
+ SplineRefigure3(spline);
+return( spline );
+}
+
+void SplinePointFree(SplinePoint *sp) {
+ // chunkfree(sp->hintmask,sizeof(HintMask));
+ free(sp->name);
+ chunkfree(sp,sizeof(SplinePoint));
+}
+
+void SplineRefigure(Spline *spline) {
+ if ( spline==NULL )
+return;
+ if ( spline->order2 )
+ SplineRefigure2(spline);
+ else
+ SplineRefigure3(spline);
+}
+
+# define RE_NearZero .00000001
+# define RE_Factor (1024.0*1024.0*1024.0*1024.0*1024.0*2.0) /* 52 bits => divide by 2^51 */
+
+int Within16RoundingErrors(bigreal v1, bigreal v2) {
+ bigreal temp=v1*v2;
+ bigreal re;
+
+ if ( temp<0 ) /* Ok, if the two values are on different sides of 0 there */
+return( false ); /* is no way they can be within a rounding error of each other */
+ else if ( temp==0 ) {
+ if ( v1==0 )
+return( v2<RE_NearZero && v2>-RE_NearZero );
+ else
+return( v1<RE_NearZero && v1>-RE_NearZero );
+ } else if ( v1>0 ) {
+ if ( v1>v2 ) { /* Rounding error from the biggest absolute value */
+ re = v1/ (RE_Factor/16);
+return( v1-v2 < re );
+ } else {
+ re = v2/ (RE_Factor/16);
+return( v2-v1 < re );
+ }
+ } else {
+ if ( v1<v2 ) {
+ re = v1/ (RE_Factor/16); /* This will be a negative number */
+return( v1-v2 > re );
+ } else {
+ re = v2/ (RE_Factor/16);
+return( v2-v1 > re );
+ }
+ }
+}
+
+/* An IEEE double has 52 bits of precision. So one unit of rounding error will be */
+/* the number divided by 2^51 */
+# define D_RE_Factor (1024.0*1024.0*1024.0*1024.0*1024.0*2.0)
+/* But that's not going to work near 0, so, since the t values we care about */
+/* are [0,1], let's use 1.0/D_RE_Factor */
+
+double CheckExtremaForSingleBitErrors(const Spline1D *sp, double t, double othert) {
+ double u1, um1;
+ double slope, slope1, slopem1;
+ int err;
+ double diff, factor;
+
+ if ( t<0 || t>1 )
+return( t );
+
+ factor = t*0x40000/D_RE_Factor;
+ if ( (diff = t-othert)<0 ) diff= -diff;
+ if ( factor>diff/4 && diff!=0 ) /* This little check is to insure we don't skip beyond the well of this extremum into the next */
+ factor = diff/4;
+
+ slope = (3*(double) sp->a*t+2*sp->b)*t+sp->c;
+ if ( slope<0 ) slope = -slope;
+
+ for ( err = 0x40000; err!=0; err>>=1 ) {
+ u1 = t+factor;
+ slope1 = (3*(double) sp->a*u1+2*sp->b)*u1+sp->c;
+ if ( slope1<0 ) slope1 = -slope1;
+
+ um1 = t-factor;
+ slopem1 = (3*(double) sp->a*um1+2*sp->b)*um1+sp->c;
+ if ( slopem1<0 ) slopem1 = -slopem1;
+
+ if ( slope1<slope && slope1<=slopem1 && u1<=1.0 ) {
+ t = u1;
+ } else if ( slopem1<slope && slopem1<=slope1 && um1>=0.0 ) {
+ t = um1;
+ }
+ factor /= 2.0;
+ }
+ /* that seems as good as it gets */
+
+return( t );
+}
+
+void SplineFindExtrema(const Spline1D *sp, extended *_t1, extended *_t2 ) {
+ extended t1= -1, t2= -1;
+ extended b2_fourac;
+
+ /* Find the extreme points on the curve */
+ /* Set to -1 if there are none or if they are outside the range [0,1] */
+ /* Order them so that t1<t2 */
+ /* If only one valid extremum it will be t1 */
+ /* (Does not check the end points unless they have derivative==0) */
+ /* (Does not check to see if d/dt==0 points are inflection points (rather than extrema) */
+ if ( sp->a!=0 ) {
+ /* cubic, possibly 2 extrema (possibly none) */
+ b2_fourac = 4*(extended) sp->b*sp->b - 12*(extended) sp->a*sp->c;
+ if ( b2_fourac>=0 ) {
+ b2_fourac = sqrt(b2_fourac);
+ t1 = (-2*sp->b - b2_fourac) / (6*sp->a);
+ t2 = (-2*sp->b + b2_fourac) / (6*sp->a);
+ t1 = CheckExtremaForSingleBitErrors(sp,t1,t2);
+ t2 = CheckExtremaForSingleBitErrors(sp,t2,t1);
+ if ( t1>t2 ) { extended temp = t1; t1 = t2; t2 = temp; }
+ else if ( t1==t2 ) t2 = -1;
+ if ( RealNear(t1,0)) t1=0; else if ( RealNear(t1,1)) t1=1;
+ if ( RealNear(t2,0)) t2=0; else if ( RealNear(t2,1)) t2=1;
+ if ( t2<=0 || t2>=1 ) t2 = -1;
+ if ( t1<=0 || t1>=1 ) { t1 = t2; t2 = -1; }
+ }
+ } else if ( sp->b!=0 ) {
+ /* Quadratic, at most one extremum */
+ t1 = -sp->c/(2.0*(extended) sp->b);
+ if ( t1<=0 || t1>=1 ) t1 = -1;
+ } else /*if ( sp->c!=0 )*/ {
+ /* linear, no extrema */
+ }
+ *_t1 = t1; *_t2 = t2;
+}
+
+int IntersectLines(BasePoint *inter,
+ BasePoint *line1_1, BasePoint *line1_2,
+ BasePoint *line2_1, BasePoint *line2_2) {
+ // A lot of functions call this with the same address as an input and the output.
+ // In order to avoid unexpected behavior, we delay writing to the output until the end.
+ bigreal s1, s2;
+ BasePoint _output;
+ BasePoint * output = &_output;
+ if ( line1_1->x == line1_2->x ) {
+ // Line 1 is vertical.
+ output->x = line1_1->x;
+ if ( line2_1->x == line2_2->x ) {
+ // Line 2 is vertical.
+ if ( line2_1->x!=line1_1->x )
+ return( false ); /* Parallel vertical lines */
+ output->y = (line1_1->y+line2_1->y)/2;
+ } else {
+ output->y = line2_1->y + (output->x-line2_1->x) * (line2_2->y - line2_1->y)/(line2_2->x - line2_1->x);
+ }
+ *inter = *output;
+ return( true );
+ } else if ( line2_1->x == line2_2->x ) {
+ // Line 2 is vertical, but we know that line 1 is not.
+ output->x = line2_1->x;
+ output->y = line1_1->y + (output->x-line1_1->x) * (line1_2->y - line1_1->y)/(line1_2->x - line1_1->x);
+ *inter = *output;
+ return( true );
+ } else {
+ // Both lines are oblique.
+ s1 = (line1_2->y - line1_1->y)/(line1_2->x - line1_1->x);
+ s2 = (line2_2->y - line2_1->y)/(line2_2->x - line2_1->x);
+ if ( RealNear(s1,s2)) {
+ if ( !RealNear(line1_1->y + (line2_1->x-line1_1->x) * s1,line2_1->y))
+ return( false );
+ output->x = (line1_2->x+line2_2->x)/2;
+ output->y = (line1_2->y+line2_2->y)/2;
+ } else {
+ output->x = (s1*line1_1->x - s2*line2_1->x - line1_1->y + line2_1->y)/(s1-s2);
+ output->y = line1_1->y + (output->x-line1_1->x) * s1;
+ }
+ *inter = *output;
+ return( true );
+ }
+}
+
+static int MinMaxWithin(Spline *spline) {
+ extended dx, dy;
+ int which;
+ extended t1, t2;
+ extended w;
+ /* We know that this "spline" is basically one dimensional. As long as its*/
+ /* extrema are between the start and end points on that line then we can */
+ /* treat it as a line. If the extrema are way outside the line segment */
+ /* then it's a line that backtracks on itself */
+
+ if ( (dx = spline->to->me.x - spline->from->me.x)<0 ) dx = -dx;
+ if ( (dy = spline->to->me.y - spline->from->me.y)<0 ) dy = -dy;
+ which = dx<dy;
+ SplineFindExtrema(&spline->splines[which],&t1,&t2);
+ if ( t1==-1 )
+return( true );
+ w = ((spline->splines[which].a*t1 + spline->splines[which].b)*t1
+ + spline->splines[which].c)*t1 + spline->splines[which].d;
+ if ( RealNear(w, (&spline->to->me.x)[which]) || RealNear(w, (&spline->from->me.x)[which]) )
+ /* Close enough */;
+ else if ( (w<(&spline->to->me.x)[which] && w<(&spline->from->me.x)[which]) ||
+ (w>(&spline->to->me.x)[which] && w>(&spline->from->me.x)[which]) )
+return( false ); /* Outside */
+
+ w = ((spline->splines[which].a*t2 + spline->splines[which].b)*t2
+ + spline->splines[which].c)*t2 + spline->splines[which].d;
+ if ( RealNear(w, (&spline->to->me.x)[which]) || RealNear(w, (&spline->from->me.x)[which]) )
+ /* Close enough */;
+ else if ( (w<(&spline->to->me.x)[which] && w<(&spline->from->me.x)[which]) ||
+ (w>(&spline->to->me.x)[which] && w>(&spline->from->me.x)[which]) )
+return( false ); /* Outside */
+
+return( true );
+}
+
+int SplineIsLinear(Spline *spline) {
+ bigreal t1,t2, t3,t4;
+ int ret;
+
+ if ( spline->knownlinear )
+return( true );
+ if ( spline->knowncurved )
+return( false );
+
+ if ( spline->splines[0].a==0 && spline->splines[0].b==0 &&
+ spline->splines[1].a==0 && spline->splines[1].b==0 )
+return( true );
+
+ /* Something is linear if the control points lie on the line between the */
+ /* two base points */
+
+ /* Vertical lines */
+ if ( RealNear(spline->from->me.x,spline->to->me.x) ) {
+ ret = RealNear(spline->from->me.x,spline->from->nextcp.x) &&
+ RealNear(spline->from->me.x,spline->to->prevcp.x);
+ if ( ret && ! ((spline->from->nextcp.y >= spline->from->me.y &&
+ spline->from->nextcp.y <= spline->to->me.y &&
+ spline->to->prevcp.y >= spline->from->me.y &&
+ spline->to->prevcp.y <= spline->to->me.y ) ||
+ (spline->from->nextcp.y <= spline->from->me.y &&
+ spline->from->nextcp.y >= spline->to->me.y &&
+ spline->to->prevcp.y <= spline->from->me.y &&
+ spline->to->prevcp.y >= spline->to->me.y )) )
+ ret = MinMaxWithin(spline);
+ /* Horizontal lines */
+ } else if ( RealNear(spline->from->me.y,spline->to->me.y) ) {
+ ret = RealNear(spline->from->me.y,spline->from->nextcp.y) &&
+ RealNear(spline->from->me.y,spline->to->prevcp.y);
+ if ( ret && ! ((spline->from->nextcp.x >= spline->from->me.x &&
+ spline->from->nextcp.x <= spline->to->me.x &&
+ spline->to->prevcp.x >= spline->from->me.x &&
+ spline->to->prevcp.x <= spline->to->me.x) ||
+ (spline->from->nextcp.x <= spline->from->me.x &&
+ spline->from->nextcp.x >= spline->to->me.x &&
+ spline->to->prevcp.x <= spline->from->me.x &&
+ spline->to->prevcp.x >= spline->to->me.x)) )
+ ret = MinMaxWithin(spline);
+ } else {
+ ret = true;
+ t1 = (spline->from->nextcp.y-spline->from->me.y)/(spline->to->me.y-spline->from->me.y);
+ t2 = (spline->from->nextcp.x-spline->from->me.x)/(spline->to->me.x-spline->from->me.x);
+ t3 = (spline->to->me.y-spline->to->prevcp.y)/(spline->to->me.y-spline->from->me.y);
+ t4 = (spline->to->me.x-spline->to->prevcp.x)/(spline->to->me.x-spline->from->me.x);
+ ret = (Within16RoundingErrors(t1,t2) || (RealApprox(t1,0) && RealApprox(t2,0))) &&
+ (Within16RoundingErrors(t3,t4) || (RealApprox(t3,0) && RealApprox(t4,0)));
+ if ( ret ) {
+ if ( t1<0 || t2<0 || t3<0 || t4<0 ||
+ t1>1 || t2>1 || t3>1 || t4>1 )
+ ret = MinMaxWithin(spline);
+ }
+ }
+ spline->knowncurved = !ret;
+ spline->knownlinear = ret;
+ if ( ret ) {
+ /* A few places that if the spline is knownlinear then its splines[?] */
+ /* are linear. So give the linear version and not that suggested by */
+ /* the control points */
+ spline->splines[0].a = spline->splines[0].b = 0;
+ spline->splines[0].d = spline->from->me.x;
+ spline->splines[0].c = spline->to->me.x-spline->from->me.x;
+ spline->splines[1].a = spline->splines[1].b = 0;
+ spline->splines[1].d = spline->from->me.y;
+ spline->splines[1].c = spline->to->me.y-spline->from->me.y;
+ }
+return( ret );
+}
+
+static bigreal FindZero5(bigreal w[7],bigreal tlow, bigreal thigh) {
+ /* Somewhere between tlow and thigh there is a value of t where w(t)==0 */
+ /* It is conceiveable that there might be 3 such ts if there are some high frequency effects */
+ /* but I ignore that for now */
+ bigreal t, test;
+ int bot_negative;
+
+ t = tlow;
+ test = ((((w[5]*t+w[4])*t+w[3])*t+w[2])*t+w[1])*t + w[0];
+ bot_negative = test<0;
+
+ for (;;) {
+ t = (thigh+tlow)/2;
+ if ( thigh==t || tlow==t )
+return( t ); /* close as we can get */
+ test = ((((w[5]*t+w[4])*t+w[3])*t+w[2])*t+w[1])*t + w[0];
+ if ( test==0 )
+return( t );
+ if ( bot_negative ) {
+ if ( test<0 )
+ tlow = t;
+ else
+ thigh = t;
+ } else {
+ if ( test<0 )
+ thigh = t;
+ else
+ tlow = t;
+ }
+ }
+}
+
+static bigreal FindZero3(bigreal w[7],bigreal tlow, bigreal thigh) {
+ /* Somewhere between tlow and thigh there is a value of t where w(t)==0 */
+ /* It is conceiveable that there might be 3 such ts if there are some high frequency effects */
+ /* but I ignore that for now */
+ bigreal t, test;
+ int bot_negative;
+
+ t = tlow;
+ test = ((w[3]*t+w[2])*t+w[1])*t + w[0];
+ bot_negative = test<0;
+
+ for (;;) {
+ t = (thigh+tlow)/2;
+ if ( thigh==t || tlow==t )
+return( t ); /* close as we can get */
+ test = ((w[3]*t+w[2])*t+w[1])*t + w[0];
+ if ( test==0 )
+return( t );
+ if ( bot_negative ) {
+ if ( test<0 )
+ tlow = t;
+ else
+ thigh = t;
+ } else {
+ if ( test<0 )
+ thigh = t;
+ else
+ tlow = t;
+ }
+ }
+}
+
+bigreal SplineMinDistanceToPoint(Spline *s, BasePoint *p) {
+ /* So to find the minimum distance we want the sqrt( (sx(t)-px)^2 + (sy(t)-py)^2 ) */
+ /* Same minima as (sx(t)-px)^2 + (sy(t)-py)^2, which is easier to deal with */
+ bigreal w[7];
+ Spline1D *x = &s->splines[0], *y = &s->splines[1];
+ bigreal off[2], best;
+
+ off[0] = (x->d-p->x); off[1] = (y->d-p->y);
+
+ w[6] = (x->a*x->a) + (y->a*y->a);
+ w[5] = 2*(x->a*x->b + y->a*y->b);
+ w[4] = (x->b*x->b) + 2*(x->a*x->c) + (y->b*y->b) + 2*(y->a*y->c);
+ w[3] = 2* (x->b*x->c + x->a*off[0] + y->b*y->c + y->a*off[1]);
+ w[2] = (x->c*x->c) + 2*(x->b*off[0]) + (y->c*y->c) + 2*y->b*off[1];
+ w[1] = 2*(x->c*off[0] + y->c*off[1]);
+ w[0] = off[0]*off[0] + off[1]*off[1];
+
+ /* Take derivative */
+ w[0] = w[1];
+ w[1] = 2*w[2];
+ w[2] = 3*w[3];
+ w[3] = 4*w[4];
+ w[4] = 5*w[5];
+ w[5] = 6*w[6];
+ w[6] = 0;
+
+ if ( w[5]!=0 ) {
+ bigreal tzeros[8], t, incr, test, lasttest, zerot;
+ int i, zcnt=0;
+ /* Well, we've got a 5th degree poly and no way to play cute tricks. */
+ /* brute force it */
+ incr = 1.0/1024;
+ lasttest = w[0];
+ for ( t = incr; t<=1.0; t += incr ) {
+ test = ((((w[5]*t+w[4])*t+w[3])*t+w[2])*t+w[1])*t + w[0];
+ if ( test==0 )
+ tzeros[zcnt++] = t;
+ else {
+ if ( lasttest!=0 && (test>0) != (lasttest>0) ) {
+ zerot = FindZero5(w,t-incr,t);
+ if ( zerot>0 )
+ tzeros[zcnt++] = zerot;
+ }
+ }
+ lasttest = test;
+ }
+ best = off[0]*off[0] + off[1]*off[1]; /* t==0 */
+ test = (x->a+x->b+x->c+off[0])*(x->a+x->b+x->c+off[0]) +
+ (y->a+y->b+y->c+off[1])*(y->a+y->b+y->c+off[1]); /* t==1 */
+ if ( best>test ) best = test;
+ for ( i=0; i<zcnt; ++i ) {
+ bigreal tx, ty;
+ tx = ((x->a*tzeros[i]+x->b)*tzeros[i]+x->c)*tzeros[i] + off[0];
+ ty = ((y->a*tzeros[i]+y->b)*tzeros[i]+y->c)*tzeros[i] + off[1];
+ test = tx*tx + ty*ty;
+ if ( best>test ) best = test;
+ }
+return( sqrt(best));
+ } else if ( w[4]==0 && w[3]!=0 ) {
+ /* Started with a quadratic -- now, find 0s of a cubic */
+ /* We could find the extrema, so we have a bunch of monotonics */
+ /* Or we could brute force it as above */
+ bigreal tzeros[8], test, zerot;
+ bigreal quad[3], disc, e[5], t1, t2;
+ int i, zcnt=0, ecnt;
+
+ quad[2] = 3*w[3]; quad[1] = 2*w[2]; quad[0] = w[1];
+ disc = (-quad[1]*quad[1] - 4*quad[2]*quad[0]);
+ e[0] = 0;
+ if ( disc<0 ) {
+ e[1] = 1.0;
+ ecnt = 2;
+ } else
+ disc = sqrt(disc);
+ t1 = (-w[1] - disc) / (2*w[2]);
+ t2 = (-w[1] + disc) / (2*w[2]);
+ if ( t1>t2 ) {
+ bigreal temp = t1;
+ t1 = t2;
+ t2 = temp;
+ }
+ ecnt=1;
+ if ( t1>0 && t1<1 )
+ e[ecnt++] = t1;
+ if ( t2>0 && t2<1 && t1!=t2 )
+ e[ecnt++] = t2;
+ e[ecnt++] = 1.0;
+ for ( i=1; i<ecnt; ++i ) {
+ zerot = FindZero3(w,e[i-1],e[i]);
+ if ( zerot>0 )
+ tzeros[zcnt++] = zerot;
+ }
+ best = off[0]*off[0] + off[1]*off[1]; /* t==0 */
+ test = (x->b+x->c+off[0])*(x->b+x->c+off[0]) +
+ (y->b+y->c+off[1])*(y->b+y->c+off[1]); /* t==1 */
+ if ( best>test ) best = test;
+ for ( i=0; i<zcnt; ++i ) {
+ bigreal tx, ty;
+ tx = (x->b*tzeros[i]+x->c)*tzeros[i] + off[0];
+ ty = (y->b*tzeros[i]+y->c)*tzeros[i] + off[1];
+ test = tx*tx + ty*ty;
+ if ( best>test ) best = test;
+ }
+return( sqrt(best));
+ } else if ( w[2]==0 && w[1]!=0 ) {
+ /* Started with a line */
+ bigreal t = -w[0]/w[1], test, best;
+ best = off[0]*off[0] + off[1]*off[1]; /* t==0 */
+ test = (x->c+off[0])*(x->c+off[0]) + (y->c+off[1])*(y->c+off[1]); /* t==1 */
+ if ( best>test ) best = test;
+ if ( t>0 && t<1 ) {
+ test = (x->c*t+off[0])*(x->c*t+off[0]) + (y->c*t+off[1])*(y->c*t+off[1]);
+ if ( best>test ) best = test;
+ }
+return(sqrt(best));
+ } else if ( w[4]!=0 && w[3]!=0 && w[2]!=0 && w[1]!=0 ) {
+ IError( "Impossible condition in SplineMinDistanceToPoint");
+ } else {
+ /* It's a point, minimum distance is the only distance */
+return( sqrt(off[0]*off[0] + off[1]*off[1]) );
+ }
+return( -1 );
+}
+
+/* This returns all real solutions, even those out of bounds */
+/* I use -999999 as an error flag, since we're really only interested in */
+/* solns near 0 and 1 that should be ok. -1 is perhaps a little too close */
+/* Sigh. When solutions are near 0, the rounding errors are appalling. */
+int _CubicSolve(const Spline1D *sp,bigreal sought, extended ts[3]) {
+ extended d, xN, yN, delta2, temp, delta, h, t2, t3, theta;
+ extended sa=sp->a, sb=sp->b, sc=sp->c, sd=sp->d-sought;
+ int i=0;
+
+ ts[0] = ts[1] = ts[2] = -999999;
+ if ( sd==0 && sa!=0 ) {
+ /* one of the roots is 0, the other two are the soln of a quadratic */
+ ts[0] = 0;
+ if ( sc==0 ) {
+ ts[1] = -sb/(extended) sa; /* two zero roots */
+ } else {
+ temp = sb*(extended) sb-4*(extended) sa*sc;
+ if ( RealNear(temp,0))
+ ts[1] = -sb/(2*(extended) sa);
+ else if ( temp>=0 ) {
+ temp = sqrt(temp);
+ ts[1] = (-sb+temp)/(2*(extended) sa);
+ ts[2] = (-sb-temp)/(2*(extended) sa);
+ }
+ }
+ } else if ( sa!=0 ) {
+ /* http://www.m-a.org.uk/eb/mg/mg077ch.pdf */
+ /* this nifty solution to the cubic neatly avoids complex arithmetic */
+ xN = -sb/(3*(extended) sa);
+ yN = ((sa*xN + sb)*xN+sc)*xN + sd;
+
+ delta2 = (sb*(extended) sb-3*(extended) sa*sc)/(9*(extended) sa*sa);
+ /*if ( RealWithin(delta2,0,.00000001) ) delta2 = 0;*/
+
+ /* the descriminant is yN^2-h^2, but delta might be <0 so avoid using h */
+ d = yN*yN - 4*sa*sa*delta2*delta2*delta2;
+ if ( ((yN>.01 || yN<-.01) && RealNear(d/yN,0)) || ((yN<=.01 && yN>=-.01) && RealNear(d,0)) )
+ d = 0;
+ if ( d>0 ) {
+ temp = sqrt(d);
+ t2 = (-yN-temp)/(2*sa);
+ t2 = (t2==0) ? 0 : (t2<0) ? -pow(-t2,1./3.) : pow(t2,1./3.);
+ t3 = (-yN+temp)/(2*sa);
+ t3 = t3==0 ? 0 : (t3<0) ? -pow(-t3,1./3.) : pow(t3,1./3.);
+ ts[0] = xN + t2 + t3;
+ } else if ( d<0 ) {
+ if ( delta2>=0 ) {
+ delta = sqrt(delta2);
+ h = 2*sa*delta2*delta;
+ temp = -yN/h;
+ if ( temp>=-1.0001 && temp<=1.0001 ) {
+ if ( temp<-1 ) temp = -1; else if ( temp>1 ) temp = 1;
+ theta = acos(temp)/3;
+ ts[i++] = xN+2*delta*cos(theta);
+ ts[i++] = xN+2*delta*cos(2.0943951+theta); /* 2*pi/3 */
+ ts[i++] = xN+2*delta*cos(4.1887902+theta); /* 4*pi/3 */
+ }
+ }
+ } else if ( /* d==0 && */ delta2!=0 ) {
+ delta = yN/(2*sa);
+ delta = delta==0 ? 0 : delta>0 ? pow(delta,1./3.) : -pow(-delta,1./3.);
+ ts[i++] = xN + delta; /* this root twice, but that's irrelevant to me */
+ ts[i++] = xN - 2*delta;
+ } else if ( /* d==0 && */ delta2==0 ) {
+ if ( xN>=-0.0001 && xN<=1.0001 ) ts[0] = xN;
+ }
+ } else if ( sb!=0 ) {
+ extended d = sc*(extended) sc-4*(extended) sb*sd;
+ if ( d<0 && RealNear(d,0)) d=0;
+ if ( d<0 )
+return(false); /* All roots imaginary */
+ d = sqrt(d);
+ ts[0] = (-sc-d)/(2*(extended) sb);
+ ts[1] = (-sc+d)/(2*(extended) sb);
+ } else if ( sc!=0 ) {
+ ts[0] = -sd/(extended) sc;
+ } else {
+ /* If it's a point then either everything is a solution, or nothing */
+ }
+return( ts[0]!=-999999 );
+}
+
+int _QuarticSolve(Quartic *q,extended ts[4]) {
+ extended extrema[5];
+ Spline1D sp;
+ int ecnt = 0, i, zcnt;
+
+ /* Two special cases */
+ if ( q->a==0 ) { /* It's really a cubic */
+ sp.a = q->b;
+ sp.b = q->c;
+ sp.c = q->d;
+ sp.d = q->e;
+ ts[3] = -999999;
+return( _CubicSolve(&sp,0,ts));
+ } else if ( q->e==0 ) { /* we can factor out a zero root */
+ sp.a = q->a;
+ sp.b = q->b;
+ sp.c = q->c;
+ sp.d = q->d;
+ ts[0] = 0;
+return( _CubicSolve(&sp,0,ts+1)+1);
+ }
+
+ sp.a = 4*q->a;
+ sp.b = 3*q->b;
+ sp.c = 2*q->c;
+ sp.d = q->d;
+ if ( _CubicSolve(&sp,0,extrema)) {
+ ecnt = 1;
+ if ( extrema[1]!=-999999 ) {
+ ecnt = 2;
+ if ( extrema[1]<extrema[0] ) {
+ extended temp = extrema[1]; extrema[1] = extrema[0]; extrema[0]=temp;
+ }
+ if ( extrema[2]!=-999999 ) {
+ ecnt = 3;
+ if ( extrema[2]<extrema[0] ) {
+ extended temp = extrema[2]; extrema[2] = extrema[0]; extrema[0]=temp;
+ }
+ if ( extrema[2]<extrema[1] ) {
+ extended temp = extrema[2]; extrema[2] = extrema[1]; extrema[1]=temp;
+ }
+ }
+ }
+ }
+ for ( i=ecnt-1; i>=0 ; --i )
+ extrema[i+1] = extrema[i];
+ /* Upper and lower bounds within which we'll search */
+ extrema[0] = -999;
+ extrema[ecnt+1] = 999;
+ ecnt += 2;
+ /* divide into monotonic sections & use binary search to find zeroes */
+ for ( i=zcnt=0; i<ecnt-1; ++i ) {
+ extended top, bottom, val;
+ extended topt, bottomt, t;
+ topt = extrema[i+1];
+ bottomt = extrema[i];
+ top = (((q->a*topt+q->b)*topt+q->c)*topt+q->d)*topt+q->e;
+ bottom = (((q->a*bottomt+q->b)*bottomt+q->c)*bottomt+q->d)*bottomt+q->e;
+ if ( top<bottom ) {
+ extended temp = top; top = bottom; bottom = temp;
+ temp = topt; topt = bottomt; bottomt = temp;
+ }
+ if ( bottom>.001 ) /* this monotonic is all above 0 */
+ continue;
+ if ( top<-.001 ) /* this monotonic is all below 0 */
+ continue;
+ if ( bottom>0 ) {
+ ts[zcnt++] = bottomt;
+ continue;
+ }
+ if ( top<0 ) {
+ ts[zcnt++] = topt;
+ continue;
+ }
+ for (;;) {
+ t = (topt+bottomt)/2;
+ if ( isnan(t) ) {
+ break;
+ } else if ( t==topt || t==bottomt ) {
+ ts[zcnt++] = t;
+ break;
+ }
+
+ val = (((q->a*t+q->b)*t+q->c)*t+q->d)*t+q->e;
+ if ( val>-.0001 && val<.0001 ) {
+ ts[zcnt++] = t;
+ break;
+ } else if ( val>0 ) {
+ top = val;
+ topt = t;
+ } else {
+ bottom = val;
+ bottomt = t;
+ }
+ }
+ }
+ for ( i=zcnt; i<4; ++i )
+ ts[i] = -999999;
+return( zcnt );
+}
+
+/* calculating the actual length of a spline is hard, this gives a very */
+/* rough (but quick) approximation */
+static bigreal SplineLenApprox(Spline *spline) {
+ bigreal len, slen, temp;
+
+ if ( (temp = spline->to->me.x-spline->from->me.x)<0 ) temp = -temp;
+ len = temp;
+ if ( (temp = spline->to->me.y-spline->from->me.y)<0 ) temp = -temp;
+ len += temp;
+ if ( !spline->to->noprevcp || !spline->from->nonextcp ) {
+ if ( (temp = spline->from->nextcp.x-spline->from->me.x)<0 ) temp = -temp;
+ slen = temp;
+ if ( (temp = spline->from->nextcp.y-spline->from->me.y)<0 ) temp = -temp;
+ slen += temp;
+ if ( (temp = spline->to->prevcp.x-spline->from->nextcp.x)<0 ) temp = -temp;
+ slen += temp;
+ if ( (temp = spline->to->prevcp.y-spline->from->nextcp.y)<0 ) temp = -temp;
+ slen += temp;
+ if ( (temp = spline->to->me.x-spline->to->prevcp.x)<0 ) temp = -temp;
+ slen += temp;
+ if ( (temp = spline->to->me.y-spline->to->prevcp.y)<0 ) temp = -temp;
+ slen += temp;
+ len = (len + slen)/2;
+ }
+return( len );
+}
+
+FitPoint *SplinesFigureFPsBetween(SplinePoint *from, SplinePoint *to,
+ int *tot) {
+ int cnt, i, j, pcnt;
+ bigreal len, slen, lbase;
+ SplinePoint *np;
+ FitPoint *fp;
+ bigreal _lens[10], *lens = _lens;
+ int _cnts[10], *cnts = _cnts;
+ /* I used just to give every spline 10 points. But that gave much more */
+ /* weight to small splines than to big ones */
+
+ cnt = 0;
+ for ( np = from->next->to; ; np = np->next->to ) {
+ ++cnt;
+ if ( np==to )
+ break;
+ }
+ if ( cnt>10 ) {
+ lens = malloc(cnt*sizeof(bigreal));
+ cnts = malloc(cnt*sizeof(int));
+ }
+ cnt = 0; len = 0;
+ for ( np = from->next->to; ; np = np->next->to ) {
+ lens[cnt] = SplineLenApprox(np->prev);
+ len += lens[cnt];
+ ++cnt;
+ if ( np==to )
+ break;
+ }
+ if ( len!=0 ) {
+ pcnt = 0;
+ for ( i=0; i<cnt; ++i ) {
+ int pnts = rint( (10*cnt*lens[i])/len );
+ if ( pnts<2 ) pnts = 2;
+ cnts[i] = pnts;
+ pcnt += pnts;
+ }
+ } else
+ pcnt = 2*cnt;
+
+ fp = malloc((pcnt+1)*sizeof(FitPoint)); i = 0;
+ if ( len==0 ) {
+ for ( i=0; i<=pcnt; ++i ) {
+ fp[i].t = i/(pcnt);
+ fp[i].p.x = from->me.x;
+ fp[i].p.y = from->me.y;
+ }
+ } else {
+ lbase = 0;
+ for ( i=cnt=0, np = from->next->to; ; np = np->next->to, ++cnt ) {
+ slen = SplineLenApprox(np->prev);
+ for ( j=0; j<cnts[cnt]; ++j ) {
+ bigreal t = j/(bigreal) cnts[cnt];
+ fp[i].t = (lbase+ t*slen)/len;
+ fp[i].p.x = ((np->prev->splines[0].a*t+np->prev->splines[0].b)*t+np->prev->splines[0].c)*t + np->prev->splines[0].d;
+ fp[i++].p.y = ((np->prev->splines[1].a*t+np->prev->splines[1].b)*t+np->prev->splines[1].c)*t + np->prev->splines[1].d;
+ }
+ lbase += slen;
+ if ( np==to )
+ break;
+ }
+ }
+ if ( cnts!=_cnts ) free(cnts);
+ if ( lens!=_lens ) free(lens);
+
+ *tot = i;
+
+return( fp );
+}
+
+static int SplinePointCategory(SplinePoint *sp) {
+ enum pointtype pt;
+
+ pt = pt_corner;
+ if ( sp->next==NULL && sp->prev==NULL )
+ ;
+ else if ( (sp->next!=NULL && sp->next->to->me.x==sp->me.x && sp->next->to->me.y==sp->me.y) ||
+ (sp->prev!=NULL && sp->prev->from->me.x==sp->me.x && sp->prev->from->me.y==sp->me.y ))
+ ;
+ else if ( sp->next==NULL ) {
+ pt = sp->noprevcp ? pt_corner : pt_curve;
+ } else if ( sp->prev==NULL ) {
+ pt = sp->nonextcp ? pt_corner : pt_curve;
+ } else if ( sp->nonextcp && sp->noprevcp ) {
+ ;
+ } else {
+ BasePoint ndir, ncdir, ncunit, pdir, pcdir, pcunit;
+ bigreal nlen, nclen, plen, pclen;
+ bigreal cross, bounds;
+
+ ncdir.x = sp->nextcp.x - sp->me.x; ncdir.y = sp->nextcp.y - sp->me.y;
+ pcdir.x = sp->prevcp.x - sp->me.x; pcdir.y = sp->prevcp.y - sp->me.y;
+ ndir.x = ndir.y = pdir.x = pdir.y = 0;
+ if ( sp->next!=NULL ) {
+ ndir.x = sp->next->to->me.x - sp->me.x; ndir.y = sp->next->to->me.y - sp->me.y;
+ }
+ if ( sp->prev!=NULL ) {
+ pdir.x = sp->prev->from->me.x - sp->me.x; pdir.y = sp->prev->from->me.y - sp->me.y;
+ }
+ nclen = sqrt(ncdir.x*ncdir.x + ncdir.y*ncdir.y);
+ pclen = sqrt(pcdir.x*pcdir.x + pcdir.y*pcdir.y);
+ nlen = sqrt(ndir.x*ndir.x + ndir.y*ndir.y);
+ plen = sqrt(pdir.x*pdir.x + pdir.y*pdir.y);
+ ncunit = ncdir; pcunit = pcdir;
+ if ( nclen!=0 ) { ncunit.x /= nclen; ncunit.y /= nclen; }
+ if ( pclen!=0 ) { pcunit.x /= pclen; pcunit.y /= pclen; }
+ if ( nlen!=0 ) { ndir.x /= nlen; ndir.y /= nlen; }
+ if ( plen!=0 ) { pdir.x /= plen; pdir.y /= plen; }
+
+ /* find out which side has the shorter control vector. Cross that vector */
+ /* with the normal of the unit vector on the other side. If the */
+ /* result is less than 1 em-unit then we've got colinear control points */
+ /* (within the resolution of the integer grid) */
+ /* Not quite... they could point in the same direction */
+ if ( sp->pointtype==pt_curve )
+ bounds = 4.0;
+ else
+ bounds = 1.0;
+ if ( nclen!=0 && pclen!=0 &&
+ ((nclen>=pclen && (cross = pcdir.x*ncunit.y - pcdir.y*ncunit.x)<bounds && cross>-bounds ) ||
+ (pclen>nclen && (cross = ncdir.x*pcunit.y - ncdir.y*pcunit.x)<bounds && cross>-bounds )) &&
+ ncdir.x*pcdir.x + ncdir.y*pcdir.y < 0 )
+ pt = pt_curve;
+ /* Cross product of control point with unit vector normal to line in */
+ /* opposite direction should be less than an em-unit for a tangent */
+ else if ( ( nclen==0 && pclen!=0
+ && (cross = pcdir.x*ndir.y-pcdir.y*ndir.x)<bounds
+ && cross>-bounds && (pcdir.x*ndir.x+pcdir.y*ndir.y)<0 )
+ ||
+ ( pclen==0 && nclen!=0
+ && (cross = ncdir.x*pdir.y-ncdir.y*pdir.x)<bounds
+ && cross>-bounds && (ncdir.x*pdir.x+ncdir.y*pdir.y)<0 ) )
+ pt = pt_tangent;
+
+ if (pt == pt_curve &&
+ ((sp->nextcp.x==sp->me.x && sp->prevcp.x==sp->me.x && sp->nextcp.y!=sp->me.y) ||
+ (sp->nextcp.y==sp->me.y && sp->prevcp.y==sp->me.y && sp->nextcp.x!=sp->me.x)))
+ pt = pt_hvcurve;
+ }
+ return pt;
+}
+
+static enum pointtype SplinePointDowngrade(int current, int geom) {
+ enum pointtype np = current;
+
+ if ( current==pt_curve && geom!=pt_curve ) {
+ if ( geom==pt_hvcurve )
+ np = pt_curve;
+ else
+ np = pt_corner;
+ } else if ( current==pt_hvcurve && geom!=pt_hvcurve ) {
+ if ( geom==pt_curve )
+ np = pt_curve;
+ else
+ np = pt_corner;
+ } else if ( current==pt_tangent && geom!=pt_tangent ) {
+ np = pt_corner;
+ }
+
+ return np;
+}
+
+// Assumes flag combinations are already verified. Only returns false
+// when called with check_compat
+int _SplinePointCategorize(SplinePoint *sp, int flags) {
+ enum pointtype geom, dg, cur;
+
+ if ( flags & pconvert_flag_none )
+ // No points selected for conversion -- keep type as is
+ return true;
+ if ( flags & pconvert_flag_smooth && sp->pointtype == pt_corner )
+ // Convert only "smooth" points, not corners
+ return true;
+
+ geom = SplinePointCategory(sp);
+ dg = SplinePointDowngrade(sp->pointtype, geom);
+
+ if ( flags & pconvert_flag_incompat && sp->pointtype == dg )
+ // Only convert points incompatible with current type
+ return true;
+
+ if ( flags & pconvert_flag_by_geom ) {
+ if ( ! ( flags & pconvert_flag_hvcurve ) && geom == pt_hvcurve )
+ sp->pointtype = pt_curve;
+ else
+ sp->pointtype = geom;
+ } else if ( flags & pconvert_flag_downgrade ) {
+ sp->pointtype = dg;
+ } else if ( flags & pconvert_flag_force_type ) {
+ if ( sp->pointtype != dg ) {
+ cur = sp->pointtype;
+ sp->pointtype = dg;
+ /* SPChangePointType(sp,cur); */
+ }
+ } else if ( flags & pconvert_flag_check_compat ) {
+ if ( sp->pointtype != dg )
+ return false;
+ }
+ return true;
+}
+
+void SplinePointCategorize(SplinePoint *sp) {
+ _SplinePointCategorize(sp, pconvert_flag_all|pconvert_flag_by_geom);
+}
+
+static void SplinePointReCategorize(SplinePoint *sp,int oldpt) {
+ SplinePointCategorize(sp);
+ if ( sp->pointtype!=oldpt ) {
+ if ( sp->pointtype==pt_curve && oldpt==pt_hvcurve &&
+ ((sp->nextcp.x == sp->me.x && sp->nextcp.y != sp->me.y ) ||
+ (sp->nextcp.y == sp->me.y && sp->nextcp.x != sp->me.x )))
+ sp->pointtype = pt_hvcurve;
+ }
+}
+
+void SplinesRemoveBetween(SplinePoint *from, SplinePoint *to, int type) {
+ int tot;
+ FitPoint *fp;
+ SplinePoint *np, oldfrom;
+ int oldfpt = from->pointtype, oldtpt = to->pointtype;
+ Spline *sp;
+ int order2 = from->next->order2;
+
+ oldfrom = *from;
+ fp = SplinesFigureFPsBetween(from,to,&tot);
+
+ if ( type==1 )
+ ApproximateSplineFromPointsSlopes(from,to,fp,tot-1,order2,mt_levien);
+ else
+ ApproximateSplineFromPoints(from,to,fp,tot-1,order2);
+
+ /* Have to do the frees after the approximation because the approx */
+ /* uses the splines to determine slopes */
+ for ( sp = oldfrom.next; ; ) {
+ np = sp->to;
+ SplineFree(sp);
+ if ( np==to )
+ break;
+ sp = np->next;
+ // SplinePointMDFree(sc,np);
+ }
+
+ free(fp);
+
+ SplinePointReCategorize(from,oldfpt);
+ SplinePointReCategorize(to,oldtpt);
+}
diff --git a/src/path/splinefit/splinefont.h b/src/path/splinefit/splinefont.h
new file mode 100644
index 0000000..83db934
--- /dev/null
+++ b/src/path/splinefit/splinefont.h
@@ -0,0 +1,191 @@
+// SPDX-License-Identifier: GPL-2.0-or-later
+
+#ifndef _SEEN_SPLINEFONT_H_
+#define _SEEN_SPLINEFONT_H_
+
+#include <glib.h>
+
+typedef double real;
+typedef double bigreal;
+typedef double extended;
+typedef int BOOL;
+
+#define chunkalloc(size) calloc(1,size)
+#define chunkfree(item,size) free(item)
+
+typedef struct basepoint {
+ real x;
+ real y;
+} BasePoint;
+
+typedef struct ipoint {
+ int x;
+ int y;
+} IPoint;
+
+enum pointtype { pt_curve, pt_corner, pt_tangent, pt_hvcurve };
+typedef struct splinepoint {
+ BasePoint me;
+ BasePoint nextcp; /* control point */
+ BasePoint prevcp; /* control point */
+ unsigned int nonextcp:1;
+ unsigned int noprevcp:1;
+ unsigned int nextcpdef:1;
+ unsigned int prevcpdef:1;
+ unsigned int selected:1; /* for UI */
+ unsigned int nextcpselected: 2; /* Is the next BCP selected */
+ unsigned int prevcpselected: 2; /* Is the prev BCP selected */
+ unsigned int pointtype:2;
+ unsigned int isintersection: 1;
+ unsigned int flexy: 1; /* When "freetype_markup" is on in charview.c:DrawPoint */
+ unsigned int flexx: 1; /* flexy means select nextcp, and flexx means draw circle around nextcp */
+ unsigned int roundx: 1; /* For true type hinting */
+ unsigned int roundy: 1; /* For true type hinting */
+ unsigned int dontinterpolate: 1; /* in ttf, don't imply point by interpolating between cps */
+ unsigned int ticked: 1;
+ unsigned int watched: 1;
+ /* 1 bits left... */
+ uint16_t ptindex; /* Temporary value used by metafont routine */
+ uint16_t ttfindex; /* Truetype point index */
+ /* Special values 0xffff => point implied by averaging control points */
+ /* 0xfffe => point created with no real number yet */
+ /* (or perhaps point in context where no number is possible as in a glyph with points & refs) */
+ uint16_t nextcpindex; /* Truetype point index */
+ struct spline *next;
+ struct spline *prev;
+ /* Inkscape: not used; HintMask *hintmask; */
+ char* name;
+} SplinePoint;
+
+
+typedef struct spline1d {
+ real a, b, c, d;
+} Spline1D;
+
+typedef struct spline {
+ unsigned int islinear: 1; /* No control points */
+ unsigned int isquadratic: 1; /* probably read in from ttf */
+ unsigned int isticked: 1;
+ unsigned int isneeded: 1; /* Used in remove overlap */
+ unsigned int isunneeded: 1; /* Used in remove overlap */
+ unsigned int exclude: 1; /* Used in remove overlap variant: exclude */
+ unsigned int ishorvert: 1;
+ unsigned int knowncurved: 1; /* We know that it curves */
+ unsigned int knownlinear: 1; /* it might have control points, but still traces out a line */
+ /* If neither knownlinear nor curved then we haven't checked */
+ unsigned int order2: 1; /* It's a bezier curve with only one cp */
+ unsigned int touched: 1;
+ unsigned int leftedge: 1;
+ unsigned int rightedge: 1;
+ unsigned int acceptableextrema: 1; /* This spline has extrema, but we don't care */
+ SplinePoint *from;
+ SplinePoint *to;
+ Spline1D splines[2]; /* splines[0] is the x spline, splines[1] is y */
+ struct linearapprox *approx;
+ /* Possible optimizations:
+ Precalculate bounding box
+ Precalculate min/max/ points of inflection
+ */
+} Spline;
+
+typedef struct splinepointlist {
+ SplinePoint *first, *last;
+ struct splinepointlist *next;
+ /* Not used: spiro_cp *spiros; */
+ uint16_t spiro_cnt, spiro_max;
+ /* These could be bit fields, but bytes are easier to access and we */
+ /* don't need the space (yet) */
+ uint8_t ticked;
+ uint8_t beziers_need_optimizer; /* If the spiros have changed in spiro mode, then reverting to bezier mode might, someday, run a simplifier */
+ uint8_t is_clip_path; /* In type3/svg fonts */
+ int start_offset; // This indicates which point is the canonical first for purposes of outputting to U. F. O..
+ char *contour_name;
+} SplinePointList, SplineSet;
+
+typedef struct dbounds {
+ real minx, maxx;
+ real miny, maxy;
+} DBounds;
+
+typedef struct quartic {
+ bigreal a,b,c,d,e;
+} Quartic;
+
+
+int RealWithin(real a,real b,real fudge);
+BOOL RealNear(real a, real b);
+
+Spline *SplineMake(SplinePoint *from, SplinePoint *to, int order2);
+Spline *SplineMake2(SplinePoint *from, SplinePoint *to);
+Spline *SplineMake3(SplinePoint *from, SplinePoint *to);
+SplinePoint *SplinePointCreate(real x, real y);
+
+void SplineRefigure3(Spline *spline);
+void SplineRefigure(Spline *spline);
+int SplineIsLinear(Spline *spline);
+void SplineFindExtrema(const Spline1D *sp, extended *_t1, extended *_t2 );
+bigreal SplineMinDistanceToPoint(Spline *s, BasePoint *p);
+
+void SplinePointFree(SplinePoint *sp);
+void SplineFree(Spline *spline);
+void SplinePointListFree(SplinePointList *spl);
+
+bigreal BPDot(BasePoint v1, BasePoint v2);
+bigreal BPCross(BasePoint v1, BasePoint v2);
+BasePoint BPRev(BasePoint v);
+
+int _CubicSolve(const Spline1D *sp,bigreal sought, extended ts[3]);
+int _QuarticSolve(Quartic *q,extended ts[4]);
+int IntersectLines(BasePoint *inter,
+ BasePoint *line1_1, BasePoint *line1_2,
+ BasePoint *line2_1, BasePoint *line2_2);
+
+#define IError(msg) g_warning(msg)
+#define TRACE g_message
+
+enum linelist_flags { cvli_onscreen=0x1, cvli_clipped=0x2 };
+
+typedef struct linelist {
+ IPoint here;
+ struct linelist *next;
+ /* The first two fields are constant for the linelist, the next ones */
+ /* refer to a particular screen. If some portion of the line from */
+ /* this point to the next one is on the screen then set cvli_onscreen */
+ /* if this point needs to be clipped then set cvli_clipped */
+ /* asend and asstart are the actual screen locations where this point */
+ /* intersects the clip edge. */
+ enum linelist_flags flags;
+ IPoint asend, asstart;
+} LineList;
+
+typedef struct linearapprox {
+ real scale;
+ unsigned int oneline: 1;
+ unsigned int onepoint: 1;
+ unsigned int any: 1; /* refers to a particular screen */
+ struct linelist *lines;
+ struct linearapprox *next;
+} LinearApprox;
+
+void LinearApproxFree(LinearApprox *la);
+
+int Within16RoundingErrors(bigreal v1, bigreal v2);
+
+enum pconvert_flags {
+ // Point selection (mutually exclusive)
+ pconvert_flag_none = 0x01,
+ pconvert_flag_all = 0x02,
+ pconvert_flag_smooth = 0x04,
+ pconvert_flag_incompat = 0x08,
+ // Conversion modes (mutually exclusive)
+ pconvert_flag_by_geom = 0x100,
+ pconvert_flag_force_type = 0x200,
+ pconvert_flag_downgrade = 0x400,
+ pconvert_flag_check_compat = 0x0800,
+ // Additional
+ pconvert_flag_hvcurve = 0x4000
+};
+
+void SplinesRemoveBetween(SplinePoint *from, SplinePoint *to, int type);
+
+#endif // _SEEN_SPLINEFONT_H_
diff --git a/src/path/splinefit/splinerefigure.c b/src/path/splinefit/splinerefigure.c
new file mode 100644
index 0000000..2f6d09b
--- /dev/null
+++ b/src/path/splinefit/splinerefigure.c
@@ -0,0 +1,117 @@
+// SPDX-License-Identifier: GPL-2.0-or-later
+/* Copyright (C) 2000-2012 by George Williams */
+/*
+ * Redistribution and use in source and binary forms, with or without
+ * modification, are permitted provided that the following conditions are met:
+
+ * Redistributions of source code must retain the above copyright notice, this
+ * list of conditions and the following disclaimer.
+
+ * Redistributions in binary form must reproduce the above copyright notice,
+ * this list of conditions and the following disclaimer in the documentation
+ * and/or other materials provided with the distribution.
+
+ * The name of the author may not be used to endorse or promote products
+ * derived from this software without specific prior written permission.
+
+ * THIS SOFTWARE IS PROVIDED BY THE AUTHOR ``AS IS'' AND ANY EXPRESS OR IMPLIED
+ * WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF
+ * MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO
+ * EVENT SHALL THE AUTHOR BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
+ * SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
+ * PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS;
+ * OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY,
+ * WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR
+ * OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF
+ * ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+ */
+
+#include <stdbool.h>
+#include <stdint.h>
+
+#include "splinerefigure.h"
+
+#include "splinefont.h"
+
+#include <math.h>
+#include <stdio.h>
+#ifdef HAVE_IEEEFP_H
+# include <ieeefp.h> /* Solaris defines isnan in ieeefp rather than math.h */
+#endif
+
+#include <stdbool.h>
+
+/* The slight errors introduced by the optimizer turn out to have nasty */
+/* side effects. An error on the order of 7e-8 in splines[1].b caused */
+/* the rasterizer to have kaniptions */
+void SplineRefigure3(Spline *spline) {
+ SplinePoint *from = spline->from, *to = spline->to;
+ Spline1D *xsp = &spline->splines[0], *ysp = &spline->splines[1];
+ Spline old;
+
+ spline->isquadratic = false;
+ if ( spline->acceptableextrema )
+ old = *spline;
+ xsp->d = from->me.x; ysp->d = from->me.y;
+ // Set noprevcp and nonextcp based on point values but then make sure both
+ // have the same value
+ from->nonextcp = from->nextcp.x==from->me.x && from->nextcp.y == from->me.y;
+ to->noprevcp = to->prevcp.x==to->me.x && to->prevcp.y == to->me.y;
+ if ( !from->nonextcp || !to->noprevcp )
+ from->nonextcp = to->noprevcp = false;
+ if ( from->nonextcp && to->noprevcp ) {
+ spline->islinear = true;
+ xsp->c = to->me.x-from->me.x;
+ ysp->c = to->me.y-from->me.y;
+ xsp->a = xsp->b = 0;
+ ysp->a = ysp->b = 0;
+ } else {
+ /* from p. 393 (Operator Details, curveto) PostScript Lang. Ref. Man. (Red book) */
+ xsp->c = 3*(from->nextcp.x-from->me.x);
+ ysp->c = 3*(from->nextcp.y-from->me.y);
+ xsp->b = 3*(to->prevcp.x-from->nextcp.x)-xsp->c;
+ ysp->b = 3*(to->prevcp.y-from->nextcp.y)-ysp->c;
+ xsp->a = to->me.x-from->me.x-xsp->c-xsp->b;
+ ysp->a = to->me.y-from->me.y-ysp->c-ysp->b;
+ if ( RealNear(xsp->c,0)) xsp->c=0;
+ if ( RealNear(ysp->c,0)) ysp->c=0;
+ if ( RealNear(xsp->b,0)) xsp->b=0;
+ if ( RealNear(ysp->b,0)) ysp->b=0;
+ if ( RealNear(xsp->a,0)) xsp->a=0;
+ if ( RealNear(ysp->a,0)) ysp->a=0;
+ if ( xsp->a!=0 && ( Within16RoundingErrors(xsp->a+from->me.x,from->me.x) ||
+ Within16RoundingErrors(xsp->a+to->me.x,to->me.x)))
+ xsp->a = 0;
+ if ( ysp->a!=0 && ( Within16RoundingErrors(ysp->a+from->me.y,from->me.y) ||
+ Within16RoundingErrors(ysp->a+to->me.y,to->me.y)))
+ ysp->a = 0;
+ SplineIsLinear(spline);
+ spline->islinear = false;
+ if ( ysp->a==0 && xsp->a==0 ) {
+ if ( ysp->b==0 && xsp->b==0 )
+ spline->islinear = true; /* This seems extremely unlikely... */
+ else
+ spline->isquadratic = true; /* Only likely if we read in a TTF */
+ }
+ }
+ if ( !isfinite(ysp->a) || !isfinite(xsp->a) || !isfinite(ysp->c) || !isfinite(xsp->c) || !isfinite(ysp->d) || !isfinite(xsp->d))
+ IError("NaN value in spline creation");
+ LinearApproxFree(spline->approx);
+ spline->approx = NULL;
+ spline->knowncurved = false;
+ spline->knownlinear = spline->islinear;
+ SplineIsLinear(spline);
+ spline->order2 = false;
+
+ if ( spline->acceptableextrema ) {
+ /* I don't check "d", because changes to that reflect simple */
+ /* translations which will not affect the shape of the spline */
+ if ( !RealNear(old.splines[0].a,spline->splines[0].a) ||
+ !RealNear(old.splines[0].b,spline->splines[0].b) ||
+ !RealNear(old.splines[0].c,spline->splines[0].c) ||
+ !RealNear(old.splines[1].a,spline->splines[1].a) ||
+ !RealNear(old.splines[1].b,spline->splines[1].b) ||
+ !RealNear(old.splines[1].c,spline->splines[1].c) )
+ spline->acceptableextrema = false;
+ }
+}
diff --git a/src/path/splinefit/splinerefigure.h b/src/path/splinefit/splinerefigure.h
new file mode 100644
index 0000000..110e2f9
--- /dev/null
+++ b/src/path/splinefit/splinerefigure.h
@@ -0,0 +1,9 @@
+// SPDX-License-Identifier: GPL-2.0-or-later
+#ifndef FONTFORGE_SPLINEREFIGURE_H
+#define FONTFORGE_SPLINEREFIGURE_H
+
+#include "splinefont.h"
+
+extern void SplineRefigure3(Spline *spline);
+
+#endif /* FONTFORGE_SPLINEREFIGURE_H */
generated by cgit v1.2.3 (git 2.39.1) at 2025-02-21 11:21:47 +0000