周波数項が複素数になったら

交流波形には、振幅・周波数・初期位相がパラメータとして存在します。

Synchro PRIMO では、リサージュ外積比の平方根の 1/2 を計算し、その値をもとに逆sin を用いて周波数Ω（＝相対角周波）を計算しますが、逆 sin の引数が虚数になった場合、どうするか。

公式ででは、　sin ( j x ) = j * sinh( x ) とされています。

Synchro PRIMO(３倍角法）では、 実は双曲線関数も、計算ができてしまいます。　ただし、3-Ls3/Ls1 < 0  となり開平ができないので、 上の変換を使いarcsin の代わりに asrsinh を使用することになります。 

ほかに、Gauss関数、振幅が指数的に増減する信号も、挙動は異なりますが、周波数（一定値）が計算されます。

 

　どの場合でも、周波数＝０に近づけるとその極限値として考えられる波形は、値が一定（電気信号ならDC)　となります。ここを共通点として、パルス、双曲線関数、三角関数、指数関数を・・・「解析接続」して、統一的に扱えないものでしょうか・・・　周波数が存在するなら、基底関数を cos, sin ではなく、cosh, sinh などで構築する「双曲線フーリエ変換」や、ガウス関数を基底とするものが可能なのかもしれません。

 

別のケースとして、 | sin Ω ] > 1 となるケース：  3-Ls3/Ls1  > 4  つまり、 Ls3/Ls1 < -1 で発生します。これも「解析接続」という方法を用いると sin z > 1 のようなケースで、複素数 z を求めることができます。 ちなみに arcsin(2) = は・・・　https://www.wolframalpha.com/input/?i=asin+%282%29&lang=ja  ←こんな答えだそうです。　sinh をベースとする曲線に、振動項が加わるようですが、今の私には解釈不能・・・です。

大学院時代の思い出３

　大学院修了まで、実は３年半かかってしまいした。この「半」はなぜこうなったか。。。

　Ｄ３になる直前の３月末にマイコプラズマ肺炎にかかってしまい、４月・５月が半潰れの状態。経験者にはわかると思いますが、とにかく体がだるいのです。修了のため授業１科目を履修をしていたのですが、予備審査の条件として落とせない科目。うちの大学院では授業はクォータ制になっていて３ヶ月で勝負が決まります。ＧＷあけの中間試験は病気真っ最中ゆえわずかに合格点に足らず。最終的には病み上がりの状態で、気力・体力を振り絞った６月上旬の期末試験なのでした。結果は後半戦での頑張りがきいて合格。

　このタイミングで国際会議の採録通知がきました。試験直後６月中旬にCamera Ready稿を提出。ジャーナルを書くための実験もほぼ終わっておりあとは論文の形にまとめるだけ・・・ところがジャーナルなんて初めてなので書き方がピント来ない。ジャーナルの締め切りは７月１２日だったと記憶しています。

　副テーマもなんとかＯＫになり、「骨子提出」という段階にはいったろころ、教務からＮＧの知らせ。なんと単位の履修要件を一部ミスしていたのです。仕方ないので夏休みに集中講義というのがありそこで２科目を一気に履修。８月はほとんどつぶれてしまいましたが、ジャーナルへの投稿も終わっており、国際会議（EUSIPCO 2013 モロッコ・マラケシュ）も後は行くだけ。９月を楽しみにしていましたが、猛暑の２科目受講と渡航準備はきつかったです。

　結局、６月上旬に期末試験、中旬にCamera Ready稿投稿、下旬が副テーマ仕上げ、７月上旬でジャーナル投稿、その１週間後に札幌で研究会、８月から集中講義、９月に国際会議、帰国して中旬に豊橋で音響学会、１０月にはジャーナルが条件付き採録・回答を作成。やめときゃいいのに、１２月の研究会にエントリー、やめときゃいいのに、予備審査でへろへろになっているのは明らかなのに、３月の音響学会にエントリー。。。　記憶が飛んでいるのですが、どうせなら限界までやってやれと思っていたのかもしれません。

　ここまでのミスで、修了は「＋３ヶ月」が確定しています。あと３ヶ月は、予備審査で野球で言うなら「派手な一発」を審査委員から食らい、論文の修正に＋３ヶ月を要したこと。この派手な一発とは、本ブログで書いた「フーリエ変換の不確定性」で当時、勉強が不十分で学位論文の記述もあまあま、ということをさらに自覚出来ていなかったという痛恨のミス・・・しかしこの３ヶ月の猶予のおかげで他のあまい箇所の加筆もでき最終的には満足のいく形になったと思っています。

　学位論文を書いてからもう５年以上過ぎましたが、しばらくは「見るのも嫌だった」自分の論文、最近になりやっと客観視できるようになりました。自分でいうのも何ですが、序論はなかなかのもので、本当に自分が書いたのか？という印象。　

