From 056f1ee39a5c0d976d3b86d306c05b9aba17f951 Mon Sep 17 00:00:00 2001 From: LyceanEM <60020395+LyceanEM@users.noreply.github.com> Date: Wed, 25 Sep 2024 23:59:16 +0100 Subject: [PATCH] Corrected type bug in lossy_propagation function to access the alpha and beta values correctly. --- .../03_frequency_domain_channel_modelling.py | 2 +- lyceanem/electromagnetics/emfunctions.py | 16 +++++++--- lyceanem/raycasting/rayfunctions.py | 29 +++++++++---------- lyceanem/utility/math_functions.py | 4 +-- 4 files changed, 28 insertions(+), 23 deletions(-) diff --git a/docs/examples/03_frequency_domain_channel_modelling.py b/docs/examples/03_frequency_domain_channel_modelling.py index c76173f..6ec5a5c 100644 --- a/docs/examples/03_frequency_domain_channel_modelling.py +++ b/docs/examples/03_frequency_domain_channel_modelling.py @@ -21,7 +21,7 @@ # Frequency and Mesh Resolution # ------------------------------ # -freq = np.asarray(26.0e9) +freq = np.asarray(22.0e9) wavelength = 3e8 / freq mesh_resolution = 0.5 * wavelength diff --git a/lyceanem/electromagnetics/emfunctions.py b/lyceanem/electromagnetics/emfunctions.py index aa5501f..ba7153e 100644 --- a/lyceanem/electromagnetics/emfunctions.py +++ b/lyceanem/electromagnetics/emfunctions.py @@ -160,6 +160,14 @@ def EthetaEphi_to_Exyz(field_data): def Exyz_to_EthetaEphi(field_data): # this function assumes a spherical field definition, will need to write a function which works based on the poynting vector/normal vector of the point + if not all( + k in field_data.point_data.keys() for k in ("theta (Radians)", "phi (Radians)") + ): + # theta and phi don't exist in the dataset + from lyceanem.geometry.geometryfunctions import theta_phi_r + + field_data = theta_phi_r(field_data) + electric_fields = extract_electric_fields(field_data) theta = field_data.point_data["theta (Radians)"] phi = field_data.point_data["phi (Radians)"] @@ -372,7 +380,7 @@ def PoyntingVector(field_data,measurement=False,aperture=None): def Directivity(field_data): - # Directivity is defined in terms of the radiation intensity in a given direction which is the power per unit solid angle, divided by the average power per unit solid angle + # calculate directivity for the given pattern if not all( k in field_data.point_data.keys() for k in ("theta (Radians)", "phi (Radians)") @@ -387,7 +395,7 @@ def Directivity(field_data): ): # E(theta) and E(phi) are not present in the data field_data = Exyz_to_EthetaEphi(field_data) - + if not all( k in field_data.point_data.keys() for k in ("Area") ): @@ -412,8 +420,8 @@ def Directivity(field_data): #Calculate Solid Angle solid_angle=field_data.point_data["Area"]/field_data.point_data["Radial Distance (m)"]**2 Utotal = Utheta + Uphi - - Uav=(np.sum(Utotal.ravel()*solid_angle)/(4*np.pi)) + + Uav=(np.nansum(Utotal.ravel()*solid_angle)/(4*np.pi)) #sin_factor = np.abs( # np.sin(field_data.point_data["theta (Radians)"]) #).reshape(-1,1) # only want positive factor diff --git a/lyceanem/raycasting/rayfunctions.py b/lyceanem/raycasting/rayfunctions.py index 5986ef5..09690f6 100644 --- a/lyceanem/raycasting/rayfunctions.py +++ b/lyceanem/raycasting/rayfunctions.py @@ -951,14 +951,14 @@ def CalculatePoyntingVectors( """ max_path_divergence = 2.0 if len(total_network) == 0: - sourcenum = np.int(np.nanmax(scattering_index[:, 0])) + 1 - sinknum = np.int(np.nanmax(scattering_index[:, 1])) + 1 + sourcenum = np.nanmax(scattering_index[:, 0]).astype(int) + 1 + sinknum = np.nanmax(scattering_index[:, 1]).astype(int) + 1 aoa_spectrum = np.zeros((sinknum, az_bins.shape[0]), dtype=np.float32) delay_spectrum = np.zeros((sinknum, az_bins.shape[0]), dtype=np.float32) else: wave_vector = (2 * np.pi) / wavelength - sourcenum = np.int(np.nanmax(scattering_index[:, 0])) + 1 - sinknum = np.int(np.nanmax(scattering_index[:, 1])) + 1 + sourcenum = np.nanmax(scattering_index[:, 0]).astype(int) + 1 + sinknum = np.nanmax(scattering_index[:, 1]).astype(int) + 1 scatterdepth = int((len(total_network[0, :]) - 3) / 3) impulse_response = np.zeros((len(time_bins), sinknum), dtype="complex") angle_response = np.zeros((len(az_bins), sinknum), dtype="complex") @@ -1485,8 +1485,8 @@ def VectorNetworkProcessv2( # """ if scattering_index.shape[0] == 0: - sourcenum = np.int(np.nanmax(scattering_index[:, 0])) + 1 - sinknum = np.int(np.nanmax(scattering_index[:, 1])) + 1 + sourcenum = np.nanmax(scattering_index[:, 0]).astype(int) + 1 + sinknum = np.nanmax(scattering_index[:, 1]).astype(int) + 1 scatter_map = np.zeros((sourcenum, sinknum, 1), dtype="complex") else: unified_model = np.append( @@ -1766,8 +1766,8 @@ def scatter_net_sortTimeEM( def BaseNetworkProcess(total_network, wavelength, scattering_index): # take the total_network matrix, and work out the path length of each connection, and hence superposition wave_vector = (2 * np.pi) / wavelength - sourcenum = np.int(np.nanmax(scattering_index[:, 0])) + 1 - sinknum = np.int(np.nanmax(scattering_index[:, 1])) + 1 + sourcenum = np.nanmax(scattering_index[:, 0]).astype(int) + 1 + sinknum = np.nanmax(scattering_index[:, 1]).astype(int) + 1 scatterdepth = int((len(total_network[1, :]) - 3) / 3) scatter_map = np.zeros((sourcenum, sinknum, scatterdepth), dtype="complex") depth_num = 0 @@ -1793,8 +1793,8 @@ def BaseNetworkProcess(total_network, wavelength, scattering_index): end_row = np.max(np.where((np.isnan(total_network[:, 6 + ((idx) * 3)])))) for idr in range(start_row, end_row + 1): - source_index = np.int(scattering_index[idr, 0]) - sink_index = np.int(scattering_index[idr, 1]) + source_index = scattering_index[idr, 0].astype(int) + sink_index = scattering_index[idr, 1].astype(int) path_distance = np.zeros((1)) if idx == 0: path_distance = ( @@ -2894,7 +2894,7 @@ def workchunkingv2( elif max_scatter == 3: full_index = np.empty((0, 4), dtype=np.int32) - chunknum = np.minimum(sources.shape[0], np.int(np.ceil(ray_estimate / 2e8))) + chunknum = np.minimum(sources.shape[0], np.ceil(ray_estimate / 2e8).astype(int)) chunknum = np.maximum(2, chunknum) source_chunking = np.linspace(0, sources.shape[0], chunknum, dtype=np.int32) for chunkindex in range(source_chunking.size - 1): @@ -2925,9 +2925,8 @@ def workchunkingv2( sources, sinks, scattering_points, filtered_index, environment, False ) if filtered_index2.shape[0] * sinks.shape[0] > 2e8: - temp_chunks = np.int( - np.ceil((filtered_index2.shape[0] * sinks.shape[0]) / 2e8) - ) + temp_chunks = np.ceil((filtered_index2.shape[0] * sinks.shape[0]) / 2e8).astype(int) + temp_chunking = np.linspace( 0, filtered_index2.shape[0], temp_chunks, dtype=np.int32 ) @@ -3111,7 +3110,7 @@ def integratedraycastersetup( ).reshape((scattering_points.shape[0] - sources - sinks), 1) # implement chunking full_index = np.empty((0, np.max(scattering_mask) + 2), dtype=np.int32) - chunknum = np.minimum(sources, np.int(np.ceil(ray_estimate / 2e8))) + chunknum = np.minimum(sources, np.ceil(ray_estimate / 2e8).astype(int)) chunknum = np.maximum(2, chunknum) source_chunking = np.linspace(0, sources, chunknum, dtype=np.int32) for chunkindex in range(source_chunking.size - 1): diff --git a/lyceanem/utility/math_functions.py b/lyceanem/utility/math_functions.py index 1bdf2fb..6edbf67 100644 --- a/lyceanem/utility/math_functions.py +++ b/lyceanem/utility/math_functions.py @@ -104,9 +104,7 @@ def discrete_transmit_power( power_at_point * power_normalisation ) / element_area.reshape(-1, 1) # calculate amplitude (V/m) - transmit_amplitude = ((transmit_power_density * impedance) ** 0.5) * ( - weights / np.linalg.norm(weights, axis=1).reshape(-1, 1) - ) + transmit_amplitude = ((transmit_power_density * impedance) ** 0.5) # transmit_excitation=transmit_amplitude_density.reshape(-1,1)*element_area.reshape(-1,1) return transmit_amplitude