Feat/pfd computation

.vscode/launch.json

@@ -9,7 +9,7 @@
             "console": "integratedTerminal",
             "args": [
                 "-p",
-                "${workspaceFolder}/sharc/input/parameters.yaml"
+                "${workspaceFolder}/sharc/campaigns/ast_mss_d2d_to_imt_cpe/input/parameter_ast_mss_d2d_to_imt_cpe_30exclusion_45beam_elev_0.0power_backoff_20load_factor_P619prop_imt-ue_imt.upto-1GHz.single-bs.urban-macro-bs_system-4.698-960MHz-block2.690km-antenna-update.yaml"
             ]
         }
     ]

sharc/antenna/antenna_f1245_fs.py

@@ -0,0 +1,196 @@
+
+# -*- coding: utf-8 -*-
+"""
+Created on Wed Apr 4 17:08:00 2018
+@author: Calil
+"""
+import matplotlib.pyplot as plt
+from sharc.antenna.antenna import Antenna
+from sharc.parameters.imt.parameters_imt import ParametersImt
+from sharc.parameters.imt.parameters_imt import ParametersImt
+from sharc.parameters.imt.parameters_antenna_imt import ParametersAntennaImt
+from sharc.parameters.parameters_antenna import ParametersAntenna
+import numpy as np
+import math
+
+
+class Antenna_f1245_fs(Antenna):
+    """Class that implements the ITU-R F.1245 antenna pattern for fixed
+    satellite service earth stations."""
+
+    def __init__(self, param: ParametersImt):
+        super().__init__()
+        self.peak_gain = param.gain
+        lmbda = 3e8 / (param.frequency * 1e6)
+        self.d_lmbda = param.diameter / lmbda
+        self.g_l = 2 + 15 * math.log10(self.d_lmbda)
+        self.phi_m = (20 / self.d_lmbda) * math.sqrt(self.peak_gain - self.g_l)
+        self.phi_r = 12.02 * math.pow(self.d_lmbda, -0.6)
+
+    def calculate_gain(self, *args, **kwargs) -> np.array:
+        """
+        Calculate the antenna gain for the given parameters.
+
+        Parameters
+        ----------
+        *args : tuple
+            Positional arguments (not used).
+        **kwargs : dict
+            Keyword arguments, expects 'phi_vec', 'theta_vec', and 'beams_l'.
+
+        Returns
+        -------
+        np.array
+            Calculated antenna gain values.
+        """
+        # phi_vec = np.absolute(kwargs["phi_vec"])
+        # theta_vec = np.absolute(kwargs["theta_vec"])
+        # beams_l = np.absolute(kwargs["beams_l"])
+        off_axis = np.absolute(kwargs["off_axis_angle_vec"])
+        if self.d_lmbda > 100:
+            gain = self.calculate_gain_greater(off_axis)
+        else:
+            gain = self.calculate_gain_less(off_axis)
+            # idx_max_gain = np.where(beams_l == -1)[0]
+            # gain = self.peak_gain
+        return gain
+
+    def calculate_gain_greater(self, phi: float) -> np.array:
+        """
+        For frequencies in the range 1 GHz to about 70 GHz, in cases where the
+        ratio between the antenna diameter and the wavelength is GREATER than
+        100, this method should be used.
+        Parameter
+        ---------
+        phi : off-axis angle [deg]
+        Returns
+        -------
+        a numpy array containing the gains in the given angles
+        """
+        gain = np.zeros(phi.shape)
+        idx_0 = np.where(phi < self.phi_m)[0]
+        gain[idx_0] = self.peak_gain - 2.5e-3 * \
+            np.power(self.d_lmbda * phi[idx_0], 2)
+        phi_thresh = max(self.phi_m, self.phi_r)
+        idx_1 = np.where((self.phi_m <= phi) & (phi < phi_thresh))[0]
+        gain[idx_1] = self.g_l
+        idx_2 = np.where((phi_thresh <= phi) & (phi < 48))[0]
+        gain[idx_2] = 29 - 25 * np.log10(phi[idx_2])
+        idx_3 = np.where((48 <= phi) & (phi <= 180))[0]
+        gain[idx_3] = -13
+        return gain
+
+    def calculate_gain_less(self, phi: float) -> np.array:
+        """
+        For frequencies in the range 1 GHz to about 70 GHz, in cases where the
+        ratio between the antenna diameter and the wavelength is LESS than
+        or equal to 100, this method should be used.
+        Parameter
+        ---------
+        phi : off-axis angle [deg]
+        Returns
+        -------
+        a numpy array containing the gains in the given angles
+        """
+        gain = np.zeros(phi.shape)
+        idx_0 = np.where(phi < self.phi_m)[0]
+        gain[idx_0] = self.peak_gain - 0.0025 * \
+            np.power(self.d_lmbda * phi[idx_0], 2)
+        idx_1 = np.where((self.phi_m <= phi) & (phi < 48))[0]
+        gain[idx_1] = 39 - 5 * \
+            math.log10(self.d_lmbda) - 25 * np.log10(phi[idx_1])
+        idx_2 = np.where((48 <= phi) & (phi < 180))[0]
+        gain[idx_2] = -3 - 5 * math.log10(self.d_lmbda)
+        return gain
+
+    def add_beam(self, phi: float, theta: float):
+        """
+        Add a new beam to the antenna.
+
+        Parameters
+        ----------
+        phi : float
+            Azimuth angle in degrees.
+        theta : float
+            Elevation angle in degrees.
+        """
+        self.beams_list.append((phi, theta))
+
+    def calculate_off_axis_angle(self, Az, b):
+        """
+        Calculate the off-axis angle between the main beam and a given direction.
+
+        Parameters
+        ----------
+        Az : float or np.array
+            Azimuth angle(s) in degrees.
+        b : float or np.array
+            Elevation angle(s) in degrees.
+
+        Returns
+        -------
+        float or np.array
+            Off-axis angle(s) in degrees.
+        """
+        Az0 = self.beams_list[0][0]
+        a = 90 - self.beams_list[0][1]
+        C = Az0 - Az
+        off_axis_rad = np.arccos(
+            np.cos(
+                np.radians(a)) *
+            np.cos(
+                np.radians(b)) +
+            np.sin(
+                np.radians(a)) *
+            np.sin(
+                np.radians(b)) *
+            np.cos(
+                np.radians(C)),
+        )
+        off_axis_deg = np.degrees(off_axis_rad)
+        return off_axis_deg
+
+
+if __name__ == '__main__':
+    off_axis_angle_vec = np.linspace(0.1, 180, num=1001)
+    # initialize antenna parameters
+    param = ParametersAntenna()
+    param.frequency = 2155
+    param_gt = ParametersAntennaImt()
+    param.gain = 33.1
+    param.diameter = 2
+    antenna_gt = Antenna_f1245_fs(param)
+    antenna_gt.add_beam(0, 0)
+    gain_gt = antenna_gt.calculate_gain(
+        off_axis_angle_vec=off_axis_angle_vec,
+    )
+    param.diameter = 3
+    antenna_gt = Antenna_f1245_fs(param)
+    gain_gt_3 = antenna_gt.calculate_gain(
+        off_axis_angle_vec=off_axis_angle_vec,
+    )
+    param.diameter = 1.8
+    antenna_gt = Antenna_f1245_fs(param)
+    gain_gt_18 = antenna_gt.calculate_gain(
+        off_axis_angle_vec=off_axis_angle_vec,
+    )
+
+    fig = plt.figure(
+        figsize=(8, 7), facecolor='w',
+        edgecolor='k',
+    )  # create a figure object
+    plt.semilogx(off_axis_angle_vec, gain_gt, "-b", label="$f = 10.7$ $GHz,$ $D = 2$ $m$")
+    plt.semilogx(off_axis_angle_vec, gain_gt_3, "-y", label="$f = 10.7$ $GHz,$ $D = 3$ $m$")
+    plt.semilogx(off_axis_angle_vec, gain_gt_18, "-g", label="$f = 10.7$ $GHz,$ $D = 1.8$ $m$")
+
+    plt.title("ITU-R F.1245 antenna radiation pattern")
+    plt.xlabel(r"Off-axis angle $\phi$ [deg]")
+    plt.ylabel("Gain relative to $G_m$ [dB]")
+    plt.legend(loc="lower left")
+    # plt.xlim((phi[0], phi[-1]))
+    plt.ylim((-20, 50))
+    # ax = plt.gca()
+    # ax.set_yticks([-30, -20, -10, 0])
+    # ax.set_xticks(np.linspace(1, 9, 9).tolist() + np.linspace(10, 100, 10).tolist())
+    plt.grid()
+    plt.show()

