fix(propagate_filter): avoid floating point issues with small negativ…

…e uncertainties via clipping
PTB-M4D · Oct 13, 2021 · bbe9d13 · bbe9d13
1 parent 20546d3
commit bbe9d13
Showing 1 changed file with 12 additions and 6 deletions.
diff --git a/src/PyDynamic/uncertainty/propagate_filter.py b/src/PyDynamic/uncertainty/propagate_filter.py
@@ -110,18 +110,20 @@ def _fir_filter(x, theta, Ux=None, Utheta=None, initial_conditions="constant"):
             )
 
         # calc subterm theta^T * Ux * theta
-        Uy += convolve(np.outer(theta, theta), Ux_extended, mode="valid")
+        Uy += convolve(np.outer(theta, theta), Ux_extended, mode="valid").clip(min=0)
 
     if Utheta is not None:
         ## extend signal Ntheta steps into the past
         x_extended = np.r_[np.full((Ntheta - 1), x0), x]
 
         # calc subterm x^T * Utheta * x
-        Uy += convolve(np.outer(x_extended, x_extended), Utheta, mode="valid")
+        Uy += convolve(np.outer(x_extended, x_extended), Utheta, mode="valid").clip(
+            min=0
+        )
 
     if (Ux is not None) and (Utheta is not None):
         # calc subterm Tr(Ux * Utheta)
-        Uy += convolve(Ux_extended, Utheta.T, mode="valid")
+        Uy += convolve(Ux_extended, Utheta.T, mode="valid").clip(min=0)
 
     return y, Uy
 
@@ -206,18 +208,22 @@ def _fir_filter_diag(
         Ux_diag_extended = np.r_[np.full((Ntheta - 1), Ux0), Ux_diag]
 
         # calc subterm theta^T * Ux * theta
-        Uy_diag += convolve(np.square(theta), Ux_diag_extended, mode="valid")
+        Uy_diag += convolve(np.square(theta), Ux_diag_extended, mode="valid").clip(
+            min=0
+        )
 
     if Utheta_diag is not None:
         ## extend signal Ntheta steps into the past
         x_extended = np.r_[np.full((Ntheta - 1), x0), x]
 
         # calc subterm x^T * Utheta * x
-        Uy_diag += convolve(np.square(x_extended), Utheta_diag, mode="valid")
+        Uy_diag += convolve(np.square(x_extended), Utheta_diag, mode="valid").clip(
+            min=0
+        )
 
     if (Ux_diag is not None) and (Utheta_diag is not None):
         # calc subterm Tr(Ux * Utheta)
-        Uy_diag += convolve(Ux_diag_extended, Utheta_diag, mode="valid")
+        Uy_diag += convolve(Ux_diag_extended, Utheta_diag, mode="valid").clip(min=0)
 
     return y, Uy_diag