Discrete-time LTI systems in Haskell II

Yesterday I wrote a post on how to implement discrete-time LTI systems in Haskell. The code I presented did not quite please my delicate aesthetic sensitivity, as it bundled the declaration of the signals with the declaration of the LTI system under study. This is not practical. Ideally, I would like to define the LTI system in a Haskell script, run it on the GHCi interpreter, then create the signals using the GHCi command line. Why? Because I can then compute the response of the system to various test input signals without reloading the script.

In this post, we will again consider the first-order causal discrete-time LTI system described by the difference equation

y (n) - \alpha y (n-1) = u (n)

where |\alpha| < 1. In my previous post, I obtained a state-space representation of this LTI system. In this post, I will propose a better implementation of the given system in Haskell.

__________

Signal and System type synonyms

What is a discrete-time signal? We can think of a discrete-time signal as a sequence of numbers or symbols. In Haskell, we often use lists of floating-point numbers to represent real-valued discrete-time signals. More generally, we have the type synonym:

type Signal a = [a]

which says that a signal of type a is a list of elements of type a. We can then have signals of various types: integer, fixed-point, floating-point, bit, binary word, etc. This should be particularly useful in case we want to simulate mixed-signal circuits.

What is a discrete-time system? It is that which maps input discrete-time signals to output discrete-time signals. We will use the type synonym:

type System a b = (Signal a -> Signal b)

which tells us that a system of type a \,\, b is a function that maps an input signal of type a to an output signal of type b. In other words, it takes a list whose elements are of type a, and it returns a list whose elements are of type b.

__________

Implementing the LTI system in Haskell

In my previous post, I created the state and output sequences using

xs = scanl f x0 us
ys = zipWith g xs us

Combining the two lines above into a single line, we obtain

ys = zipWith g (scanl f x0 us) us

which returns the output sequence when given the input sequence. Note that the state sequence is not created. Assuming zero initial condition, we can then create a system (i.e., a function that takes a list and returns a list) as follows

sys :: Floating a => System a a
sys us = zipWith g (scanl f 0.0 us) us

which takes signals of type a and returns signals of the same type, where a is a floating-point type (either Float or Double).

Putting it all together in a script, we finally obtain:

type Signal a = [a]

type System a b = (Signal a -> Signal b)

-- state-space model
f,g :: Floating a => a -> a -> a
f x u = (0.5 * x) + u  	-- state-transition
g x u = f x u		-- output

-- define LTI system
sys :: Floating a => System a a
sys us = zipWith g (scanl f 0.0 us) us

which defines the LTI system under study with \alpha = 1/2. For this choice of \alpha, the system is a lowpass filter.

__________

Example: lowpass-filtering a square wave

We run the script above on GHCi to create the desired LTI system. Let us now create an input signal, a square wave of period equal to 16 and duty cycle equal to 50%, and lowpass-filter it using the LTI system we created:

*Main> -- create square wave
*Main> let ones  = take 8 $ repeat 1.0 :: Signal Float
*Main> let zeros = take 8 $ repeat 0.0 :: Signal Float
*Main> let us = cycle (ones ++ zeros) :: Signal Float
*Main> -- create output sequence
*Main> let ys = sys us
*Main> -- check types
*Main> :type us
us :: Signal Float
*Main> :type ys
ys :: Signal Float
*Main> :type sys
sys :: Floating a => System a a

Finally, we compute the first 64 samples of the input and output signals using function take, and plot the signals using MATLAB:

In case you are wondering why the peak amplitude of the output signal is twice that of the input signal, compute the transfer function of the LTI system and note that the DC gain is equal to \frac{1}{1 - \alpha} = 2.

About these ads

Tags: , , , ,

2 Responses to “Discrete-time LTI systems in Haskell II”

  1. Dan Says:

    I was entertaining the idea of LTI modeling in haskell in my head recently and was pleased to find your post so timely for me. One item that I think warrants some more exploration (and of particular interest to me) is simulation of feedback control laws. That is, u as function of a state vector x and, more generally, t a well. The control law could be passed in as an argument function. I expect lazy evaluation would also have great advantages in simplifying the code. I haven’t yet had the opportunity to investigate this myself, but I’d be very interested in seeing a follow up post exploring these topics.

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Connecting to %s


Follow

Get every new post delivered to your Inbox.

Join 76 other followers

%d bloggers like this: