hirax.net::Keywords::「周波数解析」のブログ

■Wveletで周期解析をしてみる　　

　音声・地震などの1次元信号や、画像等の2次元信号処理の解析というのはなかなか面白そうだ。そこで、周期ムラに対してWaveletをかけて周波数解析をする練習をやってみたい。また、短時間フーリエ変換とWaveletの比較もしてみたい。音声・地震などのデータはまた別にやってみることにして、今回は画像データを扱うことにする。ただし、いきなり2次元も何なので、画像データの周期(つまり1次元的な振動)に注目して、解析を行ってみたい。

　まずは、「周期ムラのある画像」と「周期ムラのない画像」の2種類の画像を作成する。画像はいずれも数式を用いて作成した。X方向に変化する縞模様であり、表.1のような演算式になっている。一応、2次元画像ではあるが、Y方向にはなんの変化もない。2つの数式を見比べてみると判るが、いずれも2項からなり、低周波数のSinと高周波数のSinからなっている。「周波数ムラのある画像」では、その低周波数のSinの中にさらにSinがあるので、周波数がある周期で変化していることになる。一見、「周波数ムラのない画像」の方でも低周波数のSinの内部にさらにSinがあるように見えるが、0が掛けられているので、実際には存在しないのと同じである。

表.1:画像を作成するために使用した演算式

2つのSinからなり、その一方のSinの周期がムラ(一定の周期)をもっているもの。第二項目のSinの内部にさらにSinを入れることにより、周波数ムラを作っている。	2つのSin波からなり、どちらの周期も正確なもの。第二項目のSinの中のSinは0をかけてあるので、何ら影響を及ぼさない。

　そのような数式に基づいて作成した画像を図.1に示す。なお縦軸がX軸であり、横軸がY軸である。図.2(b)では周波数ムラはないが、2つの周波数成分から作成されているため、うねりが生じている。

図.1:作成した原画像

周波数ムラのある画像	周波数ムラのない画像
図.1(a):	図.1(b):

　それでは、このような画像からX軸の方向に1次元データを抽出し、周波数解析をしてみる。Y軸方向にはなんの変化もないため、無視して良い。

　まずは、Wavelet変換である。図.2がその結果である。縦軸が周波数を示している。縦軸の上方向が高周波を示し、下方向が低周波を示している。また、横軸が原画像のX方向である。白は強度が小さいことを示し、黒は強度が強いことを示している。
　いずれの画像も2つの周波数成分からなることが一目瞭然である。また、図.2(a):「周波数ムラのある画像」の方では低い周波数成分の方が、さらにある周期で周波数が変化していることがわかる。

図.2:Wavelet解析を行ったもの (Daubechiesの6次のCoifletFilterを使用)

周波数ムラのある画像	周波数ムラのない画像
図.2(a):	図.2(b):

　同じWavelet変換でも異なるFilterを用いてみると、結果は異なる。例えば、図.3がその例である。こちらの方が「周波数ムラ」がどのように生じているかを見るにはいいかもしれない。

図.3:Wavelet解析を行ったもの (Daubechiesの6次のLeastAsymmetricFilteｒを使用)

周波数ムラのある画像	周波数ムラのない画像
図.3(a):	図.3(b):

　それでは、Wavelet変換ではなくて、フーリエ変換を用いて周波数解析を行ってみる。先ほどの1次元データの全領域に対してフーリエ変換をかけてみる。その結果が図.4である。ここで、横軸が周波数を示し、右側が高周波数を示し、左側が低周波数を示している。縦軸は強度である。
　このフーリエ変換の場合も、2つの画像が2つの周波数成分からなり、図.4(a):「周波数ムラのある画像」では低周波数成分がぶれているのはわかる。しかし、その周波数ブレがどのようなものであるかまでは、わからない。

図.4:全領域にFFTをかけて、周波数解析を行ったもの

周波数ムラのある画像	周波数ムラのない画像
図.4(a):	図.4(b):

　それでは、短時間フーリエ変換をかけてみる。先ほどの1次元データに対して前の方から64点ずつ、位置をずらしながらフーリエ変換を行う。このようにすることによって、ある領域の周波数解析を行うことができる。その結果を図.5に示す。ここで、黒は強度が小さいことを示し、白は強度が大きいことを示している。横軸は原画像のX方向を示し、縦軸が周波数を示している。縦軸の上方向が高周波数を示し、下側が低い周波数を示している。結果はWaveletの解析と同様になっている。
　なお、単純のためにウィンドー処理はしていない。そのために悪影響は当然出てしまう。

