GeodeticSubdivisionHelper

GeodeticSubdivisionHelper.java

package ua.in.stasikos.worldsim;

import ua.in.stasikos.worldsim.universe.celestials.RockyBody;
import ua.in.stasikos.worldsim.util.Vector3D;

import java.util.ArrayList;
import java.util.Arrays;
import java.util.List;

public class GeoSubdivisionHelper {

    private static final int NOT_EXIST = -1;
    private static final float EPSILON = 0.0001f;
    private ArrayList<Vector3D> vertices = new ArrayList<>();
    private int[] centers = new int[12];
    private ArrayList<Triangle> faces = new ArrayList<>();

    public void setup(int subdivisions, float radius) {
        createIcosphere(subdivisions, radius);
    }

    private void createIcosphere(int subdivisions, float radius) {
        float t = (float) ((1.0f + Math.sqrt(5.0f)) / 2.0f);

        // Vertices
        vertices.add(new Vector3D(-1, t, 0).normalizeTo(radius));
        vertices.add(new Vector3D(1, t, 0).normalizeTo(radius));
        vertices.add(new Vector3D(-1, -t, 0).normalizeTo(radius));
        vertices.add(new Vector3D(1, -t, 0).normalizeTo(radius));

        vertices.add(new Vector3D(0, -1, t).normalizeTo(radius));
        vertices.add(new Vector3D(0, 1, t).normalizeTo(radius));
        vertices.add(new Vector3D(0, -1, -t).normalizeTo(radius));
        vertices.add(new Vector3D(0, 1, -t).normalizeTo(radius));

        vertices.add(new Vector3D(t, 0, -1).normalizeTo(radius));
        vertices.add(new Vector3D(t, 0, 1).normalizeTo(radius));
        vertices.add(new Vector3D(-t, 0, -1).normalizeTo(radius));
        vertices.add(new Vector3D(-t, 0, 1).normalizeTo(radius));

        // Copy vertices into separate list for centers, so they are safe from shuffling (not necessary)
        for (int i = 0; i < 12; i++) {
            centers[i] = i;
        }

        // Faces (icosahedron)
        faces.add(new Triangle(0, 11, 5));
        faces.add(new Triangle(0, 5, 1));
        faces.add(new Triangle(0, 1, 7));
        faces.add(new Triangle(0, 7, 10));
        faces.add(new Triangle(0, 10, 11));

        faces.add(new Triangle(1, 5, 9));
        faces.add(new Triangle(5, 11, 4));
        faces.add(new Triangle(11, 10, 2));
        faces.add(new Triangle(10, 7, 6));
        faces.add(new Triangle(7, 1, 8));

        faces.add(new Triangle(3, 9, 4));
        faces.add(new Triangle(3, 4, 2));
        faces.add(new Triangle(3, 2, 6));
        faces.add(new Triangle(3, 6, 8));
        faces.add(new Triangle(3, 8, 9));

        faces.add(new Triangle(4, 9, 5));
        faces.add(new Triangle(2, 4, 11));
        faces.add(new Triangle(6, 2, 10));
        faces.add(new Triangle(8, 6, 7));
        faces.add(new Triangle(9, 8, 1));

        // Triple Subdivide
        for (int i = 0; i < subdivisions; i++) {
            ArrayList<Triangle> newFaces = new ArrayList<>();
            for (Triangle face : faces) {
                int a = getOneThird(face.a, face.b);
                int b = getOneThird(face.a, face.c);
                int c = getMiddlePoint(face.b, a);
                int d = getOneThird(face.b, face.c);
                int e = getMiddlePoint(face.c, d);
                int f = getMiddlePoint(face.c, b);
                int g = getMiddlePoint(c, f);

                newFaces.add(new Triangle(face.a, a, b));
                newFaces.add(new Triangle(a, c, g));
                newFaces.add(new Triangle(a, b, g));
                newFaces.add(new Triangle(b, f, g));
                newFaces.add(new Triangle(c, face.b, d));
                newFaces.add(new Triangle(c, d, g));
                newFaces.add(new Triangle(d, g, e));
                newFaces.add(new Triangle(e, g, f));
                newFaces.add(new Triangle(f, face.c, e));
            }
            faces = newFaces;
        }
        // Find duplicate vertices in vertices list and replace them with one vertex in each triangle
        for (int i = 0; i < vertices.size(); i++) {
            Vector3D v1 = vertices.get(i);
            for (int j = i + 1; j < vertices.size(); j++) {
                Vector3D v2 = vertices.get(j);
                if (v1.equals(v2)) {
                    for (Triangle face : faces) {
                        if (face.a == j) {
                            face.a = i;
                        }
                        if (face.b == j) {
                            face.b = i;
                        }
                        if (face.c == j) {
                            face.c = i;
                        }
                    }
                }
            }
        }

        // Sanity checks, remove later
        for (Triangle face : faces) {
            if (face.a == face.b || face.a == face.c || face.b == face.c) {
                throw new RuntimeException("Face has duplicate vertices");
            }
        }

        List<Integer> duplicates = new ArrayList<>();
        for (int i = 0; i < vertices.size(); i++) {
            Vector3D v1 = vertices.get(i);
            for (int j = i + 1; j < vertices.size(); j++) {
                Vector3D v2 = vertices.get(j);
                if (v1.equals(v2)) {
                    duplicates.add(j);
                }
            }
        }

        for (int duplicate : duplicates) {
            for (Triangle face: faces) {
                if (face.a == duplicate) {
                    throw new RuntimeException("Face has duplicate vertices");
                }
                if (face.b == duplicate) {
                    throw new RuntimeException("Face has duplicate vertices");
                }
                if (face.c == duplicate) {
                    throw new RuntimeException("Face has duplicate vertices");
                }
            }
        }

        // Normalize all vertices to radius
        for (int i = 0; i < vertices.size(); i++) {
            vertices.set(i, vertices.get(i).normalizeTo(radius));
        }
    }

