Testing the FDM algorithm on experimental data
Introduction
The Filter Diagonalisation Method (FDM) can be used to find oscillator frequencies, amplitudes and phases in a recorded signal generated as a sum of the individual oscillators. One implementation of this algorithm is provided as the harminv program in Ubuntu, written by S. G. Johnson.
A similar operation can be done using the well-known Fourier transform of the signal, however, the FDM can give good results even when 1) the time record is short, 2) many oscillators are used in noisy signal, 3) weak signals are superposed on top of much stronger signals of similar frequency.
Testing on analytic data with near peaks
On the following figure, the above plot shows the rather short record of the time evolution of the input signal. The bottom plot shows: its spectrum (Fourier transform, blue line), the original oscillator frequencies (red vertical lines), peaks found in the Fourier spectrum (violet triangles) and finally, the frequencies perfectly found by running Harminv.
Obviously, Harminv gave much more accurate results than FFT, exactly finding the center frequencies of the oscillators provided, althought the central oscillator was much stronger than two on its sides.
It ought to be noted, however, that this was rather optimum setup for FDM in order to show its capabilities. Further changing the parameters (record length, record sampling, oscillator frequencies etc.) can randomly result in FDM "forgetting" or "merging" some of the oscillators.
Sources for this test
The Python script that generates the data and runs harminv can be downloaded here:python ./testing_harminv.py . Computation time is negligible.
Testing on experimental spectra of water vapour
In the following, we use a similar approach to a short time-domain terahertz spectroscopy (TDTS) record of terahertz spectra of ca. 70 cm of air in the laboratory, as is shown by the thick green line in the first subplot of the figure.
above: The time-domain record of the input data (thick green line); the exponentially decaying part used for the computation,
bottom: Fourier transform of the signal (orange line), oscillator frequencies found by harminv (blue stars) and reconstructed signal build by these oscillators (blue curve). Individual Lorentzian oscillators plot as thin blue lines.
The ringing after the pump pulse corresponds to the water vapor absorption spectra. We cut the first transient peak so that the processed signal is not dominated by broadband pumping.
Both Fourier transform and FDM were used to separate the resonant peaks. While Fourier transform gives only a very rough spectrum, FDM detected these frequencies of oscillators [THz]:
0.62562 0.75036 1.073696 1.409286 (...)
Experimental resonances from Exter M., et al.,: "Terahertz time-domain spectroscopy of water vapor":
0.558 0.753 0.989 1.099 1.115 1.165 1.209 1.412
We can see that FDM detected peaks with good accuracy even from very short experimental data. Only the very weak oscillators around 1 THz were not detected in the signal.
Sources for this test
The input file and the script. Download both files and run python airspec.py.