図.5:短時間FFTをかけたもの

周波数ムラのある画像	周波数ムラのない画像
図.5(a):	図.5(b):

　単なる全領域にわたった周波数解析と、位置と周波数が同時にわかる解析の違いは非常に大きい。使いこなすのはなかなか難しそうだが....

■オシロスコープソフトを作る　

PCを2Ch高性能オシロスコープにしたい

　オシロスコープがあると便利だが、家で使うには敷居が高いし、値段も高い。まして、FFTアナライザーがついて周波数解析も行うことができるような機械になると、遊びで買うという値段ではなくなってしまう。そこで、PCを2Ch高性能オシロスコープにするソフトを作ってみたい。

　以前、音階を調べた時に、SoundBlaster互換I/Fを使ったLabViewのサウンド入出力のViを使用してみる。目標はとにかくオシロスコープと同じ使い勝手であること、使うのが簡単であること、そして、周波数解析などが簡単に行うことができること、である。入力として、音声入力を使用しているので、たいていのPCで使うことができるし、音声入力マイクが着いているPC(たいていのノートPCは着いているだろう)なら、マイク(あるいは入力端子)を用意する必要すらない。

　というわけで、下がそのアプリケーションの画面。

作ったアプリケーションの画面 (口笛を吹いて、周波数解析を行ったところ)

　22kHzで2Ch(ステレオ)サンプリングを行っている。また、リアルタイムで任意のウィンドー処理をかけた上で、FFTをかけて周波数解析を行い、ピーク周波数の検出、表示を行うことができる。オシロスコープの表地周りは通常のオシロスコープとほぼ同じであり、使うのも簡単である。

　20KHz(ナイキスト周波数で言えば正確には10KHzか)までしか、使うことはできないが、ちょっと使いたい時には便利だ。特に、音声を解析したいならば、必要十分である。上の画面では口笛を吹いて、その音声波形を表示させ、周波数解析を同時に行っている。周波数ピークが表示され、1120Hzであるという表示がされる。

　内蔵マイクを使用すると、ただアプリケーションを走らせれば、音声周波数解析が行える。もしも、比較的高性能なマイクがPCに着いているならば、リアルタイム振動解析すら行うことができる。もちろん、マイク入力端子に何らかの入力をすれば、どんな解析もできるわけだが、何の用意もせずにできるというのは便利である。例えば、うるさいデスクトップPCの近くへこのアプリケーションを走らせたノートPCを近づけると、デスクトップPCがなぜうるさいかを簡単に調べることができる。うるさいのは、ハードディスクの周波数なのか、ファンの周波数なのかすぐにわかる。
　
　今回、作成したアプリケーションはここにおいておく。動作させたら、レンジを調整することを忘れずに。
Ocilo.lzh LZH形式　1,259kB (打ち止めです。あしからず。)
LabViewのアプリケーションライセンス上、ダウンロード数は50回までで、その数近くになったら削除することにする。

■画像ノイズ解析について考える　

考える理由

　画像ノイズ解析を目的として、2次元フーリエ変換を用いて周波数解析をすることが多い。かねがね、このやり方について疑問を感じていたので少し考えてみたい。

　その疑問とは次のようなことである。

• 通常の2D-FTでは、入力データ全領域での周波数解析を行う。従って、単発のパルスのようなノイズはバックグラウンドに埋もれてしまい、結果にはなかなか出てこない。
• 同じ理由で、2D-FTでは位置と周波数解析を同時に行うことができない。(もちろん、短時間フーリエ関数を使えば、そのような測定は行うことができる。)
• また、ホワイトノイズのようなフラットな周波数特性を持つノイズもバックグラウンドを押し上げるだけの効果しか持たないため、解析をしづらい。

　そこで、今回は単純な画像に対して、2D-FTと2D-離散Waveletの比較を行うことで、2D-FTを用いた「画像ノイズ解析」の問題について考える。

2D-FTと2D-Waveletの例

　はじめに、2D-FTと2DWaveletの例を挙げる。まずは2D-FTである。