sharc/antenna/antenna_factory.py

@@ -12,6 +12,7 @@
 from sharc.antenna.antenna_s580 import AntennaS580
 from sharc.antenna.antenna_s1528 import AntennaS1528
 from sharc.antenna.antenna_s1855 import AntennaS1855
+from sharc.antenna.antenna_f1245_fs import Antenna_f1245_fs
 from sharc.antenna.antenna_s1528 import AntennaS1528, AntennaS1528Leo, AntennaS1528Taylor
 from sharc.antenna.antenna_beamforming_imt import AntennaBeamformingImt
 
@@ -49,6 +50,8 @@ def create_antenna(
                 return AntennaS465(antenna_params.itu_r_s_465_modified)
             case "ITU-R S.1855":
                 return AntennaS1855(antenna_params.itu_r_s_1855)
+            case "ITU-R F.1245_fs":
+                return Antenna_f1245_fs(antenna_params.itu_r_f_1245_fs)
             case "ITU-R Reg. RR. Appendice 7 Annex 3":
                 return AntennaReg_RR_A7_3(antenna_params.itu_reg_rr_a7_3)
             case "MSS Adjacent":
@@ -60,6 +63,12 @@ def create_antenna(
                     azimuth,
                     elevation
                 )
+            case "Antenna System 4":
+                # Handled in station_factory.py since it requires two antennas
+                raise NotImplementedError(
+                    "Antenna System 4 requires two antenna instances and "
+                    "should be created in station_factory.py"
+                )
             case _:
                 raise ValueError(
                     f"Antenna factory does not support pattern {
@@ -88,6 +97,14 @@ def create_n_antennas(
                 antennas[i] = AntennaFactory.create_antenna(
                     antenna_params, azimuth[i], elevation[i],
                 )
+        elif antenna_params.pattern == "Antenna System 4":
+            antenna_pattern_high = AntennaS1528(antenna_params.antenna_system_4.antenna_parameters_high)
+            antenna_pattern_low = AntennaS1528(antenna_params.antenna_system_4.antenna_parameters_low)
+            for i in range(n_stations):
+                if antenna_params.antenna_system_4.beam_ground_elev_angles[i] >= 50:
+                    antennas[i] = antenna_pattern_high
+                else:
+                    antennas[i] = antenna_pattern_low
         else:
             # some antennas don't need azimuth and elevation at all
             # this makes it much faster

sharc/parameters/antenna/parameters_antenna_system4.py

@@ -0,0 +1,52 @@
+# Implementation of ParametersAntennaSystem4 class
+# The System 4 antenna was defined in WP4C Working Document 4C/356 from October 2025
+# The antenna is based on S.1528 recommends 1.2
+# It defines two sets of parameters for high elevation and low elevation beams
+from dataclasses import dataclass, field
+
+from sharc.parameters.parameters_base import ParametersBase
+from sharc.parameters.antenna.parameters_antenna_s1528 import ParametersAntennaS1528
+
+
+@dataclass
+class ParametersAntennaSystem4(ParametersBase):
+    """Dataclass containing the Antenna System 4 parameters for the simulator.
+    """
+    section_name: str = "Antenna System 4"
+
+    # Parameters for high elevation beams
+    antenna_parameters_high: ParametersAntennaS1528 = field(
+        # we don't care about frequency and bandwidth here. Just to make validation work.
+        default_factory=lambda: ParametersAntennaS1528(
+            frequency=-1,
+            bandwidth=-1
+        )
+    )
+
+    # Parameters for low elevation beams
+    antenna_parameters_low: ParametersAntennaS1528 = field(
+        # we don't care about frequency and bandwidth here. Just to make validation work.
+        default_factory=lambda: ParametersAntennaS1528(
+            frequency=-1,
+            bandwidth=-1
+        )
+    )
+
+    # Elevation angles (degrees) separating high and low elevation beams
+    # This is set during simulation runtime.
+    beam_ground_elev_angles: list[float] = None
+
+    def load_parameters_from_file(self, config_file: str):
+        """Load the parameters from file an run a sanity check.
+
+        Parameters
+        ----------
+        file_name : str
+            the path to the configuration file
+
+        Raises
+        ------
+        ValueError
+            if a parameter is not valid
+        """
+        super().load_parameters_from_file(config_file)

sharc/parameters/imt/parameters_imt_mss_dc.py

@@ -6,7 +6,11 @@
 
 from sharc.parameters.parameters_base import ParametersBase
 from sharc.parameters.parameters_orbit import ParametersOrbit
-from sharc.parameters.imt.parameters_grid import ParametersSatelliteWithServiceGrid, ParametersSelectActiveSatellite
+from sharc.parameters.imt.parameters_grid import (
+    ParametersSatelliteWithServiceGrid,
+    ParametersSelectActiveSatellite,
+    ParametersZone,
+)
 
 SHARC_ROOT_DIR = (Path(__file__) / ".." / ".." / ".." / "..").resolve()
 
@@ -189,6 +193,30 @@ def validate(self, ctx):
                 )
 
 
+@dataclass
+class ParametersPowerControlZone(ParametersBase):
+    """Dataclass for a power control zone in the IMT MSS-DC topology."""
+    geometry: ParametersZone = field(default_factory=ParametersZone)
+    power_backoff_db: float = None
+
+
+@dataclass
+class ParametersPowerControl(ParametersBase):
+    """Dataclass for power control parameters in the IMT MSS-DC topology."""
+    zones: list[ParametersPowerControlZone] = field(
+        default_factory=lambda: [ParametersPowerControlZone()])
+
+    def validate(self, ctx):
+        """
+        Validate the power control parameters.
+        """
+        super().validate(ctx)
+        for i in range(len(self.zones)):
+            self.zones[i].geometry.validate(ctx + f"zones.{i}.geometry")
+            if not isinstance(self.zones[i].power_backoff_db, float):
+                self.zones[i].power_backoff_db = 0.0  # Default to 0 dB if not set - we make sure it's not applied later.
+
+
 @dataclass
 class ParametersImtMssDc(ParametersBase):
     """Dataclass for the IMT MSS-DC topology parameters."""
@@ -211,6 +239,9 @@ class ParametersImtMssDc(ParametersBase):
     # for IMT Space Stations
     beam_radius: float = 36516.0
 
+    power_control_zones: ParametersPowerControl = field(
+        default_factory=ParametersPowerControl)
+
     sat_is_active_if: ParametersSelectActiveSatellite = field(
         default_factory=ParametersSelectActiveSatellite)