JP5566969B2

JP5566969B2 - マルチレート系の高速周波数応答同定法および高速周波数応答同定装置

Info

Publication number: JP5566969B2
Application number: JP2011164061A
Authority: JP
Inventors: 原義之石; 倉晋司高
Original assignee: Toshiba Corp
Current assignee: Toshiba Corp
Priority date: 2011-07-27
Filing date: 2011-07-27
Publication date: 2014-08-06
Anticipated expiration: 2031-07-27
Also published as: US9093114B2; CN102903370A; CN102903370B; JP2013029913A; US20130030743A1

Description

本発明の実施形態は、マルチレート系の高速周波数応答同定法および高速周波数応答同定装置に関し、たとえば制御出力のサンプリング周期が制御入力のサンプリング周期の偶数倍となる入力多重マルチレート系において、制御入力側のサンプリング周期における制御対象の周波数応答を高速に同定するアルゴリズムに関する。

磁気ディスク装置のヘッド位置決め制御系などに用いられている入力多重マルチレート系においては、制御系設計時に制御入力側のサンプリング周期における制御対象の周波数応答（伝達関数）を同定することが必要となる場合がある。入力多重マルチレート系では、制御出力（観測出力）が間引き（ダウンサンプリング）されるため、一般的な周波数応答の同定方法はそのままでは適用出来ない。従って、入力多重マルチレート系の特徴を考慮した同定アルゴリズムが必要となる。

F.Ding and T.Chen, "Identification of dual-rate systems based on finite impulse response models", International Journal of Adaptive Control and Signal Processing, 18-7, pp.589-598, (2004)

本発明の一側面は、制御出力のサンプリング周期が制御入力のそれの偶数倍である入力多重マルチレート系において、制御入力側のサンプリング周期における制御対象の周波数応答を同定することを目的とする。

本発明の一態様としてのマルチレート系の高速周波数応答同定法は、制御出力のサンプリング周期が、制御入力のそれの偶数P倍である入力多重マルチレート系の、制御入力側のサンプリング周期における制御対象の周波数応答を同定する方法であって、M系列信号生成ステップと、インパルス応答推定ステップと、周波数応答同定ステップとを備える。

前記M系列信号生成ステップは、擬似白色信号の一種であるM系列信号の周期Mpに基づき、データ取得長Mp×P−１分のM系列信号を生成する。

前記インパルス応答推定ステップは、前記生成されたM系列信号と、前記生成されたM系列信号を入力して前記制御対象から得られる出力データに基づき、前記制御対象のインパルス応答推定値を計算する。

前記周波数応答同定ステップは、前記インパルス応答推定値を離散フーリエ変換することによって、前記制御対象の周波数応答を同定する。

本発明の実施形態の同定アルゴリズムにおけるデータ取得長Lを示す図。本発明の実施形態の実行手順を示すフローチャート。計算機シミュレーションにおける入出力関係を示したブロック線図。本発明の実施形態の数式(15)を満たした場合の周波数応答同定結果を示す図。本発明の実施形態の数式(15)を満たさない場合の周波数応答同定結果を示す図。本発明の実施形態に係わる装置の機能ブロック図。

本発明の実施形態に係る同定アルゴリズムは、同定入力として擬似白色信号の一種であるM系列信号、同定モデルとしてFIR（Finite Impulse Response）モデルを用いる。また、マルチレート比は偶数P倍、すなわち、制御出力のサンプリング周期が制御入力のそれの2, 4, 6…倍であることが条件となる。

ここで、入力多重マルチレート系の伝達関数は、以下数式(1)で与えられることが知られている。

数式(1)において、Y(z)は出力の伝達関数、U(z)は入力の伝達関数、G₀(z)は入力側のサンプリング周期における制御対象の伝達関数、Pはマルチレート比である。本実施形態で同定すべき制御対象の伝達関数は、数式(1)におけるG₀(z)となる。数式(1)は、例えばマルチレート比P = 2のとき、

となる。以下のでは、マルチレート比P = 2を用いることを仮定して、本実施形態の同定アルゴリズムを説明する。

まず、本実施形態ではG₀(z) の同定モデル

を、未知パラメータ

を用いて、以下数式(3)のFIRモデルとする。

また、入力の時系列をu[0], u[1], …, u[m]とすれば、数式(2)のU(z)は数式(3)と同様に、以下数式(4)のFIRモデルで表すことが出来る。

さらに、

となる。

ここで、数式(2)のG₀(z)を

に置き換えた場合の出力、すなわちモデル出力の各時刻k = 0, 1, 2,…における時系列

は、数式(3)-(5)より以下のように単純な畳み込みによって計算することが出来る。

（以下同様）
以上より、モデル出力の時系列

は、kが奇数のとき出力が間引き（ダウンサンプリング）された応答になることが分かる。

同定モデル

の未知パラメータ

