Index: src/Tesselation/tesselation.cpp
===================================================================
--- src/Tesselation/tesselation.cpp	(revision c2f30642d960bd67dc5b7bca66d267525a8f58ab)
+++ src/Tesselation/tesselation.cpp	(revision 75d42ad5b96258a67afde9aa02cf39b04a4f6059)
@@ -35,6 +35,9 @@
 #include "CodePatterns/MemDebug.hpp"
 
+#include <algorithm>
+#include <boost/foreach.hpp>
 #include <fstream>
 #include <iomanip>
+#include <iterator>
 #include <sstream>
 
@@ -731,6 +734,7 @@
  * @param *b second endpoint
  * @param n index of Tesselation::BLS giving the line with both endpoints
- */
-void Tesselation::AddTesselationLine(const Vector * const OptCenter, const BoundaryPointSet * const candidate, class BoundaryPointSet *a, class BoundaryPointSet *b, const int n)
+ * @return true - inserted a new line, false - present line used
+ */
+bool Tesselation::AddTesselationLine(const Vector * const OptCenter, const BoundaryPointSet * const candidate, class BoundaryPointSet *a, class BoundaryPointSet *b, const int n)
 {
   bool insertNewLine = true;
@@ -784,4 +788,6 @@
     AddExistingTesselationTriangleLine(WinningLine, n);
   }
+
+  return insertNewLine;
 }
 ;
@@ -3198,4 +3204,11 @@
 ;
 
+struct CloserToPiHalf
+{
+  bool operator()(double angle, double smallestangle) {
+    return (fabs(angle - M_PI/2.) < fabs(smallestangle - M_PI/2.));
+  }
+};
+
 /** Removes a boundary point from the envelope while keeping it closed.
  * We remove the old triangles connected to the point and re-create new triangles to close the surface following this ansatz:
@@ -3262,6 +3275,4 @@
 //  TriangleMap::iterator NumberRunner = Candidates.begin();
   TesselPointList::iterator StartNode, MiddleNode, EndNode;
-  double angle;
-  double smallestangle;
   Vector Point, Reference, OrthogonalVector;
   if (count > 2) { // less than three triangles, then nothing will be created
@@ -3275,8 +3286,9 @@
       // re-create all triangles by going through connected points list
       LineList NewLines;
+      typedef std::vector<double> angles_t;
+      angles_t angles;
       for (; !connectedPath->empty();) {
         // search middle node with widest angle to next neighbours
         EndNode = connectedPath->end();
-        smallestangle = 0.;
         for (MiddleNode = connectedPath->begin(); MiddleNode != connectedPath->end(); MiddleNode++) {
           LOG(1, "INFO: MiddleNode is " << **MiddleNode << ".");
@@ -3296,16 +3308,24 @@
           OrthogonalVector = ((*MiddleNode)->getPosition()) - OldPoint;
           OrthogonalVector.MakeNormalTo(Reference);
-          angle = GetAngle(Point, Reference, OrthogonalVector);
-          //if (angle < M_PI)  // no wrong-sided triangles, please?
-          if (fabs(angle - M_PI) < fabs(smallestangle - M_PI)) { // get straightest angle (i.e. construct those triangles with smallest area first)
-            smallestangle = angle;
-            EndNode = MiddleNode;
-          }
+          angles.push_back( GetAngle(Point, Reference, OrthogonalVector) );
         }
+        const angles_t::const_iterator maxiter = std::max_element(angles.begin(), angles.end());
+        angles_t::const_iterator miniter = angles.begin();
+        // distinguish between convex and nonconvex polygon
+        if (*maxiter > M_PI) {
+          // connectedPath is not convex
+          miniter = std::min_element(angles.begin(), angles.end());
+        } else {
+          // is convex
+          miniter = std::min_element(angles.begin(), angles.end(), CloserToPiHalf());
+        }
+        EndNode = connectedPath->begin();
+        std::advance(EndNode,
+            std::distance(const_cast<const angles_t &>(angles).begin(), miniter));
+        angles.clear();
+
         MiddleNode = EndNode;
-        if (MiddleNode == connectedPath->end()) {
-          ELOG(0, "CRITICAL: Could not find a smallest angle!");
-          performCriticalExit();
-        }
+        ASSERT( MiddleNode == connectedPath->end(),
+          "Tesselation::RemovePointFromTesselatedSurface() - Could not find a smallest angle!");
         StartNode = MiddleNode;
         if (StartNode == connectedPath->begin())
@@ -3347,6 +3367,7 @@
         LOG(3, "Adding new triangle lines.");
         AddTesselationLine(NULL, NULL, TPS[0], TPS[1], 0);
-        AddTesselationLine(NULL, NULL, TPS[0], TPS[2], 1);
-        NewLines.push_back(BLS[1]);
+        // line between start and end must be new (except for very last triangle)
+        if (AddTesselationLine(NULL, NULL, TPS[0], TPS[2], 1))
+          NewLines.push_back(BLS[1]);
         AddTesselationLine(NULL, NULL, TPS[1], TPS[2], 2);
         BTS = new class BoundaryTriangleSet(BLS, TrianglesOnBoundaryCount);
Index: src/Tesselation/tesselation.hpp
===================================================================
--- src/Tesselation/tesselation.hpp	(revision c2f30642d960bd67dc5b7bca66d267525a8f58ab)
+++ src/Tesselation/tesselation.hpp	(revision 75d42ad5b96258a67afde9aa02cf39b04a4f6059)
@@ -76,5 +76,5 @@
     void AddTesselationPoint(TesselPoint* Candidate, const int n);
     void SetTesselationPoint(TesselPoint* Candidate, const int n) const;
-    void AddTesselationLine(const Vector * OptCenter, const BoundaryPointSet * const candidate, class BoundaryPointSet *a, class BoundaryPointSet *b, const int n);
+    bool AddTesselationLine(const Vector * OptCenter, const BoundaryPointSet * const candidate, class BoundaryPointSet *a, class BoundaryPointSet *b, const int n);
     void AddNewTesselationTriangleLine(class BoundaryPointSet *a, class BoundaryPointSet *b, const int n);
     void AddExistingTesselationTriangleLine(class BoundaryLineSet *FindLine, int n);