2D-FTの例(左から原画像、2D-FT結果、2D-FT結果の鳥瞰図)

　左の原画像は45度のスクリーン角のラインである。2DFTの結果にはその角度方向にピークがいくつか並んでいる。それぞれのピークの中央からの距離が周波数を示している。それはX,Y方向いずれについても言える。今回の場合はX,Y方向のスクリーンの周期が等しいため、2DFTの結果でも45度方向になっているのである。
このように、2D-FTの結果というのは周波数(X,Y両方向)と振幅がわかる。ここでのスクリーン角のような周期性を持つものの解析にはフーリエ解析というのは極めて有効である。店で見かけるインクジェットプリンターもヘッドの移動による周期ムラが激しいが、このようなムラに対してフーリエ変換を用いた周波数解析を行うのは正当であり、有効だろう。

　それでは、同じ画像に2D-Waveletをかけてみる。2D-Waveletの結果は位置と周波数強度分布情報(ホントは違うのだが)が両方出てくる。位置情報が2次元で周波数強度分布情報が1次元であるから、合わせて3次元である。そのため、表示に一工夫いる。
　第一段階として高周波成分から調べてみる。すぐにこの結果の意味がわかるだろうか？

2D-Wavelet例(左が原画像、右が一段階Waveletをかけた結果)

かなり判りづらい。この右の結果は4つの領域にわかれているが、以下の表のような意味を表している。また、いずれも灰色の部分は強度が弱く、白と黒が強度が強いことを示している。

高周波のX成分	高周波成分
低周波成分	高周波のY成分

　低周波成分が原画像と同じようであるのがわかると思う。これは2DFTと違い、Waveletでは位置情報もそのまま保持されているからである。次に、この低周波成分に対して、もう一段Waveletをかけるとこうなる。　右上から左下への対角線上のが周波数成分を示し、これで周波数成分にして3分解できたことになる。右上が一番高周波成分。その左下が次の高周波成分。右下が低周波成分である。
　もう何分割かしてみる　このようにして、画像内での位置と周波数成分が両方ともわかる。

　なお、フーリエ変換では基底関数としてSinが用いられるが、Wavelet変換では基底関数としていろいろな関数を使うことができる。今回はDaubechiesの4次のものを用いている。下がその形である。

Daubechiesの4次のフィルター

ドットのノイズを解析してみる

　それでは、今回の本題に入る。以下が原画像である。左が「2つの大きなドットからなる」画像であり、右がそれにノイズの加わった「ノイズ」画像である。ここでノイズはホワイトノイズを加えているつもりである。ドットは周期性を持つデータだが、ノイズ自体は周期性を持たない所がミソである。また、ここで言う「ノイズ」とは現実の現象とは何ら関係がない。単なる例えである。

ドット画像(左が原画像、右がノイズを加えた画像)

　まず、この2つの画像に対して、それぞれ2D-FTをかける。

２次元離散フーリエ変換を行った結果

　このグラフではXY軸とも-πからπまでの領域で示している。中央からの距離が周波数を示しており、明るいほどその周波数帯の振幅が大きいことを示している。つまり、任意の周波数帯の強度がわかる。
　右のノイズの加わった画像の2DFTの結果では、広い周波数領域で強度が上がっている。しかし、下の鳥瞰図で示した(私は立体が好きなのだ)方でもわかると思うが、バックグラウンドが持ち上がっているだけである。いずれにせよ、あまり左右の間で違いはない。今回のような64x64の画像ではなく、もっと大きい画像ではその違いははより識別不能になる。

２次元離散フーリエ変換の結果を鳥瞰図で示したもの

　さて、次に2D-Waveletで同じように計算をしてみる。下が計算結果である。どうだろうか、ノイズ(位置も周波数も)が一目で判るようには思えないだろうか？

2D-Waveletによる解析結果(左がノイズ無し、右がノイズ有り)

　今回は、自分の頭を整理するために、ただ2D-wavelet変換をかけてみた。まだまだ話しは続くのである。

■カードサイズの「画像探偵セット」!?　

お手軽線数メーターを作るのだ

