Changeset 27459a for molecuilder/src/tesselation.cpp

Timestamp:

Aug 18, 2009, 8:42:39 AM (16 years ago)

Author:

Frederik Heber <heber@…>

Children:

daf5d6

Parents:

3de1d2

Message:

Further bugfixes to three-step-procedure in convex envelope, not yet working.

• new functions StoreTrianglesinFile() and CalculateConcavityPerBoundaryPoint
• new functions Tesselation::PrintAll...() to print all points, lines and triangles currently on the boundary (simply to remove commented-out code)
• FIX: Tesselation::IsConvexRectangle(): just needs one closest point and simply checks the sign of the scalar product with both endpoints to determine in-/outside
• BUGFIX: Tesselation::PickFarthestofTwoBaselines(): NewOffset was wrong by the offset to the plane

File:

• molecuilder/src/tesselation.cpp (modified) (7 diffs)

molecuilder/src/tesselation.cpp

@@ -1751 +1751 @@
   class BoundaryPointSet *Spot = NULL;
   class BoundaryLineSet *OtherBase;
-  Vector *ClosestPoint[2];
+  Vector *ClosestPoint;
 
   int m=0;
@@ -1764 +1764 @@
 
   // get the closest point on each line to the other line
-  ClosestPoint[0] = GetClosestPointBetweenLine(out, Base, OtherBase);
-  ClosestPoint[1] = GetClosestPointBetweenLine(out, OtherBase, Base);
+  ClosestPoint = GetClosestPointBetweenLine(out, Base, OtherBase);
 
   // delete the temporary other base line
+  delete(ClosestPoint);
   delete(OtherBase);
 
   // get the distance vector from Base line to OtherBase line
-  Vector DistanceToIntersection, BaseLine;
+  Vector DistanceToIntersection[2], BaseLine;
+  double distance[2];
   BaseLine.CopyVector(Base->endpoints[1]->node->node);
   BaseLine.SubtractVector(Base->endpoints[0]->node->node);
-  DistanceToIntersection.CopyVector(ClosestPoint[0]);
-  DistanceToIntersection.SubtractVector(Base->endpoints[0]->node->node);
-  Spot = Base->endpoints[1];
-  if (BaseLine.ScalarProduct(&DistanceToIntersection) < -MYEPSILON) {
-    DistanceToIntersection.CopyVector(ClosestPoint[0]);
-    DistanceToIntersection.SubtractVector(Base->endpoints[1]->node->node);
-    Spot = Base->endpoints[0];
-    if (BaseLine.ScalarProduct(&DistanceToIntersection) < -MYEPSILON) {
-      *out << Verbose(3) << "ERROR: Something's very wrong here, both SKPs return negative?!" << endl;
-      return NULL;
-    }
-  }
-  if ((BaseLine.NormSquared() - DistanceToIntersection.NormSquared()) < -MYEPSILON) {  // distance is smaller: concave
-    *out << Verbose(3) << "INFO: Rectangle of triangles of base line " << *Base << " is concave: Base line squared length " << BaseLine.NormSquared() << " against Distance to intersection squared " << DistanceToIntersection.NormSquared() << "." << endl;
+  for (int i=0;i<2;i++) {
+    DistanceToIntersection[i].CopyVector(ClosestPoint);
+    DistanceToIntersection[i].SubtractVector(Base->endpoints[i]->node->node);
+    distance[i] = BaseLine.ScalarProduct(&DistanceToIntersection[i]);
+  }
+  if ((distance[0] * distance[1]) > MYEPSILON)  { // have same sign?
+    *out << Verbose(3) << "REJECT: Both SKPs have same sign: " << distance[0] << " and " << distance[1]  << ". " << *Base << " is concave." << endl;
+    if (distance[0] < distance[1]) {
+      Spot = Base->endpoints[0];
+    } else {
+      Spot = Base->endpoints[1];
+    }
     return Spot;
-  } else {  // base line is longer: convex
-    *out << Verbose(3) << "INFO: Rectangle of triangles of base line " << *Base << " is concave." << endl;
+  } else {  // different sign, i.e. we are in between
+    *out << Verbose(3) << "ACCEPT: Rectangle of triangles of base line " << *Base << " is convex." << endl;
     return NULL;
   }
@@ -1796 +1795 @@
 };
 
