differt.em.reflection_coefficients

reflection_coefficients(n_r, cos_theta_i)

Compute the Fresnel reflection coefficients at an interface.

Parameters:

• n_r (Inexact[ArrayLike, '*#batch']) –

  The relative refractive indices.

  This is the ratios of the refractive indices of the second media over the refractive indices of the first media.

• cos_theta_i (Float[ArrayLike, '*#batch']) – The (cosine of the) angles of incidence (or reflection).

Return type:

tuple[Complex[Array, '*batch'], Complex[Array, '*batch']]

Returns:

The reflection coefficients for s and p polarizations.

See also

fresnel_coefficients

refraction_coefficients

refractive_indices

Examples

The following example show how to compute interference patterns from line of sight and reflection on a glass ground.

>>> from differt.em import (
...     Dipole,
...     c,
...     z_0,
...     fspl,
...     poynting_vector,
...     reflection_coefficients,
...     sp_directions,
... )
>>> from differt.geometry import normalize
>>> from differt.rt import image_method

The first step is to define the antenna and the geometry of the scene. Here, we place a dipole antenna above the origin, and generate a num_positions number of positions along the horizontal line, where we will evaluate the EM fields.

>>> tx_position = jnp.array([0.0, 2.0, 0.0])
>>> rx_position = jnp.array([0.0, 2.0, 0.0])
>>> num_positions = 1000
>>> # [num_positions 3]
>>> x = jnp.logspace(0, 3, num_positions)  # From close to very far
>>> rx_positions = (
...     jnp.tile(rx_position, (num_positions, 1)).at[..., 0].add(x)
... )
>>> ant = Dipole(2.4e9)  # 2.4 GHz
>>> A_e = ant.aperture  # Effective aperture
>>> plt.xscale("symlog", linthresh=1e-1)
>>> plt.plot(
...     [tx_position[0]],
...     [tx_position[1]],
...     "o",
...     label="TX",
... )
>>> plt.plot(
...     rx_positions[::50, 0],
...     rx_positions[::50, 1],
...     "o",
...     label="RXs",
... )
>>> plt.axhline(color="k", label="Ground")
>>> plt.xlabel("x-axis (m)")
>>> plt.ylabel("y-axis (m)")
>>> plt.legend()
>>> plt.tight_layout()

Next, we compute the EM fields from the direct (line-of-sight) path. We also plot the free-space path loss (see fspl [22]) as a reference.

>>> # [num_positions 3]
>>> E_los, B_los = ant.fields(rx_positions - tx_position)
>>> # [num_positions]
>>> P_los = A_e * jnp.linalg.norm(poynting_vector(E_los, B_los), axis=-1)
>>> plt.semilogx(
...     x,
...     10 * jnp.log10(P_los / ant.reference_power),
...     label=r"$P_\text{los}$",
... )
>>> _, d = normalize(rx_positions - tx_position, keepdims=True)
>>> plt.semilogx(
...     x,
...     -fspl(d, ant.frequency, dB=True),
...     "k-.",
...     label="FSPL",
... )

After, the image_method function is used to compute the reflection points.

>>> ground_vertex = jnp.array([0.0, 0.0, 0.0])
>>> ground_normal = jnp.array([0.0, 1.0, 0.0])
>>> # [num_positions 3]
>>> reflection_points = image_method(
...     tx_position,
...     rx_positions,
...     ground_vertex[None, ...],
...     ground_normal[None, ...],
... ).squeeze(axis=-2)  # Squeeze because only one reflection
>>> # [num_positions 3], [num_positions 1]
>>> k_i, s_i = normalize(reflection_points - tx_position, keepdims=True)
>>> k_r, s_r = normalize(rx_positions - reflection_points, keepdims=True)
>>> # [num_positions 1]
>>> l = jnp.linalg.norm(rx_positions - tx_position, axis=-1, keepdims=True)
>>> tau = (s_i + s_r - l) / c  # Delay between two paths
>>> tau = tau.squeeze(axis=-1)

We then compute the EM fields at those points, and use the Fresnel reflection coefficients to compute the reflected fields.

>>> # [num_positions 3]
>>> E_i, B_i = ant.fields(reflection_points - tx_position, t=-tau)
>>> # [num_positions 1]
>>> cos_theta = jnp.sum(ground_normal * -k_i, axis=-1, keepdims=True)
>>> n_r = 1.5  # Air to glass
>>> # [num_positions 1]
>>> r_s, r_p = reflection_coefficients(n_r, cos_theta)

To apply the coefficients correctly, we must determine the polarization directions of both the incident and the reflected fields.

Important

Reflection coefficients are returned based on s and p directions. As a result, we need to first determine those local directions, and apply the corresponding reflection coefficients to the projection of the fields onto those directions [11, eq. 3.3-3.8 and 3.39, p. 70 and 77].

