Changing "normals" to "Normals" as this is vtk convention. Also adapt…

…ed convert triangles to find the triangle cells in a multicell mesh.

docs/examples/07_aperture_farfield_polarisation.py

@@ -11,6 +11,7 @@
 
 """
 import numpy as np
+import pygmsh
 
 # %%
 # Setting Farfield Resolution and Wavelength
@@ -34,9 +35,45 @@
 from lyceanem.base_classes import points,structures,antenna_structures
 
 from lyceanem.geometry.targets import meshedHorn
-structure,array_points=meshedHorn(3*wavelength, 1*wavelength, 4*wavelength, 0.05*wavelength,np.radians(10),wavelength*0.5)
-
-horn_antenna=antenna_structures(structures(solids=[structure]), points(points=[array_points]))
+structure,array_points=meshedHorn(3*wavelength, 1*wavelength, 4*wavelength, 1*wavelength,np.radians(10),wavelength*0.5)
+
+def planar_aperture(u_size,v_size,wall_thickness=2e-3,depth=2e-3,mesh_size=1.0):
+    with pygmsh.occ.Geometry() as geom:
+        geom.characteristic_length_max=mesh_size
+        u_points=np.linspace(-u_size/2,u_size/2,np.ceil(u_size/mesh_size).astype(int))
+        v_points=np.linspace(-v_size/2,v_size/2,np.ceil(v_size/mesh_size).astype(int))
+        n_points=np.linspace(0,0,1)
+        u_mesh,v_mesh,n_mesh=np.meshgrid(u_points,v_points,n_points,indexing='ij')
+        aperture_points=np.array([u_mesh.ravel(),v_mesh.ravel(),n_mesh.ravel()]).transpose()
+        source_normals = np.empty(((aperture_points.shape[0], 3)), dtype=np.float32)
+        source_normals[:] = 0.0
+        source_normals[:, 2] = 1.0
+        ground_plane=geom.add_box((-(u_size+wall_thickness*2)/2.0,-(v_size+wall_thickness*2)/2.0,-depth), extents=(u_size+wall_thickness*2,v_size+wall_thickness*2,depth),mesh_size=mesh_size)
+        point_list=[]
+        for inc in range(aperture_points.shape[0]):
+            point_list.append(geom.add_point(aperture_points[inc,:]))
+            geom.in_surface(point_list[-1],ground_plane)
+
+        #aperture=geom.add_rectangle((-u_size/2.0,-v_size/2.0,0.0), u_size, v_size,mesh_size=mesh_size)
+        #ground_plane=geom.add_box((-(u_size+wall_thickness*2)/2.0,-(v_size+wall_thickness*2)/2.0,-depth), extents=(u_size+wall_thickness*2,v_size+wall_thickness*2,depth),mesh_size=mesh_size)
+
+
+        mesh = geom.generate_mesh()
+        end_point=mesh.points.shape[0]
+        mesh.point_data['Normals']=np.zeros((mesh.points.shape[0],3))
+        mesh.point_data['Normals'][mesh.points[:,2]>-depth*0.5,2]=1.0
+        mesh.point_data['Normals'][mesh.points[:,2]>-depth*0.5,2]=-1.0
+        mesh.points=np.append(mesh.points,aperture_points,axis=0)
+        mesh.point_data['Normals']=np.append(mesh.point_data['Normals'],source_normals,axis=0)
+        mesh.cell_sets['Aperture']=[None,None,None,np.linspace(end_point,aperture_points.shape[0]+end_point-1,aperture_points.shape[0]).astype(int)]
+        mesh.cells[3].data=np.append(mesh.cells[3].data,np.linspace(end_point,aperture_points.shape[0]+end_point-1,aperture_points.shape[0]).astype(int).reshape(-1,1),axis=0)
+    return mesh
+
+mesh2=planar_aperture(3*wavelength,1*wavelength,mesh_size=wavelength*0.5)
+
+
+
+horn_antenna=antenna_structures(structures(solids=[mesh2]), points(points=[array_points]))
 
 
 from lyceanem.models.frequency_domain import calculate_farfield

lyceanem/base_classes.py

@@ -390,7 +390,7 @@ def excitation_function(
         else:
             aperture_points = self.export_all_points(point_index=point_index)
         aperture_weights = EM.calculate_conformalVectors(
-            desired_e_vector, aperture_points.point_data['normals'], self.pose[:3, :3]
+            desired_e_vector, aperture_points.point_data['Normals'], self.pose[:3, :3]
         )
         if phase_shift == "wavefront":
             source_points = np.asarray(aperture_points.points)

lyceanem/geometry/geometryfunctions.py

@@ -82,12 +82,12 @@ def mesh_transform(mesh, transform_matrix, rotate_only):
     if rotate_only:
         for i in range(mesh.points.shape[0]):
             return_mesh.points[i] = np.dot(transform_matrix, np.append(mesh.points[i], 0))[:3]
-            return_mesh.point_data['normals'][i] = np.dot(transform_matrix, np.append(mesh.point_data['normals'][i], 0))[:3]
+            return_mesh.point_data['Normals'][i] = np.dot(transform_matrix, np.append(mesh.point_data['Normals'][i], 0))[:3]
 
     else:
         for i in range(mesh.points.shape[0]):
             return_mesh.points[i] = np.dot(transform_matrix, np.append(mesh.points[i], 1))[:3]
-            return_mesh.point_data['normals'][i]= np.dot(transform_matrix, np.append(mesh.point_data['normals'][i], 0))[:3]
+            return_mesh.point_data['Normals'][i]= np.dot(transform_matrix, np.append(mesh.point_data['Normals'][i], 0))[:3]
 
 
     return return_mesh

lyceanem/geometry/targets.py

@@ -125,8 +125,8 @@ def rectReflector(majorsize, minorsize, thickness):
 
 
 
-    mesh.point_data["normals"] = pv_mesh.point_normals
-    mesh.cell_data["normals"] = pv_mesh.cell_normals
+    mesh.point_data["Normals"] = pv_mesh.point_normals
+    mesh.cell_data["Normals"] = pv_mesh.cell_normals
     red = np.zeros((8, 1), dtype=np.float32) 
     green = np.ones((8, 1), dtype=np.float32) * 0.259
     blue = np.ones((8, 1), dtype=np.float32) * 0.145
@@ -199,8 +199,8 @@ def shapeTrapezoid(x_size, y_size, length, flare_angle):
     pv_mesh = pv.PolyData( mesh_vertices, faces = triangle_list)
     pv_mesh.compute_normals(inplace=True,consistent_normals=False)
 
-    mesh.point_data["normals"] = np.asarray(pv_mesh.point_normals)
-    mesh.cell_data["normals"] = np.asarray(pv_mesh.cell_normals)
+    mesh.point_data["Normals"] = np.asarray(pv_mesh.point_normals)
+    mesh.cell_data["Normals"] = np.asarray(pv_mesh.cell_normals)
 
     red = np.zeros((8, 1), dtype=np.float32) 
     green = np.ones((8, 1), dtype=np.float32) * 0.259
@@ -476,6 +476,6 @@ def gridedReflectorPoints(
         mesh_vertices = source_coords
         mesh_normals = source_normals
 
-    mesh_points = meshio.Mesh(points=mesh_vertices, cells=[], point_data={"normals": mesh_normals})
+    mesh_points = meshio.Mesh(points=mesh_vertices, cells=[], point_data={"Normals": mesh_normals})
 
     return mesh_points

lyceanem/models/frequency_domain.py

@@ -54,7 +54,7 @@ def aperture_projection(
     )
     triangle_centroids, _ = GF.tri_centroids(aperture)
     triangle_areas = GF.tri_areas(aperture)
-    triangle_normals = np.asarray(aperture.cell_data["normals"])
+    triangle_normals = np.asarray(aperture.cell_data["Normals"])
     visible_patterns, pcd = RF.visiblespace(
         triangle_centroids,
         triangle_normals,
@@ -143,17 +143,17 @@ def calculate_farfield(
         sinks, np.zeros((len(sinks), 3), dtype=np.float32), sink_normals, lengths
     )
     sink_cloud = RF.points2pointcloud(sinks)
-    sink_cloud.point_data["normals"] = sink_normals
+    sink_cloud.point_data["Normals"] = sink_normals
     num_sources = len(np.asarray(aperture_coords.points))
     num_sinks = len(np.asarray(sink_cloud.points))
     environment_triangles=GF.mesh_conversion(antenna_solid)
 
     if project_vectors:
         conformal_E_vectors = EM.calculate_conformalVectors(
-            desired_E_axis, np.asarray(aperture_coords.point_data["normals"]), antenna_axes
+            desired_E_axis, np.asarray(aperture_coords.point_data["Normals"]), antenna_axes
         )
     else:
-        if desired_E_axis.shape[0]==np.asarray(aperture_coords.point_data["normals"]).shape[0]:
+        if desired_E_axis.shape[0]==np.asarray(aperture_coords.point_data["Normals"]).shape[0]:
             conformal_E_vectors=copy.deepcopy(desired_E_axis)
         else:
             conformal_E_vectors = np.repeat(
@@ -169,8 +169,8 @@ def calculate_farfield(
             axis=0,
         )
         unified_normals = np.append(
-            np.asarray(aperture_coords.point_data["normals"]).astype(np.float32),
-            np.asarray(sink_cloud.point_data["normals"]).astype(np.float32),
+            np.asarray(aperture_coords.point_data["Normals"]).astype(np.float32),
+            np.asarray(sink_cloud.point_data["Normals"]).astype(np.float32),
             axis=0,
         )
         unified_weights = np.ones((unified_model.shape[0], 3), dtype=np.complex64)
@@ -199,7 +199,7 @@ def calculate_farfield(
         point_informationv2[0:num_sources]["pz"] = np.asarray(
             aperture_coords.points
         ).astype(np.float32)[:, 2]
-        normals = np.asarray(aperture_coords.point_data["normals"])
+        normals = np.asarray(aperture_coords.point_data["Normals"])
         point_informationv2[0:num_sources]["nx"] = normals[:, 0]
         point_informationv2[0:num_sources]["ny"] = normals[:, 1]
         point_informationv2[0:num_sources]["nz"] = normals[:, 2]
@@ -242,11 +242,11 @@ def calculate_farfield(
         )
         unified_normals = np.append(
             np.append(
-                np.asarray(aperture_coords.point_data["normals"]).astype(np.float32),
-                np.asarray(sink_cloud.point_data["normals"]).astype(np.float32),
+                np.asarray(aperture_coords.point_data["Normals"]).astype(np.float32),
+                np.asarray(sink_cloud.point_data["Normals"]).astype(np.float32),
                 axis=0,
             ),
-            np.asarray(scatter_points.point_data["normals"]).astype(np.float32),
+            np.asarray(scatter_points.point_data["Normals"]).astype(np.float32),
             axis=0,
         )
         unified_weights = np.ones((unified_model.shape[0], 3), dtype=np.complex64)
@@ -278,13 +278,13 @@ def calculate_farfield(
             aperture_coords.points
         ).astype(np.float32)[:, 2]
         point_informationv2[0:num_sources]["nx"] = np.asarray(
-            aperture_coords.point_data["normals"]
+            aperture_coords.point_data["Normals"]
         ).astype(np.float32)[:, 0]
         point_informationv2[0:num_sources]["ny"] = np.asarray(
-            aperture_coords.point_data["normals"]
+            aperture_coords.point_data["Normals"]
         ).astype(np.float32)[:, 1]
         point_informationv2[0:num_sources]["nz"] = np.asarray(
-            aperture_coords.point_data["normals"]
+            aperture_coords.point_data["Normals"]
         ).astype(np.float32)[:, 2]
         # point_informationv2[0:num_sources]['ex']=unified_weights[0:num_sources,0]
         # point_informationv2[0:num_sources]['ey']=unified_weights[0:num_sources,1]
@@ -314,7 +314,7 @@ def calculate_farfield(
         point_informationv2[(num_sources + num_sinks) :]["pz"] = np.asarray(
             scatter_points.points
         ).astype(np.float32)[:, 2]
-        scatter_normals = np.asarray(scatter_points.point_data["normals"])
+        scatter_normals = np.asarray(scatter_points.point_data["Normals"])
         point_informationv2[(num_sources + num_sinks) :]["nx"] = scatter_normals[:, 0]
         point_informationv2[(num_sources + num_sinks) :]["ny"] = scatter_normals[:, 1]
         point_informationv2[(num_sources + num_sinks) :]["nz"] = scatter_normals[:, 2]
@@ -486,11 +486,11 @@ def calculate_scattering(
     if not multiE:
         if project_vectors:
             conformal_E_vectors = EM.calculate_conformalVectors(
-                desired_E_axis, np.asarray(aperture_coords.point_data["normals"]), antenna_axes
+                desired_E_axis, np.asarray(aperture_coords.point_data["Normals"]), antenna_axes
             )
         else:
             print("hi from here", aperture_coords.cell_data)
-            if desired_E_axis.shape[0] == np.asarray(aperture_coords.point_data["normals"]).shape[0]:
+            if desired_E_axis.shape[0] == np.asarray(aperture_coords.point_data["Normals"]).shape[0]:
                 conformal_E_vectors = copy.deepcopy(desired_E_axis)
             else:
                 conformal_E_vectors = np.repeat(
@@ -499,10 +499,10 @@ def calculate_scattering(
     else:
         if project_vectors:
             conformal_E_vectors = EM.calculate_conformalVectors(
-                desired_E_axis, np.asarray(aperture_coords.point_data["normals"]), antenna_axes
+                desired_E_axis, np.asarray(aperture_coords.point_data["Normals"]), antenna_axes
             )
         else:
-            if desired_E_axis.shape[0] == np.asarray(aperture_coords.point_data["normals"]).shape[0]:
+            if desired_E_axis.shape[0] == np.asarray(aperture_coords.point_data["Normals"]).shape[0]:
                 conformal_E_vectors = copy.deepcopy(desired_E_axis)
             else:
                 conformal_E_vectors = np.repeat(
@@ -518,8 +518,8 @@ def calculate_scattering(
             axis=0,
         )
         unified_normals = np.append(
-            np.asarray(aperture_coords.point_data["normals"]).astype(np.float32),
-            np.asarray(sink_coords.point_data["normals"]).astype(np.float32),
+            np.asarray(aperture_coords.point_data["Normals"]).astype(np.float32),
+            np.asarray(sink_coords.point_data["Normals"]).astype(np.float32),
             axis=0,
         )
         unified_weights = np.ones((unified_model.shape[0], 3), dtype=np.complex64)
@@ -545,13 +545,13 @@ def calculate_scattering(
             aperture_coords.points
         ).astype(np.float32)[:, 2]
         point_informationv2[0:num_sources]["nx"] = np.asarray(
-            aperture_coords.point_data["normals"]
+            aperture_coords.point_data["Normals"]
         ).astype(np.float32)[:, 0]
         point_informationv2[0:num_sources]["ny"] = np.asarray(
-            aperture_coords.point_data["normals"]
+            aperture_coords.point_data["Normals"]
         ).astype(np.float32)[:, 1]
         point_informationv2[0:num_sources]["nz"] = np.asarray(
-            aperture_coords.point_data["normals"]
+            aperture_coords.point_data["Normals"]
         ).astype(np.float32)[:, 2]
         # set position and velocity of sinks
         point_informationv2[num_sources : (num_sources + num_sinks)]["px"] = np.asarray(
@@ -564,13 +564,13 @@ def calculate_scattering(
             sink_coords.points
         ).astype(np.float32)[:, 2]
         point_informationv2[num_sources : (num_sources + num_sinks)]["nx"] = np.asarray(
-            sink_coords.point_data["normals"]
+            sink_coords.point_data["Normals"]
         ).astype(np.float32)[:, 0]
         point_informationv2[num_sources : (num_sources + num_sinks)]["ny"] = np.asarray(
-            sink_coords.point_data["normals"]
+            sink_coords.point_data["Normals"]
         ).astype(np.float32)[:, 1]
         point_informationv2[num_sources : (num_sources + num_sinks)]["nz"] = np.asarray(
-            sink_coords.point_data["normals"]
+            sink_coords.point_data["Normals"]
         ).astype(np.float32)[:, 2]
 
         point_informationv2[:]["ex"] = unified_weights[:, 0]
@@ -591,7 +591,7 @@ def calculate_scattering(
             if project_vectors:
                 conformal_E_vectors = EM.calculate_conformalVectors(
                     desired_E_axis[0, :].reshape(1, 3),
-                    np.asarray(aperture_coords.point_data["normals"]).astype(np.float32),
+                    np.asarray(aperture_coords.point_data["Normals"]).astype(np.float32),
                 )
             else:
                 conformal_E_vectors = np.repeat(
@@ -603,7 +603,7 @@ def calculate_scattering(
             if project_vectors:
                 conformal_E_vectors = EM.calculate_conformalVectors(
                     desired_E_axis[0, :].reshape(1, 3),
-                    np.asarray(aperture_coords.point_data["normals"]).astype(np.float32),
+                    np.asarray(aperture_coords.point_data["Normals"]).astype(np.float32),
                 )
             else:
                 if desired_E_axis.size == 3:
@@ -626,11 +626,11 @@ def calculate_scattering(
         )
         unified_normals = np.append(
             np.append(
-                np.asarray(aperture_coords.point_data["normals"]).astype(np.float32),
-                np.asarray(sink_coords.point_data["normals"]).astype(np.float32),
+                np.asarray(aperture_coords.point_data["Normals"]).astype(np.float32),
+                np.asarray(sink_coords.point_data["Normals"]).astype(np.float32),
                 axis=0,
             ),
-            np.asarray(scatter_points.point_data["normals"]).astype(np.float32),
+            np.asarray(scatter_points.point_data["Normals"]).astype(np.float32),
             axis=0,
         )
         unified_weights = np.ones((unified_model.shape[0], 3), dtype=np.complex64)
@@ -659,13 +659,13 @@ def calculate_scattering(
             aperture_coords.points
         ).astype(np.float32)[:, 2]
         point_informationv2[0:num_sources]["nx"] = np.asarray(
-            aperture_coords.point_data["normals"]
+            aperture_coords.point_data["Normals"]
         ).astype(np.float32)[:, 0]
         point_informationv2[0:num_sources]["ny"] = np.asarray(
-            aperture_coords.point_data["normals"]
+            aperture_coords.point_data["Normals"]
         ).astype(np.float32)[:, 1]
         point_informationv2[0:num_sources]["nz"] = np.asarray(
-            aperture_coords.point_data["normals"]
+            aperture_coords.point_data["Normals"]
         ).astype(np.float32)[:, 2]
         # point_informationv2[0:num_sources]['ex']=unified_weights[0:num_sources,0]
         # point_informationv2[0:num_sources]['ey']=unified_weights[0:num_sources,1]
@@ -684,13 +684,13 @@ def calculate_scattering(
         # point_informationv2[num_sources:(num_sources+num_sinks)]['vy']=0.0
         # point_informationv2[num_sources:(num_sources+num_sinks)]['vz']=0.0
         point_informationv2[num_sources : (num_sources + num_sinks)]["nx"] = np.asarray(
-            sink_coords.point_data["normals"]
+            sink_coords.point_data["Normals"]
         ).astype(np.float32)[:, 0]
         point_informationv2[num_sources : (num_sources + num_sinks)]["ny"] = np.asarray(
-            sink_coords.point_data["normals"]
+            sink_coords.point_data["Normals"]
         ).astype(np.float32)[:, 1]
         point_informationv2[num_sources : (num_sources + num_sinks)]["nz"] = np.asarray(
-            sink_coords.point_data["normals"]
+            sink_coords.point_data["Normals"]
         ).astype(np.float32)[:, 2]
         point_informationv2[(num_sources + num_sinks) :]["px"] = np.asarray(
             scatter_points.points
@@ -702,13 +702,13 @@ def calculate_scattering(
             scatter_points.points
         ).astype(np.float32)[:, 2]
         point_informationv2[(num_sources + num_sinks) :]["nx"] = np.asarray(
-            scatter_points.point_data["normals"]
+            scatter_points.point_data["Normals"]
         ).astype(np.float32)[:, 0]
         point_informationv2[(num_sources + num_sinks) :]["ny"] = np.asarray(
-            scatter_points.point_data["normals"]
+            scatter_points.point_data["Normals"]
         ).astype(np.float32)[:, 1]
         point_informationv2[(num_sources + num_sinks) :]["nz"] = np.asarray(
-            scatter_points.point_data["normals"]
+            scatter_points.point_data["Normals"]
         ).astype(np.float32)[:, 2]
         point_informationv2[:]["ex"] = unified_weights[:, 0]
         point_informationv2[:]["ey"] = unified_weights[:, 1]
@@ -741,7 +741,7 @@ def calculate_scattering(
             for e_inc in range(desired_E_axis.shape[0]):
                 conformal_E_vectors = EM.calculate_conformalVectors(
                     desired_E_axis[e_inc, :],
-                    np.asarray(aperture_coords.point_data["normals"]).astype(np.float32),
+                    np.asarray(aperture_coords.point_data["Normals"]).astype(np.float32),
                 )
                 unified_weights[0:num_sources, :] = conformal_E_vectors# / num_sources
                 point_informationv2[:]["ex"] = unified_weights[:, 0]
@@ -772,7 +772,7 @@ def calculate_scattering(
             for e_inc in range(desired_E_axis.shape[1]):
                 conformal_E_vectors = EM.calculate_conformalVectors(
                     desired_E_axis[e_inc, :],
-                    np.asarray(aperture_coords.point_data["normals"]).astype(np.float32),
+                    np.asarray(aperture_coords.point_data["Normals"]).astype(np.float32),
                 )
                 for element in range(num_sources):
                     point_informationv2[0:num_sources]["ex"] = 0.0