| | HOME | 目次 | |
2009 (C) H.Ishikawa |
７．信号処理−−雑音の中から信号を拾い出す−−
通信の分野では、伝送路や電子管などから発生するランダムな雑音により、妨害をうけて、正確に信号が受信されなくなることがあります。この章では、雑音をコンピュータで作りだし、雑音の中から信号を分離する実験をやってみましょう。これは、最近ディジタル信号処理といって、発展の著しい分野です。
７．１ ショットノイズのシミュレーション (文献（４）参照)
ショットノイズとは
ショットノイズは真空管から発生する雑音です。真空管は、図７．１．１のように、陰極と陽極があり、陰極から陽極へ電子が飛ぶことにより、電流が流れ、間にあるグリッドに電位をかけることにより、この電流を制御する構造になっています。電流は連続的に流れるように思われますが、電子という”つぶつぶなもの”が、移動するため、詳細に観察すると、陰極から電子が不規則に放出しているので、電流が I+in のように変動します。この陽極電流の動揺成分 in をショットノイズといいます。このノイズは次のようのなメカニズムにより発生すると考えられています。陰極から陽極にむかって電子１個が放出されると、その瞬間から、陽極へ靜電誘導によって、回路に電流が流れ始め、電子が陽極へ到達したあと、この電流は、消滅します。図７．１．２のように、このパルス状の電流が存続している間に（これは極めて短い時間ですが）、十分多くの電子が到達するので、その結果として、陽極電流が連続したものになります。この際、電子の発生はランダムですから、陽極電流が動揺することになるのです。
| 図７．１．１ ショットノイズ |
 |
| 図７．１．２ 拡大した陽極電流 |
 |
ショットノイズのシミュレーション
プログラムS0710.javaは、ショットノイズを、シミュレーションするプログラムです。電子１個による電流は、急竣なパルス状となります。それを、幅 P_WIDTH 、高さ１の矩形パルスとします。電子は、３．２節で学んだポアソン分布により発生すると考えられ、単位時間に P=0.1 という小さい確率とします。(1) の for ループは、時刻を１づつ進めており、そのたびに (2) のように、 P という確率で電子が飛び出します。その時幅 P_WIDTH の時間にわたって、 i[tt] という配列に１を加え、それが電子１個による電流波形になるわけです。 T_END になるまで、シミュレーションを続けます。
つぎにこの電流波形を図示しますが、その前に (3) によりヒストグラムを作成し、また平均値、分散の計算をしておきます(注)。このとき、時刻が 0 から P_WIDTH までの間は過渡状態ですから、計算からのぞきます。
ウインドウ０に電流波形を、ウインドウ１に電流の分布を示します。図７．１．３の実行結果でわかるように、ショットノイズは、歩幅１、歩く時間が不定のランダムウオークで、その分布は正規分布をしていることがわかります。
| 図７．１．３ ショットノイズのシミュレーション結果 |
 |
