Feeds:
Posts

## Matlab and Fibonacci Numbers

Matlab has long been an immensely convenient tool for processing and visualizing signals in a DSP system. The tool makes it very easy to simultaneously visualize a signal in multiple domains. The following is a simple demonstration of some of these signal-processing capabilities.

The Fibonacci sequence is the sequence of numbers formed when each number in the sequence is equal to the sum of the two Fibonacci numbers immediately preceding that number. The sequence is initialized by the numbers 0 and 1, and looks something like this: 0,1,1,2,3,5,8,13,21,34,55,89,144,233,377,etc. It is trivial to design a digital filter that produces the Fibonacci sequence as it’s impulse response. The impulse response of a discrete-time or digital filter is the response produced when the filter is stimulated by the Kronecker delta function:

$\delta[n] = \begin{cases} 1, & n = 0 \\ 0, & otherwise \end{cases}$

If the initial state of such a filter is such that all internal elements are at zero, the filter described by the following difference equation can produce an impulse response equal to the Fibonacci sequence. The input and output sequences are denoted x[n] and y[n], respectively.

$y[n] = x[n] + y[n-1] + y[n-2]$

Taking the Z-transform of this difference equation yields the following transfer function H(z) through the filter:

$\begin{array}{lcl} H(z) & = & \frac{Y(z)}{X(z)} \\ & = & \frac{1}{1 - z^{-1} - z^{-2}} \\ & = & \frac{z^{2}}{z^{2} - z - 1} \end{array}$

This transfer function has two zeros at the origin, and a pair of real-valued poles at z=-0.618 and z=1.618. A pole-zero plot can quickly be produced using the roots and zplane commands in Matlab:

zplane(roots([1,0,0]), roots([1,-1,-1]))


For the transfer function H(z), the coefficient matrices of the numerator and denominator polynomials in z-1 are given by B=[1] and A=[1,-1,-1], respectively. The frequency response of the filter can easily be computed using the freqz command in Matlab. In the following snippet of code, H is the complex-valued DFT and W is a vector of radian frequencies at which H has been computed. The following plot shows the magnitude of the frequency response.

% numerator polynomial
B = [1];
% denominator polynomial
A = [1,-1,-1];
% compute frequency response
[H,W] = freqz(B,A);
% plot the magnitude response (normalize W to pi)
plot(W/pi,abs(H));


The frequency response is computed along the unit circle in the Z-domain (z=ejw). Since the response is symmetric, the above plot only shows W/pi going from 0 to 1. As expected, the magnitude response shows peaks near 0 and 1, due to the relative proximity of the poles to these two points.

In the time-domain (or sample-domain), the response of such a system to a particular input sequence can easily be computed using the filter command in Matlab. In the following fragment of code, the above filter polynomials are used, along with initial states of zero. The input stimulus sequence x is ten samples of a Kronecker delta function, starting at n=0 and ending at n=9. The filter output y is initially 0, but then proceeds through the next ten numbers in the Fibonacci sequence.

% Kronecker delta function starting from n=0 to n=9
x = [1,zeros(1,9)];
% initial conditions of filter
init = [0,0];
% compute the filter response
y = filter(B,A,x,init);


A quicker way to compute the impulse response is through the impz command. The following line of code not only gives the first 10 values of the impulse response, it also produces a stem-line plot of the impulse response.

% compute impulse response
impz(B,A,10)


A simpler way to implement a Fibonacci sequence generator would be to discard the input x and initialize the filter to [1,0] instead of [0,0]. In any case, the Fibonacci filter example is not really practical, as the numbers increase until the internal states exceed the numerical limitations of the implemented machine. An unbounded response indicates that the filter is unstable, and this is confirmed by the location of a pole outside the unit circle in the Z-domain.

## Sampling and zero-order holds