>>> # [num_positions 3]
>>> (e_i_s, e_i_p), (e_r_s, e_r_p) = sp_directions(k_i, k_r, ground_normal)

We then transform XYZ-components into local s and p components.

>>> # [num_positions 1]
>>> E_i_s = jnp.sum(E_i * e_i_s, axis=-1, keepdims=True)
>>> E_i_p = jnp.sum(E_i * e_i_p, axis=-1, keepdims=True)
>>> B_i_s = jnp.sum(B_i * e_i_s, axis=-1, keepdims=True)
>>> B_i_p = jnp.sum(B_i * e_i_p, axis=-1, keepdims=True)

Then, we apply reflection coefficients to the local s and p components.

>>> # [num_positions 1]
>>> E_r_s = r_s * E_i_s
>>> E_r_p = r_p * E_i_p
>>> B_r_s = r_s * B_i_s
>>> B_r_p = r_p * B_i_p

And we project back to XYZ-components.

>>> E_r = E_r_s * e_r_s + E_r_p * e_r_p
>>> B_r = B_r_s * e_r_s + B_r_p * e_r_p

Finally, we apply the spreading factor and phase shift due to the propagation from the reflection points to the receiver [11, eq. 3.1, p. 63].

>>> spreading_factor = s_i / (
...     s_i + s_r
... )  # We assume that the radii of curvature are equal to 's_i'
>>> phase_shift = jnp.exp(1j * s_r * ant.wavenumber)
>>> E_r *= spreading_factor * phase_shift
>>> B_r *= spreading_factor * phase_shift
>>> P_r = A_e * jnp.linalg.norm(poynting_vector(E_r, B_r), axis=-1)
>>> plt.semilogx(
...     x,
...     10 * jnp.log10(P_r / ant.reference_power),
...     "--",
...     label=r"$P_\text{reflection}$",
... )

We also plot the total field, to better observe the interference pattern.

>>> E_tot = E_los + E_r
>>> B_tot = B_los + B_r
>>> P_tot = A_e * jnp.linalg.norm(poynting_vector(E_tot, B_tot), axis=-1)
>>> plt.semilogx(
...     x,
...     10 * jnp.log10(P_tot / ant.reference_power),
...     "-.",
...     label=r"$P_\text{total}$",
... )
>>> plt.xlabel("Distance to transmitter on x-axis (m)")
>>> plt.ylabel("Gain (dB)")
>>> plt.legend()
>>> plt.tight_layout()

From the above figure, it is clear that the ground-reflection creates an interference pattern in the received power. Moreover, we can clearly observe the Brewster angle at a distance of 6 m. This can verified by computing the Brewster angle from the relative refractive index, and matching it to the corresponding distance.

>>> brewster_angle = jnp.arctan(n_r)
>>> print(f"Brewster angle: {jnp.rad2deg(brewster_angle):.1f}°")
Brewster angle: 56.3°
>>> cos_distance = jnp.abs(jnp.cos(brewster_angle) - cos_theta)
>>> distance = x[jnp.argmin(cos_distance)]
>>> print(f"Corresponding distance: {distance:.1f} m")
Corresponding distance: 6.0 m

Two common approximations are used when estimating the path loss, see Coherent vs. Non-Coherent Radio Wave Propagation. The first one is to assume that the wavefronts are planar, and the second one is to assume that small-scale fading is negligible, or incorrectly modeled, and that the received power can be computed by non-coherent summation of each path.

>>> plt.semilogx(
...     x,
...     10 * jnp.log10(P_tot / ant.reference_power),
...     "-.",
...     color="green",
...     label="Exact",
... )
>>> P_tot_pw = A_e * jnp.sum(E_tot * E_tot, axis=-1) / z_0
>>> plt.semilogx(
...     x,
...     10 * jnp.log10(P_tot_pw / ant.reference_power),
...     ":",
...     color="red",
...     label="Plane-wave",
... )
>>> P_tot_nc = A_e * (
...     jnp.linalg.norm(  # Power from LOS path
...         poynting_vector(E_los, B_los),
...         axis=-1,
...     )
...     + jnp.linalg.norm(  # Power from reflection path
...         poynting_vector(E_r, B_r),
...         axis=-1,
...     )
... )
>>> plt.semilogx(
...     x,
...     10 * jnp.log10(P_tot_nc / ant.reference_power),
...     ":",
...     color="purple",
...     label="Non-coherent",
... )
>>> plt.xlabel("Distance to transmitter on x-axis (m)")
>>> plt.ylabel("Gain (dB)")
>>> plt.legend()
>>> plt.tight_layout()

The plane-wave approximation is relatively good for large distances, but it can be off by several dBs for short distances, as shown below.

>>> plt.xlim([1, 10])
>>> plt.ylim([-65, -35])