Appl. Comput. Math., V.11, N.1, 2012, pp.3-18

A BRIEF OVERVIEW OF SIGNAL RECONSTRUCTION VIASAMPLED-DATA H ∞ OPTIMIZATIONYUTAKA YAMAMOTO1 , MASAAKI NAGAHARA2, PRAMOD P. KHARGONEKAR3
Abstract. This paper gives a brief account of a new method for signal reconstruction by
employing the mathematical machinery from sampled-data control theory. We formulate thesignal reconstruction problem in terms of an analog performance optimization problem using astable discrete-time ﬁlter. The proposed H ∞ performance criterion naturally takes intersample
behavior into account, reﬂecting the energy distributions of the signal. We present a method
for computing optimal solutions which are guaranteed to be stable and causal. Comparisons toalternative methods are also presented.
Keywords: Sampled-Data Control Theory, Digital Signal Processing, Shannon Paradigm, Sampling Theorem, H ∞ Control Theory.
AMS Subject Classiﬁcation: 93C57, 60G35, 93B36, 94A20.

1. Introduction
Digitization of analog information such as voice, images, and video is a ubiquitous feature ofmodern communications and information technology. Therefore, reconstruction of the original
analog signal from its digital version is a problem of great interest and importance. These reconstruction procedures and used everywhere: music players, mobile telephones, cameras, television,etc. We note that there is a critical distinction between compression/transmission/recovery ofdigital data and our focus here with the reconstruction of the original analog data from whichsuch digital data are generated.
Shannon’s pioneering and much celebrated paper [27] was a landmark development in thissignal reconstruction problem. Under the assumption that the original analog signal is bandlimited (below the so-called Nyquist frequency), he showed that the original analog signal canbe exactly reconstructed by using the sampling theorem [44]. Shannon’s framework establisheda fundamental paradigm for digital signal processing. We will hereafter refer to this scheme theShannon paradigm. However, the band-limitedness assumption, necessary for perfect signal
reconstruction, is not easily satisﬁed. In many applications, the sampling rate is not highenough to allow for this assumption to hold even approximately. To address this problem, ananti-aliasing ﬁlter is often introduced to sharply cut high frequency components. This in turnleads to yet another type of distortion due to the Gibbs phenomenon (see Section 6 below).
Moreover, the sinc function, which is the impulse response of the ideal reconstruction ﬁlter,
is not causal and does not decay very fast. This slow decay rate makes it very diﬃcult toimplement. Various approximations thus become necessary. This procedure further complicatesthe digital ﬁlter design procedure, making it less transparent.
1

Department of Applied Analysis and Complex Dynamical Systems, Graduate School of Informatics, Kyoto
University, Kyoto, 606-8501, JAPAN, e-mail: yy@i.kyoto-u.ac.jp
2
Department of Applied Analysis and Complex Dynamical Systems, Graduate School of Informatics, Kyoto
University, Kyoto, 606-8501, JAPAN, e-mail: nagahara@i.kyoto-u.ac.jp
3
Department of Electrical and Computer Engineering, University of Florida, Gainesville, FL 32611-620, USAe-mail: ppk@uﬂ.edu
Manuscript received 16 December 2011.
3

4

APPL. COMPUT. MATH., V.11, N.1, 2012

Motivated by these issues, many researchers have been developing techniques that aim to solve
the signal reconstruction problem under the assumption of non-ideal signal acquisition devices.
An excellent survey of this stream of literature can be found in [30]. In this research framework,
the digital signal is assumed to be the sampled version of an analog signal processed through a
non-ideal analog ﬁlter. The reconstruction process then attempts to recover the original signal.
The idea is to place the problem into the framework of the (orthogonal or oblique) projection
theorem (in usually L2 ), and then project the signal space to the subspace generated by theshifted reconstruction functions. It is often required that the process give a consistent result,i.e., if we subject the reconstructed signal to the whole sampling process again, it should yield
the same sampled values from which it was reconstructed [28].
The principal objective of the present paper is to present a relatively brief summary of an
entirely new approach to the signal reconstruction problem. This approach takes its inspirationfrom the modern sampled-data control theory, developed in the control community since the1990’s. The fundamental accomplishment of modern sampled-data control theory is that itgives us a discrete-time controller (or ﬁlter) that optimizes the closed-loop performance withintersample behavior taken into account. In particular, it optimizes an analog (continuous-time)
performance metric. These metrics are given in terms of H ∞ or H 2 norms. This framework
gives us an ideal platform to reconstruct the original analog signals from their sampled-data
versions when original signal is not band-limited.
Chen and Francis [5] made a ﬁrst attempt to apply sampled-data control theory to signal processing (however in a discrete-time domain); see also [14]. Starting in 1995, the present authorsand our colleagues have pursued the signal reconstruction problem in the sampled-data context
to obtain an optimal analog performance via digital ﬁltering: See [19, 41, 36, 42] for generaldesign frameworks, [15, 25] for sample-rate conversion, [40] for multirate ﬁlterbank design, [1, 2]
for audio signal compression, [16] for image restoration, [26] for fractional delay ﬁlters, [17] forwavelet expansion, and [35, 39] for convergence analysis. The method has also been patented[37, 38, 12, 13] and implemented into sound processing LSI chips as a core technology by Sanyo
Semiconductors, and successfully used in mobile phones, digital voice recorders and MP3 players; their cumulative production has exceeded 30 million units as of 2011. Our objective here isto present a succinct exposition of some key aspects of this body of work.
The same philosophy of emphasizing the importance of analog performance was proposed
and pursued recently by Unser and co-workers [31, 32]. The crucial diﬀerence is however thatthey rely on L2 /H 2 type optimization and oblique projections, which are very diﬀerent from our
method here. In particular, it can raise some stability questions. The recent work of Meinsma
and Mirkin [21, 22] takes an approach that is close to ours. They give solutions for non-causal
problems and allow freedom in the choice of sample or hold devices. A detailed comparison of
our work and these related works is provided in [42]. Some other approaches (not very closely
related to our work) to extending the traditional sampling theory include: reconstruction byquasi-interpolation [6], and minimization of the worst-case regret [9].
The present paper is organized as follows: After preparing some basic notions in function
spaces in Section 2, we ﬁrst review the fundamentals in signal reconstruction using the sampling
theorem, and discuss its various drawbacks in Section 3. We will then give a fundamental setup
and formulation of our sampled-data ﬁlter design framework in Section 4, and Section 5 gives
a solution method via fast-sample/fast-hold approximation. Finally, we give some examples insignal reconstruction in Section 6.
In closing this Introduction, we note that this paper draws heavily from our very recent paper[42]. We refer the interested reader to this paper for detailed technical proofs, comparisons toother approaches, and more examples.

