The purpose of this tutorial is to drive you though the design of a digital filter that takes an analog input signal applied to the ADC, then processes incoming samples, and finally provides an analog output using a DAC.

In the first part, the tutorial provides some basic background on digital filters. You may need some prerequisite knowledge, mostly on Laplace transform and z-transform, to get a full understanding of exposed matter, but you can still grab some useful tricks even if you are new to the field of digital filtering.

The common representation of first-order Low-Pass filter transfer function in the Laplace domain is:

$$H(p)=\frac{1}{1+\tau p}$$

The cutoff frequency of such filter is given by:

$$\omega_c =\frac{1}{\tau}= 2\pi f_c$$

Let us suppose that we want a cutoff frequency $f_c=2Hz$. We therefore have:

$$\tau =\frac{1}{2\pi f_c}\approx 0.08$$

$$H(p)=\frac{1}{1+0.08\, p}$$

The Laplace variable *p* (or *s* if you live anywhere in the world but in France :-)) is useful to represent transfer functions in the continuous-time domain. When dealing with digital systems, we usually sample incoming signal at regular time intervals. The Laplace transform is no more representative of such a discrete-time system.

The Laplace-transform equivalent for sampled signals is called the z-transform.

A common approach to move a Laplace transfer function into the z-domain is to use something called the “bilinear transform”. The bilinear transform is obtained by mapping z-plane to p-plane using:

$$z=e^{\, pT}$$

where *T* is the sampling period (in seconds). In the tutorial, we will use a ** 10ms** sampling period (i.e. we will get a new input from ADC every 10ms).

The inverse mapping gives:

$$p=\frac{1}{T}ln(z)$$

Using first term of *ln* series expansion, we have:

$$p\approx \frac{2}{T}\frac{(z-1)}{(z+1)}$$

Applying the bilinear transform to our first order Low-Pass filter then gives:

$$H(p)=\frac{1}{1+\tau p}$$

$$H(z)=\frac{1}{1+\tau(\frac{2}{T}\frac{(z-1)}{(z+1)})}$$

$$H(z)=\frac{1}{\frac{T(z+1)+2\tau(z-1)}{T(z+1)}}$$

$$H(z)=\frac{Tz+T}{(T-2\tau)+z(T+2\tau))}$$

$$H(z)=\frac{\frac{T}{(T+2\tau)}+\frac{T}{(T+2\tau)}z}{\frac{(T-2\tau)}{(T+2\tau)}+z}=\frac{\alpha + \alpha z}{\beta+z}$$

With $T=0.01s$ and $\tau=0.08s$

$\alpha = \frac{T}{(T+2\tau)}=0.0588235$ and $\beta = \frac{(T-2\tau)}{(T+2\tau)}=-0.8823529$

The z-domain transfer function of our digital Low-Pass filter, sampled at 100Hz (T=10ms), and having a cutoff frequency of fc=2Hz is finally:

$$H(z) = \frac{0.0588235 + 0.0588235z}{-0.8823529+z}$$

Luckily for you, scientific computing software such as Matlab® or Scilab® are able to process linear systems transfer function from continuous-time to discrete-time (and reversely) quite easily.

The following Scilab® script does the same bilinear transform as seen above:

```
clear;
s=poly(0,'s');
sl=syslin('c',(1)/(1+0.08*s)); // Continuous-time system in transfer form
slss=tf2ss(sl); // Now in state-space form
figure(1); // Plot Bode diagram
clf(1);
bode(sl, 0.1, 10000);
t=0:0.01:0.6; // Time scale
y = csim('step',t,slss); // Simulate step response
figure(2); // Plot step response
clf(2);
plot(t,y,'r.-');
sld=cls2dls(slss,0.01); // Continuous to discrete transform
sldt=ss2tf(sld); // Converts in transfer form
disp(sldt); // Display z-transform
```

Executing the above script provides the following outputs. The obtained z transfer function confirms the one we calculated earlier.

- Bode diagram

- Step response

- z-domain transfer function

