Commit 40cf3913 authored by Daniel Boeckenhoff's avatar Daniel Boeckenhoff
Browse files

test_mesh3D works

parent b00afb4b
This source diff could not be displayed because it is too large. You can view the blob instead.
......@@ -21,7 +21,7 @@ class Sphere_Test(Base_Check, unittest.TestCase):
def test_cut_split(self):
x, y, z = sympy.symbols('x y z')
halfed = self._inst.cut(x > 0, at_intersection='split')
halfed = self._inst.cut(x + 1./100*y > 0, at_intersection='split')
if __name__ == '__main__':
unittest.main()
......@@ -1515,9 +1515,12 @@ class TensorFields(Tensors):
field = self.fields[field_index].copy()
if self.dim == field.dim:
field.transform(self.coord_sys)
if field.dim == 1:
artist = super(TensorFields, self).plot(c=field, **kwargs)
elif field.dim in (2, 3):
else:
import logging
log = logging.getLogger()
log.debug("Careful: Plotting tensors with field of"
"different dimension. No coord_sys check performed.")
if field.dim <= 3:
artist = tfields.plotting.plot_tensor_field(self, field,
**kwargs)
else:
......
......@@ -22,24 +22,30 @@ def _dist_from_plane(point, plane):
def _segment_plane_intersection(p0, p1, plane):
"""
Get the intersection between the ray spanned by p0 and p1 with the plane.
Args:
p0: array of length 3
p1: array of length 3
plane:
Returns:
points, direction
points (list of arrays of length 3): 3d points
"""
distance0 = _dist_from_plane(p0, plane)
distance1 = _dist_from_plane(p1, plane)
p0OnPlane = abs(distance0) < np.finfo(float).eps
p1OnPlane = abs(distance1) < np.finfo(float).eps
p0_on_plane = abs(distance0) < np.finfo(float).eps
p1_on_plane = abs(distance1) < np.finfo(float).eps
points = []
direction = 0
if p0OnPlane:
if p0_on_plane:
points.append(p0)
if p1OnPlane:
if p1_on_plane:
points.append(p1)
# remove duplicate points
if len(points) > 1:
points = np.unique(points, axis=0)
if p0OnPlane and p1OnPlane:
if p0_on_plane and p1_on_plane:
return points, direction
if distance0 * distance1 > np.finfo(float).eps:
......@@ -77,8 +83,8 @@ def _intersect(triangle, plane, vertices_rejected):
TODO:
align norm vectors with previous face
"""
nTrue = vertices_rejected.count(True)
lonely_bool = True if nTrue == 1 else False
n_true = vertices_rejected.count(True)
lonely_bool = True if n_true == 1 else False
index = vertices_rejected.index(lonely_bool)
s0, d0 = _segment_plane_intersection(triangle[0], triangle[1], plane)
s1, d1 = _segment_plane_intersection(triangle[1], triangle[2], plane)
......@@ -91,10 +97,6 @@ def _intersect(triangle, plane, vertices_rejected):
# TODO handle special cases. For now triangles with at least two points on plane are excluded
new_points = None
print(triangle, plane, vertices_rejected)
print(lonely_bool, index)
print(s0, s1, s2)
print(d0, d1, d2)
if len(s0) == 2:
# both points on plane
return new_points
......@@ -712,15 +714,14 @@ class Mesh3D(tfields.TensorMaps):
face_delete_indices.add(i)
if any(vertices_rejected) and not all(vertices_rejected):
# face on edge
nTrue = vertices_rejected.count(True)
lonely_bool = True if nTrue == 1 else False
n_true = vertices_rejected.count(True)
lonely_bool = True if n_true == 1 else False
triangle_points = [vertices[f] for f in face]
"""
Add the intersection points and faces
"""
intersection = _intersect(triangle_points, plane, vertices_rejected)
print(intersection)
last_idx = len(vertices) - 1
for tri_list in intersection:
new_face = []
......
......@@ -166,7 +166,6 @@ class PlotOptions(object):
vmin=vmin,
vmax=vmax,
cmap=cmap)
print(colors)
elif fmt == 'hex':
colors = [mpl.colors.to_hex(color) for color in colors]
else:
......
......@@ -326,28 +326,44 @@ def plot_mesh(vertices, faces, **kwargs):
return artist
def plot_tensor_field(points, vectors, **kwargs):
def plot_tensor_field(points, field, **kwargs):
"""
Args:
points (array_like): base vectors
vectors (array_like): direction vectors
field (): direction field
"""
po = tfields.plotting.PlotOptions(kwargs)
if points is None:
points = np.full(vectors.shape, 0.)
artists = []
xAxis, yAxis, zAxis = po.getXYZAxis()
for point, vector in zip(points, vectors):
if po.dim == 3:
artists.append(po.axis.quiver(point[xAxis], point[yAxis], point[zAxis],
vector[xAxis], vector[yAxis], vector[zAxis],
**po.plot_kwargs))
elif po.dim == 2:
artists.append(po.axis.quiver(point[xAxis], point[yAxis],
vector[xAxis], vector[yAxis],
**po.plot_kwargs))
else:
raise NotImplementedError("Dimension != 2|3")
field = np.array(field)
if len(field.shape) == 2 and field.shape[1] == 1:
# scalar
field = field.reshape(len(points))
if len(field.shape) == 1:
# scalar
colors = po.format_colors(field)
po.delNormArgs()
artists = plot_array(points, c=colors, **po.plot_kwargs)
artists.set_array(field)
elif len(field.shape) == 2:
# vector
if points is None:
points = np.full(field.shape, 0.)
xAxis, yAxis, zAxis = po.getXYZAxis()
artists = []
for point, vector in zip(points, field):
if po.dim == 3:
artists.append(
po.axis.quiver(
point[xAxis], point[yAxis], point[zAxis],
vector[xAxis], vector[yAxis], vector[zAxis],
**po.plot_kwargs))
elif po.dim == 2:
artists.append(po.axis.quiver(point[xAxis], point[yAxis],
vector[xAxis], vector[yAxis],
**po.plot_kwargs))
else:
raise NotImplementedError("Dimension != 2|3")
else:
raise NotImplementedError("Only Scalars and Vectors implemented")
return artists
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment