のろのろのろ雑記

_ [ColdFire] ColdFire×P2P地震情報(5) ノイズを除去する

　Interface 9月号付録基板を使い，P2P地震情報のサービスを作ってしまう勝手な企画です．過去の記事は，タイトルのColdFireカテゴリからどうぞ．

(2)〜(4) のおさらい

• 地震計として使うためには，「ズレ」「ノイズ」「サンプリング間隔」の3点を考える必要があるようです．
  • タイマを使うことで，サンプリング回数を保つことがそれなりにできました．
  • "0 gの電圧差""傾いた設置"の問題は，静止した状態を"地震計にとっての0 g"と考えることでクリアしました．
  • 時間の経過による傾きの変化などは，移動平均を用いることでクリアしました．

1 gあたりの感度差は，今回はそのままで

　"(2) 地震計としてのColdFire基板"において，『加速度センサによる1 gあたりの感度差は約±75 Galの差として現れることがあります．』と記しました．この差を補正するには，それぞれの軸で0 gの状態と1 gの状態を調べ，正しい値となるよう常に一定倍率を掛ける方法が考えられます．

　しかし，一般的な環境で各軸に0 gと1 gの状態を作り出すのは難しく，正確な補正は期待できません．誤った補正が掛かる可能性もあります．さらに，この感度差による誤差はそれほど大きいものではありません．

　感度差が 740 mV/g 〜 800 mV/g (Typical) 〜 860 mV/g の間で変動するということは， 92.5 % 〜 107.5 % の間で変動する，つまり誤差は±7.5 %ということになります．

　気になると言えば気になりますが，そこまで問題となるレベルではありませんので今回はそのままにします．

"移動平均" でノイズ問題もクリア

　これまで様々な問題を取り上げてきましたが，最後に『ノイズがあります．』と記した部分を解決しましょう．前回オフセットの算出に移動平均を使いましたが，これをデータにも使います．「データをならす」といったイメージで，細かなノイズを取り除きます．

　これも一種のフィルタで，「非常にすばやい変化には追従できない」という副作用があります．以下に数パターンの特性を，気象庁が計測震度の計算に用いるフィルタの特性とともに示しました．凡例は "サンプリング周波数 / 移動平均サンプル数" を表しています．

　高周波数帯（1 Hz以上）の特性が変わり，ローパスフィルタのように作用しています．気象庁フィルタと同じとはいきませんが，1 Hz〜3 Hz辺りは"20 Hz / 4 Samples"が近い値になっています．

　実は，これだけではノイズは完全に取れていません．フーリエ変換すると分かりますが，静止状態のノイズはかなり幅広く分布しています．ただ，ノイズ除去を厳しくしすぎると地震波特有の波形を見失う恐れがあるので，これ以上フィルタを掛けるのは避けます．

　それでは恒例になりましたサンプルプログラム，前回のオフセットと今回のノイズ除去を組み合わせてみました．サンプリング間隔が未調整で25 Hzくらいで出てきますが，あくまでサンプルです．

main(){
  int i,index,ind_n;
  int *sample=MemoryAlloc(120); int *nr=MemoryAlloc(6); long *offset=MemoryAlloc(12);
  // init
  InitAd(0x70); index=0; ind_n=17;
  for(i=0;i<3;i++){offset[i]=0;nr[i]=0;}
  for(i=0;i<60;i++)sample[i]=0;
  // loop
  #stop 0
  for(;;) {
    if(Getc(0)!=0)break;
    index=(index+1)%20;ind_n=(ind_n+1)%20;
    for(i=0;i<3;i++){
      offset[i]-=sample[i*20+index];nr[i]-=sample[i*20+ind_n];
      sample[i*20+index]=GetAd(i+4);
      offset[i]+=sample[i*20+index];nr[i]+=sample[i*20+index];
      PrNum(nr[i]/3-offset[i]/20);PrStr(", ");
    }
    PrStr("\r\n");SystemSleep();
  }
  MemoryFree(sample);MemoryFree(nr);MemoryFree(offset);
}

2種類の"移動平均"で，気象庁フィルタに近づく？

　サンプリング周波数20 Hzを基本とし，オフセット（基準点）に最新20サンプルの移動平均を用いるケース，データに最新3サンプルの移動平均を用いるケース，それらを組み合わせたケース，そして気象庁が計測震度の計算に用いるフィルタ，それぞれの特性を示しました．

　当初目標としていた『5HzPGA』，0.1 Hz〜5 Hzまでのバンドパスフィルタを掛けることと比べると，0.1 Hzや5 Hzで信号強度が約0.3 倍に落ち込んでいることから，『5HzPGA』として用いるには不適切な特性であることが分かります．まっとうにバンドパスしようとすると複雑な処理が必要になるので，無理があります．

