I'm learning too by answering your questions. 2010/3/2 Ian Mallett <geometrian@gmail.com>:
1 v_array #array of beautifully transformed vertices from last step 2 n_array #array of beautifully transformed normals from last step 3 f_list #list of vec4s, where every vec4 is three vertex indices and a 4 #normal index. These indices together describe a triangle-- 5 #three vertex indices into "v_array", and one normal from "n_array". 6 edges = [] 7 for v1index,v2index,v3index,nindex in f_list: 8 v1,v2,v3 = [v_array[i] for i in [vi1index,v2index,v3index]] 9 if angle_between(n_array[nindex],v1-a_constant_point)<90deg: 10 for edge in [[v1index,v2index],[v2index,v3index],[v3index,v1index]]: 11 #add "edge" to "edges" 12 #remove pairs of identical edges (also note that things like 13 #[831,326] are the same as [326,831]) 14 edges2 = [] 15 for edge in edges: 16 edges2.append(v_array[edge[0]],v_array[edge[1]]) 17 return edges2
The loop I can replace by numpy operations:
v_array array([[1, 2, 3], [4, 5, 6], [7, 8, 9]]) n_array array([[ 0.1, 0.2, 0.3], [ 0.4, 0.5, 0.6]]) f_list array([[0, 1, 2, 0], [2, 1, 0, 1]])
Retrieving the v1 vectors:
v1s = v_array[f_list[:, 0]] v1s array([[1, 2, 3], [7, 8, 9]])
Retrieving the normal vectors:
ns = n_array[f_list[:, 3]] ns array([[ 0.1, 0.2, 0.3], [ 0.4, 0.5, 0.6]])
Now how to calculate the pairwise dot product (I suppress the difference of v1 to some_point for now):
inner = numpy.inner(ns, v1s) inner array([[ 1.4, 5. ], [ 3.2, 12.2]])
This calculates *all* pairwise dot products, we have to select the diagonal of this square ndarray:
dotprods = inner[[numpy.arange(0, 2), numpy.arange(0, 2)]] dotprods array([ 1.4, 12.2])
Now we can create a boolean array saying where the dotprod is > 0 (i.e, angle < 90°), and select those triangles:
select = dotprods > 0 select array([ True, True], dtype=bool) selected = f_list[select] selected array([[0, 1, 2, 0], [2, 1, 0, 1]])
In this case it's the full list. Now build the triangles corners array:
corners = v_array[selected[:, :3]] corners array([[[1, 2, 3], [4, 5, 6], [7, 8, 9]],
[[7, 8, 9], [4, 5, 6], [1, 2, 3]]])
This has indices [triangle, vertex number (0, 1, or 2), xyz]. And compute the edges (I think you can make use of them):
edges_dist = numpy.asarray([corners[:, 1] - corners[:, 0], corners[:, 2] - corners[:, 0], corners[:, 2] - corners[:, 1]]) edges_dist array([[[ 3, 3, 3], [-3, -3, -3]],
[[ 6, 6, 6], [-6, -6, -6]], [[ 3, 3, 3], [-3, -3, -3]]]) This has indices [corner number, triangle, xyz]. I think it's easier to compare then "reversed" edges, because then edge[i, j] == -edge[k, l]? But of course:
edges = numpy.asarray([[corners[:, 0], corners[:, 1]], [corners[:, 1], corners[:, 2]], [corners[:, 2], corners[:, 0]]]) edges array([[[[1, 2, 3], [7, 8, 9]],
[[4, 5, 6], [4, 5, 6]]], [[[4, 5, 6], [4, 5, 6]], [[7, 8, 9], [1, 2, 3]]], [[[7, 8, 9], [1, 2, 3]], [[1, 2, 3], [7, 8, 9]]]]) This has indices [edge number (0, 1, or 2), corner number in edge (0 or 1), triangle]. But this may not be what you want (not flattened in triangle number). Therefore:
edges0 = numpy.asarray([corners[:, 0], corners[:, 1]]) edges1 = numpy.asarray([corners[:, 1], corners[:, 2]]) edges2 = numpy.asarray([corners[:, 2], corners[:, 0]]) edges0 array([[[1, 2, 3], [7, 8, 9]],
[[4, 5, 6], [4, 5, 6]]])
edges1 array([[[4, 5, 6], [4, 5, 6]],
[[7, 8, 9], [1, 2, 3]]])
edges2 array([[[7, 8, 9], [1, 2, 3]],
[[1, 2, 3], [7, 8, 9]]])
edges = numpy.concatenate((edges0, edges1, edges2), axis = 0) edges array([[[1, 2, 3], [7, 8, 9]],
[[4, 5, 6], [4, 5, 6]], [[4, 5, 6], [4, 5, 6]], [[7, 8, 9], [1, 2, 3]], [[7, 8, 9], [1, 2, 3]], [[1, 2, 3], [7, 8, 9]]]) This should be as intended. The indices are [flat edge number, edge side (left or right), xyz]. Now I guess you have to iterate over all pairs of them, don't know a numpy accelerated method. Maybe it's even faster to draw the edges twice than to invest O(N_edges ** 2) complexity for comparing?
I may be able to do something with lines 14-17 myself--but I'm not sure. Ok, than that's fine, and let us all know about your solution.
It may seem a bit complicated, but I hope this impression is mainly because of the many outputs ... I double-checked everything, *hope* everything is correct. So far from me, Friedrich