I am trying to identify a system by means of its differential equation (i.e., Lapace representation). I put together a rather straightforward regression algorithm (similar to Proni's method for ARMA) under the assumption that the FFT of the system response is equivalent to the Laplace transform evaluated at \$ -j\omega \$ (with \$ \omega \$ restricted to \$ (-\pi , \pi ) \$, (but this doesn't seem to change anything).
It converges quite consistently but the poles I get are both positive and negative (i.e., stable and unstable if I interpret them as I intended in the s-domain) and both inside and outside the unit circle (again unstable and stable if I interpret them as a Z-transform instead).
I am no newby to DSP, but it is clear I have some misconception somewhere. Most likely on my interpretation of the FFT as simply a sampling of the Laplace transform on the imaginary axis or in the proper scaling of this axis. I have checked the equations back and forth multiple times, and I can't see what I am doing wrong.
The gist of it is simply to solve the equation:
$$
y(s)*(1 +a_1s + a_2s^2 + a_3s^3) = x(s)*(b_0 + b_1 s + b_2 s^2)
$$
For the values of \$ a_i,b_i \$ that minimize the representation error, setting \$ s=j\omega\$ with \$ \omega \text{ in } (-\pi,\pi) \$ and making \$x(j\omega)\$ the FFT of a known input and \$ y(j\omega)\$ the FFT of the measured output (trimming the edge frequency bins as these are noisy).
I know that the FFT of a finite sampled sequence is simply sampling the Fourier transform and adding the aliased components into it (which in this case can be assumed to be zero, as the data stream is digitally filtered and subsampled).
I want to implement it in the continuous s-domain, as it is modeling an analog system and that would make the results easier to understand. So, before I throw all of this away and model it in the discrete Z-domain as an ARMA, can someone point to where my problem lies?
Best Answer
I believe that your optimization calculation has found a mathematically-correct, but non-causal solution.
To understand this, we have to dredge up some of the theory that underpins the Laplace transform. First, as you are aware, when you use the FFT, your time-domain signal is one that repeats periodically from negative infinity to positive infinity. So the Laplace Transform that we will be associating with the FFT is the bilateral transform, not the unilateral transform.
In the theory of Laplace Transforms, the Region of Convergence (ROC) matters. Given a particular expression for the Laplace Transform, the inversion is not unique until you specify the ROC. Recall that a causal system must have an ROC that includes the the portion of the s-plane to the right of all poles out to Re(s) of \$+\infty\$. A stable system must include the imaginary axis in its ROC. And, of course, as we live and breathe, a stable and causal system must have and ROC that includes both the imaginary axis and the right portion of the s-plane, and so all poles must be in the left-hand plane.
Suppose for a moment that you have a Laplace Transform with poles in both the right half plane and left half plane, but not along the imaginary axis. You can choose two ROCs to get two different inversions. If you choose the ROC to include Re(s) of \$+\infty\$, you wind up with a causal, but not stable, solution. If you choose an ROC that includes the imaginary axis, it cannot include the portion out to Re(s) of \$+\infty\$ because of the poles in the right half plane, and so you wind up with a stable, but not causal, solution.
Now, back to your calculation. You are selecting a solution such that the Fourier Transform exists, and so your ROC includes the imaginary axis. However, nothing in the calculation forces the poles to the left-hand plane....(your optimizer does not share your interest in a causal solution!), and so it is perfectly happy to put some poles in the right half plane, leading to a non-causal solution.
I don't see a way to force your poles all into the left half plane in this calculation, but I sure would like to hear about it if you think of one.