anomaly-detection-material-parameters-calibration

Sionna param calibration (research proj)
git clone https://git.ea.contact/anomaly-detection-material-parameters-calibration
Log | Files | Refs | README

ris.py (60020B)


      1 #
      2 # SPDX-FileCopyrightText: Copyright (c) 2021-2025 NVIDIA CORPORATION & AFFILIATES. All rights reserved.
      3 # SPDX-License-Identifier: Apache-2.0
      4 #
      5 """
      6 Classes and functions relating to reconfigurable intelligent surfaces
      7 """
      8 
      9 from abc import ABC
     10 from abc import abstractmethod
     11 import tensorflow as tf
     12 import matplotlib.pyplot as plt
     13 
     14 from .radio_device import RadioDevice
     15 from .scene_object import SceneObject
     16 from . import scene
     17 from .utils import rotate, normalize, outer,\
     18                    expand_to_rank
     19 
     20 class CellGrid():
     21     # pylint: disable=line-too-long
     22     r"""
     23     Class defining a cell grid that determines the physical structure of a RIS
     24 
     25     The cell grid specifies the location of unit cells within the y-z plane
     26     assuming a homogenous spacing of 0.5. The actual positions are computed by
     27     multiplying the cell positions by the wavelength and rotating them
     28     according to the RIS' orientation.
     29 
     30     A cell grid must have at least three columns and rows to ensure
     31     that discrete phase and amplitude profiles of the RIS can be interpolated.
     32 
     33     Parameters
     34     ----------
     35     num_rows : int
     36         Number of rows. Must at least be equal to three.
     37 
     38     num_cols : int
     39         Number of columns. Must at least be equal to three.
     40 
     41     dtype : tf.complex
     42         Datatype to be used in internal calculations.
     43         Defaults to `tf.complex64`.
     44 
     45     """
     46     def __init__(self,
     47                  num_rows,
     48                  num_cols,
     49                  dtype=tf.complex64):
     50 
     51         if dtype not in (tf.complex64, tf.complex128):
     52             raise ValueError("`dtype` must be tf.complex64 or tf.complex128`")
     53         self._dtype = dtype
     54         self._rdtype = dtype.real_dtype
     55 
     56         if num_rows < 3 or num_cols < 3:
     57             raise ValueError("num_rows and num_cols must be >= 3")
     58         self._num_rows = int(num_rows)
     59         self._num_cols = int(num_cols)
     60 
     61         self._cell_y_positions = tf.range(self.num_cols, dtype=self._rdtype)
     62         self._cell_y_positions -= tf.cast((self.num_cols-1.)/2., self._rdtype)
     63 
     64         self._cell_z_positions = tf.range(self.num_rows-1, -1, -1,
     65                                           dtype=self._rdtype)
     66         self._cell_z_positions -= tf.cast((self.num_rows-1.)/2., self._rdtype)
     67 
     68         z, y = tf.meshgrid(self.cell_z_positions, self.cell_y_positions)
     69         self._cell_positions = tf.stack([tf.reshape(y, [-1]),
     70                                          tf.reshape(z, [-1])], -1)
     71 
     72     @property
     73     def num_rows(self):
     74         r"""
     75         int : Number of rows
     76         """
     77         return self._num_rows
     78 
     79     @property
     80     def num_cols(self):
     81         r"""
     82         int : Number of columns
     83         """
     84         return self._num_cols
     85 
     86     @property
     87     def num_cells(self):
     88         r"""
     89         int : Number of cells
     90         """
     91         return self.num_rows * self.num_cols
     92 
     93     @property
     94     def cell_positions(self):
     95         r"""
     96         [num_cells, 2], tf.float : Cell positions ordered from
     97             top-to-bottom left-to-right
     98         """
     99         return self._cell_positions
    100 
    101     @property
    102     def cell_y_positions(self):
    103         r"""
    104         [num_cols], tf.float : y-coordinates of cells ordered
    105             from left-to-right
    106         """
    107         return self._cell_y_positions
    108 
    109     @property
    110     def cell_z_positions(self):
    111         r"""
    112         [num_rows], tf.float : z-coordinates of cells ordered
    113             from top-to-bottom
    114         """
    115         return self._cell_z_positions
    116 
    117 class Profile(ABC):
    118     # pylint: disable=line-too-long
    119     r"""Abstract class defining a phase/amplitude profile of a RIS
    120 
    121     A Profile instance is a callable that returns the profile values,
    122     gradients and Hessians at given points.
    123 
    124     Parameters
    125     ----------
    126     dtype : tf.complex
    127         Datatype to be used in internal calculations.
    128         Defaults to `tf.complex64`.
    129 
    130     Input
    131     -----
    132     points : tf.float, [num_samples, 2]
    133         Tensor of 2D coordinates defining the points on the RIS at which
    134         the profile should be evaluated.
    135         Defaults to `None`. In this case, the values for all unit cells
    136         are returned.
    137 
    138     mode : int | `None`
    139         Reradiation mode to be considered.
    140         Defaults to `None`. In this case, the values for all modes
    141         are returned.
    142 
    143     return_grads : bool
    144         If `True`, also the first- and second-order derivatives are
    145         returned.
    146         Defaults to `False`.
    147 
    148     Output
    149     ------
    150     values : [num_modes, num_samples] or [num_samples], tf.float
    151         Interpolated profile values at the sample positions
    152 
    153     grads : [num_modes, num_samples, 3] or [num_samples, 3], tf.float
    154         Gradients of the interpolated profile values
    155         at the sample positions. Only returned if `return_grads` is `True`.
    156 
    157     hessians : [num_modes, num_samples, 3, 3] or [num_samples, 3, 3] , tf.float
    158         Hessians of the interpolated profile values
    159         at the sample positions. Only returned if `return_grads` is `True`.
    160     """
    161     def __init__(self, dtype=tf.complex64):
    162         if dtype not in (tf.complex64, tf.complex128):
    163             raise ValueError("`dtype` must be tf.complex64 or tf.complex128`")
    164         self._dtype = dtype
    165         self._rdtype = dtype.real_dtype
    166 
    167     @property
    168     @abstractmethod
    169     def num_modes(self):
    170         r"""
    171         int : Number of reradiation modes
    172         """
    173         pass
    174 
    175     @abstractmethod
    176     def __call__(self, points, mode=None, return_grads=False):
    177         r"""
    178         Returns the profile values, gradients and Hessians at given points
    179 
    180         Input
    181         -----
    182         points : tf.float, [num_samples, 2]
    183             Tensor of 2D coordinates defining the points on the RIS at which
    184             the profile should be evaluated.
    185             Defaults to `None`. In this case, the values for all unit cells
    186             are returned.
    187 
    188         mode : int | `None`
    189             Reradiation mode to be considered.
    190             Defaults to `None`. In this case, the values for all modes
    191             are returned.
    192 
    193         return_grads : bool
    194             If `True`, also the first- and second-order derivatives are
    195             returned.
    196             Defaults to `False`.
    197 
    198         Output
    199         ------
    200         values : [num_modes, num_samples] or [num_samples], tf.float
    201             Interpolated profile values at the sample positions
    202 
    203         grads : [num_modes, num_samples, 3] or [num_samples, 3], tf.float
    204             Gradients of the interpolated profile values
    205             at the sample positions. Only returned if `return_grads` is `True`.
    206 
    207         hessians : [num_modes, num_samples, 3, 3] or [num_samples, 3, 3] , tf.float
    208             Hessians of the interpolated profile values
    209             at the sample positions. Only returned if `return_grads` is `True`.
    210         """
    211         pass
    212 
    213 class AmplitudeProfile(Profile):
    214     # pylint: disable=line-too-long
    215     r"""Abstract class defining an amplitude profile of a RIS
    216 
    217     An AmplitudeProfile instance is a callable that returns the profile values,
    218     gradients and Hessians at given points.
    219 
    220     Parameters
    221     ----------
    222     dtype : tf.complex
    223         Datatype to be used in internal calculations.
    224         Defaults to `tf.complex64`.
    225 
    226     Input
    227     -----
    228     points : tf.float, [num_samples, 2]
    229         Tensor of 2D coordinates defining the points on the RIS at which
    230         the profile should be evaluated.
    231         Defaults to `None`. In this case, the values for all unit cells
    232         are returned.
    233 
    234     mode : int | `None`
    235         Reradiation mode to be considered.
    236         Defaults to `None`. In this case, the values for all modes
    237         are returned.
    238 
    239     return_grads : bool
    240         If `True`, also the first- and second-order derivatives are
    241         returned.
    242         Defaults to `False`.
    243 
    244     Output
    245     ------
    246     values : [num_modes, num_samples] or [num_samples], tf.float
    247         Interpolated profile values at the sample positions
    248 
    249     grads : [num_modes, num_samples, 3] or [num_samples, 3], tf.float
    250         Gradients of the interpolated profile values
    251         at the sample positions. Only returned if `return_grads` is `True`.
    252 
    253     hessians : [num_modes, num_samples, 3, 3] or [num_samples, 3, 3] , tf.float
    254         Hessians of the interpolated profile values
    255         at the sample positions. Only returned if `return_grads` is `True`.
    256     """
    257     @property
    258     @abstractmethod
    259     def mode_powers(self):
    260         r"""
    261         [num_modes], tf.float: Relative power of reradiation modes
    262         """
    263         pass
    264 
    265 class PhaseProfile(Profile):
    266     # pylint: disable=line-too-long
    267     r"""Abstract class defining a phase profile of a RIS
    268 
    269     A PhaseProfile instance is a callable that returns the profile values,
    270     gradients and Hessians at given points.
    271 
    272     Parameters
    273     ----------
    274     dtype : tf.complex
    275         Datatype to be used in internal calculations.
    276         Defaults to `tf.complex64`.
    277 
    278     Input
    279     -----
    280     points : tf.float, [num_samples, 2]
    281         Tensor of 2D coordinates defining the points on the RIS at which
    282         the profile should be evaluated.
    283         Defaults to `None`. In this case, the values for all unit cells
    284         are returned.
    285 
    286     mode : int | `None`
    287         Reradiation mode to be considered.
    288         Defaults to `None`. In this case, the values for all modes
    289         are returned.
    290 
    291     return_grads : bool
    292         If `True`, also the first- and second-order derivatives are
    293         returned.
    294         Defaults to `False`.
    295 
    296     Output
    297     ------
    298     values : [num_modes, num_samples] or [num_samples], tf.float
    299         Interpolated profile values at the sample positions
    300 
    301     grads : [num_modes, num_samples, 3] or [num_samples, 3], tf.float
    302         Gradients of the interpolated profile values
    303         at the sample positions. Only returned if `return_grads` is `True`.
    304 
    305     hessians : [num_modes, num_samples, 3, 3] or [num_samples, 3, 3] , tf.float
    306         Hessians of the interpolated profile values
    307         at the sample positions. Only returned if `return_grads` is `True`.
    308     """
    309     pass
    310 
    311 class DiscreteProfile(Profile):
    312     # pylint: disable=line-too-long
    313     r"""Class defining a discrete phase/amplitude profile of a RIS
    314 
    315     A DiscreteProfile instance is a callable that returns the profile values,
    316     gradients and Hessians at given points.
    317 
    318     Parameters
    319     ----------
    320     cell_grid : :class:`~sionna.rt.CellGrid`
    321         Defines the physical structure of the RIS
    322 
    323     num_modes : int
    324         Number of reradiation modes.
    325         Defaults to 1.
    326 
    327     values : tf.float or tf.Variable, [num_modes, num_rows, num_cols]
    328         Values of the discrete profile for each reradiation mode
    329         and unit cell. `num_rows` and `num_cols` are defined by the
    330         `cell_grid`.
    331         Defaults to `None`.
    332 
    333     interpolator : :class:`~sionna.rt.ProfileInterpolator`
    334         Instance of a `ProfileInterpolator` that interpolates the
    335         discrete values of the profile to a continuous profile
    336         which is defined at any point on the RIS.
    337         Defaults to `None`. In this case, the
    338         :class:`~sionna.rt.LagrangeProfileInterpolator` will be used.
    339 
    340     dtype : tf.complex
    341         Datatype to be used in internal calculations.
    342         Defaults to `tf.complex64`.
    343 
    344     Input
    345     -----
    346     points : tf.float, [num_samples, 2]
    347         Tensor of 2D coordinates defining the points on the RIS at which
    348         the profile should be evaluated.
    349         Defaults to `None`. In this case, the values for all unit cells
    350         are returned.
    351 
    352     mode : int | `None`
    353         Reradiation mode to be considered.
    354         Defaults to `None`. In this case, the values for all modes
    355         are returned.
    356 
    357     return_grads : bool
    358         If `True`, also the first- and second-order derivatives are
    359         returned.
    360         Defaults to `False`.
    361 
    362     Output
    363     ------
    364     values : [num_modes, num_samples] or [num_samples], tf.float
    365         Interpolated profile values at the sample positions
    366 
    367     grads : [num_modes, num_samples, 3] or [num_samples, 3], tf.float
    368         Gradients of the interpolated profile values
    369         at the sample positions. Only returned if `return_grads` is `True`.
    370 
    371     hessians : [num_modes, num_samples, 3, 3] or [num_samples, 3, 3] , tf.float
    372         Hessians of the interpolated profile values
    373         at the sample positions. Only returned if `return_grads` is `True`.
    374     """
    375     def __init__(self,
    376                  cell_grid,
    377                  num_modes=1,
    378                  values=None,
    379                  interpolator=None,
    380                  dtype=tf.complex64):
    381 
    382         super().__init__(dtype=dtype)
    383         self._cell_grid = cell_grid
    384         self._num_modes = tf.cast(num_modes, tf.int32)
    385         if values is None:
    386             self._values = None
    387         else:
    388             self.values = values
    389         if interpolator is None:
    390             self._interpolator = LagrangeProfileInterpolator(self)
    391         else:
    392             self._interpolator = interpolator
    393 
    394     @property
    395     def shape(self):
    396         r"""
    397         tf.TensorShape : Shape of the tensor holding the values of
    398             the discrete profile
    399         """
    400         return tf.TensorShape([self.num_modes,
    401                                self.cell_grid.num_rows,
    402                                self.cell_grid.num_cols])
    403     @property
    404     def values(self):
    405         r"""
    406         [shape], tf.float : Set/get the discrete values of the profile for each
    407             reradiation mode
    408         """
    409         return self._values
    410 
    411     @values.setter
    412     def values(self, v):
    413         if not v.shape == self.shape:
    414             raise ValueError(f"`values` must have shape {self.shape}")
    415         if isinstance(v, tf.Variable):
    416             if v.dtype != self._rdtype:
    417                 msg = f"`values` must have dtype={self._rdtype}"
    418                 raise TypeError(msg)
    419             else:
    420                 self._values = v
    421         else:
    422             self._values = tf.cast(v, dtype=self._rdtype)
    423 
    424     @property
    425     def num_modes(self):
    426         r"""
    427         int : Number of reradiation modes
    428         """
    429         return self._num_modes
    430 
    431     @property
    432     def cell_grid(self):
    433         r"""
    434         :class:`~sionna.rt.CellGrid` : Defines the physical
    435             structure of the RIS
    436         """
    437         return self._cell_grid
    438 
    439     @property
    440     def spacing(self):
    441         r"""
    442         tf.float: Element spacing [m] corresponding to
    443             half a wavelength
    444         """
    445         if hasattr(scene.Scene(), "wavelength"):
    446             wavelength = tf.cast(scene.Scene().wavelength, self._rdtype)
    447             return wavelength/tf.cast(2, self._rdtype)
    448         else:
    449             # Scene is not initialized
    450             return tf.cast(0.5, self._rdtype)
    451 
    452     def show(self, mode=0):
    453         r"""Visualizes the profile as a 3D plot
    454 
    455         Input
    456         ------
    457         mode : int | `None`
    458             Reradation mode to be shown.
    459             Defaults to 0.
    460 
    461         Output
    462         ------
    463         : :class:`matplotlib.pyplot.Figure`
    464             3D plot of the profile
    465         """
    466         fig = plt.figure()
    467         ax = fig.add_subplot(111, projection='3d')
    468         y, z = tf.meshgrid(self.cell_grid.cell_y_positions*self.spacing,
    469                            self.cell_grid.cell_z_positions*self.spacing)
    470         ax.plot_surface(y, z, self.values[mode], cmap='viridis')
    471         ax.set_xlabel("y")
    472         ax.set_ylabel("z")
    473         if isinstance(self, PhaseProfile):
    474             plt.title(r"Phase profile $\chi(y, z)$")
    475         if isinstance(self, AmplitudeProfile):
    476             plt.title(r"Amplitude profile $A(y, z)$")
    477         return fig
    478 
    479     def __call__(self, points=None, mode=None, return_grads=False):
    480         r"""
    481         Returns the profile values, gradients and Hessians at given points
    482 
    483         Input
    484         -----
    485         points : tf.float, [num_samples, 2]
    486             Tensor of 2D coordinates defining the points on the RIS at which
    487             the profile should be evaluated.
    488             Defaults to `None`. In this case, the values for all unit cells
    489             are returned.
    490 
    491         mode : int | `None`
    492             Reradiation mode to be considered.
    493             Defaults to `None`. In this case, the values for all modes
    494             are returned.
    495 
    496         return_grads : bool
    497             If `True`, also the first- and second-order derivatives are
    498             returned. Only available if `points` is not `None`.
    499             Defaults to `False`.
    500 
    501         Output
    502         ------
    503         values : [num_modes, num_samples] or [num_samples], tf.float
    504             Interpolated profile values at the sample positions
    505 
    506         grads : [num_modes, num_samples, 3] or [num_samples, 3], tf.float
    507             Gradients of the interpolated profile values
    508             at the sample positions. Only returned if `return_grads` is `True`.
    509 
    510         hessians : [num_modes, num_samples, 3, 3] or [num_samples, 3, 3] , tf.float
    511             Hessians of the interpolated profile values
    512             at the sample positions. Only returned if `return_grads` is `True`.
    513         """
    514         if points is None:
    515             if mode is not None:
    516                 values = tf.transpose(self.values[mode])
    517                 values = tf.reshape(values, [-1])
    518             else:
    519                 values = tf.transpose(self.values, perm=[0,2,1])
    520                 values = tf.reshape(values, [self.num_modes, -1])
    521             return values
    522         else:
    523             return self._interpolator(points, mode, return_grads)
    524 
    525 class ProfileInterpolator(ABC):
    526     r"""
    527     Abstract class defining an interpolator of a discrete profile
    528 
    529     A ProfileInterpolator instance is a callable that interpolate
    530     the discrete profile to specified points. Optionally, the
    531     gradients and Hessians are returned.
    532 
    533     Parameters
    534     ----------
    535     discrete_profile : :class:`~sionna.rt.DiscreteProfile`
    536         Discrete profile to be interpolated
    537 
    538     Input
    539     -----
    540     points : [num_samples, 2], tf.float
    541         Positions at which to interpolate the profile
    542 
    543     mode : int | `None`
    544         Mode of the profile to interpolate. If `None`.
    545         all modes are interpolated.
    546         Defaults to `None`.
    547 
    548     return_grads : bool
    549         If `True`, gradients and Hessians are computed.
    550         Defaults to `False`.
    551 
    552     Output
    553     ------
    554     values : [num_modes, num_samples] or [num_samples], tf.float
    555         Interpolated profile values at the sample positions
    556 
    557     grads : [num_modes, num_samples, 3] or [num_samples, 3], tf.float
    558         Gradients of the interpolated profile values
    559         at the sample positions
    560 
    561     hessians : [num_modes, num_samples, 3, 3] or [num_samples,3,3], tf.float
    562         Hessians of the interpolated profile values
    563         at the sample positions
    564     """
    565     def __init__(self, discrete_profile):
    566         self._discrete_profile = discrete_profile
    567         self._dtype = discrete_profile._dtype
    568         self._rdtype = discrete_profile._rdtype
    569 
    570     @property
    571     def spacing(self):
    572         r"""
    573         tf.float: Element spacing [m] corresponding to
    574             half a wavelength
    575         """
    576         if hasattr(scene.Scene(), "wavelength"):
    577             wavelength = tf.cast(scene.Scene().wavelength, self._rdtype)
    578             return wavelength/tf.cast(2,  self._rdtype)
    579         else:
    580             # Scene is not initialized
    581             return tf.cast(0.5, self._rdtype)
    582 
    583     @property
    584     def cell_y_positions(self):
    585         r"""
    586         [num_cols], tf.float : y-coordinates of cells ordered
    587             from left-to-right
    588         """
    589         return self._discrete_profile.cell_grid.cell_y_positions*self.spacing
    590 
    591     @property
    592     def cell_z_positions(self):
    593         r"""
    594         [num_rows], tf.float : z-coordinates of cells ordered
    595             from top-to-bottom
    596         """
    597         return self._discrete_profile.cell_grid.cell_z_positions*self.spacing
    598 
    599     @property
    600     def num_rows(self):
    601         r"""
    602         int : Number of rows
    603         """
    604         return self._discrete_profile.cell_grid.num_rows
    605 
    606     @property
    607     def num_cols(self):
    608         r"""
    609         int : Number of columns
    610         """
    611         return self._discrete_profile.cell_grid.num_cols
    612 
    613     @property
    614     def values(self):
    615         r"""
    616         [shape], tf.float : Discrete values of the profile for each
    617             reradiation mode and unit cell
    618         """
    619         return self._discrete_profile.values
    620 
    621     @abstractmethod
    622     def __call__(self, points, mode=None, return_grads=False):
    623         r"""
    624         Interpolates the discrete profile to specified points
    625 
    626         Optionally, the gradients and Hessians are returned.
    627 
    628         Input
    629         -----
    630         points : [num_samples, 2], tf.float
    631             Positions at which to interpolate the profile
    632 
    633         mode : int | `None`
    634             Mode of the profile to interpolate. If `None`.
    635             all modes are interpolated.
    636             Defaults to `None`.
    637 
    638         return_grads : bool
    639             If `True`, gradients and Hessians are computed.
    640             Defaults to `False`.
    641 
    642         Output
    643         ------
    644         values : [num_modes, num_samples] or [num_samples], tf.float
    645             Interpolated profile values at the sample positions
    646 
    647         grads : [num_modes, num_samples, 3] or [num_samples, 3], tf.float
    648             Gradients of the interpolated profile values
    649             at the sample positions
    650 
    651         hessians : [num_modes, num_samples, 3, 3] or [num_samples,3,3], tf.float
    652             Hessians of the interpolated profile values
    653             at the sample positions
    654         """
    655         pass
    656 
    657 class LagrangeProfileInterpolator(ProfileInterpolator):
    658     # pylint: disable=line-too-long
    659     r"""
    660     Class defining a :class:`~sionna.rt.ProfileInterpolator` using Lagrange polynomials
    661 
    662     The class instance is a callable that interpolates a discrete profile
    663     at arbitrary positions using two-dimensional 2nd-order Lagrange interpolation.
    664 
    665     A discrete profile :math:`P(y_i,z_j)\in\mathbb{R}` defined on
    666     a grid of points :math:`y_i,z_j` for :math:`i,j \in [1,2,3]` is
    667     interpolated to position :math:`y,z` as
    668 
    669     .. math::
    670         \begin{align}
    671             P(y,z) &= \sum_{i,j} P(y_i,z_j) \ell_{i,y}(y) \ell_{j,z}(z)
    672         \end{align}
    673 
    674     where :math:`\ell_{i,y}(y)`, :math:`\ell_{j,z}(z)` are the
    675     one-dimensional 2nd-order Lagrange polynomials, defined
    676     as
    677 
    678     .. math::
    679         \begin{align}
    680             \ell_{i,y}(y) &= \prod_{j \ne i} \frac{y-y_j}{y_i-y_j} \\
    681             \ell_{j,z}(z) &= \prod_{i \ne j} \frac{z-z_i}{z_j-z_i}.
    682         \end{align}
    683 
    684     Note that the formulation above assumes for simplicity only a 3x3 grid
    685     of points. However, the implementation finds for every
    686     position the closest 3x3 grid points of the discrete profile
    687     which are used for interpolation.
    688 
    689     In order to compute spatial gradients and Hessians, we extend the the profile
    690     with a dummy :math:`x` dimension, i.e., :math:`P(x,y,z)=P(y,z)`, such that
    691 
    692     .. math::
    693         \begin{align}
    694             \nabla P(x,y,z) &= \begin{bmatrix} 0, \frac{\partial P(x,y,z)}{\partial y}, \frac{\partial P(x,y,z)}{\partial z}  \end{bmatrix}^{\textsf{T}}\\
    695             H_P(x,y,z) &= \begin{bmatrix} 0 & 0                                                 & 0 \\
    696                                             0 & \frac{\partial^2 P(x,y,z)}{\partial y^2}          & \frac{\partial^2 P(x,y,z)}{\partial y \partial z} \\
    697                                             0 & \frac{\partial^2 P(x,y,z)}{\partial z \partial y} & \frac{\partial^2 P(x,y,z)}{\partial z^2}
    698                             \end{bmatrix}
    699         \end{align}.
    700 
    701     Parameters
    702     ----------
    703     discrete_profile : :class:`~sionna.rt.DiscreteProfile`
    704         Discrete profile to be interpolated
    705 
    706     Input
    707     -----
    708     points : [num_samples, 2], tf.float
    709         Positions at which to interpolate the profile
    710 
    711     mode : int | `None`
    712         Mode of the profile to interpolate. If `None`,
    713         all modes are interpolated.
    714         Defaults to `None`.
    715 
    716     return_grads : bool
    717         If `True`, gradients and Hessians are computed.
    718         Defaults to `False`.
    719 
    720     Output
    721     ------
    722     values : [num_modes, num_samples] or [num_samples], tf.float
    723         Interpolated profile values at the sample positions
    724 
    725     grads : [num_modes, num_samples, 3] or [num_samples, 3], tf.float
    726         Gradients of the interpolated profile values
    727         at the sample positions
    728 
    729     hessians : [num_modes, num_samples, 3, 3] or [num_samples, 3, 3] , tf.float
    730         Hessians of the interpolated profile values
    731         at the sample positions
    732     """
    733 
    734     @staticmethod
    735     def lagrange_polynomials(x,
    736                              x_i,
    737                              return_derivatives=True):
    738         # pylint: disable=line-too-long
    739         r"""
    740         Compute the 2nd-order Lagrange polynomials
    741 
    742         Optionally, the first- and second-order derivatives are returned.
    743 
    744         The 2nd-order Lagrange polynomials :math:`\ell_j(x)`, :math:`j=1,2,3`,
    745         for position :math:`x\in\mathbb{R}` are computed using three distinct
    746         support positions :math:`x_i` for :math:`i=1,2,3`:
    747 
    748         .. math::
    749             \begin{align}
    750                 \ell_j(x) &= \prod_{\substack{1\leq i \leq 3 \\ i \ne j}} \frac{x-x_i}{x_j-x_i}.
    751             \end{align}
    752 
    753         Their first- and second-order derivatives are then respectively given as
    754 
    755         .. math::
    756             \begin{align}
    757                 \ell'_j(x)  &= \left(\sum_{i \ne j} x-x_i\right) \left(\prod_{i \ne j} x_j-x_i\right)^{-1} \\
    758                 \ell''_j(x) &= 2 \left(\prod_{i \ne j} x_j-x_i\right)^{-1}.
    759             \end{align}
    760 
    761         Input
    762         -----
    763         x : [batch_size], tf.float
    764             Sample positions
    765 
    766         x_i : [batch_size, 3], tf.float
    767             Support positions for every sample position
    768 
    769         return_derivatives : bool
    770             If `True`, also the first- and second-order derivatives
    771             of the Lagrange polynomials are returned.
    772             Defaults to `True`.
    773 
    774         Output
    775         ------
    776         l_i : [batch_size, 3], tf.float
    777             Lagrange polynomials for each sample position
    778 
    779         deriv_1st : [batch_size, 3], tf.float
    780             First-order derivatives for each sample position.
    781             Only returned if `return_derivatives` is `True`.
    782 
    783         deriv_2nd : [batch_size, 3], tf.float
    784             Second-order derivatives for each sample position.
    785             Only returned if `return_derivatives` is `True`.
    786         """
    787 
    788         # Compute products of differences of the sample and support points
    789         sample_diff = tf.expand_dims(x, 1) - x_i
    790         sample_prod_0 = sample_diff[:,1]*sample_diff[:,2]
    791         sample_prod_1 = sample_diff[:,0]*sample_diff[:,2]
    792         sample_prod_2 = sample_diff[:,0]*sample_diff[:,1]
    793         sample_prods = tf.stack([sample_prod_0, sample_prod_1, sample_prod_2],
    794                                 -1)
    795 
    796         # Compute products of differences of support points
    797         support_diffs = tf.expand_dims(x_i, -1) - tf.expand_dims(x_i, -2)
    798         support_diffs = tf.where(support_diffs==0, 1., support_diffs)
    799         support_prods = tf.reduce_prod(support_diffs, axis=-1)
    800 
    801         # Compute Lagrange polynomials
    802         lagrange = sample_prods/support_prods
    803 
    804         if not return_derivatives:
    805             return lagrange
    806         else:
    807             # Compute sums of differences
    808             sample_sum_0 = sample_diff[:,1] + sample_diff[:,2]
    809             sample_sum_1 = sample_diff[:,0] + sample_diff[:,2]
    810             sample_sum_2 = sample_diff[:,0] + sample_diff[:,1]
    811             sample_sums = tf.stack([sample_sum_0, sample_sum_1, sample_sum_2],
    812                                     -1)
    813             # Compute first-order derivatives
    814             deriv_1st = sample_sums/support_prods
    815 
    816             # Compute second-order derivatives
    817             deriv_2nd = tf.cast(2, support_prods.dtype)/support_prods
    818 
    819             return lagrange, deriv_1st, deriv_2nd
    820 
    821     def __call__(self, points, mode=None, return_grads=False):
    822         # pylint: disable=line-too-long
    823         r"""
    824         Interpolates a discrete profile at arbitrary position via
    825         2D 2nd-order Lagrange interpolation.
    826 
    827         A discrete profile :math:`P(y_i,z_j)\in\mathbb{R}` defined on
    828         a grid of points :math:`y_i,z_j` for :math:`i,j \in [1,2,3]` is
    829         interpolated to position :math:`y,z` as
    830 
    831         .. math::
    832             \begin{align}
    833                 P(y,z) &= \sum_{i,j} P(y_i,z_j) \ell_{i,y}(y) \ell_{j,z}(z)
    834             \end{align}
    835 
    836         where :math:`\ell_{i,y}(y)`, :math:`\ell_{j,z}(z)` are the
    837         one-dimensional 2nd-order Lagrange polynomials, defined
    838         as
    839 
    840         .. math::
    841             \begin{align}
    842                 \ell_{i,y}(y) &= \prod_{j \ne i} \frac{y-y_j}{y_i-y_j} \\
    843                 \ell_{j,z}(z) &= \prod_{i \ne j} \frac{z-z_i}{z_j-z_i}.
    844             \end{align}
    845 
    846         In order to compute spatial gradients and Hessians, we extend the the profile
    847         with a dummy :math:`x` dimension, i.e., :math:`P(x,y,z)=P(y,z)`, such that
    848 
    849         .. math::
    850             \begin{align}
    851                 \nabla P(x,y,z) &= \begin{bmatrix} 0, \frac{\partial P(x,y,z)}{\partial y}, \frac{\partial P(x,y,z)}{\partial z}  \end{bmatrix}^{\textsf{T}}\\
    852                 H_P(x,y,z) &= \begin{bmatrix} 0 & 0                                                 & 0 \\
    853                                               0 & \frac{\partial^2 P(x,y,z)}{\partial y^2}          & \frac{\partial^2 P(x,y,z)}{\partial y \partial z} \\
    854                                               0 & \frac{\partial^2 P(x,y,z)}{\partial z \partial y} & \frac{\partial^2 P(x,y,z)}{\partial z^2}
    855                              \end{bmatrix}
    856             \end{align}.
    857 
    858 
    859         Input
    860         -----
    861         points : [num_samples, 2], tf.float
    862             Positions at which to interpolate the profile
    863 
    864         mode : int | `None`
    865             Mode of the profile to interpolate. If `None`,
    866             all modes are interpolated.
    867             Defaults to `None`.
    868 
    869         return_grads : bool
    870             If `True`, gradients and Hessians are computed.
    871             Defaults to `False`.
    872 
    873         Output
    874         ------
    875         values : [num_modes, num_samples] or [num_samples], tf.float
    876             Interpolated profile values at the sample positions
    877 
    878         grads : [num_modes, num_samples, 3] or [num_samples, 3], tf.float
    879             Gradients of the interpolated profile values
    880             at the sample positions
    881 
    882         hessians : [num_modes, num_samples, 3, 3] or [num_samples, 3, 3] , tf.float
    883             Hessians of the interpolated profile values
    884             at the sample positions
    885         """
    886         num_samples = tf.shape(points)[0]
    887         # Compute absolute distances in y/z directions
    888         y_dist = tf.abs(tf.expand_dims(points[:,0], axis=1)
    889                         - tf.expand_dims(self.cell_y_positions, axis=0))
    890         z_dist = tf.abs(tf.expand_dims(points[:,1], axis=1)
    891                         - tf.expand_dims(self.cell_z_positions, axis=0))
    892 
    893         # Compute indices of three closest support points
    894         y_ind = tf.sort(tf.math.top_k(-y_dist, k=3, sorted=False)[1], -1)
    895         z_ind = tf.sort(tf.math.top_k(-z_dist, k=3, sorted=False)[1], -1)
    896 
    897         # Get support points in y and z dimensions
    898         y_i = tf.gather(self.cell_y_positions, y_ind, axis=0, batch_dims=1)
    899         z_i = tf.gather(self.cell_z_positions, z_ind, axis=0, batch_dims=1)
    900 
    901         # Compute indices of all support points
    902         support_ind = tf.reshape(tf.expand_dims(z_ind, 1)
    903                                  + tf.expand_dims(y_ind, 2)*self.num_rows,
    904                                  [num_samples, -1])
    905 
    906 
    907         # Compute support values for all modes
    908         vals = tf.transpose(self.values, perm=[2, 1, 0])
    909         if mode is not None:
    910             # Filter relevant mode
    911             vals = tf.expand_dims(vals[...,mode], -1)
    912         num_modes = tf.shape(vals)[-1]
    913         vals = tf.reshape(vals, [-1, num_modes])
    914         support_values = tf.gather(vals, support_ind, axis=0, batch_dims=1)
    915         support_values = tf.transpose(support_values, perm=[2,0,1])
    916 
    917         if not return_grads:
    918             # Compute Lagrange polynomials
    919             l_y = self.lagrange_polynomials(points[:,0], y_i, False)
    920             l_z = self.lagrange_polynomials(points[:,1], z_i, False)
    921             l_z_y = tf.reshape(tf.expand_dims(l_y, axis=-1)
    922                             * tf.expand_dims(l_z, axis=-2),
    923                             [num_samples, -1])
    924 
    925             # Compute interpolated values
    926             values = tf.reduce_sum(support_values*l_z_y, axis=-1)
    927             return tf.squeeze(values)
    928 
    929         # Compute Lagrange polynomials and derivatives
    930         l_y, d1_y, d2_y = self.lagrange_polynomials(points[:,0], y_i, True)
    931         l_z, d1_z, d2_z = self.lagrange_polynomials(points[:,1], z_i, True)
    932         l_z_y = tf.reshape(tf.expand_dims(l_y, axis=-1)
    933                            * tf.expand_dims(l_z, axis=-2),
    934                            [num_samples, -1])
    935 
    936         # Compute interpolated values
    937         values = tf.reduce_sum(support_values*l_z_y, axis=-1)
    938 
    939         # Compute gradients
    940         l_z_d_y = tf.reshape(tf.expand_dims(d1_y, axis=-1)
    941                             * tf.expand_dims(l_z, axis=-2),
    942                             [num_samples, -1])
    943 
    944         d_values_dy = tf.reduce_sum(support_values*l_z_d_y, axis=-1)
    945 
    946         l_d_z_y = tf.reshape(tf.expand_dims(l_y, axis=-1)
    947                             * tf.expand_dims(d1_z, axis=-2),
    948                             [num_samples, -1])
    949         d_values_dz = tf.reduce_sum(support_values*l_d_z_y, axis=-1)
    950 
    951         grads = tf.stack([tf.zeros_like(d_values_dy),
    952                           d_values_dy,
    953                           d_values_dz ], -1)
    954 
    955         # Compute Hessians
    956         # 1: Compute 2nd-order partial derivatives
    957         l_z_d2_y = tf.reshape(tf.expand_dims(d2_y, axis=-1)
    958                               * tf.expand_dims(l_z, axis=-2),
    959                               [num_samples, -1])
    960         d2_values_d2_y = tf.reduce_sum(support_values*l_z_d2_y, axis=-1)
    961 
    962         l_d2_z_y = tf.reshape(tf.expand_dims(l_y, axis=-1)
    963                               * tf.expand_dims(d2_z, axis=-2),
    964                               [num_samples, -1])
    965         d2_values_d2_z = tf.reduce_sum(support_values*l_d2_z_y, axis=-1)
    966 
    967         l_d_z_d_y = tf.reshape(tf.expand_dims(d1_y, axis=-1)
    968                                * tf.expand_dims(d1_z, axis=-2),
    969                                [num_samples, -1])
    970         d2_values_d_y_d_z = tf.reduce_sum(support_values*l_d_z_d_y, axis=-1)
    971 
    972         # 2: Construct rows of the Hessians
    973         row_2 = tf.stack([tf.zeros_like(d2_values_d2_y),
    974                           d2_values_d2_y,
    975                           d2_values_d_y_d_z], -1)
    976 
    977         row_3 = tf.stack([tf.zeros_like(d2_values_d2_z),
    978                           d2_values_d_y_d_z,
    979                           d2_values_d2_z], -1)
    980 
    981         row_1 = tf.zeros_like(row_2)
    982 
    983         # 3: Combine rows full Hessian matrices
    984         hessians = tf.stack([row_1, row_2, row_3], axis=2)
    985         return (values, grads, hessians)
    986 
    987 class DiscreteAmplitudeProfile(DiscreteProfile, AmplitudeProfile):
    988     # pylint: disable=line-too-long
    989     r"""Class defining a discrete amplitude profile of a RIS
    990 
    991     A discrete amplitude profile :math:`A_m` assigns to
    992     each of its units cells a possibly different amplitude value.
    993     Multiple reradiation modes can be obtained by super-positioning
    994     of profiles. The relative power of reradiation modes can
    995     be controlled via the reradiation coefficients :math:`p_m`.
    996 
    997     See :ref:`ris_primer` for more details.
    998 
    999     A class instance is a callable that returns the profile values,
   1000     gradients and Hessians at given points.
   1001 
   1002     Parameters
   1003     ----------
   1004     cell_grid : :class:`~sionna.rt.CellGrid`
   1005         Defines the physical structure of the RIS
   1006 
   1007     num_modes : int
   1008         Number of reradiation modes.
   1009         Defaults to 1.
   1010 
   1011     values : tf.float or tf.Variable, [num_modes, num_rows, num_cols]
   1012         Amplitude values for each reradiation mode
   1013         and unit cell. `num_rows` and `num_cols` are defined by the
   1014         `cell_grid`.
   1015         Defaults to `None`.
   1016 
   1017     mode_powers : tf.float, [num_modes]
   1018         Relative powers or reradition coefficients of reradiation modes.
   1019         Defaults to `None`. In this case, all reradiation modes get
   1020         an equal fraction of the total power.
   1021 
   1022     interpolator : :class:`~sionna.rt.ProfileInterpolator`
   1023         Determines how the discrete values of the profile
   1024         are interpolated to a continuous profile
   1025         which is defined at any point on the RIS.
   1026         Defaults to `None`. In this case, the
   1027         :class:`~sionna.rt.LagrangeProfileInterpolator` will be used.
   1028 
   1029     dtype : tf.complex
   1030         Datatype to be used in internal calculations.
   1031         Defaults to `tf.complex64`.
   1032 
   1033     Input
   1034     -----
   1035     points : tf.float, [num_samples, 2]
   1036         Tensor of 2D coordinates defining the points on the RIS at which
   1037         the profile should be evaluated.
   1038         Defaults to `None`. In this case, the values for all unit cells
   1039         are returned.
   1040 
   1041     mode : int | `None`
   1042         Reradiation mode to be considered.
   1043         Defaults to `None`. In this case, the values for all modes
   1044         are returned.
   1045 
   1046     return_grads : bool
   1047         If `True`, also the first- and second-order derivatives are
   1048         returned.
   1049         Defaults to `False`.
   1050 
   1051     Output
   1052     ------
   1053     values : [num_modes, num_samples] or [num_samples], tf.float
   1054         Interpolated profile values at the sample positions
   1055 
   1056     grads : [num_modes, num_samples, 3] or [num_samples, 3], tf.float
   1057         Gradients of the interpolated profile values
   1058         at the sample positions. Only returned if `return_grads` is `True`.
   1059 
   1060     hessians : [num_modes, num_samples, 3, 3] or [num_samples, 3, 3] , tf.float
   1061         Hessians of the interpolated profile values
   1062         at the sample positions. Only returned if `return_grads` is `True`.
   1063     """
   1064     def __init__(self,
   1065                  cell_grid,
   1066                  num_modes=1,
   1067                  values=None,
   1068                  mode_powers=None,
   1069                  interpolator=None,
   1070                  dtype=tf.complex64):
   1071         super().__init__(cell_grid=cell_grid,
   1072                          num_modes=num_modes,
   1073                          values=values,
   1074                          interpolator=interpolator,
   1075                          dtype=dtype)
   1076 
   1077         if values is None:
   1078             self.values = tf.ones(self.shape, self._rdtype)
   1079 
   1080         if mode_powers is None:
   1081             mode_powers = 1/tf.cast(self.num_modes, self._rdtype) * \
   1082                           tf.ones([self.num_modes], dtype=self._rdtype)
   1083         self.mode_powers = mode_powers
   1084 
   1085     @property
   1086     def mode_powers(self):
   1087         return self._mode_powers
   1088 
   1089     @mode_powers.setter
   1090     def mode_powers(self, v):
   1091         if isinstance(v, tf.Variable):
   1092             if v.dtype != self._rdtype:
   1093                 msg = f"`mode_powers` must have dtype={self._rdtype}"
   1094                 raise TypeError(msg)
   1095         else:
   1096             v = tf.cast(v, dtype=self._rdtype)
   1097 
   1098         if not v.shape==[self.num_modes]:
   1099             msg = f"`mode_powers` must have shape [{self.num_modes}]"
   1100             raise ValueError(msg)
   1101 
   1102         self._mode_powers = v
   1103 
   1104 class DiscretePhaseProfile(DiscreteProfile, PhaseProfile):
   1105     # pylint: disable=line-too-long
   1106     r"""Class defining a discrete phase profile of a RIS
   1107 
   1108     A discrete phase profile :math:`\chi_m` assigns to
   1109     each of its units cells a possibly different phase value.
   1110     Multiple reradiation modes can be created by super-positioning
   1111     of phase profiles.
   1112 
   1113     See :ref:`ris_primer` in the Primer on Electromagnetics for more details.
   1114 
   1115     A class instance is a callable that returns the profile values,
   1116     gradients and Hessians at given points.
   1117 
   1118     Parameters
   1119     ----------
   1120     cell_grid : :class:`~sionna.rt.CellGrid`
   1121         Defines the physical structure of the RIS
   1122 
   1123     num_modes : int
   1124         Number of reradiation modes.
   1125         Defaults to 1.
   1126 
   1127     values : tf.float or tf.Variable, [num_modes, num_rows, num_cols]
   1128         Phase values [rad] for each reradiation mode
   1129         and unit cell. `num_rows` and `num_cols` are defined by the
   1130         `cell_grid`.
   1131         Defaults to `None`.
   1132 
   1133     interpolator : :class:`~sionna.rt.ProfileInterpolator`
   1134         Determines how the discrete values of the profile
   1135         are interpolated to a continuous profile
   1136         which is defined at any point on the RIS.
   1137         Defaults to `None`. In this case, the
   1138         :class:`~sionna.rt.LagrangeProfileInterpolator` will be used.
   1139 
   1140     dtype : tf.complex
   1141         Datatype to be used in internal calculations.
   1142         Defaults to `tf.complex64`.
   1143 
   1144     Input
   1145     -----
   1146     points : tf.float, [num_samples, 2]
   1147         Tensor of 2D coordinates defining the points on the RIS at which
   1148         the profile should be evaluated.
   1149         Defaults to `None`. In this case, the values for all unit cells
   1150         are returned.
   1151 
   1152     mode : int | `None`
   1153         Reradiation mode to be considered.
   1154         Defaults to `None`. In this case, the values for all modes
   1155         are returned.
   1156 
   1157     return_grads : bool
   1158         If `True`, also the first- and second-order derivatives are
   1159         returned.
   1160         Defaults to `False`.
   1161 
   1162     Output
   1163     ------
   1164     values : [num_modes, num_samples] or [num_samples], tf.float
   1165         Interpolated profile values at the sample positions
   1166 
   1167     grads : [num_modes, num_samples, 3] or [num_samples, 3], tf.float
   1168         Gradients of the interpolated profile values
   1169         at the sample positions. Only returned if `return_grads` is `True`.
   1170 
   1171     hessians : [num_modes, num_samples, 3, 3] or [num_samples, 3, 3] , tf.float
   1172         Hessians of the interpolated profile values
   1173         at the sample positions. Only returned if `return_grads` is `True`.
   1174     """
   1175     def __init__(self,
   1176                  cell_grid,
   1177                  num_modes=1,
   1178                  values=None,
   1179                  interpolator=None,
   1180                  dtype=tf.complex64):
   1181         super().__init__(cell_grid=cell_grid,
   1182                          num_modes=num_modes,
   1183                          values=values,
   1184                          interpolator=interpolator,
   1185                          dtype=dtype)
   1186 
   1187         if values is None:
   1188             self.values = tf.zeros(self.shape, self._rdtype)
   1189 
   1190 class RIS(RadioDevice, SceneObject):
   1191     # pylint: disable=line-too-long
   1192     r"""
   1193     Class defining a reconfigurable intelligent surface (RIS)
   1194 
   1195     A RIS consists of a planar arrangement of unit cells
   1196     with :math:`\lambda/2` spacing.
   1197     It's :class:`~sionna.rt.PhaseProfile` :math:`\chi_m` and
   1198     :class:`~sionna.rt.AmplitudeProfile` :math:`A_m` can be
   1199     configured after the RIS is instantiated. Both together
   1200     define the spatial modulation coefficient :math:`\Gamma` which
   1201     determines how the RIS reflects electro-magnetic waves.
   1202 
   1203     See :ref:`ris_primer` in the Primer on Electromagnetics for
   1204     more details or have a look at the `tutorial notebook <https://nvlabs.github.io/sionna/examples/Sionna_Ray_Tracing_RIS.html>`_.
   1205 
   1206     An RIS instance is a callable that computes the spatial modulation coefficient
   1207     and gradients/Hessians of the underlying phase profile for provided
   1208     points on the RIS' surface.
   1209 
   1210     Parameters
   1211     ----------
   1212     name : str
   1213         Name
   1214 
   1215     position : [3], float
   1216         Position :math:`(x,y,z)` as three-dimensional vector
   1217 
   1218     num_rows : int
   1219         Number of rows. Must at least be equal to three.
   1220 
   1221     num_cols : int
   1222         Number of columns. Must at least be equal to three.
   1223 
   1224     num_modes : int
   1225         Number of reradiation modes.
   1226         Defaults to 1.
   1227 
   1228     orientation : [3], float
   1229         Orientation :math:`(\alpha, \beta, \gamma)` specified
   1230         through three angles corresponding to a 3D rotation
   1231         as defined in :eq:`rotation`.
   1232         This parameter is ignored if ``look_at`` is not `None`.
   1233         Defaults to [0,0,0]. In this case, the normal vector of
   1234         the RIS points towards the positive x-axis.
   1235 
   1236     velocity : [3], float
   1237         Velocity vector [m/s]. Used for the computation of
   1238         path-specific Doppler shifts.
   1239 
   1240     look_at : [3], float | :class:`~sionna.rt.Transmitter` | :class:`~sionna.rt.Receiver` | :class:`~sionna.rt.RIS` | :class:`~sionna.rt.Camera` | `None`
   1241         A position or the instance of a :class:`~sionna.rt.Transmitter`,
   1242         :class:`~sionna.rt.Receiver`, :class:`~sionna.rt.RIS`, or :class:`~sionna.rt.Camera` to look at.
   1243         If set to `None`, then ``orientation`` is used to orientate the device.
   1244 
   1245     color : [3], float
   1246         Defines the RGB (red, green, blue) ``color`` parameter for the device as displayed in the previewer and renderer.
   1247         Each RGB component must have a value within the range :math:`\in [0,1]`.
   1248         Defaults to `[0.862,0.078,0.235]`.
   1249 
   1250     dtype : tf.complex
   1251         Datatype to be used in internal calculations.
   1252         Defaults to `tf.complex64`.
   1253 
   1254     Input
   1255     -----
   1256     points : tf.float, [num_samples, 2]
   1257         Tensor of 2D coordinates defining the points on the RIS at which
   1258         the spatial modulation profile should be evaluated.
   1259         Defaults to `None`. In this case, the values for all unit cells
   1260         are returned.
   1261 
   1262     mode : int | `None`
   1263         Reradiation mode to be considered.
   1264         Defaults to `None`. In this case, the values for all modes
   1265         are returned.
   1266 
   1267     return_grads : bool
   1268         If `True`, also the first- and second-order derivatives are
   1269         returned.
   1270         Defaults to `False`.
   1271 
   1272     Output
   1273     ------
   1274     gamma : [num_modes, num_samples] or [num_samples], tf.complex
   1275         Spatial modulation coefficient at the sample positions
   1276 
   1277     grads : [num_modes, num_samples, 3] or [num_samples, 3], tf.float
   1278         Gradients of the interpolated phase profile values
   1279         at the sample positions. Only returned if `return_grads` is `True`.
   1280 
   1281     hessians : [num_modes, num_samples, 3, 3] or [num_samples, 3, 3] , tf.float
   1282         Hessians of the interpolated phase profile values
   1283         at the sample positions. Only returned if `return_grads` is `True`.
   1284     """
   1285     def __init__(self,
   1286                 name,
   1287                 position,
   1288                 num_rows,
   1289                 num_cols,
   1290                 num_modes=1,
   1291                 orientation=(0.,0.,0.),
   1292                 velocity=(0.,0.,0.),
   1293                 look_at=None,
   1294                 color=(0.862,0.078,0.235),
   1295                 dtype=tf.complex64):
   1296 
   1297         # Initialize the parent classes
   1298         # RadioDevice and SceneObject inherit from Object
   1299         # Python will initialize in the following order:
   1300         # RadioDevice->SceneObject->Object
   1301         super().__init__(name=name,
   1302                          position=position,
   1303                          orientation=orientation,
   1304                          look_at=look_at,
   1305                          radio_material=None,
   1306                          color=color,
   1307                          dtype=dtype)
   1308 
   1309         # Set velocity vector
   1310         self.velocity = tf.cast(velocity, dtype=dtype.real_dtype)
   1311 
   1312         if num_rows < 3 or num_cols < 3:
   1313             raise ValueError("num_rows and num_cols must be >= 3")
   1314 
   1315         # Set immutable properties
   1316         self._num_modes = int(num_modes)
   1317         self._cell_grid = CellGrid(num_rows, num_cols, self._dtype)
   1318 
   1319         # Init amplitude profile
   1320         self.amplitude_profile = DiscreteAmplitudeProfile(self.cell_grid,
   1321                                                      num_modes=self.num_modes,
   1322                                                      dtype=self._dtype)
   1323 
   1324         # Init phase profile
   1325         self.phase_profile = DiscretePhaseProfile(self.cell_grid,
   1326                                              num_modes=self.num_modes,
   1327                                              dtype=self._dtype)
   1328 
   1329     @property
   1330     def cell_grid(self):
   1331         r"""
   1332         :class:`~sionna.rt.CellGrid` : Defines the physical
   1333             structure of the RIS
   1334         """
   1335         return self._cell_grid
   1336 
   1337     @property
   1338     def cell_positions(self):
   1339         r"""
   1340         [num_cells, 2], tf.float : Cell positions in the
   1341             local coordinate system (LCS) of the RIS, ordered
   1342             from top-to-bottom left-to-right.
   1343         """
   1344         return self.cell_grid.cell_positions*self.spacing
   1345 
   1346     @property
   1347     def cell_world_positions(self):
   1348         r"""
   1349         [num_cells, 3], tf.float : Cell positions in the
   1350             global coordinate system (GCS) of the RIS, ordered
   1351             from top-to-bottom left-to-right.
   1352         """
   1353         x_coord = tf.zeros([self.num_cells, 1], self._rdtype)
   1354         pos = tf.concat([x_coord, self.cell_positions], axis=-1)
   1355         pos = rotate(pos, self.orientation)
   1356         pos += tf.expand_dims(self.position, 0)
   1357         return pos
   1358 
   1359     @property
   1360     def world_normal(self):
   1361         r"""
   1362         [3], tf.float : Normal vector of the RIS in the
   1363             global coordinate system (GCS)
   1364         """
   1365         n_hat = tf.constant([1,0,0], self._rdtype)
   1366         return rotate(n_hat, self.orientation)
   1367 
   1368     @property
   1369     def num_rows(self):
   1370         r"""
   1371         int : Number of rows
   1372         """
   1373         return self.cell_grid.num_rows
   1374 
   1375     @property
   1376     def num_cols(self):
   1377         r"""
   1378         int : Number of columns
   1379         """
   1380         return self.cell_grid.num_cols
   1381 
   1382     @property
   1383     def num_cells(self):
   1384         r"""
   1385         int : Number of cells
   1386         """
   1387         return self.num_rows*self.num_cols
   1388 
   1389     @property
   1390     def num_modes(self):
   1391         r"""
   1392         int : Number of reradiation modes
   1393         """
   1394         return self._num_modes
   1395 
   1396     @property
   1397     def spacing(self):
   1398         r"""
   1399         tf.float: Element spacing [m] corresponding to
   1400             half a wavelength
   1401         """
   1402         if hasattr(scene.Scene(), "wavelength"):
   1403             wavelength = tf.cast(scene.Scene().wavelength, self._rdtype)
   1404             return wavelength/tf.cast(2, self._rdtype)
   1405         else:
   1406             # Scene is not initialized
   1407             return tf.cast(0.5, self._rdtype)
   1408 
   1409     @property
   1410     def size(self):
   1411         """
   1412         [2], tf.float : Size of the RIS (width, height) [m]
   1413         """
   1414         return tf.stack([self.spacing * self.num_cols,
   1415                           self.spacing * self.num_rows], axis=0)
   1416 
   1417     @property
   1418     def velocity(self):
   1419         """
   1420         [3], tf.float : Get/set the velocity vector [m/s]
   1421         """
   1422         return self._velocity
   1423 
   1424     @velocity.setter
   1425     def velocity(self, v):
   1426         if not tf.shape(v)==3:
   1427             raise ValueError("`velocity` must have shape [3]")
   1428         self._velocity = tf.cast(v, self._dtype.real_dtype)
   1429 
   1430     @property
   1431     def amplitude_profile(self):
   1432         r"""
   1433         :class:`~sionna.rt.AmplitudeProfile` : Set/get amplitude profile
   1434         """
   1435         return self._amplitude_profile
   1436 
   1437     @amplitude_profile.setter
   1438     def amplitude_profile(self, v):
   1439         if not isinstance(v, AmplitudeProfile):
   1440             raise ValueError("Not a valid AmplitudeProfile")
   1441         self._amplitude_profile = v
   1442 
   1443     @property
   1444     def phase_profile(self):
   1445         r"""
   1446         :class:`~sionna.rt.PhaseProfile` : Set/get phase profile
   1447         """
   1448         return self._phase_profile
   1449 
   1450     @phase_profile.setter
   1451     def phase_profile(self, v):
   1452         if not isinstance(v, PhaseProfile):
   1453             raise ValueError("Not a valid PhaseProfile")
   1454         self._phase_profile = v
   1455 
   1456     def phase_gradient_reflector(self, sources, targets):
   1457         # pylint: disable=line-too-long
   1458         r"""
   1459         Configures the RIS as ideal phase gradient reflector
   1460 
   1461         For an incoming direction :math:`\hat{\mathbf{k}}_i`
   1462         and desired outgoing direction :math:`\hat{\mathbf{k}}_r`,
   1463         the necessary phase gradient along the RIS with normal
   1464         :math:`\hat{\mathbf{n}}` can be computed as
   1465         (e.g., Eq.(12) [Vitucci24]_):
   1466 
   1467         .. math::
   1468             \nabla\chi_m = k_0\left( \mathbf{I}- \hat{\mathbf{n}}\hat{\mathbf{n}}^\textsf{T} \right) \left(\hat{\mathbf{k}}_i - \hat{\mathbf{k}}_r  \right).
   1469 
   1470         The phase profile is obtained by assigning zero phase to the first
   1471         unit cell and evolving the other phases linearly according to the gradient
   1472         across the entire RIS.
   1473 
   1474         Multiple reradiation modes can be configured.
   1475 
   1476         The amplitude profile is set to one everywhere with a uniform relative
   1477         power allocation across modes.
   1478 
   1479         Input
   1480         -----
   1481         sources : tf.float, [3] or [num_modes, 3]
   1482             Tensor defining for every reradiation mode
   1483             a source from which the incoming wave originates.
   1484 
   1485         targets : tf.float, [3] or [num_modes, 3]
   1486             Tensor defining for every reradiation mode
   1487             a target towards which the incoming wave should be
   1488             reflected.
   1489         """
   1490         # Convert inputs to tensors
   1491         sources = tf.cast(sources, self._rdtype)
   1492         targets = tf.cast(targets, self._rdtype)
   1493         sources = expand_to_rank(sources, 2, 0)
   1494         targets = expand_to_rank(targets, 2, 0)
   1495         shape = [self.num_modes, 3]
   1496 
   1497         # Ensure the desired shape [num_modes, 3]
   1498         for i, x in enumerate([sources, targets]):
   1499             if not (tf.shape(x)==shape).numpy().all():
   1500                 msg = f"Wrong shape of input {i+1}. " + \
   1501                       f"Expected {shape}, got {x.shape}"
   1502                 raise ValueError(msg)
   1503 
   1504         # Compute incoming and outgoing directions
   1505         # [num_modes, 3]
   1506         k_i, _ = normalize(self.position[tf.newaxis] - sources)
   1507         k_r, _ = normalize(targets - self.position[tf.newaxis])
   1508 
   1509         # Tangent projection operator - Eq.(10)
   1510         # [1, 3]
   1511         normal = self.world_normal[tf.newaxis]
   1512         # [1, 3, 3]
   1513         p = tf.eye(3, dtype=self._rdtype) - outer(normal,normal)
   1514 
   1515         # Compute phase gradient - Eq.(12)
   1516         # [num_modes, 3]
   1517         grad = self.scene.wavenumber * tf.linalg.matvec(p, k_i-k_r)
   1518         # Rotate phase gradient to LCS of the RIS and keep y/z components
   1519         # [num_modes, 1, 1, 2]
   1520         grad = rotate(grad, self.orientation, inverse=True)[:,1:]
   1521         grad = tf.reshape(grad, [self.num_modes, 1, 1, 2])
   1522 
   1523         # Using the top-left cell as reference, compute the offsets
   1524         # [1, num_rows, num_cols, 2]
   1525         offsets = self.cell_positions - self.cell_positions[:1]
   1526         offsets = tf.reshape(offsets, [self.num_cols, self.num_rows, 2])
   1527         offsets = tf.transpose(offsets, perm=[1,0,2])
   1528         offsets = tf.expand_dims(offsets, 0)
   1529 
   1530         # Compute phase profile based on the constant gradient assumption
   1531         # [num_modes, num_rows, num_cols]
   1532         phases = tf.reduce_sum(offsets*grad, axis=-1)
   1533         self.phase_profile.values = phases
   1534 
   1535         # Set a neutral amplitude profile
   1536         self.amplitude_profile.values = tf.ones_like(phases)
   1537         mode_powers = 1/tf.cast(self.num_modes, self._rdtype) * \
   1538                           tf.ones([self.num_modes], dtype=self._rdtype)
   1539         self.amplitude_profile.mode_powers = mode_powers
   1540 
   1541     def focusing_lens(self, sources, targets):
   1542         # pylint: disable=line-too-long
   1543         r"""
   1544         Configures the RIS as focusing lens
   1545 
   1546         The phase profile is configured in such a way that
   1547         the fields of all rays add up coherently at a specific
   1548         point. In other words, the phase profile undoes the
   1549         distance-based phase shift of every ray connecting a
   1550         source to a target via a specific unit cell.
   1551 
   1552         For a source and target at positions
   1553         :math:`\mathbf{s}` and :math:`\mathbf{t}`, the phase
   1554         :math:`\chi_m(\mathbf{x})` of a unit cell located at :math:`\mathbf{x}`
   1555         is computed as (e.g., Sec. IV-2 [Degli-Esposti22]_)
   1556 
   1557         .. math::
   1558             \chi_m(\mathbf{x}) = k_0 \left(\lVert\mathbf{s}-\mathbf{x}\rVert + \lVert\mathbf{s}-\mathbf{t}\rVert\right).
   1559 
   1560         Multiple reradiation modes can be configured.
   1561 
   1562         The amplitude profile is set to one everywhere with a uniform relative
   1563         power allocation across modes.
   1564 
   1565         Input
   1566         -----
   1567         sources : tf.float, [3] or [num_modes, 3]
   1568             Tensor defining for every reradiation mode
   1569             a source from which the incoming wave originates.
   1570 
   1571         targets : tf.float, [3] or [num_modes, 3]
   1572             Tensor defining for every reradiation mode
   1573             a target towards which the incoming wave should be
   1574             reflected.
   1575         """
   1576         # Convert inputs to tensors
   1577         sources = tf.cast(sources, self._rdtype)
   1578         targets = tf.cast(targets, self._rdtype)
   1579         sources = expand_to_rank(sources, 2, 0)
   1580         targets = expand_to_rank(targets, 2, 0)
   1581         shape = [self.num_modes, 3]
   1582 
   1583         # Ensure the desired shape [num_modes, 3]
   1584         for i, x in enumerate([sources, targets]):
   1585             if not (tf.shape(x)==shape).numpy().all():
   1586                 msg = f"Wrong shape of input {i+1}. " + \
   1587                       f"Expected {shape}, got {x.shape}"
   1588                 raise ValueError(msg)
   1589 
   1590         # Compute incoming and outgoing distances
   1591         # [num_modes, num_cells]
   1592         d_i = normalize(self.cell_world_positions[tf.newaxis] - sources[:,tf.newaxis])[1]
   1593         d_o = normalize(self.cell_world_positions[tf.newaxis] - targets[:,tf.newaxis])[1]
   1594 
   1595         # Compute phases such that the total phase shifts for all cells
   1596         # are equal
   1597         phases = self.scene.wavenumber * (d_i+d_o)
   1598         phases = tf.reshape(phases, [self.num_modes, self.num_cols, self.num_rows])
   1599         phases = tf.transpose(phases, perm=[0,2,1])
   1600         self.phase_profile.values = phases
   1601 
   1602         # Set a neutral amplitude profile
   1603         self.amplitude_profile.values = tf.ones_like(phases)
   1604         mode_powers = 1/tf.cast(self.num_modes, self._rdtype) * \
   1605                           tf.ones([self.num_modes], dtype=self._rdtype)
   1606         self.amplitude_profile.mode_powers = mode_powers
   1607 
   1608     def __call__(self, points=None, mode=None, return_grads=False):
   1609         # pylint: disable=line-too-long
   1610         r"""
   1611         Computes the spatial modulation coefficient and gradients/Hessians of phase profile
   1612 
   1613         Input
   1614         -----
   1615         points : tf.float, [num_samples, 2]
   1616             Tensor of 2D coordinates defining the points on the RIS at which
   1617             the spatial modulation profile should be evaluated.
   1618             Defaults to `None`. In this case, the values for all unit cells
   1619             are returned.
   1620 
   1621         mode : int | `None`
   1622             Reradiation mode to be considered.
   1623             Defaults to `None`. In this case, the values for all modes
   1624             are returned.
   1625 
   1626         return_grads : bool
   1627             If `True`, also the first- and second-order derivatives are
   1628             returned.
   1629             Defaults to `False`.
   1630 
   1631         Output
   1632         ------
   1633         gamma : [num_modes, num_samples] or [num_samples], tf.complex
   1634             Spatial modulation coefficient at the sample positions
   1635 
   1636         grads : [num_modes, num_samples, 3] or [num_samples, 3], tf.float
   1637             Gradients of the interpolated phase profile values
   1638             at the sample positions. Only returned if `return_grads` is `True`.
   1639 
   1640         hessians : [num_modes, num_samples, 3, 3] or [num_samples, 3, 3] , tf.float
   1641             Hessians of the interpolated phase profile values
   1642             at the sample positions. Only returned if `return_grads` is `True`.
   1643         """
   1644         # Get amplitudes
   1645         a = self.amplitude_profile(points, mode)
   1646 
   1647         # Get mode powers
   1648         p = self.amplitude_profile.mode_powers
   1649 
   1650         # Get phases and (optionally) phase gradients and Hessians
   1651         if return_grads and points is not None:
   1652             chi, grads, hessians = self.phase_profile(points, mode, True)
   1653         else:
   1654             chi = self.phase_profile(points, mode, False)
   1655 
   1656         # Compute spatial modulation coefficient
   1657         zero = tf.cast(0, self._rdtype)
   1658         gamma = tf.complex(a, zero)
   1659         chi = tf.complex(zero, chi)
   1660         p = tf.complex(tf.sqrt(p), zero)
   1661         gamma *= tf.exp(chi)
   1662         if mode is None:
   1663             gamma*= tf.reshape(p, [-1, 1])
   1664         else:
   1665             gamma *= p[mode]
   1666 
   1667         if return_grads and points is not None:
   1668             return gamma, grads, hessians
   1669         else:
   1670             return gamma