Y. YAMAMOTO et al.: A BRIEF OVERVIEW OF SIGNAL RECONSTRUCTION ...

5

2. Preliminaries
2.1. H ∞ optimization. Let us introduce some basic function spaces and performance measures. Let L2 (a, b) (or L2 [a, b), L2 [a, b], etc.) be the space of Lebesgue square integrable functions on the interval (a, b), a < b. For a function f valued in Rn or Cn , its L2 -norm is denoted
by
 b
1/2
2
|f (t)| dt,(1)
f 2 =
anwhere | · | denotes the Euclidean norm in C . Let H 2 denote the space of Cn -valued functions

f that are analytic on the open right half plane C+ := {s : Res > 0} and satisfy
 ∞|f (x + jy)|2 dy < ∞.
supx>0 −∞
2
The H -norm of a function f ∈ H 2 is deﬁned byf 2 := sup
x>0

1
2π

 ∞
−∞

2

|f (x + jy)| dy

1/2
.

(2)

It is well known that Laplace transform gives an isometry between L2 [0, ∞) and H 2 .
The space H ∞ denotes the Hardy space of functions analytic on C+ and bounded there. It
is a Banach space with normf ∞ := sup |f (s)|.

(3)

s∈C+

An element f of H ∞ admits nontangential limit to the imaginary axis almost everywhere,
which we denote by f (jω), ω ∈ R. Then the H ∞ -norm of f ∈ H ∞ is equal tof ∞ = esssup |f (jω)|.
−∞<ω<∞(4)

Now let G be the transfer function of a ﬁnite-dimensional, asymptotically stable linear
continuous-time system. Then G belongs to H ∞ , and its “size” is measured by the H ∞ -norm,i.e., the supremum (or maximum) of the Bode magnitude plot as in (4).
The steady-state response of G against a sinusoid ejωt is given by G(jω)ejωt , and its magnitude
is bounded as|G(jω)ejωt | ≤sup−∞<ω<∞|G(jω)| · |ejωt | = G∞ .

In general, for u ∈ H 2 , it is known that(5)
Gu2 ≤ G∞ · u2 ,∞2and this bound is tight. Hence the H norm gives the L energy induced-gain, and minimizing
it yields a system that works uniformly well for the whole frequency range1.
For this reason, it is recognized that the H ∞ -norm criterion is often superior to the H 2 -norm
criterion, where the H 2 -norm for a stable matrix transfer function is deﬁned as
  ∞
 1/2
1
∗
trace {G (jω)G(jω)}dt.
G2 :=
2π −∞The H ∞ norm has been used successfully in the control literature [7, 8, 11].
1However, it is to be noted that it is not possible to uniformly attenuate |G(jω)|. If we attenuate G(jω) for acertain frequency range, it will yield an ampliﬁcation at another range. Due to this eﬀect, one usually introduce
a frequency weighting W (s), and minimize |W (jω)G(jω)|.

6

APPL. COMPUT. MATH., V.11, N.1, 2012

2.2. Lifting, transfer functions, and frequency responses. The H ∞ norm criterion isnaturally extended to sampled-data systems. The problem here is that such systems have twotime sets: continuous and discrete. Hence the overall system is not time-invariant in the classicalsense. This diﬃculty can be remedied by the now-standard technique called lifting [3, 4, 33, 34],
which converts a linear time-invariant continuous-time system to an inﬁnite-dimensional discretetime system. It is then possible to naturally extend the notion of the H ∞ -norm to sampled-datasystems. To be more precise, let G denote the input/output operator of such a system. Thenits H ∞ -norm is deﬁned to be the induced norm against all L2 inputs:G∞ :=Gu2.
u∈L2 ,u=0 u2
sup

(6)

Via lifting, this norm is indeed equivalent to the maximum gain of the frequency response
operator of G as in (3).
We start by placing a continuous-time signal in a discrete-time framework. Take a continuoustime signal w(t), and consider the following mapping L (with a suitable domain and codomain)that maps w into a sequence of functions as(Lw)[k] := w[k] := {w(kh + θ)}θ∈[0,h), k = 0, 1, 2, . . . .

(7)

See Fig. 1. The operator L is called lifting.

w(t)

0
0

h

2h

3h

h 0
0

4h

h 0
1

h0
2

h
3

k

Figure 1. Lifting: a continuous-time signal w(t) (left)
is converted to a function-valued discrete time signal (right)

