SICE Annual Conference 2008
August 20-22, 2008, The University Electro-Communications, Japan

Interpolation of Nonuniformly Decimated Signals
via Sampled-data H ∞ Optimization
M. Nagahara†1 , M. Ogura†2 , and Y. Yamamoto†31†

Graduate School of Informatics, Kyoto University, Kyoto, Japannagahara@ieee.org, 2 ogura@acs.i.kyoto-u.ac.jp, 3 yy@i.kyoto-u.ac.jp

Abstract: In this paper, we consider signal interpolation of discrete-time signals which are decimated nonuniformly.
A conventional interpolation method is based on the sampling theorem, and the resulting system consists of an ideal
ﬁlter with complex coefﬁcients. On the other hand we adopt sampled-data H ∞ optimization, which can take account ofintersample behavior, and the optimal ﬁlter with real coefﬁcients is obtained. An example shows the effectiveness of ourmethod. By examples, we also show that there is the optimal decimation pattern.
Keywords: nonuniform decimation, interpolation, sampled-data control

1. INTRODUCTION
Interpolation is a fundamental operation in digital signal processing, and has many applications such as signal
reconstruction, signal compression/expansion, and resizing/rotating digital images, see [11], [3]. If digital data
to be interpolated are located uniformly on the time axis,
the uniform interpolation is executed by an expander anda digital ﬁlter (called an interpolation ﬁlter) [11], whichis conventionally designed via the sampling theorem.
Periodic nonuniform interpolation (or decimation)
also plays an important role in signal processing, such as
signal compression by nonuniform ﬁlterbanks [8], superresolution image processing [9], and time-interleaved AD
converters [10]. The design has been studied by many
researchers [12], [8], [13], [4], [5], in which the design
methods are based on the generalized sampling theorem.
The optimal ﬁlter (or the perfect reconstruction ﬁlter) is
an ideal lowpass ﬁlter with complex coefﬁcients [12],[11]. Since the ideal ﬁlter cannot be realized, approximation methods are also proposed, see in particular [12],[13].
On the other hand, real signals such as audio signal (esp. orchestral music) breaks the band-limiting assumption in the sampling theorem, that is, they havesome frequency components beyond the Nyquist frequency. In view of this, we have to take account of
the whole frequency range in designing interpolation systems. Sampled-data H ∞ optimization [1], [7] is very adequate for this purpose.
In this article, we ﬁrst deﬁne nonuniform decimation/interpolation. This deﬁnition includes the block decimation introduced in [8]. Then we formulate the interpolation problem as a sampled-data H ∞ optimization. The
optimal ﬁlter is given by a periodic system, which can berealized by a multirate ﬁlterbank. Design examples showthe effectiveness of our method. Moreover, we consider
by examples what is the optimal decimation pattern. We
show that although the decimation rate is the same, theoptimal value can differ if the pattern differs. That is,
the performance depends on the decimation pattern. This
property can be used in designing signal compression.

2. NONUNIFORM DECIMATION AND
INTERPOLATION
Consider a discrete-time signal x := {x0 , x1 , x2 , . . .}shown in Fig. 1. Then nonuniform decimation by M :=[1, 1, 0] (we call this a decimation pattern) is deﬁned asfollows (see Fig. 2).
(↓ M )x := {x0 , x1 , x3 , x4 , x6 , . . .}.

(1)

That is, we ﬁrst divide the time axis into segments oflength three (the number of the elements of M ), then,in each segment, retain the samples corresponding to 1
in M and discard the one corresponding to 0. This decimation includes so-called block decimation, in which theﬁrst R1 samples of each segment of R2 samples are retained while the rest are discarded [8]. By using our notation, the block decimation R2 : R1 is represented as ↓ Mwith
M = [1, . . . , 1, 0, . . . , 0].
� �� � � �� �R1R2 −R1Then we consider interpolation. First, we deﬁne thenonuniform expander ↑ M with M = [1, 1, 0] by(↑ M )x := {x0 , x1 , 0, x2 , x3 , 0, x4 , . . .}.
That is, we ﬁrst divide the time axis into segment oflength two (the number of the elements 1 of M ), theninsert 0 into the portion corresponding to 0 in M . Applying this to the decimated sequence (1), we havev := (↑ M )(↓ M )x = {x0 , x1 , 0, x3 , x4 , 0, x6 , . . .}.
Then, the interpolation is completed by ﬁltering v by adigital ﬁlter K (see Fig. 3 (a)).

