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

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

Info

Publication number
JP5566969B2
JP5566969B2 JP2011164061A JP2011164061A JP5566969B2 JP 5566969 B2 JP5566969 B2 JP 5566969B2 JP 2011164061 A JP2011164061 A JP 2011164061A JP 2011164061 A JP2011164061 A JP 2011164061A JP 5566969 B2 JP5566969 B2 JP 5566969B2
Authority
JP
Japan
Prior art keywords
sequence signal
frequency response
impulse response
control
input
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
JP2011164061A
Other languages
English (en)
Other versions
JP2013029913A (ja
Inventor
原 義 之 石
倉 晋 司 高
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Toshiba Corp
Original Assignee
Toshiba Corp
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Toshiba Corp filed Critical Toshiba Corp
Priority to JP2011164061A priority Critical patent/JP5566969B2/ja
Priority to US13/418,460 priority patent/US9093114B2/en
Priority to CN201210068435.0A priority patent/CN102903370B/zh
Publication of JP2013029913A publication Critical patent/JP2013029913A/ja
Application granted granted Critical
Publication of JP5566969B2 publication Critical patent/JP5566969B2/ja
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G11—INFORMATION STORAGE
    • G11B—INFORMATION STORAGE BASED ON RELATIVE MOVEMENT BETWEEN RECORD CARRIER AND TRANSDUCER
    • G11B20/00—Signal processing not specific to the method of recording or reproducing; Circuits therefor
    • G11B20/10—Digital recording or reproducing
    • G11B20/10009—Improvement or modification of read or write signals
    • G—PHYSICS
    • G11—INFORMATION STORAGE
    • G11B—INFORMATION STORAGE BASED ON RELATIVE MOVEMENT BETWEEN RECORD CARRIER AND TRANSDUCER
    • G11B5/00—Recording by magnetisation or demagnetisation of a record carrier; Reproducing by magnetic means; Record carriers therefor
    • G11B5/48—Disposition or mounting of heads or head supports relative to record carriers ; arrangements of heads, e.g. for scanning the record carrier to increase the relative speed
    • G11B5/54—Disposition or mounting of heads or head supports relative to record carriers ; arrangements of heads, e.g. for scanning the record carrier to increase the relative speed with provision for moving the head into or out of its operative position or across tracks
    • G11B5/55—Track change, selection or acquisition by displacement of the head
    • G11B5/5521—Track change, selection or acquisition by displacement of the head across disk tracks

Landscapes

  • Engineering & Computer Science (AREA)
  • Signal Processing (AREA)
  • Feedback Control In General (AREA)
  • Moving Of The Head To Find And Align With The Track (AREA)
  • Testing And Monitoring For Control Systems (AREA)

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)で与えられることが知られている。
Figure 0005566969
 数式(1)において、Y(z)は出力の伝達関数、U(z)は入力の伝達関数、G0(z)は入力側のサンプリング周期における制御対象の伝達関数、Pはマルチレート比である。本実施形態で同定すべき制御対象の伝達関数は、数式(1)におけるG0(z)となる。数式(1)は、例えばマルチレート比P = 2のとき、
Figure 0005566969
 となる。以下のでは、マルチレート比P = 2を用いることを仮定して、本実施形態の同定アルゴリズムを説明する。
まず、本実施形態ではG0(z) の同定モデル
Figure 0005566969
 を、未知パラメータ
Figure 0005566969
 を用いて、以下数式(3)のFIRモデルとする。
Figure 0005566969
 また、入力の時系列をu[0], u[1], …, u[m]とすれば、数式(2)のU(z)は数式(3)と同様に、以下数式(4)のFIRモデルで表すことが出来る。
Figure 0005566969
 さらに、
Figure 0005566969
 となる。
ここで、数式(2)のG0(z)を
Figure 0005566969
 に置き換えた場合の出力、すなわちモデル出力の各時刻k = 0, 1, 2,…における時系列
Figure 0005566969
 は、数式(3)-(5)より以下のように単純な畳み込みによって計算することが出来る。
Figure 0005566969
Figure 0005566969
Figure 0005566969
Figure 0005566969
Figure 0005566969
 （以下同様）
以上より、モデル出力の時系列
Figure 0005566969
 は、kが奇数のとき出力が間引き（ダウンサンプリング）された応答になることが分かる。
同定モデル
Figure 0005566969
 の未知パラメータ