This idea makes it possible to view time-invariant, or even periodically time-varying continuoustime systems as linear, time-invariant discrete-time systems.
Using this operator, one can describe a linear, time-invariant continuous-time system with alinear, time-invariant discrete-time system. Consider the following linear system:ẋ(t) = Ax(t) + Bu(t),y(t) = Cx(t),(8)

where x(t) ∈ Rn , u(t) ∈ Rm and y(t) ∈ Rp are the state, input and output of this system,
respectively. Let us assume, e.g., u ∈ L2loc [0, ∞), the set of locally square-integrable functions
on [0, ∞). The idea is that we view the continuous-time system (8) as one with discrete-timingt = kh, k = 0, 1, 2, . . . such that it receives function-valued inputs at these instants and produces
function-valued outputs at these times also. Suppose that (8) is at state x(kh) at time t = kh.
ThenAhx(kh + h) = e hx(kh) +eA(h−τ ) Bu(kh + τ )dτ,0Aθy(kh + θ) = Ce θx(kh) +0eA(θ−τ ) Bu(kh + τ )dτ.

Y. YAMAMOTO et al.: A BRIEF OVERVIEW OF SIGNAL RECONSTRUCTION ...

7

where 0 ≤ θ < h denotes the intersample parameter. Lifting the input u(t) and the output y(t)as per (7), we can rewrite these formulas as a lifted discrete-time system [4, 34]:x[k + 1] = Ax[k] + Bu[k],y[k] = Cx[k] + Du[k],k = 0, 1, 2, . . . ,where x[k] = x(kh), u[k] = (Lu)[k], y[k] = (Ly)[k], andA : Rn → Rn : x → eAh x, heA(h−τ ) Bu(τ )dτ
B : L2 [0, h) → Rn : u →0

C : Rn → L2 [0, h) : x → CeAθ x, θ
2
2
CeA(θ−τ ) Bu(τ )dτ,
D : L [0, h) → L [0, h) : u →

(9)0

where θ ∈ [0, h) describes the intersample parameter. Observe that the operators A, B, C, D
above do not depend on time k, and hence system (9) is a time-invariant discrete-time system,
albeit with inﬁnite-dimensional input and output spaces. Hence it is straightforward to connect
this system with a discrete-time controller (or a ﬁlter), and the obtained sampled-data system
is again a linear, time-invariant discrete-time system without sacriﬁcing any intersampling
information. The resulting system can also be described by a 4-tuple of operators A, B, C, D,
and its transfer function (operator) of the lifted system is deﬁned as
G(z) = D + C(zI − A)−1 B
with such A, B, C, D. Note that for each ﬁxed z ∈ C \ σ(eAh ), (σ(eAh ):= the set of eigenvaluesof eAh ), G(z) is a linear operator acting on L2 [0, h) into itself. The frequency response operator
is then deﬁned as G(ejωh ), and the gain at frequency ω is deﬁned as  jωh 
G e
v
jωh
.
G(e ) = supv
v∈L2 [0,h)v=0The H ∞ norm of G then becomesG∞ =

G(ejωh ),supω∈[0,2π/h)

which is known to be identical to the L2 -induced norm given by (6) [4].
3. Signal reconstruction and sampling theorem
Consider the block diagram depicted in Fig. 2.
wc

F

yc

yd

K(z)

yK,d

Φ

y

Figure 2. Signal Reconstruction System

In this diagram, the signal wc ∈ L2 denotes the external analog signal to be reconstructed.
It is ﬁltered by an analog ﬁlter (acquisition device) F , and then sampled by the sampler withsampling period h. If f (t) denotes the impulse response of the analog ﬁlter F , then the discretetime signal yd [k] is easily seen to be given by

8

APPL. COMPUT. MATH., V.11, N.1, 2012

yd [k] = (f ∗ wc )[k] = fˇ(· − kh)wc ,(10)
ˇ
where f (t) := f (−t) is the mirror image (with respect to time) of f , and f, g denotes the innerproduct in L2 . The obtained signal yd is then processed by a discrete-time ﬁlter K and thenthe ﬁltered discrete-time signal yK,d is converted back to an analog signal y via a reconstruction
device Φ. Denoting by φ the impulse response of Φ, the reconstructed y is given by∞

yK,d[k]φ(t − kh).

y(t) =(11)k=−∞

In the Shannon paradigm, the analog ﬁlter F is taken to be the ideal ﬁlter, and φ above isthe sinc function [44, 30]. As mentioned in the Introduction, this has several limitations. To
take care of this, one often employs an approximation of the ideal ﬁlter with respect to H 2 norm[10], and this unfortunately yields a sharp ringing eﬀect in the frequency domain.
Unser and co-workers published series of papers of generalized sampling theorems where theacquisition device F is not the ideal ﬁlter [28, 29, 30]. First deﬁne the subspace
∞

Vf :=

α[k]f (t − kh) : {α[k]} ∈ 2

(12)

k=−∞

generated by the translates of the impulse response of the acquisition ﬁlter, and the reconstruction space
∞

Vφ :=

β[k]φ(t − kh) : {β[k]} ∈ 2

(13)

k=−∞

generated by the translates of the reconstruction function φ. From the consistency requirement[28], a key step in their procedure is the oblique projection of L2 onto Vφ perpendicular to Vf .
A precise comparison of this approach with our work is given in [42].
4. H ∞ signal reconstruction problem
We are now ready to precisely state our signal reconstruction problem. The basic features
are the following:• We allow a ﬁnite step preview for reconstruction.
• The acquisition device, sampling and hold elements are ﬁxed.
Consider the block diagram Fig. 3.
yc (· − mh)

