diff --git a/docs/examples/01_aperture_projection.py b/docs/examples/01_aperture_projection.py index 561c512..22429bf 100644 --- a/docs/examples/01_aperture_projection.py +++ b/docs/examples/01_aperture_projection.py @@ -53,7 +53,11 @@ body, array, _ = data.exampleUAV(10e9) - +import pyvista as pv +pl=pv.Plotter() +pl.add_mesh(pv.from_meshio(body),color="green") +pl.add_mesh(pv.from_meshio(array),color="aqua") +pl.show() # %% ## .. image:: ../_static/open3d_structure.png @@ -64,6 +68,10 @@ surface_array.cells[0].data = np.asarray(array.cells[0].data)[: (array.cells[0].data).shape[0] // 2, :] surface_array.cell_data["Normals"] = np.array(array.cell_data["Normals"])[: (array.cells[0].data).shape[0] // 2] +#Recalculate Normals +from lyceanem.geometry.geometryfunctions import compute_normals +surface_array=compute_normals(surface_array) + # %% # Structures diff --git a/docs/examples/02_coherently_polarised_array.py b/docs/examples/02_coherently_polarised_array.py index adad35a..5579152 100644 --- a/docs/examples/02_coherently_polarised_array.py +++ b/docs/examples/02_coherently_polarised_array.py @@ -61,6 +61,11 @@ import pyvista as pv +pl=pv.Plotter() +pl.add_mesh(pv.from_meshio(body),color="green") +pl.add_mesh(pv.from_meshio(array),color="aqua") +pl.add_axes() +pl.show() source_points = surface_array.points @@ -106,9 +111,9 @@ azimuth_resolution=az_res, elevation_resolution=elev_res ) UAV_Static_Pattern.pattern[:, :, 0] = Etheta -UAV_Static_Pattern.pattern[:, :, 0] = Ephi +UAV_Static_Pattern.pattern[:, :, 1] = Ephi -UAV_Static_Pattern.display_pattern() +UAV_Static_Pattern.display_pattern(desired_pattern='Power') # %% # .. image:: ../_static/sphx_glr_02_coherently_polarised_array_001.png diff --git a/lyceanem/geometry/geometryfunctions.py b/lyceanem/geometry/geometryfunctions.py index e83b6f6..a6380b6 100644 --- a/lyceanem/geometry/geometryfunctions.py +++ b/lyceanem/geometry/geometryfunctions.py @@ -11,25 +11,6 @@ from ..electromagnetics.emfunctions import transform_em -def tri_areas(triangle_mesh): - """ - Calculate the area of each triangle in a triangle mesh. - - Args: - triangle_mesh (np.ndarray): A (N, 3, 3) array of N triangles, each with 3 vertices. - - Returns: - np.ndarray: A (N,) array of the area of each triangle. - """ - assert ( - triangle_mesh.cells[0].type == "triangle" - ), "Only triangle meshes are supported." - triangle_indices = triangle_mesh.cells[0].data - v0 = triangle_mesh.points[triangle_indices[:, 0]] - v1 = triangle_mesh.points[triangle_indices[:, 1]] - v2 = triangle_mesh.points[triangle_indices[:, 2]] - return 0.5 * np.linalg.norm(np.cross(v1 - v0, v2 - v0), axis=1) - def cell_centroids(field_data): """ @@ -177,7 +158,13 @@ def mesh_transform(mesh, transform_matrix, rotate_only): return return_mesh - +def locate_cell_index(field_data,cell_type='triangle'): + for inc, cell in enumerate(field_data.cells): + if cell.type==cell_type: + desired_index=inc + + return desired_index + def compute_areas(field_data): cell_areas = [] for inc, cell in enumerate(field_data.cells): diff --git a/lyceanem/models/frequency_domain.py b/lyceanem/models/frequency_domain.py index 44c70fb..1005cdd 100644 --- a/lyceanem/models/frequency_domain.py +++ b/lyceanem/models/frequency_domain.py @@ -52,14 +52,24 @@ def aperture_projection( directivity_envelope = np.zeros( (elev_range.shape[0], az_range.shape[0]), dtype=np.float32 ) - triangle_centroids, _ = GF.tri_centroids(aperture) - triangle_areas = GF.tri_areas(aperture) - triangle_normals = np.asarray(aperture.cell_data["Normals"]) + triangle_centroids = GF.cell_centroids(aperture) + aperture = GF.compute_areas(aperture) + aperture=GF.compute_normals(aperture) + triangle_cell_index=GF.locate_cell_index(aperture) + triangle_normals = aperture.cell_data["Normals"][triangle_cell_index] + # Ensure no clash with triangles in raycaster + triangle_centroids.points+=1e-6*triangle_normals + import pyvista as pv + pl=pv.Plotter() + pl.add_mesh(pv.from_meshio(aperture),scalars="Area") + pl.add_mesh(pv.from_meshio(environment.solids[0]),color="red") + pl.add_mesh(pv.from_meshio(triangle_centroids),color="green") + pl.show() visible_patterns, pcd = RF.visiblespace( triangle_centroids, triangle_normals, blocking_triangles, - vertex_area=triangle_areas, + vertex_area=aperture.point_data['Area'], az_range=az_range, elev_range=elev_range, shell_range=farfield_distance, diff --git a/lyceanem/raycasting/rayfunctions.py b/lyceanem/raycasting/rayfunctions.py index 5a78a84..4096948 100644 --- a/lyceanem/raycasting/rayfunctions.py +++ b/lyceanem/raycasting/rayfunctions.py @@ -348,7 +348,7 @@ def visiblespace( """ azaz, elel = np.meshgrid(az_range, elev_range) - sourcenum = len(source_coords) + sourcenum = source_coords.points.shape[0] sinks = np.zeros((len(np.ravel(azaz)), 3), dtype=np.float32) sinks[:, 0], sinks[:, 1], sinks[:, 2] = azeltocart( @@ -363,10 +363,10 @@ def visiblespace( # need to create farfield sinks in az,elev coordinates, then convert to xyz sink coordinates, and generate index # missed_points,hit_points,missed_index,hit_index,shadow_rays=chunkingRaycaster1D(source_coords,sinks,np.zeros((1,3),dtype=np.float32),initial_index,environment,1,terminate_flag=True) hit_index, _ = workchunkingv2( - source_coords, sinks, np.empty((0, 3), dtype=np.float32), environment, 1 + source_coords.points, sinks, np.empty((0, 3), dtype=np.float32), environment, 1 ) unified_model = np.append( - source_coords.astype(np.float32), sinks.astype(np.float32), axis=0 + source_coords.points.astype(np.float32), sinks.astype(np.float32), axis=0 ) # filtered_network2,final_network2,filtered_index2,final_index2,shadow_rays directions = np.zeros((len(hit_index), 3), dtype=np.float32) @@ -386,6 +386,7 @@ def visiblespace( ) # angles[np.isnan(angles)]=0 + vertex_area[np.isnan(vertex_area)]=0 # visible_patterns=quickpatterncreator(az_range,elev_range,source_coords,angles,vertex_area,hit_index) if len(vertex_area) == 1: if vertex_area == 0: diff --git a/lyceanem/tests/reflectordata.py b/lyceanem/tests/reflectordata.py index 3b0aa1e..7ebe156 100644 --- a/lyceanem/tests/reflectordata.py +++ b/lyceanem/tests/reflectordata.py @@ -54,22 +54,20 @@ def exampleUAV(frequency): array, np.array([0.25, 0, 0]) + np.array([-0.18, 0, 0.0125]) ) - def structure_cells(array): + #def structure_cells(array): ## add collumn of 3s to beggining of each row - array = np.append( - np.ones((array.shape[0], 1), dtype=np.int32) * 3, array, axis=1 - ) - return array - - pyvista_array = pv.PolyData(array.points, structure_cells(array.cells[0].data)) - pyvista_body = pv.PolyData(body.points, structure_cells(body.cells[0].data)) - pyvista_array.compute_normals(inplace=True) - pyvista_body.compute_normals(inplace=True) - - array.point_data["Normals"] = pyvista_array.point_normals - body.point_data["Normals"] = pyvista_body.point_normals - array.cell_data["Normals"] = pyvista_array.cell_normals - body.cell_data["Normals"] = pyvista_body.cell_normals + # array = np.append( + # np.ones((array.shape[0], 1), dtype=np.int32) * 3, array, axis=1 + # ) + # return array + + #pyvista_array = pv.PolyData(array.points, structure_cells(array.cells[0].data)) + #pyvista_body = pv.PolyData(body.points, structure_cells(body.cells[0].data)) + #pyvista_array.compute_normals(inplace=True) + #pyvista_body.compute_normals(inplace=True) + + array = GF.compute_normals(array) + body = GF.compute_normals(body) wavelength = 3e8 / frequency mesh_sep = wavelength * 0.5