| S0710.java ショットノイズのシミュレーション | Download | 実行 | |
/* S0710.java
* ショットノイズのシミュレーション
* (C) H.Ishikawa 2008
*/
package simulation;
import java.applet.*;
import java.awt.*;
import java.awt.event.*;
import window.Window;
public class S0710 extends Applet implements ActionListener {
Button button0;
Button button1;
boolean sw = false;
public void init() {
button0 = new Button(" 実行 ");
button1 = new Button("クリア");
add(button0);
add(button1);
button0.addActionListener(this);
button1.addActionListener(this);
}
public void actionPerformed(ActionEvent e) {
String label = e.getActionCommand();
if (label.equals(" 実行 "))
sw = true;
else
sw = false;
repaint();
}
public void paint(Graphics g) {
Window w;
w = new Window();
int height;
int width;
int SPACE = 30;
int T_END = 1000000;
double P = 0.01;
int P_WIDTH = 10000; /* パルス幅 */
int t、 tt;
int i[] = new int[T_END]; /* 電流 */
int f[] = new int[150]; /* ヒストグラム */
double x = 0.0;
double x2 = 0.0;
int n = 0;
int ii、 j;
/* スクリーンの大きさを読み取る */
width = getSize().width;
height = getSize().height;
/* グラフィックの準備 */
w.setWindow(0、 0.0、 0.0、 (double) T_END、 150.0、 SPACE、 height - SPACE、
width / 2、 SPACE);
w.axis(0、 "t"、 T_END / 10、 "i"、 50.0、 g);
w.moveTo(0、 0.0、 0.0、 g);
w.setWindow(1、 0.0、 0.0、 0.1、 150.0、 width / 2 + SPACE、 height - SPACE、
width - SPACE、 SPACE);
w.axis(1、 "freq"、 0.01、 "i"、 50、 g);
/* メイン */
if (sw == true) {
w.moveTo(1、 0.0、 0.0、 g);
for (t = 0; t < T_END; t++) { //(1)
if (P > Math.random()) { //(2)
for (tt = t; tt <= t + P_WIDTH; tt++) {
if (tt < T_END) {
i[tt] = i[tt] + 1;
}
}
}
}
for (t = P_WIDTH; t < T_END; t++) {
f[i[t]] = f[i[t]] + 1; //(3)
x = x + (double) i[t];
x2 = x2 + (double) i[t] * (double) i[t];
n = n + 1;
}
/* 表示 */
for (t = 0; t < T_END; t++) {
g.setColor(Color.green);
w.lineTo(0、 (double) t、 (double) (i[t])、 g);
}
for (ii = 0; ii < 150; ii++) {
g.setColor(Color.blue);
w.line(1、 0.0、 (double) ii、 (double) f[ii] / (double) n、
(double) ii、 g);
}
g.setColor(Color.black);
g.drawString("mean = " + (x / (double) n) + " var = "
+ (x2 - x * x / (double) n) / (double) n、 10、 height);
stop();
}
}
}
|
| この章始め |
７．２ ＣＲフィルタによる信号の分離 (文献（６）参照)
ガウス雑音
前節でのべた真空管の陽極電流は、直流分に正規分布の雑音（これをガウス雑音という）が重量したものでした。このほかに雑音には、熱雑音というものがあります。熱雑音は抵抗から発する雑音で、これも正規分布となります。それは、電流は非常に多数の電子が作る微小電流を合成したもので、それぞれが抵抗の中で、不規則な熱運動をするため、抵抗の両端に、正規分布の雑音電圧が発生するのです。
今、正弦波とこのガウス雑音が重量された波形を、図７．２．１のようにコンデンサCと抵抗RのＣＲローパスフィルタに加えると、どのような波形が得られるか、シミュレーションしてみることにしましょう。
ガウス雑音は、４．５節でのべた正規乱数の関数を用いて作ることができます。
| 図７．２．１ CRフィルタ |
 |
ＣＲフィルタ
図７．２．１において次の微分方程式が成立します。
 |
・・・・・ |
式７．２．１ |
すなわち
 |
・・・・・ |
式７．２．２ |
この微分方程式を解くには、ルンゲクッタ法を用います(注)。
プログラムS0720.javaにおいて、メソッド
double noise(double en);
は、熱雑音電圧 en のガウス雑音を発生させます。また、
double cr(double e、 double v);
は式７．２．２を表わしています。(1) から (2) までで、ルンゲクッタ法により、微分方程式を解いています。
アウトプット
２つのウィンドを用い、ウインドウ０には、正弦波と雑音波の合計を表示し、ウインドウ１にフィルタの出力 v を表示しています。図７．２．２のように、ＣＲフィルタを通った波形は、こまかい雑音が少なく、スムースになっているのがわかります。ＣＲフィルタは典型的なローパスフィルタで、雑音を減らす効果があります。
| 図７．２．２ ＣＲフィルタのシミュレーション結果 |
 |
