Fixes to Newell algorithm implementation

Signed-off-by: Antonio Ospite <ospite@studenti.unina.it>

diff --git a/vrm.py b/vrm.py
index 9095307..82922f8 100755 (executable)
--- a/vrm.py
+++ b/vrm.py
@@ -92,7 +92,7 @@ class config:
     polygons = dict()
     polygons['SHOW'] = True
     polygons['SHADING'] = 'FLAT'
     polygons = dict()
     polygons['SHOW'] = True
     polygons['SHADING'] = 'FLAT'

-    polygons['HSR'] = 'PAINTER' # 'PAINTER' or 'NEWELL'
+    #polygons['HSR'] = 'PAINTER' # 'PAINTER' or 'NEWELL'

     polygons['HSR'] = 'NEWELL'
     # Hidden to the user for now
     polygons['EXPANSION_TRICK'] = True
     polygons['HSR'] = 'NEWELL'
     # Hidden to the user for now
     polygons['EXPANSION_TRICK'] = True


@@ -113,12 +113,20 @@ class config:
 
 
 
 
 
 

-# Debug utility function
+# Utility functions

 print_debug = True
 def debug(msg):
     if print_debug:
         sys.stderr.write(msg)
 
 print_debug = True
 def debug(msg):
     if print_debug:
         sys.stderr.write(msg)
 

+def sign(x):
+    if x == 0:
+        return 0
+    elif x < 0:
+        return -1
+    else:
+        return 1
+

 
 # ---------------------------------------------------------------------
 #
 
 # ---------------------------------------------------------------------
 #


@@ -1025,6 +1033,13 @@ class Renderer:
             self._doModelingTransformation(mesh, obj.matrix)
 
             self._doBackFaceCulling(mesh)
             self._doModelingTransformation(mesh, obj.matrix)
 
             self._doBackFaceCulling(mesh)

+            if True:
+                for f in mesh.faces:
+                    f.sel = 1-f.sel
+                mesh.flipNormals()
+                for f in mesh.faces:
+                    f.sel = 1
+

 
             self._doLighting(mesh)
 
 
             self._doLighting(mesh)
 


@@ -1466,49 +1481,55 @@ class Renderer:
                 cmp(max([v.co[2] for v in f1]), max([v.co[2] for v in f2]))
                 )
 
                 cmp(max([v.co[2] for v in f1]), max([v.co[2] for v in f2]))
                 )
 

+        mesh.quadToTriangle(0)

 
 

-        def isOnSegment(v1, v2, p):
-
-            # when p is at extreme points
-            if p == v1 or p == v2:
-                return False
+        from split import Distance, isOnSegment

 
 

+        def projectionsOverlap(P, Q):

 
 

-            EPS = 10e-7
+            for i in range(0, len(P.v)):

 
 

-            l1 = (v1-p).length
-            l2 = (v2-p).length
-            l = (v1-v2).length
+                v1 = Vector(P.v[i-1])
+                v1[2] = 0
+                v2 = Vector(P.v[i])
+                v2[2] = 0

 
 

-            print "l: ", l, " l1: ", l1, " l2: ", l2, "diff: %.9f" % (l - (l1+l2) )
-
-            if abs(l - (l1+l2)) < EPS:
-                return True
-            else:
-                return False
+                EPS = 10e-7

 
 

+                for j in range(0, len(Q.v)):
+                    v3 = Vector(Q.v[j-1])
+                    v3[2] = 0
+                    v4 = Vector(Q.v[j])
+                    v4[2] = 0
+                    
+                    ret = LineIntersect(v1, v2, v3, v4)
+                    # if line v1-v2 and v3-v4 intersect both return
+                    # values are the same.
+                    if ret and ret[0] == ret[1]  and isOnSegment(v1, v2,
+                            ret[0], True) and isOnSegment(v3, v4, ret[1], True):

 
 
 
 

-        def Distance(point, face):
-            """ Calculate the distance between a point and a face.
+                        l1 = (ret[0] - v1).length
+                        l2 = (ret[0] - v2).length

 
 

-            An alternative but more expensive method can be:
+                        l3 = (ret[1] - v3).length
+                        l4 = (ret[1] - v4).length

 
 

-                ip = Intersect(Vector(face[0]), Vector(face[1]), Vector(face[2]),
-                        Vector(face.no), Vector(point), 0)
-
-                d = Vector(ip - point).length
-            """
+                        if  (l1 < EPS or l2 < EPS) and (l3 < EPS or l4 < EPS):
+                            continue

 
 

