fix Raman gain estimation during design

- Replaced multiple calls to the span_loss function
  with recording the span loss result in the fiber elements,
  reducing computation time.
- Updated Raman gain estimation based on design target powers to ensure
  accurate Edfa gain calculation or gain reduction during design.
- display the computed and design Raman gain in RamanFiber string output
- do not add padding on Raman fibers
- Added to_json function to preserve user input SimParams values,
  which were previously overwritten by initializing SimParams
  with fake values during design.

Next step is to allow users to balance computation time and
target accuracy of the design by inputing their own SimParams
and ref channels design values.

Signed-off-by: EstherLerouzic <[email protected]>
Change-Id: I1ca4954d0621858cefb3d776a538131992cae3e3

gnpy/core/elements.py

@@ -242,7 +242,6 @@ def __init__(self, *args, params=None, **kwargs):
         # on the path, since it depends on the equalization definition on the degree.
         self.ref_pch_out_dbm = None
         self.loss = 0  # auto-design interest
-
         # Optical power of carriers are equalized by the ROADM, so that the experienced loss is not the same for
         # different carriers. The ref_effective_loss records the loss for a reference carrier.
         self.ref_effective_loss = None
@@ -707,11 +706,16 @@ def __init__(self, *args, params=None, **kwargs):
     def to_json(self):
         return dict(super().to_json, operational=self.operational)
 
+    def __str__(self):
+        return super().__str__() + f'\n  reference gain (dB):         {round(self.estimated_gain, 2)}' \
+            + f'\n  actual gain (dB):            {round(self.actual_raman_gain, 2)}'
+
     def propagate(self, spectral_info: SpectralInformation):
         """Modifies the spectral information computing the attenuation, the non-linear interference generation,
         the CD and PMD accumulation.
         """
         # apply the attenuation due to the input connector loss
+        pin = watt2dbm(sum(spectral_info.signal))
         attenuation_in_db = self.params.con_in + self.params.att_in
         spectral_info.apply_attenuation_db(attenuation_in_db)
 
@@ -741,6 +745,8 @@ def propagate(self, spectral_info: SpectralInformation):
         spectral_info.apply_attenuation_db(attenuation_out_db)
         self.pch_out_dbm = watt2dbm(spectral_info.signal + spectral_info.nli + spectral_info.ase)
         self.propagated_labels = spectral_info.label
+        pout = watt2dbm(sum(spectral_info.signal))
+        self.actual_raman_gain = self.loss + pout - pin
 
 
 class Edfa(_Node):

gnpy/core/network.py

@@ -19,6 +19,7 @@
 from gnpy.core.info import ReferenceCarrier, create_input_spectral_information
 from gnpy.tools import json_io
 from gnpy.core.parameters import SimParams
+from gnpy.core.science_utils import RamanSolver
 
 
 logger = getLogger(__name__)
@@ -141,7 +142,6 @@ def target_power(network, node, equipment):  # get_fiber_dp
     POWER_SLOPE = 0.3
     dp_range = list(equipment['Span']['default'].delta_power_range_db)
     node_loss = span_loss(network, node, equipment)
-
     try:
         dp = round2float((node_loss - SPAN_LOSS_REF) * POWER_SLOPE, dp_range[2])
         dp = max(dp_range[0], dp)
@@ -187,63 +187,74 @@ def next_node_generator(network, node):
         yield from next_node_generator(network, next_node)
 
 
-def estimate_raman_gain(node, equipment):
+def estimate_raman_gain(node, equipment, power_dbm):
     """If node is RamanFiber, then estimate the possible Raman gain if any
-    for this purpose propagate a fake signal in a copy.
-    to be accurate the nb of channel should be the same as in SI, but this increases computation time
+    for this purpose computes stimulated_raman_scattering loss_profile. This may be time consuming.
     """
