3 Analyzing Signals
One of the major tasks of signal processing is the determination of characteristics of a signal, in both the time and frequency domains. Signals and Systems provides a number of functions to assist in this process. Two general functions to report on a variety of signal characteristics are included (one for continuous and one for discrete signals). (Note that all capabilities of these omnibus functions can be accessed individually in separate functions throughout the packages.) Also, a variety of plotting functions can be used to generate standard graphic representations of signals.
3.1 Reporting Signal Characteristics
The signal analysis functions ASPAnalyze and DSPAnalyze create reports of various signal characteristics. For ASPAnalyze, the components include a plot with respect to the independent variable, the Laplace transform with its region of convergence, a polezero plot, and the frequency response with a magnitudephase plot. DSPAnalyze generates the Z transform in place of the Laplace transform, and uses the discretetime Fourier transform in place of a Fourier transform. You can choose which of these elements are to be displayed, and the order in which they are returned.
Analyzing signals.
The routines handle only one and twodimensional signals, primarily because visualization techniques for higherdimensional signals are not well defined. The analysis proceeds with the named independent variable(s), while any other parameters in the signal expression are given the value while graphing. Explicit alternate values to be used for graphing can be specified by the ParameterValues option. This option takes a list of transformation rules to be applied to the signal before it is plotted.
First, verify that the signal processing functions are available.
In[1]:=Needs["SignalProcessing`"]
Signals are analyzed with a single command. Note that a value is specified for the parameter a during graphing.
In[2]:=ASPAnalyze[ t Exp[a t] Cos[3 Pi t/16] UnitStep[t], {t, 0, 3 Pi}, ParameterValues > {a > 1} ]
Out[2]=
Options specific to the analysis function.
Important options include AnalysisReport and AnalysisOutput. The AnalysisReport option accepts a list of tags that determine which components of the report will be generated, though the order in which they are generated cannot be specified, due to dependencies among the computations. These include SignalPlot, Transform, PoleZeroPlot, FrequencyResponse, and MagnitudePhasePlot. Some of these may not be generated when specified if the computation cannot be performed (e.g., when the transform cannot be found). The AnalysisOutput option determines which objects will be returned. It accepts the same tags as AnalysisReport. Null objects will be returned for items that could not be computed.
First, a straightforward analysis of a discretetime pulse is generated.
In[3]:=DSPAnalyze[DiscretePulse[10, n], {n, 5, 15}]
Out[3]=
The transform routines will attempt to use the Global` variables z and w in the returned Z and Laplace transforms (z1, z2, w1, and w2 in twodimensional transforms). If the user has already assigned values to these symbols, unique symbols of the form z$n, w$n, etc., will be used instead.
Set a value for w (a variable used by the analysis routines).
In[4]:=w = 3
Out[4]=
These options cause only the polezero plot to be displayed, and the discretetime Fourier transform to be returned with the plot. An option is passed to the plotting function to suppress the printing of the list of poles and zeros.
In[5]:=DSPAnalyze[DiscretePulse[10, n], {n, 5, 15}, AnalysisReport > {PoleZeroPlot}, AnalysisOutput > {FrequencyResponse, PoleZeroPlot}, PoleZeroPlotOptions > {Justification > False} ]
Out[5]=
Clear the value that was set before the last computation.
In[6]:=Clear[w]
This sets up a signal by first determining a passband in the frequency domain, and then transforming to the time domain.
In[7]:=(lowpass2Dtile = ContinuousPulse[Pi, w1 + Pi/2]* ContinuousPulse[Pi, w2 + Pi/2]; lowpass2D = InverseDiscreteTimeFourierTransform[ lowpass2Dtile, {w1, w2}, {n1, n2} ] )
Out[7]=
Here is an analysis of a downsampled version of the twodimensional signal.
In[8]:=DSPAnalyze[ Downsample[{{2,1},{1,3}}, {n1, n2}][lowpass2D], {n1, 5, 5}, {n2, 5, 5}]
Out[8]=
3.2 Plotting Signals
In signal processing, it is useful to get a feel for signal behavior by plotting the signal. A variety of visualization techniques can be found in Signals and Systems. Extended versions of a number of standard Mathematica graphics functions have been created by adding the prefix Signal to the function names and are used for basic plots of time (or spatial) response. Magnitudephase plots (and the variant called the Bode plot) are used for frequency response analysis, while polezero, rootlocus, and Nyquist plots allow various signal characteristics, such as stability, to be examined.
Signal plotting functions modeled after standard Mathematica functions.
The common functions for displaying the time (or spatial) response of a signal are enhanced versions of traditional Mathematica plotting functions. The functions for continuous signals add the ability to handle complex values and impulse functions, while those for discrete signals generate a "fencepost" plot (also known as "lollipop" or "stem" plots) in one dimension, and a density plot in two dimensions.
A complexvalued signal with a singularity function can be displayed by SignalPlot.
In[9]:=SignalPlot[ Exp[I t] + DiracDelta[t  2], {t, 0, 2 Pi} ]
Out[9]=
Here is a plot of a signal defined over a twodimensional support.
In[10]:=SignalPlot3D[ DirichletSinc[6, Sqrt[t1^2 + t2^2]], {t1, 3, 3}, {t2, 3, 3}, PlotPoints > 30 ]
Out[10]=
For SignalPlot (and only SignalPlot) the PlotStyle option has been redefined somewhat. Because there are two curves (the real part and the imaginary part) associated with each function, it is expected that a style will be specified as a pair, the first part of which applies to the real curve and the second to the imaginary curve. If only a single style is given, it will apply to both curves.
Here are three functions displayed with different styles.
In[11]:=SignalPlot[ {x^2/2, Sin[x] + I Cos[x], UnitStep[x] + I/2 UnitStep[x]}, {x, 0, 2}, PlotStyle > {GrayLevel[0], {Hue[0], Hue[.6]}, {{Dashing[{0.15, 0.15}], GrayLevel[0]}, {Dashing[{0.05, 0.05}], GrayLevel[0.5]}} } ]
Out[11]=
Option specific to continuous twodimensional plotting functions.
Dirac delta functions are graphed in more than one way in the signals and systems literature. One technique is to let the height of the impulse represent the area under the delta function, while another is to make them of uniform height, since they have infinite value at a point. Both of these can be generated with the DiracDeltaScaling option to SignalPlot and MagnitudePhasePlot. When set to True, the heights are scaled to the areas under the deltas. When False, the impulses are drawn as "infinite" pulses, running from the horizontal axis to the top of the graph.
Here, the delta functions are drawn scaled to the areas under the deltas.
In[12]:=SignalPlot[ DiracDelta[t + 1]  2 DiracDelta[t  1] + 3 DiracDelta[t  3], {t, 2, 4}, DiracDeltaScaling > True ]
Out[12]=
In this graph, the deltas are drawn as uniform impulses.
In[13]:=SignalPlot[ DiracDelta[t + 1]  2 DiracDelta[t  1] + 3 DiracDelta[t  3], {t, 2, 4}, DiracDeltaScaling > False ]
Out[13]=
A plot of a discrete signal can also be generated. Note that the independent variable can only take integer values for a discrete signal.
In[14]:=DiscreteSignalPlot[ Sin[n/3] DiscretePulse[5, n  2], {n, 0, 7} ]
Out[14]=
A convenient representation for a discrete signal over a twodimensional support is a density plot. The usual DensityPlot options can be passed through DiscreteSignalPlot3D.
In[15]:=DiscreteSignalDensityPlot[ t1 t2 DiscretePulse[{5, 7}, {t1, t2}], {t1, 2, 6}, {t2, 2, 9}, Mesh > False ]
Out[15]=
Here is the same function plotted as a surface. Note that the samples are located at the polygon intersections, not in the center of each polygon as they are for the density plot.
In[16]:=DiscreteSignalPlot3D[ t1 t2 DiscretePulse[{5, 7}, {t1, t2}], {t1, 2, 6}, {t2, 2, 9}, Mesh > False ]
Out[16]=
These functions accept essentially the same options as their equivalent standard functions. The default behavior of some options is modified; for instance, axes labels are generated by default, as is a plot label. For all continuous signal plots, the default value of the MaxBend option has been set to , to allow the creation of smoother plots.
Functions for plotting the frequency response.
A common technique for visualizing a frequency response is to separate the magnitude and phase, and plot each. The implementation of MagnitudePhasePlot generates two graphics, and returns a list of graphics objects. BodePlot is essentially a call to the magnitudephase plot with axes scales set for logarithmic points. These functions assume that a frequency response has been computed, and is being passed as the argument to the plotting function.
Here is a plot of the frequency response of an eighthorder comb filter.
In[17]:=MagnitudePhasePlot[ (1  z^(8))/(1  z^(8)/8)/. z > Exp[I w], {w, 0, 2 Pi} ];
The magnitudephase plot is the customary tool for frequency analysis of signals. Here is the discretetime Fourier transform of a twodimensional signal.
In[18]:=MagnitudePhasePlot3D[ DiscreteTimeFourierTransform[ DiscretePulse[3, n1  1] * DiscretePulse[3, n2  1], {n1, n2}, {w1, w2} ], {w1, 2 Pi, 2 Pi}, {w2, 2 Pi, 2 Pi} ];
The Bode plot is equivalent to MagnitudePhasePlot with certain options set to different default values. Here is a Bode plot of the frequency response of a particular analog transfer function.
In[19]:=BodePlot[ 0.01 (0.01 s + s^2 + 1)/ (10 s^2 + 1.1 s^3 + 0.635 s^4 + s^5/16)/. s > I w, {w, 0.1, 10} ];
Here is a discrete transform of a pulse of width starting at the origin.
In[20]:=trans = DiscreteFourierTransform[ DiscretePulse[10, n], 20, n, w ]
Out[20]=
Here is a magnitudephase plot of the result. Note that the transform consisted of an infinite summation of shifted functions (the effect of the Periodic operator). Since the exponential and sinc functions have infinite domain, the plotting function wasn't able to compute a closed form for the behavior of Periodic, and instead only computed a single cycle. Hence, the resulting image ignores aliasing effects of the shifted copies of the exponential and sinc terms in the transform.
In[21]:=DiscreteMagnitudePhasePlot[ trans, {w, 20, 20} ]
Out[21]=
Options specific to MagnitudePhasePlot.
MagnitudePhasePlot and MagnitudePhasePlot3D accept options standard to Plot and Plot3D, respectively. They also accept certain other options dictating the scaling to be used in the plots, and specifying whether the plots are over continuous or discrete frequency domains.
The DomainScale and MagnitudeScale options can take the values Log or Linear, specifying logarithmic or linear scaling on the horizontal axis of both graphs or the vertical axis of the first graph, respectively. If logarithmic scaling is used, the domain in question must be positive. The MagnitudeScale option can also be set to None, which will prevent the magnitude plot from being displayed.
PhaseScale determines whether the phase angle is plotted in degrees or radians when given the value Degree or Radian, respectively. This option affects only the vertical axis of the second (phase) graph. The PhaseScale option can also be set to None, which will prevent the phase plot from being displayed.
The Domain option will cause the plot to be generated in the continuous frequency domain (if set to Continuous) or the discrete domain (if given as Discrete). The range will take on only integer values in the discrete domain.
Of these specific options, BodePlot accepts only PhaseRangeScale. MagnitudePhasePlot also accepts the DiracDeltaScaling option described earlier.
Here is the magnitude frequency response of a secondorder Butterworth filter plotted on a logarithmic scale. Note that as the frequency grows much greater than the cutoff frequency, the rolloff approaches decibels per octave, where is the order of the Butterworth filter.
In[22]:=MagnitudePhasePlot[ FourierTransform[ AnalogFilter[ GetRootList[s^2 + Sqrt[2] s + 1, s], {}, t ], t, w ], {w, 1, 10}, DomainScale > Log, MagnitudeScale > Log, PhaseScale > None, GridLines > {{2,4,8},{0, 12, 24, 36, 48}}, Frame > True, Axes > False ]
Out[22]=
Polezero plots in Signals and Systems.
The basic polezero plot is designed for plotting poles and zeros of rational transfer functions, along with appropriate regions of convergence. Often, these are taken from the results of ZTransform and LaplaceTransform. If so, most of the information necessary to generate the plot is taken from the transform object, and does not need to be specified again. Alternate forms of the syntax for PoleZeroPlot allow arbitrary functions to be used, although interpretation of the region of convergence in the Z or Laplace domain then needs to be specified by the Domain option. The plot range will usually be determined automatically (by the locations of the poles and zeros), but can be overridden by use of the PlotRange option.
Here is a polezero plot of a Laplace transform.
In[23]:=PoleZeroPlot[ LaplaceTransform[ t Exp[t] UnitStep[t], t, s ] ]
Out[23]=
This is a polezero plot of a Z transform.
In[24]:=PoleZeroPlot[ ZTransform[ 0.5^n DiscreteStep[n], n, z ] ]
Out[24]=
Here is a plot of the combfilter transfer function used earlier. Note that with the ZTransform default domain, the unit circle is drawn to clarify the plot.
In[25]:=PoleZeroPlot[ (1  z^(8))/(1  z^(8)/8), z ]
Out[25]=
If a region of convergence is specified, the appropriate area is shaded.
In[26]:=PoleZeroPlot[ z^2/(z^2 + 0.25), {z, .5, Infinity} ]
Out[26]=
Here is a polezero plot of a signal over a twodimensional support. This particular system is unstable; note that poles and zeros can be found outside of the region of convergence. The lower bound of the region of convergence is marked by the dashed line. The zeros and poles are drawn in separate graphs, in which the first two plots correspond to the first variable, with the second variable mapped by Exp[I w1], while the second two graphs reverse this order, mapping the first variable to Exp[I w2].
In[27]:=PoleZeroPlot[ ZTransform[ ((n1 + n2)! (4/5)^n1 * (2/7)^n2 DiscreteStep[n1, n2])/ (n1! n2!), {n1, n2}, {z1, z2} ] ]
Out[27]=
Options for PoleZeroPlot.
The Justification option causes the locations of the poles and zeros to be printed, while the PoleZeroPhasors option causes lines to be drawn from the poles and zeros to the specified point.
In[28]:=PoleZeroPlot[ (1  z^(8))/(1  z^(8)/8), z, Domain > ZTransform, Justification > All, PoleZeroPhasor > {.5, .5} ]
Out[28]=
The rootlocus plot.
The rootlocus plot is used to investigate the evolution of the positions of roots of some function as a parameter changes over a series of values. The RootLocusPlot function will attempt to symbolically factor the function so that separate root curves can be plotted parametrically. The initial locations of the poles and zeros are displayed in red, while the final locations are in black.
Here is a rootlocus plot.
In[29]:=RootLocusPlot[ 1/(s ((s + 4)^2 + 16) + k), s, {k, 0, 70} ]
Out[29]=
Standard graphics options can be used with RootLocusPlot, and will be applied to all of the images.
The Nyquist plot.
The Nyquist plot is useful for investigating the frequency response of a transfer function in the complex plane.
Here is the Nyquist plot of a transfer function.
In[30]:=NyquistPlot[1/(I w + 1)^2, {w, 10, 10}, Frame > True, Axes > False ]
Out[30]=
NyquistPlot will accept standard Graphics options, as well as the option PlotPoints. The PlotPoints option is set to a default value of , since curves plotted by this method tend to have considerable variability.