Figure 0005566969
 を推定するには、実際のハードウェア（例えば磁気ディスク装置のヘッド位置決め制御系）に入力信号u[0], u[1], …, u[m]を印加した場合の出力の時系列y[k]と、モデル出力の時系列
Figure 0005566969
 との誤差e[k]の二乗和を最小化する問題（最小二乗法）を解けば良い。マルチレート比P = 2の入力多重マルチレート系では、kが奇数のときは
Figure 0005566969
 とも間引きされ0になるのでe[k] (k = 0, 1, 2,…) は、以下数式(6)のようになる。
Figure 0005566969
 従って、解くべき最小二乗法は、
Figure 0005566969
 となる。数式(7)の解は、未知パラメータ
Figure 0005566969
 と入力u[0], u[1], u[2],…から成る以下数式(8)のベクトル、
Figure 0005566969
 を用いて以下数式(9)で与えられる。
Figure 0005566969
次に、数式(9)を実際に計算することを考える。本実施形態ではu[k]に擬似白色信号の一種であるM系列を用いるので、数式(9)の総和回数Nを、M系列の周期Mpより、
Figure 0005566969
 と選ぶと、数式(9)右辺の1番目の括弧内、
Figure 0005566969
 は、
Figure 0005566969
 となる。ここで、M系列の周期性を考慮すると数式(12)は、
Figure 0005566969
 と変形出来、M系列1周期分の自己相関行列となることが分かる。従って、白色信号の自己相関の性質より数式(13)の非対角項は、以下数式(14)のように近似的に0とすることが出来る。
Figure 0005566969
数式(14)のσu 2はu[k]の分散である。以上より、数式(9)右辺の逆行列演算は単純な除算に置き換えることが出来、計算負荷を大幅に低減することが出来る。すなわち、計算負荷の大きい逆行列演算を行う必要はない。数式(14)を成立させるために必要最小限の同定実験時の取得データ長Lは、
Figure 0005566969
 となる。同定実験時の取得データ長の概念を図１に示す。図１では、マルチレート比P = 2、M系列の周期Mp=7の場合のデータ取得長が示される。
以上より、マルチレート比P = 2の場合の未知パラメータ推定アルゴリズムは、
Figure 0005566969
 となる。一般的な偶数倍マルチレート比（P = 2, 4, 6, …）の場合は、
Figure 0005566969
 となる。
以上の数式(15),(17)が本実施形態の基本的なアルゴリズムである。また、数式(17)によって計算される
Figure 0005566969
 はG0(z)のインパルス応答推定値（FIRモデルの分子多項式係数）であるため、
