diff --git a/sharc/antenna/antenna_s1528.py b/sharc/antenna/antenna_s1528.py index b8ee9076..de9a966b 100644 --- a/sharc/antenna/antenna_s1528.py +++ b/sharc/antenna/antenna_s1528.py @@ -175,12 +175,18 @@ def __init__(self, param: ParametersAntennaS1528): # the system design self.l_s = param.antenna_l_s + # far-out side-lobe level [dBi] + if param.far_out_side_lobe is None: + self.l_f = 0 + else: + self.l_f = param.far_out_side_lobe + # for elliptical antennas, this is the ratio major axis/minor axis # we assume circular antennas, so z = 1 - self.z = 1 - - # far-out side-lobe level [dBi] - self.l_f = 0 + if param.major_minor_axis_ratio is None: + self.z = 1 + else: + self.z = param.major_minor_axis_ratio # back-lobe level self.l_b = np.maximum( @@ -190,19 +196,22 @@ def __init__(self, param: ParametersAntennaS1528): # one-half the 3 dB beamwidth in the plane of interest self.psi_b = param.antenna_3_dB_bw / 2 - if self.l_s == -15: - self.a = 2.58 * math.sqrt(1 - 1.4 * math.log10(self.z)) - elif self.l_s == -20: - self.a = 2.58 * math.sqrt(1 - 1.0 * math.log10(self.z)) - elif self.l_s == -25: - self.a = 2.58 * math.sqrt(1 - 0.6 * math.log10(self.z)) - elif self.l_s == -30: - self.a = 2.58 * math.sqrt(1 - 0.4 * math.log10(self.z)) + if math.isclose(self.z, 1.0): + self.a = 2.58 else: - sys.stderr.write( - "ERROR\nInvalid AntennaS1528 L_s parameter: " + str(self.l_s), - ) - sys.exit(1) + if self.l_s == -15: + self.a = 2.58 * math.sqrt(1 - 1.4 * math.log10(self.z)) + elif self.l_s == -20: + self.a = 2.58 * math.sqrt(1 - 1.0 * math.log10(self.z)) + elif self.l_s == -25: + self.a = 2.58 * math.sqrt(1 - 0.6 * math.log10(self.z)) + elif self.l_s == -30: + self.a = 2.58 * math.sqrt(1 - 0.4 * math.log10(self.z)) + else: + sys.stderr.write( + f"ERROR\nInvalid AntennaS1528 L_s parameter {self.l_s} for z={self.z}" + ) + sys.exit(1) self.b = 6.32 self.alpha = 1.5 diff --git a/sharc/parameters/antenna/parameters_antenna_s1528.py b/sharc/parameters/antenna/parameters_antenna_s1528.py index 4cd5ea35..ed2f204f 100644 --- a/sharc/parameters/antenna/parameters_antenna_s1528.py +++ b/sharc/parameters/antenna/parameters_antenna_s1528.py @@ -37,6 +37,14 @@ class ParametersAntennaS1528(ParametersBase): l_r: float = None l_t: float = None + # Recommends 1.2 only + # For elliptical antennas, this is the ratio major axis/minor axis + # we assume circular antennas, so z = 1 + major_minor_axis_ratio: float = None + + # Far-out side-lobe level + far_out_side_lobe: float = None + def load_parameters_from_file(self, config_file: str): """Load the parameters from file an run a sanity check. diff --git a/sharc/parameters/imt/parameters_grid.py b/sharc/parameters/imt/parameters_grid.py index bd6e0144..551a5d2d 100644 --- a/sharc/parameters/imt/parameters_grid.py +++ b/sharc/parameters/imt/parameters_grid.py @@ -339,6 +339,9 @@ class ParametersSatelliteWithServiceGrid(ParametersTerrestrialGrid): beam_radius: float = None + # [deg] + minimum_service_angle: float = 5.0 + def validate(self, ctx): """Validates instance parameters. Parameters @@ -353,6 +356,9 @@ def validate(self, ctx): ) self.cell_radius = self.beam_radius + if self.minimum_service_angle < 0. or self.minimum_service_angle > 90: + raise ValueError(f"{ctx}.minimum_service_angle should be in [0, 90]") + if not isinstance( self.eligible_sats_margin_from_border, float) and not isinstance( diff --git a/sharc/propagation/atmosphere.py b/sharc/propagation/atmosphere.py index 846ce151..9d520803 100644 --- a/sharc/propagation/atmosphere.py +++ b/sharc/propagation/atmosphere.py @@ -385,6 +385,9 @@ def get_reference_atmosphere_p835( """ h_km = altitude / 1000 + season = season.lower() + if season not in ['winter', 'summer']: + raise ValueError(f"Invalid season name {season}.") if latitude <= 22: # low latitude diff --git a/sharc/topology/topology_imt_mss_dc.py b/sharc/topology/topology_imt_mss_dc.py index bd0f8cbb..337cbe43 100644 --- a/sharc/topology/topology_imt_mss_dc.py +++ b/sharc/topology/topology_imt_mss_dc.py @@ -75,117 +75,101 @@ def get_coordinates( ): """Compute the coordinates of the visible space stations.""" orbit_params.sat_is_active_if.validate("orbit_params.sat_is_active_if") - # Calculate the total number of satellites across all orbits - total_satellites = sum( - orbit.n_planes * - orbit.sats_per_plane for orbit in orbit_params.orbits) - if any([ - not hasattr(orbit_params, attr) - for attr in ["sat_is_active_if", "orbits", "beam_radius", "num_beams", "beam_positioning"] - ]): - raise ValueError( - "Parameter passed to TopologyImtMssDc needs to contain all of the attributes:\n" - '["sat_is_active_if", "orbits", "beam_radius", "num_beams", "beam_positioning"]') - - idx_orbit = np.zeros( - total_satellites, - dtype=int) # Add orbit index array - - # List to store indices of active satellites - active_satellite_idxs = [] + num_base_stations = 0 MAX_ITER = 10000 # Maximum iterations to find at least one visible satellite - i = 0 # Iteration counter for ensuring satellite visibility - while len(active_satellite_idxs) == 0: - # Initialize arrays to store satellite positions, angles and - # distance from center of earth - all_positions = { - "R": [], - "lat": [], - "lon": [], - "sx": [], - "sy": [], - "sz": [], - "alt": []} - all_elevations = [] # Store satellite elevations - all_azimuths = [] # Store satellite azimuths - - current_sat_idx = 0 # Index tracker for satellites across all orbits - - # Iterate through each orbit defined in the parameters - for orbit_idx, param in enumerate(orbit_params.orbits): - orbit = OrbitModel( - Nsp=param.sats_per_plane, # Satellites per plane - Np=param.n_planes, # Number of orbital planes - phasing=param.phasing_deg, # Phasing angle in degrees - long_asc=param.long_asc_deg, # Longitude of ascending node in degrees - omega=param.omega_deg, # Argument of perigee in degrees - delta=param.inclination_deg, # Orbital inclination in degrees - hp=param.perigee_alt_km, # Perigee altitude in kilometers - ha=param.apogee_alt_km, # Apogee altitude in kilometers - Mo=param.initial_mean_anomaly, # Initial mean anomaly in degrees - # whether to use only time as random variable - model_time_as_random_variable=param.model_time_as_random_variable, - t_min=param.t_min, - t_max=param.t_max, - ) - # Generate random positions for satellites in this orbit - pos_vec = orbit.get_orbit_positions_random( - rng=random_number_gen) - - # Determine the number of satellites in this orbit - num_satellites = len(pos_vec["sx"]) - - # Assign orbit index to satellites - idx_orbit[current_sat_idx:current_sat_idx + - num_satellites] = orbit_idx - - # Extract satellite positions and calculate distances - sx, sy, sz = pos_vec['sx'], pos_vec['sy'], pos_vec['sz'] - # Distance from Earth's center - r = np.sqrt(sx**2 + sy**2 + sz**2) - - # When getting azimuth and elevation, we need to consider sx, sy and sz points - # from the center of earth to the satellite, and we need to point the satellite - # towards the center of earth - # Calculate elevation angles - elevations = np.degrees(np.arcsin(-sz / r)) - # Calculate azimuth angles - azimuths = np.degrees(np.arctan2(-sy, -sx)) - - # Append satellite positions and angles to global lists - all_positions['lat'].extend(pos_vec['lat']) # Latitudes - all_positions['lon'].extend(pos_vec['lon']) # Longitudes - all_positions['sx'].extend(sx) # X-coordinates - all_positions['sy'].extend(sy) # Y-coordinates - all_positions['sz'].extend(sz) # Z-coordinates - all_positions["R"].extend(r) - all_positions["alt"].extend(pos_vec['alt']) - all_elevations.extend(elevations) # Elevation angles - all_azimuths.extend(azimuths) # Azimuth angles - - active_sats_mask = np.ones(len(pos_vec['lat']), dtype=bool) - - if "MINIMUM_ELEVATION_FROM_ES" in orbit_params.sat_is_active_if.conditions: - # Calculate satellite visibility from base stations - elev_from_bs = calc_elevation( - coordinate_system.ref_lat, # Latitude of base station - pos_vec['lat'], # Latitude of satellites - coordinate_system.ref_long, # Longitude of base station - pos_vec['lon'], # Longitude of satellites - # Perigee altitude in kilometers - sat_height=pos_vec['alt'] * 1e3, - es_height=coordinate_system.ref_alt, + while num_base_stations == 0: + # Calculate the total number of satellites across all orbits + total_satellites = sum( + orbit.n_planes * + orbit.sats_per_plane for orbit in orbit_params.orbits) + if any([ + not hasattr(orbit_params, attr) + for attr in ["sat_is_active_if", "orbits", "beam_radius", "num_beams", "beam_positioning"] + ]): + raise ValueError( + "Parameter passed to TopologyImtMssDc needs to contain all of the attributes:\n" + '["sat_is_active_if", "orbits", "beam_radius", "num_beams", "beam_positioning"]') + + idx_orbit = np.zeros( + total_satellites, + dtype=int) # Add orbit index array + + # List to store indices of active satellites + active_satellite_idxs = [] + + i = 0 # Iteration counter for ensuring satellite visibility + while len(active_satellite_idxs) == 0: + # Initialize arrays to store satellite positions, angles and + # distance from center of earth + all_positions = { + "R": [], + "lat": [], + "lon": [], + "sx": [], + "sy": [], + "sz": [], + "alt": []} + all_elevations = [] # Store satellite elevations + all_azimuths = [] # Store satellite azimuths + + current_sat_idx = 0 # Index tracker for satellites across all orbits + + # Iterate through each orbit defined in the parameters + for orbit_idx, param in enumerate(orbit_params.orbits): + orbit = OrbitModel( + Nsp=param.sats_per_plane, # Satellites per plane + Np=param.n_planes, # Number of orbital planes + phasing=param.phasing_deg, # Phasing angle in degrees + long_asc=param.long_asc_deg, # Longitude of ascending node in degrees + omega=param.omega_deg, # Argument of perigee in degrees + delta=param.inclination_deg, # Orbital inclination in degrees + hp=param.perigee_alt_km, # Perigee altitude in kilometers + ha=param.apogee_alt_km, # Apogee altitude in kilometers + Mo=param.initial_mean_anomaly, # Initial mean anomaly in degrees + # whether to use only time as random variable + model_time_as_random_variable=param.model_time_as_random_variable, + t_min=param.t_min, + t_max=param.t_max, ) - - # Determine visible satellites based on minimum elevation - # angle - active_sats_mask = active_sats_mask & (elev_from_bs.flatten( - ) >= orbit_params.sat_is_active_if.minimum_elevation_from_es) - - if "MAXIMUM_ELEVATION_FROM_ES" in orbit_params.sat_is_active_if.conditions: - # no need to recalculate if already calculated above - if "MINIMUM_ELEVATION_FROM_ES" not in orbit_params.sat_is_active_if.conditions: + # Generate random positions for satellites in this orbit + pos_vec = orbit.get_orbit_positions_random( + rng=random_number_gen) + + # Determine the number of satellites in this orbit + num_satellites = len(pos_vec["sx"]) + + # Assign orbit index to satellites + idx_orbit[current_sat_idx:current_sat_idx + + num_satellites] = orbit_idx + + # Extract satellite positions and calculate distances + sx, sy, sz = pos_vec['sx'], pos_vec['sy'], pos_vec['sz'] + # Distance from Earth's center + r = np.sqrt(sx**2 + sy**2 + sz**2) + + # When getting azimuth and elevation, we need to consider sx, sy and sz points + # from the center of earth to the satellite, and we need to point the satellite + # towards the center of earth + # Calculate elevation angles + elevations = np.degrees(np.arcsin(-sz / r)) + # Calculate azimuth angles + azimuths = np.degrees(np.arctan2(-sy, -sx)) + + # Append satellite positions and angles to global lists + all_positions['lat'].extend(pos_vec['lat']) # Latitudes + all_positions['lon'].extend(pos_vec['lon']) # Longitudes + all_positions['sx'].extend(sx) # X-coordinates + all_positions['sy'].extend(sy) # Y-coordinates + all_positions['sz'].extend(sz) # Z-coordinates + all_positions["R"].extend(r) + all_positions["alt"].extend(pos_vec['alt']) + all_elevations.extend(elevations) # Elevation angles + all_azimuths.extend(azimuths) # Azimuth angles + + active_sats_mask = np.ones(len(pos_vec['lat']), dtype=bool) + + if "MINIMUM_ELEVATION_FROM_ES" in orbit_params.sat_is_active_if.conditions: # Calculate satellite visibility from base stations elev_from_bs = calc_elevation( coordinate_system.ref_lat, # Latitude of base station @@ -197,122 +181,141 @@ def get_coordinates( es_height=coordinate_system.ref_alt, ) - # Determine visible satellites based on minimum elevation - # angle - active_sats_mask = active_sats_mask & (elev_from_bs.flatten( - ) <= orbit_params.sat_is_active_if.maximum_elevation_from_es) - - if "LAT_LONG_INSIDE_COUNTRY" in orbit_params.sat_is_active_if.conditions: - flat_active_lon = pos_vec["lon"].flatten()[ - active_sats_mask] - flat_active_lat = pos_vec["lat"].flatten()[ - active_sats_mask] - - # create points(lon, lat) to compare to country - sats_points = gpd.points_from_xy( - flat_active_lon, flat_active_lat, crs=EARTH_DEFAULT_CRS) - - # Check if the satellite is inside the country polygon - polygon_mask = np.zeros_like(active_sats_mask) - polygon_mask[active_sats_mask] = sats_points.within( - orbit_params.sat_is_active_if.lat_long_inside_country.filter_polygon) - - active_sats_mask = active_sats_mask & polygon_mask - - visible_sat_idxs = np.arange( - current_sat_idx, current_sat_idx + len(pos_vec['lat']), dtype=int - )[active_sats_mask] - active_satellite_idxs.extend(visible_sat_idxs) - - # Update the index tracker for the next orbit - current_sat_idx += len(sx) - - i += 1 # Increment iteration counter - if i >= MAX_ITER: # Check if maximum iterations reached - raise RuntimeError( - "Maximum iterations reached, and no satellite was selected within the minimum elevation criteria." - ) - # We have the list of visible satellites, now create a Topolgy of this subset and move the coordinate system - # reference. - # Convert X-coordinates to meters - all_space_station_x = np.ravel(np.array(all_positions['sx'])) * 1e3 - # Convert Y-coordinates to meters - all_space_station_y = np.ravel(np.array(all_positions['sy'])) * 1e3 - # Convert Z-coordinates to meters - all_space_station_z = np.ravel(np.array(all_positions['sz'])) * 1e3 - all_elevation = np.ravel(np.array(all_elevations)) # Elevation angles - all_azimuth = np.ravel(np.array(all_azimuths)) # Azimuth angles - all_lat = np.ravel(np.array(all_positions['lat'])) - all_lon = np.ravel(np.array(all_positions['lon'])) - all_sat_altitude = np.ravel(np.array(all_positions['alt'])) * 1e3 - - total_active_satellites = len(active_satellite_idxs) - space_station_x = all_space_station_x[active_satellite_idxs] - space_station_y = all_space_station_y[active_satellite_idxs] - space_station_z = all_space_station_z[active_satellite_idxs] - elevation = all_elevation[active_satellite_idxs] - azimuth = all_azimuth[active_satellite_idxs] - lat = all_lat[active_satellite_idxs] - lon = all_lon[active_satellite_idxs] - sat_altitude = all_sat_altitude[active_satellite_idxs] - sat_altitude = all_sat_altitude[active_satellite_idxs] + # Determine visible satellites based on minimum elevation + # angle + active_sats_mask = active_sats_mask & (elev_from_bs.flatten( + ) >= orbit_params.sat_is_active_if.minimum_elevation_from_es) + + if "MAXIMUM_ELEVATION_FROM_ES" in orbit_params.sat_is_active_if.conditions: + # no need to recalculate if already calculated above + if "MINIMUM_ELEVATION_FROM_ES" not in orbit_params.sat_is_active_if.conditions: + # Calculate satellite visibility from base stations + elev_from_bs = calc_elevation( + coordinate_system.ref_lat, # Latitude of base station + pos_vec['lat'], # Latitude of satellites + coordinate_system.ref_long, # Longitude of base station + pos_vec['lon'], # Longitude of satellites + # Perigee altitude in kilometers + sat_height=pos_vec['alt'] * 1e3, + es_height=coordinate_system.ref_alt, + ) + + # Determine visible satellites based on minimum elevation + # angle + active_sats_mask = active_sats_mask & (elev_from_bs.flatten( + ) <= orbit_params.sat_is_active_if.maximum_elevation_from_es) + + if "LAT_LONG_INSIDE_COUNTRY" in orbit_params.sat_is_active_if.conditions: + flat_active_lon = pos_vec["lon"].flatten()[ + active_sats_mask] + flat_active_lat = pos_vec["lat"].flatten()[ + active_sats_mask] + + # create points(lon, lat) to compare to country + sats_points = gpd.points_from_xy( + flat_active_lon, flat_active_lat, crs=EARTH_DEFAULT_CRS) + + # Check if the satellite is inside the country polygon + polygon_mask = np.zeros_like(active_sats_mask) + polygon_mask[active_sats_mask] = sats_points.within( + orbit_params.sat_is_active_if.lat_long_inside_country.filter_polygon) + + active_sats_mask = active_sats_mask & polygon_mask + + visible_sat_idxs = np.arange( + current_sat_idx, current_sat_idx + len(pos_vec['lat']), dtype=int + )[active_sats_mask] + active_satellite_idxs.extend(visible_sat_idxs) + + # Update the index tracker for the next orbit + current_sat_idx += len(sx) + + i += 1 # Increment iteration counter + if i >= MAX_ITER: # Check if maximum iterations reached + raise RuntimeError( + "Maximum iterations reached, and no satellite was selected within the minimum elevation criteria." + ) + # We have the list of visible satellites, now create a Topolgy of this subset and move the coordinate system + # reference. + # Convert X-coordinates to meters + all_space_station_x = np.ravel(np.array(all_positions['sx'])) * 1e3 + # Convert Y-coordinates to meters + all_space_station_y = np.ravel(np.array(all_positions['sy'])) * 1e3 + # Convert Z-coordinates to meters + all_space_station_z = np.ravel(np.array(all_positions['sz'])) * 1e3 + all_elevation = np.ravel(np.array(all_elevations)) # Elevation angles + all_azimuth = np.ravel(np.array(all_azimuths)) # Azimuth angles + all_lat = np.ravel(np.array(all_positions['lat'])) + all_lon = np.ravel(np.array(all_positions['lon'])) + all_sat_altitude = np.ravel(np.array(all_positions['alt'])) * 1e3 + + total_active_satellites = len(active_satellite_idxs) + space_station_x = all_space_station_x[active_satellite_idxs] + space_station_y = all_space_station_y[active_satellite_idxs] + space_station_z = all_space_station_z[active_satellite_idxs] + elevation = all_elevation[active_satellite_idxs] + azimuth = all_azimuth[active_satellite_idxs] + lat = all_lat[active_satellite_idxs] + lon = all_lon[active_satellite_idxs] + sat_altitude = all_sat_altitude[active_satellite_idxs] + sat_altitude = all_sat_altitude[active_satellite_idxs] + + # Convert the ECEF coordinates to the transformed cartesian coordinates and set the Space Station positions + # used to generetate the IMT Base Stations + space_station_x, space_station_y, space_station_z = \ + coordinate_system.ecef2enu(space_station_x, space_station_y, space_station_z) + + # Rotate the azimuth and elevation angles off the center beam the new + # transformed cartesian coordinates + r = 1 + # transform pointing vectors, without considering geodesical earth + # coord system + pointing_vec_x, pointing_vec_y, pointing_vec_z = polar_to_cartesian( + r, all_azimuth, all_elevation) + pointing_vec_x, pointing_vec_y, pointing_vec_z = \ + coordinate_system.ecef2enu( + pointing_vec_x, pointing_vec_y, pointing_vec_z, translate=0) + _, all_azimuth, all_elevation = cartesian_to_polar( + pointing_vec_x, pointing_vec_y, pointing_vec_z) - # Convert the ECEF coordinates to the transformed cartesian coordinates and set the Space Station positions - # used to generetate the IMT Base Stations - space_station_x, space_station_y, space_station_z = \ - coordinate_system.ecef2enu(space_station_x, space_station_y, space_station_z) - - # Rotate the azimuth and elevation angles off the center beam the new - # transformed cartesian coordinates - r = 1 - # transform pointing vectors, without considering geodesical earth - # coord system - pointing_vec_x, pointing_vec_y, pointing_vec_z = polar_to_cartesian( - r, all_azimuth, all_elevation) - pointing_vec_x, pointing_vec_y, pointing_vec_z = \ - coordinate_system.ecef2enu( - pointing_vec_x, pointing_vec_y, pointing_vec_z, translate=0) - _, all_azimuth, all_elevation = cartesian_to_polar( - pointing_vec_x, pointing_vec_y, pointing_vec_z) - - beams_elev, beams_azim, sx, sy = TopologyImtMssDc.get_satellite_pointing( - random_number_gen, - coordinate_system, - orbit_params, - total_active_satellites, - all_space_station_x, all_space_station_y, all_space_station_z, - all_azimuth, - all_elevation, - all_lat, all_lon, all_sat_altitude, - active_satellite_idxs - ) + beams_elev, beams_azim, sx, sy = TopologyImtMssDc.get_satellite_pointing( + random_number_gen, + coordinate_system, + orbit_params, + total_active_satellites, + all_space_station_x, all_space_station_y, all_space_station_z, + all_azimuth, + all_elevation, + all_lat, all_lon, all_sat_altitude, + active_satellite_idxs + ) - # In SHARC each sector is treated as a separate base station, so we need to repeat the satellite positions - # for each sector. - sat_ocurr = [len(x) for x in beams_elev] - - elevation = np.array( - functools.reduce( - lambda x, - y: list(x) + - list(y), - beams_elev)) - azimuth = np.array( - functools.reduce( - lambda x, - y: list(x) + - list(y), - beams_azim)) - - space_station_x = np.repeat(space_station_x, sat_ocurr) - space_station_y = np.repeat(space_station_y, sat_ocurr) - space_station_z = np.repeat(space_station_z, sat_ocurr) - - num_base_stations = np.sum(sat_ocurr) - lat = np.repeat(lat, sat_ocurr) - lon = np.repeat(lon, sat_ocurr) - - altitudes = np.repeat(sat_altitude, sat_ocurr) + # In SHARC each sector is treated as a separate base station, so we need to repeat the satellite positions + # for each sector. + sat_ocurr = [len(x) for x in beams_elev] + + elevation = np.array( + functools.reduce( + lambda x, + y: list(x) + + list(y), + beams_elev)) + azimuth = np.array( + functools.reduce( + lambda x, + y: list(x) + + list(y), + beams_azim)) + + space_station_x = np.repeat(space_station_x, sat_ocurr) + space_station_y = np.repeat(space_station_y, sat_ocurr) + space_station_z = np.repeat(space_station_z, sat_ocurr) + + num_base_stations = np.sum(sat_ocurr) + lat = np.repeat(lat, sat_ocurr) + lon = np.repeat(lon, sat_ocurr) + + altitudes = np.repeat(sat_altitude, sat_ocurr) assert (space_station_x.shape == (num_base_stations,)) assert (space_station_y.shape == (num_base_stations,)) @@ -408,7 +411,7 @@ def get_satellite_pointing( eligible_sats_msk &= polygon_mask eligible_sats_idx = np.where(eligible_sats_msk)[0] - elev = calc_elevation( + all_elevations = calc_elevation( grid_lat[:, np.newaxis], all_sat_lat[eligible_sats_msk][np.newaxis, :], grid_lon[:, np.newaxis], @@ -417,7 +420,7 @@ def get_satellite_pointing( es_height=0, ) - best_sats = elev.argmax(axis=-1) + best_sats = all_elevations.argmax(axis=-1) best_sats_true = eligible_sats_idx[best_sats] @@ -450,8 +453,12 @@ def get_satellite_pointing( sat_points_towards = defaultdict(list) - for i, sat in enumerate(best_sats_true): - sat_points_towards[sat].append(i) + min_service_angle = orbit_params.beam_positioning.service_grid.minimum_service_angle + for i, sat in enumerate(best_sats): + if all_elevations[i][sat] < min_service_angle: + continue + actual_best_sat_idx = eligible_sats_idx[sat] + sat_points_towards[actual_best_sat_idx].append(i) # now only return the angles that # the caller asked with the active_sat_idxs parameter diff --git a/tests/parameters/parameters_for_testing.yaml b/tests/parameters/parameters_for_testing.yaml index 1cd0d9fc..86631172 100644 --- a/tests/parameters/parameters_for_testing.yaml +++ b/tests/parameters/parameters_for_testing.yaml @@ -264,6 +264,7 @@ imt: # if positive, makes border smaller by x km # if negative, makes border bigger by x km eligible_sats_margin_from_border: -2.1 + minimum_service_angle: 10.1 sat_is_active_if: # for a satellite to be active, it needs to respect ALL conditions conditions: @@ -998,6 +999,9 @@ single_space_station: # Antenna diameter [m] # if no diameter is passed, a diameter will be assumed according to document diameter: 2.12 + itu_r_s_1528: + major_minor_axis_ratio: 1.0 + far_out_side_lobe: -25 ########################################################################### # Selected channel model channel_model: "P619" @@ -1096,6 +1100,7 @@ mss_d2d: # if positive, makes border smaller by x km # if negative, makes border bigger by x km eligible_sats_margin_from_border: -2.1 + minimum_service_angle: 12.1 sat_is_active_if: # for a satellite to be active, it needs to respect ALL conditions conditions: diff --git a/tests/parameters/test_parameters.py b/tests/parameters/test_parameters.py index e174b53d..e5fd967d 100644 --- a/tests/parameters/test_parameters.py +++ b/tests/parameters/test_parameters.py @@ -313,6 +313,8 @@ def test_parameters_imt(self): 0.11) self.assertEqual( self.parameters.imt.topology.mss_dc.beam_positioning.service_grid.eligible_sats_margin_from_border, -2.1) + self.assertEqual( + self.parameters.imt.topology.mss_dc.beam_positioning.service_grid.minimum_service_angle, 10.1) self.assertEqual(len( self.parameters.imt.topology.mss_dc.beam_positioning.service_grid.grid_in_zone.from_countries.country_names), 2) self.assertEqual( @@ -709,6 +711,8 @@ def test_parametes_mss_d2d(self): 0.11) self.assertEqual( self.parameters.mss_d2d.beam_positioning.service_grid.eligible_sats_margin_from_border, -2.1) + self.assertEqual( + self.parameters.mss_d2d.beam_positioning.service_grid.minimum_service_angle, 12.1) self.assertEqual( len(self.parameters.mss_d2d.beam_positioning.service_grid.grid_in_zone.from_countries.country_names), 2) self.assertEqual( @@ -871,6 +875,20 @@ def test_parameters_single_space_station(self): self.parameters.single_space_station.antenna.gain, ) + self.assertEqual( + self.parameters.single_space_station.antenna.itu_r_s_1528.antenna_gain, + self.parameters.single_space_station.antenna.gain, + ) + self.assertEqual( + self.parameters.single_space_station.antenna.itu_r_s_1528.major_minor_axis_ratio, + 1.0 + ) + + self.assertEqual( + self.parameters.single_space_station.antenna.itu_r_s_1528.far_out_side_lobe, + -25 + ) + self.assertEqual( self.parameters.single_space_station.param_p619.earth_station_alt_m, self.parameters.single_space_station.geometry.es_altitude, diff --git a/tests/test_antenna_s1528.py b/tests/test_antenna_s1528.py index 0d2aae6c..b025506b 100644 --- a/tests/test_antenna_s1528.py +++ b/tests/test_antenna_s1528.py @@ -23,6 +23,8 @@ def setUp(self): param.antenna_gain = 39 param.antenna_pattern = "ITU-R S.1528-0" param.antenna_3_dB_bw = 2 + param.major_minor_axis_ratio = 1 + param.far_out_side_lobe = 0 param.antenna_l_s = -20 self.antenna20 = AntennaS1528(param) diff --git a/tests/test_topology_imt_mss_dc.py b/tests/test_topology_imt_mss_dc.py index 481e7ac8..69a715ca 100644 --- a/tests/test_topology_imt_mss_dc.py +++ b/tests/test_topology_imt_mss_dc.py @@ -6,6 +6,7 @@ from sharc.station_manager import StationManager from sharc.parameters.parameters_orbit import ParametersOrbit from sharc.support.sharc_geom import CoordinateSystem, lla2ecef +from sharc.satellite.utils.sat_utils import calc_elevation class TestTopologyImtMssDc(unittest.TestCase): @@ -156,6 +157,51 @@ def test_visible_satellites(self): # oblateness npt.assert_array_less(min_elevation_angle, xy_plane_elevations) + def test_minimum_service_angle(self): + """Test minimum visibility angle for service grid service.""" + orbit = ParametersOrbit( + n_planes=1, + sats_per_plane=1, + phasing_deg=3.9, + long_asc_deg=18.0, + inclination_deg=54.5, + perigee_alt_km=525, + apogee_alt_km=525, + ) + self.coordinate_system.set_reference(0.0, 0.0, 0) + + self.params.orbits = [orbit] + self.params.beam_positioning.type = "SERVICE_GRID" + self.params.beam_positioning.service_grid.grid_in_zone.type = "CIRCLE" + self.params.beam_positioning.service_grid.grid_in_zone.circle.center_lat = 0.0 + self.params.beam_positioning.service_grid.grid_in_zone.circle.center_lon = 0.0 + self.params.beam_positioning.service_grid.grid_in_zone.circle.radius_km = 10 * 111.0 + self.params.validate("") + + self.imt_mss_dc_topology = TopologyImtMssDc( + self.params, self.coordinate_system) + + n_previous_selected = np.inf + for a in [5, 50, 80]: + self.imt_mss_dc_topology.orbit_params.beam_positioning.service_grid.minimum_service_angle = a + self.imt_mss_dc_topology.calculate_coordinates( + random_number_gen=np.random.RandomState(8)) + + lon_lat_grid = self.params.beam_positioning.service_grid.lon_lat_grid + elev_from_bs = calc_elevation( + lon_lat_grid[1], + self.imt_mss_dc_topology.lat[0], + lon_lat_grid[0], + self.imt_mss_dc_topology.lon[0], + sat_height=self.imt_mss_dc_topology.height[0], + es_height=0.0, + ) + n_selected = np.sum(elev_from_bs >= a) + self.assertLess(n_selected, n_previous_selected) + self.assertLess(n_selected, len(elev_from_bs)) + self.assertEqual(n_selected, self.imt_mss_dc_topology.num_base_stations) + n_previous_selected = n_selected + if __name__ == '__main__': unittest.main()