From d90ed713f34ed6ef1ca41c8c3eb314bd92121a4d Mon Sep 17 00:00:00 2001
From: Ardavan Oskooi <ardavano@google.com>
Date: Fri, 8 Dec 2023 05:13:57 +0000
Subject: [PATCH] fix bugs in script

---
 .../Near_to_Far_Field_Spectra.md              | 20 ++++++++++---------
 python/examples/disc_extraction_efficiency.py | 20 ++++++++++---------
 2 files changed, 22 insertions(+), 18 deletions(-)

diff --git a/doc/docs/Python_Tutorials/Near_to_Far_Field_Spectra.md b/doc/docs/Python_Tutorials/Near_to_Far_Field_Spectra.md
index a5fd9c041..36fec9ea3 100644
--- a/doc/docs/Python_Tutorials/Near_to_Far_Field_Spectra.md
+++ b/doc/docs/Python_Tutorials/Near_to_Far_Field_Spectra.md
@@ -671,6 +671,8 @@ WAVELENGTH_UM = 1.0  # wavelength (in vacuum)
 # radius of quarter circle for computing flux in far field
 FARFIELD_RADIUS_UM = 1000 * WAVELENGTH_UM
 
+RESOLUTION = 50  # pixels/μm
+
 
 def plot_radiation_pattern_polar(theta_rad: np.ndarray, radial_flux: np.ndarray):
     """Plots the radiation pattern in polar coordinates.
@@ -806,8 +808,6 @@ def dipole_in_disc(
     Returns:
       A 2-tuple of the total flux and the radiation pattern.
     """
-    resolution = 50  # pixels/μm
-
     t_pml_um = 0.5  # thickness of PML
     t_air_um = 1.0  # thickness of air padding above disc
     length_r_um = 5.0  # length of cell in r
@@ -826,11 +826,6 @@ def dipole_in_disc(
         mp.PML(t_pml_um, direction=mp.Z, side=mp.High),
     ]
 
-    # An Er source at r=0 needs to be slighty offset.
-    # https://github.com/NanoComp/meep/issues/2704
-    if rpos == 0:
-        rpos = 1.5 / resolution
-
     src_cmpt = mp.Er
     src_pt = mp.Vector3(rpos, 0, -0.5 * size_z + h_disc * t_disc_um)
     sources = [
@@ -850,7 +845,7 @@ def dipole_in_disc(
     ]
 
     sim = mp.Simulation(
-        resolution=resolution,
+        resolution=RESOLUTION,
         cell_size=cell_size,
         dimensions=mp.CYLINDRICAL,
         m=m,
@@ -885,7 +880,7 @@ def dipole_in_disc(
         ),
     )
 
-    delta_vol = 2 * np.pi * rpos / (resolution**2)
+    delta_vol = 2 * np.pi * rpos / (RESOLUTION**2)
     dipole_flux = -np.real(sim.ldos_Fdata[0] * np.conj(sim.ldos_Jdata[0])) * delta_vol
 
     dipole_radiation_pattern = radiation_pattern(theta_rad, sim, n2f_mon)
@@ -899,6 +894,10 @@ if __name__ == "__main__":
     num_dipoles = 11
     dipole_rpos = np.linspace(0, DISC_RADIUS_UM, num_dipoles)
 
+    # An Er source at r=0 needs to be slighty offset.
+    # https://github.com/NanoComp/meep/issues/2704
+    dipole_rpos[0] = 1.5 / RESOLUTION
+
     delta_rpos = dipole_rpos[2] - dipole_rpos[1]
 
     # number of angular grid points in [0, π/2]
@@ -967,6 +966,9 @@ if __name__ == "__main__":
             f"dipole:, {rpos:.4f}, {m}, " f"{dipole_radiation_pattern_total_flux:.6f}"
         )
 
+    flux_total /= num_dipoles
+    radiation_pattern_total /= num_dipoles
+
     radiation_pattern_total_flux = radiation_pattern_flux(
         theta_rad, radiation_pattern_total
     )
diff --git a/python/examples/disc_extraction_efficiency.py b/python/examples/disc_extraction_efficiency.py
index 5f7979886..b6e5446d4 100644
--- a/python/examples/disc_extraction_efficiency.py
+++ b/python/examples/disc_extraction_efficiency.py
@@ -21,6 +21,8 @@
 # radius of quarter circle for computing flux in far field
 FARFIELD_RADIUS_UM = 1000 * WAVELENGTH_UM
 
+RESOLUTION = 50  # pixels/μm
+
 
 def plot_radiation_pattern_polar(theta_rad: np.ndarray, radial_flux: np.ndarray):
     """Plots the radiation pattern in polar coordinates.
@@ -156,8 +158,6 @@ def dipole_in_disc(
     Returns:
       A 2-tuple of the total flux and the radiation pattern.
     """
-    resolution = 50  # pixels/μm
-
     t_pml_um = 0.5  # thickness of PML
     t_air_um = 1.0  # thickness of air padding above disc
     length_r_um = 5.0  # length of cell in r
@@ -176,11 +176,6 @@ def dipole_in_disc(
         mp.PML(t_pml_um, direction=mp.Z, side=mp.High),
     ]
 
-    # An Er source at r=0 needs to be slighty offset.
-    # https://github.com/NanoComp/meep/issues/2704
-    if rpos == 0:
-        rpos = 1.5 / resolution
-
     src_cmpt = mp.Er
     src_pt = mp.Vector3(rpos, 0, -0.5 * size_z + h_disc * t_disc_um)
     sources = [
@@ -200,7 +195,7 @@ def dipole_in_disc(
     ]
 
     sim = mp.Simulation(
-        resolution=resolution,
+        resolution=RESOLUTION,
         cell_size=cell_size,
         dimensions=mp.CYLINDRICAL,
         m=m,
@@ -235,7 +230,7 @@ def dipole_in_disc(
         ),
     )
 
-    delta_vol = 2 * np.pi * rpos / (resolution**2)
+    delta_vol = 2 * np.pi * rpos / (RESOLUTION**2)
     dipole_flux = -np.real(sim.ldos_Fdata[0] * np.conj(sim.ldos_Jdata[0])) * delta_vol
 
     dipole_radiation_pattern = radiation_pattern(theta_rad, sim, n2f_mon)
@@ -249,6 +244,10 @@ def dipole_in_disc(
     num_dipoles = 11
     dipole_rpos = np.linspace(0, DISC_RADIUS_UM, num_dipoles)
 
+    # An Er source at r=0 needs to be slighty offset.
+    # https://github.com/NanoComp/meep/issues/2704
+    dipole_rpos[0] = 1.5 / RESOLUTION
+
     delta_rpos = dipole_rpos[2] - dipole_rpos[1]
 
     # number of angular grid points in [0, π/2]
@@ -317,6 +316,9 @@ def dipole_in_disc(
             f"dipole:, {rpos:.4f}, {m}, " f"{dipole_radiation_pattern_total_flux:.6f}"
         )
 
+    flux_total /= num_dipoles
+    radiation_pattern_total /= num_dipoles
+
     radiation_pattern_total_flux = radiation_pattern_flux(
         theta_rad, radiation_pattern_total
     )