64-4-2
不等間隔にサンプルされた信号からの H ∞ 最適復元
○和田敏裕 永原正章 山本裕 （京都大学）

Digital Signal Reconstruction from Nonuniform Samples via H ∞ Optimization∗T. Wada, M. Nagahara and Y. Yamamoto (Kyoto University)
Abstract– This paper considers a signal reconstruction problem from non-uniformly sampled signals. The
reconstruction system is constructed of two systems. One is a resampler which converts non-uniform samples
to uniform ones. The other is a digital filter which is aimed at reconstructing the original analog signals,and reducing noise caused by the non-uniform sampler. We design these systems by sampled-data H-infinityoptimization. A design example shows the result.
Key Words: Nonuniform sampling, Signal reconstruction, Sampled-data H ∞ control

1

はじめに離散時間システムの設計問題の応用として有用なものの 1 つに, サンプリング周期の変換が挙げられる. 近
年のデジタル信号処理の発達でその応用も多岐にわたり, 場合によってはリアルタイム性が求められたり, 巨大なデータ量を持つ離散時間信号を短時間のうちに処理しなければならないこともある. これに対して, 処理時の計算量を軽減するためにあらかじめオフラインで変換器を設計しておくことが有効である. 特に音声や画像といったメディアを表す信号の場合は, あらかじめ信号の周波数モデルを仮定することができるため, より良い変換器が設計できる可能性がある.
そのようなアプローチの 1 つとして, 変換前と変換後のサンプリング周期が整数比になっているものに関しては, サンプル値 H ∞ 制御の枠組みを用いた設計法が既に提案されている 8, 9) . しかし, サンプリングに関してより複雑な状況を呈する問題も多くある. 本研究では, そのうち不等間隔にサンプルされた離散時間信号を等間隔にリサンプルする問題を考える.
この種のリサンプル処理については, これまでも様々な手法が提案されてきた. しかしそれらは, 元の曲線に完全帯域制限 5, 10, 11) や区分的多項式 16) といった強い仮定を置いている. さもなくば, サンプリング自体に周期性などの条件を置いている 1, 17) . また, そういった制約が無い手法の場合は, 信号の長さに対して非常に大きな処理時間を必要とする 3) .
本稿では, 不等間隔にサンプルされた信号のリサンプリングに関して, 代数演算と離散時間フィルタを用いた処理系を提案し, その設計法を示す. この処理系は前述のような強い仮定を必要とせず, また信号長に対して処理時間が比較的短いという特徴を持つ.

2

不等間隔なサンプリング

定義 2 (不等間隔サンプラ) サンプル点集合 T に対して, 不等間隔サンプラ ST を(ST uc )[k] := (uc (tk ), tk ),tk ∈ T ,k∈Zと定義する.
不等間隔にサンプルされた信号は, サンプリングの時刻とその時刻での信号のサンプル値の両方の情報を持つ.
本手法では, リサンプリングの処理系 R を次のように構成する.

R = KH.
ただし, K は連続時間信号から離散時間信号への線形システムであり, H は次で定義するホールドである.
定義 3 (不等間隔ゼロ次ホールド) 不等間隔サンプラ
ST の出力 ug [k] = (uk , tk ) に対して，システム H を(Hug [k])(t) := uk , ∀t ∈ [tk , tk+1 )と定義する. これを不等間隔ゼロ次ホールドと呼ぶ.
つぎに, 考慮する不等間隔信号を限定する. まず, 不等間隔サンプリングされる前の連続時間信号 ws を考える.
仮定 4 信号 ws の周波数モデルを, 厳密にプロパーな有理伝達関数 W1 とする.
これを理想サンプラでサンプルし, ゼロ次ホールドした
信号 wr0 と, T で不等間隔サンプルし, 不等間隔ゼロ次ホールドした信号 wr との差を wn とおく．すなわち，

wn := wr0 − wr = (HhZO Shidl − HST )ws .

不等間隔なサンプリングを扱うために, いくつかの概念を導入する.

信号 wr0 , wr および誤差 wn の例をそれぞれ, Fig. 1,Fig. 2 および Fig. 3 に示す．この誤差をノイズとみなし, 次のような仮定を置く.

定義 1 (サンプル点集合) 任意の k ∈ Z に対し時刻tk ∈ R (ただし tk < tk+1 ) が与えられたとき, 集合

仮定 5 信号 wn の周波数モデルを, プロパーな有理伝達関数 W2 とする.

