Skip to content

Implementation of convert_spectrum_to_IIR() function of PropagationParameter::CFrequencySpectrum class

Issue

FrequencySpectrum class is missing the conversion function from spectrum to IIR. This procedure is more complicated than the other way round, as a curve has to be fitted to the spectrum. There are different algorithms one can use but my research showed that there is no "simple" one.

Consequence

Right now, Propgation Modules can set the IIR inside the FrequencySpectrum class but calling getFrequencySpectrum() throws a NotImplemented exception. Therefore, only magnitude and phase spectras are allowed/work atm.

Questions to be anserwed

  • Which algorithm to use to convert from spectrum to IIR?
  • For which input is the algorithm stable?
  • There is a tradeoff between number of coefficients in the IIR and accuracy of the fitted curve (if I understood correctly.)
    • How many coefficients do we want/need? Maybe make it variable?
    • We only use Biquads. Therefore, our IIR is either only a biquad or needs to be compoused out of multiple Biquads. The second option must be taken into account when choosing the algorithm for the conversion.

My research

  • Matlab has the method we need called invfreqz, however it's implementation is rather old and I've read somewhere that it may not always be stable.
  • This reddit post seems to have the excat same issue. I am unsure how helpful the answers are though.
    • One answered with the 'Textbook method' which seems intersting. The algorithm described is the one used in Matlab as well, just without iterations.
    • Steiglitz-McBride algorithm mentioned
  • Stackexchange post regarding the accuracy of matlabs invfreqz function in regards to audio filters. May give some insight into accuracy and stability? It also mentions a few algorithms:
    • Berchins FDLS: While googling I found this github repo which implements the FDLS technique in python. Maybe it is helpful.
    • Fixed-pole parallel filter design: Haven't looked into it.
    • Deep learning: We probably do not want this method.
  • ITADSP has some kind of conversion function, however a conversion from FIR to IIR, see here. Additionally, The Yulewalk function seems to be not fully implemented or contains a bug. The Burg algorithm is the only one actually setting Filter Coefficients. One may have a look at the implementation but it's not documented.