Signal Processing Packages and Notebooks
for the Mathematica Environment
Brian Evans and James McClellan
Digital Signal Processing Group
School of Electrical Engineering
Georgia Institute of Technology
Atlanta, GA 30332-0250
EMAIL: evans@eedsp.gatech.edu
Since the version of the signal processing packages released by
The Mathematica Journal, I have made changes to just about every file.
The major changes implemented by the new version of the signal
processing packages are discussed in the SignalProcessingExamples
notebook, but they are summarized here:
-------------------------- Version 2.0 ----------------------------------
(1) added object SignalPlot that plots the real and imaginary components
of a 1-D function simultaneously on the same plot, as well as
Dirac delta functions. It also plots 2-D signals. (The signal
processing packages already plot 1-D discrete-time functions in
lollipop style via DiscreteGraphics.)
(2) the option Dialogue -> All now works for the forward bilateral Laplace
transform (which means that all transform rule bases can completely
justify their answers)
(3) can plot upsampled sequences again
(4) unified piecewise convolution (either function can be in list
or expression form)
(5) added many transform pairs to the continuous-time Fourier transform
rule bases CTFTransform and InvCTFTransform
(6) the Fourier transform rule bases now handle the convolution operator:
CTFTransform[ Convolve[t][CPulse[1, t + 1/2], CPulse[1, t + 1/2]], t, w ]
(7) RootLocus plots a root locus for one varying parameter
(8) CFIR and CIIR structures have an optional third argument: Roots -> roots.
Change is reflected in AnalogFilters notebook.
(9) Computational objects for IIR structures, called IIRFunction, now exist.
Also, FIR structures can be rewritten as a formula. Their values
can now be plotted. Values are calculated as they are needed, and
IIRFunction remembers (caches) previous calculations.
(10) SequencePlot now provides a standard way to plot discrete-time expressions
(11) Mathematica 1.2 has a problem performing partial fractions decomposition
on real-valued polynomials even if you do all the factoring yourself.
So, we have added a new option to the InvLaPlace object called Apart.
If set to All, then we kludge around the Apart deficiency so that
a partial fractions decomposition is formed. Otherwise, we inverse
transform the function as an continuous-time IIR filter whose input
is the inverse transform of the numerator. This is made possible
by a new function called MyApart (see below).
(12) The z-transform of IIR structures returns a valid region of convergence.
(13) RealQ[z] returns True if and only if the head (tag) of z is Real.
We have added RealValuedQ[z] which returns True if the imaginary
component of z is 0 or 0.0 (this behaves the way RealQ used to).
-------------------------- Version 2.1 ----------------------------------
(14) ZSolve, a difference equation solver, is more robust.
ZSolve can also justify its answers in the same way as the transforms do.
(15) DTFT of Step[n] is now correct (was missing a pi term).
(16) LSolve can also justify its answers in the same way as the transforms do.
(17) Convolution is now represented by more than one operator.
Linear convolution in the discrete domain is Convolve[n] and
linear convolution in the continuous domain is CConvolve[t].
-------------------------- Version 2.2 ----------------------------------
(18) Discrete-time convolution is now implemented.
(19) New notebook called "SignalProcessingIntroduction" serves as an
introduction to Mathematica as well as DSP.
-------------------------- Version 2.21 ---------------------------------
(20) SequencePlot for 2-D signals now samples functions at integer
indices.
-------------------------- Version 2.22 ---------------------------------
(21) Stable will resolve stability of non-separable signals more often
(e.g., Stable will return False when applied to the z-transform of
(1/2)^n1 (4/5)^n2 Multinomial[n1,n2] Step[n1,n2]
instead of an unresolved expression.)
-------------------------- Version 2.23 ---------------------------------
(22) Cleaned up descriptions and implementation of parameterized operators
like Shift, Upsample, and Z.
-------------------------- Version 2.3 ----------------------------------
(23) The TransformLookup option for the transform rule bases allows one
to specify additional transform pairs. This is useful when a trans-
form rule base does not contain a pair that you need. It is also
useful for transforming abstract functions (i.e., saying that x[n]
becomes X[z]). For example,
In:= ZTransform[a^n x[n], n, z, TransformLookup -> { x[n] :> X[z] }]
Out= ZTransData[ X[a z], Rminus[0], Rplus[Infinity], ZVariables[z]]
The N-dimensional case is not as straightforward because each
transform rule base applies a one-dimensional rule base N times
in the order that the N time variables are given. Therefore,
users must specify intermediate transform pairs:
In:= LaPlace[ x[t1,t2], {t1,t2}, {s1,s2},
TransformLookup -> { x[t1,t2] :> X1[s1,t2],
X1[s1,t2] :> X[s1,s2 } ]
Out= LTransData[ X[s1, s2], Rminus[{0, 0}],
Rplus[{Infinity, Infinity}], LVariables[{s1, s2}] ]
Here, X1[s1, t2] denotes the Laplace transform of x[t1,t2] with
respect to the first variable t1. It turns out that you could use
X instead of X1 and the transform will still work out correctly.
-------------------------- Version 2.31 ---------------------------------
(24) Debugged the TeX forms of the signals (functions) and systems
(operators) introduced by the signal processing packages, esp.
their multidimensional forms.
(25) Added AllSubsets function which returns all of the possible subsets
of a set. For example, AllSubsets[ {v1,v2,v3} ] returns
{ v1, v2, v3, {v1, v2}, {v1, v3}, {v2, v3}, {v1, v2, v3} }.
Similarly, AllSubsets[ {v1, v2, v3}, Plus] returns
{ v1, v2, v3, v1 + v2, v1 + v3, v2 + v3, v1 + v2 + v3 }.
------------------------- Version 2.32 --------------------------------
(26) The Laplace transform of (t - 1) CStep[t] was incorrect. It returned
the same transform as (t - 1) CStep[t - 1] did. One rule was patched
in "LaPlace.m".
(27) The inverse Laplace transform cannot handle complex-valued roots
of a rational polynomial because of limitations of 1.2 Mathematica's
Apart function. We have encoded a new function MyApart which
works around the problem but it is very, very slow. Now, inverse
transforms of functions like "s / ( s^2 + 4 )^2" are possible.
------------------------- Version 2.33 --------------------------------
(28) The inverse Laplace transform was not always correctly inverting
terms that corresponded to integration in the time domain. Patched
one rule and added another to the post-processing rules.
(29) The primitive MyApart did not normalize the coefficient of the
highest power term in the denominator before rooting it.
------------------------- Version 2.34 --------------------------------
(30) The inverse Laplace transform was still not always correctly inverting
terms that corresponded to integration in the time domain. Patched
one rule and added another to the post-processing rules.
------------------------- Version 2.4 --------------------------------
(31) Integrating Delta functions over a finite interval did not work
when the coefficient on the variable was one. Fixed.
(36) The routines supporting the rules for multirate multidimensional
signal processing (i.e., lattice and integer matrix theory) have
been gathered in "Multirate.m" in the Support directory, including
"DistinctCosetVectors", "SmithNormalForm", and "SmithReducedForm".
(33) Inverse z-transforms of affine log functions are fixed.
(34) WORKS UNDER ALL VERSIONS OF MATHEMATICA NOW!!!!!
In making sure that the signal processing packages work for BOTH
Mathematica 1.2 and 2.0, the following changes have taken place:
(a) aliases have been discontinued so the former aliased objects
AliasedSinc, AliasSinc, ASinc, InverseLaPlaceTransform,
LaPlaceTransform, and MagnitudePhasePlot are now set equal to
Dirichlet, Dirichlet, Dirichlet, InvLaPlace, LaPlace, and
MagPhasePlot, respectively. Pretty much a transparent change.
Added new alias: RootLocusPlot for RootLocus.
(b) Mathematica 2.0 has a new undocumented function called ListQ
which replaces the one provide by us (our definition for ListQ
is still used for Mathematica 1.2 and lower).
(c) Mathematica 2.0 represents collections (headless expressions
like "(0,2,1,0)") using the head (data tag) Sequence.
So, ToCollection has been modified to return a collection
with the proper syntax under both Mathematica 1.2 and 2.0.
As an indirect result, the functions GetAllExponents and
GetAllFactors have been rewritten to return lists instead
of collections. We have removed the function ExtractVariables,
but GetVariables was kept so use it instead.
(d) The evaluation strategy for the Which construct has changed.
Now, each conditional statement in the construct MUST evaluate
to True or False. We now force all of our Which conditional
statements to evaluate to True or False by using TrueQ.
(35) Mathematica 2.0 offers the DeclarePackage which allows one to specify
the objects defined by a package. That way, when a user invokes
a new primitive, Mathematica 2.0 will automatically load the
package that defines it and then continue evaluating the expression.
This works for most objects but does not handle the case where
additional packages build up a definition incrementally.
The initialization files "Analog.m", "Digital.m", and "Support.m"
now only define the ``stubs'' for the objects in the analog, digital,
and support signal processing packages, respectively. This means
that you can place Needs[ "SignalProcessing`SignalProcessing`" ]
in the initialization file without noticeably slowing down the
initialization process. That way, students can go about their
business without having to know which packages to load in beforehand.
(36) One warning about Mathematica 2.0: it now assumes that variables
represent complex (instead of real) numbers, so Sqrt[x^2] remains
unsimplified. Although this is usually beneficial, this assumption
sometimes clashes when you want to treat a variable as real-valued.
For example, Plot[ x^(1/3), {x, -5, 5} ] cannot handle x in [-5, 0)
because it takes the wrong branch cut at the origin which can be
seen by plotting the real part of x^(1/3). It does plot the
correct curve over x in [0, 5].
(37) For the most part, the notebook interface in 2.0 is upwardly
compatible with that of Mathematica 1.2. However, 2.0 includes the
preamble with all PostScript graphics now (an extra 20 kb of
information per graphics cell) so some of our notebooks, like the
piecewise convolution, will double in size. This change caused us
to break up the z-transform notebook into three notebooks.
------------------------- Future Plans --------------------------------
(38) By evaluating Needs[ "SignalProcessing`ObjectOriented`" ], one
will be able to attach properties to (parameterized) operators using
the routine DefSystem. The default properties (ASSOCIATIVE, LINEAR,
SHIFTINVARIANT, etc.) are maintained in the variable SPproperties.
They can be attached to Mathematica's built-in operators (Re, Im,
Conjugate, etc.) or to the new operators (Shift, Upsample, etc.).
(39) Syntax checking and simplification rules have already been
developed for the new signals (functions) and systems (operators).
Next, we are developing a comprehensive set of rearrangement
rules for signals and systems based on their properties.
(40) New tutorial notebook on the discrete-time Fourier transform.
Because of changes (1), (5), and (6), we can recreate (for the
most part) the pictorial dictionary of Fourier transforms in
chapter 19 of Bracewell's book entitled "The Fourier Transform and
Its Applications".
Because of changes (9) and (12) and because of subtle modifications to
the ExtractVariables object, DSPAnalyze now works correctly for
IIR and FIR structures.
I have included a release of all of the files that compose the
signal processing packages. We have updated all of the notebooks
because of the above changes. We have included all eight notebooks.
The signal processing packages alter these standard Mathematica functions:
(a) Integrate so that it can handle generalized functions like the
Dirac delta (Delta), step (CStep), and pulse (Pulse) functions.
(b) Simplify applied to an And expression removes redundant conditions
(c) Det of a number is that number
(d) Dot product of two numbers is the product of the two numbers
(e) The TeX forms of Re and Im primitives are Re and Im instead of
R and I.