Optimize ivas_filter_process_fx part 1.

    Word16 state_e[IVAS_FILTER_MAX_STAGES][IVAS_BIQUAD_FILT_LEN];

#endif

#ifdef OPT_2239_IVAS_FILTER_PROCESS

    Word16 q_diff; // q_diff = q_out - q_in -> q_out = q_diff + q_in

#endif

} ivas_filters_process_state_t;

import numpy as np

from scipy import signal

import matplotlib.pyplot as plt

def generate_biquad_butterworth_coeffs(cutoff_freq, fs, order=2, ftype='lowpass', output='ba'):

    Wn = cutoff_freq

    if output == 'sos':

        sos_coeffs = signal.butter(order, Wn, btype=ftype, analog=False, output='sos', fs=fs)

        return sos_coeffs

    elif output == 'ba':

        b_coeffs, a_coeffs = signal.butter(order, Wn, btype=ftype, analog=False, output='ba', fs=fs)

        return b_coeffs, a_coeffs

    else:

        raise ValueError("Output format must be 'sos' or 'ba'")

def plot_biquad_frequency_response(b, a, fs):

    w, h = signal.freqz(b, a, fs=fs)

    magnitude_db = 20 * np.log10(np.abs(h))

    phase_degrees = np.unwrap(np.angle(h)) * 180 / np.pi

    fig, (ax1, ax2) = plt.subplots(2, 1, sharex=True, tight_layout=True)

    ax1.plot(w, magnitude_db, 'b')

    ax1.set_ylabel('Magnitude (dB)', color='b')

    ax1.set_xlabel('Frequency (Hz)')

    ax1.grid()

    ax2.plot(w, phase_degrees, 'g')

    ax2.set_ylabel('Phase (degrees)', color='g')

    ax2.set_xlabel('Frequency (Hz)')

    ax2.grid()

    plt.show()

def max_abs_output_biquad(b, a, n_samples=10000):

    if a[0] != 1.0:

        b = np.array(b) / a[0]

        a = np.array(a) / a[0]

        a[0] = 1.0

    impulse_input = np.zeros(n_samples)

    impulse_input[0] = 1.0

    impulse_response = signal.lfilter(b, a, impulse_input)

    max_gain = np.sum(np.abs(impulse_response))

    return max_gain

# sampling_rate = 44100.0 # Hz

# cutoff_frequency = 5000.0 # Hz

# b_coeffs,a_coeffs = generate_biquad_butterworth_coeffs(cutoff_frequency, sampling_rate)

# plot_biquad_frequency_response(b_coeffs, a_coeffs, sampling_rate)

# # Calculate the maximum gain

# gain = max_abs_output_biquad(b_coeffs, a_coeffs)

# print(f"Numerator coefficients (b): {b_coeffs}")

# print(f"Denominator coefficients (a): {a_coeffs}")

# print(f"Maximum possible time-domain gain (L1 norm): {gain:.4f}")

IVAS_C_HPF_48k = 0.675231906655777

IVAS_C_HPF_32k = 0.554854910159853

IVAS_C_HPF_16k = 0.307863971328499

IVAS_C_FAST_48k = 0.995842001845110

IVAS_C_FAST_32k = 0.993769490623395

IVAS_C_FAST_16k = 0.987577800493881

IVAS_C_SLOW_48k = 0.999739617238810

IVAS_C_SLOW_32k = 0.999609451284012

IVAS_C_SLOW_16k = 0.999219055096324

#env_hpf

sampling_rate = 48000.0 # Hz

b_coeffs = [ +IVAS_C_HPF_48k, -IVAS_C_HPF_48k, 0.0 ]

a_coeffs = [ 1.0, -IVAS_C_HPF_48k, 0.0 ]

plot_biquad_frequency_response(b_coeffs, a_coeffs, sampling_rate)

gain = max_abs_output_biquad(b_coeffs, a_coeffs)

print(f"Maximum possible time-domain gain (L1 norm): {gain:.4f}")

#env_fast

sampling_rate = 48000.0 # Hz

b_coeffs = [ 1.0 - IVAS_C_FAST_48k, 0.0, 0.0 ]

a_coeffs = [ 1.0, -IVAS_C_FAST_48k, 0.0 ]

plot_biquad_frequency_response(b_coeffs, a_coeffs, sampling_rate)

gain = max_abs_output_biquad(b_coeffs, a_coeffs)

print(f"Maximum possible time-domain gain (L1 norm): {gain:.4f}")

#env_slow

sampling_rate = 48000.0 # Hz

b_coeffs = [ 1.0 - IVAS_C_SLOW_48k, 0.0, 0.0 ]

a_coeffs = [ 1.0, -IVAS_C_SLOW_48k, 0.0 ]

plot_biquad_frequency_response(b_coeffs, a_coeffs, sampling_rate)

gain = max_abs_output_biquad(b_coeffs, a_coeffs)

print(f"Maximum possible time-domain gain (L1 norm): {gain:.4f}")

/* #################### End BASOP optimization switches ############################ */

/* #################### Start BASOP optimization switches ############################ */

#define OPT_2239_IVAS_FILTER_PROCESS                   /* Dolby: Issue 2239, optimize ivas_filter_process_fx. */

/* #################### End BASOP optimization switches ############################ */

#endif