The z-transform is a discrete-time representation of the filter that relates output samples to input samples:

$$H(z) = \frac{\alpha + \alpha z}{\beta+z}=\frac{Y(z)}{X(z)}$$

Where $X(z)$ represents the input samples and $Y(z)$ represents the output samples.

Let us turn this into something more practical from a programming perspective:

$$\frac{Y(z)}{X(z)} = \frac{\alpha + \alpha z}{\beta+z}$$

$$Y(z)\beta + Y(z)z = X(z)\alpha + X(z)\alpha z$$

Above equation is processed each time a new sample is available. While X(z) represents the current incoming signal sample, X(z)z represents the sample that will come (in the future) the next sampling period. More generally, multiplying X(z) or Y(z) by z^{n} is the same as putting an index *n* on the sample being processed.

n = 0 → current sample

n = 1 → next sample

n = -1 → previous sample

Because no filter can predict what the future will be, we can divide all the parts of the above equation by *z* to have:

$$Y(z)\beta z^{-1} + Y(z) = X(z)\alpha z^{-1} + X(z)\alpha$$

$$Y(z) = X(z)\alpha z^{-1} + X(z)\alpha-Y(z)\beta z^{-1}$$

The *z ^{-1}* factor represents a delay of one sample. The above equation can therefore be written:

$$y_n = \alpha x_{n-1} + \alpha x_n -\beta y_{n-1}$$

$$y_n = \alpha(x_{n-1} + x_n) -\beta y_{n-1}$$

where:

- y
_{n}is the output sample computed at the n^{th}iteration - y
_{n-1}is the output sample computed at the (n-1)^{th}iteration - x
_{n}is the input sample read at the n^{th}iteration - x
_{n-1}is the input sample read at the (n-1)^{th}iteration

The above equation is called the recursion equation, and can be directly applied to the filter computation within a program.

For instance, you can experiment the above equation with Excel (or any other tool you like):

- Create a column of sample indexes
**n** - Create a column of
**x**input samples (some ‘0’ and then some ‘1’ for a step response) - Create a column for
**y**outputs by using the above equation - Plot both
**x**and**y**data with**n**as x-axis.

Looking at the above step response, you may recognize a first order low-pass filter. Moreover, you'll notice that time is not there since each output sample is calculated after both the previous one and the input sample. This should make clear in your mind that actual filter time constant $\tau$ now depends on the time interval between each sample calculation, i.e. the sampling period $T$ so that x-axis can be changed from $n$ the sample index into time $t$.

In the previous equation, we observe that each output sample calculation requires:

- One addition → $x_{n-1} {\color{Red} \oplus} x_n$
- Two products → $\alpha {\color{Red} \otimes} (x_{n-1} + x_n)$ and $\beta {\color{Red} \otimes} y_{n-1}$
- Another addition → $\alpha(x_{n-1}+x_n) {\color{Red} \oplus} -\beta y_{n-1}$

If we assume that the sampling frequency if high compared to input signal frequency (which is usually the case), we can further simplify our equation by assuming that x_{n-1} and x_{n} are not very different.

$$y_n=\alpha(x_{n-1}+x_n)-\beta y_{n-1}$$

with $x_{n-1}\approx x_n$