e−mhs
wc

F (s)

yc

ec
−

yK,d
Fa (s)

yd

Sh

↑L

K(z)

Hh/L

P (s)

Figure 3. Error system of a sampled-data design ﬁlter

The external continuous-time signal wc ∈ L2 is ﬁrst ﬁltered, or band-limited (mildly but not
perfectly) by going through the analog low-pass ﬁlter F (s), which is linear and time-invariant,and ﬁnite-dimensional. This F (s) is a rational function of s which is strictly proper (i.e., thedegree of the numerator polynomial is less than that of the denominator). As is well known, it
is represented by a linear, time-invariant system

Y. YAMAMOTO et al.: A BRIEF OVERVIEW OF SIGNAL RECONSTRUCTION ...

9dx(t) = Ax(t) + Bu(t),dty(t) = Cx(t),where A, B, C are constant matrices of appropriate sizes, and F (s) = C(sI − A)−1 B. Hence F is
not an ideal ﬁlter unlike the case of the Shannon paradigm, and is physically realizable throughthe above state space model. The signal wc is the external signal that drives F and produces
the actual signal yc to be processed. That is, we assume that the original analog signals to besampled are in the following subspace of L2 :
F L2 := yc ∈ L2 : yc = F wc , wc ∈ L2 .
It is proved in [24] that the band-limited signal subspace
BL := yc ∈ L2 : supp ŷc ⊂ (−π/h, π/h) ,is a proper subset of F L2 , that is, BL  F L2 . The ﬁlter F (s) is chosen based on the following
guidelines:• a frequency distribution of input analog signals obtained by averaging or enclosing gainsof their Fourier transforms.
• a dynamical model of signal generator such as musical instruments.
The example in Section 6 gives a brief guideline on how to choose F (s) based on the envelopeof energy distributions of the signal. Note that when F is ideal, then the class we are dealing
with agrees with the ideal sampling theorem.
The produced signal yc is then sampled by ideal sampler Sh , ﬁltered by an anti-aliasingﬁlter Fa (s), and becomes a discrete-time signal yd with sampling period h. This signal is thenupsampled by ↑ L to allow for processing (interpolation) between the original sampling period
h. The digital ﬁlter K(z) processes this upsampled signal to produce yK,d. The signal yK,d thengoes through the zero-order hold Hh/L and becomes a continuous-time signal. It is then further
processed by an analog low-pass ﬁlter P (s) to become the ﬁnal analog output yp .
In the upper part of the diagram, we allow m steps of delay for the analog signal yc andobtain yc (t − mh). This is a setup for allowing a “preview” of yc for m samples by the proper
ﬁlter transfer function K(z). It is very eﬀective compared to reconstruction without a preview.
This also takes care of certain processing delays caused by the processing ﬁlter. The integer m
is a design parameter that can be chosen by the designer. This is in marked contrast to theconventional design methodologies: These methods usually allow a non-causal impulse response
for reconstruction, e.g., [32, 21, 22]. But in real implementation, one has to truncate it, and it
is often unclear how many steps one would need to obtain a desired accuracy. In the present
setup, one can prespecify an allowable step of delays (preview), and obtain an optimal design
under such a constraint.
Finally, the processed signal yp is compared with this delayed yc (t − mh) and subtracted fromit to obtain the error signal ec . The design objective is to make the error as small as possible.
Observe also that this design framework is formulated in the continuous-time domain in contrast
to the usual discrete-time setups.
We must specify a performance index to give a precise meaning to this problem. The following
L2 induced norm from wc to ec (or the sampled-data H ∞ norm) is the one we take:J :=ec 2.
wc 2
wc ∈L2 ,wc =0
sup

We thus arrive at the following design problem:

(14)

10

APPL. COMPUT. MATH., V.11, N.1, 2012

Problem 1. Let Tew denote the input/output operator from wc to ec (·) := yc (·) − uc (· − mh)in Fig. 3. Given an attenuation level γ > 0, ﬁnd, if one exists, a stable digital (discrete-time)ﬁlter K(z) such that
J := Tew ∞ =supwc ∈L2 [0,∞)Tew wc 2< γ.
wc 2(15)

The performance index (15) intends to minimize the maximum error induced by an (unknown)
input wc that gives rise to the largest norm of ec among all inputs. This is made possible by the
H ∞ design methodology. Note that the actual error is not known, but due to the min-max natureof the problem, we can minimize the worst transmission error. Observe also that this setup allows
for a capability of minimizing continuous-time phase errors due to the continuous-time natureof the performance index, as opposed to the conventional gain-phase design principles.
This min-max problem diﬀers sharply from the orthogonal projection based methods. Also,
due to sampling, Tew is not even time-invariant (in continuous-time).
It is now known however that this problem is reducible to a linear time-invariant problem via
lifting; see Subsection 2.2; the problem is now solvable via now-standard H ∞ control theory, see,e.g., [4, 3] (see also [8] for standard treatments of H ∞ control in the continuous-time setting).
The existence of e−mhs makes this an inﬁnite-dimensional H ∞ problem; see [19, 35, 39, 41],
etc. The simplest solution is to employ the so-called fast-sampling/fast-hold approximation,
which we will outline in the next section.

