#=============================================================================================== # ray2edge_dist() # # Note: this routine is fairly specific to the task at hand, so the name of it shouldn't be taken # as a general problem solver. # # What this routine does is compute the distance from a vertex to the opposing edge of a triangle, # along the vector described by the line from the vertex to the intersection point that was # calculated in the linecast_loop() code. # # The methodology I'll use to come up with this is as follows: # # 1.create a 'plane' from the edge being tested, perpendicular to the plane described # by the intersection point and the two points of the edge . # # 2.calculate the intersection between a ray cast from the vertex in question (along # the vertex-->intersection vector), with the plane created in step one, above. # # ...note that once we have the new plane from step #1, this is basically the same methodology # we used back in linecast_loop() to compute the original intersection point, just twisted # around in space to deal with the new ray direction. # #=============================================================================================== def ray2edge_dist(vert,e1,e2,pnorm,isect): """ Computes the distance from one vertex of a triangle to the opposing edge. vert - one vertex of a triangle e1,e2 - the vertices that make up the opposing edge pnorm - the normal of the triangle isect - the point of intersection (some point within the triangle) """ # Create plane for e1..e2 edge, perpendicular to triangle evec = vecsub(e2,e1) # Get edge vector perp = vector_normalize(vector_crossproduct(evec,pnorm)) # Get vector perpendicular to polygon plane = vector_dotproduct(perp,e1) # Convert to plane rvec = vector_normalize(vecsub(isect,vert)) # Get ray vector from vert-->intersect ndota = vector_dotproduct(perp,rvec) ndotv = vector_dotproduct(perp,vert) if ndota <= 0.0: # divide-by-zero sanity check - due to some math rounding error that if pcage: testing_boxes(e1,0.003,"edge_MISS_",0,1) return FLT_MAX # was adjusted for the initial point_in_triangle test, it's possible # that the point is at the vertex being tested, which could cause the # ray vector to be wacky, so just return a 'maximum' distance. # Compute intersection distance from vert to edge, along vector from vert-->intersect edist = (plane - ndotv) / ndota return edist #=============================================================================================== # get_raycast_weights() # # For each vertex of the triangle, it's weighting value is computed as: # # weight = 1.0 - (ISectDist divided by EdgeDist) # # [...if interested, see additional notes and diagrams here: # http:#www.renderosity.com/mod/forumpro/showthread.php?thread_id=2677445&page=2 ] # # So far, we have an ISectDist for each vertex-->intersection (idv0..idv2). # What we don't have yet is the distance from each vertex to the opposing edge, # along the vector described by the line from the vertex to the intersection. # # Optimization: The nature of these weights is that the sum of them add up to 1.0, # so we really only need to calculate the first 2 weights and then we can deduce # the 3rd based on those results. #=============================================================================================== def get_raycast_weights(ixPoint, polyNorm, v0, v1, v2) """ """ # compute weight for v0, using v1..v2 edge edist = ray2edge_dist(v0,v1,v2,polyNorm,ixPoint) if edist >= (FLT_MAX-0.0000001): # FLT_MAX is the signal that the intersection was _at_ this vertex weights[0] = 1.0; weights[1] = 0.0; weights[2] = 0.0; # (or at least so close that the error-correcting code let it through) return weights if edist <= 0.0: # don't divide by zero (this would be a degenerate triangle) weights[0] = 0.0; else: idist = point_distance(v0, ixPoint) # distance from vert of the tripoly to the intersection point weights[0] = 1.0 - (idist / edist) # solve for weight = 1.0 - (ISectDist divided by EdgeDist) if weights[0] < 0.0: weights[0] = 0.0; # repeat the above for v1, using v2..v0 edge edist = ray2edge_dist(v1,v2,v0,polyNorm,ixPoint); if edist >= (FLT_MAX-0.0000001): # FLT_MAX is the signal that the intersection was _at_ this vertex weights[0] = 0.0; weights[1] = 1.0; weights[2] = 0.0; # (or at least so close that the error-correcting code let it through) return weights if edist <= 0.0: # don't divide by zero (this would be a degenerate triangle) weights[1] = 0.0; else: idist = point_distance(v1, ixPoint) # distance from vert of the tripoly to the intersection point weights[1] = 1.0 - (idist / edist) # solve for weight = 1.0 - (ISectDist divided by EdgeDist) if weights[1] < 0.0: weights[1] = 0.0; # wgt2 gets whatever is left over... weights[2] = 1.0 - (weights[0] + weights[1]) if weights[2] < 0.0: weights[2] = 0.0 return weights #=============================================================================================== # correlateToNearVertList() # # This is the Python version of the mesh.CorrelateToNearVertList() routine. # # NOTE: I was about half-way through porting this over to Python code when I realized that I was # actually just porting it to a Python version of the C++ code, but still requiring/relying # on the .pyd. I went ahead and started pythonizing more of it, but there are still some # dependencies left that will need to be addressed. # # With that said, see the TODO section, below # #----------------------------------------------------------------------------------------------- #*************************************** TODO: ************************************************* #----------------------------------------------------------------------------------------------- # 1. # Some variables/lists needed by this routine (currently fetched from .pyd Mesh class, below, # but would otherwise likely need to be passed in)... # # sourceTriPolys -> list of tripolys that make up the source mesh # sourceVertList -> list of vertices that make up the source mesh # srcVertTriPolys -> list of tripoly indices per vertex of source mesh # sourceNormList -> list of source mesh Normals # # ...probably also going to need a list of plane equations for the source mesh tripolys, as # well as a list of normals for those tripolys (the code below is still set up to work with # the .pyd TriPoly class, where the plane and normal are already available) # #----------------------------------------------------------------------------------------------- # 2. # Error-checking has not been ported to Python (it's commented out, below)... this should be # implemented, noting that the C++ callError() routine would just print out an error message # then return (Null or whatever this normally returns as an error/failure condition). # #----------------------------------------------------------------------------------------------- # 3. # Return value(s)... the current Hybrid code is expecting a filled in HitPoint class returned, # but currently only uses the weights and the tripoly hit (in order to get the vertex indices # that those weights are associated with - so it could just also return those indices). Either # way, I mostly just left that unfinished, so I'll leave that to you :). # #----------------------------------------------------------------------------------------------- # 4. # The 2 routines above here were either grabbed directly from out old scripts, or modified from # them, so they should be good-to-go, but they might be calling routines from some older code # as well (which you may no longer have implemented in the new script). #----------------------------------------------------------------------------------------------- #=============================================================================================== def correlateToNearVertList(targetVert, targetNorm, sourceMesh, chkVertList, avgNorms=1): """ """ sourceTriPolys = sourceMesh.GetTriPolys() # complete list of tripolys that make up the source mesh #if( !sourceTriPolys ) return NULL; sourceVertList = sourceMesh.GetWorldVertices() # complete list of vertices that make up the above tripolys #if( !sourceVertList ) return NULL; srcVertTriPolys = sourceMesh.GetVertTriPolys() # complete list of tripoly indices per vertex of source mesh #if( !srcVertTriPolys ) return NULL; numSourceVerts = len(sourceVertList) #if( !numSourceVerts ) callError(funcname, "NumSourceVerts is 0"); numCheckVerts = len(chkVertList) #if( !numCheckVerts ) callError(funcname, "NumCorrelateVerts is 0"); rayNorm = targetNorm # copy the normal (so we can safely modify it) rayVert = targetVert # copy the vertex, just for naming consistency #------------------------------------------------------------------------------------ # experimental normal averaging... currently, it averages the target vertex normal # with the normals of the closest source mesh vertices #------------------------------------------------------------------------------------ if( avgNorms ): sourceNormList = sourceMesh.GetWorldNormals() # complete list of source mesh Normals #if( !sourceNormList ) return NULL; normCnt = 1 for ci in range(numCheckVerts): cvNdx = chkVertList[ci] #if( cvNdx < 0 || cvNdx >= numSourceVerts ) #{ # char erbuf[256]; # sprintf(erbuf, "cvNdx [%d of %d] is out of range", cvNdx, numSourceVerts); # callError(funcname, erbuf); #} chkNorm = sourceNormList[cvNdx] if not chkNorm: continue rayNorm[0] += chkNorm[0] rayNorm[1] += chkNorm[1] rayNorm[2] += chkNorm[2] normCnt++ scale = 1.0 / normCnt; rayNorm = ~(rayNorm.Scale(scale)) # average the normal(s) and normalize #------------------------------------------------------------------------------------ # set up a few more variables before looping... #------------------------------------------------------------------------------------ matchdist = 100000 # signed distance to intersection point (just initialized to a large value) hitPolyNdx = -1 Weights = [0.0, 0.0, 0.0] for ci in range(numCheckVerts): cvNdx = chkVertList[ci] #if( cvNdx < 0 || cvNdx >= numSourceVerts ) #{ # char erbuf[256]; # sprintf(erbuf, "cvNdx [%d of %d] is out of range", cvNdx, numSourceVerts); # callError(funcname, erbuf); #} chkTriPolyIndices = srcVertTriPolys[cvNdx] # list of tripolys associated with this vertex if not chkTriPolyIndices: continue polyCount = len(chkTriPolyIndices) if not polyCount: continue for pi in range(polyCount): pNdx = chkTriPolyIndices[pi] checkTriPoly = sourceTriPolys[pNdx] # get TriPoly polyNorm = checkTriPoly.normal ndota = vector_dotproduct(polyNorm, rayNorm) if ndota <= 0.0: continue # if angle reversed, move on #================================================================================= # If we got past the above test, then we need to see if the ray # actually intersects the tri-poly.. #================================================================================= # determine ray intersection distance with plane of tripoly ndotv = vector_dotproduct(polyNorm, rayVert) isect = (checkTriPoly.plane - ndotv) / ndota; # create the intersection point testPoint = [rayVert[0] + (rayNorm[0] * isect), rayVert[1] + (rayNorm[1] * isect), rayVert[2] + (rayNorm[2] * isect)] #ray_segment(rayVert, rayNorm, isect) v0 = sourceVertList[checkTriPoly.v0] v1 = sourceVertList[checkTriPoly.v1] v2 = sourceVertList[checkTriPoly.v2] #------ if( !VECTOR_point_in_triangle(testPoint, v0, v1, v2) ) continue; ----------------- eu = [v1[0]-v0[0], v1[1]-v0[1], v1[2]-v0[2]] #vecsub(v1,v0) ev = [v2[0]-v0[0], v2[1]-v0[1], v2[2]-v0[2]] #vecsub(v2,v0) uu = eu[0] * eu[0] + eu[1] * eu[1] + eu[2] * eu[2] #vector_dotproduct(eu,eu) vv = ev[0] * ev[0] + ev[1] * ev[1] + ev[2] * ev[2] #vector_dotproduct(ev,ev) uv = eu[0] * ev[0] + eu[1] * ev[1] + eu[2] * ev[2] #vector_dotproduct(eu,ev) wx = [testPoint[0]-v0[0], testPoint[1]-v0[1], testPoint[2]-v0[2]] #vecsub(point,v0) wu = wx[0] * eu[0] + wx[1] * eu[1] + wx[2] * eu[2] #vector_dotproduct(wx,eu) wv = wx[0] * ev[0] + wx[1] * ev[1] + wx[2] * ev[2] #vector_dotproduct(wx,ev) D = (uv * uv) - (uu * vv) # Get and test parametric coords s = (uv * wv - vv * wu) / D if (s < -0.00000022) or (s-1.0) > 0.00000006: continue # point is outside T t = ((uv * wu) - (uu * wv)) / D if (t < -0.00000022) or ((s+t) - 1.0 > 0.00000006): continue # point is outside T #-------------------------------------------------------------------------------------------- #----------------------------------------------------------------------------------- # if we got here, then the intersection point is in fact within the triangle so # check the distance to see if it's closer than earlier hits, but keep looping # in case there's a better one. #----------------------------------------------------------------------------------- if abs(isect) < abs(matchdist): #================================================================================= # This code could be executed multiple times for each vertex/hit as better matches # are found, that's rarely more than a few matches for a relatively small number # of vertices. Since I've already got pointers to all the needed data handy, we # can afford to eat a few cycles now and save a few later. #================================================================================= matchdist = isect; # update matchdist for future tests hitDist = isect hitPolyNdx = pNdx; ixPoint = testPoint ixDelta = [ixPoint[0]-rayVert[0], ixPoint[1]-rayVert[1], ixPoint[2]-rayVert[2]] if hitPolyNdx != -1: #================================================================================= # if we got a hit, we compute the weights and break out of the loop #================================================================================= checkTriPoly = sourceTriPolys[hitPolyNdx] # need to get the right poly... polyNorm = checkTriPoly.normal; # it's normal... v0 = sourceVertList[checkTriPoly.v0] # and the verts for it again... v1 = sourceVertList[checkTriPoly.v1] v2 = sourceVertList[checkTriPoly.v2] Weights = get_raycast_weights(ixPoint, polyNorm, v0, v1, v2) # get the weights... return Weights, v0, v1, v2; # ???? return Weights, 0, 0, 0; # ????