    private int getOneThird(int p1, int p2) {
        Vector3D v1 = vertices.get(p1);
        Vector3D v2 = vertices.get(p2);
        Vector3D third = new Vector3D(
                (2 * v1.getX() + v2.getX()) / 3,
                (2 * v1.getY() + v2.getY()) / 3,
                (2 * v1.getZ() + v2.getZ()) / 3
        );

        int index = vertices.size();
        vertices.add(third);
        return index;
    }

    private int getMiddlePoint(int p1, int p2) {
        Vector3D v1 = vertices.get(p1);
        Vector3D v2 = vertices.get(p2);
        Vector3D middle = new Vector3D(
                (v1.getX() + v2.getX()) / 2,
                (v1.getY() + v2.getY()) / 2,
                (v1.getZ() + v2.getZ()) / 2
        );

        int index = vertices.size();
        vertices.add(middle);
        return index;
    }


    public List<Polygon> getRegions() {
        List<Polygon> regions = new ArrayList<>();

        System.out.println("Triangles: " + faces.size());

        // Start with icosahedron vertices
        int[] adjacent;

        for (int center : centers) {
            Polygon pentagon = new Pentagon();
            pentagon.center = center;
            adjacent = new int[5];
            for (int i = 0; i < adjacent.length; i++) {
                adjacent[i] = NOT_EXIST;
            }

            // Look for triangles with vertices in center of pentagon and remove them from face list
            int remaining = 5;
            for (int i = 0; i < faces.size(); i++) {
                Triangle face = faces.get(i);
                if (center == face.a) {
                    tryAdd(adjacent, face.b);
                    tryAdd(adjacent, face.c);
                    faces.remove(i);
                    remaining--;
                    i--;
                } else if (center == face.b) {
                    tryAdd(adjacent, face.a);
                    tryAdd(adjacent, face.c);
                    faces.remove(i);
                    remaining--;
                    i--;
                } else if (center == face.c) {
                    tryAdd(adjacent, face.a);
                    tryAdd(adjacent, face.b);
                    faces.remove(i);
                    remaining--;
                    i--;
                }
                if (remaining == 0) {
                    break;
                }
            }

            for (int i = 0; i < adjacent.length; i++) {
                pentagon.polyVertices[i] = adjacent[i];
                float[] latLong = getLatLong(vertices.get(adjacent[i]), 6371);
                System.out.println("Vertex " + i + ": " + vertices.get(adjacent[i]) + " index: " + adjacent[i] + " lat: " + latLong[0] + " long: " + latLong[1]);
            }

            regions.add(pentagon);
            System.out.println("Pentagon: " + pentagon);
        }

        System.out.println("Triangles left: " + faces.size());

        // do same for hexagons, using the remaining triangles and known vertices of all pentagons
        while (faces.size() > 0) {
            for (int i = 0; i < faces.size(); i++) {
                Triangle face = faces.get(i);

                // Does it have two adjacent vertices with any existing region?
                for (int j = 0; j < regions.size(); j++) {
                    int[] adjacentVertices = regions.get(j).getPolyVertices();
                    int center = findNotAdjacent(face, adjacentVertices);
                    if (center != NOT_EXIST) {
                        Polygon hexagon = new Hexagon();
                        hexagon.center = center;
                        adjacent = new int[6];
                        for (int k = 0; k < adjacent.length; k++) {
                            adjacent[k] = NOT_EXIST;
                        }

                        // Look for triangles with vertices in center of hexagon and remove them from face list
                        int remaining = 6;
                        for (int l = 0; l < faces.size(); l++) {
                            Triangle t = faces.get(l);
                            if (center == t.a) {
                                tryAdd(adjacent, t.b);
                                tryAdd(adjacent, t.c);
                                faces.remove(l);
                                l--;
                                remaining--;
                            } else if (center == t.b) {
                                tryAdd(adjacent, t.a);
                                tryAdd(adjacent, t.c);
                                faces.remove(l);
                                l--;
                                remaining--;
                            } else if (center == t.c) {
                                tryAdd(adjacent, t.a);
                                tryAdd(adjacent, t.b);
                                faces.remove(l);
                                remaining--;
                                l--;
                            }
                            if (remaining == 0) {
                                break;
                            }
                        }

                        for (int k = 0; k < adjacent.length; k++) {
                            hexagon.polyVertices[k] = adjacent[k];
                            if (adjacent[k] != NOT_EXIST) {
                                System.out.println("Vertex " + k + ": " + vertices.get(adjacent[k]));
                            } else {
                                System.out.println("Vertex " + k + ": " + "null");
//                                throw new RuntimeException("Vertex is null, should not happen ever");
                            }
                        }
                        // These are supposed to have same index for both lists.
                        regions.add(hexagon);
                        regions.get(j).addNeighbor(hexagon);
                        hexagon.addNeighbor(regions.get(j));
                        System.out.println("Hexagon: " + hexagon);
                        break; // If we repeat, we will find more neighbors, but triangle search will fail as we took them off list
                    }
                }
            }
            System.out.println("Triangles left: " + faces.size());
//            faces.remove(0);
        }

        // find remaining neighbors for all regions
        for (int i = 0; i < regions.size(); i++) {
            Polygon region = regions.get(i);
            for (int j = 0; j < regions.size(); j++) {
                if (i != j) {
                    Polygon neighbor = regions.get(j);
                    int[] vertices = region.getPolyVertices();
                    int[] neighborVertices = neighbor.getPolyVertices();
                    for (int k = 0; k < vertices.length; k++) {
                        for (int l = 0; l < neighborVertices.length; l++) {
                            if (vertices[k] == neighborVertices[l]) {
                                region.addNeighbor(neighbor);
                                neighbor.addNeighbor(region);
                                break;
                            }
                        }
                    }
                }
            }
        }

        return regions;
    }

