4
\$\begingroup\$

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?

\$\endgroup\$
11
  • \$\begingroup\$ Your assumption of mapping the DFT directly to the CTFT implies that your continuous time signals have a sampling rate of 1Hz. Is this intentional and expected? \$\endgroup\$ – Ste Kulov Nov 6 '20 at 6:07
  • \$\begingroup\$ @SteKulov i am aware of that. I can address that later simply via the time/frequency scaling property of Laplace & Fourier transforms. So I am not too concerned about that aspect yet. My main concern is the flipping signs in the poles of the estimates, which I have now verified with artificial data. \$\endgroup\$ – Edgar Brown Nov 6 '20 at 7:06
  • \$\begingroup\$ "poles I get are both positive and negative". Are you taking phase into account ? Or are you matching only the *gain". A pole at say, \$-1+2i\$ has the same effect as a pole at \$+1+2i\$ on the magnitude response. The contribution to the phase is what allows to tell them apart. Same situation for discrete-time stuff, but inside and outside unit circle. It is not clear to me if you are not utilising the phase info. Can you add a simple calculation example which demonstrates the problem ? \$\endgroup\$ – AJN Nov 6 '20 at 13:01
  • \$\begingroup\$ Are you using the result of the fft as complex numbers or only the magnitude of the fft result ? \$\endgroup\$ – AJN Nov 6 '20 at 13:02
  • \$\begingroup\$ @EdgarBrown Please show a sample calculation for a simple system which replicates the problem. a.k.a a minimum non-working example. \$\endgroup\$ – AJN Nov 6 '20 at 13:05
3
\$\begingroup\$

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.

\$\endgroup\$
2
  • \$\begingroup\$ Yes, this goes along the lines of what I realized with these discussions. I am not constraining the problem enough. I am looking on how to do that, I think I could do it with some inequalities and Lyapunov exponents if I can figure out a way to express those constraints. \$\endgroup\$ – Edgar Brown Nov 8 '20 at 12:31
  • \$\begingroup\$ If you find a good way to do it, I think your innovation would be a pretty big deal. It would be a great technique for strongly band-limited systems. Good luck! \$\endgroup\$ – rpm2718 Nov 8 '20 at 12:45
1
\$\begingroup\$

DFT is not the same as a CFT which, in turn, means that you can't treat the DFT as being a slice on the Laplace domain. You are applying a DFT on an \$s\$ operator, which can't be done, unless you transform \$s\rightarrow z\$, and that comes in one of the several ways to do it: impulse invariance, step invariance, bilinear transform, other methods.

So the DFT and the CFT are similar, but not the same, because one deals with the continuous and infinite \$j\omega\$, the other with the discrete and finite \$j\Omega\$, and the Laplace and the Z-transform are similar, but one deals with an infinite plain defined in \$s=\alpha+j\omega\$, and the other with the finite and very well defined domain of \$z=\text{e}^{j\Omega}\$.

Which means that your problem can only be solved in two ways: either work on the continuous transfer function in the Laplace domain with CFT, or apply a Z-transform and work with DFT/FFT.

\$\endgroup\$
9
  • \$\begingroup\$ AFAIK this is not strictly true. The sampled (in time and frequency) transforms are directly related to their continuous counterparts. You can have the continuous transform of a sampled sequence and the sampled transform of a continuous sequence. Sampling in one domain is simply periodicity in the other. A DFT simply combines these two aspects. The way I see it, a z-transform should not even be involved in this. \$\endgroup\$ – Edgar Brown Nov 6 '20 at 7:13
  • 1
    \$\begingroup\$ @EdgarBrown You can't apply a continuous transform to a discrete sequence, or vice-versa. The whole point of the Z-transform is to map the infinite domain of \$s\$ onto the finite region of \$z\$. If you do not do that, all the infinity of frequencies outside will be aliased into \$z\$. The fact that they are related does not mean they are the same. One is periodic, the other is not. You cannot mix the two. \$\endgroup\$ – a concerned citizen Nov 6 '20 at 7:22
  • \$\begingroup\$ of course you can. As long as you know what you are doing. Yes, the infinite frequencies will be aliased in, you are right about that, but this is not a problem if you account for that and in this case the signal has been defined with a finite spectrum response. The z-transform can be seen as a convenient representation tool, to avoid accounting for all of that periodicity in the Laplace domain. But these should be directly equivalent. \$\endgroup\$ – Edgar Brown Nov 6 '20 at 8:01
  • \$\begingroup\$ @EdgarBrown The Z-transform is a mapping tool and the Laplace domain is not periodic. They are not similar simply because you think they should be. Just because you plot your transfer function within a finite spectrum, defined by yourself, it doesn't mean that its response is not infinite. And, to be clear, you can apply a DFT to a continuous time transfer function, but it doesn't mean that the results will be meaningful other than for non-discrete sequences, and vice-versa. But if you want to think otherwise, so be it. \$\endgroup\$ – a concerned citizen Nov 6 '20 at 8:11
  • 1
    \$\begingroup\$ If you insist on performing an FFT, then you have to consider a very high sampling rate compared to the highest bandwidth frequency, to avoid (considerable) leakage. However, you'd be dealing with a virtual impulse invariance method, which means Z-transform, volens-nolens. \$\endgroup\$ – a concerned citizen Nov 6 '20 at 16:26

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.