Figure 0005566969
 を離散フーリエ変換すればG0(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分の信号および出力データを取得したら、制御対象のインパルス応答推定値を計算する。具体的に、制御対象の同定モデルであるＦＩＲモデルのインパルス応答推定値
Figure 0005566969
 を、数式(17) により計算する（S6）。そして、
Figure 0005566969
 を離散フーリエ変換する（S7）。離散フーリエ変換の結果が、求めるべき同定対象の周波数応答となる。
図６に、本実施形態に係わる装置の機能ブロック図を示す。
M系列信号生成部１１は、擬似白色信号の一種であるM系列信号を、データ取得長L＝Mp×P−１分、生成する。
生成したM系列信号の各時刻の信号がそれぞれ同定対象１２に入力され、出力データが制御対象（同定対象）１２から出力される。
記憶部１５は、生成部１１で生成されたM系列信号と、各時刻の信号に対応して同定対象１２から得られた出力データを記憶する。
インパルス応答推定部１３は、記憶部１５に記憶されたM系列信号および出力データに基づき、制御対象１２のインパルス応答推定値を計算する。この計算は、数式（A）により行う。
Figure 0005566969
数式(A)は、前述した数式(17)と同じものである。前述した通り、u[k]は、時刻ｋ（k番目）で取得したM系列信号の成分、y[k]は、u[k]を同定対象１２の入力としたときに得られる出力データ、σu 2はu[k]の分散である。
フーリエ変換部１４は、
Figure 0005566969
 を離散フーリエ変換することによって、制御対象の周波数応答を同定する。
以上、本実施形態によれば、入力多重マルチレート系において、制御入力側のサンプリング周期における制御対象の周波数応答を、高速かつ簡易な演算によって同定することが可能となる。
以下、本発明の実施例を計算機シミュレーションによって示す。
まず、複数の共振モードを持つ機械システムを想定した、以下数式(18)の離散時間系伝達関数を、計算機シミュレーションにおける同定対象とする。
Figure 0005566969
サンプリング周期は50[μs]とする。この同定対象にサンプリング周期50[μs]でM系列信号を入力し、マルチレート比P = 2でダウンサンプリングされた出力（サンプリング周期100[μs]）を計算機シミュレーションによって得る。このときの同定対象の入出力関係を示したブロック線図を図３に示す。M系列のシフトレジスタ段数は13段（周期Mp = 8191）とした。また、数式(17)で推定するインパルス応答
Figure 0005566969
 の長さは1024とした。これは、数式(8)においてn = 1023、すなわち、
Figure 0005566969
 とすることに相当する。最後に、数式(17)を実行し得られたインパルス応答推定値
Figure 0005566969
 に対して1024点のFFTを実行し同定対象の周波数応答を得る。ここで、入出力の取得データ長を本実施形態(15)式に従い、
Figure 0005566969
 とした場合と、数式(15)式に従わず任意に設定した値、
Figure 0005566969
 とした場合とで数式(17)の計算を行い、取得データ長が同定精度に与える影響を調べる。
以上の設定の下で計算機シミュレーションを実行した結果（同定対象の周波数応答）を図４および図５に示す。図４および図５において、実線が同定された周波数応答、点線が同定対象の真の周波数応答である。また、図中5000[Hz]に表示されている破線はサンプリング周期100[μs]時のナイキスト周波数を示している。本実施形態の同定手法と異なり、マルチレート系の特性を陽に考慮しない同定手法では、ナイキスト周波数以上（破線右側）の周波数応答を同定することは困難である。
図４より、本実施形態の数式(15)を満たす取得データ長Lを設定した場合は、同定された周波数応答が真の周波数応答と良く一致しており、ナイキスト周波数以上の帯域においても正確に同定出来ていることが分かる。一方、図５より、数式(15)を満たさない任意の取得データ長Lを設定した場合、周波数応答の大まかな一致は見られるものの、数式(15)を満たす場合に比べて、真の周波数応答との誤差が大きくなっていることが分かる。以上の計算機シミュレーションによる実施例により、本実施形態(15),(17)式の効果が検証された。
なお、磁気ディスク装置のヘッド位置決め制御系のように、ナイキスト周波数以上に共振モードを有する制御対象に対して、本実施形態のアルゴリズムによって得られる制御入力側のサンプリング周期における周波数応答は、共振モードを安定化するために用いられるマルチレートノッチフィルターの設計に利用することが出来る。

Claims (4)

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

Priority Applications (3)

Application Number Priority Date Filing Date Title
JP2011164061A JP5566969B2 (ja) 2011-07-27 2011-07-27 マルチレート系の高速周波数応答同定法および高速周波数応答同定装置
US13/418,460 US9093114B2 (en) 2011-07-27 2012-03-13 Method for identifying high-speed frequency response of multirate system and apparatus therefor
CN201210068435.0A CN102903370B (zh) 2011-07-27 2012-03-15 多速率系统的高速频率响应识别法和装置

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
JP2011164061A JP5566969B2 (ja) 2011-07-27 2011-07-27 マルチレート系の高速周波数応答同定法および高速周波数応答同定装置

Publications (2)

Publication Number Publication Date
JP2013029913A JP2013029913A (ja) 2013-02-07
JP5566969B2 true JP5566969B2 (ja) 2014-08-06

Family

ID=47575572

Family Applications (1)

Application Number Title Priority Date Filing Date
JP2011164061A Active JP5566969B2 (ja) 2011-07-27 2011-07-27 マルチレート系の高速周波数応答同定法および高速周波数応答同定装置

Country Status (3)

Country Link
US (1) US9093114B2 (ja)
JP (1) JP5566969B2 (ja)
CN (1) CN102903370B (ja)

Families Citing this family (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP5813151B2 (ja) 2014-02-21 2015-11-17 ファナック株式会社 制御ループの周波数特性を算出する機能を有する数値制御装置
CN105242111B (zh) * 2015-09-17 2018-02-27 清华大学 一种采用类脉冲激励的频响函数测量方法
US10862717B2 (en) * 2019-02-05 2020-12-08 Tektronix, Inc. Multirate data for S-parameter extraction
JP2024039769A (ja) 2022-09-12 2024-03-25 株式会社東芝 ディスク装置の製造方法、および、プログラム

Family Cites Families (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US5721648A (en) * 1992-04-10 1998-02-24 Seagate Technology, Inc. Multirate digital control system for use with a system having a linear transfer function, such as a head positioning system in a magnetic disc drive
JP2001249702A (ja) 1999-12-27 2001-09-14 Hitachi Ltd 位置決め制御方法
US20030090275A1 (en) * 2001-09-20 2003-05-15 Acorn Technologies, Inc. Determining a continuous-time transfer function of a system from sampled data measurements with application to disk drives
CN1926768B (zh) 2004-03-03 2010-07-14 独立行政法人科学技术振兴机构 信号处理装置和方法
JP2006195543A (ja) * 2005-01-11 2006-07-27 Fuji Electric Holdings Co Ltd モデル同定装置およびモデル同定プログラム
US8050160B2 (en) * 2008-08-29 2011-11-01 Kabushiki Kaisha Toshiba Characterizing frequency response of a multirate system

Also Published As

Publication number Publication date
US9093114B2 (en) 2015-07-28
US20130030743A1 (en) 2013-01-31
JP2013029913A (ja) 2013-02-07
CN102903370A (zh) 2013-01-30
CN102903370B (zh) 2015-06-10

Similar Documents

Publication Publication Date Title
JP4994448B2 (ja) モーダルパラメータの推定方法および装置
JP5566969B2 (ja) マルチレート系の高速周波数応答同定法および高速周波数応答同定装置
US20140254805A1 (en) Systems and methods for protecting a speaker
JP6143943B2 (ja) 二つのサブシステムの協調シミュレーション方法及び装置
KR20130054236A (ko) 시변 필터들을 이용하는 전체 파동장 역산
Muhamad et al. On the orthogonalised reverse path method for nonlinear system identification
JP2008107334A (ja) 計測震度概算装置、それを用いた計測震度概算システム及び計測震度概算方法
JP2007141227A (ja) ランダムサンプルタイミング法及びこれを使用するシステム
CN109374966A (zh) 一种电网频率估计方法
JP2011047953A (ja) 振動センサ及びその動作方法
WO2016131184A1 (zh) 用于测距的方法和装置
CN106575118A (zh) 过程工厂中一致的稳定状态行为的检测
JP2014041565A (ja) 時系列データ解析装置、方法、及びプログラム
Ding et al. Parameter identification of multi-input, single-output systems based on FIR models and least squares principle
JP2018151589A5 (ja)
Mendrok et al. Operational modal filter and its applications
JP2006238664A (ja) 電力系統の動特性定数推定方法と装置、電力系統解析方法と装置、電力系統安定化方法と装置、およびプログラム
JP5738778B2 (ja) 最適モデル推定装置、方法、及びプログラム
US8326907B2 (en) Method for establishing a simulating signal suitable for estimating a complex exponential signal
JP5649654B2 (ja) 導体長計測装置及び導体長計測方法
JP5852935B2 (ja) 伝達関数推定装置、伝達関数推定方法、および、伝達関数推定プログラム
JP6154777B2 (ja) 高速畳込近似装置、高速畳込近似方法、プログラム
JP2004521338A5 (ja)
JP6931830B2 (ja) インダクタンスの計測方法、計測システム、及びプログラム
JP6543298B2 (ja) 振動波形データ間の同期ずれ量計算方法、及び、振動波形データ間の同期補正方法

Legal Events

Date Code Title Description
A621 Written request for application examination

Free format text: JAPANESE INTERMEDIATE CODE: A621

Effective date: 20130906

A977 Report on retrieval

Free format text: JAPANESE INTERMEDIATE CODE: A971007

Effective date: 20140514

TRDD Decision of grant or rejection written
A01 Written decision to grant a patent or to grant a registration (utility model)

Free format text: JAPANESE INTERMEDIATE CODE: A01

Effective date: 20140520

A61 First payment of annual fees (during grant procedure)

Free format text: JAPANESE INTERMEDIATE CODE: A61

Effective date: 20140618

R151 Written notification of patent or utility model registration

Ref document number: 5566969

Country of ref document: JP

Free format text: JAPANESE INTERMEDIATE CODE: R151