Consider the design of a Savitzky-Golay filter of window length 101 and (high) polynomial order 20. Using scipy version 1.10.1, the filter coefficients can be obtained in python as:
from scipy import signal h = signal.savgol_coeffs(window_length=101, polyorder=20)
For reference, in Matlab r2022b
order=20; framelen=101; b = sgolay(order,framelen); h = b((framelen+1)/2,:);
which looks more reasonable. For lower order polynomials, the results coincide much more closely.
It looks like scipy's implementation solves a least squares problem where the design matrix is computed as:
order = np.arange(polyorder + 1).reshape(-1, 1) A = x ** order
leading to a condition number > 1e29!
Is there some other way to get these coefficients in python?