+void Tesselation::PrintAllBoundaryPoints(ofstream *out)
+{
+  // print all lines
+  *out << Verbose(1) << "Printing all boundary points for debugging:" << endl;
+  for (PointMap::iterator PointRunner = PointsOnBoundary.begin();PointRunner != PointsOnBoundary.end(); PointRunner++)
+    *out << Verbose(2) << *(PointRunner->second) << endl;
+};
+
+void Tesselation::PrintAllBoundaryLines(ofstream *out)
+{
+  // print all lines
+  *out << Verbose(1) << "Printing all boundary lines for debugging:" << endl;
+  for (LineMap::iterator LineRunner = LinesOnBoundary.begin(); LineRunner != LinesOnBoundary.end(); LineRunner++)
+    *out << Verbose(2) << *(LineRunner->second) << endl;
+};
+
+void Tesselation::PrintAllBoundaryTriangles(ofstream *out)
+{
+  // print all triangles
+  *out << Verbose(1) << "Printing all boundary triangles for debugging:" << endl;
+  for (TriangleMap::iterator TriangleRunner = TrianglesOnBoundary.begin(); TriangleRunner != TrianglesOnBoundary.end(); TriangleRunner++)
+    *out << Verbose(2) << *(TriangleRunner->second) << endl;
+};
 
 /** For a given boundary line \a *Base and its two triangles, picks the central baseline that is "higher".
@@ -1826 +1848 @@
   Distance.SubtractVector(ClosestPoint[0]);
 
-  // delete the temporary other base line
+  // delete the temporary other base line and the closest points
+  delete(ClosestPoint[0]);
+  delete(ClosestPoint[1]);
   delete(OtherBase);
 
@@ -1843 +1867 @@
       BaseLineNormal.AddVector(&(runner->second->NormalVector));
     }
-    BaseLineNormal.Scale(-1./2.); // has to point inside for BoundaryTriangleSet::GetNormalVector()
+    BaseLineNormal.Scale(1./2.);
 
     if (Distance.ScalarProduct(&BaseLineNormal) > MYEPSILON) { // Distance points outwards, hence OtherBase higher than Base -> flip
@@ -1866 +1890 @@
   // construct the plane of the two baselines (i.e. take both their directional vectors)
   Vector Normal;
-  Vector OtherBaseline;
-  Normal.CopyVector(Base->endpoints[1]->node->node);
-  Normal.SubtractVector(Base->endpoints[0]->node->node);
+  Vector Baseline, OtherBaseline;
+  Baseline.CopyVector(Base->endpoints[1]->node->node);
+  Baseline.SubtractVector(Base->endpoints[0]->node->node);
   OtherBaseline.CopyVector(OtherBase->endpoints[1]->node->node);
   OtherBaseline.SubtractVector(OtherBase->endpoints[0]->node->node);
+  Normal.CopyVector(&Baseline);
   Normal.VectorProduct(&OtherBaseline);
   Normal.Normalize();
-
-  // project one offset point of OtherBase onto this plane
+  *out << Verbose(4) << "First direction is " << Baseline << ", second direction is " << OtherBaseline << ", normal of intersection plane is " << Normal << "." << endl;
+
+  // project one offset point of OtherBase onto this plane (and add plane offset vector)
   Vector NewOffset;
   NewOffset.CopyVector(OtherBase->endpoints[0]->node->node);
   NewOffset.ProjectOntoPlane(&Normal);
-  Vector NewDirection;
-  NewDirection.CopyVector(OtherBase->endpoints[1]->node->node);
-  NewDirection.ProjectOntoPlane(&Normal);
+  NewOffset.AddVector(Base->endpoints[0]->node->node);
 
   // calculate the intersection between this projected baseline and Base
   Vector *Intersection = new Vector;
-  Intersection->GetIntersectionOfTwoLinesOnPlane(out, Base->endpoints[0]->node->node, Base->endpoints[1]->node->node, &NewOffset, &NewDirection, &Normal);
+  Intersection->GetIntersectionOfTwoLinesOnPlane(out, Base->endpoints[0]->node->node, Base->endpoints[1]->node->node, &NewOffset, &OtherBaseline, &Normal);
+  Normal.CopyVector(Intersection);
+  Normal.SubtractVector(Base->endpoints[0]->node->node);
+  *out << Verbose(3) << "Found closest point on " << *Base << " at " << *Intersection << ", factor in line is " << fabs(Normal.ScalarProduct(&Baseline)/Baseline.NormSquared()) << "." << endl;
 
   return Intersection;
@@ -2938 +2965 @@
 double getAngle(const Vector &point, const Vector &reference, const Vector OrthogonalVector)
 {
-  if (reference.IsNull())
+  if (reference.IsZero())
     return M_PI;
 
   // calculate both angles and correct with in-plane vector
-  if (point.IsNull())
+  if (point.IsZero())
     return M_PI;
   double phi = point.Angle(&reference);