| S0720.java ＣＲフィルタのシミュレーション | Download | 実行 | |
/* S0720.java
* ＣＲフィルタのシミュレーション
* (C) H.Ishikawa 2008
*/
package simulation;
import java.applet.*;
import java.awt.*;
import java.awt.event.*;
import window.Window;
public class S0720 extends Applet implements ActionListener {
double C = 1.0;
double R = 0.1;
double E0 = 1.0; /* 信号の振幅 */
double OMEGA = 2.0 * Math.PI * 0.1; /* 信号の角周波数 */
double E_NOISE = 0.1; /* 雑音の標準偏差 */
double T_END = 10.0; /* シミュレーションの終了 */
double DT = 0.01; /* シミュレーションのきざみ */
Button button0;
Button button1;
boolean sw = false;
public void init() {
button0 = new Button(" 実行 ");
button1 = new Button("クリア");
add(button0);
add(button1);
button0.addActionListener(this);
button1.addActionListener(this);
}
public void actionPerformed(ActionEvent e) {
String label = e.getActionCommand();
if (label.equals(" 実行 "))
sw = true;
else
sw = false;
repaint();
}
public void paint(Graphics g) {
Window w;
w = new Window();
int height;
int width;
int SPACE = 30;
double e、 v;
double t;
double k1、 k2、 k3、 k4;
/* スクリーンの大きさを読み取る */
width = getSize().width;
height = getSize().height;
/* グラフィックの準備 */
w.setWindow(0、 0.0、 -2.0、 T_END、 2.0、 SPACE、 height / 2 - SPACE / 2、
width - SPACE、 SPACE);
w.axis(0、 "t"、 1.0、 "e"、 0.5、 g);
w.setWindow(1、 0.0、 -2.0、 T_END、 2.0、 SPACE、 height - SPACE、 width
- SPACE、 height / 2 + SPACE / 2);
w.axis(1、 "t"、 1.0、 "v"、 0.5、 g);
w.moveTo(0、 0.0、 0.0、 g);
w.moveTo(1、 0.0、 0.0、 g);
/* メイン */
if (sw == true) {
v = 0.0;
for (t = 0.0; t <= T_END; t = t + DT) {
e = noise(E_NOISE) + E0 * Math.sin(OMEGA * t);
k1 = DT * cr(e、 v); //(1)ここから
k2 = DT * cr(e、 v + k1 / 2.0);
k3 = DT * cr(e、 v + k2 / 2.0);
k4 = DT * cr(e、 v + k3);
v = v + (k1 + 2.0 * k2 + 2.0 * k3 + k4) / 6.0; //(2)ここまでルンゲクッタ
g.setColor(Color.green);
w.lineTo(0、 t、 e、 g);
w.lineTo(1、 t、 v、 g);
}
}
stop();
}
public double cr(double e、 double v) { //回路方程式
return ((e - v) / (C * R));
}
public double noise(double en) { //ノイズ発生
double xw = 0.0;
for (int i = 1; i <= 12; i++) {
xw = xw + Math.random();
}
return (en * (xw - 6.0));
}
}
|
| この章始め |
７．３ 移動平均による信号の分離
(文献（６）参照)
移動平均
前節と同様の正弦波と熱雑音が重量された波形を、こんどは移動平均という方法で処理し、信号を分離することを試みましょう。
移動平均は、統計手法としてよく用いられる方法で、振幅の変化の激しい時系列について、その中の細かい変化を無視し、その大勢をわかりやすくするものです。具体的には、ある測定値を中心として、その前後の一定の個数（ M ）を平均して、その平均化されたものを移動平均といいます。すなわち、図７．３．１のように、x1
、 x2、 ･･･xi の移動平均 ui は、
 |
・・・・・ |
式７．３．１ |
となります。
| 図７．３．１ 移動平均 |
 |
たとえば、 M=3 の場合、移動平均 ui は
 |
・・・・・ |
式７．３．２ |
となり、下線部分に同じ計算が入っています。平均化点数が小さいときには、問題がないのですが、この点数が大きくなると時間ロスが大きくなるので、
 |
・・・・・ |
式７．３．３ |
として、つぎつぎに計算していきます。
プログラムの説明
プログラムS0730.javaについて説明します。(1)式７．３．３の Si-1 をプログラムでは sum と表わしています。(2)また、 xi-M すなわち、 i の時点から M 前の x の値が必要になりますので、長さ M の配列 mem[M] を用意します。(3)この配列に
mem[i % M] = x;
として書き込み
mem[(i - M) % M]
として取り出せば、 xi-M が取り出せます。(4)が式７．３．３のかっこの中の計算をしているわけです。プログラムのそのほかの部分はプログラムS0720.javaと同様です。出力は、図７．３．２のように、ウインドウ０は入力の波形を、ウインドウ１で、移動平均後の波形を示します。
| 図７．３．２ 移動平均による信号の分離 |
 |
