From f0c2574ff43104e348c17ec4c54c134d24ffc7dd Mon Sep 17 00:00:00 2001
From: zymon <s@zymon.org>
Date: Tue, 11 Jun 2024 14:45:29 +0200
Subject: [PATCH] ISM rectangular post-filter refactor

---
 src/ISM.jl   | 14 ++------------
 src/utils.jl | 13 +++++--------
 2 files changed, 7 insertions(+), 20 deletions(-)

diff --git a/src/ISM.jl b/src/ISM.jl
index 6c7cf34..04b4c41 100644
--- a/src/ISM.jl
+++ b/src/ISM.jl
@@ -39,16 +39,6 @@ function ISM(
     h = zeros(config.N)
     # TODO: h = zeros(T, config.N)
     ISM!(h, tx, rx, room, config)
-
-    if config.hp
-        return AllenBerkley_highpass100(h, config.fs)
-        # TODO: Tego nie ma w ISM!
-        # TODO: Zmień to na funkcie operującą na zaalokowanym już h
-        # AllenBerkley_highpass100!(h, config.fs)
-        # return h
-    else
-        return h
-    end
 end
 
 
@@ -158,16 +148,16 @@ function ISM!(
     room::RectangularRoom,
     config::ISMConfig;
 ) where {T<:AbstractFloat}
-
     # Impulse parameters
     Δt = 1 / config.fs
     ω = 2π / config.Wd
     twid = π * config.fs
-
     for image_source ∈ image_source_generator(tx, rx, room, config)
         image_source |> isnothing && continue
         A, τ, _, _ = image_source
         insert_impuse!(h, A, τ, config.Wd, config.fs, Δt, twid, ω, config.N)
     end
+    config.hp && AllenBerkley_highpass100!(h, config.fs)
+    return h
 end
 
diff --git a/src/utils.jl b/src/utils.jl
index 3bcd509..c953f82 100644
--- a/src/utils.jl
+++ b/src/utils.jl
@@ -35,29 +35,26 @@ function Sabine_RT60(T60, L::Tuple, c)
     sqrt(1-α)
 end
 
-
-
 """
 b = [1, -B1, -B2]
 a = [1, A1, R1] 
 """
-function AllenBerkley_highpass100(x, fs)
-    o = x .* 0
+function AllenBerkley_highpass100!(x, fs)
     Y = zeros(3)
 
-    W = 2π*100/fs
+    W = 2π / fs * 100
 
     R1 = exp(-W)
     B1 = 2*R1*cos(W)
-    B2 = -R1 * R1   
+    B2 = -R1 * R1
     A1 = -(1+R1)
 
     for i = eachindex(x)
         Y[3] = Y[2]
         Y[2] = Y[1]
         Y[1] = B2*Y[3] + B1*Y[2] + x[i]
-        o[i] = Y[1] + A1*Y[2] + R1*Y[3]
+        x[i] = Y[1] + A1*Y[2] + R1*Y[3]
     end
 
-    return o
+    return x
 end