From 0304aa3ab9f9ce651d0f9bb3f171a6c8dc901f00 Mon Sep 17 00:00:00 2001
From: Radek Oslejsek <oslejsek@fi.muni.cz>
Date: Sun, 21 Mar 2021 08:45:37 +0100
Subject: [PATCH] Refactored and improved the code for symmetry estimation

---
 .../cz/fidentis/analyst/EfficiencyTests.java  |  41 ++-
 .../analyst/symmetry/ApproxSymmetryPlane.java | 152 ++++++++-
 .../cz/fidentis/analyst/symmetry/Plane.java   |  42 +++
 .../analyst/symmetry/SignificantPoints.java   |  96 ++++++
 .../analyst/symmetry/SymmetryEstimator.java   | 294 ++++--------------
 .../visitors/mesh/CurvatureVisitor.java       |  55 ++--
 .../mesh/GaussianCurvatureVisitor.java        |  14 +-
 .../visitors/mesh/MaxCurvatureVisitor.java    |  13 +-
 .../analyst/mesh/core/MeshFacetImpl.java      |   4 +-
 9 files changed, 425 insertions(+), 286 deletions(-)
 create mode 100644 Comparison/src/main/java/cz/fidentis/analyst/symmetry/SignificantPoints.java

diff --git a/Comparison/src/main/java/cz/fidentis/analyst/EfficiencyTests.java b/Comparison/src/main/java/cz/fidentis/analyst/EfficiencyTests.java
index d0531829..c4de9f35 100644
--- a/Comparison/src/main/java/cz/fidentis/analyst/EfficiencyTests.java
+++ b/Comparison/src/main/java/cz/fidentis/analyst/EfficiencyTests.java
@@ -6,9 +6,10 @@ import cz.fidentis.analyst.mesh.core.MeshFacet;
 import cz.fidentis.analyst.symmetry.Config;
 import cz.fidentis.analyst.symmetry.Plane;
 import cz.fidentis.analyst.symmetry.SymmetryEstimator;
-import cz.fidentis.analyst.visitors.mesh.BoundingBoxVisitor;
+import cz.fidentis.analyst.visitors.mesh.GaussianCurvatureVisitor;
 import cz.fidentis.analyst.visitors.mesh.HausdorffDistance;
 import cz.fidentis.analyst.visitors.mesh.HausdorffDistance.Strategy;
+import cz.fidentis.analyst.visitors.mesh.MaxCurvatureVisitor;
 import java.io.File;
 import java.io.IOException;
 import java.util.List;