-    f_min = equipment['SI']['default'].f_min
-    f_max = equipment['SI']['default'].f_max
-    roll_off = equipment['SI']['default'].roll_off
-    baud_rate = equipment['SI']['default'].baud_rate
-    power_dbm = equipment['SI']['default'].power_dbm
-    power = dbm2watt(equipment['SI']['default'].power_dbm)
-    spacing = equipment['SI']['default'].spacing
-    tx_osnr = equipment['SI']['default'].tx_osnr
-
-    sim_params = {
-        "raman_params": {
-            "flag": True,
-            "result_spatial_resolution": 10e3,
-            "solver_spatial_resolution": 50
-        },
-        "nli_params": {
-            "method": "ggn_spectrally_separated",
-            "dispersion_tolerance": 1,
-            "phase_shift_tolerance": 0.1,
-            "computed_channels": [1, 18, 37, 56, 75]
-        }
-    }
     if isinstance(node, elements.RamanFiber):
+        if hasattr(node, "estimated_gain"):
+            return node.estimated_gain
+        f_min = equipment['SI']['default'].f_min
+        f_max = equipment['SI']['default'].f_max
+        roll_off = equipment['SI']['default'].roll_off
+        baud_rate = equipment['SI']['default'].baud_rate
+        power = dbm2watt(power_dbm)
+        spacing = equipment['SI']['default'].spacing
+        tx_osnr = equipment['SI']['default'].tx_osnr
+
+        # reduce the nb of channels to speed up
+        spacing = spacing * 3
+        power = power * 3
+
+        sim_params = {
+            "raman_params": {
+                "flag": True,
+                "result_spatial_resolution": 50e3,
+                "solver_spatial_resolution": 100
+            }
+        }
+
         # in order to take into account gain generated in RamanFiber, propagate in the RamanFiber with
-        # SI reference channel.
-        spectral_info_input = create_input_spectral_information(f_min=f_min, f_max=f_max, roll_off=roll_off,
-                                                                baud_rate=baud_rate, power=power, spacing=spacing,
-                                                                tx_osnr=tx_osnr)
-        n_copy = deepcopy(node)
-        # need to set ref_pch_in_dbm in order to correctly run propagate of the element, because this
-        # setting has not yet been done by autodesign
-        n_copy.ref_pch_in_dbm = power_dbm
+        if hasattr(node, "estimated_gain"):
+            # do not compute twice to save on time
+            return node.estimated_gain
+        spectral_info = create_input_spectral_information(f_min=f_min, f_max=f_max, roll_off=roll_off,
+                                                          baud_rate=baud_rate, power=power, spacing=spacing,
+                                                          tx_osnr=tx_osnr)
+        pin = watt2dbm(sum(spectral_info.signal))
+        attenuation_in_db = node.params.con_in + node.params.att_in
+        spectral_info.apply_attenuation_db(attenuation_in_db)
+        save_sim_params = {"raman_params": SimParams._shared_dict['raman_params'].to_json(),
+                           "nli_params": SimParams._shared_dict['nli_params'].to_json()}
         SimParams.set_params(sim_params)
-        pin = watt2dbm(sum(spectral_info_input.signal))
-        spectral_info_out = n_copy(spectral_info_input)
-        pout = watt2dbm(sum(spectral_info_out.signal))
-        estimated_gain = pout - pin + node.loss
+        stimulated_raman_scattering = RamanSolver.calculate_stimulated_raman_scattering(spectral_info, node)
+        attenuation_fiber = stimulated_raman_scattering.loss_profile[:spectral_info.number_of_channels, -1]
+        spectral_info.apply_attenuation_lin(attenuation_fiber)
+        attenuation_out_db = node.params.con_out
+        spectral_info.apply_attenuation_db(attenuation_out_db)
+        pout = watt2dbm(sum(spectral_info.signal))
+        estimated_loss = pin - pout
+        estimated_gain = node.loss - estimated_loss
+        node.estimated_gain = estimated_gain
+        SimParams.set_params(save_sim_params)
         return round(estimated_gain, 2)
     else:
         return 0.0
 
 
