refactor signal generation to support single-channel mode and adjust birefringence parameters

src/single-core-data-gen/add_pypho.py

@@ -26,7 +26,7 @@ while not (__parent_dir / "pypho" / "pypho").exists() and __parent_dir != Path("
     __parent_dir = __parent_dir.parent
 
 if __parent_dir != Path("/"):
-    sys.path.append(str(__parent_dir / "pypho"))
+    sys.path.insert(0, str(__parent_dir / "pypho"))
     __log.append(f"Added '{__parent_dir/ "pypho"}' to 'PATH'")
 else:
     __log.append('pypho not found')

src/single-core-data-gen/generate_signal.py

@@ -38,7 +38,8 @@ gamma = 1.14
 alpha = 0.2
 D = 17
 S = 0
-birefsteps = 1
+birefsteps = 0
+max_delta_beta = 0.4
 ; birefseed = 0xC0FFEE
 
 [signal]
@@ -97,6 +98,7 @@ class pam_generator:
         pulse_shape="gauss",
         fwhm=0.33,
         seed=None,
+        single_channel=False
     ) -> None:
         self.glova = glova
         self.pulse_shape = pulse_shape
@@ -104,46 +106,51 @@ class pam_generator:
         self.mod_order = mod_order
         self.fwhm = fwhm
         self.seed = seed
+        self.single_channel = single_channel
 
     def __call__(self, E, symbols, max_jitter=0):
-        symbols_x = symbols[0] / (self.mod_order or np.max(symbols[0]))
-        symbols_y = symbols[1] / (self.mod_order or np.max(symbols[1]))
-
-        diffs_x = np.diff(symbols_x, prepend=symbols_x[0])
-        diffs_y = np.diff(symbols_y, prepend=symbols_y[0])
-
         max_jitter = int(round(max_jitter * self.glova.sps))