　一方，気象庁フィルタと比べると多くの部分が近似しています．気象庁の計測震度「っぽいもの」を目指す可能性が急浮上してきました．がんばってみましょう．

今回のまとめ

• ノイズについては，移動平均を使うことで緩和しました．
• 何故か気象庁フィルタに近い特性を示したので，その方向でがんばることにしました

次回以降の予定

　これで全問題をクリアし，"地震計としてのColdFire基板"が一応出来上がりました．次回以降は，以下のような点を考えて"ColdFire基板 地震計"を制作していきます．

• 震度の計算が必要ですが，気象庁の方法を実装するのは無理があります．また，『5HzPGA』と別物になるため，5HzPGAのデータに基づく算出も出来ません．
  • 計測震度と，その元になった波形データをフィルタに通した値との相関を見出せるかもしれません．
  • ここから最大加速度だけでがんばれるかもしれません．
• 通信仕様を全く考えていません．

参考

• http://www.seisvol.kishou.go.jp/eq/kyoshin/kaisetsu/calc_sindo.htm
• http://www.jma.go.jp/jma/kishou/know/shindo/keisoku.html

_ [ColdFire] ColdFire×P2P地震情報(6) 計測震度との対応を見る

　Interface 9月号付録基板を使い，P2P地震情報のサービスを作ってしまう勝手な企画です．過去の記事は，タイトルのColdFireカテゴリからどうぞ．

前回までのまとめ

• 「ズレ」「ノイズ」「サンプリング間隔」の問題をクリアし，"地震計としてのColdFire基板"が出来上がりました．
• これから，いよいよ"ColdFire基板 地震計"を制作していきます。

計測震度との対応を見る

　揺れの大きさとして「震度」が一般的であることや，震度階級の元となる計測震度が揺れの周期などを十分考慮して決定されていることを考えれば，ColdFire基板でもおおよそ相当する震度を出したいものです．

　そこで，これまでズレやノイズを解消するために作った『フィルタ』を通した値が，気象庁の計測震度とどう関係しているかを見ていきます．

　幸い，気象庁には100 Hzサンプリングの加速度データと計測震度が一部公開されています．20 Hzサンプリングにするため5個に1個の間隔でデータを取得してフィルタを通し，ベクトルを用いて加速度の大きさを求めます．

　フィルタを通さない場合の3方向，通した場合の3方向と2方向，3ケースを比べました．フィルタを通すことにより，計測震度との対応が良くなっていることが分かります．±0.5くらい誤差がありますが，ひとますこれでやってみましょう．

まずはシミュレーション

　まずは，表計算ソフトを用いてシミュレーションしてみます．ExcelやOpenOffice.org Calcを例にすると，フィルタはAVERAGE，ベクトルは各軸の2乗を足したものをSQRTすればよいので，プログラムと全く同等のシミュレーションができます．

　震度の算出については，Excelが "y = 0.7968 Ln(x) + 0.9717" と言っているのでこれを用います．最初に，静止した状態のノイズレベルがどの程度の震度になるのかシミュレーションしてみました．

　なんということでしょう！　静止しているのに震度3だなんて…

　つづく…

参考

• http://www.seisvol.kishou.go.jp/eq/kyoshin/jishin/index.html

_ [ColdFire] ColdFire×P2P地震情報(7) "33.3 Hzの壁"を超える

　Interface 9月号付録基板を使い，P2P地震情報のサービスを作ってしまう勝手な企画です．過去の記事は，タイトルのColdFireカテゴリからどうぞ．

前回のまとめ

• 計測震度とフィルタとの対応を見ました．
• 静止状態は震度3相当でした…

"33.3 Hzの壁"

　静止状態で震度3になってしまったということは，すなわち震度4以上でないと働かないということであって，それを震度計・地震計とは呼べないということです．

　何とかしてノイジーな状態を解消しなければなりません．ひとまず，"オーバーサンプリング"を考えてみます．オーバーサンプリングの効果は私も詳しく知りませんが，とりあえず震度3という表示にビックリしたので何でもやってみます．

loop_axis(){
 char t1; int i,cr;
 t1=CreateTimer(0,30000,0);
 InitAd(0x70);
 for(i=0;i<1000;i++){
  GetAd(4);GetAd(5);GetAd(6);
 }
 cr=GetTimerCount(t1);
 PrNum(cr);PrStr("\r\n");
 CreateTimer(t1,0,0);
}

Ax_Timer::loop_axis
26999

Ax_Timer::loop_axis
27000

　上記のコードは，3軸を1000回サンプリングした時の所要時間を計測するものです．タイマは10 ms単位ですから， "(30000 - RESULT) / 100" と計算すれば所要秒数が分かります．

　計算してみると，ちょうど30 s掛かっていることが分かります．サンプリングのためにGetAd関数を3軸×1000回＝3000回呼び出していますから，「もしかして，GetAd関数はキッチリ10 ms掛かるようになっているんじゃないか？」，つまり3軸取得時には『33.3 Hzの壁』があるのではないか，という疑問が出てきます．

