Choosing the Filter for a Weighing Scale

A kitchen scale samples its load cell fs times per second. The weight of what is put on it changes slowly (below fpass hertz), but the readings also pick up vibration from the worktop and hum from the mains, all above fstop hertz. The firmware will smooth the readings with an FIR filter, and the designer has several candidate weight lists candidates. A longer filter costs more computing and delays the reading, so the designer wants the shortest one that meets the specification:

• passband: at every checked frequency f <= fpass, the gain must be within ripple_db decibels of 0 dB: abs(gain_db) <= ripple_db (the weight is shown faithfully),
• stopband: at every checked frequency f >= fstop, the gain must be at most -atten_db decibels: gain_db <= -atten_db (vibration and hum are suppressed).

The frequencies checked are the grid f_i = i * fs / 400 for i = 0, 1, ..., 200, from 0 Hz up to fs / 2. The gain in decibels is computed as in "What the Drone's Sensor Filter Lets Through": w = 2 * pi * f / fs, H = sum of b[k] * exp(-1j * w * k), gain_db = 20 * log10(|H|), and -120.0 if |H| is below 1e-6.

Write pick_filter(candidates, fs, fpass, fstop, ripple_db, atten_db) that returns a pair:

1. the index in candidates of the shortest filter that meets both limits (the lowest index among equally short ones), or -1 if none does,
2. a list with one pair (deviation, stop_level) per candidate, in order: the largest abs(gain_db) over the passband frequencies of the grid, and the largest gain_db over its stopband frequencies.

The tests build some candidates with the helper windowed_sinc(taps, cutoff, fs), which returns a standard low-pass design with taps weights; you may use it with Run to make your own. To see why a candidate passes or fails, press Run with plot(grid, levels, title="response"), where levels holds its gain in dB at each grid frequency.

Examples

Input:  candidates = [[1], [0.5, 0.5]], fs = 1000, fpass = 101, fstop = 401,
        ripple_db = 1, atten_db = 10
Output: (1, [(0.0, 0.0), (0.4358734890997429, -10.413160154277307)])
Explanation: [1] passes everything unchanged, so it suppresses nothing. The average of two
neighbours loses at most 0.44 dB below 101 Hz and is at least 10.4 dB down from 402.5 Hz
(the first grid frequency at or above 401 Hz), so it meets both limits.

Input:  candidates = [[0.125] * 8, windowed_sinc(31, 40, 1000), windowed_sinc(61, 40, 1000),
                      windowed_sinc(41, 60, 1000)], fs = 1000, fpass = 21, fstop = 101,
        ripple_db = 1, atten_db = 40
Output: (3, [(0.3631805585590642, -12.79751814357273), (1.2218600203142347, -43.42691016340525),
             (0.4865088802248915, -61.600743633573806), (0.09053541206330945, -43.82932795363644)])
Explanation: the moving average is far from 40 dB down; the 31-weight design cuts too close
and loses 1.2 dB of the passband; both 41 and 61 weights pass, and 41 is shorter.

Constraints

• 0 <= fpass < fstop <= fs / 2, so the passband and the stopband each contain at least one grid frequency
• at most 10 candidates of at most 101 weights each
• fpass and fstop never fall exactly on a grid frequency, and no deviation or stop level lies within 1e-6 of its limit
• answers are compared with a tolerance of 1e-6

Goals

• Judge a filter by its frequency response on a grid of frequencies
• Check a passband limit and a stopband limit, both in decibels
• Pick the cheapest filter (fewest weights) that meets the specification