　今でも思い出すのは、１２月末に学位論文の「白紙のテンプレート」を見たときの絶望感（目次を書いたあと、この分量は無理・・・とマジで思った）。予備審査は２月末or３月上旬。執筆に与えられた時間は約６０日。とにかく頭をクリーンにしたく、宿泊客がほとんどゼロの冬の「白峰温泉」に３日間こもったのは懐かしい思い出です。

　苦しかったＤ３の後半で、「やけくそ」で発表した２件の研究会での発表（執筆開始直前の１２月と予備審査直後の翌３月）、これが現在のSynchro PRIMO を決定づけたものでした。２つのリサージュ外積を割り算するアイディアは９月のマラケシュで着想したのでした。たぶんＤ３（すでに５１歳）の１年は今までの人生でもっとも密度の濃かった一年だと思います。

　うちの大学院は別名「トラの穴」とも呼ばれていることを、どこからか聞きました。たしかにいくつかの勝負を勝たないと出られない、そんな感じです。本ブログではあえて大学院名は書きませんが、私のジャーナルをみればすぐに分かるかと思います。

 

位相差計算の「最下限」

　簡易式の位相差計算方法で、どこまで「小さな位相角」がはかれるか試してみました。

つい先日、発表されたブラックホールの視角は 1マイクロ秒角。 1度の 1/3600 さらに 100万分の1というオーダーです。

本Blogで紹介した方法で計算してみると、 「１マイクロ秒角」で有効数字３桁程度で計算できていました。

２つの信号の「位相差」を極小の状態でも計算できるので、そのうち・・・だれか効果的なアプリケｰションを開発してくれることを夢見ています。現実問題は、いかにPureな正弦波をとりだすか。または、含まれるノイズの影響をどこまで取り除けるかだろうと思います。

瞬時周波数と成分分離

　フーリエ変換使用すると、「成分」を分離することができます。 Synchro PRIMOでは、少ないサンプル数で「単一成分」を精密に計算できますが、「複数成分」の同時計測は可能なのでしょうか。

　以前、２成分でシミュレーションを行い理論式も導出してみましたが、意味がさっぱりわかりません。しかし最近になってひとつ重要なことに気づきました。　

　本ブログで紹介しましたが、「２成分の混合波形」は、それぞれの周波数、振幅によって定義される「瞬時周波数」の振る舞いを呈します。この計算式は L.O.コーエンの文献にも解説されていることを紹介しました。　逆に言うと、２成分ということがわかっているなら、瞬時周波数の振る舞いから、それぞれの成分の周波数・振幅を計算できることになります。それぞれの成分の振幅比が大きい場合（＝どちらかの成分が優位）瞬時周波数は脈動し、脈動する周期の逆数が「周波数差」として、中心周波数からの変動量で、おそらく・・・振幅が計算できるはずです。

　では３成分以上になるとどうか・・・数学的には、未知変数の数＋１のデータが要求されます。成分数が無限になると・・・必要なデータ量も増えていきます。実現のためには「時間窓」を大きくし、瞬時性を犠牲にするか、さらに高精度のサンプリングを行うしか方法はありません。　

　ここで現実的な実装を考えてみます。サンプリング周波数を高くしてみましょう。隣接サンプルとの値の差がどんどんなくなり、振幅方向のさらに高い分解能が必要になります。振幅方向の分解能が高いと、今度はノイズが問題になってきます。

　以上をまとめると、成分数（＝帯域？）、瞬時性、ノイズの３つの制約が存在し、「超えられない壁」が存在するような予感がします。　

　実は「時間・周波数・ノイズ」の制約は、なにかシャノンの定理に絡んでいるな・・・と博士課程の段階で気づいていました。シャノンの「通信路容量」は、C = B* log(1+SNR) となります（＝シャノン・ハートレーの定理）　。

　Synchro PRIMO でもシャノンの定理が見え隠れしているのです。（余談：この関連について博士論文の審査の際、別の先生からコメントがあったらしいです）

Synchro PRIMO 計算式が示唆する「周波数」

周波数は普通は０以上の実数です。しかし周波数変換・変復調をおこなうと計算上「負の周波数」が登場します。

では、複素数になったらどうなるか？　一応「複素周波数」というものはあり z = ρ + j*ωt  のように表現できます。実際の波形には１）振動項、２）減衰項が存在します。　

　波形を、f(t) = 1/2 { exp(z) - exp(z) } と考えると振動項・減衰項もあらわれ、なんとなく見えてきます。ただし、この式のzをどう呼べばいいのでしょうか。とりあえず「複素位相角」としておきます。

 Synchro PRIMO では、指数的に減衰しながら振動する波形の周波数は正しく求めることができます。上記の複素周波数がうまく処理できているからです。しかし、一般の振幅変調波形ですと、一定の値が計算できません。さあここからどうするか。