@@ -38,7 +39,11 @@ public class EfficiencyTests {
         face1 = new HumanFace(faceFile2);
         face2 = new HumanFace(faceFile4);
         
-        System.out.println(measureSymmetryPlane(face1, true) + "\tmsec:\tSymmetry plane computation");
+        face1.getMeshModel().compute(new GaussianCurvatureVisitor(false)); // initialize everything, then measure
+        
+        System.out.println(measureCurvature(face1, new GaussianCurvatureVisitor(false), false) + "\tmsec:\tGaussian curvature");
+        System.out.println(measureCurvature(face1, new MaxCurvatureVisitor(false), false) + "\tmsec:\tMax curvature");
+        System.out.println(measureSymmetryPlane(face1, true) + "\tmsec:\tSymmetry plane");
         
         System.out.println(measureKdTreeCreation(face1) + "\tmsec:\tKd-tree creation of first face");
         System.out.println(measureKdTreeCreation(face2) + "\tmsec:\tKd-tree creation of second face");
@@ -46,17 +51,17 @@ public class EfficiencyTests {
         boolean relativeDist = true;
         boolean printDetails = false;
         
-        //System.out.println();
-        //System.out.println("Hausdorff distance sequentially:");
-        //testAndPrint(face1, new HausdorffDistance(face2.getMeshModel(), Strategy.POINT_TO_POINT_ABSOLUTE, false, false, false), printDetails);
-        //testAndPrint(face1, new HausdorffDistance(face2.getMeshModel(), Strategy.POINT_TO_POINT, relativeDist, false, false), printDetails);
-        //testAndPrint(face1, new HausdorffDistance(face2.getMeshModel(), Strategy.POINT_TO_TRIANGLE_APPROXIMATE, relativeDist, false, false), printDetails);
+        System.out.println();
+        System.out.println("Hausdorff distance sequentially:");
+        testAndPrint(face1, new HausdorffDistance(face2.getMeshModel(), Strategy.POINT_TO_POINT_ABSOLUTE, false, false, false), printDetails);
+        testAndPrint(face1, new HausdorffDistance(face2.getMeshModel(), Strategy.POINT_TO_POINT, relativeDist, false, false), printDetails);
+        testAndPrint(face1, new HausdorffDistance(face2.getMeshModel(), Strategy.POINT_TO_TRIANGLE_APPROXIMATE, relativeDist, false, false), printDetails);
 
-        //System.out.println();
-        //System.out.println("Hausdorff distance consurrently:");
-        //testAndPrint(face1, new HausdorffDistance(face2.getMeshModel(), Strategy.POINT_TO_POINT_ABSOLUTE, false, false, true), printDetails);
-        //testAndPrint(face1, new HausdorffDistance(face2.getMeshModel(), Strategy.POINT_TO_POINT, relativeDist, false, true), printDetails);
-        //testAndPrint(face1, new HausdorffDistance(face2.getMeshModel(), Strategy.POINT_TO_TRIANGLE_APPROXIMATE, relativeDist, false, true), printDetails);
+        System.out.println();
+        System.out.println("Hausdorff distance consurrently:");
+        testAndPrint(face1, new HausdorffDistance(face2.getMeshModel(), Strategy.POINT_TO_POINT_ABSOLUTE, false, false, true), printDetails);
+        testAndPrint(face1, new HausdorffDistance(face2.getMeshModel(), Strategy.POINT_TO_POINT, relativeDist, false, true), printDetails);
+        testAndPrint(face1, new HausdorffDistance(face2.getMeshModel(), Strategy.POINT_TO_TRIANGLE_APPROXIMATE, relativeDist, false, true), printDetails);
     }
     
     protected static void testAndPrint(HumanFace face, HausdorffDistance vis, boolean printDetails) {
@@ -98,4 +103,16 @@ public class EfficiencyTests {
         
         return retTime;
     }
+    
+    private static long measureCurvature(HumanFace face, GaussianCurvatureVisitor vis, boolean printDetails) {
+        long startTime = System.currentTimeMillis();
+        face.getMeshModel().compute(vis);
+        long retTime =  System.currentTimeMillis() - startTime;
+        
+        if (printDetails) {
+            
+        }
+        
+        return retTime;
+    }
 }
diff --git a/Comparison/src/main/java/cz/fidentis/analyst/symmetry/ApproxSymmetryPlane.java b/Comparison/src/main/java/cz/fidentis/analyst/symmetry/ApproxSymmetryPlane.java
index 7bc26d63..c645abcf 100644
--- a/Comparison/src/main/java/cz/fidentis/analyst/symmetry/ApproxSymmetryPlane.java
+++ b/Comparison/src/main/java/cz/fidentis/analyst/symmetry/ApproxSymmetryPlane.java
@@ -1,7 +1,10 @@
 package cz.fidentis.analyst.symmetry;
 
+import cz.fidentis.analyst.mesh.core.MeshFacet;
+import cz.fidentis.analyst.mesh.core.MeshPoint;
+import javax.vecmath.Vector3d;
+
 /**
- * 
  * Symmetry plane with votes used for the decision about symmetry estimate of 3D models.
  *
  * @author Natalia Bebjakova
@@ -11,6 +14,29 @@ public class ApproxSymmetryPlane extends Plane implements Comparable<ApproxSymme
     
     private int votes;
 
+    /**
+     * Constructor.
+     * 
+     * @param facet Mesh facet
+     * @param sigPoints Mesh vertices with the most significant curvature
+     * @param config Symmetry plane configuration
+     * @param i index of the first significant point used for the plane construction
+     * @param j index of the second significant point used for the plane construction
+     * @param maxDistance Distance limit
+     * @throws UnsupportedOperationException if the symmetry plane cannot be created
+     */
+    public ApproxSymmetryPlane(
+            MeshFacet facet, SignificantPoints sigPoints, Config config, int i, int j, double maxDistance
+    ) throws UnsupportedOperationException {
+        
+        if (i == j) {
+            throw new UnsupportedOperationException();
+        }
+        
+        setNormAndDist(facet, sigPoints, config, i, j);
+        computeVotes(facet, sigPoints, config, maxDistance);
+    }
+
     /**
      * returns number of votes that were given to plane while computing the symmetry 
      * 
@@ -19,16 +45,121 @@ public class ApproxSymmetryPlane extends Plane implements Comparable<ApproxSymme
     public int getVotes() {
         return votes;
     }
+    
+    private void setNormAndDist(MeshFacet facet, SignificantPoints sigPoints, Config config, int i, int j) {
+        double sigCurvI = sigPoints.getCurvature(i);
+        double sigCurvJ = sigPoints.getCurvature(j);
+        int sigPointI = sigPoints.getVertexIndex(i);
+        int sigPointJ = sigPoints.getVertexIndex(j);
+        
+        double minRatio = config.getMinCurvRatio();
+        double maxRatio = 1.0 / minRatio;
+        double ratioIJ = sigCurvI / sigCurvJ; 
+        if (ratioIJ < minRatio || ratioIJ > maxRatio) {
+            throw new UnsupportedOperationException();
+        }
+                
+        MeshPoint meshPointI = facet.getVertex(sigPointI);
+        MeshPoint meshPointJ = facet.getVertex(sigPointJ);
+                
+        Vector3d p1 = new Vector3d(meshPointI.getPosition());
+        Vector3d p2 = new Vector3d(meshPointJ.getPosition());
+                        
+        Vector3d avrg = new Vector3d(p1);
+        avrg.add(p2);
+        avrg.scale(0.5);
+                        
+        Vector3d normal = new Vector3d(p1);
+        normal.sub(p2);
+        normal.normalize(); 
+                        
+        double d = -(normal.x * avrg.x) - (normal.y * avrg.y) - (normal.z * avrg.z);
+                       
+        Vector3d ni = new Vector3d(meshPointI.getNormal());
+        Vector3d nj = new Vector3d(meshPointI.getNormal());
+        ni.normalize();
+        nj.normalize();
+
+        Vector3d normVec = new Vector3d(ni);
+        normVec.sub(nj);
+        normVec.normalize();
+        double normCos = normVec.dot(normal);
+
+        if (Math.abs(normCos) < config.getMinNormAngleCos()) {
+            throw new UnsupportedOperationException();
+        }
+
+        setNormal(normal);
+        setDistance(d);
+    }
 
     /**
-     * Constructor.
-     * @param plane Original plane without votes 
-     * @param votes number of votes given to the plane
-     * @throws IllegalArgumentExpcption if the @code{plane} argument is null
+     * Computes votes for given plane 
+     *
+     * @param facet Mesh facet
+     * @param sigPoints Mesh vertices with the most significant curvature
+     * @param config Symmetry plane configuration
+     * @param maxDist Distance limit
      */
-    public ApproxSymmetryPlane(Plane plane, int votes) {
-        super(plane);
-        this.votes = votes;
+    private void computeVotes(MeshFacet facet, SignificantPoints sigPoints, Config config, double maxDist) {
+        
+        normalize();
+        
+        Vector3d normal = getNormal();
+        double d = getDistance();
+        double maxCurvRatio = 1.0 / config.getMinCurvRatio();
+        
+        if (!facet.hasVertexNormals()) {
+            facet.calculateVertexNormals();
+        }
+
+        for (int i = 0; i < sigPoints.size(); i++) {
+            for (int j = 0; j < sigPoints.size(); j++) {
+                if (i == j) {
+                    continue;
+                }
+                
+                double ratioIJ = sigPoints.getCurvature(i) / sigPoints.getCurvature(j);
+                if (ratioIJ < config.getMinCurvRatio() || ratioIJ > maxCurvRatio) {
+                    continue;
+                }
+                
+                int sigPointI = sigPoints.getVertexIndex(i);
+                int sigPointJ = sigPoints.getVertexIndex(j);
+                    
+                Vector3d vec = new Vector3d(facet.getVertices().get(sigPointI).getPosition());
+                vec.sub(facet.getVertices().get(sigPointJ).getPosition());
+                vec.normalize();
+                double cos = vec.dot(normal);
+                    
+                Vector3d ni = new Vector3d(facet.getVertex(sigPointI).getNormal());
+                Vector3d nj = new Vector3d(facet.getVertex(sigPointJ).getNormal());
+                ni.normalize();
+                nj.normalize();
+                    
+                Vector3d normVec = ni;
+                normVec.sub(nj);
+                normVec.normalize();
+                double normCos = normVec.dot(normal);
+ 
+                if (Math.abs(cos) < config.getMinAngleCos() || Math.abs(normCos) < config.getMinNormAngleCos()) {
+                    continue;
+                }
+                
+                Vector3d avrg = new Vector3d(facet.getVertices().get(sigPointI).getPosition());
+                Vector3d aux = new Vector3d(facet.getVertices().get(sigPointJ).getPosition());
+                avrg.add(aux);
+                avrg.scale(0.5);
+                        
+                double dist = normal.x * avrg.x + normal.y * avrg.y + normal.z * avrg.z + d;
+                        
+                dist = Math.abs(dist);
+                        
+                if (dist <= maxDist) {                 
+                    votes++;  
+                }   
+            }
+        }
     }
     
     /**
@@ -41,4 +172,9 @@ public class ApproxSymmetryPlane extends Plane implements Comparable<ApproxSymme
     public int compareTo(ApproxSymmetryPlane other) {
         return Integer.compare(votes, other.votes);
     }
+    
+    @Override
+    public String toString() {
+        return this.getNormal() + " " + getDistance() + " " + votes;
+    }
 }
diff --git a/Comparison/src/main/java/cz/fidentis/analyst/symmetry/Plane.java b/Comparison/src/main/java/cz/fidentis/analyst/symmetry/Plane.java
index b2827693..f27a4d69 100644
--- a/Comparison/src/main/java/cz/fidentis/analyst/symmetry/Plane.java
+++ b/Comparison/src/main/java/cz/fidentis/analyst/symmetry/Plane.java
@@ -1,5 +1,6 @@
 package cz.fidentis.analyst.symmetry;
 
+import java.util.List;
 import javax.vecmath.Vector3d;
 
 /**
@@ -44,6 +45,39 @@ public class Plane {
         this(plane.getNormal(), plane.getDistance());
     }
     
+    /**
+     * Creates average plane from existing planes.
+     * 
+     * @param planes Source planes
+     * @throws {@code IllegalArgumentException} the {@code planes} list is {@code null} or empty
+     */
+    public Plane(List<Plane> planes) {
+        if (planes == null || planes.isEmpty()) {
+            throw new IllegalArgumentException("planes");
+        }
+        
+        Vector3d n = new Vector3d();
+        double d = 0;
+        Vector3d refDir = planes.get(0).getNormal();
+        for (int i = 0; i < planes.size(); i++) {
+            Vector3d normDir = planes.get(i).getNormal();
+            if (normDir.dot(refDir) < 0) {
+                n.sub(normDir);
+                d -= planes.get(i).getDistance();
+            } else {
+                n.add(normDir);
+                d += planes.get(i).getDistance();
+            }
+        }
+        
+        setNormal(n);
+        setDistance(d);
+        normalize();
+    }
+    
+    protected Plane() {
+    }
+    
     /**
      * Normalize the plane
      */
@@ -74,4 +108,12 @@ public class Plane {
     public double getDistance() {
         return distance;
     }
+    
+    protected void setNormal(Vector3d normal) {
+        this.normal = normal;
+    }
+    
+    protected void setDistance(double dist) {
+        this.distance = dist;
+    }
 }
diff --git a/Comparison/src/main/java/cz/fidentis/analyst/symmetry/SignificantPoints.java b/Comparison/src/main/java/cz/fidentis/analyst/symmetry/SignificantPoints.java
new file mode 100644
index 00000000..7c3e3eb7
--- /dev/null
+++ b/Comparison/src/main/java/cz/fidentis/analyst/symmetry/SignificantPoints.java
@@ -0,0 +1,96 @@
+package cz.fidentis.analyst.symmetry;
+
+import java.util.ArrayList;
+import java.util.Collections;
+import java.util.List;
+
+/**
+ * Finds and stores the top X points with the most significant curvarure.
+ * 
+ * @author Radek Oslejsek
+ * @author Natalia Bebjakova
+ */
+public class SignificantPoints {
+
+    private List<VertCurvature> significantPoints;
+    
+    /**
+     * Constructor.
+     * @param curvatures List of points' curvatures
+     * @param max Maximal number of points
+     */
+    public SignificantPoints(List<Double> curvatures, int max) {
+        significantPoints = new ArrayList<>(curvatures.size());
+        for (int i = 0; i < curvatures.size(); i++) { 
+            significantPoints.add(new VertCurvature(i, curvatures.get(i)));
+        }
+        
+        Collections.sort(significantPoints);
+        if (max < significantPoints.size()) {
+            significantPoints = new ArrayList<>(significantPoints.subList(0, max));
+        }
+    }
+    
+    /**
+     * Returns index of the i-th most significant mesh vertex.
+     * 
+     * @param i The "i" (index of the ordering)
+     * @return index of the i-th most significant mesh vertex.
+     * @throws IllegalArgumentException if the {@code i} paramter is out of rage
+     */
+    public int getVertexIndex(int i) {
+        if (i >= 0 && i < significantPoints.size()) {
+            return significantPoints.get(i).vertIndex;
+        } else {
+            throw new IllegalArgumentException("i");
+        }
+    }
+    
+    /**
+     * Returns curvature of the i-th most significant mesh vertex.
+     * 
+     * @param i The "i" (index of the ordering)
+     * @return curvature of the i-th most significant mesh vertex.
+     * @throws IllegalArgumentException if the {@code i} paramter is out of rage
+     */
+    public double getCurvature(int i) {
+        if (i >= 0 && i < significantPoints.size()) {
+            return significantPoints.get(i).curvature;
+        } else {
+            throw new IllegalArgumentException("i");
+        }
+    }
+    
+    /**
+     * Returns numnber of the most significant points.
+     * @return numnber of the most significant points.
+     */
+    public int size() {
+        return significantPoints.size();
+    }
+    
+    /**
+     * Helper class for sorting points with respect to their curvature.
+     * @author Radek Oslejsek
+     */
+    private class VertCurvature implements Comparable<VertCurvature> {
+        
+        private final int vertIndex;
+        private final double curvature;
+        
+        VertCurvature(int vertIndex, double curvature) {
+            this.vertIndex = vertIndex;
+            this.curvature = curvature;
+        }
+
+        @Override
+        public int compareTo(VertCurvature arg) {
+            return -Double.compare(this.curvature, arg.curvature);
+        }
+        
+        public String toString() {
+            return ""+curvature;
+        }
+    }
+    
+}
diff --git a/Comparison/src/main/java/cz/fidentis/analyst/symmetry/SymmetryEstimator.java b/Comparison/src/main/java/cz/fidentis/analyst/symmetry/SymmetryEstimator.java
index d4025e8e..c11b9552 100644
--- a/Comparison/src/main/java/cz/fidentis/analyst/symmetry/SymmetryEstimator.java
+++ b/Comparison/src/main/java/cz/fidentis/analyst/symmetry/SymmetryEstimator.java
@@ -2,7 +2,6 @@ package cz.fidentis.analyst.symmetry;
 
 import cz.fidentis.analyst.mesh.core.CornerTableRow;
 import cz.fidentis.analyst.mesh.core.MeshFacet;
-import cz.fidentis.analyst.mesh.core.MeshPoint;
 import cz.fidentis.analyst.visitors.mesh.BoundingBox;
 import cz.fidentis.analyst.mesh.core.MeshFacetImpl;
 import cz.fidentis.analyst.mesh.core.MeshPointImpl;
@@ -13,13 +12,15 @@ import cz.fidentis.analyst.visitors.mesh.GaussianCurvatureVisitor;
 import cz.fidentis.analyst.visitors.mesh.MaxCurvatureVisitor;
 import cz.fidentis.analyst.visitors.mesh.TriangleListVisitor;
 import java.util.ArrayList;
-import java.util.Collections;
 import java.util.List;
-import javax.swing.ImageIcon;
-import javax.swing.JOptionPane;
+import java.util.concurrent.Callable;
+import java.util.concurrent.ExecutionException;
+import java.util.concurrent.ExecutorService;
+import java.util.concurrent.Executors;
+import java.util.concurrent.Future;
+import java.util.logging.Level;
+import java.util.logging.Logger;
 import javax.swing.JPanel;
-import javax.swing.ProgressMonitor;
-import javax.swing.UIManager;
 import javax.vecmath.Vector3d;
 
 /**
@@ -29,6 +30,7 @@ import javax.vecmath.Vector3d;
  * On many different 3D models, it exists other values of config that will have better impact on result. 
  * 
  * @author Natalia Bebjakova
+ * @author Radek Oslejsek
  */
 public class SymmetryEstimator {
     
@@ -38,8 +40,6 @@ public class SymmetryEstimator {
     
     private BoundingBox boundingBox; // Represent min-max box. It is automatically maintained by given point array.
     
-    private List<MeshTriangle> triangles; // Helping array of triangles computed from corner table 
-    
     private final Config config;
     
     /**
@@ -60,7 +60,10 @@ public class SymmetryEstimator {
         
         TriangleListVisitor vis = new TriangleListVisitor(false);
         f.accept(vis);
-        triangles = vis.getTriangles();
+        
+        BoundingBoxVisitor visitor = new BoundingBoxVisitor(false);
+        facet.accept(visitor);
+        boundingBox = visitor.getBoundingBox();
     }
     
     /**
@@ -90,20 +93,6 @@ public class SymmetryEstimator {
         return facet;
     }
 
-    /**
-     * If bounding box is not created yet, it's created now. 
-     * 
-     * @return Represent min-max box of the boundries of the model. 
-     */
-    public BoundingBox getBoundingBox() {
-        if (boundingBox == null) {
-            BoundingBoxVisitor visitor = new BoundingBoxVisitor(false);
-            facet.accept(visitor);
-            boundingBox = visitor.getBoundingBox();
-        }
-        return boundingBox; 
-    }
-    
     /**
      * Computes the approximate plane of symmetry.
      * TO DO: ZRUSIT VSTUPNI PARAMETR!!!
@@ -112,113 +101,75 @@ public class SymmetryEstimator {
      * @return approximate plane of symmtetry
      */
     public Plane getApproxSymmetryPlane(JPanel panel) {
-        
         if (!facet.hasVertexNormals()) {
             facet.calculateVertexNormals();
         }
         
-        List<ApproxSymmetryPlane> planes = new ArrayList<>();
+        List<Double> curvatures = calculateCurvatures(facet);
+        final SignificantPoints sigPoints = new SignificantPoints(curvatures, config.getSignificantPointCount());
         
-        List<Double> curvatures = initCurvatures(facet);
-        List<Double> sortedCurvatures = new ArrayList<>(curvatures);
-        Collections.sort(sortedCurvatures);
-           
-        if (config.getSignificantPointCount() > facet.getNumberOfVertices()) {
-            config.setSignificantPointCount((facet.getNumberOfVertices()) - 1);
-        }
-        double bottomCurvature = sortedCurvatures.get(sortedCurvatures.size() - 1 - config.getSignificantPointCount());   
+        // Initiate structure for concurrent computation:
+        ExecutorService executor = Executors.newFixedThreadPool(Runtime.getRuntime().availableProcessors());
+        final List<Future<ApproxSymmetryPlane>> results = new ArrayList<>(sigPoints.size()*sigPoints.size());
         
-        List<Integer> significantPoints = new ArrayList<>();
-        List<Double> significantCurvatures = new ArrayList<>();
-        for (int i = 0; i < facet.getNumberOfVertices(); i++) {
-            if (curvatures.get(i) >= bottomCurvature) {
-                significantCurvatures.add(curvatures.get(i));
-                significantPoints.add(i);
+        // Compute candidate planes concurrently:
+        final double maxDistance = boundingBox.getMaxDiag() * config.getMaxRelDistance();
+        for (int i = 0; i < sigPoints.size(); i++) {
+            for (int j = 0; j < sigPoints.size(); j++) {
+                int finalI = i;
+                int finalJ = j;
+                
+                Callable<ApproxSymmetryPlane> worker = new  Callable<ApproxSymmetryPlane>() {
+                    @Override
+                    public ApproxSymmetryPlane call() throws Exception {
+                        try {
+                            return new ApproxSymmetryPlane(facet, sigPoints, config, finalI, finalJ, maxDistance);
+                        } catch(UnsupportedOperationException ex) {
+                            return null;
+                        }
+                    }
+                };
+                
+                results.add(executor.submit(worker));
             }
         }
         
-        Plane plane = new Plane(0, 0, 0, 0);
-        int lastVotes = 0;
+        // Wait until all symmetry planes are computed:
+        executor.shutdown();
+        while (!executor.isTerminated()){
+        }
         
-        for (int i = 0; i < significantPoints.size(); i++) {
-            for (int j = 0; j < significantPoints.size(); j++) {
-                if (i == j) {
+        // Get the best-fitting planes:
+        List<Plane> planes = new ArrayList<>();
+        try {
+            int maxVotes = 0;
+            for (Future<ApproxSymmetryPlane> res: results) {
+                ApproxSymmetryPlane newPlane = res.get();
+                if (newPlane == null) {
                     continue;
                 }
-                
-                double minRatio = config.getMinCurvRatio();
-                double maxRatio = 1.0 / minRatio;
-                if (significantCurvatures.get(i) / significantCurvatures.get(j) < minRatio ||
-                        significantCurvatures.get(i) / significantCurvatures.get(j) > maxRatio) {
-                    continue;
-                }
-                
-                MeshPoint meshPointI = facet.getVertex(significantPoints.get(i));
-                MeshPoint meshPointJ = facet.getVertex(significantPoints.get(j));
-                
-                Vector3d p1 = new Vector3d(meshPointI.getPosition());
-                Vector3d p2 = new Vector3d(meshPointJ.getPosition());
-                        
-                Vector3d avrg = new Vector3d(p1);
-                avrg.add(p2);
-                avrg.scale(0.5);
-                        
-                Vector3d normal = new Vector3d(p1);
-                normal.sub(p2);
-                normal.normalize(); 
-                        
-                double d = -(normal.x * avrg.x) - (normal.y * avrg.y) - (normal.z * avrg.z);
-                       
-                Vector3d ni = new Vector3d(meshPointI.getNormal());
-                Vector3d nj = new Vector3d(meshPointI.getNormal());
-                ni.normalize();
-                nj.normalize();
-
-                Vector3d normVec = new Vector3d(ni);
-                normVec.sub(nj);
-                normVec.normalize();
-                double normCos = normVec.dot(normal);
-
-                if (Math.abs(normCos) < config.getMinNormAngleCos()) {
-                    continue;
-                }
-                
-                Plane newPlane = new Plane(normal.x, normal.y, normal.z, d);
-                int currentVotes = getVotes(
-                        newPlane, 
-                        significantCurvatures, 
-                        significantPoints, 
-                        config.getMinCurvRatio(), 
-                        config.getMinAngleCos(), 
-                        config.getMinNormAngleCos(), 
-                        getBoundingBox().getMaxDiag() * config.getMaxRelDistance()
-                );
-
-                planes.add(new ApproxSymmetryPlane(newPlane, currentVotes));
-                
-                if (currentVotes > lastVotes) {
-                    lastVotes = currentVotes;
-                    plane = newPlane;
+                if (newPlane.getVotes() > maxVotes) {
+                    maxVotes = newPlane.getVotes();
+                    planes.clear();
+                    planes.add(newPlane);
+                } else if (newPlane.getVotes() == maxVotes) {
+                    planes.add(newPlane);
                 }
             }
+        } catch (final InterruptedException | ExecutionException ex) {
+            Logger.getLogger(SymmetryEstimator.class.getName()).log(Level.SEVERE, null, ex);
         }
         
-        Collections.sort(planes);
-        List<ApproxSymmetryPlane> finalPlanes = new ArrayList<>();
-        for (int i = 0; i < planes.size(); i++) {
-            if (planes.get(i).getVotes() == lastVotes){
-                finalPlanes.add(planes.get(i));
-            }
-        }
-        Plane finalPlane = computeNewPlane(finalPlanes);
-        if (config.isAveraging()){
-            plane = finalPlane;
+        // Compute the average plane, if required...
+        if (config.isAveraging() && !planes.isEmpty()) {
+            return new Plane(planes);
         }
-
-        return plane;
+        
+        // ... or return the random best-fitting plane:
+        return planes.isEmpty() ? null : planes.get(0);
     }
     
-    private List<Double> initCurvatures(MeshFacet facet) {
+    private List<Double> calculateCurvatures(MeshFacet facet) {
         CurvatureVisitor vis;
         if (facet.getNumberOfVertices() < 2500) {
             // searching for Maximum curvature in vertex using Gaussian curvature 
@@ -227,7 +178,7 @@ public class SymmetryEstimator {
             vis = new MaxCurvatureVisitor(false);
         } else {
             //curvatures.add(getGaussianCurvature(i));
-            vis = new GaussianCurvatureVisitor(false);            
+            vis = new GaussianCurvatureVisitor(false);
         }
         facet.accept(vis);
         List<Double> curvatures = new ArrayList<>(vis.getCurvatures().get(facet));
@@ -237,49 +188,6 @@ public class SymmetryEstimator {
             }
         }
         return curvatures;
-        
-        /*
-        for (int i = 0; i < facet.getNumberOfVertices(); i++) {
-            if (facet.getNumberOfVertices() < 2500) {
-                // searching for Maximum curvature in vertex using Gaussian curvature 
-                // and Mean curvature leads to longer calculation but better results for real face models
-                //curvatures.add(getMaxCurvature(i));
-                vis = new MaxCurvatureVisitor(false);
-            } else {
-                //curvatures.add(getGaussianCurvature(i));
-                vis = new GaussianCurvatureVisitor(false);
-            }
-            if (Double.isNaN(curvatures.get(i))){
-                curvatures.set(i, Double.MIN_VALUE);
-            }
-        }
-        return curvatures;
-        */
-    }
-    
-    private Plane computeNewPlane(List<ApproxSymmetryPlane> finalPlanes) {
-        if (finalPlanes.isEmpty()) {
-            return null;
-        }
-        double newA = 0, newB = 0, newC = 0, newD = 0;
-        Vector3d refDir = finalPlanes.get(0).getNormal();
-        for (int i = 0; i < finalPlanes.size(); i++) {
-            Vector3d normDir = finalPlanes.get(i).getNormal();
-            if (normDir.dot(refDir) < 0) {
-                newA -= normDir.x;
-                newB -= normDir.y;
-                newC -= normDir.z;
-                newD -= finalPlanes.get(i).getDistance();
-            } else {
-                newA += normDir.x;
-                newB += normDir.y;
-                newC += normDir.z;
-                newD += finalPlanes.get(i).getDistance();
-            }
-        }
-        Plane finalPlane = new Plane(newA, newB, newC, newD);
-        finalPlane.normalize();
-        return finalPlane;
     }
     
     /**
@@ -289,7 +197,7 @@ public class SymmetryEstimator {
      */
     public SymmetryEstimator mergeWithPlane(Plane plane) {
         Vector3d normal = plane.getNormal();
-        Vector3d midPoint = getBoundingBox().getMidPoint().getPosition();
+        Vector3d midPoint = boundingBox.getMidPoint().getPosition();
 
         double alpha = -((normal.x * midPoint.x) + 
                 (normal.y * midPoint.y) + (normal.z * midPoint.z) +
@@ -317,7 +225,7 @@ public class SymmetryEstimator {
         b.normalize();
         
         SymmetryEstimator planeMesh = new SymmetryEstimator(midPointOnPlane, a, b,
-                (getBoundingBox().getMaxPoint().subtractPosition(getBoundingBox().getMinPoint())).getPosition().x);
+                (boundingBox.getMaxPoint().subtractPosition(boundingBox.getMinPoint())).getPosition().x);
        
         return mergeMeshWith(planeMesh);
     }
@@ -381,76 +289,4 @@ public class SymmetryEstimator {
         return facet;
     }
 
-    /**
-     * Computes votes for given plane 
-     * 
-     * @param plane Plane for which votes are computed
-     * @param curvatures significant curvatures chosen for computing
-     * @param points significant points chosen for computing
-     * @param minCurvRatio optional parameter from configuration
-     * @param minAngleCos optional parameter from configuration
-     * @param minNormAngleCos optional parameter from configuration
-     * @param maxDist optional parameter from configuration
-     * @return total votes given to plane while computing the symmetry
-     */
-    private int getVotes(Plane plane,
-            List<Double> curvatures,
-            List<Integer> points,
-            double minCurvRatio,
-            double minAngleCos,
-            double minNormAngleCos,
-            double maxDist) {
-        
-        plane.normalize();
-        
-        Vector3d normal = plane.getNormal();
-        double d = plane.getDistance();
-        double maxCurvRatio = 1.0 / minCurvRatio;
-        int votes = 0;
-        //List<Vector3d> normals = calculateNormals();
-        
-        //if (!facet.hasVertexNormals()) {
-        //    facet.calculateVertexNormals();
-        //}
-
-        for (int i = 0; i < curvatures.size(); i++) {
-            for (int j = 0; j < curvatures.size(); j++) {
-                if (i != j && (curvatures.get(i) / curvatures.get(j) >= minCurvRatio
-                        && curvatures.get(i) / curvatures.get(j) <= maxCurvRatio)) {
-                    
-                    Vector3d vec = new Vector3d(facet.getVertices().get(points.get(i)).getPosition());
-                    vec.sub(facet.getVertices().get(points.get(j)).getPosition());
-                    vec.normalize();
-                    double cos = vec.dot(normal);
-                    
-                    Vector3d ni = new Vector3d(facet.getVertex(points.get(i)).getNormal());
-                    Vector3d nj = new Vector3d(facet.getVertex(points.get(j)).getNormal());
-                    ni.normalize();
-                    nj.normalize();
-                    
-                    Vector3d normVec = ni;
-                    normVec.sub(nj);
-                    normVec.normalize();
-                    double normCos = normVec.dot(normal);
- 
-                    if (Math.abs(cos) >= minAngleCos && Math.abs(normCos) >= minNormAngleCos) {
-                        Vector3d avrg = new Vector3d(facet.getVertices().get(points.get(i)).getPosition());
-                        Vector3d aux = new Vector3d(facet.getVertices().get(points.get(i)).getPosition());
-                        avrg.add(aux);
-                        avrg.scale(0.5);
-                        
-                        double dist = normal.x * avrg.x + normal.y * avrg.y + normal.z * avrg.z + d;
-                        
-                        dist = Math.abs(dist);
-                        
-                        if (dist <= maxDist) {                 
-                            votes++;  
-                        }   
-                    }
-                }
-            }
-        }
-        return votes;
-    }
-    
 }
diff --git a/Comparison/src/main/java/cz/fidentis/analyst/visitors/mesh/CurvatureVisitor.java b/Comparison/src/main/java/cz/fidentis/analyst/visitors/mesh/CurvatureVisitor.java
index a6f2c06d..94adbc7f 100644
--- a/Comparison/src/main/java/cz/fidentis/analyst/visitors/mesh/CurvatureVisitor.java
+++ b/Comparison/src/main/java/cz/fidentis/analyst/visitors/mesh/CurvatureVisitor.java
@@ -1,57 +1,66 @@
-/*
- * To change this license header, choose License Headers in Project Properties.
- * To change this template file, choose Tools | Templates
- * and open the template in the editor.
- */
 package cz.fidentis.analyst.visitors.mesh;
 
 import cz.fidentis.analyst.mesh.MeshVisitor;
 import cz.fidentis.analyst.mesh.core.MeshFacet;
-import cz.fidentis.analyst.mesh.core.MeshPoint;
 import cz.fidentis.analyst.mesh.core.MeshTriangle;
-import cz.fidentis.analyst.symmetry.SymmetryEstimator;
 import java.util.ArrayList;
 import java.util.Collections;
 import java.util.HashMap;
 import java.util.List;
 import java.util.Map;
-import javax.vecmath.Vector3d;
+import java.util.concurrent.ExecutorService;
+import java.util.concurrent.Executors;
+import java.util.concurrent.Future;
 
 /**
- *
- * @author oslejsek
+ * Abstract class for algorithms calculting curvatures of mesh vertices.
+ * 
+ * @author Natalia Bebjakova
+ * @author Radek Oslejsek
  */
 public abstract class CurvatureVisitor extends MeshVisitor {
     
-    private Map<MeshFacet, List<Double>> curvatures = new HashMap<>();
+    private final Map<MeshFacet, List<Double>> curvatures = new HashMap<>();
     
     /**
      * Constructor.
      * 
-     * @param asynchronous If {@code true}, then asynchronous visitor is created.
+     * @param asynchronous If {@code true}, then asynchronous visitor is created, 
+     * i.e., has to be invoked by {@code Executor}
      */
     public CurvatureVisitor(boolean asynchronous) {
         super(asynchronous);
     }
- 
+    
+    /**
+     * Returns curvatures for all inspected mesh facets. The order corresponds to
+     * the order of vertices, i.e., i-th value represents the curvature of i-th mesh vertex.
+     * 
+     * @return curvatures of inspected mesh facets.
+     */
+    public Map<MeshFacet, List<Double>> getCurvatures() {
+        return Collections.unmodifiableMap(curvatures);
+    }
+    
     @Override
-    protected synchronized void visitMeshFacet(MeshFacet facet) {
-        if (curvatures.containsKey(facet)) {
-            return; // already visited facet
+    protected void visitMeshFacet(final MeshFacet facet) {
+        synchronized (this) {
+            if (curvatures.containsKey(facet)) {
+                return; // already visited facet
+            }
+            curvatures.put(facet, new ArrayList<>());
+            facet.calculateVoronoiPoints();
         }
         
-        curvatures.put(facet, new ArrayList<>());
+        ExecutorService executor = Executors.newFixedThreadPool(Runtime.getRuntime().availableProcessors());
+        final List<Future<Double>> results = new ArrayList<>(facet.getNumberOfVertices());
         
-        List<MeshTriangle> triangles = facet.getTriangles();
+        final List<MeshTriangle> triangles = facet.getTriangles();
         for (int i = 0; i < facet.getNumberOfVertices(); i++) {
-           curvatures.get(facet).add(calculateCurvature(facet, triangles, i));
+            curvatures.get(facet).add(calculateCurvature(facet, triangles, i));
         }
     }
     
-    public Map<MeshFacet, List<Double>> getCurvatures() {
-        return Collections.unmodifiableMap(curvatures);
-    }
-    
     protected abstract double calculateCurvature(MeshFacet facet, List<MeshTriangle> triangles, int centerIndex);
     
 }
diff --git a/Comparison/src/main/java/cz/fidentis/analyst/visitors/mesh/GaussianCurvatureVisitor.java b/Comparison/src/main/java/cz/fidentis/analyst/visitors/mesh/GaussianCurvatureVisitor.java
index f6dd5051..64b2280b 100644
--- a/Comparison/src/main/java/cz/fidentis/analyst/visitors/mesh/GaussianCurvatureVisitor.java
+++ b/Comparison/src/main/java/cz/fidentis/analyst/visitors/mesh/GaussianCurvatureVisitor.java
@@ -1,8 +1,3 @@
-/*
- * To change this license header, choose License Headers in Project Properties.
- * To change this template file, choose Tools | Templates
- * and open the template in the editor.
- */
 package cz.fidentis.analyst.visitors.mesh;
 
 import cz.fidentis.analyst.mesh.core.MeshFacet;
@@ -11,15 +6,18 @@ import java.util.List;
 import javax.vecmath.Vector3d;
 
 /**
- *
- * @author oslejsek
+ * Calculates Gaussian curvatures of mesh vertices.
+ * 
+ * @author Natalia Bebjakova
+ * @author Radek Oslejsek
  */
 public class GaussianCurvatureVisitor extends CurvatureVisitor {
     
     /**
      * Constructor.
      * 
-     * @param asynchronous If {@code true}, then asynchronous visitor is created.
+     * @param asynchronous If {@code true}, then asynchronous visitor is created, 
+     * i.e., has to be invoked by {@code Executor}
      */
     public GaussianCurvatureVisitor(boolean asynchronous) {
         super(asynchronous);
diff --git a/Comparison/src/main/java/cz/fidentis/analyst/visitors/mesh/MaxCurvatureVisitor.java b/Comparison/src/main/java/cz/fidentis/analyst/visitors/mesh/MaxCurvatureVisitor.java
index d5018e08..f63e40c7 100644
--- a/Comparison/src/main/java/cz/fidentis/analyst/visitors/mesh/MaxCurvatureVisitor.java
+++ b/Comparison/src/main/java/cz/fidentis/analyst/visitors/mesh/MaxCurvatureVisitor.java
@@ -7,15 +7,20 @@ import java.util.List;
 import javax.vecmath.Vector3d;
 
 /**
- *
- * @author oslejsek
+ * Calculates curvatures of mesh vertices more precisely using Laplacian operator.
+ * This algorithm is roughly 5-times slower than the pure Gaussian curvature algorithm.
+ * 
+ * @author Natalia Bebjakova
+ * @author Radek Oslejsek
  */
-public class MaxCurvatureVisitor extends GaussianCurvatureVisitor{
+public class MaxCurvatureVisitor extends GaussianCurvatureVisitor {
     
     /**
      * Constructor.
      * 
-     * @param asynchronous If {@code true}, then asynchronous visitor is created.
+     * @param asynchronous If {@code true}, then asynchronous visitor is created, 
+     * i.e., has to be invoked by {@code Executor}
+     * @param parallel If {@code true}, then the algorithm computation runs concurrently utilizing all CPU cores
      */
     public MaxCurvatureVisitor(boolean asynchronous) {
         super(asynchronous);
diff --git a/MeshModel/src/main/java/cz/fidentis/analyst/mesh/core/MeshFacetImpl.java b/MeshModel/src/main/java/cz/fidentis/analyst/mesh/core/MeshFacetImpl.java
index 5b813c82..32fadb52 100644
--- a/MeshModel/src/main/java/cz/fidentis/analyst/mesh/core/MeshFacetImpl.java
+++ b/MeshModel/src/main/java/cz/fidentis/analyst/mesh/core/MeshFacetImpl.java
@@ -85,7 +85,7 @@ public class MeshFacetImpl implements MeshFacet {
      * @throws RuntimeException if there are duplicate meth points in the mesh facet
      */
     @Override
-    public void calculateVertexNormals() {
+    public synchronized void calculateVertexNormals() {
         Map<Vector3d, Vector3d> normalMap = new HashMap<>(); // key = mesh point, value = normal
         
         // init normals:
@@ -213,7 +213,7 @@ public class MeshFacetImpl implements MeshFacet {
     }
     
     @Override
-    public List<Vector3d> calculateVoronoiPoints() {
+    public synchronized List<Vector3d> calculateVoronoiPoints() {
         if (voronoiPoints == null) {
             voronoiPoints = new ArrayList<>(getNumTriangles());
             for (MeshTriangle tri: this) {
-- 
GitLab