5
$\begingroup$

Given the transfer function $H(z) = \frac{1}{4}(1+z^{-1}+z^{-2}+z^{-3})$ and the input signal $x[n] = \cos(\frac{\pi}{50}n)+\cos(\frac{4\pi}{5}n)$, how can I find the output signal $y[n]$ in mathematica?


What I tried

H[z_] := 1/4 (1 + 1/z + 1/z^2 + 1/z^3)
x[n_] := Cos[π/50 n] + Cos[(4 π)/5 n]


zt = ZTransform[x[n], n, z, Assumptions -> {n ∈ Integers, n >= 0, z >= 0}];

InverseZTransform[Abs[zt * H[z]], z, n, 
                  Assumptions -> {n ∈ Integers, n >= 0, z >= 0}]

However, I get a bunch of error messages from this.

Assumptions::mepreal :  "In attempting to decide whether a solution <<...>.
is real, $MaxExtraPrecision  50.` was encountered. The solution was assumed real. 
Increasing the value of $MaxExtraPrecision may help resolve the uncertainty.>>

I'm really not sure how to filter this signal in Mathematica.

$\endgroup$
1
  • $\begingroup$ I don't understand what filtering is in this context, but if you drop the Abs, the InverseZTransform gives a reasonable result. $\endgroup$ Oct 17, 2015 at 15:47

2 Answers 2

10
$\begingroup$

Mathematica considers the unilateral (or one-sided) z-transform of a sequence $x[n]$, that is, $X(z)=\sum_{n=0}^\infty x[n] z^{-n}$, which implies a causal sequence $x[n]$ (defined only for non-negative $n$). In signal processing applications, where non-causal sequences are often considered, it is more convenient to use the two-sided version of the z-transform, i.e., $X(z)=\sum_{n=-\infty}^\infty x[n] z^{-n}$.

However, even if Mathematica had an implementation of the two-sided (inverse) z-transform, in the majority of cases (w.r.t. signal processing applications at least) there is no need to employ the formal definition of the transform as "standard" sequences are typically involved which can be simply handled by inspection.

In your case, it is easy to see that the inverse z-transform of $H(z)$ is $h[n]=\frac{1}{4}(\delta[n]+\delta[n-1]+\delta[n-2]+\delta[n-3])$. It follows at once that the output of the filter will be $y[n]=\frac{1}{4}\sum_{k=0}^3 x[n-k]$.

Plotting the input and output sequences is trivial. (The discrete-time signals are plotted as piecewise linear for visualization purposes.)

x[n_] := Cos[\[Pi]/50 n] + Cos[(4 \[Pi])/5 n];
y[n_] := Sum[x[n - k], {k, 0, 3}]/4.;

Grid[{{   
 DiscretePlot[x[n], {n, -50, 50}, Frame -> True, Axes -> False, 
 FrameLabel -> {"n", "x[n]"}, ImageSize -> 250, Joined -> True, 
 Filling -> None],
 DiscretePlot[y[n], {n, -50, 50}, Frame -> True, Axes -> False, 
 FrameLabel -> {"n", "y[n]"}, ImageSize -> 250, Joined -> True, 
 Filling -> None]    
}}]

enter image description here

Note that the output of the filter has suppressed the higher frequency cosine term of the input signal

$\endgroup$
5
$\begingroup$

You could use the built-in functions OutputResponse or RecurrenceFilter.

tfm = TransferFunctionModel[1/4 (1 + z^-1 + z^-2 + z^-3), z^-1, 
  SamplingPeriod -> 1]

OutputResponse[N@tfm, 
Cos[(π n)/50] + Cos[(4 π n)/5], {n, 0, 200}] // ListLinePlot

RecurrenceFilter[N@tfm, Table[Cos[(π n)/50] + Cos[(4 π n)/5],
{n, 0, 200}]] // ListLinePlot

enter image description here

$\endgroup$
6
  • $\begingroup$ I like your response to the question but why have you used z^-1 in the transfer function rather than z? $\endgroup$ Jan 11, 2019 at 6:46
  • $\begingroup$ Simply because that is the transfer function in the question. $\endgroup$ Jan 11, 2019 at 14:23
  • $\begingroup$ I don’t think I was clear - my fault. You wrote TransferFunctionModel[1/4(1+z^-1+z^-2+z^-3), z^-1, SamplingPeriod->1] I was expecting TransferFunctionModel[1/4(1+z^-1+z^-2+z^-3), z, SamplingPeriod->1]. The difference is the 2nd argument where you have z^-1 and I have z. I noticed this because I tried a BodePlot of your transfer function and it comes up with a bunch of errors, but if I change to z then BodePlot works but OutputResponse doesn’t. $\endgroup$ Jan 13, 2019 at 2:24
  • $\begingroup$ With $z^{-1}$ the transfer function is considered as a polynomial in $z^{-1}$, whereas with $z$ it is a polynomial in $z$, which has the effect of adding additional poles and zeros. In either case I don't get any error with OutputResponse but with $z^{-1}$ BodePlot returns an empty plot. What version are you using? $\endgroup$ Jan 14, 2019 at 0:11
  • $\begingroup$ Thanks for your response. I am using 11 on a Mac. For z: OutputResponse hangs (CPU usage maxed out - I have to abort), but BodePlot looks correct. For 1/z: OutputResponse looks good, but BodePlot has a bunch of errors and gives empty plots. $\endgroup$ Jan 14, 2019 at 1:08

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

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