During a recent interview, the great musician Neil Young expressed his desire for high-quality formats for music downloads. In terms of popular music, high-quality refers to file formats which preserve high-resolution data (e.g. 24-bit) sampled at rates much higher than necessary (e.g. 192 kHz). Whereas recording the original sources at higher sampling rates may provide some benefits with respect to the particular equipment being used, preserving these sample rates for distribution of the final mixed music makes no sense. An excellent discussion of why this format is unnecessary can be found here (xiph.org).

Over the years, many “audiophiles” have insisted on creating new file formats, distributing the audio files at absurd sampling rates, for the sake of “remaining faithful to the original audio waveform.” While many people know about the sampling theorem, there is a common misconception present in the minds of the lay-person when they look at the image of a sampled waveform and try to apply their intuition: they see the output of a the sampling process as a disjointed and distorted-looking stair-step response.

It has become common practice to represent the sampled waveform through an analog-to-digital converter (ADC) as a stair-step response (including on this blog). This representation is not strictly correct because it presumes that signals produced at the output of an ADC have a continuous-time representation. What actually emerges from an ADC is a signal in the discrete-time domain, where the waveform discontinuous and only exists at the sampling instants. This may seem like a trivial point, but there are ramifications for the untrained eye. When someone who is not well-versed in signal processing theory views an image showing the classic stair-step sampled waveform, their mind intuitively views this as a grossly degraded version of the original waveform. This leads to scores of “audiophiles” to incorrectly assume that an audio signal sampled at 192 kHz is inherently “more accurate” than more traditional (and sufficient) rates of 44.1 kHz (compact disc) or 48 kHz.

In reality, the output of an ADC looks more like a discontinuous sequence of points (“dots”) which when interpolated recreate the original signal. When such an image is shown to the human eye, the sampled waveform does not appear as distorted as the stair-step representation. The digital (discrete-time) circuitry that follows the ADC has no concept of what the signal might look like in-between the samples. The signal only exists at the active clock edges, and as long as Nyquist is satisfied the samples accurately represent the input waveform (assuming all setup-time and hold-time constraints also remain satisfied).

In order to understand how the stair-step response comes about, we need to consider the operation of a digital-to-analog converter (DAC). When converting a discrete-time signal to a continuous-time waveform, something known as a reconstruction filter is required. This reconstruction filter is specially designed to produce a continuous-time output when provided with a discrete-time input. A common type of DAC reconstruction filter is the zero-order hold, which is implemented by simply holding constant the previous sample until the next sample is encountered. The zero-order hold reconstruction filter is what leads to the aforementioned stair-step representation of the input signal. The sight of this repulsive-looking waveform leads to further questions. What do those “stair-steps” represent? Are they harmful? How do we remove these effects to recreate the original smooth signal? In order to answer these questions, we must dig deeper.

The filtering operation is basically a time-domain convolution of the input signal with the filter’s impulse response. This corresponds to a multiplication in the frequency domain. The impulse response of a zero-order hold reconstruction filter is a single square pulse, with a width equal to the sampling period. Its frequency-domain representation looks like a sinc function, which continues forever in both positive and negative frequency directions, with nulls at multiples of the sampling rate. Any discrete-time signal has a frequency-domain representation which contains an infinite number of copies of the input signal band, spaced at multiples of the sampling rate. The time-domain convolution of this signal with the reconstruction filter is equivalent to the multiplication of their frequency-domain representations.

As a result, the reconstructed signal still contains an infinite number of copies of the original waveform, albeit attenuated as we move further and further away from the origin in the frequency domain. The presence of these higher-frequency copies is what leads to the stair-step shape of the signal waveform. As long as the repeated copies can be removed without harming the primary signal band, the original signal can be perfectly reproduced without any loss. These copies need to be filtered out in order to leave us with a clean single spectral copy of the original waveform. This is usually achieved using a low-pass filter at the output. Throughout this signal-chain, there are practical issues that need to be dealt with, such as correcting the pass-band droop in both the reconstruction and the low-pass filters, as well as compensating for any phase non-linearities.