$$y_{n}^{'}\approx 2\alpha x_n-\beta y_{n-1}$$

Now, we can notice that:

$$2\alpha =\frac{2T}{(T+2\tau)}\; and\; \beta=\frac{(T-2\tau)}{(T+2\tau)}$$

$$\beta+1=\frac{(T-2\tau)}{(T+2\tau)}+1=\frac{(T-2\tau)+(T+2\tau)}{(T+2\tau)}=\frac{2T}{(T+2\tau)}=2\alpha$$

$$2\alpha = \beta+1$$

The above equation rewrites:

$$y_{n}^{'}\approx (\beta+1) x_n-\beta y_{n-1}$$

Let us write

$$\gamma = (\beta+1)$$

$$y_{n}^{'}\approx \gamma x_n+(1-\gamma) y_{n-1}$$

With $\gamma = \frac{2T}{(T+2\tau)}$

We have removed one addition. Using again:

T=0.01s τ=0.08s

We have finally:

$$y_{n}^{'}\approx 0.11765 x_n+0.88235 y_{n-1}$$

The plot below compares both step responses of first (y) and simplified (y’) version of our filter. Based on that result, we can accept the assumption that *x _{n-1}*≈

Before going into the filter implementation, let us prepare an application that simply copies the input voltage from the ADC channel to the DAC output. The next sections review the necessary setup for involved peripherals.

Below is the setup for **ADC**. It is a one channel conversion therefore **DMA **is not involved:

```
/*
* ADC_Init()
* Initialize ADC for single channel conversion
* - Channel 11 -> pin PC1
*/
void BSP_ADC_Init()
{
// Enable GPIOC clock
RCC->AHBENR |= RCC_AHBENR_GPIOCEN;
// Configure pin PC1 as analog
GPIOC->MODER &= ~(GPIO_MODER_MODER1_Msk);
GPIOC->MODER |= (0x03 <<GPIO_MODER_MODER1_Pos);
// Enable ADC clock
RCC->APB2ENR |= RCC_APB2ENR_ADC1EN;
// Reset ADC configuration
ADC1->CR = 0x00000000;
ADC1->CFGR1 = 0x00000000;
ADC1->CFGR2 = 0x00000000;
ADC1->CHSELR = 0x00000000;
// Enable 12-bit continuous conversion mode
ADC1->CFGR1 |= (0x00 <<ADC_CFGR1_RES_Pos) | ADC_CFGR1_CONT;
// Select PCLK/2 as ADC clock
ADC1->CFGR2 |= (0x01 <<ADC_CFGR2_CKMODE_Pos);
// Set sampling time to 28.5 ADC clock cycles
ADC1->SMPR = 0x03;
// Select channel 11
ADC1->CHSELR |= ADC_CHSELR_CHSEL11;
// Enable ADC
ADC1->CR |= ADC_CR_ADEN;
// Start conversion
ADC1->CR |= ADC_CR_ADSTART;
}
```

3.2. DAC configuration

Same thing for the **DAC **setup. We have a basic configuration:

```
/*
* DAC_Init()
* Initialize DAC for a single output
* on channel 1 -> pin PA4
*/
void BSP_DAC_Init()
{
// Enable GPIOA clock
RCC->AHBENR |= RCC_AHBENR_GPIOAEN;
// Configure pin PA4 as analog
GPIOA->MODER &= ~GPIO_MODER_MODER4_Msk;
GPIOA->MODER |= (0x03 <<GPIO_MODER_MODER4_Pos);
// Enable DAC clock
RCC->APB1ENR |= RCC_APB1ENR_DACEN;
// Reset DAC configuration
DAC->CR = 0x00000000;
// Enable Channel 1
DAC->CR |= DAC_CR_EN1;
}
```

**TIM6 **timebase is tuned for 10ms period between update interrupts:

```
/*
* BSP_TIMER_Timebase_Init()
* TIM6 at 48MHz
* Prescaler = 48000 -> Counting frequency is 1kHz
* Auto-reload = 10 -> Update period is 10ms
* Enable Update Interrupt
*/
void BSP_TIMER_Timebase_Init()
{
// Enable TIM6 clock
RCC->APB1ENR |= RCC_APB1ENR_TIM6EN;
// Reset TIM6 configuration
TIM6->CR1 = 0x0000;
TIM6->CR2 = 0x0000;
// Set TIM6 prescaler
// Fck = 48MHz -> /48000 = 1kHz counting frequency
TIM6->PSC = (uint16_t) 48000 -1;
// Set TIM6 auto-reload register for 10ms
TIM6->ARR = (uint16_t) 10 -1;
// Enable auto-reload preload
TIM6->CR1 |= TIM_CR1_ARPE;
// Enable Interrupt upon Update Event
TIM6->DIER |= TIM_DIER_UIE;
// Start TIM6 counter
TIM6->CR1 |= TIM_CR1_CEN;
}
```

And **NVIC **is configure to allow **TIM6 **interrupts with priority 1:

```
/*
* BSP_NVIC_Init()
* Setup NVIC controller for desired interrupts
*/
void BSP_NVIC_Init()
{
// Set priority level 1 for TIM6 interrupt
NVIC_SetPriority(TIM6_DAC_IRQn, 1);
// Enable TIM6 interrupts
NVIC_EnableIRQ(TIM6_DAC_IRQn);
}
```

3.4. Follower stage

Finally, the **main()** function implementing follower stage is given below. in and out variables have been set as global in order to allow STM-Studio® monitoring:

```
// Global variables
...
uint8_t timebase_irq = 0;
uint16_t in, out;
// Main program
void main()
{
// Configure System Clock
SystemClock_Config();
// Initialize 10ms timebase
BSP_TIMER_Timebase_Init();
BSP_NVIC_Init();
// Initialize ADC
BSP_ADC_Init();
// Initialize DAC
BSP_DAC_Init();
// Initialize LED pin
BSP_LED_Init();
// Main loop
while(1)
{
// Do every 10ms
if (timebase_irq == 1)
{
// Start measure
BSP_LED_On();
// Input signal <-- ADC
while ( (ADC1->ISR & ADC_ISR_EOC) != ADC_ISR_EOC );
in = ADC1->DR;
// Follower stage
out = in;
// Output signal --> DAC
DAC->DHR12R1 = out;
// Stop measure
BSP_LED_Off();
timebase_irq = 0;
}
}
}
```

In order to test the follower stage, let us connect the user button pin **PC13 **to the **ADC **input pin **PC1**. This is a convenient way to generate full-scale input steps (0V → 3.3V) to the filter input pin **PC1**. Make sure that you have no code processing the push-button event. You don’t even need to call the **BSP_PB_Init()** function as all pins are input upon reset and we’re not going to read button state anyway.

Save all, build and run the application.

If you have an oscilloscope at hand, you can make some verifications:

- Probing
**PA5**, you should see a short pulse that occurs every 10ms, according to the settings of the timebase:

- Zooming horizontal, the width of the
**PA5**pulse tells you about the time it takes to execute the filter loop. At this moment, the process takes**1.5µs**to complete. Remember that we are doing nothing more than output = input (follower stage). It has to be short…

- Probing both
**PC1**and**PA4**and playing with the user push-button, you can make sure that the follower stage actually follows (**DAC**output on**PA4**is a copy of**ADC**input on**PC1**):

Alternatively, you can monitor both in and out variables with STM-Studio® while playing with the user button:

When you have checked that everything is working as expected, you are ready to start the filter implementation.

Commit name "Digital filter - Follower"Push onto Gitlab |

Remember that our filter is described by the following recursion equation:

$$y_{n}= 0.11765 x_n+0.88235 y_{n-1}$$

$$y_{n}= c_1\times x_n+c_2\times y_{n-1}$$

In a first time, we will implement this equation directly using *floats *for the two coefficients $c_1$ and $c_2$ .

Edit the main() function as shown below.

```
...
void main()
{
float c1, c2; // Filter coefs
float x = 0; // Filter input
float y = 0; // FIlter ouput
...
// Setup filter coefs
c1 = 0.11765f;
c2 = 0.88235f;
// Main loop
while(1)
{
// Do every 10ms
if (timebase_irq == 1)
{
// Start measure
BSP_LED_On();
// Input signal <-- ADC
while ( (ADC1->ISR & ADC_ISR_EOC) != ADC_ISR_EOC );
in = ADC1->DR;
// Filter stage
x = (float)in;
y = (c1*x) + (c2*y);
out = (uint16_t)y;
// Output signal --> DAC
DAC->DHR12R1 = out;
// Stop measure
BSP_LED_Off();
timebase_irq = 0;
}
}
}
```

Note that in *y = (c1*x) + (c2*y)* we already have y_{n} on the left-hand side of '=' and y_{n-1} on the right-hand side.

Build the project and run.

You can use STM-Studio® to get a quick view on the result:

For precise measurements, it’s better using an oscilloscope. Probe both input (**PC1**) and output (**PA4**) pins:

Remember that we want *τ *=80ms. About 95% of the final value should be reached within 3*τ* = 240ms. The figure above is a detailed view of the step (push-button) response. We can see that the filter behaves perfectly as expected. The filter final value is almost reached after 5*τ*.

Next, probe **PA5 **(filter processing time). In the figure below, we first note that the filter is processed quite fast compared to the sampling period (10ms), meaning that CPU is spending most of the time doing nothing, waiting for the next time base interrupt.

We can measure an execution time of about 16µs:

If you play with the user Push-Button, you may notice that processing time is not the same when input signal is high or low. It is actually a little bit faster (13.4µs) when input is 0 (button down) because arithmetic operations are likely simplified with null operands. Therefore, one can see that computation time is not deterministic, and depends on operands values. That’s a consequence of conditional branches in the assembly code that result from floating point arithmetic synthesis.

Commit name "Digital Filter - Floating point arithmetics"Push onto Gitlab |

Using **float **type to represent real numbers such as our coefficients *c1* and *c2* is a natural thing, but you may know that real numbers may also be represented by integers, assuming a virtual decimal point. Such representation is called **fixed-point representation** (‘float’ being for floating-point representation).

We are aware that **Cortex-M0** CPU does not have any dedicated hardware to perform floating point arithmetic (**Floating Point Unit, FPU**). So it is likely that using floats is slowing down code execution (i.e. floating point arithmetic produces many assembly lines to complete).

Try the code below where filter variables are represented by a 12-bit decimal place fixed-point:

```
...
// Main program
void main()
{
uint16_t c1, c2; // Filter coefs
uint16_t x = 0; // Filter input
uint32_t y = 0; // Filter ouput
...
// Setup filter coefs
c1 = 482; // 0.11765f *2^12
c2 = 3614; // 0.88235f *2^12
// Main loop
while(1)
{
// Do every 10ms
if (timebase_irq == 1)
{
// Start measure
BSP_LED_On();
// Input signal <-- ADC
while ( (ADC1->ISR & ADC_ISR_EOC) != ADC_ISR_EOC );
in = ADC1->DR;
// Filter stage
x = in;
y = (c1*x) + (c2*y);
y = (y >>12U); // Remove decimal part
out = y;
// Output signal --> DAC
DAC->DHR12R1 = out;
// Stop measure
BSP_LED_Off();
timebase_irq = 0;
}
}
}
```

Build the project and flash the target.

Make sure that the filter behaves exactly as before (i.e. when using floats):

If you are using STM-Studio®, you can also see the filter behavior. Just set the sampling period to its maximum (sequential loop) and display all the data points:

When this is done, use the **PA5 **pin to measure the filter processing time. You’ll find that this implementation is about **8x faster** (2µs) than the one making use of floats. For the very same behavior! Play with the user button and you will also notice that processing time is now deterministic (i.e. computation time does not depend on operands). That's very good indeed.

Commit name "Digital filter - Fixed point arithmetics"Push onto Gitlab |

In this tutorial, you have implemented a digital low-pass filter. With a little help from Scilab® you should be able to synthesize in the digital domain any filter as long as you have its transfer function.

The tutorial took you into 2 coding approaches, the second being 8x faster that the first one for the very same result.

Processing efficiency depends on your coding skills. It has consequences on application performance and more generally on the business (money).

A code that is more efficient allows you:

- to buy a less powerful, less expensive processor for the same task
- to reduce MCU clock speed and therefore reduce power consumption (i.e. increase battery life in the case of mobile application, reduce heat dissipation requirements…)
- spend more time in low-power modes

These benefits are not just “icing on the cake”. It is fundamental that you produce efficient code every time you can.

- Log in to post comments