-            plNormal = Vector(face.no)
-            plVert0 = Vector(face[0])
+                        debug("Projections OVERLAP!!\n")
+                        debug("line1:"+
+                                " M "+ str(v1[0])+','+str(v1[1]) + ' L ' + str(v2[0])+','+str(v2[1]) + '\n' +
+                                " M "+ str(v3[0])+','+str(v3[1]) + ' L ' + str(v4[0])+','+str(v4[1]) + '\n' +
+                                "\n")
+                        debug("return: "+ str(ret)+"\n")
+                        return True

 
 

-            #d = abs( (point * plNormal ) - (plVert0 * plNormal) )
-            d = (point * plNormal ) - (plVert0 * plNormal)
-            debug("d: %.10f - sel: %d, %s\n" % (d, face.sel, str(point)) )
+            return False

 
 

-            return d

 
 

+        from facesplit import facesplit

 
         # FIXME: using NMesh to sort faces. We should avoid that!
         nmesh = NMesh.GetRaw(mesh.name)
 
         # FIXME: using NMesh to sort faces. We should avoid that!
         nmesh = NMesh.GetRaw(mesh.name)


@@ -1527,31 +1548,58 @@ class Renderer:
         facelist = nmesh.faces[:]
         maplist = []
 
         facelist = nmesh.faces[:]
         maplist = []
 

-        #EPS = 10e-7
-        EPS = 0
+        EPS = 10e-8
+        #EPS = 0

 
         global progress
 
         global progress

+        # The steps are _at_least_ equal to len(facelist), we do not count the
+        # feces coming out from plitting!!

         progress.setActivity("HSR: Newell", len(facelist))
         progress.setQuiet(True)
 
         progress.setActivity("HSR: Newell", len(facelist))
         progress.setQuiet(True)
 

-        #while len(facelist)-1:
+        
+        steps = -1
+        split_done = 0
+        marked_face = 0
+

         while len(facelist):
         while len(facelist):

+            print "\n----------------------"

             P = facelist[0]
             P = facelist[0]

+            
+            #steps += 1
+            #if steps == 3:
+            #    maplist = facelist
+            #    break
+            print len(facelist)
+            if len(facelist) == 33:
+                #maplist = facelist
+                break
+
+
+            #pSign = 1
+            #if P.normal[2] < 0:
+            #    pSign = -1
+            pSign = sign(P.normal[2])
+
+            # We can discard faces thar are perpendicular to the view
+            if pSign == 0:
+                facelist.remove(P)
+                continue
+

 
 

-            pSign = 1
-            if P.sel == 0:
-                pSign = -1
+            split_done = 0
+            face_marked = 0

 
 

-            #while False:

             for Q in facelist[1:]:
 
                 debug("P.smooth: " + str(P.smooth) + "\n")
                 debug("Q.smooth: " + str(Q.smooth) + "\n")
                 debug("\n")
 
             for Q in facelist[1:]:
 
                 debug("P.smooth: " + str(P.smooth) + "\n")
                 debug("Q.smooth: " + str(Q.smooth) + "\n")
                 debug("\n")
 

-                qSign = 1
-                if Q.sel == 0:
-                    qSign = -1
+                #qSign = 1
+                #if Q.normal[2] < 0:
+                #    qSign = -1
+                qSign = sign(Q.normal[2])

  
                 # We need to test only those Qs whose furthest vertex
                 # is closer to the observer than the closest vertex of P.
  
                 # We need to test only those Qs whose furthest vertex
                 # is closer to the observer than the closest vertex of P.


@@ -1563,16 +1611,17 @@ class Renderer:
                 if not ZOverlap:
                     debug("\nTest 0\n")
                     debug("NOT Z OVERLAP!\n")
                 if not ZOverlap:
                     debug("\nTest 0\n")
                     debug("NOT Z OVERLAP!\n")

-                    if not Q.smooth:
-                        # We can safely print P
+                    if Q.smooth == 0:
+                        # If Q is not marked then we can safely print P

                         break
                     else:
                         break
                     else:

+                        debug("met a marked face\n")

                         continue
                 
                 # Test 1: X extent overlapping
                 xP = [v.co[0] for v in P.v]
                 xQ = [v.co[0] for v in Q.v]
                         continue
                 
                 # Test 1: X extent overlapping
                 xP = [v.co[0] for v in P.v]
                 xQ = [v.co[0] for v in Q.v]

-                notXOverlap = (max(xP) < min(xQ)) or (max(xQ) < min(xP))
+                notXOverlap = (max(xP) <= min(xQ)) or (max(xQ) <= min(xP))

 
                 if notXOverlap:
                     debug("\nTest 1\n")
 
                 if notXOverlap:
                     debug("\nTest 1\n")


@@ -1582,7 +1631,7 @@ class Renderer:
                 # Test 2: Y extent Overlapping
                 yP = [v.co[1] for v in P.v]
                 yQ = [v.co[1] for v in Q.v]
                 # Test 2: Y extent Overlapping
                 yP = [v.co[1] for v in P.v]
                 yQ = [v.co[1] for v in Q.v]