"33.3 Hzの壁"を超える

　この壁を超えます．新適当マイコン電子工作研究所あたりを見ていると，どうやらメモリマップドI/Oのようです．ということは，A/Dコンバータ(ADC)とてメモリアクセスで簡単に出来るはず…

　リファレンスマニュアルを調べると，"Chapter 28 Analog-to-Digital Converter (ADC)"，493ページ以降に書かれています．

• 8チャネルの入力と，2つのA/Dコンバータ
  • AN0〜3はコンバータAに，AN4〜7はコンバータBに接続
• サンプリングするチャネルを指定し，シーケンシャルまたはパラレルサンプリング
  • 1回のサンプリング命令で最大8サンプル（チャネル）をサンプリングしデータを格納
• 2チャネルの差分値，オフセットなどの設定が可能

　GetAd関数からは想像もつかない作りになっています（？）．この辺りを自由自在に操れると，SilentCでもまだまだ頑張れる気がしてきます．

　実際にこれを操作する方法も，リファレンスマニュアルに載っています．ここでは簡単に示します．

1. ADLST1とADLST2の設定(28.4.4/500ページ)．リザルトレジスタが8サンプル分用意されているので，どのリザルトレジスタにどのチャネルからのサンプリングデータを注ぎ込むか指定する．
   • パラレルサンプリング実行時にはチャネル指定に制限があるので注意する．
2. 必要であれば，ADOFSnでオフセットを設定(28.4.11/509ページ)．
3. 必要であれば，CTRL1でサンプリング方法を指定．
4. CTRL1を設定してサンプリングし，ADRSLTnを読み取ることの繰り返し．
   • ADRSLTnはシフトされていることに注意する．

　これを使って全チャネル(AN0〜7)を読み取ってみましょう．