| S0730.java 移動平均による信号の分離 | Download | 実行 | |
/* S0730.java
* 移動平均による信号の分離
* (C) H.Ishikawa 2008
*/
package simulation;
import java.applet.*;
import java.awt.*;
import java.awt.event.*;
import window.Window;
public class S0730 extends Applet implements ActionListener {
double X0 = 1.0; /* 信号の振幅 */
double OMEGA = 2.0 * Math.PI * 0.1; /* 信号の角周波数 */
double X_NOISE = 0.1; /* 雑音の標準偏差 */
double T_END = 10.0; /* シミュレーションの終了 */
double DT = 0.01; /* シミュレーションのきざみ */
int M = 11; /* 平均化点数 */
Button button0;
Button button1;
boolean sw = false;
public void init() {
button0 = new Button(" 実行 ");
button1 = new Button("クリア");
add(button0);
add(button1);
button0.addActionListener(this);
button1.addActionListener(this);
}
public void actionPerformed(ActionEvent e) {
String label = e.getActionCommand();
if (label.equals(" 実行 "))
sw = true;
else
sw = false;
repaint();
}
public void paint(Graphics g) {
Window w;
w = new Window();
int height;
int width;
int SPACE = 30;
double x = 0.0; /* 雑音＋信号 */
double u = 0.0; /* 出力 */
double sum = 0.0; /* (1)本文参照 */
double mem[] = new double[M]; /* (2)本文参照 */
int i; /* シミュレーションのクロック */
/* スクリーンの大きさを読み取る */
width = getSize().width;
height = getSize().height;
/* グラフィックの準備 */
w.setWindow(0、 0.0、 -2.0、 T_END、 2.0、 SPACE、 height / 2 - SPACE / 2、
width - SPACE、 SPACE);
w.axis(0、 "t"、 1.0、 "e"、 0.5、 g);
w.setWindow(1、 0.0、 -2.0、 T_END、 2.0、 SPACE、 height - SPACE、 width
- SPACE、 height / 2 + SPACE / 2);
w.axis(1、 "t"、 1.0、 "v"、 0.5、 g);
w.moveTo(0、 0.0、 0.0、 g);
w.moveTo(1、 0.0、 0.0、 g);
/* メイン */
if (sw == true) {
for (i = 0; i <= T_END/DT; i ++) {
x = noise(X_NOISE) + X0 * Math.sin(OMEGA * i * DT);
if (i < M) { /* 0からM-1までの始めの部分 */
sum = sum + x;
mem[i] = x;
} else { /* M以上の場合 */
sum = sum + x - mem[(i - M) % M]; /* (4)本文参照 */
mem[i % M] = x; /* (3)本文参照 */
}
u = sum / M;
g.setColor(Color.green);
w.lineTo(0、 i * DT、 x、 g);
w.lineTo(1、 i * DT、 u、 g);
}
stop();
}
}
double noise(double xn) {
double xw = 0.0;
int i;
for (i = 1; i <= 12; i ++){
xw = xw + Math.random();
}
return(xn * (xw - 6.0));
}
}
|
| この章始め |
７．４ 同期加算
雑音にうもれた信号を取り出す
雑音の方が信号より大きいという悪条件の中で、信号を取り出すにはどうしたらいいのでしょうか。もしも同じ信号が、くり返し発せられるということが、わかっていれば、図７．４．１のように波形を加算し、平均することにより、雑音の影響を少なくすることができます。ただし信号を時系列上で同じ位置に固定して加算しなけれはいけません。これを信号と同期をとるといい、同期加算をよばれるゆえんです。アンサンブル平均ともいいます。
| 図７．４．１ 同期加算の原理 |
 |