-                notYOverlap = (max(yP) < min(yQ)) or (max(yQ) < min(yP))
+                notYOverlap = (max(yP) <= min(yQ)) or (max(yQ) <= min(yP))

 
                 if notYOverlap:
                     debug("\nTest 2\n")
 
                 if notYOverlap:
                     debug("\nTest 2\n")


@@ -1593,9 +1642,8 @@ class Renderer:
                 # Test 3: P vertices are all behind the plane of Q
                 n = 0
                 for Pi in P:
                 # Test 3: P vertices are all behind the plane of Q
                 n = 0
                 for Pi in P:

-                    print P.col[0]

                     d = qSign * Distance(Vector(Pi), Q)
                     d = qSign * Distance(Vector(Pi), Q)

-                    if d > EPS:
+                    if d >= -EPS:

                         n += 1
                 pVerticesBehindPlaneQ = (n == len(P))
 
                         n += 1
                 pVerticesBehindPlaneQ = (n == len(P))
 


@@ -1608,7 +1656,6 @@ class Renderer:
                 # Test 4: Q vertices in front of the plane of P
                 n = 0
                 for Qi in Q:
                 # Test 4: Q vertices in front of the plane of P
                 n = 0
                 for Qi in Q:

-                    print Q.col[0]

                     d = pSign * Distance(Vector(Qi), P)
                     if d <= EPS:
                         n += 1
                     d = pSign * Distance(Vector(Qi), P)
                     if d <= EPS:
                         n += 1


@@ -1621,52 +1668,28 @@ class Renderer:
 
                 # Test 5: Line Intersections... TODO
                 # Check if polygons effectively overlap each other, not only
 
                 # Test 5: Line Intersections... TODO
                 # Check if polygons effectively overlap each other, not only

-                # boundig boxes as dome before.
+                # boundig boxes as done before.

                 # Since we We are working in normalized projection coordinates
                 # we kust check if polygons intersect.
 
                 # Since we We are working in normalized projection coordinates
                 # we kust check if polygons intersect.
 

-                def projectionsOverlap(P, Q):
-
-                    for i in range(0, len(P.v)):
-
-                        v1 = Vector(P.v[i-1])
-                        v1[2] = 0
-                        v2 = Vector(P.v[i])
-                        v2[2] = 0
-
-                        for j in range(0, len(Q.v)):
-                            v3 = Vector(Q.v[j-1])
-                            v3[2] = 0
-                            v4 = Vector(Q.v[j])
-                            v4[2] = 0
-                            
-                            ret = LineIntersect(v1, v2, v3, v4)
-                            # if line v1-v2 and v3-v4 intersect both return
-                            # values are the same.
-                            if ret and ret[0] == ret[1]  and isOnSegment(v1, v2,
-                                    ret[0]) and isOnSegment(v3, v4, ret[1]):
-                                debug("Projections OVERLAP!!\n")
-                                debug("line1:"+
-                                        " M "+ str(v1[0])+','+str(v1[1]) + ' L ' + str(v2[0])+','+str(v2[1]) + '\n' +
-                                        " M "+ str(v3[0])+','+str(v3[1]) + ' L ' + str(v4[0])+','+str(v4[1]) + '\n' +
-                                        "\n")
-                                debug("return: "+ str(ret)+"\n")
-                                return True
-
-                    return False
-

                 if not projectionsOverlap(P, Q):
                     debug("\nTest 5\n")
                     debug("Projections do not overlap!\n")
                     continue
 
 
                 if not projectionsOverlap(P, Q):
                     debug("\nTest 5\n")
                     debug("Projections do not overlap!\n")
                     continue
 
 

-                # We do not know if P obscures Q.
+                # We still do not know if P obscures Q.
+
+                # But if Q is marked we do a split trying to resolve a
+                # difficulty (maybe a visibility cycle).

                 if Q.smooth == 1:
                 if Q.smooth == 1:

-                    # Split P or Q, TODO
-                    debug("Cycle detected!\n")
+                    # Split P or Q
+                    debug("Possibly a cycle detected!\n")

                     debug("Split here!!\n")
                     debug("Split here!!\n")

-                    continue
+                    old_facelist = facelist[:]
+                    facelist = facesplit(P, Q, facelist, nmesh)
+                    split_done = 1
+                    break 

 
 
                 # The question now is: Does Q obscure P?
 
 
                 # The question now is: Does Q obscure P?


@@ -1674,9 +1697,8 @@ class Renderer:
                 # Test 3bis: Q vertices are all behind the plane of P
                 n = 0
                 for Qi in Q:
                 # Test 3bis: Q vertices are all behind the plane of P
                 n = 0
                 for Qi in Q:

-                    print Q.col[0]

                     d = pSign * Distance(Vector(Qi), P)
                     d = pSign * Distance(Vector(Qi), P)

-                    if d > EPS:
+                    if d >= -EPS:

                         n += 1
                 qVerticesBehindPlaneP = (n == len(Q))
 
                         n += 1
                 qVerticesBehindPlaneP = (n == len(Q))
 


@@ -1688,7 +1710,6 @@ class Renderer:
                 # Test 4bis: P vertices in front of the plane of Q
                 n = 0
                 for Pi in P:
                 # Test 4bis: P vertices in front of the plane of Q
                 n = 0
                 for Pi in P:

-                    print P.col[0]

                     d = qSign * Distance(Vector(Pi), Q)
                     if d <= EPS:
                         n += 1
                     d = qSign * Distance(Vector(Pi), Q)
                     if d <= EPS:
                         n += 1


@@ -1698,58 +1719,84 @@ class Renderer:
                     debug("\nTest 4bis\n")
                     debug("P IN FRONT OF Q!\n")
 
                     debug("\nTest 4bis\n")
                     debug("P IN FRONT OF Q!\n")
 

-
-                import intersection
-
+                
+                # We don't even know if Q does obscure P, so they should
+                # intersect each other, split one of them in two parts.

                 if not qVerticesBehindPlaneP and not pVerticesInFrontPlaneQ:
                     debug("\nSimple Intersection?\n")
                 if not qVerticesBehindPlaneP and not pVerticesInFrontPlaneQ:
                     debug("\nSimple Intersection?\n")

-                    # Split P or Q, TODO
-                    print "Test 3bis or 4bis failed"
-                    print "Split here!!2\n"
-
-                    """newfaces = intersection.splitOn(P, Q, 0)
-                    print newfaces
-                    facelist.remove(Q)
-                    for nf in newfaces:
-                        if nf:
-                            nf.col = Q.col
-                            facelist.append(nf)
-                    """
-
-                    break
-
-                # We do not know
-                if Q.smooth:
-                    # split P or Q
-                    print "Split here!!\n"
-                    """
-                    newfaces = intersection.splitOn(P, Q, 0)
-                    facelist.remove(Q)
-                    for nf in newfaces:
-                        if nf:
-                            nf.col = Q.col
-                            facelist.append(nf)
-
-                    """ 
-                    break
+                    debug("Test 3bis or 4bis failed\n")
+                    debug("Split here!!2\n")

 
 

-                Q.smooth = 1
+                    old_facelist = facelist[:]
+                    facelist = facesplit(P, Q, facelist, nmesh)
+
+                    steps += 1
+                    if steps == 2:
+                        maplist = [P, Q]
+                        print P, Q
+                    split_done = 1
+                    break 
+
+                    

                 facelist.remove(Q)
                 facelist.insert(0, Q)
                 facelist.remove(Q)
                 facelist.insert(0, Q)

+                Q.smooth = 1
+                face_marked = 1
+
+                # Make merked faces BLUE. so to see them
+                #for c in Q.col:
+                #    c.r = 0
+                #    c.g = 0
+                #    c.b = 255
+                #    c.a = 255
+                
+                debug("Q marked!\n")
+                print [f.smooth for f in facelist]
+                break

            
             # Write P!                     
            
             # Write P!                     

-            P = facelist[0]
-            facelist.remove(P)
-            maplist.append(P)
+            if split_done == 0 and face_marked == 0:
+                P = facelist[0]
+                facelist.remove(P)
+                maplist.append(P)
+
+                progress.update()
+                #if progress.progressModel.getProgress() == 100:
+                #    break
+            if steps == 2:
+                """
+                for c in Q.col:
+                    c.r = 0
+                    c.g = 0
+                    c.b = 255
+                    c.a = 255
+                for c in P.col:
+                    c.r = 0
+                    c.g = 0
+                    c.b = 255
+                    c.a = 255
+                """
+                print steps
+                #maplist.append(P)
+                #maplist.append(Q)
+
+               # for f in facelist:
+               #     if f not in old_facelist:
+               #         print "splitted?"
+               #         maplist.append(f)
+
+                break
+            """
+            """

 
 

-            progress .update()
+         

 
 

-        

         nmesh.faces = maplist
 
         for f in nmesh.faces:
             f.sel = 1
         nmesh.update()
         nmesh.faces = maplist
 
         for f in nmesh.faces:
             f.sel = 1
         nmesh.update()

+        print nmesh.faces

 
     def _doHiddenSurfaceRemoval(self, mesh):
         """Do HSR for the given mesh.
 
     def _doHiddenSurfaceRemoval(self, mesh):
         """Do HSR for the given mesh.