This example shows how to improve the frequency-domain accuracy of a system with a time delay that is a fractional multiple of the sample time.

For systems with time delays that are not integer multiples of the sample time, the `Tustin`

and `Matched`

methods by default round the time delays to the nearest multiple of the sample time. To improve the accuracy of these methods for such systems, `c2d`

can optionally approximate the fractional portion of the time delay by a discrete-time all-pass filter (a Thiran filter). In this example, discretize the system both without and with an approximation of the fractional portion of the delay and compare the results.

Create a continuous-time transfer function with a transfer delay of 2.5 s.

```
G = tf(1,[1,0.2,4],'ioDelay',2.5);
```

Discretize `G`

using a sample time of 1 s. `G`

has a sharp resonance at 2 rad/s. At a sample time of 1 s, that peak is close to the Nyquist frequency. For a frequency-domain match that preserves dynamics near the peak, use the Tustin method with prewarp frequency 2 rad/s.

discopts = c2dOptions('Method','tustin','PrewarpFrequency',2); Gt = c2d(G,1,discopts)

Warning: Rounding delays to the nearest multiple of the sampling period. For more accuracy in the time domain, use the ZOH or FOH methods. For more accuracy in the frequency domain, use Thiran filters to approximate the fractional delays (type "help c2dOptions" for more details). Gt = 0.1693 z^2 + 0.3386 z + 0.1693 z^(-3) * ------------------------------ z^2 + 0.7961 z + 0.913 Sample time: 1 seconds Discrete-time transfer function.

The software warns you that it rounds the fractional time delay to the nearest multiple of the sample time. In this example, the time delay of 2.5 times the sample time (2.5 s) converts to an additional factor of `z^(-3)`

in `Gt`

.

Compare `Gt`

to the continuous-time system `G`

.

plotopts = bodeoptions; plotopts.Ylim = {[-100,20],[-1080,0]}; bodeplot(G,Gt,plotopts); legend('G','Gt')

There is a phase lag between the discretized system `Gt`

and the continuous-time system `G`

, which grows as the frequency approaches the Nyquist frequency. This phase lag is largely due to the rounding of the fractional time delay. In this example, the fractional time delay is half the sample time.

Discretize `G`

again using a third-order discrete-time all-pass filter (Thiran filter) to approximate the half-period portion of the delay.

discopts.FractDelayApproxOrder = 3; Gtf = c2d(G,1,discopts);

The `FractDelayApproxOrder`

option specifies the order of the Thiran filter that approximates the fractional portion of the delay. The other options in `discopts`

are unchanged. Thus `Gtf`

is a Tustin discretization of `G`

with prewarp at 2 rad/s.

Compare `Gtf`

to `G`

and `Gt`

.

plotopts.PhaseMatching = 'on'; bodeplot(G,Gt,Gtf,plotopts); legend('G','Gt','Gtf','Location','SouthWest')

The magnitudes of `Gt`

and `Gtf`

are identical. However, the phase of `Gtf`

provides a better match to the phase of the continuous-time system through the resonance. As the frequency approaches the Nyquist frequency, this phase match deteriorates. A higher-order approximation of the fractional delay would improve the phase matching closer to the Nyquist frequencies. However, each additional order of approximation adds an additional order (or state) to the discretized system.

If your application requires accurate frequency-matching near the Nyquist frequency, use `c2dOptions`

to make `c2d`

approximate the fractional portion of the time delay as a Thiran filter.

`c2d`

| `c2dOptions`

| `thiran`

Was this topic helpful?