Load Cell Probing Algorithm Testing

[garethky]

How to summarize this…

Investigation

This python notebook is a deep dive into Prusa’s use of Butterworth Bandpass filters in their endstop implementation. I’ve established that mkfilter, Matlab and SciPy all work the same and can generate the same filter coefficient as seen in their code. They also filter data in the same way, the outputs from the mkfilter code and SciPy’s lfilter are identical.

Prusa couldn’t come up with 1 filter design to run all their printers. They have 3, one for the MK4, from for the XL when probing in Z and one for the XL when probing in X & Y.

Their filters are more narrow than I expected:

Printer	Low Frequency Cutoff	High Frequency Cutoff
MK4	0.8Hz	12Hz
XL	11.2Hz	17.6Hz

As a plot of frequency response:

I understand the reason for the low frequency cutoff, that’s about rejecting drift. But I don’t understand the reasoning behind such a narrow passband.

When applied to a tap I see 3 potential problems with these filters:

1. The lag. The filtered lines are shifted to the right from the real force plot. This means the real force is significantly higher than the reported filtered force when they trigger. So 40g isn’t really 40g, its more like 200g.

2. The XL graph doesn’t even break the 40g trigger threshold (red horizontal lines).

3. Lastly all that wiggling at the start of the graph is the filter settling. That’s a long time to wait.

I cant pick a single filter design for all printers, I want faster triggering and I want that settling time gone.

So after much experimentation I have this:

1. The “Proposed Filter” is a highpass filter with the same low frequency cutoff as the corresponding Prusa filter, combined with a notch filter at 60Hz. Using a highpass filter keeps all the high frequency signal components and this results in the filters tracking very closely with the blue/true force line.

2. The filters are immediately flat. Using sosfilt_zi, the filter steady state initialization vector is calculated and then multiplied it by the first sample. This results in no settling time.

All of these filters perform the same at rejecting low frequency drift in the test data:

That is “continuous taring”. The proposed filters are just a little bit more noisy than Prusa’s, but still less than the raw force line. I very much like that tradeoff.

Implementation

I settled on the Second Order Sections (sos) filter format. Its not quite as memory compact as the b, a format but the code to execute it is very simple (right out of the SciPy source code) and its more resistant to 32 bit floating point errors:

def filter_step(self, x: float) -> float:
    x_cur: float = x
    for section in range(len(self.sos)):
        x_new = self.sos[section][0] * x_cur + self.zi[section][0]
        self.zi[section][0] = (self.sos[section][1]) * x_cur
                              - (self.sos[section][4]) * x_new
                              + (self.zi[section][1])
        self.zi[section][1] = (self.sos[section][2]) * x_cur
                              - (self.sos[section][5]) * x_new
        x_cur = x_new
    return x_cur

Translating that to pure c was not too difficult.

I know that this is a lot of math and theory. I’ve committed a filter_workbench notebook to the branch that would be part of a PR. Hopefully this helps printer designers work out a good filter for their printer.

I’d like to offer regular users an “easy mode” config to express filters like this:

[load_cell_probe]
    endstop_filter: {
              'highpass': 0.8,
              'notch': 60.0
       }

This would allow a lot of flexibility and readability. But using this config would require SciPy be installed to do all the math to convert this into a filter and initialization vector.

SciPy can be installed just like NumPy:

# for python 2.7
sudo apt-get install python-scipy

# for python 3
sudo apt-get install python3-scipy

For manufacturers and printer designers that need more control, I could take in the raw filter design:

[load_cell_probe]
    raw_endstop_filter: {
        'sos': [[ 0.92487988,  1.84975976,  0.92487988,  1.,  1.86554535,  0.87403728]
                [ 1.,         -2.,          1.,          1., -1.97845941,  0.97868912]
                [ 0.60876834, -0.50578302,  0.60876834,  1., -0.50578302,  0.21753668]],
        'zi': [[ 0.06440686  0.06020638]
               [-0.98928674  0.98928674]
               [ 0.          0.        ]]
       }

This is very flexible and avoids a dependency on SciPy, but its a potential footgun because who knows what that block of numbers does? (You can find out via SciPy sosfreqz, but who will do that?) Users blindly copy config files they don’t understand and that scares me.

I’d like to hear what you all think of this idea.