同期加算の原理
３．５節の中心極限定理において、いかなる分布の乱数でも、それを m 回加算したあとの分布を見ると、もとの乱数の分散を σ2 とすると、分散が mσ2 の正規乱数に近づくということを知りました。すなわち、標準偏差で言うと、
になるということです。このことは加算したあと加算回数で割ると（すなわち平均をとると）、標準偏差は
となり m 回の加算で、標準偏差が
になることになります。ガウス雑音の平均振幅は標準偏差であらわされますから、ガウス雑音をつぎつぎ m 回加算して平均をとると、その大きさは
になるはずです。一方、信号の方は、タイミングをあわせて、（同期して）加算されるので、その平均振幅は変わりません。したがって、信号と雑音の比（Ｓ／Ｎ比とよばれる）は、加算回数にしたがって、改善されていきます。これが同期加算の原理です。
シミュレーションプログラム
プログラムS0740.javaにより、同期加算をシミュレーションしてみましょう。信号は、メソッド
doube signal(double t、 double peak、 double p_width);
により発生させます。これは高さ peak 、幅 p_width の単パルスです。雑音は例によってガウス雑音を
double noise (double sd);
により、正規乱数として発生させます。信号と雑音の合計を x とし、各クロック i ごとにそれを配列 mem[i] に加えて平均を計算していきます。
ディスプレイ上には、４つのウインドウ（ウインドウ番号 j＝0、 1、 2、 3 ）を用意し、それぞれに加算回数 ０、３、１５、６３に対応させ、その回数に達したそきに表示させることとします。
| S0740.java 同期加算による信号の分離 | Download | 実行 | |
/* S0740.java
* 同期加算による信号の分離
* (C) H.Ishikawa 2008
*/
package simulation;
import java.applet.*;
import java.awt.*;
import java.awt.event.*;
import window.Window;
public class S0740 extends Applet implements ActionListener {
double X0 = 0.5; /* 信号のピーク値 */
double W = 0.1; /* 信号の幅 */
double NOISE = 0.50; /* 雑音の標準偏差 */
double T_END = 10.0; /* シミュレーションの終了時刻 */
double DT = 0.01; /* シミュレーションのキザミ */
int M = 1000; /* memの大きさ T_END/DT */
int N = 4; /* 加算回数が0、3、15、63のときに表示 */
int N_END = 64 ; /* 加算回数の最大 */
Button button0;
Button button1;
boolean sw = false;
public void init() {
button0 = new Button(" 実行 ");
button1 = new Button("クリア");
add(button0);
add(button1);
button0.addActionListener(this);
button1.addActionListener(this);
}
public void actionPerformed(ActionEvent e) {
String label = e.getActionCommand();
if (label.equals(" 実行 "))
sw = true;
else
sw = false;
repaint();
}
public void paint(Graphics g) {
Window w;
w = new Window();
int height;
int width;
int SPACE = 30;
int m; /* 加算回数 */
int mm = 1; /* 表示する加算回数 */
int j = 0; /* 表示するウィンドの番号 */
double x = 0.0; /* 信号+雑音 */
double mem[] = new double [M]; /* 加算メモリ */
int i;
/* スクリーンの大きさを読み取る */
width = getSize().width;
height = getSize().height;
/* グラフィックの準備 */
w.setWindow(0、 0.0、-2.0、T_END、2.0、
SPACE、height/4-SPACE、width-SPACE、0);
w.axis(0、 ""、 1.0、 "m=0"、 2.0、 g);
w.moveTo(0、 0.0、 0.0、 g);
w.setWindow(1、 0.0、-2.0、T_END、2.0、
SPACE、height/2-SPACE、width-SPACE、height/4);
w.axis(1、 ""、 1.0、 "m=3"、 2.0、 g);
w.moveTo(1、 0.0、 0.0、 g);
w.setWindow(2、 0.0、-2.0、T_END、2.0、
SPACE、height/4*3-SPACE、width-SPACE、height/2);
w.axis(2、 ""、 1.0、 "m=15"、 2.0、 g);
w.moveTo(2、 0.0、 0.0、 g);
w.setWindow(3、 0.0、-2.0、T_END、2.0、
SPACE、height-SPACE、width-SPACE、height/4*3);
w.axis(3、 ""、 1.0、 "m=63"、 2.0、 g);
w.moveTo(3、 0.0、 0.0、 g);
/* メイン */
if (sw == true) {
for (m = 0; m <= N_END; m ++) {
for (i = 0; i < T_END / DT; i ++) {
x = noise(NOISE) + signal(i * DT、 X0、 W);
mem[i] = (mem[i] * m + x) / (m + 1);
}
if (m == mm - 1) {
for (i = 0; i < T_END / DT; i ++) {
g.setColor(Color.blue);
w.lineTo(j、 i * DT、 mem[i]、 g); /* jのウィンドウに表示 */
}
mm = mm * N;
j = j + 1;
}
}
stop();
}
}
/* 信号パルス 大きさpeak パルス幅p_width */
double signal(double t、 double peak、 double p_width) {
double t0 = 1.7; /* パルスの始まり時刻 */
if (t0 <= t && t <= t0 + p_width) {
return(peak);
} else {
return(0.0);
}
}
/* ガウス雑音 */
double noise(double sd) {
double xw = 0.0;
long i;
for (i = 1; i <= 12; i ++){
xw = xw + Math.random();
}
return(sd * (xw - 6.0));
}
}
|
実行結果
プログラムS0740.javaを実行してみましょう。信号の振幅 0.5 に対して、雑音の大きさが同じ標準偏差 0.5 として計算しています。原波形（ウインドウ０）において、信号はどこにあるのかまったくわかりません。それが、加算回数を増やして行くと、信号がうかびあがってくるのがわかります．加算回数６３回では、十分に信号を認識することができます。
| 図７．４．２ 同期加算による信号の分離 |
 |
| この章始め |
演習
１． ショットノイズは、歩幅１、歩く時間が不定のランダムウオークであることに、注意しなさい。（６章演習２．参照）
（解答省略）
２． S0730.javaにおいては、平均化点数を１１としていた。さまざまな平均化点数を用いてシミュレーションをおこない、出力波形を観測しなさい。
（解答省略）
３． 雑音にうもれた正絃波を、同期加算により、分離しなさい。
（解答） A0730.java | 実行 |
| この章始め |
| | HOME | 目次 | |
2009 (C) H.Ishikawa |