    private int findNotAdjacent(Triangle face, int[] adjacentVertices) {
        for (int i = 0; i < adjacentVertices.length; i++) {
            if (face.a == adjacentVertices[i]) {
                for (int j = 0; j < adjacentVertices.length; j++) {
                    if (face.b == adjacentVertices[j]) {
                        return face.c;
                    } else if (face.c == adjacentVertices[j]) {
                        return face.b;
                    }
                }
            } else if (face.b == adjacentVertices[i]) {
                for (int j = 0; j < adjacentVertices.length; j++) {
                    if (face.a == adjacentVertices[j]) {
                        return face.c;
                    } else if (face.c == adjacentVertices[j]) {
                        return face.a;
                    }
                }
            } else if (face.c == adjacentVertices[i]) {
                for (int j = 0; j < adjacentVertices.length; j++) {
                    if (face.a == adjacentVertices[j]) {
                        return face.b;
                    } else if (face.b == adjacentVertices[j]) {
                        return face.a;
                    }

                }
            }
        }
        return NOT_EXIST;
    }

    public Vector3D getVertex(int index) {
        return vertices.get(index);
    }

    private void tryAdd(int[] adjacent, int vertex) {
        for (int i = 0; i < adjacent.length; i++) {
            if (adjacent[i] == vertex) { // Already added
                return;
            }

            // Sanity check
            if (adjacent[i] != NOT_EXIST) {
                if (vertices.get(adjacent[i]).equals(vertices.get(vertex))) {
                    System.out.println("Duplicate vertices: " + vertices.get(adjacent[i]) + " and " + vertices.get(vertex));
                     throw new RuntimeException("Duplicate vertices");
                }
            }

            if (adjacent[i] == NOT_EXIST) {
                adjacent[i] = vertex;
                break;
            }
        }
    }