を推定するには、実際のハードウェア（例えば磁気ディスク装置のヘッド位置決め制御系）に入力信号u[0], u[1], …, u[m]を印加した場合の出力の時系列y[k]と、モデル出力の時系列

との誤差e[k]の二乗和を最小化する問題（最小二乗法）を解けば良い。マルチレート比P = 2の入力多重マルチレート系では、kが奇数のときは

とも間引きされ0になるのでe[k] (k = 0, 1, 2,…) は、以下数式(6)のようになる。

従って、解くべき最小二乗法は、

となる。数式(7)の解は、未知パラメータ

と入力u[0], u[1], u[2],…から成る以下数式(8)のベクトル、

を用いて以下数式(9)で与えられる。

次に、数式(9)を実際に計算することを考える。本実施形態ではu[k]に擬似白色信号の一種であるM系列を用いるので、数式(9)の総和回数Nを、M系列の周期M_pより、

と選ぶと、数式(9)右辺の1番目の括弧内、

は、

となる。ここで、M系列の周期性を考慮すると数式(12)は、

と変形出来、M系列1周期分の自己相関行列となることが分かる。従って、白色信号の自己相関の性質より数式(13)の非対角項は、以下数式(14)のように近似的に0とすることが出来る。

数式(14)のσ_u ²はu[k]の分散である。以上より、数式(9)右辺の逆行列演算は単純な除算に置き換えることが出来、計算負荷を大幅に低減することが出来る。すなわち、計算負荷の大きい逆行列演算を行う必要はない。数式(14)を成立させるために必要最小限の同定実験時の取得データ長Lは、

となる。同定実験時の取得データ長の概念を図１に示す。図１では、マルチレート比P = 2、M系列の周期Mp=7の場合のデータ取得長が示される。

以上より、マルチレート比P = 2の場合の未知パラメータ推定アルゴリズムは、

となる。一般的な偶数倍マルチレート比（P = 2, 4, 6, …）の場合は、

となる。

以上の数式(15),(17)が本実施形態の基本的なアルゴリズムである。また、数式(17)によって計算される

はG₀(z)のインパルス応答推定値（FIRモデルの分子多項式係数）であるため、

を離散フーリエ変換すればG₀(z)の周波数応答を求めることが出来る。離散フーリエ変換はFFT（Fast Fourier Transform）を用いて容易に実行可能である。

以上示した本実施形態のアルゴリズムは、同定入力信号をM系列信号と規定しているため、一般的に正弦波掃引法などと比較して同定実験に要する時間を短縮することが出来る。また、数式(17)は単純な加算及び乗除算で構成されているため、本実施形態を適用するハードウェアに搭載されるマイクロプロセッサや外部の計算機に、簡易な記述で実装することが出来る。

本実施形態の実行手順を示したフローチャートは図２となる。この実行手順に従った処理は、ハードウェアに搭載されるマイクロプロセッサや、外部の計算機で実行される。

ステップS1において、取得したデータ長Lが、Mp×P−１に達したかを判定する。

Mp×P−１に達していないときは、次の時刻のM系列信号を生成し（S2）、同定対象（例：磁気ディスク装置）にM系列信号を入力する（S3）。なお、M系列信号の生成アルゴリズムは事前に与えておく。

制御対象（同定対象）の出力を観測し（S4）、M系列と出力データを、内部メモリまたは外部記憶装置等の記憶部に保存する（S5）。データ長Lが、Mp×P−1に達するまで、以上のステップS1〜S4を繰り返す。

データ長LがMp×P-1に達したら、すなわち、データ長L分の信号および出力データを取得したら、制御対象のインパルス応答推定値を計算する。具体的に、制御対象の同定モデルであるＦＩＲモデルのインパルス応答推定値

を、数式(17) により計算する（S6）。そして、

を離散フーリエ変換する（S7）。離散フーリエ変換の結果が、求めるべき同定対象の周波数応答となる。

図６に、本実施形態に係わる装置の機能ブロック図を示す。

M系列信号生成部１１は、擬似白色信号の一種であるM系列信号を、データ取得長L＝Mp×P−１分、生成する。

生成したM系列信号の各時刻の信号がそれぞれ同定対象１２に入力され、出力データが制御対象（同定対象）１２から出力される。

記憶部１５は、生成部１１で生成されたM系列信号と、各時刻の信号に対応して同定対象１２から得られた出力データを記憶する。

インパルス応答推定部１３は、記憶部１５に記憶されたM系列信号および出力データに基づき、制御対象１２のインパルス応答推定値を計算する。この計算は、数式（A）により行う。

数式(A)は、前述した数式(17)と同じものである。前述した通り、u[k]は、時刻ｋ（k番目）で取得したM系列信号の成分、y[k]は、u[k]を同定対象１２の入力としたときに得られる出力データ、σ_u ²はu[k]の分散である。

