makeweight: always use geometric mean for default gain_mid

src/weights.jl

@@ -6,7 +6,7 @@ Create a weighting function that goes from gain `low` at zero frequency, through
 
 # Arguments:
 - `low`: A number specifying the DC gain 
-- `mid`: A number specifying the frequency at which the gain is 1, or a tuple `(freq, gain)`. If `gain_mid` is not specified, the geometric mean of `high` and `low` is used.
+- `mid`: A number specifying the frequency at which the gain is 1, or a tuple `(freq, gain)`. If `gain_mid` is not specified, it defaults to `1` whenever `low` and `high` straddle 1 (the usual case), and to the geometric mean `√(low*high)` otherwise.
 - `high`: A number specifying the gain at ∞
 
 ```@example
@@ -18,7 +18,14 @@ vline!([5], l=(:black, :dash), primary=false)
 ```
 """
 function makeweight(low, mid::Number, high)
-    makeweight(low, (mid, high < 1 ? √(high*low) : 1), high)
+    # Default `gain_mid` must lie strictly between `low` and `high` for the
+    # downstream formula to produce real coefficients. When `low` and `high`
+    # straddle 1 the historical default `gain_mid = 1` satisfies that and is
+    # preserved (back-compat). When both are on the same side of 1, fall back
+    # to the geometric mean `√(low*high)` so the previously-broken case where
+    # both `low > 1` and `high > 1` no longer produces complex/NaN poles.
+    gain_mid = (low <= 1 && high >= 1) ? one(promote_type(typeof(low), typeof(high))) : √(high*low)
+    makeweight(low, (mid, gain_mid), high)
 end
 
 

test/test_weights.jl

@@ -12,6 +12,25 @@ w = gain_and_delay_uncertainty(1, 2, 1)
 w = makeweight(0.1, 1, 2)
 @test dcgain(w)[] ≈ 0.1
 @test evalfr(w, 10000im)[] ≈ 2 atol=1e-3
+# `low` and `high` straddle 1, so the default gain_mid is 1 (preserved back-compat).
+@test abs(evalfr(w, 1im)[]) ≈ 1 atol=1e-6
+
+# Regression: previously, when both `low > 1` and `high > 1`, the default gain_mid was
+# hard-coded to 1 which violates the formula's `low < mag < high` precondition and
+# produced complex/NaN poles. Now falls back to √(low*high), which lies between them.
+w_both_gt_1 = makeweight(2, 1, 3)
+@test all(isreal, denvec(tf(w_both_gt_1))[1])
+@test all(isreal, numvec(tf(w_both_gt_1))[1])
+@test dcgain(w_both_gt_1)[] ≈ 2
+@test evalfr(w_both_gt_1, 10000im)[] ≈ 3 atol=1e-3
+@test abs(evalfr(w_both_gt_1, 1im)[]) ≈ √(2*3) atol=1e-6
+
+# Symmetric same-side case: both below 1. Default gain_mid = √(low*high).
+w_both_lt_1 = makeweight(0.5, 1, 0.1)
+@test all(isreal, denvec(tf(w_both_lt_1))[1])
+@test dcgain(w_both_lt_1)[] ≈ 0.5
+@test evalfr(w_both_lt_1, 10000im)[] ≈ 0.1 atol=1e-3
+@test abs(evalfr(w_both_lt_1, 1im)[]) ≈ √(0.5*0.1) atol=1e-6
 
 
 w = neglected_lag(1)