T := {tk ∈ R; k ∈ Z}
をサンプル点集合と呼ぶ.

第7回制御部門大会（2007年3月6日〜8日・調布）

仮定 4, 仮定 5 の下で, K の設計は文献 9) のような,ノイズ付加時の信号復元問題とみなすことができ, 次のように定式化できる.

SY0003/07/0000-6442 (c) 2007 SICE

Aw1`
I

¸ ·Aw1`
I,0Cw10

Wj = Lh Wj L−1
h ,(j = 1, 2),-1

idl −1Sidl
h := Sh Lh ,·Ŵp (z) :=

1

00.1

0.2

0.3

0.4 0.5 0.6
time [sec]

0.7

0.8

0.9

110-1
0.1

0.2

0.3

0.4 0.5 0.6 0.7 0.8 0.9time [sec]Fig. 2: ws (実線) と wr0 (破線) の例

10.1

10-1
0.2

0.3

0.4 0.5 0.6 0.7
time [sec]Fig. 3: wr − wr0 の例

0.80.9問題 6 出力のサンプル周期 h, 許す遅れ N , 及び W1 ,W2 , γr > 0 が与えられたとき,k(DhN − HhZO KHhZO Shidl )W1 rk2< γr ,krk2r∈L2supkHhZO KW2 nk2sup< γnknk2n∈L2の下で γn を最小化する K を求めよ. ただし, Dh は時間 h の遅れ作用素, W1 , W2 は伝達関数 W1 , W2 によってあらわされる入出力作用素とする．

γn に関する最小化問題は, ある γn に関する可解問題の反復によって解くことができる. また, γr , γn はともに W1 , W2 に含めることができるので, 以下では問題 6 として γr = γn = 1 に関する可解問題を考える.

3

主要結果

問題 6 は Periodic な系であり, サンプル値制御理論の枠組みを用いて取り扱うことができる. システムǦu を#"∆Ǧr110∆HhZO
Ǧu :=
Ǧr21 G̃n2104)

とおく. ただし, Lh をリフティング作用素 19) , U を離散時間の 1 ステップの遅れ作用素として,Ǧr11 := −Cw1 U N Wp Bw1 − Dw1 ΠL U N ΠR ,Ǧr21 := HhZO Sidl
h Cw1 Wp Bw1 ,° ZO ° °°−1G̃n21 := °Hh ° °∆HhZO ° W2 ,∆ := (I − Dw1 (I − ΠL ΠR )D∗w1 )−1/2 ,ΠL := Bw1 (Bw1 B∗w1 )† ,

ΠR := B∗w1 ,:= Ŵ1 (z),HhZO := Lh HhZO

とおく.
さて H ∞ 制御問題の正則性に関する仮定 13)
仮定 7 任意の θ ∈ [0, 2π) において·¸Au` − ejθ I0は左可逆. かつ
Cu1
Du12·
¸
Au` − ejθ I Bu1は右可逆.
Cu2Du2110¸とする1 . Ǧu は有限次元の状態空間の実現を持ち,Au` Bu10
ˆ (z) Cu10Du12  := Ǧu
Cu2 Du210Fig. 1: ws (実線) と wr (破線) の例

0Bw1Dw1の下で, 問題 6 が可解ならば K は離散時間システム Kd
及び積分型サンプラ Shint によりK := Kd Shintと構成できる 18) . 仮定 7 は, 個々の θ については文献
12)の手法を用いることにより判定できる.
さて, ここで次のようなサンプラを定義しよう.
step

定義 8 任意の θ ∈ R+ に対して, 作用素 Sθ をZ θ
step22Sθ : L 3 uc 7→ vd ∈ ` , vd [k] =u(kθ +σ)dσ0と定義する.
step
この Sθstep を用いて Shint を Mφ ULL Sh/Lと近似する.
Mφ は定数行列であり, LL は L 倍の離散時間リフティ
ング作用素 4) である. これは積分型サンプラの区分求積法による近似に相当し, 問題 6 が可解ならば十分大きな L に対してこの近似は収束する 18) .
このとき, 離散時間システム K 及び行列 Mφ は, Ke :=KMφ として同時に求めることができる. すなわち, 問題 6 は問題 9 に帰着する.

問題 9 Fig. 4 において,

supr̆∈`2ky̆k2< 1,kr̆k2

supn̆∈`2ky̆k2<1kn̆k2を満たす Ke を求めよ. ただし Fig. 4 中
step −1
Ssteph/L := Sh/L Lh/L ,W2f := Lh/L W2 L−1h/Lである.