　「どんなものでも、自分の目で眺めてみた〜い」と、ワタシはいつでも思う。世の中スベテのものを、自分の目で眺めてみた〜いと思う。しょんべん小僧が空中に描き出す放物線、巨乳ギャルにロックオンするオッパイ星人の目の動き、ビデオにかかるモザイクの向こう、はたまた田代まさしが恋い焦がれるミニスカートの結界の秘密、とにかく世の中のものスベテを何でもかんでも眺めてみたい、覗いてみたい、とワタシはいつも思っているのである。（とはいえ、誤解されると困るので念のために書いておくが、もちろんミニスカートの中を覗いたりはしないのだ）

　そんなわけで、ワタシのケータイのストラップには「ちっちゃなちっちゃな虫メガネ」がついている。この虫メガネを武器にして、ワタシは色んなモノを覗くのがクセになっている。スーツ姿で出張している時だって、おもむろにこの虫眼鏡を取り出して、色んなものを覗いてみたりしているのである。
 

ワタシのケータイのストラップには「ちっちゃなちっちゃな虫メガネ」がついている

　だから、毎朝届けられる新聞に折り込まれているチラシやカタログを眺めるときだって、そんなカタログに「ちっちゃなちっちゃな虫メガネ」を向けてみて、その「虫メガネ」を通して、カタログがどんな風に印刷されているかをよく眺めてみる。下の左のようなカラーの綺麗なカタログだって、「虫メガネ」を通して眺めてみると、右の拡大写真みたいに、四色（シアン、マゼンダ、イエロー、ブラック）が規則正しく並んでいるようすが見えてくるのである。離れてみればキレイな写真が「虫メガネ」を通して眺めてみるだけで、こんな風に様子が変わるなんてとても不思議な気分になったりするのである。まるで、女性の化粧のように不思議で、こんな四色の手品はとても面白いのである。
 

「虫メガネ」を通して眺めてみると、
四色（シアン、マゼンダ、イエロー、ブラック）が規則正しく並んでいる

カタログの一部

左の青い矩形部の拡大写真 （左上と右上部分はマゼンダとシアンだけを見やすくしてみた）

　で、こんなカタログの拡大図を眺めていると、四色に分けられた色がどんな風に並べられ形作られているかを、知りたくてたまらなくなったりする。つまり、「どんな間隔でこの色は配置されているのだろう？」とか「一体、どんな角度でこの色は並べられているのだろう？」とか思うわけである。もちろん、この画像に二次元フーリエ変換などをかけさえすれば（周波数解析をすれば）、「どんな角度で・どんな間隔で色が並べられているか」ということは知ることができるけれど、まさかワタシの頭の中でそんな作業ができるわけはない。かといって、このカタログを読み込むための「画像読みとり装置」や「解析をするためのコンピューター」を毎日持ち歩くなんてこともできるわけもない。

　そこで、「どんな間隔でこの色は配置されているのだろう？」とか「一体、どんな角度でこの色は並べられているのだろう？」とかいう疑問の答えがすぐ判るように、先日こんなカード、ぴったりクレジットカードサイズの透明シートに規則的なパターンを印刷したカード、を作ってみた。名付けて、Peco-Chartなのである。
　これはもちろん、判る人には判るだろうが、ハンディ「線数メーター」というモノである。一言で言えば、画像の周波数解析をとっても簡単にすることができる手品の小道具のようなカードなのだ。
 

Peco-Chart

　例えば、さっきのカタログの上にこのPeco-Chartを重ねて置いてみると、アラ不思議、何やらヘンな模様、不思議なモアレ模様が浮かび上がってくる。下の左の写真、あるいは右の拡大写真を眺めてみればマゼンダとシアンの同心円がハッキリと浮かび上がっているのが見えるだろう。例えば、シアン色の場合は75°の角度で175線（175線/inch）位の位置、そしてマゼンダ色の場合は45°の角度で同じく175線位の位置を中心として、同心円状のモアレ模様が浮かび上がっている。

　つまり、「このカタログはシアンは75°の角度方向に1インチ辺り175個のドットが並べられていて、マゼンダのドットは45°の角度方向に並べられている」、ということを、このPeco-Chartを重ねて置いてみさえすればたちどころに知ることができるのである。
 

さっきのカタログにPeco-Chartを重ねてみると…!?

何やらモアレが見える	左の写真の拡大図
	カタログの拡大写真