-def span_loss(network, node, equipment):
-    """Total loss of a span (Fiber and Fused nodes) which contains the given node"""
+def span_loss(network, node, equipment, input_power=None):
+    """Total loss of a span (Fiber and Fused nodes) which contains the given node
+    Do not recompute, if it was already computed: records it in design_span_loss"""
+    if hasattr(node, "design_span_loss"):
+        return node.design_span_loss
     loss = node.loss if node.passive else 0
     loss += sum(n.loss for n in prev_node_generator(network, node))
     loss += sum(n.loss for n in next_node_generator(network, node))
     # add the possible Raman gain
-    gain = estimate_raman_gain(node, equipment)
-    gain += sum(estimate_raman_gain(n, equipment) for n in prev_node_generator(network, node))
-    gain += sum(estimate_raman_gain(n, equipment) for n in next_node_generator(network, node))
-
+    gain = estimate_raman_gain(node, equipment, input_power)
+    gain += sum(estimate_raman_gain(n, equipment, input_power) for n in prev_node_generator(network, node))
+    gain += sum(estimate_raman_gain(n, equipment, input_power) for n in next_node_generator(network, node))
+    node.design_span_loss = loss - gain
     return loss - gain
 
 
@@ -399,7 +410,8 @@ def set_egress_amplifier(network, this_node, equipment, pref_ch_db, pref_total_d
                     node.target_pch_out_dbm = round(node.delta_p + pref_ch_db, 2)
                 elif node.delta_p is None:
                     node.target_pch_out_dbm = None
-
+            elif isinstance(node, elements.RamanFiber):
+                _ = span_loss(network, node, equipment, input_power=pref_ch_db + dp)
             prev_dp = dp
             prev_voa = voa
             prev_node = node
@@ -701,6 +713,9 @@ def add_fiber_padding(network, fibers, padding, equipment):
         next_node = get_next_node(fiber, network)
         if isinstance(next_node, elements.Fused):
             continue
+        # do not pad if this is a Raman Fiber
+        if isinstance(fiber, elements.RamanFiber):
+            continue
         this_span_loss = span_loss(network, fiber, equipment)
         if this_span_loss < padding:
             # add a padding att_in at the input of the 1st fiber:
@@ -710,6 +725,7 @@ def add_fiber_padding(network, fibers, padding, equipment):
             # just after a roadm: need to check that first_fiber is really a fiber
             if isinstance(first_fiber, elements.Fiber):
                 first_fiber.params.att_in = first_fiber.params.att_in + padding - this_span_loss
+                first_fiber.design_span_loss += first_fiber.params.att_in
 
 
 def add_missing_elements_in_network(network, equipment):

gnpy/core/parameters.py

@@ -46,6 +46,11 @@ def __init__(self, flag=False, result_spatial_resolution=10e3, solver_spatial_re
         self.result_spatial_resolution = result_spatial_resolution  # [m]
         self.solver_spatial_resolution = solver_spatial_resolution  # [m]
 
+    def to_json(self):
+        return {"flag": self.flag,
+                "result_spatial_resolution": self.result_spatial_resolution,
+                "solver_spatial_resolution": self.solver_spatial_resolution}
+
 
 class NLIParams(Parameters):
     def __init__(self, method='gn_model_analytic', dispersion_tolerance=1, phase_shift_tolerance=0.1,
@@ -62,6 +67,12 @@ def __init__(self, method='gn_model_analytic', dispersion_tolerance=1, phase_shi
         self.phase_shift_tolerance = phase_shift_tolerance
         self.computed_channels = computed_channels
 
+    def to_json(self):
+        return {"method": self.method,
+                "dispersion_tolerance": self.dispersion_tolerance,
+                "phase_shift_tolerance": self.phase_shift_tolerance,
+                "computed_channels": self.computed_channels}
+
 
 class SimParams(Parameters):
     _shared_dict = {'nli_params': NLIParams(), 'raman_params': RamanParams()}