5. Solution via fast-sampling/fast-hold approximation
While Problem 1 is known to be reducible to a ﬁnite-dimensional problem [19, 23], it is notnecessarily appealing computationally. It is often more convenient to resort to an approximation
method. We employ the fast sample/hold approximation [18, 4, 35, 39]. This method approximates continuous-time inputs and outputs via a sampler and hold that operate in the periodh/N for some positive integer N . The method usually works fairly well for N ∼ 5L, where
L is the upsampling ratio given in Section 4, and the convergence of such an approximationis shown in [35, 39]. We show here the design procedure of K(z) by the fast sampling/holdapproximation.
The error system in Fig. 3 is a multirate system due to the upsampler ↑ L. We ﬁrst reduce
this system to a single-rate one. Introduce the discrete-time lifting, also known as the polyphasedecomposition [43], LL and its inverse L−1
L as

T
,
LL := (↓ L) 1 z · · · z L−1
1 z −1 · · · z −L+1 (↑ L).
L−1
L :=

(16)

Then K(z)(↑ L) can be rewritten by a lifted system asK(z)(↑ L) = L−1
L K(z),K(z):= LL K(z)L−1L

1 0 ···

0

T

.


The ﬁlter K(z)is an LTI (linear and time-invariant), single-input/L-output system. Deﬁne−1Hh := Hh/L LL , and we obtain the following equality
 h K(z)S
Hh/L K(z)(↑ L)Sh = H
h.
Hence the multirate system in Fig. 3 is reduced to the single-rate system shown in Fig. 4.

Y. YAMAMOTO et al.: A BRIEF OVERVIEW OF SIGNAL RECONSTRUCTION ...

11

e−mhs
wc

yc

F (s)

ec
−
Fa (s)

yd

Sh


K(z)

h
H

P (s)

Figure 4. Reduced single-rate error system Tew

We then employ the fast sample/hold approximation for the error system Tew in Fig. 4. We
connect fast sample and hold devices Sh/N , Hh/N , and the discrete-time lifting LN with theerror system Tew as shown in Fig. 5.
w̃d

wd

L−1
N

Hh/N

wc

ec

Tew

Sh/N

ed

LN

ẽd

Figure 5. Fast-sampling/Fast-hold (FSFH) transformation

We call this fast-sampling/fast-hold (FSFH) transformation and denote it by FSFH(Tew , h, N ).
Before we give the design formula of FSFH(Tew , h, N ), we introduce FSFH transformation forcontinuous-time LTI systems. For brevity of notation, let us adopt the following shorthand notation for continuous-time transfer function D + C(sI − A)−1 B or discrete-time transfer function
D + C(zI − A)−1 B:


A B
C D


:=

D + C(sI − A)−1 B, for continuous-time systems,
D + C(zI − A)−1 B, for discrete-time systems.

Let c2d denote the step-invariant transformation [4], that is,

 




h
A B
A B
eAh 0 eAt Bdt.
, h := Sh
Hh =
c2d
C D
C D
C
D
Also, let lift denote the discrete-time lifting transformation [4], that is,

lift

A B
C D




,N



A B
C D



:= LN
L−1
N =
⎤
⎡
AN −1 BAN −2 B . . . B
AN
⎢
C
D
0
... 0 ⎥
⎥
⎢
⎢
.. ⎥
..
⎢ CA
. . ⎥
CB
D
⎥.
⎢
⎢
..
..
..
.
.. 0 ⎥
⎦
⎣
.
.
.
N
−1
N
−2
N
−3
CA
B CA
B ... DCA

Then, for continuous-time LTI system F , the FSFH transformation of F is given by
FSFH(F, h, N ) := lift c2d(F, h/N ), N .
By the FSFH transformation, the sampled-data error system Tew can be approximated by adiscrete-time LTI system as in the following theorem:
Theorem 1. Let N = Ll, where l is a positive integer, and deﬁne the discrete-time LTI
system TN as follows:FN (z),
TN (z) = z −m FN (z) − PN (z)H K(z)S
where

(17)

12

APPL. COMPUT. MATH., V.11, N.1, 2012

FN := FSFH(F, h, N ), PN := FSFH(P, h, N ), FN := FSFH(F Fa , h, N ),H := diag {Il } ∈ RN ×L , Il := [1, 1, . . . , 1]T ∈ Rl , S := [1, 0, . . . , 0] ∈ R1×N .
 and for each ω ∈ [0, 2π/h), the frequency response
Then, for each ﬁxed K
TN (ejωh ) → Tew (ejωh ),(18)

as N → ∞, and this convergence is uniform with respect to ω ∈ [0, 2π/h). Furthermore, this if K ranges over a compact set of ﬁlters.
convergence is also uniform in K
The proof is almost the same as in [42, Theorem 1].
 our design problem (15) can be apIn view of the uniformity of convergence TN ∞ in K,
proximated byTN ∞ < γ.
satisfying the above
This is a discrete-time H ∞ optimization problem. To obtain a ﬁlter K(z)inequality, we can adopt numerical softwares as MATLAB with robust control toolbox [20], by
the generalized plant representation depicted in Fig. 6, where w̃d = LN wd and ẽd = LN ed .
w̃dz −m FN (z)S FN (z)ud

−PN (z)H
0



ẽd
yd


K(z)
Figure 6. Discrete-time system for H ∞ optimization

Once the optimal ﬁlter K(z)is obtained, one can obtain the interpolation ﬁlter K(z) by thefollowing formula:K(z) =

1 z −1 · · ·z −L+1

 L ).