Fig. 4 中, 網掛け部分は全体として時不変かつ因果的な系になっていることに注意すると, 文献 14) の手法により等価な有限次元離散時間系を得ることができる.
結果, 問題は次の問題に帰着する.
1 上付き † は Moore-Penrose 擬逆行列を表す. また, G をシステムとして, Ĝ は G に対応する伝達関数行列を表す.

n̆r̆

L−1LW1W2fSstep
h/L

Sidl
hHhZO


0 0
0 Aw
Φb := 
0 00 000>
Bw1 Bw10,>−Aw10
00µ·¸¶An BneĜne (z) := dliftL,Cn Dne

LL

An := Ψc11 ,>Bne Bne:= Ψc11 Ψ0c12 ,UNy̆−Ssteph/L

LL

KeUHhZO

Cne := Ψd21 ,0Dne := Ψd21 Ψd14 + Ψd22 Ψ0d24 ,Ψc := e(h/L)Φc , Ψ0c := e−(h/L)Φc ,·¸>
Aw2 −Bw2 Bw2
Φc :=,0−A>w2
Ψd := e(h/L)Φd , Ψ0d := e−(h/L)Φd ,
>>Aw2 0 −Bw2 Bw2−Bw2 Dw2>
> Cw2 0 −Dw2 Bw2−Dw2 Dw2,Φd := >> 00−Aw2−Cw2 +Fig. 4: 設計問題
問題 10 システム Ge を"
U N −1 Gre110Ge :=Gre21G̃neとして°°°°° F` (Ĝe , K̂e ) °Dre12#·0·<1∞を満たす因果的かつ時不変な線形システム Ke を求めよ. ただし
ArBre1
 Cre1
Dre11 ·¸Ĝre11
 (h/L)Cw10 (z) := ,Ĝre21..
.. .
. (h/L)Cw1−1G̃ne := hkDre12 kGne ,0Ar := Ψb22 ,>Dre11 := Dra Bre1 , Bre1 Bre1:= Ψb23 Ψ>b22 , > Cre1 £¤>  DraCre1 Dra Dre12>
Dre12


Σ11 Σ12 Σ13
Σ22 Σ23  ,:= Σ>12>Σ13 Σ>
Σ3323†
Σ23 := (Ψa13 − Ψb13 Ψ>
b22 Ψa33 )Ψa23 ,Σ33 := Σ23 (Ψb24 − Ψa24 ) + Ψa14 − Ψb13 Ψ>
b22 Ψa34 ,Ψa := ehΦa , Ψb := ehΦb ,0−Cw10>
0AwBw1 Bw1
Φa := >>
0 −Cw1 Cw1−Aw1000I0 ,
Cw1 0

Aw1
Cw1

Bw1
Dw1

Aw2
Cw2

Bw2
Dw2

¸000:= W1 (s),¸:= W2 (s)である. dliftL は伝達関数に対する離散時間リフティング 4) を表す.
問題 10 は, 一般的な離散時間 H ∞ 制御問題であり, 様々な解法が知られている 6) . また, MATLAB 等を用いると数値的に解くことができる 2, 7) . 式中の行列の指数の計算は, 効率の良い計算方法が知られている 15) .

4

サンプラの実装

不等間隔にサンプルされた信号に対する処理系は
step
stepK̃H := Ke ULL Sh/LH であるが, Sh/L
H は T を入力信号のサンプリング集合として,¶·noh hk, (k + 1) ∩ Ttl̀k , . . . , tĺk :=L Lとしたとき
step(Sh/LHug )[k] =+Σ11 := Σ12 (Ψb22 − Ψa22 ) − Ψ>b22 Ψa32 ,†>Σ12 := (I − Ψb22 Ψb33 )Ψa23 ,Σ13 := Σ12 (Ψb24 − Ψa24 ) − Ψ>b22 Ψa34 ,Σ22 := Ψb33 Ψ†b23 − Ψa33 Ψ†a23 ,0ĺX
k −1¡µ¶htl̀k − k ug [`lk − 1]L

(tl+1 − tl )ug [l]l=l̀kµ+h(k + 1) − tĺkL¢¶ug [´lk ]のように計算できる.
stepSh/L
H は 1 ステップ分の非因果性を含んでいるが,その直後に接続された遅れにより, システム全体は因果的になっている.
Ke は通常の離散時間のシステムであり, 実装可能である.

5数値例