フーリエ変換部１４は、

を離散フーリエ変換することによって、制御対象の周波数応答を同定する。

以上、本実施形態によれば、入力多重マルチレート系において、制御入力側のサンプリング周期における制御対象の周波数応答を、高速かつ簡易な演算によって同定することが可能となる。

以下、本発明の実施例を計算機シミュレーションによって示す。

まず、複数の共振モードを持つ機械システムを想定した、以下数式(18)の離散時間系伝達関数を、計算機シミュレーションにおける同定対象とする。

サンプリング周期は50[μs]とする。この同定対象にサンプリング周期50[μs]でM系列信号を入力し、マルチレート比P = 2でダウンサンプリングされた出力（サンプリング周期100[μs]）を計算機シミュレーションによって得る。このときの同定対象の入出力関係を示したブロック線図を図３に示す。M系列のシフトレジスタ段数は13段（周期M_p = 8191）とした。また、数式(17)で推定するインパルス応答

の長さは1024とした。これは、数式(8)においてn = 1023、すなわち、

とすることに相当する。最後に、数式(17)を実行し得られたインパルス応答推定値

に対して1024点のFFTを実行し同定対象の周波数応答を得る。ここで、入出力の取得データ長を本実施形態(15)式に従い、

とした場合と、数式(15)式に従わず任意に設定した値、

とした場合とで数式(17)の計算を行い、取得データ長が同定精度に与える影響を調べる。

以上の設定の下で計算機シミュレーションを実行した結果（同定対象の周波数応答）を図４および図５に示す。図４および図５において、実線が同定された周波数応答、点線が同定対象の真の周波数応答である。また、図中5000[Hz]に表示されている破線はサンプリング周期100[μs]時のナイキスト周波数を示している。本実施形態の同定手法と異なり、マルチレート系の特性を陽に考慮しない同定手法では、ナイキスト周波数以上（破線右側）の周波数応答を同定することは困難である。

図４より、本実施形態の数式(15)を満たす取得データ長Lを設定した場合は、同定された周波数応答が真の周波数応答と良く一致しており、ナイキスト周波数以上の帯域においても正確に同定出来ていることが分かる。一方、図５より、数式(15)を満たさない任意の取得データ長Lを設定した場合、周波数応答の大まかな一致は見られるものの、数式(15)を満たす場合に比べて、真の周波数応答との誤差が大きくなっていることが分かる。以上の計算機シミュレーションによる実施例により、本実施形態(15),(17)式の効果が検証された。

なお、磁気ディスク装置のヘッド位置決め制御系のように、ナイキスト周波数以上に共振モードを有する制御対象に対して、本実施形態のアルゴリズムによって得られる制御入力側のサンプリング周期における周波数応答は、共振モードを安定化するために用いられるマルチレートノッチフィルターの設計に利用することが出来る。

Claims

制御出力のサンプリング周期が、制御入力のそれの偶数P倍である入力多重マルチレート系の、制御入力側のサンプリング周期における制御対象の周波数応答を同定する方法であって、
擬似白色信号の一種であるM系列信号の周期Mpに基づき、データ取得長Mp×P−１分のM系列信号を生成するM系列信号生成ステップと、
前記生成されたM系列信号と、前記生成されたM系列信号を入力して前記制御対象から得られる出力データに基づき、前記制御対象のインパルス応答推定値を計算するインパルス応答推定ステップと、
前記インパルス応答推定値を離散フーリエ変換することによって、前記制御対象の周波数応答を同定する周波数応答同定ステップと
を備えた方法。
前記インパルス応答推定ステップは、
u[k]を、前記M系列信号における時刻kの信号、y[k]を、u[k]を入力したときの前記制御対象の出力データ、σ_u ²をu[k]の分散としたとき、

により、前記インパルス応答推定値を計算する、
ことを特徴とする請求項１に記載の方法。
制御出力のサンプリング周期が、制御入力のそれの偶数P倍である入力多重マルチレート系の、制御入力側のサンプリング周期における制御対象の周波数応答を同定する装置であって、
擬似白色信号の一種であるM系列信号の周期Mpに基づき、データ取得長Mp×P−１分のM系列信号を生成するM系列信号生成部と、
前記生成されたM系列信号と、前記生成されたM系列信号を入力して前記制御対象から得られる出力データ基づき、前記制御対象のインパルス応答推定値を計算するインパルス応答推定部と、
前記インパルス応答推定値を離散フーリエ変換することによって、前記制御対象の周波数応答を同定する周波数応答同定部と
を備えた装置。
前記インパルス応答推定部は、
u[k]を、前記M系列信号における時刻kの信号、y[k]を、u[k]を入力したときの前記制御対象の出力データ、σ_u ²をu[k]の分散としたとき、

により、前記インパルス応答推定値を計算する、
ことを特徴とする請求項３に記載の装置。