K(z6. Design examples
Let us make a comparison with a usual linear phase ﬁlter—the equi-ripple FIR ﬁlter obtained
by Parks-McClellan method [43, 45] with 64 taps. Parks-McClellan method is widely used fordesigning FIR ﬁlters for interpolation. We design the proposed ﬁlter K(z) with interpolation
ratio L = 4, sampling period h = 1, and delay step m = 4. The analog ﬁlters F (s), Fa (s), andP (s) are given by1, T = 7.0187,(T s + 1)(0.1T s + 1)
1
1
, P (s) =.
Fa (s) =0.01s + 1
0.05s + 1
Reﬂecting a typical energy distribution of orchestral music, the time constant T = 7.0187 istaken to be equivalent to 1 kHz with sampling frequency 44.1 kHz. It therefore corresponds toan energy distribution that decays by −20 dB per decade from 1 kHz and −40 dB per decadefrom 10 kHz.
F (s) =

Y. YAMAMOTO et al.: A BRIEF OVERVIEW OF SIGNAL RECONSTRUCTION ...

13

20

0

Gain (dB)

−20

−40

−60

−80

−100
0

0.5

1

1.5
2
Freqeuncy (rad/sec)

2.5

3Figure 7. Bode gain plot of the proposed ﬁlter (solid) and the equi-ripple ﬁlter (dash)Fig. 7 shows the Bode gain plots of the proposed ﬁlter and the equi-ripple FIR ﬁlter with
64 taps. We can see that the equi-ripple ﬁlter has a sharp decay around the cutoﬀ frequency
ω = π/4, while the ﬁlter obtained by the proposed method shows a rather mild decay.

1.5

1Amplitude

0.5

0

−0.5

−1

−1.5
0

5

Time (sec)10

15

Figure 8. Response of the equi-ripple ﬁlter against a rectangular wave

Fig. 8 shows the response of the equi-ripple ﬁlter against a rectangular wave. It exhibits a very
sharp ringing eﬀect. This is because the ﬁlter has a sharp cut-oﬀ characteristic, and inevitably
introduces the well-known Gibbs phenomenon due to the fact that the frequency components
beyond the pass-band are sharply truncated. In contrast, Fig. 9 shows the response of the ﬁlter
designed by the present method. It shows virtually no ringing.

14

APPL. COMPUT. MATH., V.11, N.1, 2012

1.5

1

Amplitude

0.5

0

−0.5

−1

−1.5
0

5

Time (sec)10

15

Figure 9. Response of a sampled-data design ﬁlter against a rectangular wave

Fig. 10 shows the frequency response gain of the sampled-data error system Tew . We designﬁlters by the proposed method for reconstruction delay m = 1, 2, 3, 4, . . . , and compare themwith the 64-tap equi-ripple ﬁlter. For all m, the H ∞ norm of Tew does not change, but thegain reduces for almost all frequencies as m increases. For m ≥ 4, the frequency response gain
remains unchanged. This indicates that m = 4 is optimal from the point of view that smaller
delay is better for signal processing. The equi-ripple ﬁlter exhibits large errors in the whole
frequency range as compared with the proposed ﬁlter when m ≥ 4. Note that the equi-ripple
ﬁlter has the reconstruction delay m = 16. These errors give an explanation of the ringing eﬀect
in Fig. 8.

conventional
proposed (m=1)
proposed (m=2)
proposed (m=3)
proposed (m≥4)

−34
−35

Gain (dB)

−36
−37
−38
−39
−40
−41
−2

10

−1

10
Frequency (rad/sec)

10

0

Figure 10. Frequency response of error system Tew with sampled-data designed
ﬁlter (solid) and the equi-ripple ﬁlter (dash)Fig. 11 shows the relation between the upsampling ratio L and the achieved performanceTew ∞ . As illustrated in this ﬁgure, the performance improves as the upsampling ratio L

Y. YAMAMOTO et al.: A BRIEF OVERVIEW OF SIGNAL RECONSTRUCTION ...

15

increases, but for L ≥ 4, the performance is nearly unchanged. This ﬁgure suggests that L = 4
is nearly optimal since for large L leads to high sampling frequency Lπ (rad/sec).

0.04

0.035

||Tew||

∞

0.03

0.025

0.02

0.015
1

2

3

4

5

L

6

7

8

9

10

Figure 11. Upsampling ratio L versus performance Tew ∞

7. Concluding remarks
We have presented an overview on a new framework for digital signal processing. The fundamental philosophy is the emphasis on analog (continuous-time) performance with discrete-time
signal processing. This naturally leads to a technical diﬃculty because of the two diﬀerent timesets involved: continuous and discrete. Leveraging sampled-data H ∞ control theory, we havepresented computable procedures for designing optimal, stable, causal ﬁlters. These ﬁlters are
optimal with respect to a uniform analog performance measure. Our methodology is applicable
to a wide variety of theoretical and application problems in digital signal processing. We hopethat this article helps the readers familiarize themselves to this new approach and that it will
be be more widely used in the future.
Acknowledgments
This research is supported in part by the JSPS Grant-in-Aid for Scientiﬁc Research (B) No.
2136020318360203, Grant-in-Aid for Exploratory Research No. 22656095, and the MEXT Grantin-Aid for Young Scientists (B) No. 22760317. Pramod Khargonekar is also supported in partby Eckis Professor endowment.

References[1] Ashida, S., Nagahara, M., Yamamoto, Y. Audio signal compression via sampled-data control theory, Proc.
SICE Annual Conference, 2003, pp.1182-1185.
[2] Ashida, S., Kakemizu, H., Nagahara, M., Yamamoto, Y. Sampled-data audio signal compression with Huﬀman coding, Proc. SICE Annual Conference, 2004, pp.972-976.
[3] Bamieh, B., Pearson, J.B. A general framework for linear periodic systems with applications to H∞ sampleddata control, IEEE Trans. Autom. Control, AC-37, 1992, pp.418-435[4] Chen, T., Francis, B.A. Optimal Sampled-Data Control Systems, Springer, 1995.
[5] Chen, T., Francis, B.A. Design of multirate ﬁlter banks by H∞ optimization, IEEE Trans. Signal Processing,SP-43, 1995, pp.2822-2830.

16

APPL. COMPUT. MATH., V.11, N.1, 2012

[6] Condat, L., Blu, T., Unser, M. Beyond interpolation: optimal reconstruction by quasi-interpolation, IEEE
International Conference on Image Processing, Vol.1, 2005, pp.33-36.
[7] Doyle, J.C., Francis, B.A., Tannenbaum, A.R. Feedback Control Theory, Macmillan, New York, 1992.
[8] Doyle, J.C., Glover, K., Khargonekar, P.P., Francis, B.A. State-space solutions to standard H∞ and H2
control problems, IEEE Trans. Autom. Contr., AC-34: 1989, pp.831-847[9] Eldar, Y.C., Dvorkind, T.G. A minimum squared-error framework for generalized sampling, IEEE Trans.
Signal Processing, SP-54, N.6, 2006, pp.2155-2167.
[10] Fliege, N.J. Multirate Digital Signal Processing, John Wiley, New York, 1994.
[11] Francis, B.A. A course in h∞ control theory, Springer, 1987.
[12] Fujiyama, K., Iwasaki, N., Hirasawa Y., Yamamoto, Y. High frequency compensator and reproducing device,
US patent No. 07324024B2, 2008.
[13] Fujiyama, K., Iwasaki, N., Hirasawa, Y., Yamamoto, Y. High frequency compensator and reproducing device,
Chinese patent No. 648701, 2010.
[14] Hassibi, B., Erdogan, A. T., Kailath, T. MIMO linear equalization with an H ∞ criterion, IEEE Trans. Signal
Processing, SP-54, N.2, , 2006, pp.499-511.
[15] Ishii, H., Yamamoto, Y., Francis, B.A. Sample-rate conversion via sampled-data H ∞ control, Proc. 38th
IEEE CDC, 1999, pp.3440-3445[16] Kakemizu, H., Nagahara, M., Kobayashi, A., Yamamoto, Y. Noise reduction of JPEG images by sampleddata H ∞ optimal epsilon ﬁlters, Proc. SICE Annual Conference, 2005, pp.1080-1085.
[17] Kashima, K., Yamamoto, Y., Nagahara, M. Optimal wavelet expansion via sampled-data control theory,
IEEE Signal Processing Lett., V.11. 2004, pp.79-82.
[18] Keller, J.P., Anderson, B.D.O. A new approach to the discretization of continuous-time systems, IEEE Trans.
Autom. Control, AC-37, N.2, 1992, pp.214-223.
[19] Khargonekar, P.P., Yamamoto, Y. Delayed signal reconstruction using sampled-data control, Proc. 35th IEEE
CDC, 1996, pp.1259-1263.
[20] Mathworks, robust control toolbox, http://www.mathworks.com/products/robust/[21] Meinsma, G., Mirkin, L. Sampling from a system-theoretic viewpoint: Part I-concepts and tools, IEEE Trans.
Signal Processing, SP-58, N.7, 2010, pp.578-590.
[22] Meinsma, G., Mirkin, L. Sampling from a system-theoretic viewpoint: Part II-noncausal solutions, IEEE
Trans. Signal Processing, SP-58, N.7, 2010, pp.3591-3606.
[23] Mirkin, L., Tadmor, G. Yet another H ∞ discretization, IEEE Trans. Autom. Contr., AC-48, 2003, pp.891894.
[24] Nagahara, M. Ogura, M., Yamamoto, Y. H ∞ design of periodically nonuniform interpolation and decimation
for non-band-limited signals, SICE Journal of Control, Measurement, and System Integration, V.4-5, 2011,
pp.341-348.
[25] Nagahara, M., Yamamoto,Y. A new design for sample-rate converters, Proc. 39th IEEE CDC2000, Sydney,
Australia, 2000, pp.4296-4301.
[26] Nagahara, M., Yamamoto, Y. Optimal design of fractional delay FIR ﬁlters without band-limiting assumption, Proc. IEEE International Conference on Acoustics, Speech, and Signal Processing, V.4, 2005, pp.221224.
[27] Shannon, C.E. Communication in the presence of noise, Proc. IRE, 1949, pp.10-21.
[28] Unser, M., Aldroubi, A. A general sampling theory for nonideal acquisition devices, IEEE Trans. on Signal
Processing, SP-42, N.11, pp.2915-2925, 1994.
[29] Unser, M., Zerubia, J. Generalized sampling: stability and performance analysis, IEEE Trans. Signal Processing, SP-45, N.12, 1997, pp.2941-2950.
[30] Unser, M. Sampling - 50 years after Shannon, Proc. of the IEEE, V.88-4, 2000.
[31] Unser, M., Blu, T. Cardinal exponential splines: Part I - theory and ﬁltering algorithms, IEEE Trans. Signal
Processing, SP-53, N.4, 2005, pp.1425-1438.
[32] Unser, M. Cardinal exponential splines: Part II – Think analog, act digital, IEEE Trans. Signal Processing,
SP-53, N.4, pp.1439-1449, 2005.
[33] Yamamoto, Y.A function space approach to sampled-data control systems and tracking problems, IEEE
Trans. Autom. Control, AC-39, 1994, pp.703-713.
[34] Yamamoto, Y. Digital control, Wiley Encyclopedia of Electrical and Electronics Engineering, V.5, 1999,
pp.445-457.
[35] Yamamoto, Y., Madievski, A.G., Anderson, B.D.O. Approximation of frequency response for sampled-data
control systems, Automatica, V.35-4, 1999, pp.729-734.
[36] Yamamoto,Y. A new approach to signal processing via sampled-data control theory, Australian J. Elec. &Electronics Eng., V.2, 2005, pp.141-148.
[37] Yamamoto, Y. Digital analog converters and a design method for the pertinent ﬁlters, Japanese patent,
N.3820331, 2006.

Y. YAMAMOTO et al.: A BRIEF OVERVIEW OF SIGNAL RECONSTRUCTION ...

17[38] Yamamoto, Y., Nagahara, M. Sample-rate converters, Japanese patent, N.3851757, 2006.
[39] Yamamoto, Y., Anderson B.D.O., Nagahara, M. Approximating sampled-data systems with applications todigital redesign, Proc. 41st IEEE CDC, 2002, pp.3724-3729.
[40] Yamamoto, Y., Fujioka, H., Khargonekar, P.P. Signal reconstruction via sampled-data control with multirate
ﬁlter banks, Proc. 36th IEEE CDC. 1997, pp.3395-3400.
[41] Yamamoto, Y., Nagahara, M., Fujioka, H. Multirate signal reconstruction and ﬁlter design via sampled-data
H ∞ control, Proc. MTNS 2000, Perpignan, France, 2000.
[42] Yamamoto,Y., Nagahara, M., Khargonekar, P.P. Signal reconstruction via H ∞ sampled-data control theoryBeyond the Shannon paradigm, IEEE Trans. Signal Processing (to be published).
[43] Vaidyanathan, P.P. Multirate Systems and Filter Banks, Prentice Hall, Englewood Cliﬀs, 1993.
[44] Zayed, A.I. Advances in Shannon’s Sampling Theory, Boca Raton, FL, CRC Press, 1996.
[45] Zelniker, G., Taylor, F.J. Advanced Digital Signal Processing, Marcel Dekker, 1994.

Yutaka Yamamoto - received the Ph.D. degree inMathematics from the University of Florida, Gainesville,in 1978, under the guidance of Professor Rudolf E.
Kalman. He is currently a Professor in the Department
of Applied Analysis and Complex Dynamical Systems,
Graduate School of Informatics, Kyoto University, Kyoto, Japan. His current research interests include thetheory of sampled data control systems, its application
to digital signal processing, realization and robust control
of distributed parameter systems and repetitive control.
Dr. Yamamoto was the recipient of G. S. Axelby Outstanding Paper Award of the IEEE Control Systems
Society in 1996, the Commendation for Science and Technology by the Minister of Education, Culture, Sports,
Science and Technology, Prize for Science and Technology, in Research Category in 2007, Distinguished MemberAward of the Control Systems Society of the IEEE in 2009, and several other awards from the Society of Instrument
and Control Engineers (SICE) and the Institute of Systems, Control and Information Engineers (ISCIE). He hasserved as Senior Editor of the Trans. Automatic Control (TAC) for 2010-2011, and as an AE for several journals
including TAC, Automatica, MCSS, SCL. He was a chair of the Steering Committee of MTNS for 2006-2008, and
the General Chair of the MTNS 2006 in Kyoto.
He is currently President-Elect of the Control Systems Society (CSS) of the IEEE, and was a Vice President
for 2005-2008 of CSS, and President of ISCIE of Japan for 2008-2009. He is a fellow of the IEEE and SICE,
Japan.

Masaaki Nagahara - received the Bachelor’s degree in Engineering from Kobe University in 1998,
the Master’s degree and the Doctoral degree in informatics from Kyoto University in 2000 and 2003. He
is currently an Assistant Professor at the Graduate
School of Informatics, Kyoto University. His researchinterests include digital signal processing and digital
control systems. He was a recipient of Outstanding
Paper Awards from the Society of Instrument andControl Engineers (SICE) in 1999. He is a memberof IEEE, ISCIE, SICE, and IEICE.

18

APPL. COMPUT. MATH., V.11, N.1, 2012

Pramod P. Khargonekar - received B. Tech. Degree in Electrical Engineering from the Indian Instituteof Technology, Bombay, India and M.S. degree in Mathematics and Ph.D. degree in Electrical Engineering fromthe University of Florida. He has held faculty positions at
the University of Minnesota and The University of Michigan. He was Chairman of the Department of Electrical Engineering and Computer Science at Michigan from1997 to 2001, and also held the title Claude E. Shannon
Professor of Engineering Science there.
From 2001 to 2009, he was Dean of the College of Engineering and is now Eckis Professor Electrical andComputer Engineering at the University of Florida. Khargonekar’s current research interests are focused on
renewable energy and electric grid, neural engineering, and systems theory.
Dr. Khargonekar is a recipient of the NSF Presidential Young Investigator Award, the American AutomaticControl Council’s Donald Eckman Award, the Japan Society for Promotion of Science Fellowship(s), and aDistinguished Alumnus Award from the Indian Institute of Technology, Bombay. He is a co-recipient of the IEEE
W. R. G. Baker Prize Award, the IEEE CSS George S. Axelby Best Paper Award, and the AACC Hugo Schuck
Best Paper Award. He was a Springer Professor the University of California, Berkeley in 2010. He is a Fellowof IEEE and is on the list of Highly Cited Researchers from Thomson-Reuters (ISI Web of Knowledge). At the
University of Michigan, he received a teaching excellence award from the EECS department, a research excellence
award from the College of Engineering, and the Arthur F. Thurnau Professorship. At the University of Minnesota,
he received the George Taylor Distinguished Research Award.