adc(){
  int i;
  char *ctr; ctr=0x40190000; //ADC CTRL1
  int *dat; dat=0x40190012; //ADC ADRSLT
  InitAd(0x70);
  #stop 0
  for(;;) {
    if(Getc(0)!=0)break;
    ctr[0]=32; //サンプリング
    for(i=0;i<8;i++){
      PrNum(dat[i]>>3);PrStr(", "); //AN0-AN7
    }
  PrStr("\r\n");
  SystemSleep();
}

2548, 2352, 2256, 2211, 1921, 2238, 2771, 2275,
2550, 2335, 2158, 2138, 1923, 2238, 2789, 2195,
2544, 2306, 1972, 1990, 1924, 2250, 2802, 2272,
2548, 2304, 2029, 2181, 1930, 2249, 2801, 2470,
2551, 2305, 2288, 2388, 1910, 2210, 2790, 2306,

OK

　AN4,5,6がX,Y,Z軸です．ちゃんと読み取れました．サンプリングにはある程度時間が掛かりますが，SilentCがコードを読む間に終わっているようなのでウェイトは入れません．さて，これでどれだけ早くなるでしょうか．

loop_axis(){
  char t1, *k; int *a,i,cr;
  k=0x40190000; a=0x40190012;
  t1=CreateTimer(0,30000,0);
  InitAd(0x70);
  for(i=0;i<1000;i++){
    k[0]=32;
    cr=a[4] >> 3;cr=a[5] >> 3;cr=a[6] >> 3;
  }
  cr=GetTimerCount(t1);
  PrNum(cr);PrStr("\r\n");
  CreateTimer(t1,0,0);
}

Ax_Timer::loop_axis
29710

Ax_Timer::loop_axis
29705

　ビットシフトして値を出すところまで行いました．およそ3 s，GetAd関数を呼ぶ時よりも10倍早くなっています！　移動平均などの計算時間があるのでストレートに10倍のサンプリング周波数は難しいですが，GetAd関数での待ち時間がなくなるのは大きいです．

参考

• http://noritan-micon.blog.so-net.ne.jp/
• http://www.freescale.com/files/32bit/doc/ref_manual/MCF52235RM.pdf
• http://www.cqpub.co.jp/interface/download/2008/09/kairozu.pdf

_ [ColdFire] ColdFire×P2P地震情報(A) PWMと圧電ブザーで音を鳴らす

　Interface 9月号付録基板を使い，P2P地震情報のサービスを作ってしまう勝手な企画です．が，番号がアルファベットのものはおまけです．

二番煎じに関する注意事項

　これは新適当マイコン電子工作研究所で行われていることの二番煎じです．が，PWMの出力をクリスタルイヤホンで…　などと見ていると，「これだ」と思ってそのまま突っ走ってしまったわけです．お付き合いください．

経緯

1. 「PWMの出力を…　クリスタルイヤホン？」というわけで，検索する．
2. 「ふむふむ，圧電効果とかいうのか…」
3. 「圧電？　圧電ブザーってあるね？　つーか同じじゃね？」
4. 「じゃあ圧電ブザー繋げて同じことすれば音鳴るんじゃね？」

材料

• ColdFire基板 (Interface 9月号付録，RJ-45とDCジャックを実装済み)
• 何に使ったか覚えていない導線
• 古いマシンから切り取った圧電ブザー
  • "HYCOM HY-05"と書かれていて，どうやら3.0 V〜8.0 Vで動作するらしい．
• 私の時間

手順1: 本当に鳴るのか

　まずは，MCF52233付録基板 - SilentCが使っていないモジュールを探せ(4)に記載されているとおりにやってみます．

　(PDF)リファレンスマニュアルや(PDF)回路図を見ると，PNQ[1]＝IRQ1のようなので，次のように接続します．

　すると，確かに「プツッ，プツッ，プツッ…」と3.6 Hzくらいで鳴っているのが聞こえます．いけるようです．

手順2: プリスケーラ？　スケーラ？　デューティ？　周波数を設定するには…

　圧電ブザーは高速でON/OFFを繰り返すパルス信号によって，その信号の周波数の音を出すようです．ですから，この周波数さえ設定できれば音の高低を自由自在に操れるわけです．

　が，WindowsのBeep関数のように「周波数」という設定項目はありません．あるのは，プリスケーラ，スケーラ，デューティ，意味の分からない項目ばかりです．じっくり向き合う必要があるようです．

手順3: ColdFireの"Pulse Width Modulation(PWM)"って？

　と言いつつ，PWMの解説は他所に譲ります．ここでは，ColdFireで用意されている設定項目について見ていきます．リファレンスマニュアルをご用意ください．

PWM Clock Select Register (PWMCLK) - 29.2.3/534ページ

PWMのクロック源を選ぶ．AやBを指定するとプリスケーラ(PWMPRCLK)のみ，SAやSBを指定するとプリスケーラとスケーラ(PWMSCLA/PWMSLCB)の両方を通る．

PWM Prescale Clock Select Register (PWMPRCLK) - 29.2.4/535ページ

プリスケーラの分周比を選ぶ．簡単に言うと，内部バスクロックの何分の1のクロックで動くようにするかを指定する．ただし，分母は2の0乗,1乗,2乗,…,7乗．

PWM Scale A/B Register (PWMSCLA/B) - 29.2.7〜/538ページ〜

プリスケーラクロックの何分の1で動くようにするかを指定する．ただし，分母は2,4,6,…,512．

PWM Channel Period Registers (PWMPERn) - 29.2.10/540ページ

PWMCLKなどで指定したクロックを元にカウントを行うが，そのカウンタの上限値，1サイクルに掛かるクロック数．

PWM Channel Duty Registers (PWMDTYn) - 29.2.11/541ページ

1サイクルに掛かるクロック数に対し，信号がONになるクロック数．PWMPERnと同じにすれば常にONだし，半分ならON/OFFがピッタリ交互だし，0ならOFFのまま．

　プリスケーラ，スケーラ，ピリオド（デューティ）を組み合わせることで，やっと周波数を指定出来るというわけです．これはなかなかしんどい．

　スケーラを通る（PWMCLKにSA/SBを指定した）場合の周波数の計算式を置いておきます．ただし，デューティ比（信号ONの割合）は50 %とします．

=   60*1000*1000   / ( 2 ^ PWMPRCLK ) / ( 2 * PWMSCLA/B ) / PWMPERn
  内部バスクロック     プリスケーラ          スケーラ       ピリオド

手順4: 音階に対応する組み合わせ

　早速組み合わせてみました．平均律を用いています．0.6 Hzくらいズレているところもありますが，お遊びということで．

　ポイントは，「プリスケーラを1減らすと，周波数が2倍＝1オクターブ上になる」こと，そして逆もまたあるということです．上記のような1オクターブ分の音階表を持っておけば，あとはプリスケーラの値だけで対応できます．

　プリスケーラは0〜7まで指定できるので，この場合だとG#2〜G 8までは自由に出せることになります．この範囲を外れてしまうと所々音が鳴らなくなります．

手順5: 実装

mel(){
  // char *meloは可読性のため改行を入れたので，実行時には削除すること
  char *melo="5E-85D#85E-85D#85E-84B-85D-85C-84A-83E-83A-84C-84E-84A-84B-83E-8
              3G#84E-84G#84B-85C-83E-83A-84E-85E-85D#85E-85D#85E-84B-85D-85C-8
              4A-83E-83A-84C-84E-84A-84B-83E-83G#84D-85C-84B-84A-83E-83A-8SS-8
              EE";
  char *pt=melo; int bt,slp;
  bt=10000/(16000/60); // テンポを指定し，4分音符の長さを10 ms単位で計算
                       // "10000 / (テンポ*10 / 60)"
                       // 小数計算が出来ないため簡易精度確保

  int *pnq = 0x40100068;  // PNQ[1]
  char *pwm = 0x401b0000; // PWM

  *pnq |= 12; pwm[2]=2;   // PNQ[1]とPWM1を接続，PWM1を有効に

  #stop 0
  for(;;){
    if(Getc(0)!=0)break;
    if(pt[0]==69)break;

    if(pt[1]==67&&pt[2]==45){pwm[3]=4;pwm[8]=256;pwm[21]=14;pwm[29]=7;} //C
    if(pt[1]==67&&pt[2]==35){pwm[3]=3;pwm[8]=199;pwm[21]=34;pwm[29]=17;} //C#
    if(pt[1]==68&&pt[2]==45){pwm[3]=3;pwm[8]=103;pwm[21]=62;pwm[29]=31;} //D
    if(pt[1]==68&&pt[2]==35){pwm[3]=3;pwm[8]=137;pwm[21]=44;pwm[29]=22;} //D#
    if(pt[1]==69&&pt[2]==45){pwm[3]=3;pwm[8]=39;pwm[21]=146;pwm[29]=73;} //E
    if(pt[1]==70&&pt[2]==45){pwm[3]=3;pwm[8]=179;pwm[21]=30;pwm[29]=15;} //F
    if(pt[1]==70&&pt[2]==35){pwm[3]=3;pwm[8]=181;pwm[21]=28;pwm[29]=14;} //F#
    if(pt[1]==71&&pt[2]==45){pwm[3]=5;pwm[8]=46;pwm[21]=26;pwm[29]=13;} //G
    if(pt[1]==71&&pt[2]==35){pwm[3]=2;pwm[8]=215;pwm[21]=42;pwm[29]=21;} //G#
    if(pt[1]==65&&pt[2]==45){pwm[3]=2;pwm[8]=147;pwm[21]=58;pwm[29]=29;} //A
    if(pt[1]==65&&pt[2]==35){pwm[3]=2;pwm[8]=149;pwm[21]=54;pwm[29]=27;} //A#
    if(pt[1]==66&&pt[2]==45){pwm[3]=2;pwm[8]=211;pwm[21]=36;pwm[29]=18;} //B
    if(pt[1]==83){*pwm=0;}else{*pwm=2;}
    pwm[3]-=pt[0]-51; //所詮ビープなので，指定より2オクターブ高く出して聞こえやすく

    if(pt[3]==52){slp=bt;}
    if(pt[3]==56){slp=bt/2;}
    if(pt[3]==57){slp=bt/3;}
    if(pt[3]==50){slp=bt*2;}
    if(pt[3]==49){slp=bt*4;}

    Sleep(slp);

    pt+=4;
  }
  *pwm=0;
  PrStr("E.\r\n");
}

　char *meloでメロディを作ります．4バイトで1個の音を奏でる独自形式です．

[オクターブ][音階][シャープ][長さ]
 ※音階はA〜GまたはS（休み）で，シャープは#または-（ありまたはなし）
 ※長さは 1（全音符），2（2分音符），4，8，9（3連8分音符）のみ

　forループ内で次のような処理を繰り返します．

1. オクターブが E であれば，ループを抜ける．
2. 音階とシャープを見て，プリスケーラなどを設定
3. 休みであればPWMオフ，そうでなければPWMオン
4. オクターブを見て，プリスケーラを加減
5. 長さを見て，演奏時間を計算してスリープ
6. 次の位置へポインタをセット

　若干順番が乱れていますが，気にしないことにします．かなり強引にやってしまいましたが，これを実行すると「エリーゼのために」が流れます．

　制限として，約256バイトに収める必要があるようです．それを超えたり，終了命令（オクターブに E）を忘れたりすると，範囲外にまで演奏を続行して変な音になります．

手順6: 実演

• ニコニコ動画: ColdFire基板と圧電ブザーで音を鳴らしてみた(YouTube版)

　SilentCからここまで操作できちゃうとは．プログラミングガイドを読んで「こんなものか」と思っていたのがバカみたいです．スバラシイ．

参考

• http://noritan-micon.blog.so-net.ne.jp/
• http://www.freescale.com/files/32bit/doc/ref_manual/MCF52235RM.pdf
• http://www.cqpub.co.jp/interface/download/2008/09/kairozu.pdf
• http://www.apache-band.com/apchhp/music_theory/frequency_temperament.html

_ [ColdFire] ColdFire×P2P地震情報(8) 実際どこまで出せるのか

　Interface 9月号付録基板を使い，P2P地震情報のサービスを作ってしまう勝手な企画です．過去の記事は，タイトルのColdFireカテゴリからどうぞ．

前回のまとめ

• A/Dコンバータへ直接アクセスして，"33.3 Hzの壁"を超えました．

実際どこまで出せるのか

　A/Dコンバータへ直接アクセスすることで，どれくらいまでオーバーサンプリング出来るでしょうか．

直接送信

main()
{
char soc;
char *ctr;ctr=0x40190000;
int *dat; dat=0x40190012;
int *sd=MemoryAlloc(8);
soc=CreateSocket(0);if(soc<0){PrStr("Soc!");return;}
Bind(soc,6909,0);sd[3]=3338;
PrStr("START!\r\n");
for(;;) {
  if(Getc(0)!=0)break;
  *ctr=32;
  sd[0]=dat[4]>>3;sd[1]=dat[5]>>3;sd[2]=dat[6]>>3;
  SendTo(soc,0xc0a80b69,6909,sd,8);
}
PrStr("GAME OVER\r\n");
CloseSocket(soc);
}

　取得値をそのままUDP送信した場合です．表示を見ると，258.56 Hz出ていることが分かります．

フィルタ送信

main(){
  char soc; int i,j,index,ind_n;
  int *sample=MemoryAlloc(600); long *nr=MemoryAlloc(12); long *offset=MemoryAlloc(12);
  int *data=MemoryAlloc(8); char *k; int *a;
  k=0x40190000; a=0x40190012;
  // init
  soc=CreateSocket(0);if(soc<0){PrStr("Soc!");return;}
  Bind(soc,6909,0); data[3]=3338;
  InitAd(0x70); index=0; ind_n=85;
  for(i=0;i<3;i++){offset[i]=0;nr[i]=0;}
  for(i=0;i<300;i++)sample[i]=0;
  // loop
  #stop 0
  for(;;) {
    if(Getc(0)!=0)break;
    *k=32;
    index=(index+1)%100;ind_n=(ind_n+1)%100;
    for(i=0;i<3;i++){
      j=i*100+index;
      offset[i]-=sample[j];nr[i]-=sample[i*100+ind_n];
      sample[j]=a[i+4] >> 3;
      offset[i]+=sample[j];nr[i]+=sample[j];
      data[i]=nr[i]/15-offset[i]/100;
    }
    SendTo(soc,0xc0a80b69,6909,data,8);
  }
  MemoryFree(sample);MemoryFree(nr);MemoryFree(offset);MemoryFree(data);
  CloseSocket(soc);
}

　移動平均によるフィルタを掛けてUDP送信した場合です．Unsignedで受け取ったので値が大変なことになっていますが，それよりも問題はサンプリング周波数．41.15 Hzと，先ほどの6分の1くらいしか出ていません．

　やはり，インタプリタで行うには限度があるようです．

どうしようか

　"あくまでSilentCで"ということで考えると，次のような方法があります．

• 集計サーバに取得値をそのまま送りつけて，震度を計算してもらう
  • データ量は 6bytes * 100 Hz = 600bytes/s = 4800bps
  • インターネットごしにUDP送信するのか？　TCP送信だとこんな速度出るのか？
• LAN内のパソコンに送りつけて，震度を計算．ないときは低いサンプリング周波数で単独動作
  • 切り替えをうまくできるだろうか？
  • 単独動作時の低いサンプリング周波数によるノイジーな状態を許容できるだろうか？
• 30 Hzや40 Hzくらいで，今までよりはマシにする
  • そんなに変わるだろうか？

　しばらく悩むことになりそうです．

_ [ColdFire] ColdFire×P2P地震情報(9) 熱にご注意，SilentCベンチ

　Interface 9月号付録基板を使い，P2P地震情報のサービスを作ってしまう勝手な企画です．過去の記事は，タイトルのColdFireカテゴリからどうぞ．

　今回は，加速度センサやSilentCの気になるところを調べてみました．

熱にご注意

　P2P地震情報の中の人から「加速度センサの値って熱の影響受けやすいみたい．実験したからデータあげる」と言われたので，データを見てみました．ドライヤーで加熱したそうな．

　って，変わりすぎ！　データシート，データシート…　加速度センサのデータシートでは，0 gの感度差は±2.0 mg/℃，感度差は±0.03 %/℃が仕様とされているようです．いや，緑色のX軸は明らかにその範囲外．温度センサを付けていなかったので，温めすぎた可能性も否定できませんが．

　しかしながら，このような緩やかな変化は移動平均を使ってオフセットを出せば概ね解決します．問題になるのは，熱でノイズが増すかどうか．X軸のFFTを出してみました．

　ドライヤーを入れる前が平常時，入れた後が加熱時，入れてしばらく経ったのが過熱時です．平常時と比べると，加熱時，過熱時はノイズレベルが1段，2段ほど増えていることがわかります．

　ドライヤーの風による影響も含まれていると思いますが，同じ風量を当てた「加熱時」と「過熱時」とでノイズレベルに差があることから，ドライヤーの風だけではなく，熱（温度）が関係していると言ってもよいでしょう．

　つまり，ノイジーな状態を緩和するには，基板の熱対策も重要になってくる…　かもしれない，ということです．データシートきっちり読める人，人柱実験が大好きな人を募集しています（？）．

SilentCベンチ

　SilentCはインタプリタです．"filename::funcname"という風に別ファイルの関数を呼べるのは大変便利ですが，パフォーマンスが心配になります．せっかくなので，調べてみましょう．

main(){
 char c,soc; int i;
 //事前処理が必要な場合はここ
 c=CreateTimer(0,30000,0);
 for(i=0;i<10000;i++){
   //計測したい処理
 }
 PrNum(GetTimerCount(c));
 CreateTimer(c,0,0);
 //事後処理が必要な場合はここ
}

　"(30000-出力値)*10/10000"と計算すると，処理の所要時間がms単位で出ます．ただし，forループの処理時間を考慮しなくてはなりません．予めループ内に何もない状態で計測を行い，それをループ処理そのものに掛かった時間と仮定して一律差し引くようにしました（0.4543 ms）．

　なかなか面白いところがありますね．

• 関数を呼び出す場合，同じファイルの末尾に書くよりも別ファイルの先頭に書いたほうが早い．
  • 上から 別ファイル1行目の関数呼び出し，同じファイル末尾の関数呼び出し，同じファイル1行目の関数呼び出し です．
• メモリマップドI/Oなどでポインタを使って読み書きする場合，"pointer[0]=0;"を"*pointer=0;"とすると約2倍早い（箇所は限定される）．
• 途中で終了するための入力チェックは結構重たい．UDPの受信チェックなんかはさらに重たい．
• 値を10進数で書いても16進数で書いてもあまり変わらない．バイト数の関係か，10進数で書いたほうが気持ち早い．
• 除算はビットシフトのほうがやや早いが，数値で書いても遅くはない．

参考

• http://www.freescale.co.jp/pdf/MMA7360Lr1.pdf
• http://www.freescale.com/files/32bit/doc/data_sheet/MCF52235DS.pdf

_ [ColdFire] ColdFire×P2P地震情報(10) よし，サーバに送りつけよう

　Interface 9月号付録基板を使い，P2P地震情報のサービスを作ってしまう勝手な企画です．過去の記事は，タイトルのColdFireカテゴリからどうぞ．

オーバーサンプリングで効果的にノイズ除去できるのか？

　そういえば，オーバーサンプリングをすることでノイズ除去できるかどうか調べていませんでした．ということで調べました．

20 Hzサンプリング

20 Hzサンプリング，20サンプル移動平均をオフセットに，3サンプル移動平均をデータに用いた．

100 Hzサンプリング

100 Hzサンプリング，100サンプル移動平均をオフセットに，15サンプル移動平均をデータに用いた．

　もとになるデータは気象庁の加速度データ．フィルタを通した時の最大加速度と計測震度との関係を近似式として算出して用います．

20 Hzサンプリング時

0.7968 * Ln(最大加速度) + 0.9717

100 Hzサンプリング時

0.8055 * Ln(最大加速度) + 0.9568

　劇的な効果とは言えませんが，震度換算で0.6ほどノイズが減りました．うーん，もっと減ってくれてもいいんですが…

　ノイズで悩んでばかりでは進まないので，これ以上は後回しにしてさっさと作りましょう．

よし，サーバに送りつけよう

　ColdFire×P2P地震情報(8)で実装方法をあれこれ悩んでいましたが，結局『集計サーバに取得値をそのまま送りつけて，震度を計算してもらう』ことにしました．

• より小さな地震を検出できる
• ネットワーク化して情報を活用できる
• P2P地震情報の各情報もやり取りできる

　こうしたメリットは大きいです．「ネットワークに依存する」というデメリットがありますが，切断時には単独動作するようなモードを設ければ良いでしょう．

通信仕様を決めよう

　実装する前に，通信仕様を決めましょう．

• 加速度センサの値を送信して，震度換算値を取得する
  • データ記録，集計，公表はサーバ側で行う
• 地震情報や津波予報の概要，緊急地震速報の発表情報を取得する
  • 都道府県などによって送信条件を指定できるとよい
• ColdFire基板（SilentC）専用にならない程度に，実装が楽になるような何か

　ということで，暫定仕様をでっち上げました．SilentCで実装する以上，仕様を隠す意味はどこにもありませんので，この辺に置いておきましょう．

EPSP ColdFire Gateway（EPSP-CFG） プロトコル仕様 Ver.0x01

• クライアント-サーバ型，サーバはP2P地震情報（EPSP）ネットワークとのゲートウェイも兼ねる
• データ量の小さいバイナリプロトコル
• 加速度センサの活用と諸情報の受信機能

　次回以降，この仕様をもとに実装していきます．

参考

• http://www.seisvol.kishou.go.jp/eq/kyoshin/jishin/index.html

_ [地震][PC] ウェザーニューズ(WNI)の「The Last 10-Second」が大変なことになっている件（→なっていた件）

　※2008年9月29日追記：9月29日17時30分ごろ、WNIのサーバが約2倍に増強されました。全サーバに余裕が生まれ、スムーズに接続できるようになっています：ウェザーニューズ(WNI)の「The Last 10-Second」が大変なことになって"いた"件

　※2008年9月21日時点の不確かなデータに基づくものです。一切の正確性・最新性その他を保証しません。

　ウェザーニューズの緊急地震速報サービス「The Last 10-Second（L10SまたはTLTSとも）」が大変なことになっているらしい件。

まとめ

• 頻繁に（20分おきに）切れて、なかなか繋がらない。
  • 頻繁に切れる原因はよく分からない。
  • なかなか繋がらないのは「WNIのサーバがいっぱいいっぱいだから」。

• 20分おきに2分切れる場合、「緊急地震速報を受信出来ない時間」が9.1％発生！
  • 防災目的でマジメに使う場合は、改善を求めるか他サービスの利用を考えるべし。

現象

　私はサービスを使っていないのですが、サービスを使う人曰く、

• 20分おきに切れる
• ログを見ると『Try connect』が複数回続いて『Communicating』になる

　とのこと。つまり、頻繁に切れた上になかなか繋がらない（接続が回復しない）ようです。

サービスの仕組みと原因

　詳しい人に調べてもらいました。私自身は未検証です。

1. HTTPでサーバリストを取得、リスト内のサーバへHTTP接続・認証
2. 認証後は、定期的に接続確認のための通信と、緊急地震速報の通信が行われる

　という仕組みだそうです。ダウンロードページのシステム要件にも"HTTPアクセスが可能なこと"とありますが、おそらくプロキシのある環境での利用を考慮しているのでしょう。

　で、定期的に接続確認しているんだったら切れないだろうと思うのですが、「突然切断処理が行われる」らしい。

　さらに、その後繋がりにくい現象については、次のように言われました。

1. サーバへHTTP接続・認証するとき、「セッションが限界だ」と言われて、
2. 取得したサーバリスト内の別のサーバを選択し、やり直す

• この「やり直す」のが何度も起きている

　つまりは、「たくさんのサーバの中で、接続できるサーバを何度も何度も探しているから時間が掛かる」のを「繋がりにくい」という風に感じていることになります。

で、どれくらい繋がりにくいの？

　…と言ったら、しばらくしてアドレスが返ってきました。

• WNIのLastなんちゃら監視所

　「(・∀・)イイ」が接続できるサーバの数、「（´・ω・`）ショボーン」が接続できない（セッション限界）サーバの数で、グラフの青と赤にそれぞれ対応しているとのこと。ここの1週間分のデータをグラフにしたものが最初の画像です。

　夜間帯を見ると「(・∀・)イイ」が0になっていて、ほとんど接続できない状態になっているようです…　接続リトライが行われるので厳密には「時間が経つと繋がる」ようですが、これはほぼ限界に近いのではないでしょうか…

対策

諦める

WNIに改善を求める

　免責事項や利用規約に「こういうことが起きても知らないよ」と書いてあるかもしれませんが。

複数のマシンで使わない

　※そういう仕様になっていたら読み飛ばしてください。

　1アカウントにつき1台でのみ使用し、サーバにやさしい会員になっておく。そして、WNIがサーバ増強を予定していることを期待するというもの。

他のサービスを使う

　「OCN緊急地震速報」や「緊急地震速報 フレッツタイプ」を使う。どちらもIPv6マルチキャストを使用しているので、少々コストが掛かるものの配信遅延はきわめて小さい。サーバの障害耐性は不明。

P2P地震情報の緊急地震速報 配信試験（オープンβ）を使う

　P2P地震情報の「緊急地震速報 配信試験（オープンβ）」を使う。無償で使える一方、「一般向け緊急地震速報の『発表』」しか分からないし、6秒くらい遅延するしで他サービスより大きく劣る。

　というわけで、エロイ人ありがとうございました。

_ [地震][PC] ウェザーニューズ(WNI)の「The Last 10-Second」が大変なことになって"いた"件

　先日 ウェザーニューズ(WNI)の「The Last 10-Second」が大変なことになっている件 と題してサービスのサーバがいっぱいいっぱいであることをお伝えしましたが、どうやら改善されたようです。

• WNIのLastなんちゃら監視所

　9月29日17時30分を過ぎたあたりで、サーバ数が18から39と約2倍になっています。そして、利用者の接続が各サーバに分散したのか、全サーバで "(・∀・)イイ" 表示になりました。

　これだけ増えれば、利用者の多い時間帯でも大丈夫でしょう。ああ、でも20分おきに切れる現象は直ってるんでしょうか？