-        digital_x = self.generate_digital_signal(diffs_x, max_jitter)
-        digital_y = self.generate_digital_signal(diffs_y, max_jitter)
-
-        digital_x = np.pad(
-            digital_x, (0, self.glova.sps // 2), "constant", constant_values=(0, 0)
-        )
-        digital_y = np.pad(
-            digital_y, (0, self.glova.sps // 2), "constant", constant_values=(0, 0)
-        )
-
         if self.pulse_shape == "gauss":
             wavelet = self.gauss(oversampling=6)
         else:
             raise ValueError(f"Unknown pulse shape: {self.pulse_shape}")
-
+        
+        # prepare symbols
+        symbols_x = symbols[0] / (self.mod_order or np.max(symbols[0]))
+        diffs_x = np.diff(symbols_x, prepend=symbols_x[0])
+        digital_x = self.generate_digital_signal(diffs_x, max_jitter)
+        digital_x = np.pad(
+            digital_x, (0, self.glova.sps // 2), "constant", constant_values=(0, 0)
+        )
+        
         # create analog signal of diff of symbols
         E_x = np.convolve(digital_x, wavelet)
-        E_y = np.convolve(digital_y, wavelet)
-
+        
         # convert to pam and set modulation depth (scale and move up such that 1 stays at 1)
         E_x = np.cumsum(E_x) * self.modulation_depth + (1 - self.modulation_depth)
-        E_y = np.cumsum(E_y) * self.modulation_depth + (1 - self.modulation_depth)
-
+        
         # cut off the wavelet tails
         E_x = E_x[self.glova.sps // 2 + len(wavelet) // 2 - 1 : -len(wavelet) // 2]
-        E_y = E_y[self.glova.sps // 2 + len(wavelet) // 2 - 1 : -len(wavelet) // 2]
-
+        
         # modulate the laser
         E[0]["E"][0] = np.sqrt(np.multiply(np.square(E[0]["E"][0]), E_x))
-        E[0]["E"][1] = np.sqrt(np.multiply(np.square(E[0]["E"][1]), E_y))
+        
+        if not self.single_channel:
+            symbols_y = symbols[1] / (self.mod_order or np.max(symbols[1]))
+            diffs_y = np.diff(symbols_y, prepend=symbols_y[0])
+            digital_y = self.generate_digital_signal(diffs_y, max_jitter)
+            digital_y = np.pad(
+                digital_y, (0, self.glova.sps // 2), "constant", constant_values=(0, 0)
+            )
+            E_y = np.convolve(digital_y, wavelet)
 
+            E_y = np.cumsum(E_y) * self.modulation_depth + (1 - self.modulation_depth)
+
+            E_y = E_y[self.glova.sps // 2 + len(wavelet) // 2 - 1 : -len(wavelet) // 2]
+
+            E[0]["E"][1] = np.sqrt(np.multiply(np.square(E[0]["E"][1]), E_y))
+        else:
+            E[0]["E"][1] = np.zeros_like(E[0]["E"][0], dtype=E[0]["E"][0].dtype)
         return E
 
     def generate_digital_signal(self, symbols, max_jitter=0):
@@ -212,6 +219,7 @@ def initialize_fiber_and_data(config, input_data_override=None):
             pulse_shape=config["signal"]["pulse_shape"],
             fwhm=config["signal"]["fwhm"],
             seed=config["signal"]["jitter_seed"],
+            single_channel=False
         )
 
         symbols_x = symbolsrc(pattern="random")
@@ -243,8 +251,9 @@ def initialize_fiber_and_data(config, input_data_override=None):
         py_fiber.birefarray = pypho.birefringence_segment.create_pmd_fibre(
             py_fiber.l,
             py_fiber.l / config["fiber"]["birefsteps"],
-            0,
-            config["fiber"]["birefseed"],
+            # maxDeltaD=config["fiber"]["d"]/5,
+            maxDeltaBeta = config["fiber"].get("max_delta_beta", 0),
+            seed=config["fiber"]["birefseed"],
         )
     c_params = pypho.cfiber.ParamsWrapper.from_fiber(
         py_fiber, max_step=1e3 if py_fiber.gamma == 0 else 200
@@ -288,6 +297,7 @@ def save_data(data, config):
         f"D = {config['fiber']['d']}",
         f"S = {config['fiber']['s']}",
         f"birefsteps = {config['fiber']['birefsteps']}",
+        f"max_delta_beta = {config['fiber'].get('max_delta_beta', 0)}",
         f"birefseed = {hex(birefseed)}" if birefseed else "; birefseed = not set",
         "",
         "[signal]",
@@ -327,6 +337,7 @@ def save_data(data, config):
         config["fiber"]["s"],
         f"{config['signal']['modulation'].upper()}{config['signal']['mod_order']}",
         config["fiber"]["birefsteps"],
+        config["fiber"].get("max_delta_beta", 0),
     )
 
     lookup_file = "-".join(map(str, filename_components)) + ".ini"
@@ -339,9 +350,10 @@ def save_data(data, config):
     print("Saved data to", save_dir / save_file)
 
 
-def length_loop(config, lengths, incremental=False):
+def length_loop(config, lengths, incremental=False, bireflength=None, save=True):
     lengths = sorted(lengths)
     input_override = None
+    birefsteps_running = 0
     for lind, length in enumerate(lengths):
         # print(f"\nGenerating data for fiber length {length}")
         if lind > 0 and incremental:
@@ -354,6 +366,9 @@ def length_loop(config, lengths, incremental=False):
         else:
             print(f"\nGenerating data for fiber length {length}m")
         config["fiber"]["length"] = length
+        if bireflength is not None:
+            config["fiber"]["birefsteps"] = length // bireflength
+            birefsteps_running += config["fiber"]["birefsteps"]
         # set the input data to the output data of the previous run
         cfiber, cdata, noise, edfa = initialize_fiber_and_data(
             config, input_data_override=input_override
@@ -379,11 +394,14 @@ def length_loop(config, lengths, incremental=False):
             input_override = (cdata.E_out, noise)
             cdata.E_in = cdata_orig.E_in
             config["fiber"]["length"] = lengths[lind]
+            if bireflength is not None:
+                config["fiber"]["birefsteps"] = birefsteps_running
         
         E_tmp = [{'E': cdata.E_out, 'noise': noise*(-cfiber.params.l*cfiber.params.alpha)}]
         E_tmp = edfa(E=E_tmp)
         cdata.E_out = E_tmp[0]['E']
-        save_data(cdata, config)
+        if save:
+            save_data(cdata, config)
 
     in_out_eyes(cfiber, cdata)
 
@@ -502,8 +520,10 @@ if __name__ == "__main__":
     add_pypho.show_log()
     config = get_config()
 
-    lengths = np.arange(50000, 100000+6000, 1000)
-    length_loop(config, lengths, incremental=False)
+    # lengths = np.arange(50000, 150000+6000, 10000)
+    lengths = [10000]
+    length_loop(config, lengths, incremental=False, bireflength=1000, save=True)
+    # birefringence is constant over coupling length -> several 100m -> bireflength=1000 (m)
 
-    single_run_with_plot(config, save=False)
+    # single_run_with_plot(config, save=False)