| 前のページ | English Here 目次 | 次のページ |
| Copyright Maeda Yutaka |
ニュース！ 本ＨＰ 工業数学が工学社より出版されました。  ２．線形システム ２．１ 畳み込み積分(Convolution) 畳み込み積分は  と定義されています。しょっちゅう見かけるのですが、この式の直感的な意味が理解しづらく思えます。 ｆとｇを合成して作ったｈは（ｔ）の関数です。 （τ）は何の役目が有るのでしょうか？ 具体的な例を基に、もう少し掘り下げて意味を考えてみましょう。 線形システム（RLC回路等）への入力として単位インパルス関数を考えます。ｇ（ｔ）はその応答関数とします。  この例はRC回路のもっともポピュラーな応答の例で、ｇは指数関数となります。 入力インパルスとｇの間には線形関係があり、インパルスをｆ倍するとgもｆ倍されます。  ここで、入力インパルスの大きさｆが、時間軸上の位置によって変わっているものとしましょう。 その場合以下の様に、時間軸上のインパルス位置をτとすると、ｆ（τ）と関数の形で表現できます。  一方、応答関数ｇの任意の時刻ｔにおける電圧は、図のとおりｇ（ｔ-τ）となります。  さて、以下の図はインパルスが3個連続して入力された場合です。その応答は、それぞれの応答の合計値、つまりそれぞれの応答を重ね合わせた結果となります。  そこで、入力信号の波形をインパルスの連続体と考え、Σによる合計を積分に置き換えると、 という、先ほどの畳み込み積分の式が出来上がりました。 ここで最初の疑問、『τは何の役目があるのか？』ですが、この例の場合は『入力波形の形状をスキャンするための時間軸上のパラメーター』と言うことになります。 そしてｈ（ｔ）は、その入力波形に対する応答波形を表しています。 一般に畳み込み積分の演算子を以下の様にアスタリスク’＊’で表します。  また、説明は省略しますが、以下が成り立ちます。  ところでこの例でお分かりの様に、畳み込みが出来るのは、重ね合わせが出来る、つまり『線形性』と、時間がそのままシフトできる、つまり『時不変性』によるもので、線形時不変システムの特徴と言えます。 逆に線形時不変システムでは単位インパル応答さえ分かれば、どんな入力に対する応答も、畳み込み積分により求められると言う事です。 尚、離散信号においては畳み込み和が同様の演算となります。 フーリエ変換（FT）では関数の畳み込み積分は、変換により関数の積となります。 FT f(t)*g(t) ←→ F(ω）G（ω） 変換は対を成しているので、以下のように逆も成り立ちます。 FT f(t)g(t) ←→ F(ω）*G（ω） これは、ここで紹介したその他の変換（離散フーリエ変換、ラプラス変換等）でも同様に成り立つ普遍的で、非常に重要な原理といえます。 ２．２ 直交変換と固有値解析 固有値問題とは物理、工学等の様々な分野で共通に遭遇する課題パターンの一つであると言われています。固有値解析とは、その場合に利用される解法であり、応用範囲の広い普遍的な技法と言えます。 今後いろいろな場面で必要となる事が予想されるため、直交変換における事例を考察しておきます。 DCTやDFTなどの直交変換を考えます。線形変換の式として以下の様に現されます。  変換行列aにより、ｘの時間ベクトルからｙの周波数ベクトルに変換されました。 さらに以下の様にｙの各要素に対しスカラーの係数列λを掛けて大きさを変えてみます。Λの様に対角線上以外が０である行列を対角行列（Diagonal matrix)と呼びます。  行列の式では、以下の様になります。  ,  周波数領域で各要素の大きさを変える事は、フィルターの操作と言えます。身近なところではステレオセットに付いているグラフィック・イコライザーアンプのスライドボリュームを操作する様なものでしょうか。 このｙ’を逆変換により元の時間領域に戻してみます。  当然、時間領域に戻されたと言っても、ｘ’は元のｘとは異なってしまいます。 この式は新たな変換行列Wを定義することにより、以下の様な線形変換の式に整理できます。  （  ）このWがイコライザーアンプに相当し、ｘより入力された音を変形してｘ’に出力するというイメージになります。 では問題です。『もしWが未知であった場合、ｘとｘ’からWを推定し、さらにそれが変換行列aとフィルターλより構成されているという事実を知る事が出来るのでしょうか？』 これがいわゆる固有値解析の一例です。Scilabを使って確認してみましょう。 DCTで周波数領域に変換後、フィルタリングした上で、再度IDCTにより時間領域に戻します。 ソースプログラム SciLab_src_EVAL1 フィルター無し（λを全て１）の状態での、インパルスｘに対する応答出力ｘｘです。   応答出力：ｘｘ 入力と同じインパルスが再現されてます。 フィルターに低域通過型の係数を設定してみます。   応答出力：ｘｘ 指数関数的な応答が得られました。フィルターの効果です。この場合の変換行列ｗです。  複雑な数字が並んでいます。しかしよく見ると縦横対称になっています。対称行列あるいはエルミート行列であることが推測されます。 実際にはこの変換行列は未知であり、実験等により推定する事になります。入力インパルスの位置を一つ一つずらしながら、対応する列を埋めて行く事により求められます。 今、未知のWが実験で求まったものとします。このデータよりフィルターの係数λを求めて見ましょう。この係数は固有値として求められます。Scilabでは固有値はspec(W)で計算できます。以下がその結果です。  順序は別にして、ちゃんとフィルターの係数が求められました。 どういう理屈でフィルターの係数が固有値なるものに相当するのか調べてみましょう。 変換行列は元々以下の通り、フィルター係数による対角行列Λを直交行列aと、その逆行列とで挟んだものでした。  この式の両辺にaの逆行列を右からかけて変形します。  今回、aは直交行列であったので逆行列は転置行列に置き換えることが出来ます。  ここで、変換行列aを以下の様に行ベクトルから構成されていると考えます。  , ただし は ，i=1～４の行ベクトルそうすると は上記、行ベクトルが転置されて列ベクトルに対応するため、以下の4組の式に分解されます。 ， i=1～４λはいわゆる固有値であり、aiはその固有値に対応する固有ベクトルです。
以上の事は、適切なベクトルを選ぶことによって、行列Wが固有値というスカラー量に変換できる事を示しています。そして、この様な関係を見つけることを固有値問題と言います。以下の式は量子力学における有名なシュレディンガーの波動方程式です。  この式も固有値問題となっています。観測量Hに対し、Ψは固有ベクトルに相当する波動関数、Eは固有値でありエネルギーを表します。つまりこれを解くと周波数スペクトルが求まることのなります。 今回の例では固有値はフィルターの係数でした。 ↓ そして固有ベクトルはDCT変換の行ベクトルに対応します。何らかの方法で行列Wが観測出来た場合に、それを変換行列と対角行列（この例ではフィルターの係数）に分解することは、工学的にも有用な手法です。この様な技法を対角化と言います。 未知の行列Ｗより固有値、固有ベクトルを求める方法は線形代数の教科書を参考にしてください。 明示的にフーリエ変換等の直行行列と逆行列でスカラー列を挟んだ正方行列の固有値はスカラー列と等しく、物理的にはフィルターです。また固有ベクトルは挟んだときの直行行列です。（一般的に順番は再現されませんが） ここで紹介したフィルター（固有値）は、この後解説するラプラス変換で重要な役割をはたすこと となります。 2．3 ラプラス変換・ｚ変換 ラプラス変換 フーリエ変換では虚数軸を周波数を示す係数軸として利用しているというイメージです。  なぜならe-iωtはtをパラメータとして複素平面でクルクルと回転するわけであり、ωはその回転速度係数となります。虚数軸上の回転速度係数ωの位置と、波形を対比させてみると以下の図の様になります。  次にeαt の様に指数部が実数の場合を考えてみます。α＝０では当然オール１の直流です。これを境にα＜０では１から０に収束します。α＞０では１から発散しています。  ラプラス変換はこの2つの軸をを組み合わせて線型システムの過渡動作を解析しようと言った処でしょうか。 変換式は以下の様にフーリエ変換の回転因子に振幅の発散・収束を表す実数による指数関数のe-αtが掛かってきます。  実数部と虚数部をまとめて複素数ｓと表記しします。積分が0から始まっているのは、ラプラス変換では因果関数と言って時刻0以降の信号のみを扱うためです。 ｓの係数による違いを複素平面上のそれぞれの位置に表示してみましょう。  虚数軸が周波数を、実数軸が発散・収束を表している事が良く分かります。 -------------------------------------------------------- ニュース！ 本ＨＰ 工業数学が工学社より出版されました。 -------------------------------------------------------- 以下本書増強部のダイジェストです。詳細は本書にてご確認ください。事例も豊富に記載されております。 尚、書籍には数式の動作確認をＳｉｌａｂによりステップバイステップで学習可能な工夫をしております。 ラプラス変換の謎 ラプラス変換は過渡現象の解析に利用されます。理由は過渡解析に必要となる微分方程式を代数的に扱えるからです。これは一般的に演算子法と呼ばれる数学分野となります。 ところでラプラス変換を勉強し始めて最初に感じる謎は 「なぜラプラス変換で微分方程式が解けるのか？フーリエ変換に実数係数αを足しただけなのに。αって何なんだ？」 ではないでしょうか？ 以下この謎を解くことを目的にラプラス変換を学習してみましょう。まずはステップ関数のラプラス変換からです。 ラプラス変換とは 以下のように時間ゼロから立ち上がりレベル’１’を継続する関数をユニットステップ関数：u(t) といいます。  ラプラス変換すると  大きさがｈの場合は  遅延時間τのステップ関数のラプラス変換は 
 以上を組み合わせて任意のパルスが作れます。 たとえば、以下のような大きさｈで幅τの場合  ゼロ時間立ち上がりのパルスに対し、τ時間遅らせて反転したパルスを足し合わせて  となります。 このような調子で色々な時間関数のラプラス変換をまとめると ラプラス変換－逆変換表
これ等を電気回路の過渡解析に応用してみましょう。 以下のＣＲ回路があり、ｘにレベルｈとなるパルス信号が加えられた場合の応答出力ｙを求めます。  キルヒホッフの法則より  ↓  微分方程式を解かずに、これをそのままラプラス変換します。  左辺のｘはパルス信号でありラプラス変換Ｘは  右辺は積分演算を変換表から拾って  左右両辺を結びＹについて解くと、  これがラプラス領域で求めた代数方程式の解となります。一切、微分方程式を解いていません。 後は時間関数を変換表やその組み合わせより逆引きします。  それらしい時間領域での指数関数が現れました。 そうです。これがラプラス変換マジックです！ いやー、なんとも驚きですね。 今回の例はごく単純なものでしたが。もっと複雑な微分方程式となるものも、同様の方法、つまりラプラス領域の代数方程式で解くことができます。 フーリエ変換に＋αしただけでどうしてこんなことができるのでしょうか？ 以下その秘密を探ってみましょう。 微分演算の性質 関数  を微分すると となります。 つまり三つの各項に微分演算を分配して足すことができます。 もちろん以下の通り、各項が三角関数でも問題ありません。 微分演算を各項に適用して加算することができます。  とこで、ろもうお気づきでしょうか？ 元の式はフーリエ逆変換と同じ形式をしています。 つまり時間信号をフーリエ変換により周波数領域に展開して、微積分演算を各項（周波数別に並んだ三角関数）に適用することが可能であり、その場合は単純な三角関数の微分で済んでしまうということです。 そして三角関数の微積分は複素指数関数に置き換えて表現すると
であり、iωが微積分演算子の役割をしている、いたってシンプルな表現となります。 以上をイメージで示すと、 三角関数の微積分  なお、ここではフーリエ変換・逆変換やラプラス変換・逆変換の間にある１/（2π）などの係数は本質的なものでは無いので省略して解説します。 さて、あらためてラプラス変換表で微積分を確認すると
このs=iω+αでα=0と置いた場合がフーリエ変換となるので、同じことを示しています。 つまりラプラス変換で微積分を扱う能力はもともとフーリエ変換に備わっていたと言う事です。 ところでラプラス変換の場合はどうして[iω]の代わりに[s]になるのでしょうか？ フーリエ変換の基底は、 [eiωt] でした。 複素表現された三角関数波群ですね。これを微分すると iω[eiωt] この頭に下りてきた[iω]を微分演算子だと騒いでいたわけです。 周波数領域に変換された時間信号は、各周波数別の基底に対する成分に分解されています。各基底を微分した係数を、周波数領域の信号にかけることは、周波数領域の信号を微分することと等価になるわけです。 一方ラプラス変換の基底は [est] です。微分すると s[est] です。つまり[s]がラプラス変換における微分演算子になる分けです。 拡張されたαの役割 次の謎は「αって何なんだ？」でしたね。 フーリエ変換が存在するためにはフーリエ積分が有限値に収束（絶対収束）する必要が有ります。フーリエ変換で様々な関数に演算子法を使おうとすると、少々適用範囲が狭くなります。 たとえば関数：etのフーリエ積分は収束しません。しかし以下の図の様にe-αtをかけると収束します。この例ではα=2としました。   強制圧縮・・・少々強引ですね。 しかし一旦積分が収束し周波数領域での関数が確定する（つまりラプラス変換が求まる）と、周波数領域（ラプラス領域）で先ほどの微積分等の演算を施すことができます。そしてその後、逆変換により時間関数に戻した時点で、e+2tを掛けることでつじつまが合ってしまいます。 ラプラス変換で微分方程式が解ける理由（暫定版） 以下の微分方程式でｘが与えられた場合、ｙについて解くことを考察しましょう。 [d/dt] ｙ(t) = x(t) 時間領域の信号をフーリエ変換したスペクトルをそのまま逆変換すると元の時間信号にもどります。図にするとこんな感じでしょうか。 [ ・ ]は行列、ベクタとの類似性を暗示しています。 逆変換（基底） スペクトル 変換 時間信号 x(t) ← [ eiωt ] ← X(ω) ← [ e－iωt ] ← x(t) 変換はＤＦＴで言うとフーリエ変換行列を時間信号ベクトルに掛けることに相当し、逆変換は基底行列に周波数信号（スペクトル）ベクトルを掛けることに相当します。この場合い入力と出力は同じものになります。 前出の微分用フィルタ[iω]を挿入してみます。 逆変換（基底） 微分フィルタ スペクトル 変換 時間信号 x’(t) ← [ eiωt ] ← [iω] X(ω) ← [ e-iωt ] ← x(t) この場合 x’(t)は x(t)の微分出力となります。 数式で言うと、 x’(t) = [dx/dt] x(t) です。 ここで、微分用フィルタを逆特性を持ったフィルタに取替えてみます。スカラーなので単純に逆数でＯＫです。 逆変換（基底） 逆特性フィルタ 変換 時間信号 ｘf(t) ← [ eiωt ] ← [1/iω] Ｘ(ω) ← [ e-iωt ] ← ｘ(t) このｘf(t)をもう一度、微分演算の操作に戻すと x(t) ← [ eiωt ] ← [iω] ＜[1/iω] X(ω)＞ ← [ e-iωt ] ← xｆ(t) ↑ ↑ キャンセル ｘf(t)は逆特性が与えらているのでスペクトルは＜・＞の中のようになります。このため微分フィルタ[iω]をキャンセルしx(t)を出力します。つまりｘf(t)が求めるｙでした。 ｘf(t) ↓ [d/dt] ｙ(t) = x(t) 以上がフーリエ変換による簡単な１階の微分方程式を解くメカニズムです。 周波数領域では信号（スペクトル）はスカラー化されているため、操作は代数演算で記述されます。ある演算をキャンセルするためには、単純に逆数を取れば良いだけですが、いま一つ実感がわきません。 微分演算をフィルタとして捉えると、微分方程式を解く行為は、逆特性フィルタで補償する操作と考えられます。この様に考えると少しは実感がわくのではないでしょうか。 この例に限らず、より高階の微分を含む複雑な微分方程式も、同様の仕組みで解けます。 そして、そこには固有値と共通の、面白い仕組みが存在するのですが・・・ 固有値と伝達関数 以下、前出の微分演算フィルターの信号フローです。 直感による理解しやすさのためにラプラスの周波数領域変数： s=iω+α で、実数部αを0として解説します。 ～～～周波数領域～～～ x’(t) ←[ eiωt ] ←｛ Y(ω) ←[iω] ←X(ω) ｝←[ e-iωt ] ←x(t) この中で、以下の様に変換・フィルタ・逆変換をまとめて変換Wと置いた場合、 [W]=[eiωt]・[iω]・[e-iωt] [iω]が[W]の固有値となり、[e-iωt]が固有ベクトルとなります。さらに周波数領域内で[iω]の入り口と出口に注目すると Y（ω）= iω・X(ω) の関係が見出せます。 周波数領域でのフィルタ係数は入力と出力の比を表し、周波数の代数関数としてH（ω）と表現されます。 H（ω） = Y（ω）/X(ω) = iω この関数を一般的に伝達関数（実際はH(s)です・・・）と呼びます。また変換Wを行列とした場合、伝達関数はWの固有値とみなすことができます。 回路内にR,L,Cをいろいろ持ち込んで、微分動作（iω）や積分動作（1/(iω)）を掛けたり足したりと組み合わせていくと、伝達関数は分数の形をした式になります。 一般形（以下もう一度、変数をωからｓに戻して説明します）で書くと、 多項式型 Σn=0～N ans-n H(s) = --------------- Σm=0～M bms-m 因数分解型 (s-z1) (s-z2) ・・・(s-zn) H(s) = H0-------------------- (s-p1) (s-p2) ・・・(s-pm) 因数分解型で分母のs=pとなる点は極と呼ばれ、その（複素平面上の）点の位置は伝達関数の特性（収束発散・振動）を解析する上で重要です。 ＜続く＞ 以上、[特講]工業数学で追加された「ラプラス変換の謎」に関して増強部です。今回ＨＰでは書籍と同様、微積分演算が可能な理由や微分方程式が解ける理由を解説しましたが、少々手抜きダイジェスト版です。詳細にご興味のある方はぜひ下記のご購入をご検討ください。 豊富な図や事例、その他見どころ、沢山あります。 -------------------------------------------------------- ニュース！ 本ＨＰ 工業数学が工学社より出版されました。 -------------------------------------------------------- 雑談：ラプラス変換表の眺め方 ラプラス変換の表を丸暗記するのもいいのですが、各項目の関連性に注目すると面白いことに気づきます。その分覚えやすいかもしれません・・・ δとu(t)の関係 δのラプラス変換は[１]です。δは積分するとユニットステップ関数u(t)になります。 そこでラプラス領域におけるユニットステップ関数は、δのラプラス変換[１]に対し積分演算子[1/s]が掛かったものと考えると納得できます。 時間領域 ラプラス（周波数）領域 δ －＜ラプラス変換＞→ [１] ｜ ｜ 〈積分：∫〉 〈 積分：1/s 〉 ↓ ↓ u(t) －＜ラプラス変換＞→ [１][1/s] u(t)とe-atとの関係 フーリエ変換の周波数軸は[iω]で虚数軸上に固定されていました。ラプラス変換では以下の様に、虚数軸からα離れた縦軸を周波数軸として取ることができます。 虚数軸[iω] α=2 ↑ | ： | ： | ： ---------------------→ 実数軸[α] | ： | ： | ： ここで、たとえばα=2と固定してみます。そうすると周波数軸αは虚数軸から2はなれた直線上となります。この場合u(t)のラプラス変換[1/s]はこの直線α上で以下の様に定義されたことになります。 α=2でのu(t)のラプラス変換： [1/（2+iω）] もともとαは積分収束用に掛ける指数関数の圧縮係数であり、逆変換時に同じ値で符号反転させて掛けると、伸張して元の時間関数に戻せます。（αを固定して変換・逆変換する時の例えです） これに対して、u(t)のラプラス変換後、つまり[1/（2+iω）]に、さらに実数部2に[＋a]とすると、[1/（2+a+iω）]。 これを逆変換+伸張しても余計なaの分だけ元に戻らずu(t)が指数減衰（あるいは指数発散）することになります。 逆に言うと、指数関数のラプラス変換はu(t)のラプラス変換[1/s]において実数部にオフセットの[+a]を加えたものと考えられます。これを図にすると、 u(t) －＜ラプラス変換＞→ [1/s] ｜ ｜ 〈指数減衰・発散〉 〈 実数部オフセットa 〉 ↓ ↓ [ e-at ] －＜ラプラス変換＞→ [1/(s+a)] e-atとsin，cosとの関係 これは以下の様にオイラーの公式により、互いに逆回転する複素指数関数の和と差であらわされます。   割と素直なのe-atの応用例でした。 数学が”１”から始まった？様に、ラプラス変換はδから始まった？と言えるのかもしれません。 ＜続く＞ ｚ変換 -------------------------------------------------------- ニュース！ 本ＨＰ 工業数学が工学社より出版されました。 -------------------------------------------------------- 以下、z変換の増強部の概要を順次追加します。 離散フーリエ変換（DFT）により、デジタル信号処理の基盤が整ったといえます。 離散時間信号や離散周波数信号が数列間の変換として扱えるようになりました。 フーリエ変換やラプラス変換で議論してきた微分方程式は数列上どのように扱えるのでしょうか？ 差分方程式 以下の差分方程式が連続システムの微分方程式に対応します。 差分方程の例： xの入力に対してa倍のｙが出力される線形システム｛ｈ｝を想定すると x[n]→｛h｝→y[n] 方程式では、 y[n]=a・x[n] この場合、システムの解析は簡単です。入力と出力は単純なスカラー倍です。 それでは、このシステムに対し1時刻前の出力のｂ倍を加えてみましょう。 方程式では、 y[n]=a・x[n]+b・y[n-1] a=0.7, b=0.3の場合を手計算（表計算に入れると簡単です）でシミュレーションしてみましょう。 ユニットステップ関数入力例 x[n] y[n-1] y[n] a ・1 + b ・0 = 0.7 a ・1 + b ・0.7 = 0.91 a ・1 + b ・0.91 = 0.973 ： ： 単位インパルス（δ）入力例 x[n] y[n-1] y[n] a ・1 + b ・0 = 0.7 a ・0 + b ・0.7 = 0.21 a ・0 + b ・0.21 = 0.063 ： ： a=0.7, b=-0.3の場合 単位インパルス（δ）入力例 x[n] y[n-1] y[n] a ・1 + b ・0 = 0.7 a ・0 + b ・0.7 = -0.21 a ・0 + b・-0.21 = 0.063 ： ： ぱっと見、過渡応答そのもの！・・・に見えませんか？ 微分方程式と同様、変化の度合いを表す係数を含んでいます。 この過渡応答もどきとアナログ回路の過渡応答の類似性を調べるために、以下のようなRCによる線形システム、つまりアナログフィルターを考えます。  入力電圧をEｘ、出力電圧をEｙとすると、以下の関係式となります。  この式の中の抵抗Rによるドロップ電圧ERは以下により求まります。  よって、以下の線形微分方程式が求められました。  この式において、係数を適当に取ると、以下の様に離散時間（⊿ｔ）の線形差分方程式に近似できます。  以上の考察により、「時間シフト足しこみ係数」の様なもので表現された差分方程式は、離散時間信号におけるフィルターを表現していると、捉えることができそうです。 アナログフィルターの解析に、微分方程式を解くことが必要となり、ラプラス変換が強力なツールとなりました。同様に「差分方程式用ラプラス変換」が有ると嬉しくなりますね！ ・・・と、この様な動機？により、ラプラス変換の式も離散化を行ってみましょう。 サンプリング周期Ｔsとして複素指数部を以下の様にｚとおきます。  そうすると、ラプラス変換の離散変換版として以下の式を得ます。  この変換をｚ変換とよびます。 早速、ユニットインパルスδ[n]のｚ変換を求めて見ましょう。  数式 ZT-1 あっさりと Z[δ[n]]=1 が求まりました。 パルスを遅らせた場合のZ[δ(n-1)]ではどうでしょうか？答えは Z[δ(n-1)]=z-1 です。 上記数式 ZT-1の変換式パルスの位置を右にk個ずらすとz-kとなることは明白です。 という事で、少なくともδ[n]に関しては以下が成り立ちます。 Z[δ(n-k)]=Z[δ[n]]・z-k 数式 ZT-2 当然逆方向（パルス位置を進める）は、 Z[δ(n＋k)]=Z[δ[n]]・z-k が成り立ちます。ただし、詳細は省略しますが、母体となったラプラス変換の定義と同様により因果関数に限定する必要があり、時刻0より先に食い込んだ信号は0になるような補正が必要となります。 次にユニットステップ関数のｚ変換を求めて見ましょう。  数式 ZT-3 と、求まりました。 早速数式 ZT-2がこの変換にも成り立つか試してみましょう。時間領域の任意の関数は連続するインパルスδの集合体（数式表現では畳み込みu(n)*δ(n))と考えられます。δに成り立つこの時間推移（遅れ・進み）の公式はその集合体uにも成り立つはずです。 数式 ZT-3において、１時刻遅れのユニットステップ関数Z[u[n-1]]では、第１項がゼロとなります。したがって Z[u[n-1]]=Z[u[n]]-1=z/(z-1)-1=1/(z-1) 一方、z-1を掛けた場合は Z[u[n]]・z-1=｛z/(z-1)｝・z-1=1/(z-1) と、一致しました。関数u(t)でも時間推移則が成り立っています。 ここでは簡単な例示のみとしましたが、このz-1は時間遅延演算子と言ってよいでしょう。そしてこの演算子は非常に重要でありデジタルフィルターの基礎となるものです。なぜならこれにより差分方程式を表現できるからです。 ２．4 フィルターの基礎 ラプラス変換の伝達関数の多項式型表記を少し修正してｚ変換用とします。  分母の第一項を１としました。伝達関数は所詮比例係数なので、このように定義しても特に問題ありません。 これを使った周波数領域での表現は  であり、伝達関数Hの定義をを代入し、変形すると  ↓  ここでz-kが遅延演算子であることに注目つつ、全体を逆変換します。  あれ！ なんと、これは差分方程式では？ と、この様に周波数領域で部品を組み合わせて機能設計し、その係数を以下のように時間領域に展開することが可能となります。  aはフィードフォアードの係数、bはフードバック係数、そしてz-1の箱は１時刻遅延を表します。 早速上記の構造をScilabで実装してみましょう。 簡単な例として前出の差分方程式で手計算した例：第０次フィードフォアード係数を0.7、第1次フィードバックを-0.3とし単位インパルスを入力します。 dfilter.sci ソースプログラム（右クリックで保存） 以下インパルスに対する応答をシミュレーションしています。  アナログ回路の応答でお馴染みの出力結果となりました。もちろん手計算の結果とも一致しています。尚、差分方程式例でのｂの定義はプラスマイナス逆になっているのでご注意ください。 フィルタと畳み込みの関係 ある線形システムにおいて、単位インパルス入力に対するインパルス応答が以下の様（でたらめな数値ですが・・・）に求まっていたとします。  入力関数を以下の通りとします。  出力は以下の通り、畳み積分の離散信号版である畳み込み和で求められます。  畳み込み和の筆算方法です。意外に簡単にですので是非、試してみてください。  結局、出力ｙ（ｔ）は以下の通りとなりました。  確認のためscilabで畳み込み和の関数：lconv(g,x)を作ってみました。 ソースプログラム SciLab_src_lconv 以下、実行結果です。  ｇ（ｔ）、x(t)それぞれの周波数領域表現(この場合はｚ変換）をG（z)、X（z）とすると、 Z-1[ G（z)・X（z）] = ｇ（ｔ）*x(t) G（z)がフィルタの特性を表す伝達関数とみなせます。このことはｚ変換に限らず、フーリエ変換やラプラス変換でも同様のことが言えます。 周波数領域表現（Frequency-domain representation) 次はフィルターの周波数特性をｚ変換を利用して分析をしてみましょう。  ここで、α=0とおくと、ｚ変換の伝達関数H(z)は離散時間フーリエ変換（DTFT)による伝達関数F（ω）と同様となります。したがって、 振幅特性は  位相特性は  データ列の入出力と、複数個の係数ですべてが決まります。サンプリング時間Tsは必要なときに意識すればよいのでTs=１としωを正規化角周波数（-π～π[rad/s]）として扱います。 したがって時間領域でのフィルタの係数が指定されて、その特性を知りたい場合は、   を計算すれば良いことになります。 具体的には上記の式に各係数を割り当て、その複素指数関の組み合わせを解析的に解く（ツールとしてはMaximaなどの数式解析ソフトを使う？）なり、ωをScilab等で離散的にスイープすれば特性を知ることができます。 例１） 時間領域での1次帰還型フィルターの例と同様、以下の通りとしました。 第0次フィードフォアード係数： b0= 0.7 第1次フィードバック係数： a1=-0.3 Maximaによる数式ベースの解析 ・・・[特講]工業数学未掲載（補足です） 以下の様に、ほぼ数式の通り記述するだけで簡単に特性が求まります。
青が振幅特性、赤が位相特性です。横軸は正規化角周波数です。 Scilabによる数値解析 respfreq Scilab ソースプログラム（右クリックで保存）  青が振幅特性、緑が位相特性です。横軸は正規化角周波数です。 Scilabによる数値解析でもMaxmaと同様の結果を得ました。 例２） 3次無帰還型フィルターの例です。 フィードフォアード係数：［1，３，３，１］/８ フィードバック係数：０ Scilabによる数値解析 ＊上記例１)のScilabのコード上、コメントアウトしている行を有効にして実行してください。  -------------------------------------------------------- 以上、 [特講]工業数学（工学社出版）での増強部の概要・補足を順次追加中！-------------------------------------------------------- |
| 前のページ | 目次 | 次のページ |
,
,
,
