Podwójne konturowanie - Znalezienie punktu charakterystycznego, normalne wyłączone

Wykonuję ten samouczek, aby wdrożyć Dual Contouring http://www.sandboxie.com/misc/isosurf/isosurfaces.html

Moje źródło danych to siatka 16 x 16 x 16; Przesuwam tę siatkę od dołu do góry, od lewej do prawej, blisko do daleko.

Dla każdego indeksu mojej siatki tworzę strukturę kostki:

public Cube(int x, int y, int z, Func<int, int, int, IsoData> d, float isoLevel) {
            this.pos = new Vector3(x,y,z);
            //only create vertices need for edges
            Vector3[] v = new Vector3[4];
            v[0] = new Vector3 (x + 1, y + 1, z);
            v[1] = new Vector3 (x + 1, y, z + 1);
            v[2] = new Vector3 (x + 1, y + 1, z + 1);
            v[3] = new Vector3 (x, y + 1, z + 1);
            //create edges from vertices
            this.edges = new Edge[3];
            edges[0] = new Edge (v[1], v[2], d, isoLevel);
            edges[1] = new Edge (v[2], v[3], d, isoLevel);
            edges[2] = new Edge (v[0], v[2], d, isoLevel);
        }

Ze względu na to, jak przemierzam siatkę, muszę patrzeć tylko na 4 wierzchołki i 3 krawędzie. W tym obrazie wierzchołki 2, 5, 6, 7 odpowiadają moim wierzchołkom 0, 1, 2, 3, a krawędzie 5, 6, 10 odpowiadają moim krawędziom 0, 1, 2.

Krawędź wygląda następująco:

    public Edge(Vector3 p0, Vector3 p1, Func<int, int, int, IsoData> d, float isoLevel) {
        //get density values for edge vertices, save in vector , d = density function, data.z = isolevel 
        this.data = new Vector3(d ((int)p0.x, (int)p0.y, (int)p0.z).Value, d ((int)p1.x, (int)p1.y, (int)p1.z).Value, isoLevel);
        //get intersection point
        this.mid = LerpByDensity(p0,p1,data);
        //calculate normals by gradient of surface
        Vector3 n0 = new Vector3(d((int)(p0.x+1),   (int)p0.y,      (int)p0.z       ).Value - data.x,
                                 d((int)p0.x,       (int)(p0.y+1),  (int)p0.z       ).Value - data.x,
                                 d((int)p0.x,       (int)p0.y,      (int)(p0.z+1)   ).Value - data.x);

        Vector3 n1 = new Vector3(d((int)(p1.x+1),   (int)p1.y,      (int)p1.z       ).Value - data.y,
                                 d((int)p1.x,       (int)(p1.y+1),  (int)p1.z       ).Value - data.y,
                                 d((int)p1.x,       (int)p1.y,      (int)(p1.z+1)   ).Value - data.y);
        //calculate normal by averaging normal of edge vertices
        this.normal = LerpByDensity(n0,n1,data);
    }

Następnie sprawdzam wszystkie krawędzie pod kątem zmiany znaku, jeśli jest taka, znajduję otaczające kostki i uzyskuję punkt funkcji tych kostek.

Teraz to działa, jeśli ustawię punkt funkcji na środek sześcianu, a następnie uzyskam blokowy wygląd Minecrafta. Ale nie tego chcę.

Aby znaleźć punkt funkcji, chciałem to zrobić jak w tym poście: https://gamedev.stackexchange.com/a/83757/49583

Zasadniczo zaczynasz wierzchołek w środku komórki. Następnie uśredniasz wszystkie wektory pobrane z wierzchołka do każdej płaszczyzny i przesuwasz wierzchołek wzdłuż tej wypadkowej, i powtarzasz ten krok określoną liczbę razy. Przekonałem się, że przesunięcie o ~ 70% wzdłuż wypadkowej ustabilizuje się w najmniejszej liczbie iteracji.

Mam więc klasę samolotu:

private class Plane {

        public Vector3 normal;
        public float distance;

        public Plane(Vector3 point, Vector3 normal) {
            this.normal = Vector3.Normalize(normal);
            this.distance = -Vector3.Dot(normal,point);
        }

        public float Distance(Vector3 point) {
            return Vector3.Dot(this.normal, point) + this.distance;
        }

        public Vector3 ShortestDistanceVector(Vector3 point) {
            return this.normal * Distance(point);
        }
 }

oraz funkcję umożliwiającą uzyskanie punktu charakterystycznego, w którym tworzę 3 płaszczyzny, po jednej dla każdej krawędzi i uśredniając odległość do środka:

 public Vector3 FeaturePoint {
            get {
                Vector3 c = Center;
 //                 return c; //minecraft style

                Plane p0 = new Plane(edges[0].mid,edges[0].normal);
                Plane p1 = new Plane(edges[1].mid,edges[1].normal);
                Plane p2 = new Plane(edges[2].mid,edges[2].normal);

                int iterations = 5;
                for(int i = 0; i < iterations; i++) {
                    Vector3 v0 = p0.ShortestDistanceVector(c);
                    Vector3 v1 = p1.ShortestDistanceVector(c);
                    Vector3 v2 = p2.ShortestDistanceVector(c);
                    Vector3 avg = (v0+v1+v2)/3;
                    c += avg * 0.7f;
                }

                return c;
            }
        }