    public float[] getLatLong(Vector3D v, float radius) {
        // Convert x, y, z points on sphere in 3D space to latitude and longitude
        float longitude = (float) Math.toDegrees(Math.atan2(v.getY(), v.getX()));
        float latitude = (float) Math.toDegrees((Math.asin(v.getZ() / radius)));
        return new float[]{latitude, longitude};
    }

    public class Polygon {
        int[] polyVertices;
        Polygon[] neighbors;
        int center;

        Polygon(int vertices) {
            this.polyVertices = new int[vertices];
            this.neighbors = new Polygon[vertices];
        }

        public int[] getPolyVertices() {
            return polyVertices;
        }

        public void addNeighbor(Polygon polygon) {
            for (int i = 0; i < neighbors.length; i++) {
                if (neighbors[i] == null) {
                    neighbors[i] = polygon;
                    break;
                }
            }
        }

        public int getCenter() {
            return center;
        }

        public int getVerticesCount() {
            return polyVertices.length;
        }

        @Override
        public String toString() {
            return "Polygon{" +
                    "vertices=" + Arrays.toString(polyVertices) +
                    ", center=" + center +
                    '}';
        }

        public int getVertex(int i) {
            return polyVertices[i];
        }

        public void sortVertices() {
            // Use x, y, z to arrange vertices in clockwise order
            for (int i = 0; i < polyVertices.length; i++) {
                if (polyVertices[i] == NOT_EXIST) {
                    break;
                }
                float angle = (float) Math.atan2(vertices.get(polyVertices[i]).getY(), vertices.get(polyVertices[i]).getX());
                for (int j = i + 1; j < polyVertices.length; j++) {
                    if (polyVertices[j] == NOT_EXIST) {
                        break;
                    }
                    float angle2 = (float) Math.atan2(vertices.get(polyVertices[j]).getY(), vertices.get(polyVertices[j]).getX());
                    if (angle2 < angle) {
                        int temp = polyVertices[i];
                        polyVertices[i] = polyVertices[j];
                        polyVertices[j] = temp;
                    }
                }
            }
        }
    }

    public class Pentagon extends Polygon {
        Pentagon() {
            super(5);
        }
    }

    public class Hexagon extends Polygon {
        Hexagon() {
            super(6);
        }
    }

    public class Triangle {
        @SuppressWarnings("checkstyle:VisibilityModifier")
        int a;
        @SuppressWarnings("checkstyle:VisibilityModifier")
        int b;
        @SuppressWarnings("checkstyle:VisibilityModifier")
        int c;

        Triangle(int a, int b, int c) {
            this.a = a;
            this.b = b;
            this.c = c;
        }

    }

    public static void main(String[] args) {
        GeoSubdivisionHelper gsh = new GeoSubdivisionHelper();
        gsh.setup(4, 6371);
        RockyBody body = new RockyBody(null, "Test", 6371, 6371, 0);
        List<Polygon> polygons = gsh.getRegions();

        for (Polygon region : polygons) {
            int[] polyVertices = region.getPolyVertices();
            for (int i = 0; i < polyVertices.length; i++) {
                for (int j = i + 1; j < polyVertices.length; j++) {
                    if (polyVertices[i] != NOT_EXIST && polyVertices[i] == polyVertices[j]) {
                        throw new RuntimeException("Duplicate vertices in region " + region);
                    } else if (polyVertices[i] != NOT_EXIST && polyVertices[j] != NOT_EXIST) {
                        Vector3D v1 = gsh.getVertex(polyVertices[i]);
                        Vector3D v2 = gsh.getVertex(polyVertices[j]);
                        if (v1.equals(v2)) {
                            System.out.printf("Duplicate vertices: %s, %s\n", v1, v2);
                            throw new RuntimeException("Duplicate vertices in region " + region + " at " + i + " and " + j);
                        }
                    }

                }
            }
            System.out.println(region);
        }
    }

}