3. DESIGN OF INTERPOLATIONFILTER
3.1 Design Problem
In this section, we consider a general case where thedecimation pattern is deﬁned by- 1151 -M := [b1 , b2 , . . . , bM ],bi ∈ {0, 1}.

PR0001/08/0000-1151 ¥400 © 2008 SICEx↓Mv↑MyK(a)xLM

�KEL−1Ny(b)Fig. 3 (a) Nonuniform decimation and interpolation, (b)lifted system

Fig. 1 Discrete-time signal xe−Lsw

Let i1 , i2 , . . . , iN (i1 < i2 < . . . < iN ) be the indices of bi ’s such that bi1 = · · · = biN = 1 (N is thenumber of ones in M ). Then the interpolation processK(↑ M )(↓ M ) is periodically time-varying, more precisely, it is (N, M )-periodic [6] or (N, M )-shift-invariant[2]. By this fact, the multirate system can be equivalently
represented as a time-invariant system via discrete-timelifting [1], which is also called polyphase decomposition[11]. Let LM be the discrete-time lifting operator, that is, 
x0xM




 x1   xM +1  LM : {x0 , x1 , . . .} �→  .  ,  .  , . . . .
 ..   .. xM −1x2M −1
By this, the interpolation process is equivalently represented asNShLME�K

Fig. 4 Error system

1 1 0 1 1 0 1 1 0 1 1 0 1Fig. 2 Nonuniform decimation (M = [1, 1, 0])�K(↑ M )(↓ M ) = L−1
N KELM ,� := LN KL−1 ,KF(2)(3)

where E = [Eij ] is an N × M matrix whose elements
are deﬁned as follows:�
1, (i, j) = (1, i1 ), (2, i2 ), . . . , (N, iN ),Eij =0, otherwise.
For example, if M = [1, 1, 0] (M = 3, N = 2, i1 = 1,i2 = 2) then�
�
1 0 0
E=.
0 1 0
� is a time-invariant system withNote that the lifted ﬁlter KN -inputs and N -outputs.
� consider
To design the interpolation ﬁlter K (or K),
the error system shown in Fig. 4. In this ﬁgure, F is a linear time-invariant continuous-time system whose transfer
function is ﬁnite-dimensional and strictly proper, which
is a model of the original analog signal. The block Sh
represents the ideal sampler with sampling period h, andL−1N

Hh−eHh the zero-order hold with the same sampling period.
The delay e−Ls is a design parameter which control thereconstruction delay and the performance Then our problem is formulated as a sampled-data H ∞ optimization.
� which minimizes
Problem 1: Find the optimal ﬁlter K
the H ∞ norm (the L2 -induced norm) of the error system
Tew from the continuous-time signal w to the error e. The
L2 -induced norm of Tew is deﬁned by
�Tew w�2.
�w�2w∈L�Tew � := sup2(4)w�=0

Note that the norm (4) equals the H ∞ norm of thesampled-data system Tew [1]. The problem is therefore
called a sampled-data H ∞ optimization problem.
3.2 Filter design and implementation
The problem formulated above is a standard sampleddata signal reconstruction problem. The H ∞ optimal ﬁl� can be obtained by the fast-sampling method, seeter K[7].
The ﬁlter K is obtained by the inverse discrete-time
lifting of (3), that is�
K = L−1
N KLN .

Then the interpolation system is K(↑ M ), see Fig.3 (a)and the equation (2).
There is however another simpler way to implement
the interpolation system, by using a multirate ﬁlterbank,see Fig. 5. In this ﬁlterbank, Fi1 (z), Fi2 (z), . . ., FiN (z)are obtained by the following equation:1Fi1 (z)−1 Fi2 (z)  � N � z.
 = K(z ) 
