.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_examples\08_array_validation.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code. .. rst-class:: sphx-glr-example-title .. _sphx_glr_auto_examples_08_array_validation.py: Antenna Array Validation ========================== This example uses some custom functions to import antenna array measurements from a 32 element antenna array panel, and simulate the same antenna array within LyceanEM. Some basic functions are included to import the hdf5 file, and demostrate a suitable method for generating an equivalent model antenna array using gmsh. The measurement file can be found at https://osf.io/64yaj/files/osfstorage/6867f249514c2235de49f20e , and should be placed in the same folder as this script. .. GENERATED FROM PYTHON SOURCE LINES 11-159 .. code-block:: Python import pyvista as pv import meshio import numpy as np from scipy.constants import speed_of_light def Beamform(Antennas, weights): """ Assuming all the sample fields in the Antennas list are the same points with only the element differing, the weights array will be used to calculate a resultant array pattern, """ field_shape = Antennas[0].point_data["Ex-Real"].shape beamformed_map = np.zeros((field_shape[0], 3), dtype=np.complex64) for element in range(len(Antennas)): beamformed_map[:, 0] += ( ( Antennas[element].point_data["Ex-Real"] + 1j * Antennas[element].point_data["Ex-Imag"] ) * weights[element].reshape(1, -1) ).ravel() beamformed_map[:, 1] += ( ( Antennas[element].point_data["Ey-Real"] + 1j * Antennas[element].point_data["Ey-Imag"] ) * weights[element].reshape(1, -1) ).ravel() beamformed_map[:, 2] += ( ( Antennas[element].point_data["Ez-Real"] + 1j * Antennas[element].point_data["Ez-Imag"] ) * weights[element].reshape(1, -1) ).ravel() beamformed_pattern = meshio.Mesh( points=Antennas[0].points, cells=[("triangle", Antennas[0].cells[0].data)], point_data={ "Normals": Antennas[0].point_data["Normals"], "Ex-Real": np.real(beamformed_map[:, 0]), "Ex-Imag": np.imag(beamformed_map[:, 0]), "Ey-Real": np.real(beamformed_map[:, 1]), "Ey-Imag": np.imag(beamformed_map[:, 1]), "Ez-Real": np.real(beamformed_map[:, 2]), "Ez-Imag": np.imag(beamformed_map[:, 2]), }, ) return beamformed_pattern def Panel_Array_data(): from pathlib import Path from lyceanem.geometry.geometryfunctions import compute_areas, compute_normals from lyceanem.electromagnetics.emfunctions import update_electric_fields dataset_location = Path("Panel_Array_Measurements.hdf5") import h5py with h5py.File(dataset_location, "r") as f1: # print(f1.name) # print(f1['Fields'].keys()) antennas = [] for field_inc in range(f1["Metadata"]["Element Positions"][:].shape[0]): pattern_keys = "Element {}".format(field_inc) temp_mesh = meshio.Mesh( points=f1["Sample Positions"]["Field Positions"][:, :], cells=[("triangle", f1["Sample Positions"]["Triangles"][:, :])], ) temp_exeyez = f1["Fields"][pattern_keys][:, 2:] temp_mesh = update_electric_fields( temp_mesh, temp_exeyez[:, 0], temp_exeyez[:, 1], temp_exeyez[:, 2] ) antennas.append(temp_mesh) frequency = f1["Metadata"]["Frequency (Hz)"][()] element_positions = f1["Metadata"]["Element Positions"][:, :] from lyceanem.electromagnetics.emfunctions import ( PoyntingVector, Directivity, field_magnitude_phase, ) new_list = [] for pattern in antennas: pattern = compute_areas(pattern) pattern.point_data["Normals"] = ( pattern.points - np.mean(pattern.points, axis=0) ) / np.linalg.norm( (pattern.points - np.mean(pattern.points, axis=0)), axis=1 ).reshape( -1, 1 ) pattern = PoyntingVector(pattern) pattern = Directivity(pattern) pattern = field_magnitude_phase(pattern) new_list.append(pattern) return new_list, element_positions, frequency def Array_Panel(mesh_size): from lyceanem.utility.mesh_functions import pyvista_to_meshio major_size = 345.5e-3 minor_size = 170.5e-3 element_sep = 43e-3 element_point_major = 301e-3 element_point_minor = 129e-3 columns = 8 rows = 4 element_area = element_sep**2 x, y, z = np.meshgrid( np.linspace(-element_point_major / 2, element_point_major / 2, columns), np.linspace(-element_point_minor / 2, element_point_minor / 2, rows), 1e-6, ) element_points = np.array([x.ravel(), y.ravel(), z.ravel()]).transpose() element_mesh = pv.PolyData(np.array([x.ravel(), y.ravel(), z.ravel()]).transpose()) element_mesh.point_data["Area"] = ( np.ones((element_mesh.points.shape[0])) * element_area ) element_mesh.point_data["Normals"] = np.repeat( np.array([0, 0, 1]).reshape(1, 3), element_mesh.points.shape[0], axis=0 ) from lyceanem.geometry.geometryfunctions import compute_areas, compute_normals import gmsh gmsh.initialize() gmsh.option.setNumber("Mesh.CharacteristicLengthMax", mesh_size) box = gmsh.model.occ.addBox( -major_size / 2, -minor_size / 2, -5e-3, major_size, minor_size, 5e-3 ) gmsh.model.occ.synchronize() gmsh.model.mesh.generate(dim=2) file_name = "temp_mesh.stl" gmsh.write(file_name) gmsh.finalize() panel_mesh = compute_normals(compute_areas(meshio.read(file_name))) return element_mesh, panel_mesh .. GENERATED FROM PYTHON SOURCE LINES 160-166 Importing measured Validation Data ------------------------------------------- The function Panel_Array_data imports the measured validation data for the Panel antenna array, including sample positions within the University of Bristol Anechoic chamber, and the measured antenna patterns for each element in the array. This is then used to define the wavelength of operation, and set the desired maximum mesh length for the model antenna array. The function Array_Panel is defined in such a way that the elements are generated with the correct ordering for the antenna array as it was measured and recorded in the chamber and validation data. .. GENERATED FROM PYTHON SOURCE LINES 166-171 .. code-block:: Python antennas, element_positions, frequency = Panel_Array_data() wavelength = speed_of_light / frequency elements, panel = Array_Panel(wavelength * 0.5) .. rst-class:: sphx-glr-script-out .. code-block:: none C:\Users\lycea\miniconda3\envs\CudaDevelopment\Lib\site-packages\meshio\stl\_stl.py:40: RuntimeWarning: overflow encountered in scalar multiply if 84 + num_triangles * 50 == filesize_bytes: .. GENERATED FROM PYTHON SOURCE LINES 172-176 Visualisation of the Modelled Array ---------------------------------------- Once the antenna array has been generated for simulation, it can then be visualised within the model space, and a rendering exported. .. GENERATED FROM PYTHON SOURCE LINES 176-183 .. code-block:: Python pl = pv.Plotter() pl.add_mesh(elements) pl.add_mesh(pv.from_meshio(panel), show_edges=True) pl.add_axes() pl.show(screenshot="ArrayasSimulated.png") .. tab-set:: .. tab-item:: Static Scene .. image-sg:: /auto_examples/images/sphx_glr_08_array_validation_001.png :alt: 08 array validation :srcset: /auto_examples/images/sphx_glr_08_array_validation_001.png :class: sphx-glr-single-img .. tab-item:: Interactive Scene .. offlineviewer:: C:\Users\lycea\PycharmProjects\LyceanEM-Python\docs\source\auto_examples\images\sphx_glr_08_array_validation_001.vtksz .. GENERATED FROM PYTHON SOURCE LINES 184-188 Definition of the sample field and simulation ---------------------------------------------- In order to make comparisons straightforward, the sample field positions are copied from the measured data while removing the actual field data itself. The excitation_function is called to allow the definition of electric field vectors for each antenna with the same polarization as produced by the patch antennas of the validation antenna array. .. GENERATED FROM PYTHON SOURCE LINES 188-223 .. code-block:: Python import lyceanem.models.frequency_domain as FD from lyceanem.electromagnetics.emfunctions import excitation_function desired_e_vector = np.array([1.0, 0.0, 0]) weights = excitation_function( elements, desired_e_vector, wavelength=wavelength, phase_shift="wavefront", steering_vector=np.array([0, 0, 1]), ) sample_field = antennas[0].copy() sample_field.point_data.clear() sample_field.point_data["Normals"] = antennas[0].point_data["Normals"] from lyceanem.base_classes import structures blockers = structures(solids=[panel]) Ex, Ey, Ez = FD.calculate_scattering( aperture_coords=elements, sink_coords=sample_field, scatter_points=panel, antenna_solid=blockers, desired_E_axis=weights, wavelength=wavelength, scattering=0, project_vectors=False, beta=(2 * np.pi) / wavelength, ) .. rst-class:: sphx-glr-script-out .. code-block:: none C:\Users\lycea\miniconda3\envs\CudaDevelopment\Lib\site-packages\lyceanem\electromagnetics\empropagation.py:3719: ComplexWarning: Casting complex values to real discards the imaginary part uvn_axes[2, :] = point_vector C:\Users\lycea\miniconda3\envs\CudaDevelopment\Lib\site-packages\lyceanem\electromagnetics\empropagation.py:3731: ComplexWarning: Casting complex values to real discards the imaginary part uvn_axes[0, :] = np.cross(local_axes[0, :], point_vector) / np.linalg.norm( C:\Users\lycea\miniconda3\envs\CudaDevelopment\Lib\site-packages\lyceanem\electromagnetics\empropagation.py:3758: ComplexWarning: Casting complex values to real discards the imaginary part uvn_axes[1, :] = np.cross(point_vector, uvn_axes[0, :]) / np.linalg.norm( .. GENERATED FROM PYTHON SOURCE LINES 224-228 Processing of the Modelled results --------------------------------------- The function update_electric_fields is used to populate the sample field with the simulated farfield pattern of the antenna array with the calculated weights. A wide variety of different quantities can then be calculated to allow consistent analysis of the measured and simulated antenna patterns. .. GENERATED FROM PYTHON SOURCE LINES 228-289 .. code-block:: Python from lyceanem.electromagnetics.emfunctions import ( update_electric_fields, PoyntingVector, Directivity, Exyz_to_EthetaEphi, field_magnitude_phase, ) from lyceanem.electromagnetics.beamforming import create_display_mesh measured_array_pattern = field_magnitude_phase( Directivity(PoyntingVector(Beamform(antennas, weights[:, 1]))) ) dynamic_range = 60 sample_field = field_magnitude_phase( Directivity(PoyntingVector(update_electric_fields(sample_field, Ex, Ey, Ez))) ) display = create_display_mesh( sample_field, label="D(Total)", field_radius=1.0, dynamic_range=dynamic_range ) display.point_data["D(Total)_(dBi)"] = 10 * np.log10(display.point_data["D(Total)"]) measured_display = create_display_mesh( measured_array_pattern, label="D(Total)", field_radius=1.0, dynamic_range=dynamic_range, ) measured_display.point_data["D(Total)_(dBi)"] = 10 * np.log10( measured_display.point_data["D(Total)"] ) print("Maximum Directivity ", np.max(display.point_data["D(Total)_(dBi)"])) print( "Measured Maximum Directivity ", np.max(measured_display.point_data["D(Total)_(dBi)"]), ) pl = pv.Plotter(shape=(1, 2), border=False) pl.subplot(0, 0) pl.add_mesh(elements, color="aqua") pl.add_mesh(panel, color="grey") pl.add_mesh(display, scalars="D(Total)_(dBi)", clim=[22 - dynamic_range, 22]) pl.add_axes() pl.add_title("LyceanEM") pl.subplot(0, 1) pl.add_mesh(elements, color="aqua") pl.add_mesh(panel, color="grey") pl.add_mesh(measured_display, scalars="D(Total)_(dBi)", clim=[22 - dynamic_range, 22]) pl.add_title("Measured") pl.link_views() pl.view_isometric() pl.show(screenshot="MeasurementtoSimulationComparison.png") simulated = display.point_data["D(Total)"] measured = measured_display.point_data["D(Total)"] print( "Correlation {}".format(np.corrcoef(simulated, measured)[0, 1]), " between measurement and simulation", ) .. tab-set:: .. tab-item:: Static Scene .. image-sg:: /auto_examples/images/sphx_glr_08_array_validation_002.png :alt: 08 array validation :srcset: /auto_examples/images/sphx_glr_08_array_validation_002.png :class: sphx-glr-single-img .. tab-item:: Interactive Scene .. offlineviewer:: C:\Users\lycea\PycharmProjects\LyceanEM-Python\docs\source\auto_examples\images\sphx_glr_08_array_validation_002.vtksz .. rst-class:: sphx-glr-script-out .. code-block:: none C:\Users\lycea\miniconda3\envs\CudaDevelopment\Lib\site-packages\lyceanem\electromagnetics\emfunctions.py:539: RuntimeWarning: divide by zero encountered in log10 field_data.point_data["Poynting_Vector_(Magnitude_(dBW/m2))"] = 10 * np.log10( C:\Users\lycea\miniconda3\envs\CudaDevelopment\Lib\site-packages\lyceanem\electromagnetics\beamforming.py:1617: RuntimeWarning: divide by zero encountered in log10 logdata = log_multiplier * np.log10(pattern_mesh.point_data[label]) C:\Users\lycea\PycharmProjects\LyceanEM-Python\docs\source\examples\08_array_validation.py:249: RuntimeWarning: divide by zero encountered in log10 display.point_data["D(Total)_(dBi)"] = 10 * np.log10(display.point_data["D(Total)"]) Maximum Directivity 20.07648569434074 Measured Maximum Directivity 21.082011910586985 Correlation 0.9966437988952569 between measurement and simulation .. rst-class:: sphx-glr-timing **Total running time of the script:** (1 minutes 8.238 seconds) .. _sphx_glr_download_auto_examples_08_array_validation.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: 08_array_validation.ipynb <08_array_validation.ipynb>` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: 08_array_validation.py <08_array_validation.py>` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: 08_array_validation.zip <08_array_validation.zip>` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_