diff --git a/.pylintrc b/.pylintrc index 2482c1ed8d..5533de9d3c 100644 --- a/.pylintrc +++ b/.pylintrc @@ -1,7 +1,7 @@ [MASTER] extension-pkg-allow-list=pydantic ignore=material_library.py, plugins -good-names=ax, im, Lx, Ly, Lz, x0, y0, z0, x, y, z, f, t, y1, y2, x1, x2, xs, ys, zs, Ax, Nx, Ny, Nz, dl, rr, E, H, xx, yy, zz, dx, dy, Jx, Jy, Jz, Ex, Ey, Ez, Mx, My, Mz, Hx, Hy, Hz, dz, e, fp, dt, a, c, +good-names=ax, im, Lx, Ly, Lz, x0, y0, z0, x, y, z, f, t, y1, y2, x1, x2, xs, ys, zs, Ax, Nx, Ny, Nz, dl, rr, E, H, xx, yy, zz, dx, dy, Jx, Jy, Jz, Ex, Ey, Ez, Mx, My, Mz, Hx, Hy, Hz, dz, e, fp, dt, a, c, k0, Jx_k, Jy_k, My_k, Mx_k, N_theta, N_phi, L_theta, L_phi [BASIC] diff --git a/tidy3d/convert.py b/tidy3d/convert.py index 708af10e5e..1f91c4ce0e 100644 --- a/tidy3d/convert.py +++ b/tidy3d/convert.py @@ -359,7 +359,7 @@ def old_json_monitors(sim: Simulation) -> Dict: { "type": "FrequencyMonitor", "store": store, - "interpolate": True, + "interpolate": False, } ) elif isinstance(monitor, FieldTimeMonitor): diff --git a/tidy3d/plugins/near2far/near2far.py b/tidy3d/plugins/near2far/near2far.py index 96339c5aa2..99ad82f8e5 100644 --- a/tidy3d/plugins/near2far/near2far.py +++ b/tidy3d/plugins/near2far/near2far.py @@ -3,63 +3,72 @@ import numpy as np from ...constants import C_0, ETA_0 -from ...components.data import FieldData +from ...components.data import SimulationData +from ...log import SetupError + +# TODO: implement new version with simulation.discretize class Near2Far: """Near field to far field transformation tool.""" - def __init__(self, field_data: FieldData): + def __init__(self, sim_data: SimulationData, mon_name: str, frequency: float): """Constructs near field to far field transformation object from monitor data. Parameters ---------- - field_data : FieldData - Description + sim_data : :class:`.SimulationData` + Container for simulation data containing a near field monitor. + mon_name : str + Name of the :class:`.FieldMonitor` to use as source of near field. + Must be a :class:`.FieldMonitor` and stored in ``sim_data``. + frequency : float + Frequency to select from the :class:`.FieldMonitor` to use for projection. + Must be a frequency stored in the :class:`FieldMonitor`. """ - # get frequency info - self.f = float(field_data.data.Ez.f[0].values) - self.k0 = 2 * np.pi * self.f / C_0 - - # get normal axis (ignore components) - xs = field_data.data.Ez.x.values - ys = field_data.data.Ez.y.values - zs = field_data.data.Ez.z.values - mon_size = [xs.shape[-1], ys.shape[-1], zs.shape[-1]] - self.axis = mon_size.index(1) - assert self.axis == 2, "Currently only works for z normal." - - # get normal and planar coordinates - # zs, (xs, ys) = field_data.geometry.pop_axis((xs, ys, zs), axis=self.axis) - self.z0 = zs[0] - self.xs = xs - self.ys = ys - self.xx, self.yy = np.meshgrid(self.xs, self.ys, indexing="ij") - self.dx = np.mean(np.diff(xs)) - self.dy = np.mean(np.diff(ys)) - - # get tangential near fields - # for em_field in "EH": - # for component in "xyz": - # field_name = em_field + component - # assert field_name in field_data.field, f"missing field: {field_name}" - - self.Ex = field_data.data.Ex.values - self.Ey = field_data.data.Ey.values - self.Ez = field_data.data.Ez.values - self.Hx = field_data.data.Hx.values - self.Hy = field_data.data.Hy.values - self.Hz = field_data.data.Hz.values - - # _, (self.Ex, self.Ey) = field_data.geometry.pop_axis(E, axis=self.axis) - # _, (self.Hx, self.Hy) = field_data.geometry.pop_axis(H, axis=self.axis) + try: + field_data = sim_data[mon_name] + except Exception as e: + raise SetupError( + f"No data for monitor named '{mon_name}' " "found in supplied sim_data." + ) from e + + if any(field_name not in field_data for field_name in ("Ex", "Ey", "Hx", "Hy", "Ez")): + raise SetupError(f"Monitor named '{mon_name}' doesn't store all field values") + + try: + Ex = field_data["Ex"].sel(f=frequency) + Ey = field_data["Ey"].sel(f=frequency) + # self.Ez = field_data['Ez'].sel(f=frequency) + Hx = field_data["Hx"].sel(f=frequency) + Hy = field_data["Hy"].sel(f=frequency) + Hz = field_data["Hz"].sel(f=frequency) + except Exception as e: + raise SetupError( + f"Frequency {frequency} not found in all fields " f"from monitor '{mon_name}'." + ) from e + + self.k0 = 2 * np.pi * frequency / C_0 + + # grid sizes + self.dx = np.mean(np.diff(Hz.x.values)) + self.dy = np.mean(np.diff(Hz.y.values)) + + # interpolate to centers + x_centers = Hz.x.values + y_centers = Hz.y.values + Ex = Ex.interp(x=x_centers, y=y_centers) + Ey = Ey.interp(x=x_centers, y=y_centers) + Hx = Hx.interp(x=x_centers, y=y_centers) + Hy = Hy.interp(x=x_centers, y=y_centers) + self.xx, self.yy = np.meshgrid(x_centers, y_centers, indexing="ij") # compute equivalent sources - self.Jx = -self.Hy - self.Jy = self.Hx - self.Mx = self.Ey - self.My = -self.Ex + self.Jx = -np.squeeze(Hy.values) + self.Jy = np.squeeze(Hx.values) + self.Mx = np.squeeze(Ey.values) + self.My = -np.squeeze(Ex.values) def _radiation_vectors(self, theta, phi): """Compute radiation vectors at an angle in spherical coordinates @@ -78,49 +87,32 @@ def _radiation_vectors(self, theta, phi): """ # precompute trig functions and add extra dimensions - theta = np.array(theta) - phi = np.array(phi) - sin_theta = np.sin(theta).reshape((-1, 1, 1)) - cos_theta = np.cos(theta).reshape((-1, 1, 1)) - sin_phi = np.sin(phi).reshape((-1, 1, 1)) - cos_phi = np.cos(phi).reshape((-1, 1, 1)) + sin_theta = np.sin(theta) + cos_theta = np.cos(theta) + sin_phi = np.sin(phi) + cos_phi = np.cos(phi) # precompute fourier transform phase term {dx dy e^(ikrcos(psi))} - FT_phase_x = np.exp(1j * self.k0 * self.xx * sin_theta * cos_phi) - FT_phase_y = np.exp(1j * self.k0 * self.yy * sin_theta * sin_phi) - FT_phase = self.dx * self.dy * FT_phase_x * FT_phase_y - - # multiply the phase terms with the current sources - Jx_phased = np.sum(self.Jx * FT_phase, axis=(-2, -1)) - Jy_phased = np.sum(self.Jy * FT_phase, axis=(-2, -1)) - Mx_phased = np.sum(self.Mx * FT_phase, axis=(-2, -1)) - My_phased = np.sum(self.My * FT_phase, axis=(-2, -1)) - - # get rid of extra dimensions - cos_phi = np.squeeze(cos_phi) - sin_phi = np.squeeze(sin_phi) - cos_theta = np.squeeze(cos_theta) - sin_theta = np.squeeze(sin_theta) + phase_x = np.exp(1j * self.k0 * self.xx * sin_theta * cos_phi) + phase_y = np.exp(1j * self.k0 * self.yy * sin_theta * cos_phi) + phase = self.dx * self.dy * phase_x * phase_y + + Jx_k = np.sum(self.Jx * phase) + Jy_k = np.sum(self.Jy * phase) + Mx_k = np.sum(self.Mx * phase) + My_k = np.sum(self.My * phase) # N_theta (8.33a) - integrand_x = +Jx_phased * cos_theta * cos_phi - integrand_y = +Jy_phased * cos_theta * sin_phi - N_theta = integrand_x + integrand_y + N_theta = Jx_k * cos_theta * cos_phi + Jy_k * cos_theta * sin_phi # N_phi (8.33b) - integrand_x = -Jx_phased * sin_phi - integrand_y = +Jy_phased * cos_phi - N_phi = integrand_x + integrand_y + N_phi = -Jx_k * sin_phi + Jy_k * cos_phi # L_theta (8.34a) - integrand_x = +Mx_phased * cos_theta * cos_phi - integrand_y = +My_phased * cos_theta * sin_phi - L_theta = integrand_x + integrand_y + L_theta = Mx_k * cos_theta * cos_phi + My_k * cos_theta * sin_phi # L_phi (8.34b) - integrand_x = -Mx_phased * sin_phi - integrand_y = +My_phased * cos_phi - L_phi = integrand_x + integrand_y + L_phi = -Mx_k * sin_phi + My_k * cos_phi return N_theta, N_phi, L_theta, L_phi @@ -143,16 +135,23 @@ def fields_spherical(self, r, theta, phi): (Er, Etheta, Ephi), (Hr, Htheta, Hphi), fields in polar coordinates. """ + + # project radiation vectors to distance r away for given angles N_theta, N_phi, L_theta, L_phi = self._radiation_vectors(theta, phi) scalar_proj_r = 1j * self.k0 * np.exp(-1j * self.k0 * r) / (4 * np.pi * r) + + # assemble E felds E_theta = -scalar_proj_r * (L_phi + ETA_0 * N_theta) E_phi = scalar_proj_r * (L_theta - ETA_0 * N_phi) E_r = np.zeros_like(E_phi) E = np.stack((E_r, E_theta, E_phi)) + + # assemble H fields H_theta = -E_phi / ETA_0 H_phi = E_theta / ETA_0 H_r = np.zeros_like(H_phi) H = np.stack((H_r, H_theta, H_phi)) + return E, H def fields_cartesian(self, x, y, z):