1# /// script
2# dependencies = [
3#   "matplotlib",
4#   "scipy",
5# ]
6# ///
7
8"""
9Generate a FIR low-pass filter using the Remez exchange algorithm, and plot the frequency response.
10"""
11
12import argparse
13from math import pi
14
15import matplotlib.pyplot as plt
16import numpy as np
17from scipy.signal import freqz, remez
18
19
20def main():
21    args = parse_args()
22
23    weight_pb = 1 - 10 ** (-args.pbr / 20)
24    weight_sb = 10 ** (-args.sba / 20)
25    coefficients = remez(
26        args.n,
27        [0, args.fp, args.fe, args.fs / 2],
28        [1, 0],
29        weight=[weight_pb, weight_sb],
30        fs=args.fs,
31    )
32    print([float(coeff) for coeff in coefficients])
33
34    w, h = freqz(coefficients, worN=1 << 20)
35
36    figure = plt.figure()
37    axes = figure.add_subplot(
38        title="Filter Frequency Response",
39        xlabel="Frequency (Hz)",
40        ylabel="Gain (dB)",
41        xlim=(0, 40000),
42        ylim=(-80, 20),
43    )
44    axes.grid(visible=True)
45
46    axes.plot(w * (args.fs / 2) / pi, 20 * np.log10(abs(h)))
47    axes.plot(2 * [24000], [-100, 100])
48    axes.plot([0, args.fs], 2 * [-3])
49
50    plt.show(block=True)
51
52
53def parse_args():
54    arg_parser = argparse.ArgumentParser(
55        description="Generate a FIR low-pass filter using the Remez exchange algorithm"
56    )
57    arg_parser.add_argument(
58        "-fs", type=float, required=True, help="Source frequency (Hz)"
59    )
60    arg_parser.add_argument("-n", type=int, required=True, help="Number of taps")
61    arg_parser.add_argument(
62        "-fp",
63        default=16000,
64        type=float,
65        required=False,
66        help="Transition band start frequency (Hz) (default=16000)",
67    )
68    arg_parser.add_argument(
69        "-fe",
70        default=24000,
71        type=float,
72        required=False,
73        help="Stopband edge frequency (Hz) (default=24000)",
74    )
75    arg_parser.add_argument(
76        "-pbr",
77        default=0.1,
78        type=float,
79        required=False,
80        help="Desired passband ripple (dB) (default=0.1)",
81    )
82    arg_parser.add_argument(
83        "-sba",
84        default=60,
85        type=float,
86        required=False,
87        help="Desired stopband attenuation (dB) (default=60)",
88    )
89    return arg_parser.parse_args()
90
91
92if __name__ == "__main__":
93    main()