..
..



.
.
−N +1
FiN (z)zThe uniform decimator ↓ M and expander ↑ M where M
is a positive integer is deﬁned by- 1152 -y = (↓ M )x,
w = (↑ M )v,

yn = xM n , n = 0, 1, 2, . . .
�xn/M , if n = 0, M, 2M, . . .,wn =0, otherwise.

z i1 −1↓M

↑M

Fi1 (z)+z i2 −1↓M

↑M

Fi2 (z)+−21−21.5Gain (dB)−22z iN −1↓M

↑MFiN (z)−23Fig. 5 Nonuniform ﬁlterbankz↓3↑3F1 (z)↓3↑3F2 (z)

−22.5

−23.5
−24 −310+10−2−110
Frequency (rad/s)

10010

1Fig. 7 Frequency response: proposed (solid) and conventional (dots)

Fig. 6 Nonuniform ﬁlterbank (M = [1, 1, 0])2
1.5

By using our deﬁnition of nonuniform decimator and expander, the uniform ones are given by↓ M =↓ M ,↑ M =↑ M ,10.5

M = [1, 0, 0, . . . , 0].
� �� �0M −1−0.5−1Fig. 6 shows an example of a nonuniform ﬁlterbank whenM = [1, 1, 0].

−1.5

4. DESIGN EXAMPLES−2150

In this section, we show design examples.
4.1 Optimal ﬁlter design
Here we design the optimal ﬁlter K (or Fi1 ,. . . FiN inFig. 5). The design parameters are as follows: the decimation vector M = [1, 1, 0], the sampling period h = 1,the time delay L = 12. The analog characteristic of thecontinuous-time input signals is modeled byF (s) =1.
10s + 1(5)For comparison, we adopt the method of Hilbert transformer [12], [11] as a conventional one. Note that thismethod is based on the sampling theorem, assuming thatthe original analog signal is fully band-limited up to thefrequency ω = 2π/3 (2/3 of the Nyquist frequency π).
Note also that the conventional ﬁlter requires very large
delay (L = 61.5). Fig. 7 shows the frequency responsesof the error system Tew shown in Fig. 4. The conventional
interpolation shows a large error in high frequency, while
the sampled-data H ∞ optimal interpolation shows a ﬂatresponse. To illustrate the difference between these responses, we simulate interpolation of a rectangular wave.
Fig. 8 shows the time response. The conventional interpolation causes large ripples, while our interpolation showsa better response. This is because the rectangular wave
has high frequency components around the edges, andour interpolation takes account of such frequency components.

155

160

165

170

175

Fig. 8 Time response: proposed (solid), conventional(dash), and input signal(dots)
4.2 Decimation pattern analysis
Consider M = 3 and N = 2. Then there are threepatterns of decimation: M1 = [1, 1, 0], M2 = [1, 0, 1],and M3 = [0, 1, 1]. These are essentially the same except
for delays, that is,z −1 (↑ M2 )(↓ M2 ) = (↑ M1 )(↓ M1 )z −1 ,(↑ M3 )(↓ M3 )z −1 = z −1 (↑ M1 )(↓ M1 ).
The optimal values of our H ∞ optimization for thesedecimations are therefore the same. However, when M =4 and N = 2, there can be difference. In this case, the essential patterns are M = [1, 1, 0, 0] or M = [1, 0, 1, 0].
Set the design parameters F (s) as (5), h = 1, and L = M(the length of the segment). Table 1 shows the optimalvalues γ = min �Tew � (see (4)). This result shows thatthe pattern M = [1, 0, 1, 0] (or M = [0, 1, 0, 1]) is thebetter, which is equals to the uniform decimation ↓ 2.
We then design when the segment length M = 5. Table 2 shows the result. By this, the optimal value γ depends on the position of the zeros in M , not depends on
the number of the ones in M . For example, although thepattern (F) retains more samples than the pattern (E), theoptimal values are the same. Table 3 shows the optimalγ when M = 7 and N = 4. In this case, there are 5 the- 1153 -Table 1 Optimal value γ for M = 4 and N = 2

