Stabilizing the inverse transform of a system

Question

Currently I am trying to understand an unknown mechatronic system. I used white noise as an input signal and measured the output of the system. Then I estimated a second order SISO transfer function (Output/Input) in the frequency domain using Matlab.

sys=Out/In
h = idfrd(sys, freq*2*pi, Ts);
np = 2;
nz = 2;
Gs = tfest(h, np, nz)
opts = c2dOptions('Method', 'tustin', 'PrewarpFrequency', 2*pi*120);
Gz = c2d(Gs, 1/4000, opts);

$$ Gs=\frac{3.765\cdot s^2 - 528.7 s - 5.113 \cdot 10^4}{s^2+183 \cdot s + 4.117 \cdot 10^5} $$

A bode Plot of the derived second order system $Gs$ (blue) and the measured data (orange) shows good agreement in amplitude and phase above $60\;Hz$.

Further, a pole-zero plot of $Gz$ shows that all poles are stable (inside unit circle of z-plane), but that there is one zero outside of the unit circle of the z-plane.

Now, I want to counteract the behaviour of the mechatronic system in a way that I achieve a unity output. So by modifying the input to the system I want to basically get rid of the system characteristics and achieve a pure unit transfer behaviour.

IGz=1/Gz; % so that IGz * Gz = 1

Mathematically, I need to multiply the input signal with the inverse transform of the system. Here a bode plot of the system $Gs$ (blue) and its inverse transform $IGz$ (orange) basically shows what I want to do.

But because of the zero outside the unit circle, the inverse transfer function now has a pole outside the unit circle and thus is unstable. Therefore I wonder how I can counteract the behaviour of the system while not getting an unstable system?

---

Attempt

So far, I have tried to replace the unstable pole by reflecting the pole inside the unit circle as advised here: $(z-p)$ by $(z-1/p)$ (as the unstable pole is purely real I don't need to do the complex conjugate). This keeps the amplitude response, but changes the phase response (blue: unstable $IGz$, orange: stabilized $IGz$):

However, I would like to keep the phase response as close as possible to my inverse transfer function as well. But, as I don't care so much about what happens below $60\; Hz$ I wonder if there is an approximation of my inverse system which is good above $~60\;Hz$ in amplitude and in phase, but can deviate below $~60\;Hz$ and is stable. I wonder if there is an established way how one would derive such a system? Here is what I have tried so far.

pIG=pole(IGz);
zIG=zero(IGz);

z=tf('z',1/4000);
Test=((z-zIG(1))*(z-zIG(2)))/((z-1/pIG(1))*(z-pIG(2)));
Fac=evalfr(IGz,100)/evalfr(Test,100);
Test1=Fac*Test;
bode(IGz,Test1)

The result looks good at higher frequencies, however the phase could be better below $100\;Hz$ or so (blue: unstable $IGz$, orange: stabilized $IGz$ with phase kept roughly the same above $100\;Hz$):

Edit

I just looked at the impulse responses. The system $Gs$ of course is oscillating, but decaying.

However, the impulse response of the inverse system is, of course, unstable.

The impulse Response of the second stabilized system is still overshooting with high amplitude. I wonder if this works, since I, of course, am limited in output w.r.t. the controller.

With the input $u$ to the controller and the output $y$, I want to do the following $u \rightarrow IGz (controller) \rightarrow y \rightarrow Gz (mechatronic\; system) \rightarrow y_2$. During all this, the output $y$ of the controller can't get too big and it should obey $\frac{y_2}{u}=1$.

How can I best come up with an $IGz$ that achieves this?

Answer 1

Like you mentionned, you cannot cancel a right-half-plane zero (or a zero outside the unit circle) by placing a pole on it. A unstable pole in your compensator will make the command of your controller unbounded (i.e. it will reach infinity).

There are no ways to cancel a right-half-plane zero. It's sometimes possible to "remove" a right-half-plane zero by chansing a sensor position.

From Murray's book http://www.cds.caltech.edu/~murray/books/AM08/pdf/fbs-limits_18Aug2019.pdf

The existence of the right half-plane zero can be removed if we choose to measure the location of the vehicle by the position of the center of the rear wheels instead of the center of mass.

Edit : A right-half plane zero can also be caused by a time delay. Delays to computation can be reduced, but some other delays cannot be avoided.

Edit 2 : I will try to answer your other questions later.

Answer 2

As I understood, you're trying to implement a deconvolution by directly inverting the estimated (or known) wavelet of the system response. The process fails due to the wavelet being non-minimum phase thus uninvertible. You also say that you don't care about the particular band being uninvertible since you're interested in frequencies above 60Hz.

This is a well-known problem in applied signal processing used in geophysics. There are two possible solutions for this situation:

1. Find a minimum-phase equivalent of your system response wavelet, $mpG$, find a stable inverse of it, $ImpG$, and a Wiener filter $WGmpG$ converting $G$ to $mpG$ and it's inverse $IWGmpG$. Your desired input to the mechatronical system would then be $y = u \star IWgmpG \star ImpG$
2. Band-limited deconvolution with pre-whitening: invert only the required hi-pass version $hpG$ of your wavelet $G$ adding some white-noise obtaining a stable but less accurate inversion $IwhpG = \frac {1} {hpG + \epsilon \cdot ||hpG ||} $. You'd then have $y=u \star IwhpG$.

You can also combine band-limiting, pre-whitening and min-phasing approaches to achieve the best composite inversion.