Coverage for /usr/local/lib/python3.11/dist-packages/pyrocko/himesh.py: 0%
67 statements
« prev ^ index » next coverage.py v6.5.0, created at 2023-10-06 06:59 +0000
« prev ^ index » next coverage.py v6.5.0, created at 2023-10-06 06:59 +0000
1# https://pyrocko.org - GPLv3
2#
3# The Pyrocko Developers, 21st Century
4# ---|P------/S----------~Lg----------
6'''
7Hierarchical mesh with points-to-faces support, based on
8:py:mod:`pyrocko.icosphere`.
9'''
11import numpy as num
13from pyrocko.guts import Object, Int
15from pyrocko import icosphere, geometry as g
17guts_prefix = 'sparrow'
19and_ = num.logical_and
22class Level(object):
23 def __init__(self, vertices, faces):
24 self.vertices = vertices
25 self.faces = faces
26 self.next = None
27 self.nvertices = vertices.shape[0]
28 self.nfaces = faces.shape[0]
29 self.prepare_projections()
31 def points_in_face(self, iface, points):
32 uv = g.triangle_barycentric_coordinates(
33 self.bary_trans[iface, :, :],
34 self.vertices[self.faces[iface, 0]],
35 points)
37 return and_(
38 and_(uv[:, 0] >= 0., uv[:, 1] >= 0.), uv[:, 0] + uv[:, 1] <= 1.)
40 def points_to_face(self, points, icandidates=None):
42 if icandidates is None:
43 icandidates = num.arange(self.nfaces)
45 npoints = points.shape[0]
46 ifaces = num.zeros(npoints, dtype=int)
47 points_flat_all = num.zeros((npoints, 3))
48 for icandidate in icandidates:
49 points_flat, mask = self.project(icandidate, points)
50 ipoints = num.where(num.logical_and(
51 self.points_in_face(icandidate, points_flat),
52 mask))[0]
54 if self.next:
55 icandidates_child = num.where(
56 self.next.parents == icandidate)[0]
57 ifaces[ipoints], points_flat_all[ipoints, :] \
58 = self.next.points_to_face(
59 points[ipoints], icandidates_child)
61 else:
62 ifaces[ipoints] = icandidate
63 points_flat_all[ipoints, :] = points_flat[ipoints, :]
65 return ifaces, points_flat_all
67 def project(self, iface, points):
68 denom = g.vdot(points, self.aplanes[iface, :])
69 mask = denom > 0.0
70 points = points.copy()
71 points[mask, :] /= denom[mask, num.newaxis]
72 points[num.logical_not(mask), :] = 0.0
73 return points, mask
75 def prepare_projections(self):
76 planes = g.triangle_planes(
77 self.vertices, self.faces)
79 planes /= g.vdot(
80 planes,
81 self.vertices[self.faces[:, 0]])[:, num.newaxis]
83 self.aplanes = planes
85 self.bary_trans = g.triangle_barycentric_transforms(
86 self.vertices, self.faces)
89class HiMesh(Object):
91 order = Int.T(default=0)
93 def __init__(self, **kwargs):
94 Object.__init__(self, **kwargs)
95 self.levels = self.make_levels()
97 def make_levels(self):
98 # root = None
99 level = None
100 levels = []
101 for vertices, faces in icosphere.iter_icospheres(
102 self.order, inflate=False):
104 new = Level(vertices, faces)
105 if not levels:
106 new.parents = None
107 else:
108 new.parents = g.refine_triangle_parents(faces.shape[0])
109 level.next = new
111 level = new
112 levels.append(new)
114 return levels
116 def get_vertices(self):
117 return g.normalize(self.levels[-1].vertices)
119 def get_faces(self):
120 return self.levels[-1].faces
122 def points_to_faces(self, points):
123 return self.levels[0].points_to_face(points)[0]