The point of all this is that what actually emerges at the output of an ADC is a series of instantaneous sample dots, floating in time and space, and ready to be consumed by the next discrete-time (digital) processing circuit. The human brain finds it much easier to spatially interpolate these points and imagine these to be a reasonably accurate representation. However, a stair-step depiction of the waveform is not only rejected by our intuition, but strictly speaking, it is also not what actually emerges from the ADC as digital samples. The stair-step waveform more closely represents an intermediate signal within a DAC that happens to use a zero-order hold reconstruction filter, and this is the wrong waveform to which the lay-person’s intuition should be applied.

## Journal publication

Today, my paper on Tunable Mismatch Shaping for Quadrature Bandpass Delta-Sigma Data Converters was published in the Journal of Signal Processing Systems, Volume 65, Number 2, pp 199-210, November 2011.

From the abstract:

This paper presents an architecture for quadrature bandpass mismatch shaping that allows the center frequency of the mismatch suppression band to be tunable over the entire Nyquist range. The approach is based on the previously reported complex-valued tree-based mismatch shaper, and extends this to allow tunable operation. The proposed design has been implemented using VHDL and synthesized to logic gates. The hardware complexity and mismatch shaping performance of the proposed architecture are compared to that of a reference architecture, which uses separate tunable mismatch shapers for each complex component path. Simulation results show consistent mismatch shaping performance across the entire tuning range.

## ICASSP paper accepted

My paper entitled A novel technique for Tunable Mismatch Shaping in Oversampled Digital-to-Analog Converters has been accepted for presentation at the 2010 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP). The conference takes place from March 14-19, 2010, at the Sheraton Dallas Hotel in Dallas, Texas. ICASSP is a technical conference focusing on signal processing theory and its applications. In its 35th year, the conference is expected to feature over 120 lecture and poster sessions.

## Quantization and SQNR

This is a continuation of a discussion about quantization and analog-to-digital converters. In that discussion, the normalized quantization step through an N-bit ADC was denoted q, where q = 1/2N. The ADC encoder transfer function yielded a quantization error range over the interval [-q/2,+q/2].

Quantization is a highly nonlinear process. Denoting the input and output of a quantizer as u[n] and uq[n], respectively, the error from quantization uq[n]-u[n] can be re-arranged to yield the additive noise model of quantization error: uq[n]=u[n]+e, where e is the quantization error.

The figure below shows the quantization error for a full-scale sine wave over a single period. Also shown is the quantization error for a full-scale sawtooth ramp signal.

Although the quantization error from the sinusoid is signal-dependent and nonlinear, the commonly used additive noise model assumes a stochastic process in order to simplify the analysis. In particular, the error is treated as an independent and identically distributed (i.i.d.) random variable.

If the quantization error is modeled as a random variable with a uniform distribution, the probability density function is given by:

$p(e) = \begin{cases} 0, & x \mbox{ \textless } -q/2 \\ 1/q, & -q/2 \leq x \leq +q/2 \\ 0, & x \mbox{ \textgreater } +q/2 \end{cases}$

The root mean-square (RMS) quantization error with such a distribution can thus be derived:

$\begin{array}{lcl} e_{RMS}^2 & = & E(e^2) \\ & = & \int\limits_{-q/2}^{+q/2} e^2 \cdot p(e) \,d e \\ & = & \frac{1}{q} \int\limits_{-q/2}^{+q/2} e^2 \,d e \\ & = & \frac{q^2}{12} \\ e_{RMS} & = & \frac{q}{\sqrt{12}} \end{array}$

The RMS value of a full-scale sinusoid whose peak-to-peak swing has been normalized to unity is given by:

$sig_{RMS} = \frac{1}{2\sqrt{2}}$

The signal-to-quantization-noise ratio (SQNR) through the ADC can then be computed and expressed in decibels (dB) as:

