Code

update 2geom to rev 1114 and fix error.
authorjohanengelen <johanengelen@users.sourceforge.net>
Thu, 30 Aug 2007 19:47:20 +0000 (19:47 +0000)
committerjohanengelen <johanengelen@users.sourceforge.net>
Thu, 30 Aug 2007 19:47:20 +0000 (19:47 +0000)
src/2geom/Makefile_insert
src/2geom/d2-sbasis.cpp
src/2geom/d2-sbasis.h
src/2geom/d2.cpp [deleted file]

index 4a87c1a93d3a21a69bce3f0115d71516ec149c3d..4e8fb13d8f57ef0aba768afafb023d756e6187b8 100644 (file)
@@ -24,7 +24,6 @@
        2geom/crossing.cpp      \\r
        2geom/d2-sbasis.cpp     \\r
        2geom/d2-sbasis.h       \\r
-       2geom/d2.cpp    \\r
        2geom/d2.h      \\r
        2geom/geom.cpp  \\r
        2geom/geom.h    \\r
index b68e3e6c558e61d5cce8e9432efc64626e681f22..c5c005c24b711f152d62622e9851245e5738eb00 100644 (file)
  * This is due to the trickinesses of template submatching. */
 
 namespace Geom {
-\r
-SBasis L2(D2<SBasis> const & a, unsigned k) { return sqrt(dot(a, a), k); }\r
-\r
-D2<SBasis> multiply(Linear const & a, D2<SBasis> const & b) {\r
-    return D2<SBasis>(multiply(a, b[X]), multiply(a, b[Y]));\r
-}\r
-\r
-D2<SBasis> multiply(SBasis const & a, D2<SBasis> const & b) {\r
-    return D2<SBasis>(multiply(a, b[X]), multiply(a, b[Y]));\r
-}\r
-\r
-D2<SBasis> truncate(D2<SBasis> const & a, unsigned terms) {\r
-    return D2<SBasis>(truncate(a[X], terms), truncate(a[Y], terms));\r
-}\r
-\r
-unsigned sbasis_size(D2<SBasis> const & a) {\r
-    return std::max((unsigned) a[0].size(), (unsigned) a[1].size());\r
-}\r
-\r
-//TODO: Is this sensical? shouldn't it be like pythagorean or something?\r
-double tail_error(D2<SBasis> const & a, unsigned tail) {\r
-    return std::max(a[0].tailError(tail), a[1].tailError(tail));\r
-}\r
-\r
-Piecewise<D2<SBasis> > sectionize(D2<Piecewise<SBasis> > const &a) {\r
-    Piecewise<SBasis> x = partition(a[0], a[1].cuts), y = partition(a[1], a[0].cuts);\r
-    assert(x.size() == y.size());\r
-    Piecewise<D2<SBasis> > ret;\r
-    for(unsigned i = 0; i < x.size(); i++)\r
-        ret.push_seg(D2<SBasis>(x[i], y[i]));\r
-    ret.cuts.insert(ret.cuts.end(), x.cuts.begin(), x.cuts.end());\r
-    return ret;\r
-}\r
-\r
-D2<Piecewise<SBasis> > make_cuts_independant(Piecewise<D2<SBasis> > const &a) {\r
-    D2<Piecewise<SBasis> > ret;\r
-    for(unsigned d = 0; d < 2; d++) {\r
-        for(unsigned i = 0; i < a.size(); i++)\r
-            ret[d].push_seg(a[i][d]);\r
-        ret[d].cuts.insert(ret[d].cuts.end(), a.cuts.begin(), a.cuts.end());\r
-    }\r
-    return ret;\r
-}\r
-\r
-Piecewise<D2<SBasis> > rot90(Piecewise<D2<SBasis> > const &M){\r
-  Piecewise<D2<SBasis> > result;\r
-  if (M.empty()) return M;\r
-  result.push_cut(M.cuts[0]);\r
-  for (unsigned i=0; i<M.size(); i++){\r
-    result.push(rot90(M[i]),M.cuts[i+1]);\r
-  }\r
-  return result;\r
-}\r
-\r
-Piecewise<SBasis> dot(Piecewise<D2<SBasis> > const &a, \r
-                     Piecewise<D2<SBasis> > const &b){\r
-  Piecewise<SBasis > result;\r
-  if (a.empty() || b.empty()) return result;\r
-  Piecewise<D2<SBasis> > aa = partition(a,b.cuts);\r
-  Piecewise<D2<SBasis> > bb = partition(b,a.cuts);\r
-\r
-  result.push_cut(aa.cuts.front());\r
-  for (unsigned i=0; i<aa.size(); i++){\r
-    result.push(dot(aa.segs[i],bb.segs[i]),aa.cuts[i+1]);\r
-  }\r
-  return result;\r
-}\r
-\r
-Piecewise<SBasis> cross(Piecewise<D2<SBasis> > const &a, \r
-                       Piecewise<D2<SBasis> > const &b){\r
-  Piecewise<SBasis > result;\r
-  if (a.empty() || b.empty()) return result;\r
-  Piecewise<D2<SBasis> > aa = partition(a,b.cuts);\r
-  Piecewise<D2<SBasis> > bb = partition(b,a.cuts);\r
-\r
-  result.push_cut(aa.cuts.front());\r
-  for (unsigned i=0; i<a.size(); i++){\r
-    result.push(cross(aa.segs[i],bb.segs[i]),aa.cuts[i+1]);\r
-  }\r
-  return result;\r
-}\r
-\r
+
+SBasis L2(D2<SBasis> const & a, unsigned k) { return sqrt(dot(a, a), k); }
+
+D2<SBasis> multiply(Linear const & a, D2<SBasis> const & b) {
+    return D2<SBasis>(multiply(a, b[X]), multiply(a, b[Y]));
+}
+
+D2<SBasis> multiply(SBasis const & a, D2<SBasis> const & b) {
+    return D2<SBasis>(multiply(a, b[X]), multiply(a, b[Y]));
+}
+
+D2<SBasis> truncate(D2<SBasis> const & a, unsigned terms) {
+    return D2<SBasis>(truncate(a[X], terms), truncate(a[Y], terms));
+}
+
+unsigned sbasis_size(D2<SBasis> const & a) {
+    return std::max((unsigned) a[0].size(), (unsigned) a[1].size());
+}
+
+//TODO: Is this sensical? shouldn't it be like pythagorean or something?
+double tail_error(D2<SBasis> const & a, unsigned tail) {
+    return std::max(a[0].tailError(tail), a[1].tailError(tail));
+}
+
+Piecewise<D2<SBasis> > sectionize(D2<Piecewise<SBasis> > const &a) {
+    Piecewise<SBasis> x = partition(a[0], a[1].cuts), y = partition(a[1], a[0].cuts);
+    assert(x.size() == y.size());
+    Piecewise<D2<SBasis> > ret;
+    for(unsigned i = 0; i < x.size(); i++)
+        ret.push_seg(D2<SBasis>(x[i], y[i]));
+    ret.cuts.insert(ret.cuts.end(), x.cuts.begin(), x.cuts.end());
+    return ret;
+}
+
+D2<Piecewise<SBasis> > make_cuts_independant(Piecewise<D2<SBasis> > const &a) {
+    D2<Piecewise<SBasis> > ret;
+    for(unsigned d = 0; d < 2; d++) {
+        for(unsigned i = 0; i < a.size(); i++)
+            ret[d].push_seg(a[i][d]);
+        ret[d].cuts.insert(ret[d].cuts.end(), a.cuts.begin(), a.cuts.end());
+    }
+    return ret;
+}
+
+Piecewise<D2<SBasis> > rot90(Piecewise<D2<SBasis> > const &M){
+  Piecewise<D2<SBasis> > result;
+  if (M.empty()) return M;
+  result.push_cut(M.cuts[0]);
+  for (unsigned i=0; i<M.size(); i++){
+    result.push(rot90(M[i]),M.cuts[i+1]);
+  }
+  return result;
+}
+
+Piecewise<SBasis> dot(Piecewise<D2<SBasis> > const &a, 
+                     Piecewise<D2<SBasis> > const &b){
+  Piecewise<SBasis > result;
+  if (a.empty() || b.empty()) return result;
+  Piecewise<D2<SBasis> > aa = partition(a,b.cuts);
+  Piecewise<D2<SBasis> > bb = partition(b,a.cuts);
+
+  result.push_cut(aa.cuts.front());
+  for (unsigned i=0; i<aa.size(); i++){
+    result.push(dot(aa.segs[i],bb.segs[i]),aa.cuts[i+1]);
+  }
+  return result;
+}
+
+Piecewise<SBasis> cross(Piecewise<D2<SBasis> > const &a, 
+                       Piecewise<D2<SBasis> > const &b){
+  Piecewise<SBasis > result;
+  if (a.empty() || b.empty()) return result;
+  Piecewise<D2<SBasis> > aa = partition(a,b.cuts);
+  Piecewise<D2<SBasis> > bb = partition(b,a.cuts);
+
+  result.push_cut(aa.cuts.front());
+  for (unsigned i=0; i<a.size(); i++){
+    result.push(cross(aa.segs[i],bb.segs[i]),aa.cuts[i+1]);
+  }
+  return result;
+}
+
+Piecewise<D2<SBasis> > operator*(Piecewise<D2<SBasis> > const &a, Matrix const &m) {
+  Piecewise<D2<SBasis> > result;
+  if(a.empty()) return result;
+  result.push_cut(a.cuts[0]);
+  for (unsigned i = 0; i < a.size(); i++) {
+    result.push(a[i] * m, a.cuts[i+1]);
+  }
+  return result;
+}
+
 /* Replaced by remove_short_cuts in piecewise.h
-//this recursively removes the shortest cut interval until none is shorter than tol.\r
-//TODO: code this in a more efficient way!\r
-Piecewise<D2<SBasis> > remove_short_cuts(Piecewise<D2<SBasis> > const &f, double tol){\r
-    double min = tol;\r
-    unsigned idx = f.size();\r
-    for(unsigned i=0; i<f.size(); i++){\r
-        if (min > f.cuts[i+1]-f.cuts[i]){\r
-            min = f.cuts[i+1]-f.cuts[i];\r
-            idx = int(i);\r
-        }\r
-    }\r
-    if (idx==f.size()){\r
-        return f;\r
-    }\r
-    if (f.size()==1) {\r
-        //removing this seg would result in an empty pw<d2<sb>>...\r
-        return f;\r
-    }\r
-    Piecewise<D2<SBasis> > new_f=f;\r
-    for (int dim=0; dim<2; dim++){\r
-        double v = Hat(f.segs.at(idx)[dim][0]);\r
-        //TODO: what about closed curves?\r
-        if (idx>0 && f.segs.at(idx-1).at1()==f.segs.at(idx).at0()) \r
-            new_f.segs.at(idx-1)[dim][0][1] = v;\r
-        if (idx<f.size() && f.segs.at(idx+1).at0()==f.segs.at(idx).at1()) \r
-            new_f.segs.at(idx+1)[dim][0][0] = v;\r
-    }\r
-    double t = (f.cuts.at(idx)+f.cuts.at(idx+1))/2;\r
-    new_f.cuts.at(idx+1) = t;    \r
-    \r
-    new_f.segs.erase(new_f.segs.begin()+idx);\r
-    new_f.cuts.erase(new_f.cuts.begin()+idx);        \r
-    return remove_short_cuts(new_f, tol);\r
-}\r
+//this recursively removes the shortest cut interval until none is shorter than tol.
+//TODO: code this in a more efficient way!
+Piecewise<D2<SBasis> > remove_short_cuts(Piecewise<D2<SBasis> > const &f, double tol){
+    double min = tol;
+    unsigned idx = f.size();
+    for(unsigned i=0; i<f.size(); i++){
+        if (min > f.cuts[i+1]-f.cuts[i]){
+            min = f.cuts[i+1]-f.cuts[i];
+            idx = int(i);
+        }
+    }
+    if (idx==f.size()){
+        return f;
+    }
+    if (f.size()==1) {
+        //removing this seg would result in an empty pw<d2<sb>>...
+        return f;
+    }
+    Piecewise<D2<SBasis> > new_f=f;
+    for (int dim=0; dim<2; dim++){
+        double v = Hat(f.segs.at(idx)[dim][0]);
+        //TODO: what about closed curves?
+        if (idx>0 && f.segs.at(idx-1).at1()==f.segs.at(idx).at0()) 
+            new_f.segs.at(idx-1)[dim][0][1] = v;
+        if (idx<f.size() && f.segs.at(idx+1).at0()==f.segs.at(idx).at1()) 
+            new_f.segs.at(idx+1)[dim][0][0] = v;
+    }
+    double t = (f.cuts.at(idx)+f.cuts.at(idx+1))/2;
+    new_f.cuts.at(idx+1) = t;    
+    
+    new_f.segs.erase(new_f.segs.begin()+idx);
+    new_f.cuts.erase(new_f.cuts.begin()+idx);        
+    return remove_short_cuts(new_f, tol);
+}
 */
-\r
-//if tol>0, only force continuity where the jump is smaller than tol.\r
-Piecewise<D2<SBasis> > force_continuity(Piecewise<D2<SBasis> > const &f, \r
-                                        double tol,\r
-                                        bool closed){\r
-    if (f.size()==0) return f;\r
-    Piecewise<D2<SBasis> > result=f;\r
-    unsigned cur   = (closed)? 0:1;\r
-    unsigned prev  = (closed)? f.size()-1:0;\r
-    while(cur<f.size()){\r
-        Point pt0 = f.segs[prev].at1();\r
-        Point pt1 = f.segs[cur ].at0();\r
-        if (tol<=0 || L2sq(pt0-pt1)<tol*tol){\r
-            pt0 = (pt0+pt1)/2;\r
-            for (unsigned dim=0; dim<2; dim++){\r
-                result.segs[prev][dim][0][1]=pt0[dim];\r
-                result.segs[cur ][dim][0][0]=pt0[dim];\r
-            }\r
-        }\r
-        prev = cur++;\r
-    }\r
-    return result;\r
-}\r\r
+
+//if tol>0, only force continuity where the jump is smaller than tol.
+Piecewise<D2<SBasis> > force_continuity(Piecewise<D2<SBasis> > const &f, 
+                                        double tol,
+                                        bool closed){
+    if (f.size()==0) return f;
+    Piecewise<D2<SBasis> > result=f;
+    unsigned cur   = (closed)? 0:1;
+    unsigned prev  = (closed)? f.size()-1:0;
+    while(cur<f.size()){
+        Point pt0 = f.segs[prev].at1();
+        Point pt1 = f.segs[cur ].at0();
+        if (tol<=0 || L2sq(pt0-pt1)<tol*tol){
+            pt0 = (pt0+pt1)/2;
+            for (unsigned dim=0; dim<2; dim++){
+                result.segs[prev][dim][0][1]=pt0[dim];
+                result.segs[cur ][dim][0][0]=pt0[dim];
+            }
+        }
+        prev = cur++;
+    }
+    return result;
+}
 }