本稿では, 例として倍率 3/π の解像度変換器を設計した. こういった, 無理数比や公倍数の大きな比の解像度変換は, 既存の手法による因果的な処理が困難である.

4030Gain [dB]20100-10-20-30-40-500.0001

0.001

0.01
0.1
Frequency [Hz]110Fig. 5: ゲイン線図. 実線: W1 , 破線: W2Fig. 8: 3/π 倍解像度変換画像

0.5
0.45
0.4
0.35
0.3
0.25
0.2
0.15
0.1
0.05
0-0.050510Time [sec]1520Fig. 6: K̃ のインパルス応答
Ke の設計には, h = 1, N = 6, L = 10 とし, W1 とW2 には Fig. 5 のような伝達関数を用いた. 得られたK̃ のインパルス応答を Fig. 6 に示す.
さて, この処理系を用いて Fig. 7 を 3/π 倍解像度変換したものが Fig. 8 である.

6終わりに不等間隔にサンプルされた信号からの H ∞ 最適な復元法を提案した. 問題は, 離散時間の H ∞ 制御問題に帰着する. さらに, 提案法を用いれば従来の手法では困難な解像度変換なども, 高品位かつ高速に実行できる.

参考文献
1) A. Aldroubi and K. Grŏchenig: Nonuniform sampling
and reconstruction in shift-invariant spaces, SIAM
Review, 43-4, 585/620 (2001)
2) G. Balas and R. Chiang and A. Packard and M. Safonov: Robust Control Toolbox Version 3, The MathWorks (2005)

Fig. 7: 原画像 (lena)

3) D. Capel and A. Zisserman: Computer vision applied
to super resolution, IEEE Signal Processing Magazine, 20-3, 75/86 (2003)
4) T. Chen and B. Francis: Optimal Sampled-Data Control Systems, Springer (1995)
5) Y. C. Eldar and A. V. Oppenheim: Filterbank reconstruction of bandlimited signals from nonuniform andgeneralized samples, IEEE Trans. Signal Processing,48-10, 2864/2875 (2000)
6) P. Gahinet and P. Apkarian: A linear matrix inequality approach to H ∞ control, Int. J. Robust and Nonlinear Contr., 4, 421/448 (1994)
7) P. Gahinet and A. Nemirovski and A. J. Laub and M.
Chilali: LMI Control Toolbox, The MathWorks (1995)
8) H. Ishii, Y. Yamamoto, and B. A. Francis: Samplerate conversion via sampled-data H ∞ control, In
Proc. of the 38th Conf. on Dicision and Control,3440/3445 (1999)9) 掛水: サンプル値 H ∞ 制御理論によるディジタル画像復元, 修士論文, 京都大学 (2005)
10) J. Kybic, T. Blu, and M. Unser: Generalized sampling: A variational approach—Part I: Theory, IEEETrans. Signal Processing, 50-8, 1965/1976 (2002)11) P. Marziliano and M. Vetterli: Reconstruction of irregularly sampled discrete-time bandlimited signals
with unknown sampling locations, IEEE Trans. Signal Processing, 48-12, 3462/3471 (2000)
12) L. Mirkin and Z. J. Palmor: A new representationof the parameters of lifted system, IEEE Trans. Automat. Contr., 44-4, 833/840 (1999)
13) L. Mirkin and H. P. Rotstein and Z. J. Palmor: H 2and H ∞ design of sampled-data system using lifting.
Part I: General framework and solutions, SIAM J.
Control Optim., 38-1, 175/196 (1999)
14) L. Mirkin and G. Tadmor: Yet Another Discretization, IEEE Trans. Automat. Contr., 48-5, 891/894(2003)
15) C. Moler and C. V. Loan: Nineteen dubious ways
to compute the exponential of a matrix, twenty-fiveyears later, SIAM Review, 45-1, 3/49 (2003)16) M. Unser, A. Aldroubi, and M. Eden: B-spline signal processing: Part I—Theory, IEEE Trans. Signal
Processing, 41-2, 821/833 (1993)17) M. Unser and J. Zerubia: A generalized sampling theory without band-limiting constraints, IEEE Trans.
Curcuits and Systems—II: Analog and Digital Signal
Processing, 45-8, 959/969 (1998)18) 和田: 不等間隔にサンプルされた信号からの H ∞ 最適復元, 修士論文, 京都大学 (2007)
19) Y. Yamamoto: A function space approach tosampled-data control systems and tracking problems,IEEE Trans. Automat. Contr., 39-4, 703/713 (1994)