$\begin{array}{lcl} SQNR_{sig} & = & 10log_{10} \left( \frac{sig_{RMS}}{e_{RMS}} \right)^2 \\ & = & 10log_{10} \left( \frac{1}{2\sqrt{2}} / \frac{q}{\sqrt{12}} \right)^2 \end{array}$

Substituting q=1/2N gives:

$\begin{array}{lcl} SQNR_{sig} & = & 10log_{10} \left( 2^N \cdot \sqrt{3/2} \right)^2 \\ & = & 20log_{10} \left( 2^N \right) + 20log_{10} \left( \sqrt{3/2} \right) \\ & = & (6.02N + 1.76) dB \end{array}$

This is the well-known equation for SNR or dynamic range through an N-bit ADC using the additive noise model of quantization error, and in the absence of all other noise sources like thermal noise in the analog circuitry, dither and sampling jitter. Note that no over-sampling is assumed here.

This analysis assumes that quantization errors are uniformly distributed over the quantization interval. In reality, the errors are not uniformly distributed for a sinusoidal input. For example, referring back to the time-domain quantization error from a sinusoid and a sawtooth ramp shown in the figure above, the respective error distributions are shown in the figure below.

The quantization error of the sawtooth wave appears to be uniformly distributed, but that of the sinusoid is clearly not. This is due to the signal-dependence of the sinusoid’s quantization error. Since the sawtooth actually produces uniformly distributed quantization errors, it is instructive to compute the SQNR from quantizing such a signal.

The RMS value of a full-scale sawtooth whose peak-to-peak swing has been normalized to unity is given by:

$saw_{RMS} = \frac{1}{2\sqrt{3}}$

Using the RMS quantization error derived above for a uniformly distributed quantization error, the SQNR of a sawtooth wave applied to an ADC can be expressed as:

$\begin{array}{lcl} SQNR_{saw} & = & 10log_{10} \left( \frac{saw_{RMS}}{e_{RMS}} \right)^2 \\ & = & 10log_{10} \left( \frac{1}{2\sqrt{3}} / \frac{q}{\sqrt{12}} \right)^2 \\ & = & 10log_{10} \left( 2^N \right)^2 \\ & = & (6.02N) dB \end{array}$

In general, the computed SQNR depends on the signal source and the model used for the quantization error. For sinusoidal inputs, the approximation of uniformly distributed quantization error improves as the ADC precision increases.

The figure below compares the error distribution of the sawtooth with that of four ADC resolutions (3 bits, 6 bits, 9 bits, and 12 bits). Clearly, the distribution approaches the quantization model of a sawtooth as the ADC resolution is increased.

Modeling the SQNR as 6dB per bit of ADC precision is a good approximation, especially as the ADC precision asymptotically increases. For many signal processing applications, the usefulness of approximating the quantization error as an i.i.d. noise source, far exceeds the inaccuracy of the model.

An Analog-to-Digital Converter (ADC) does exactly what the name implies: it converts an analog electrical signal to a digital representation. Specifically, the analog signal is a continuous-time continuous-amplitude signal, and the digital signal representation produced by the typical ADC is a sequence of discrete-time discrete-amplitude samples. The process of conversion from a high-resolution signal to a low-resolution signal is also known as quantization.

The two main types of ADCs are oversampling converters and Nyquist-rate converters, and there are several architectures for these. In most cases, there is some form of uniform quantization being performed on a high-resolution signal, in order to represent it in terms of a finite set of quantization levels. The error that results from quantization is referred to as quantization error. The spectral representation of random quantization error is known as quantization noise.

An ADC requires a clock signal to synchronize the instances when the analog signal is sampled. The clock frequency is referred to as the sampling rate, i.e. the rate at which samples are taken, and can be denoted fs. It is important that the clock have little or no clock jitter, which creates uncertainty in the sampling instant, and hence increases the quantization error.

