勉強しないとな~blog

ちゃんと勉強せねば…な電気設計エンジニアです。

マイコンと加速度センサでタコメータを作る - 4. 実データ解析2

RXタコメータをバイクに乗せて比較的長距離走ってみたので、これを解析します。

変更点

いくつか変更したことがあります。

加速度ログデータフォーマットの変更

加速度ログデータの先頭に設定値を残すよう、プログラムの変更を行っています。これで、データから自動で測定間隔、測定分解能を取得してグラフの軸に反映させることができます。

データファイルフォーマット :

アドレス(byte) サイズ(byte) 内容
0x0 0x4 ヘッダ識別子 'LOGH' (0x4C4F4748)
0x4 0x4 ヘッダバージョン 0x1
0x8 0x4 ヘッダデータサイズ(ヘッダ識別子、バージョンは含まない) 0x8
0xC 0x4 測定レート ADXL345レジスタに設定する4bitコード
0x10 0x4 測定レンジ ADXL345レジスタに設定する2bitコード
0x14以降 加速度データ
  • ヘッダ識別子を入れることで、この情報領域を持っていないバージョンのデータでも同じ解析ツールで解析することができます。
  • ヘッダバージョンを入れることで、またあとでフォーマットを変えたくなった場合でも後方互換の解析ツールを用意することができます。
  • ヘッダデータサイズを入れることで、ヘッダの読み飛ばしが簡単にできるようになります。

グラフ表示

横軸を秒単位で、縦軸をG単位(重力加速度単位)にするようにしました。先ほどのヘッダから情報を取得します。

これを反映して、ログファイル読み込みのOctaveスクリプトを変更しました。