　こんなペラペラのカードで周波数解析ができるなんてとても不思議に思えたりもするけれど、ちょっと考えてみればこれはごく当たり前の話である。

• モアレがタネの科学のふろく- うれし、懐かし、学研の科学 - (2000.05.02)

でも似たようなことをしたように、モアレというものは「二種類以上の何らかの模様（パターン）が干渉して発生する」ものである。つまり、ある意味「二つのパターンの相関をとる」ということである。そしてまた、少し考えてみれば「画像の周波数解析」というものは「対象となる画像」と「基準関数（三角関数etc.）」の間で相関を調べることと同じである。だから、「基準となるパターン」を「対象となる画像」の上に重ねてみた時に見えるモアレのパターンは「対象となる画像」の周波数解析結果を実は示していると考えてみても良いのである。だから、このPeco-Chartはクレジットカードサイズのペラペラなちっぽけなヤツではあるのだけれど、実は色んな画像の周波数解析をしてくれるスゴイヤツだったのである。そして、こんなポケットに入るほど小さい線数メーターははなかなか無いので、、カード入れからコイツを華麗に取り出してみせたりすると、うらやましがる人もとても多く（仕事柄、画像出力に関わる人達が多いから）、なかなかに気持ちが良いのである。

　とはいえ、自慢してばかりでは何なので、さらに大量に配布すべく新たなPeco-Chart二号機をデザインしてみた。それが、下の名付けてPecochart^proである。
 

Pecochart^pro

このPecochart^proの謳い文句はその名の通り「プロ仕様」というわけで、画像出力に関わっている何人ものベータテスター達（自分も含めて）の感想をもとにして、

1. 70線から350線までの線数とスクリーン角度の測定ができる「線数・角度メーター」（分解能を2線単位から1線単位へと二倍向上、スクリーン角度のガイドは2.5度刻み）
2. 90線から410線までの線数を高精度（0.25線刻み）に計測することができる「線数メーター」（新機能）
3. 線・文字の太さを計測できる「線幅スケール」（新機能）
4. 〜30級、〜30ポイントまでの文字サイズ（級数、ポイント）測定ができる「文字スケール」（Peco-Chartと同じ）
5. 「8cm定規」（Peco-Chartと同じものを使い勝手はそのままにコンパクト化）
6. そして、便利な「画像に関する換算表」（内容を従来比75%増量）

をこれまでと同じクレジットカードサイズに満載してみた自慢の品なのである。とはいえ、こちらのPecochart^proの方はまだデザインしてみただけで、フィルムへの焼き付けなどはまだしていない。だから、もしこのPecochart^proを欲しい〜とか、こんな機能が欲しい〜とか、印刷屋さんのノベルティグッズに使いたい〜とか思うような人がもしもいるならば、ワタシ(jun@hirax.net)までぜひ連絡して欲しいと思う。とりあえず、画像出力機器に携わる人には必携の自慢ツールになること間違いなし、と勝手に思っていたりするのである。

　それにしても、何か自分に役に立つツールを作るというのは本当に楽しい作業だった。手作りツールぎゃらりい脇色彩研究所ではないけれど、こんな「hirax.netオリジナルグッズ」をいっぱい作って「探偵セット」ならぬ「できるかな?セット」として、面白メールをくれた方にプレゼントとかしてみたら楽しいんだろうなぁ、と思うのである。というわけで、そんなツールをせっせと作るのだぁ、なんてことを実は計画中なのでした、ハイ。
 

■アダルトビデオの科学 今日のトリビア 編　

　「帰ってきたオッパイ星人の力学 アダルトビデオ鑑賞 編(仮題)」改め、「アダルトビデオの科学 今日のトリビア編 - 横揺れ縦揺れ2.5ヘルツ -」です。

アダルトビデオの振動方向はジャンルによって違うが、一秒間に2.5往復ほどである。

　ブロックマッチング、モルレーのウェーブレットで時系列周波数解析、そんなこんなでエロエロを解析。えっ、こんなテクノロジーを使っても解析する対象を間違えているんでしょうか？しかし、そんな科学や新しいテクノロジーを牽引している原動力も（時には）エロエロだったりするわけです。

　気づけば今年初の「できるかな?」、しかも聖バレンタインな日にこともあろうに「アダルトビデオの科学」　うぅ。それでも…「それもまた良し」ということにしておきましょうか。