An ADC also typically requires a reference voltage, denoted by VREF, which determines the valid voltage range over which the analog input signal can be converted. The input range of an ADC that only operates on positive voltages would go from zero volts, or circuit ground, to VREF. If the analog signal takes on values outside this voltage range, a well-designed input circuit will non-catastrophically limit the ADC to either minimum or maximum voltage, depending on the input signal. As expected, this would produce either a minimum or maximum digital value at the output.

The most common representation used for the digital samples produced by an ADC is a string of binary digits (or bits), where 00..0 represents the smallest analog input, and 11..1 represents the largest. These are sometimes referred to as ADC output codes. An N-bit binary number can represent at most 2N unique levels, and therefore, an N-bit ADC can produce 2N unique codes.

The quantization step or width of each ADC code can be denoted q, where q = VREF/2N. The nominal ADC code width is expected to be equal to a single LSB (least-significant bit), which is the right-most bit in a binary word representation. When the code width is normalized to VREF, q = 1/2N.

In the figure above, an example of a 3-bit ADC encoder transfer function is shown on the left, relating the digital output to the analog input. The encoder transfer function is arranged so that any input signal less than q/2 produces the smallest digital code, 000, input signals between q/2 and 3q/2 produce the next digital code 001, and so on. Alternate arrangements are possible, depending on the specific application requirements of the ADC.

The quantization error resulting from using this encoder transfer function is shown in the figure above on the right, and in this case, it takes values over the interval [-q/2,+q/2]. This assumes that the ADC input is appropriately limited and the digital output code is saturated when the input signal goes outside the operating range of the ADC.

The figure below shows the result when a full-scale sine wave is provided at the input to an ADC with this encoder transfer function.

It should be apparent that quantization is a highly non-linear process, and this makes it very difficult to perform an exact analysis of an otherwise linear system. In order to use classical linear analysis, it is necessary to derive a suitable linearized model of the quantizer, and this will be covered in a future post.

## Tricks for dB scale conversion

The use of the decibel scale is ubiquitous in electronic systems. For example, the dynamic range through a system is measured by the signal-to-noise ratio in decibels, or dB. When working in the lab, or viewing the output of any measurement device, it is useful know how to make quick mental conversions between the linear and dB scales.

The amplitude or power of a signal is typically quantified relative to a fixed reference. A full-scale signal has a ratio of 1:1, and is expressed as 20.log(1/1) = 0 dB. Remember that multiplication (division) in the linear domain is equivalent to addition (subtraction) in the log domain. One needs to only remember a few values in order to compute most conversions:

1. A ten-fold (10x) increase or decrease in linear signal amplitude results in a +20 dB or -20 dB change on the decibel scale, respectively.
2. A doubling or halving (2x or ½x) of linear signal amplitude results in [approximately] a +6 dB or -6 dB change on the decibel scale, respectively.
3. A tripling (3x) of linear signal amplitude can be approximated by using 3 ≈ √10. The square-root is equivalent to a power of half, and in the log domain, this simply halves the dB value. This results in ½ x 20 dB = 10 dB.

Using these basic rules, it is easy to quickly compute the linear ratios corresponding the dB value. Note that in the list above, the 3x ↔ 10 dB conversion is the greatest source of error in the final approximation.

Some examples:

• Convert 54 dB to the linear scale: note that 54 dB = (60 – 6) dB, which is equivalent to 1000 x ½ = 500 in the linear domain (this is a good approximation to the actual value of 501.2)
• Convert 7x to the dB scale: note that 7 = √49 ≈ √50 = √(100 x ½), which is equivalent to ½ x (40 – 6) = 17 dB (the answer should be 16.9)
• Convert 30 dB to the linear scale: note that 30 dB = (20 + 10) dB, which is equivalent to (10 x 3) = 30 in the linear domain (the answer should be 31.6). Alternatively, we can use 30 dB = (40 – 10) dB, which converts to (100 ÷ 3) = 33.3 (the magnitude of the approximation error is about the same as that of the first answer)

In the computerized world of today (or when no one brings a calculator to the lab), these mental shortcuts can be very useful as a quick sanity check.