　振幅変調された信号の基本周波数は一定なのでしょうか？　ここは直感的な理解と複素解析から得られる結果が一致しないことも覚悟して検討しなくてはなりません。　上記の「複素位相位相角 z」 の微分が「複素周波数」であるとすればよさそうです。そう仮定すると、σ＝一定の減衰をともなう振動では、瞬時周波数の計算には影響しないことがなんとなく分かります。

　振幅変調波形の瞬時周波数を複素解析的にもとめる前に複素平面で位相角の動きをイメージしてみます。　複素位相 z は普通は虚数軸方向に「速度ωで走って」いますが、振幅変調がはいると、実数軸方向に「ぶれる」ことになり、瞬時的な速度＝瞬時周波数は「一定でない」と考えられます。ここが盲点で、我々が無線工学でイメージしている周波数と変調された信号の瞬時周波数は厳密には異なるということが示唆されているのです。　

　Synchro PRIMOの計測（計算）では、指数減衰する三角関数の周波数は厳密に求めることができます。（指数形式で入力信号を与えると、計算式のLsを展開することで容易に証明できます）

　では正弦波で振幅変調された波形の瞬時周波数はどうなるか。log (1+m*sin(ωt)) (mは変調指数）に比例した変動ではないかと考えています。

 

大学院時代の思い出２

体力的・気力的に一番きつかったこと

学位審査にはいるためには、１）授業の単位取得、２）国際会議での発表（少なくとも採択がきまっている）、３）論文誌掲載　が必要でした。

この１）〜３）の勝負所がＤ３の３月〜9月に一気にきてしまって、当時のスケジュールを見ると今でも目眩がします。博士課程にいて最も大変とされているのが論文誌への掲載。査読付きでないといけません。国際会議も同様で、酷い内容だと一発でRejectを食らいます。私のいた大学院では、授業の単位の一つに「副テーマ」というのがあってこれも担当教官に提出し合格しないといけません。あと３月の音響学会で発表した内容で、ちょっとした賞をいただいてしまったため、9月の研究発表会のエントリーも余儀なくされました。

　物事、調子よく進んでいるときは少々の困難はなんとかなるものです。 

すべてを、　I MUST DO otherwise....  と考えるのではなく、　I WILL DO, then !  と考えれば少々のミス・想定外なことは前向きに乗り切れました。

　私の場合、メンタルな意味で苦しんだ記憶はないのですが、研究の最終アウトプットをどういう形で示すのかは悩みました。D2になるころ、「C++ソフトウェアで実現可能、PCで実行可能な水準の精密周波数解析方法」という枠で進めることを決心しましたが、当時採用していた計算式には近似が含まれており、ここがみっともないと指摘され、ぐぬぬ・・・　約半年後、近似を用いない計算式が判明しジャーナル掲載へ一歩近づきました。しかし、この時点では、現在のSynchro PRIMO法のヒントとなる「ある数式」には気づいていなかったのです。

 

大学院時代の思い出

　本ブログで解説している方法は大学院で研究したものを修了後、さらに発展させたものです。しかし、その着想はどこからきたのかと質問されたなら、その説明は少々ややこしくなります。

　実は大学院に行ったのは２０１１年〜２０１４年。私の生まれは昭和３０年代ですので・・・年齢的な状況はお察しのとおりかと思います。よく聞かれる質問は、なんでそんなに勉強したいのか？　なのですが、これも学部卒(1985)以降のことを説明しないとなかなか理解してもらえません。

　実は、学部４年の時 PC-9801で稼働するFFT (8086アセンブラー）を自作し卒業研究をしたのですが、このときにフーリエ変換の不確定性を知りました。（当時つかった教科書は手元にあります）当時、卒論のテーマの裏で今でいう「楽音分離」をやってみたく、PC-9801 のBASIC言語をつかって遊んでいました。不確定性に気付いたときショックだった記憶があります。最終的の卒論は「１／ｆノイズのスペクトル」を自作FFTでシミュレーションしたのですが、結局、修士課程には進学せず、この疑問はずっとしまい込むことになってしまいました。

 

 

雑音の周波数

　２つの、異なる周波数成分を足し合わせても、「平均周波数」は振幅の大きい側に支配されることは以前述べました。今度は、成分が無限に含まれているとされる「雑音」について、周波数はどう表現できるかを考えてみたいと思います。

　私がいま考えている方法は、現波形と１サンプルずらした波形をリサージュ図に描画したとし、その「回転数」をみる方法です。厳密には確率的信号処理の考えで検討しなければならないのですが、ノイズのスペクトルの包絡線の違いでさまざまな計算値がでると思います。

　複数の周波数成分が「入り乱れている」信号の周波数領域の解析において伝統的にF0が重要視されますが、調波の構造が比較的単純でも「ミッシング・ファンフダメンタル」のような現象がおこると、高調波をすべて計測し、その構造からF0を推定するという厄介な作業となります。

　今回提唱している「平均周波数」（名前を工夫して、代表周波数みたいなものも可能）の定義と数理が明確になれば、いちいちフーリエ変換しなくても、かなりのことがわかるのではないかと思います。

　フーリエ変換は便利な手法であるため、広い領域で使われています。時間の区間が±無限であれば）任意の波形は成分に分けられる（＝順変換）ということはたぶん成立するのでしょうが、逆に、任意の成分を用意すれば任意の波形は合成できる（逆変換）は　微分不能な波形を生成することがあります（リーマン関数など）。連続波形なおのに導関数が定まらないというケースです。　どう解釈してよいのか今の私には分かりません。

　私見ではありますが、フーリエ変換の原点にかえって「正しい使い方」をもっと知ってもらいたいなと思っています。ＦＦＴが考案されて以降、濫用されているような気がしてなりません。

連続スペクトルなるものの考察

「連続スペクトル」「スペクトル包絡」などという言葉がありますが、さらに意味を考えてみます。

　連続系の信号で、明らかに連続スペクトルになるものには、一例として、「インパルス」「ガウス関数（ガウシアン）」があります。インパルスの時間幅はゼロ・振幅は無限大・積分すると１という定義はともあれ、現実的には、一定の時間幅をもった矩形波で考えてもいいでしょう。

　では、インパルスは平坦な連続スペクトルなので、周波数がすこしずつ異なるサイン波を無限に重ねればインパルスを合成できそうです。。次の思考実験ですが、パルスをつかって三角関数は合成できるでしょうか。離散系では、sin関数であっても、δ関数の和として表現できます。ここで注意ですが、δ関数のスペクトルは平坦。そんなδ関数をつかって合成すると、１本の線スペクトルになります。この理由は、線スペクトルといっても、各周波数で位相角はきまっていて、位相まで含めて合計すると１本の線になってしまうということです。δ関数の連続スペクトルを合成するとさまざななスペクトルが表現可能、こういう見方も可能です。

　今度はガウス関数。フーリエ変換した結果もガウス関数という、面白い性質をもっています。　ガウス関数をつかって任意の波形が合成できるなら、フーリエ変換も同じような重ね合わせとなってしまうでしょう。この方法は・・・考え中です。

  現時点で悟っていることは、「厳密な連続スペクトルはパルスでしか得られない」、「線スペクトルはいくつ足し合わせても線スペクトルの塊」・・・有限和と無限和の境目がむずかしいです。

 

 

無限和による表現と無限積による表現

　ちょっとした思考実験です。　フーリエ級数（変換）のように、任意の波形が「基底波形」の無限和で表現できるのはご存じのとおりですが、無限積ではきないか。

　まずは有限和で一例を。　f(t) = sin(ω1 t),  g(t) = sin(ω2 t )  として　 y = f + g を考えてみます。三角関数の和積変換を用いると、 y = 2 * sin ( (ω1+ω2)/2) * cos ((ω1-ω2)/2) と表現できます。 ２成分で振幅が等しい場合は、和を積の形に変形できます。

では、振幅が異なる場合；f(t) = sin(ω1 t),  g(t) = 2*sin(ω2 t ) ならどうするか？

y = (f + (1/2 g)) + (1/2)g 　とすると、 y = 2 * sin ( (ω1+ω2)/2) * cos ((ω1-ω2)/2) +sin(ω2 t ) となります。振幅変調のかかった成分と、ω２成分の残骸。　さらにこの２つを統合できないものでしょうか。

　次のケース。振幅はすべて等しいとして、３成分の合成ではどうなるか。 y = sin A + sin B + sin C のようなケース。高校でならった和積の３成分への拡張です。実は数学オリンピッククラスでは当たり前の公式だそうです。　https://mathtrain.jp/waseki　ただし A+B+C = π　の条件があります。この制約を取っ払うのも一苦労。３成分ができれば、Ｎ成分への拡張も・・・根性いれればできるでしょう。

　あとは、振幅が違い周波数が同じ場合ですが、合成公式：https://mathtrain.jp/asinbcos　を使えば１つに集約可能。振幅が違っても、うまく合成公式を活用して、すべての成分の振幅をそろえて、Ｎ成分の「和積」をつかえば目的のものは完成です。　その中ででてくる各項の周波数は何を意味しているのか。　我々が知ってる級数の和と同じことを表しているわけですから、その解釈に、いままで築かなかった「意味」があるのかもしれません。

今時点ではっきりしているのは、振幅の等しい２つの成分の和は、べつの２成分の積で表すことができること。ここをとっかかりに、なにかをブレークスルーしたいと思っています。