Ale to nie działa, wierzchołki są wszędzie. Gdzie jest błąd? Czy faktycznie mogę obliczyć normalną krawędź, uśredniając normalną wierzchołek krawędzi? Nie mogę uzyskać gęstości w środkowym punkcie krawędzi, ponieważ mam tylko siatkę liczb całkowitych jako źródło danych ...

Edycja: Znalazłem również tutaj http://www.mathsisfun.com/algebra/systems-linear-equations-matrices.html , że mogę użyć macierzy do obliczenia przecięcia 3 płaszczyzn, przynajmniej tak to rozumiałem, więc Stworzyłem tę metodę

 public static Vector3 GetIntersection(Plane p0, Plane p1, Plane p2) {              
            Vector3 b = new Vector3(-p0.distance, -p1.distance, -p2.distance);

            Matrix4x4 A = new Matrix4x4 ();
            A.SetRow (0, new Vector4 (p0.normal.x, p0.normal.y, p0.normal.z, 0));
            A.SetRow (1, new Vector4 (p1.normal.x, p1.normal.y, p1.normal.z, 0));
            A.SetRow (2, new Vector4 (p2.normal.x, p2.normal.y, p2.normal.z, 0));
            A.SetRow (3, new Vector4 (0, 0, 0, 1));

            Matrix4x4 Ainv = Matrix4x4.Inverse(A);

            Vector3 result = Ainv * b;
            return result;
        }

które z tymi danymi

        Plane p0 = new Plane (new Vector3 (2, 0, 0), new Vector3 (1, 0, 0));
        Plane p1 = new Plane (new Vector3 (0, 2, 0), new Vector3 (0, 1, 0));
        Plane p2 = new Plane (new Vector3 (0, 0, 2), new Vector3 (0, 0, 1));

        Vector3 cq = Plane.GetIntersection (p0, p1, p2);

oblicza przecięcie na (2.0, 2.0, 2.0), więc zakładam, że działa poprawnie. Mimo to nieprawidłowe wierzchołki. Naprawdę myślę, że to moje normalne.

Przede wszystkim twoje normalne powinny być w porządku, jeśli są obliczane na podstawie różnic do przodu / do tyłu / centralnie. Problem polega na tym, że przesunąłeś punkt środkowy w kierunku przeciwnym do niewłaściwego kierunku w funkcji FeaturePoint, co skutkuje oddaleniem się od minimum.

Vector3 c = Center;
Plane p0 = new Plane(edges[0].mid,edges[0].normal);
Plane p1 = new Plane(edges[1].mid,edges[1].normal);
Plane p2 = new Plane(edges[2].mid,edges[2].normal);

int iterations = 5;
for(int i = 0; i < iterations; i++) {
    Vector3 v0 = p0.GetDistanceToPoint(c) * edges[0].normal;
    Vector3 v1 = p1.GetDistanceToPoint(c) * edges[1].normal;
    Vector3 v2 = p2.GetDistanceToPoint(c) * edges[2].normal;
    Vector3 avg = (v0+v1+v2)/3;
    c -= avg * 0.7f; // Error was here!
}
return c;

Dzieje się tak, ponieważ Twój kod nie jest zbieżny z punktem i dlatego wyskakuje z okna wokseli. Nie wiem, czy kod z Czy ktoś może wyjaśnić podwójne konturowanie? miał zastosować rzutowanie, w którym punkt jest rzutowany na płaszczyznę poprzez:

distance = Vector3.Dot(point - origin, normal);
projectedPoint = point - distance * normal;

ale to ta sama metoda. Jeśli przepiszesz projekcję do oryginalnego kodu, spowoduje to:

    Vector3 v0 = c - p0.GetDistanceToPoint(c) * edges[0].normal;
    Vector3 v1 = c - p1.GetDistanceToPoint(c) * edges[1].normal;
    Vector3 v2 = c - p2.GetDistanceToPoint(c) * edges[2].normal;
    c = (v0+v1+v2)/3;

które można przepisać na:

    Vector3 v0 = p0.GetDistanceToPoint(c) * edges[0].normal;
    Vector3 v1 = p1.GetDistanceToPoint(c) * edges[1].normal;
    Vector3 v2 = p2.GetDistanceToPoint(c) * edges[2].normal;
    c = c - (v0+v1+v2)/3;

i dlatego skutkuje pierwszym kodem. Rzutując punkt na trzech płaszczyznach nieplanarnych, powoli zbiega się w kierunku minimum, ponieważ minimalizujesz odległość od każdej płaszczyzny do punktu, jak pokazano na zdjęciu.

Czerwone kropki oznaczają punkt charakterystyczny, niebieskie linie normalne, a fioletowy punkt rzutowany na płaszczyznę. Nie musisz także używać współczynnika 0,7, ponieważ bez niego powinien on zbiegać się szybciej. Jeśli używasz tej metody, strzeż się, że algorytm może nie działać, jeśli masz nie przecinające się płaszczyzny.