index dab2559ca56d29af52b7fb626bb64de271b2adea..95f1ca0dd3562f58150d3da9731ed91c35198695 100644 (file)
@@ -4,82 +4,85 @@
 #ifndef __2GEOM_SBASIS_CURVE_H
 #define __2GEOM_SBASIS_CURVE_H
 
-#include "sbasis.h"\r
-#include "sbasis-2d.h"\r
+#include "sbasis.h"
+#include "sbasis-2d.h"
 #include "piecewise.h"
+#include "matrix.h"
 
 //TODO: implement intersect
 
 namespace Geom {
-\r
-inline D2<SBasis> compose(D2<SBasis> const & a, SBasis const & b) {\r
-    return D2<SBasis>(compose(a[X], b), compose(a[Y], b));\r
-}\r
-\r
-SBasis L2(D2<SBasis> const & a, unsigned k);\r
-double L2(D2<double> const & a);\r
-\r
-D2<SBasis> multiply(Linear const & a, D2<SBasis> const & b);\r
-inline D2<SBasis> operator*(Linear const & a, D2<SBasis> const & b) { return multiply(a, b); }\r
-D2<SBasis> multiply(SBasis const & a, D2<SBasis> const & b);\r
-inline D2<SBasis> operator*(SBasis const & a, D2<SBasis> const & b) { return multiply(a, b); }\r
-D2<SBasis> truncate(D2<SBasis> const & a, unsigned terms);\r
-\r
-unsigned sbasis_size(D2<SBasis> const & a);\r
-double tail_error(D2<SBasis> const & a, unsigned tail);\r
-\r
-//Piecewise<D2<SBasis> > specific decls:\r
-\r
-Piecewise<D2<SBasis> > sectionize(D2<Piecewise<SBasis> > const &a);\r
-D2<Piecewise<SBasis> > make_cuts_independant(Piecewise<D2<SBasis> > const &a);\r
-Piecewise<D2<SBasis> > rot90(Piecewise<D2<SBasis> > const &a);\r
-Piecewise<SBasis> dot(Piecewise<D2<SBasis> > const &a, Piecewise<D2<SBasis> > const &b);\r
-Piecewise<SBasis> cross(Piecewise<D2<SBasis> > const &a, Piecewise<D2<SBasis> > const &b);\r
-\r
-Piecewise<D2<SBasis> > force_continuity(Piecewise<D2<SBasis> > const &f, \r
-                                        double tol=0,\r
-                                        bool closed=false);\r
-\r
-class CoordIterator\r
-: public std::iterator<std::input_iterator_tag, SBasis const>\r
-{\r
-public:\r
-  CoordIterator(std::vector<D2<SBasis> >::const_iterator const &iter, unsigned d) : impl_(iter), ix_(d) {}\r
-\r
-  inline bool operator==(CoordIterator const &other) { return other.impl_ == impl_; }\r
-  inline bool operator!=(CoordIterator const &other) { return other.impl_ != impl_; }\r
-\r
-  inline SBasis operator*() const {\r
-        return (*impl_)[ix_];\r
-  }\r
-\r
-  inline CoordIterator &operator++() {\r
-    ++impl_;\r
-    return *this;\r
-  }\r
-  inline CoordIterator operator++(int) {\r
-    CoordIterator old=*this;\r
-    ++(*this);\r
-    return old;\r
-  }\r
-\r
-private:\r
-  std::vector<D2<SBasis> >::const_iterator impl_;\r
-  unsigned ix_;\r
-};\r
-\r
-inline CoordIterator iterateCoord(Piecewise<D2<SBasis> > const &a, unsigned d) {\r
-    return CoordIterator(a.segs.begin(), d);\r
-}\r
-\r
-//bounds specializations with order\r
-inline Rect bounds_fast(D2<SBasis> const & s, unsigned order=0) {\r
-    return Rect(bounds_fast(s[X], order),\r
-                bounds_fast(s[Y], order));\r
-}\r
-inline Rect bounds_local(D2<SBasis> const & s, Interval i, unsigned order=0) {\r
-    return Rect(bounds_local(s[X], i, order),\r
-                bounds_local(s[Y], i, order));\r
+
+inline D2<SBasis> compose(D2<SBasis> const & a, SBasis const & b) {
+    return D2<SBasis>(compose(a[X], b), compose(a[Y], b));
+}
+
+SBasis L2(D2<SBasis> const & a, unsigned k);
+double L2(D2<double> const & a);
+
+D2<SBasis> multiply(Linear const & a, D2<SBasis> const & b);
+inline D2<SBasis> operator*(Linear const & a, D2<SBasis> const & b) { return multiply(a, b); }
+D2<SBasis> multiply(SBasis const & a, D2<SBasis> const & b);
+inline D2<SBasis> operator*(SBasis const & a, D2<SBasis> const & b) { return multiply(a, b); }
+D2<SBasis> truncate(D2<SBasis> const & a, unsigned terms);
+
+unsigned sbasis_size(D2<SBasis> const & a);
+double tail_error(D2<SBasis> const & a, unsigned tail);
+
+//Piecewise<D2<SBasis> > specific decls:
+
+Piecewise<D2<SBasis> > sectionize(D2<Piecewise<SBasis> > const &a);
+D2<Piecewise<SBasis> > make_cuts_independant(Piecewise<D2<SBasis> > const &a);
+Piecewise<D2<SBasis> > rot90(Piecewise<D2<SBasis> > const &a);
+Piecewise<SBasis> dot(Piecewise<D2<SBasis> > const &a, Piecewise<D2<SBasis> > const &b);
+Piecewise<SBasis> cross(Piecewise<D2<SBasis> > const &a, Piecewise<D2<SBasis> > const &b);
+
+Piecewise<D2<SBasis> > operator*(Piecewise<D2<SBasis> > const &a, Matrix const &m);
+
+Piecewise<D2<SBasis> > force_continuity(Piecewise<D2<SBasis> > const &f, 
+                                        double tol=0,
+                                        bool closed=false);
+
+class CoordIterator
+: public std::iterator<std::input_iterator_tag, SBasis const>
+{
+public:
+  CoordIterator(std::vector<D2<SBasis> >::const_iterator const &iter, unsigned d) : impl_(iter), ix_(d) {}
+
+  inline bool operator==(CoordIterator const &other) { return other.impl_ == impl_; }
+  inline bool operator!=(CoordIterator const &other) { return other.impl_ != impl_; }
+
+  inline SBasis operator*() const {
+        return (*impl_)[ix_];
+  }
+
+  inline CoordIterator &operator++() {
+    ++impl_;
+    return *this;
+  }
+  inline CoordIterator operator++(int) {
+    CoordIterator old=*this;
+    ++(*this);
+    return old;
+  }
+
+private:
+  std::vector<D2<SBasis> >::const_iterator impl_;
+  unsigned ix_;
+};
+
+inline CoordIterator iterateCoord(Piecewise<D2<SBasis> > const &a, unsigned d) {
+    return CoordIterator(a.segs.begin(), d);
+}
+
+//bounds specializations with order
+inline Rect bounds_fast(D2<SBasis> const & s, unsigned order=0) {
+    return Rect(bounds_fast(s[X], order),
+                bounds_fast(s[Y], order));
+}
+inline Rect bounds_local(D2<SBasis> const & s, Interval i, unsigned order=0) {
+    return Rect(bounds_local(s[X], i, order),
+                bounds_local(s[Y], i, order));
 }
 
 }
diff --git a/src/2geom/d2.cpp b/src/2geom/d2.cpp
deleted file mode 100644 (file)
index 8653806..0000000
+++ /dev/null
@@ -1,177 +0,0 @@
-/*\r
- * d2.cpp - Lifts one dimensional objects into 2d \r
- *\r
- * Copyright 2007 Michael Sloan <mgsloan@gmail.com>\r
- *\r
- * This library is free software; you can redistribute it and/or\r
- * modify it either under the terms of the GNU Lesser General Public\r
- * License version 2.1 as published by the Free Software Foundation\r
- * (the "LGPL") or, at your option, under the terms of the Mozilla\r
- * Public License Version 1.1 (the "MPL"). If you do not alter this\r
- * notice, a recipient may use your version of this file under either\r
- * the MPL or the LGPL.\r
- *\r
- * You should have received a copy of the LGPL along with this library\r
- * in the file COPYING-LGPL-2.1; if not, output to the Free Software\r
- * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA\r
- * You should have received a copy of the MPL along with this library\r
- * in the file COPYING-MPL-1.1\r
- *\r
- * The contents of this file are subject to the Mozilla Public License\r
- * Version 1.1 (the "License"); you may not use this file except in\r
- * compliance with the License. You may obtain a copy of the License at\r
- * http://www.mozilla.org/MPL/\r
- *\r
- * This software is distributed on an "AS IS" basis, WITHOUT WARRANTY\r
- * OF ANY KIND, either express or implied. See the LGPL or the MPL for\r
- * the specific language governing rights and limitations.\r
- *\r
- */\r
-\r
-#include "d2.h"\r
-\r
-namespace Geom {\r
-\r
-SBasis L2(D2<SBasis> const & a, unsigned k) { return sqrt(dot(a, a), k); }\r
-double L2(D2<double> const & a) { return hypot(a[0], a[1]); }\r
-\r
-D2<SBasis> multiply(Linear const & a, D2<SBasis> const & b) {\r
-    return D2<SBasis>(multiply(a, b[X]), multiply(a, b[Y]));\r
-}\r
-\r
-D2<SBasis> multiply(SBasis const & a, D2<SBasis> const & b) {\r
-    return D2<SBasis>(multiply(a, b[X]), multiply(a, b[Y]));\r
-}\r
-\r
-D2<SBasis> truncate(D2<SBasis> const & a, unsigned terms) {\r
-    return D2<SBasis>(truncate(a[X], terms), truncate(a[Y], terms));\r
-}\r
-\r
-unsigned sbasis_size(D2<SBasis> const & a) {\r
-    return std::max((unsigned) a[0].size(), (unsigned) a[1].size());\r
-}\r
-\r
-//TODO: Is this sensical? shouldn't it be like pythagorean or something?\r
-double tail_error(D2<SBasis> const & a, unsigned tail) {\r
-    return std::max(a[0].tailError(tail), a[1].tailError(tail));\r
-}\r
-\r
-Piecewise<D2<SBasis> > sectionize(D2<Piecewise<SBasis> > const &a) {\r
-    Piecewise<SBasis> x = partition(a[0], a[1].cuts), y = partition(a[1], a[0].cuts);\r
-    assert(x.size() == y.size());\r
-    Piecewise<D2<SBasis> > ret;\r
-    for(unsigned i = 0; i < x.size(); i++)\r
-        ret.push_seg(D2<SBasis>(x[i], y[i]));\r
-    ret.cuts.insert(ret.cuts.end(), x.cuts.begin(), x.cuts.end());\r
-    return ret;\r
-}\r
-\r
-D2<Piecewise<SBasis> > make_cuts_independant(Piecewise<D2<SBasis> > const &a) {\r
-    D2<Piecewise<SBasis> > ret;\r
-    for(unsigned d = 0; d < 2; d++) {\r
-        for(unsigned i = 0; i < a.size(); i++)\r
-            ret[d].push_seg(a[i][d]);\r
-        ret[d].cuts.insert(ret[d].cuts.end(), a.cuts.begin(), a.cuts.end());\r
-    }\r
-    return ret;\r
-}\r
-\r
-Piecewise<D2<SBasis> > rot90(Piecewise<D2<SBasis> > const &M){\r
-  Piecewise<D2<SBasis> > result;\r
-  if (M.empty()) return M;\r
-  result.push_cut(M.cuts[0]);\r
-  for (unsigned i=0; i<M.size(); i++){\r
-    result.push(rot90(M[i]),M.cuts[i+1]);\r
-  }\r
-  return result;\r
-}\r
-\r
-Piecewise<SBasis> dot(Piecewise<D2<SBasis> > const &a, \r
-                     Piecewise<D2<SBasis> > const &b){\r
-  Piecewise<SBasis > result;\r
-  if (a.empty() || b.empty()) return result;\r
-  Piecewise<D2<SBasis> > aa = partition(a,b.cuts);\r
-  Piecewise<D2<SBasis> > bb = partition(b,a.cuts);\r
-\r
-  result.push_cut(aa.cuts.front());\r
-  for (unsigned i=0; i<aa.size(); i++){\r
-    result.push(dot(aa.segs[i],bb.segs[i]),aa.cuts[i+1]);\r
-  }\r
-  return result;\r
-}\r
-\r
-Piecewise<SBasis> cross(Piecewise<D2<SBasis> > const &a, \r
-                       Piecewise<D2<SBasis> > const &b){\r
-  Piecewise<SBasis > result;\r
-  if (a.empty() || b.empty()) return result;\r
-  Piecewise<D2<SBasis> > aa = partition(a,b.cuts);\r
-  Piecewise<D2<SBasis> > bb = partition(b,a.cuts);\r
-\r
-  result.push_cut(aa.cuts.front());\r
-  for (unsigned i=0; i<a.size(); i++){\r
-    result.push(cross(aa.segs[i],bb.segs[i]),aa.cuts[i+1]);\r
-  }\r
-  return result;\r
-}\r
-\r
-/* Replaced by remove_short_cuts in piecewise.h
-//this recursively removes the shortest cut interval until none is shorter than tol.\r
-//TODO: code this in a more efficient way!\r
-Piecewise<D2<SBasis> > remove_short_cuts(Piecewise<D2<SBasis> > const &f, double tol){\r
-    double min = tol;\r
-    unsigned idx = f.size();\r
-    for(unsigned i=0; i<f.size(); i++){\r
-        if (min > f.cuts[i+1]-f.cuts[i]){\r
-            min = f.cuts[i+1]-f.cuts[i];\r
-            idx = int(i);\r
-        }\r
-    }\r
-    if (idx==f.size()){\r
-        return f;\r
-    }\r
-    if (f.size()==1) {\r
-        //removing this seg would result in an empty pw<d2<sb>>...\r
-        return f;\r
-    }\r
-    Piecewise<D2<SBasis> > new_f=f;\r
-    for (int dim=0; dim<2; dim++){\r
-        double v = Hat(f.segs.at(idx)[dim][0]);\r
-        //TODO: what about closed curves?\r
-        if (idx>0 && f.segs.at(idx-1).at1()==f.segs.at(idx).at0()) \r
-            new_f.segs.at(idx-1)[dim][0][1] = v;\r
-        if (idx<f.size() && f.segs.at(idx+1).at0()==f.segs.at(idx).at1()) \r
-            new_f.segs.at(idx+1)[dim][0][0] = v;\r
-    }\r
-    double t = (f.cuts.at(idx)+f.cuts.at(idx+1))/2;\r
-    new_f.cuts.at(idx+1) = t;    \r
-    \r
-    new_f.segs.erase(new_f.segs.begin()+idx);\r
-    new_f.cuts.erase(new_f.cuts.begin()+idx);        \r
-    return remove_short_cuts(new_f, tol);\r
-}\r
-*/
-\r
-//if tol>0, only force continuity where the jump is smaller than tol.\r
-Piecewise<D2<SBasis> > force_continuity(Piecewise<D2<SBasis> > const &f, \r
-                                        double tol,\r
-                                        bool closed){\r
-    if (f.size()==0) return f;\r
-    Piecewise<D2<SBasis> > result=f;\r
-    unsigned cur   = (closed)? 0:1;\r
-    unsigned prev  = (closed)? f.size()-1:0;\r
-    while(cur<f.size()){\r
-        Point pt0 = f.segs[prev].at1();\r
-        Point pt1 = f.segs[cur ].at0();\r
-        if (tol<=0 || L2sq(pt0-pt1)<tol*tol){\r
-            pt0 = (pt0+pt1)/2;\r
-            for (unsigned dim=0; dim<2; dim++){\r
-                result.segs[prev][dim][0][1]=pt0[dim];\r
-                result.segs[cur ][dim][0][0]=pt0[dim];\r
-            }\r
-        }\r
-        prev = cur++;\r
-    }\r
-    return result;\r
-}\r
-  \r
-};\r