Fullscreen HTML5Med Res (??MB)HTML4

Circuit Surgery Regular Regular clinic clinic by by Ian Ian Bell Bell Topics in digital signal processing – sinc filters W e are looking at various topics related to digital signal processing (DSP). DSP covers a wide range of electronics applications where signals are manipulated, analysed, generated, stored or displayed as digital data but originate from and/or are converted to real-world signals for interaction with humans or other parts of the physical world. Fig.1 shows the key elements of a generic DSP system with a signal path from an analog input via digital processing to an analog output. This does not necessarily represent every DSP system (not all have all the parts shown), but it serves as a reference for the various subsystems we will consider. Last month, we discussed the general structure of basic digital filters, which allows us implement a wide range of filters by using different sets of coefficients. This month, we will start to look at how those coefficients are calculated to achieve the desired filter. Digital filter design Last month, we also continued our discussion on implementing digital filters in LTspice The LTspice electronics simulation software package enables us to obtain transient and frequency response analyses based on the maths of a digital filter (we use LTspice behavioural sources to implement the filter maths). However, additional circuitry, such as antialiasing and reconstruction filters, could be implemented using more realistic circuit models to help investigate overall system performance. The additional circuitry can be also be implemented in idealised behavioural form to investigate the general principles of DSP systems. Analogue In Antialiasing filter Sample and hold The example filter used over the last couple of articles was a moving average filter. This was chosen as an example because it is easy to understand the function it is performing and the calculations behind it. Moving average filters are very effective at removing random noise from signals, particularly where the timing needs to be preserved as much as possible through the filtering process. However, they are not very useful for frequency separation, ie, rejecting one range of frequencies while passing another. Sinc filters There are many approaches to designing digital filters, so we will not provide comprehensive coverage of the types and techniques in this series. However, it is useful to go beyond the moving average filter and look at a slightly more complex, but still reasonably straightforward, approach called the ‘windowed sinc’ filter. This type of filter can provide good frequency separation, but it is not as good as a moving average in terms of timing and pulse response behaviour. We will discuss the sinc aspect of the filter this month, and leave windowing for next month. The term ‘sinc’ in the filter’s name indicates the mathematical function used. We introduced the sinc function in the October 2024 Circuit Surgery column in the context of DAC frequency response. sinc ( x )= To recap, the sinc function is defined as: ( πx ) When x is zero, wesin sinc ( x )= get 0/0, which is mathematicallyπ problematic πx (due to division by zero). However, for consistency 2 f B sincπ ( 2 f B t ) Digital ADC Analogue 2 f cDAC sincπ (Reconstruction 2 f ck) Digital processing filter i c π Out N−1 i− ( ( 2 )) Fig.1: a generic digital signal processing (DSP) system a =2 f structure. sinc 2 f Practical Electronics | February | 2025 sin ( x ) x and utility, sinc(0) is sin defined ( x ) to be 1. A sincfunction ( x )= is the normalised variant of this x sinc function, denoted by a subscript π: sinc π ( x )= sin ( πx ) πx This function occurs frequently in 2 f B sincπ ( 2 f B t ) DSP and telecommunications. The importance of the sinc function comes 2 f c sincπ ( 2between f c k ) the time from the relationship and frequency domain representations N−1 in of signals and frequency responses; a i =2 f c sinc 2 f c i− of ideal freπ contexts particular, in the 2 quency responses and stepped signals. ( ( )) Spectra Converting a signal waveform (time domain) to its frequency domain representation gives us the spectrum of the signal – a plot of the frequency content. We have shown the theoretical signal spectra several times in this series on DSP, and have used LTspice to plot the spectra of simulated waveforms. In the lab, many oscilloscopes and all spectrum analysers can display the spectra of real signals. Spectra provide useful information on signal quality, show the need for or effect of filtering, and help us understand phenomena such as aliasing. As we have mentioned before, mathematically, a function in the time domain is converted to the frequency domain using the Fourier transformation. LTspice uses an algorithm called the Fast Fourier Transform (FFT) to calculate waveform spectra. We can also convert from the frequency domain to the time domain using the inverse Fourier transform. Simulation files Most months, LTSpice is used to support descriptions and analysis in Circuit Surgery. The examples and files are available for download from the PE website: https://bit.ly/pe-downloads c 19 x sinc π ( x )= The mathematics are different for continuous and sampled signals. For sampled signals, the discrete Fourier transform (DFT) and its inverse are used (the FFT is an efficient implementation of the DFT). Fortunately, for basic design and analysis, we do not have to perform detailed mathematics. It is either built into software tools or measurement instruments, or forms the basis of techniques which do not require heavy maths. Impulse response Converting from a frequency response (the gain of a circuit such as a filter or amplifier versus frequency) to the time domain may be less familiar than obtaining signal spectra. However, it is another important concept and the result is the impulse response of the circuit. We introduced the idea of impulse responses last month. To recap, for continuous time systems, an impulse is an infinitely short pulse of unit area, but for sampled signal processing, it is a single sample of value 1. The impulse response is the output produced when the input to a circuit is an impulse at time zero and the input is otherwise zero for all time. We discussed the structure of basic digital filters last month. For the Finite Impulse Response (FIR) filter shown in Fig.2, the input and the input delayed by multiples of the sampling period are multiplied by the N coefficients a0 to aN-1 and summed to obtain the output. We showed that with an impulse input, the filter will output the values of its coefficients in sequence. To put it another way, the coefficients are equal to the impulse response of the filter, so if we know the impulse response of the desired filter we have the coefficient values required to implement it. a0 x(n) × a0x(n) Σ a1 Δt x(n–1) × a1x(n–1) Σ a2 Δt x(n–2) × a2x(n–2) y(n) It follows that if we want to try to create an ideal filter, it helps to find the impulse response via the inverse Fourier transform of an ideal frequency response (brick wall filter). The result for an ideal low-pass filter is that the impulse response is a sinc function. This is illustrated in Fig.3 for a continuous time filter. The upper plot shows the ideal frequency sin ( x ) response sinc filter ( x )=with a bandwidth for a low-pass x (cutoff frequency) of fB and gain of one (unity gain) in the pass ( πx ) The imsin band. sinc π (of x )= pulse function the filter (where t is πx time) is written as: 2 f B sincπ ( 2 f B t ) The frequency response plot in Fig.3 2 fpositive f c knegative ) c sincπ ( 2 shows both and frequencies, as we have done before, and as is commonly the case N−1 when disa i =2 fsignal c sincprocessing π 2 f c i− theory. The cussing 2 impulse response time axis is scaled by 2fB, so the graph is generally applicable and clearly shows the integral zero-crossing points of the normalised sinc function. ( ( )) Coefficients Fig.3 and the sinc function above apply to continuous time functions. For DSP with sampled data, things are slightly different. The cutoff frequency fc is expressed as a fraction of the sampling frequency fs, so for a bandwidth as shown in Fig.3, we would use f c = f B / f s. Given that the Nyquist limit constrains signals to be below half the sampling frequency to prevent aliasing, fc is in the range 0 to 0.5. The sinc function for the impulse response above is continuous in time t, but for a digital filter, we have a discrete set of coefficients that we can refer to using an index k, which references the sample times (these occur at times kT where T is the sample period). For the Gain Frequency response 1.2 1.0 0.8 0.6 0.4 0.2 0 –0.2 –fB f fB Memory × Fig.2: the structure of a finite impulse response (FIR) filter. 20 In the ideal case, the index (k) ranges N−1 from negative infinity to positive infinity a i =2 f c sinc π 2 f c i− because the sinc function has 2non-zero values over an infinite range. In reality, of course, we cannot have an infinite number of coefficients. The simplest solution is to use a reduced range that is practical. This means that the filter response will be non-ideal. As we will discuss in more detail next month, simply truncating the set of coefficients like this works, but does not produce the best filter response. The sinc function is symmetrical about zero, so the total number of coefficients (N) should be odd to have the same number of coefficients either side of zero. Therefore, the index (k) ranges from –(N – 1)/2 through zero to +(N – 1)/2. The positive part of the index range implies future values of the signal (current time is zero), so as discussed in December 2024 Circuit Surgery, such filters are referred to as non-causal and cannot be implemented. The solution is straightforward: the coefficients are time shifted so that the –(N – 1)/2 to +(N – 1)/2 index range maps to 0 to N – 1. For for 25 sinexample, (x) sinc x )= coefficients, the( theoretical index range x is –12 to +12, but in the practical implementation, the shiftedsin coefficient ( πx ) index sinc )=Fig.2) is 0 to 24. range (the i inπ a( x i in πx for the filter To calculate coefficients structure shown in Fig.3 with an index 2 f sinc ( 2–f 1,B twe ) use k = i i in the range Bof 0 toπN – (N – 1)/2 as the index in the impulse 2 f c sincabove. ) gives: π (2 f c k response function This ( ( ( ( a i =2 f c sincπ 2 f c i− )) N−1 2 )) The non-causal theoretical filter has zero delay, but the coefficient time-shift for the causal realisation results in a signal delay from input to output. This delay is equal to the coefficient shift, which is (N – 1)/2. So, for example, a 25-coefficient filter will have a delay of 12 sample periods. ∆fT 1.0 2fB 0.5 aNx(n–N–1) Processing 2 f c sincπ ( 2 f c k ) Gain 2fBsinc(2fBt) Σ Impulse response x(n–N–1) digital sinc2filter, the impulse response f sinc π (2 f B t ) is given by: B Inverse Fourier transform aN–1 Δt sin ( πx ) πx –10 –5 0 5 10 2fBt Fig.3: the impulse response of an ideal low-pass filter. 0.0 0 fC 0.5 f/fS Fig.4: the cutoff frequency and transition bandwidth of a low-pass filter. Practical Electronics | February | 2025 If it is necessary to combine or compare the filtered and non-filtered signals, the non-filtered signal will have to be delayed by the number of sample periods equal to the filter’s delay. Transition bandwidth The number of coefficients determines how sharp the cutoff of the filter is; the more coefficients, the sharper the cutoff. This can be expressed numerically as the transition bandwidth (fT) or roll-off. The transition bandwidth is the fraction of the sampling frequency (fs) over which the filter gain changes from the pass band to the stop band. Exactly where the transition band starts and ends is a somewhat arbitrary definition but could be, for example, from 99% to 1% of the gain in the pass band. Fig.4 shows an example where the transition bandwidth is about 0.1. A commonly used approximation for the width of the transition band is 4 / N where there are N coefficients, but this applies to the improved (windowed) versions, which we have not discussed yet. Even then, it varies with the details of implementation. For simple truncation of the coefficients, as described above, the transition bandwidth is better estimated using 0.9 / N. The filter coefficients (ai) can be calculated as described above. Before doing this, the number of coefficients has to be decided. This could be simply a decision based on what is considered achievable. Alternatively, the desired transition bandwidth could be used to calculate the number of coefficients required – assuming the relationship such as fT = 0.9 / N is known for the filter being designed. And as can be seen from the reciprocal relationship to N, sharper cutoffs require more coefficients. If the filter is implemented in digital hardware, this implies a larger circuit. For software implementations, more calculations are required, which requires more time and therefore reduces the maximum sampling rate which can be used. After all, the calculation process has to be completed within each sampling period for the filter to keep up with the incoming signal data. Increasing the number of coefficients also increases the delay of the signal through the filter. One disadvantage of sinc filters is that they are relatively inefficient in terms of processing time. However, faster approaches tend to be more complex and difficult to work with. Calculating the coefficients If we apply a constant value of 1 to Practical Electronics | February | 2025 Although there are software tools and online calculators to find digital filter coefficients, it is instructive to perform some of these calculations to see how things work, so we will go through an example here. Filter coefficients can be calculated manually, but it is more efficient to use a spreadsheet such as Microsoft Excel and LibreOffice Calc. Alternatively, if you can write code (for example in C or Python), it only takes a few lines to perform the calculations to find low-pass sinc filter coefficients. a low-pass sinc filter, the output will settle to a value equal to the sum of all the coefficients. This should be clear from Fig.2 – if the input and all the stored samples are 1, the output is equal to the sum of all the coefficients multiplied by 1, that is the sum of the coefficients. This is the filter’s gain at DC or very low frequencies, well into the pass band. If a filter of unity gain is required, the coefficients obtained from the impulse response can be normalised by dividing each of the values by the sum of the coefficients. Fig.5: this spreadsheet can calculate the coefficients for a sinc low-pass filter. Amplitude 0.2 0.15 0.1 0.05 –6 –5 –4 –3 –2 –1 0 –0.05 0 1 2 3 4 5 6 Sample index Fig.6: the impulse response of the filter obtained from the spreadsheet in Fig.5. 21 The spreadsheet shown in Fig.5 calculates the low-pass sinc filter coefficients from the sampling frequency (f s), required bandwidth (f B) and number of coefficients (N) using the equation for ai given above. N must be odd. These values are entered on the first three rows. In row 5, the cutoff relative to the sample frequency is calculated as fc = fB / fs using the formula =C2/C1. This value occurs in the functions used to calculate the coefficients and is referenced using $C$5. For readers not very familiar with spreadsheets, if the content of a cell (box in the grid) starts with an equals sign (=), it is interpreted as a formula and the spreadsheet performs the given calculation. The formulae can use basic arithmetic and built-in maths operations (such as SIN() for sine). The values of cells are referred to by their row and column coordinates. For example, the value of 10000 for fs in the example is obtained using C1. It is common for the same formula to be used on many rows, using data that varies down a column (or vice versa). Therefore, Excel automatic adjusts the cell references when performing copy and paste operations on formulae. This makes it easy to create many rows of similar calculation by pasting. This goes wrong where a constant value is required in all the functions over multiple rows (eg, f c is used in all the coefficient calculations). To use a constant cell reference, a $ character is included in cell coordinate, eg, $C$5 rather than C5 for the cell containing f c. In row 6, the index shift is calculated as (N – 1) / 2 using the formula =(C3-1)/2. This value is used to find the index value k to use in the coefficient calculation from the index i used in the filter. It is referred to using $C$6. In this example, we use a low value of 13 for N to keep the spreadsheet small. The values of i can be seen ranging from 0 to 12 (N – 1) in column A in rows 10 to 22. The value of i in rows 11 to 22 is found by adding one to the value in the previous row. For example, cell A11 has the formula =A10+1. Columns B to D break the coefficient calculation into steps. This is not strictly necessary, but makes the spreadsheet formulae easier to work with, as large formulae can be difficult to read. Column B is used to calculate the index k from index i using the index shift value from cell C6. For example, in row 10, the formula is =A10-$C$6. Column C is used to calculate the value that is going into the sinc function (x) as x = 2fck using the formula (again in row 10) =2*$C$5*B10. In column D, the coefficient ai is calculated as 2fc. sincπ(x) using the formula (in row 10): =2*$C$5*(SIN(PI()*C10)/(PI()*C10)) The function PI() evaluates to the value of π. This formula is not applicable when x = 0, which occurs when i is in the middle of its range (6 in the example). Here, the value of the sinc function is set to 1, so the spreadsheet formula is different in that row. In this case, cell D16 has the formula =2*$C$5*1 to calculate 2fc × 1. The 1 is included as a reminder of sinc(0) being used, and the cell is highlighted. Fig.7: an LTspice schematic implementing a low-pass sinc filter. 22 Practical Electronics | February | 2025 Spreadsheet cell D24 uses the formula =SUM(D10:D22) to calculate the sum of all the coefficients. In column E, the value of each coefficient from column D is divided by the sum value using (for row 10) =D10/$D$24. This produces the set of normalised coefficients, resulting in a filter with a gain of 1 in the pass band. Fig.6 shows the impulse response of the filter obtained from the spreadsheet in Fig.5. This is a plot of the normalised coefficient values against sample index k. Simulation example The LTspice schematic in Fig.7 implements the digital filter calculations for a 25-coefficient, 1kHz low-pass sinc filter with a sampling rate of 48kHz (a 20.83μs sample period). This schematic is just the digital filter part, but uses the same approach as previous examples. If required, sampling, antialiasing, reconstruction filters and other system or analysis components can be added in line with previous discussions. As explained last month, we have placed the coefficients and filter calculation function in an include file, rather than directly on the schematic. This code is shown in Fig.8. The reason we do this is to make it easier to see and edit the coefficients and to avoid cluttering the schematic up with too much extraneous text. The coefficients used here were calculated using a spreadsheet similar to the example discussed above. The results from simulating the filter are shown in Fig.9. You will find similar plots in many DSP chip data sheets as such PE digital filters are widely used. .func a(i) {table(i, + 0,0.02951210551, + 1,0.03191959113, + 2,0.03420780599, + 3,0.03635417297, + 4,0.03833734989, + 5,0.04013747349, + 6,0.04173642024, + 7,0.04311799631, + 8,0.04426815733, + 9,0.04517517611, + 10,0.04582976922, + 11,0.04622523114, + 12,0.04635750502, + 13,0.04622523114, + 14,0.04582976922, + 15,0.04517517611, + 16,0.04426815733, + 17,0.04311799631, + 18,0.04173642024, + 19,0.04013747349, + 20,0.03833734989, + 21,0.03635417297, + 22,0.03420780599, + 23,0.03191959113, + 24,0.02951210551)} .func filter(x0,x1,x2,x3,x4,x5,x6,x7,x8,x9,x10,x11,x12, + x13,x14,x15,x16,x17,x18,x19,x20,x21,x22,x23,x24) { + a(0)*x24+a(1)*x23+a(2)*x22+a(3)*x21+a(4)*x20+ + a(5)*x19+a(6)*x18+a(7)*x17+a(8)*x16+a(9)*x15+ + a(10)*x14+a(11)*x13+a(12)*x12+a(13)*x11+a(14)*x10+ + a(15)*x9+a(16)*x8+a(17)*x7+a(18)*x6+a(19)*x5+ + a(20)*x4+a(21)*x3+a(22)*x2+a(23)*x1+a(24)*x0} Fig.8: the code in the include file (Coeffs25.txt) used in the LTspice schematic in Fig.6. Fig.9: the simulated frequency response of the filter in Fig.6. Practical Electronics | February | 2025 23