Table 3 Optimal value γ for M = 7, N = 4

Decimation Pattern

Decimation Patternγ
0.2293γ0.3089

M = [1100], [1001], [0110], [0011](A)
0.1529M = [1010], [0101](B)

Table 2 Optimal value γ for M = 5 and N = 1, 2, 3, 4
Patternγ

Pattern

0.3813(A)

0.2320γ(C)

0.2303(D)

0.1536(E)(D)
0.15470.3062(B)(E)
0.2303(C)

0.1536(F)

essential patterns (A) to (E). We can see that γ depends
on the maximal number of the consecutive zeros in M
(we here call this the consecutive number).
The results shows that the reconstruction performance
depends on the consecutive number and not on the number of retained samples. By this observation, we can conjecture that the optimal decimation pattern M is the pattern in which the zeros are least consecutive. In otherwords, most uniformly distributed pattern is the best. Inview of this, the block decimation introduced in [8] cannot be optimal.

5. CONCLUSION
We have proposed an interpolation method of nonuniform decimation via sampled-data H ∞ optimization. By
a design example we have shown the effectiveness of ourmethod. We have also suggested by examples that there
can be the optimal decimation pattern.

REFERENCES[1]T. Chen and B. Francis, Optimal Sampled-Data
Control Systems, Springer, 1995.
[2] T. Chen, L. Qiu, and E.-W. Bai, General multirate
building structures with application to nonuniform
ﬁlter banks, IEEE Trans. Circuits and Systems — II,
Vol. 45, No. 8, pp. 948–958, 1998.
[3] R. Klette and P. Zamperoni, Handbook of Image
Processing Operators, John Wiley & Sons, 1996.
[4] Y.-P. Lin and P. P. Vaidyanathan, Periodicallynonuniform sampling of bandpass signals, IEEE

Trans. Circuits and Systems — II, Vol. 45, No. 3,
1998.
[5] R. J. Marks II and S. Narayanan, Interpolation of
discrete periodic nonuniform decimation using alias
unraveling, Proc. of IEEE ISCAS, pp. 281–284,
2002.
[6] D. G. Meyer, A parametrization of stabilizing controllers for multirate sampled-data systems, IEEE
Trans. Automatic Control, Vol. 35, pp. 233–236,
1990.
[7] M. Nagahara and Y. Yamamoto, A new design forsample-rate converters, Proc. 39th IEEE CDC, pp.
4296–4301, 2000.
[8] K. Nayebi, T. P. Barnwell, III, and M. J. T. Smith
Nonuniform ﬁlter banks: A reconstruction and design theory, IEEE Trans. Signal Processing, Vol. 41,
No. 3, pp. 1114–1127, 1993.
[9] S. C. Park, M. K. Park, and M. G. Kang,Super-resolution image reconstruction: A technicaloverview, IEEE Signal Processing Magazine, Vol.
20, pp. 21–36, 2003.
[10] T. Strohmer and J. Tanner, Fast reconstruction algorithms for periodic nonuniform sampling withapplications to time-interleaved ADCs, Proc. ofICASSP 2007, Vol. 3, pp. III-881–III-884, 2007.
[11] P. P. Vaidyanathan, Multirate Systems and Filter
Banks, Prentice-Hall, 1993.
[12] P. P. Vaidyanathan and V. C. Liu, Efﬁcient reconstruction of band-limited sequences from nonuniformly decimated versions by use of polyphase ﬁlter banks, IEEE Trans. on ASSP, Vol. 38, pp. 1927–1936, 1990.
[13] L. Vandendorpe, P. Delogne, B. Maison, and L. Cuvelier, MMSE design of interpolation and downsampling FIR ﬁlters in the context of periodic
nonuniform sampling, IEEE Trans. Signal Processing, Vol. 45, No. 5, 1997.

- 1154 -