function [gData ovrn] = read_adxl345_log_bin_2byte(fname)

  fid = fopen(fname, "r");

  # Read header
  hdrId = fread(fid, 1, "uint32");
  
  if (hdrId == 0x4C4F4748) # 'LOGH', little endian
    hdrVer = fread(fid, 1, "uint32");
    hdrSize = fread(fid, 1, "uint32");
    
    printf("Header detected, Ver:%d, Size:%d\n", hdrVer, hdrSize);

    if(hdrVer == 1)
      hdrRate = fread(fid, 1, "uint32");
      hdrRange = fread(fid, 1, "uint32");
      
      rateCodes = [0x0, 0x1, 0x2, 0x3, 0x4, 0x5, 0x6, 0x7, ...
                   0x8, 0x9, 0xA, 0xB, 0xC, 0xD, 0xE, 0xF];
      rateVals = [0.1, 0.2, 0.39, 0.78, 1.56, 3.13, 6.5, 12.5, ...
                  25, 50, 100, 200, 400, 800, 1600, 3200]; # Measure frequency in Hz
      
      rangeCodes = [0x0, 0x1, 0x2, 0x3];
      rangeVals = [2.0, 4.0, 8.0, 16.0]; # Measure range in G

      rate = rateVals(find(rateCodes == hdrRate)(1));
      range = rangeVals(find(rangeCodes == hdrRange)(1));
      
      dt = 1 / rate;
      gUnit = range * 2 / 1024.0; # 10bit ADC, -range ~ +range
      
      printf("Rate:%d[Hz], Range:%d[G]\n", rate, range);
      
    endif
    
    dataTop = 12 + hdrSize;
    
  else # Header was not detected
    dt = 1;
    gUnit = 0.004; # 4mg/LSB
    fseek(fid, 0, SEEK_SET);
    dataTop = 0;
    
    printf("Header was not detected\n");
  endif

  # Data body
  gData = fread(fid, Inf, "int16", 2);
  gData = [dt*[1:size(gData)(1)].', gUnit*gData];
    # gData(:,1) : Time
    # gData(:,2) : Acceleration
  fseek(fid, dataTop + 2, SEEK_SET);
  ovrn = fread(fid, Inf, "int16", 2);
  fclose(fid);
endfunction

今までのスクリプトでは、読み出し結果として返すgDataは単に加速度データ列でしたが、このスクリプトでは秒単位の時刻列も追加しています。
また、加速度の単位もGに直しています。

# gData(:,1) : Time
# gData(:,2) : Acceleration

実データでグラフ表示

Octaveグラフ表示の手順は以下の通りです。

まずバイナリファイル読み込み。
読み取ったヘッダ情報等も表示されます。

>> [gData ovrn] = read_adxl345_log_bin_2byte("adxl345_log_003.bin");
Header detected, Ver:1, Size:8
Rate:1600[Hz], Range:2[G]

以下のようにしてやると、書き込みバッファのオーバーランが起こったかどうか、どこで起こったかが分かります。

 >> find(ovrn == 1)
ans =

    3196
    3197
    6396
    6397
   17532
   17533
   28668
   28669
   39804
   39805
   41468
   41469
   52604
   52605
   63740
   63741

やっぱりオーバランしているようです。
それほど頻繁にデータが抜けているわけではないので、これで進めます。

次に、前に作ったスクリプトFFT時系列データを作ります。

>> gFFT = fft_slice(gData(1:end,2),512);

OctaveGUIを使っていると、変数の状態を見ることができます。

f:id:nokixa:20200906145240p:plain

"ワークスペースウィンドウ"を見ると、gFFTが512x214のサイズで確保されたことが分かります。
時間方向に214点取れたようです。

前と同様にpcolorを使用してグラフ表示しますが、今回は軸の目盛りで縦軸に周波数の値、横軸に秒単位の時間と、希望の値を表示するようにしたいと思います。

まず、横軸の目盛りを示す配列を用意します。

>> x = gData(512:512:end,1);
>> size(x)
ans =

   124     1

512点のFFTを取っているので、gDataの時間データを1/512に間引いたものとしています。サイズも合っています。

次に縦軸、これは時間刻みの逆数を最大周波数として、これを512分割したものになります。

>> y = [1:512]/512/gData(1,1);

これらを使って、pcolorでグラフ表示します。

>> pcolor(x, y, abs(gFFT));

f:id:nokixa:20200906162555p:plain

横軸、縦軸の目盛りが前回と変わりました。
ラベルも付けておきます。

>> xlabel("t[s]");
>> ylabel("Freq[Hz]");

f:id:nokixa:20200906162932p:plain

周波数範囲後半の半分は見る意味がないのと、DC成分は省きたい、前回見たところではそれほど広い周波数範囲は見ないで良さそうだった、というところから、表示範囲を変えます。
pcolorをやり直すと、ラベル等がリセットされるので、もう一度出します。ついでに表示色凡例も出します。

>> pcolor(x, y(3:128), abs(gFFT(3:128,:)));
>> xlabel("t[s]");
>> ylabel("Freq[Hz]");
>> colorbar

f:id:nokixa:20200906163547p:plain

それらしいものが表示されました。
40秒間くらいのデータですが、確か停止状態→発進→停止のログだったかと思います。
ドラッグスター250は5速までありますが、このログを見るとギアチェンジしたと思われる箇所(回転数が急激に落ちるところ)が4か所あるので、それっぽいです。
いい感じのグラフ表示ができました。

このグラフを見ると、はっきりした周波数のピークが見られるので、このピーク周波数が回転数に対応すると思ってよさそうです。

別データ

前のデータは40秒程度のデータですが、より長時間取ったデータをグラフ表示しました。

f:id:nokixa:20200906170541p:plain

6分ぐらいのログデータになっています。
結構汚くなっているのと、なぜか値の範囲を変更しないと波形が見えませんでした…
うまくグラフが見えなかったので、凡例も消しました。

>> caxis([0 20]);
>> colorbar("off");

これはまた後で表示を見直してみようかな。

まとめ

グラフ表示ができ、どんな加速度が測定できているかという様子が見えてきました。

ここから回転数の推定をしていこうと思いますが、1つ課題として、これが正しいのかどうかを評価する、ということがあります。
この評価方法も考えていこうと思います。