SPRESENSE で超音波風速計を作ってみようとした話

SPRESENSE で超音波を扱うのも慣れてきたので、風速計が作れるんじゃね？と思いトライしてみました。しかし、風速計の世界はそんなに甘いものではないということを思い知る結果となりました。

# 超音波風速計の原理

超音波風速計の原理は次の通りです。相対する超音波トランスデューサとレシーバで超音波の到達時間の差を測ることで風速が分かります。
<br/>

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/598183/cdd2c757-161d-6a84-60b6-670b155ed6d1.png)

これを十字方向に組むことで、風の方向と速度がわかるというわけです。

<img src ="https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/598183/56438cf9-35e7-357f-9357-c1fe31740b85.png" width=350/>


# SPRESENSE の構成
SPRESENSEでは２つの超音波トランスデューサとレシーバを接続するので次のような構成となります。超音波トランスデューサの出力は交流にするためコンデンサでACに変えています。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/598183/b0af0fff-d3d1-9e22-65f8-9d8033f77610.png)

# SPRESENSEの超音波計測用プログラム
SPRESENSEのプログラムは[超音波測距のプログラム](https://qiita.com/TaroYoshino/items/9b3cf5c5b52a2cab5eeb)を拡張したものです。しきい値の判定は最大値に設定しました。半分にすると誤差成分が大きいため、値が大きくばらつくためです。

```c++
#include <FrontEnd.h>
#include <MemoryUtil.h>
#include <arch/board/board.h>

#define SAMPLE_SIZE (768)

//*****
#include <sys/ioctl.h>
#include <stdio.h>
#include <fcntl.h>   
#include <nuttx/timers/pwm.h>
int fd0;
struct pwm_info_s info0;
int fd1;
struct pwm_info_s info1;//*****

//***
uint32_t fire_time = 0;
uint32_t last_time = 0;
uint32_t current_time = 0;
//***

FrontEnd *theFrontEnd;

const int32_t channel_num = AS_CHANNEL_STEREO;
const int32_t bit_length  = AS_BITLENGTH_16;
const int32_t sample_size = SAMPLE_SIZE;
const int32_t frame_size  = sample_size * (bit_length / 8) * channel_num;
bool isErr = false;

void frontend_attention_cb(const ErrorAttentionParam *param) {
  Serial.println("ERROR: Attention! Something happened at FrontEnd");
  if (param->error_code >= AS_ATTENTION_CODE_WARNING) isErr = true;
}

static bool frontend_done_cb(AsMicFrontendEvent ev, uint32_t result, uint32_t detail){
  UNUSED(ev);  UNUSED(result);  UNUSED(detail);  
  return true;
}


static void frontend_pcm_cb(AsPcmDataParam pcm) {
  static uint8_t input[frame_size];

  frontend_signal_input(pcm, input, frame_size);
  signal_process((int16_t*)input, (int16_t*)NULL, sample_size);
  return;
}

void frontend_signal_input(AsPcmDataParam pcm, uint8_t* input, const uint32_t frame_size) {
  /* clean up the input buffer */
  memset(input, 0, frame_size);

  if (!pcm.is_valid) {
    Serial.println("WARNING: Invalid data! @frontend_signal_input");
    return;
  }

  //***
  digitalWrite(4, LOW);
  last_time = current_time;
  current_time = micros();
  //***

  if (pcm.size > frame_size) {
    Serial.print("WARNING: Captured size is too big! -");
    Serial.print(String(pcm.size));
    Serial.println("- @frontend_signal_input");
    pcm.size = frame_size;
  } 
  
  /* copy the signal to signal_input buffer */
  if (pcm.size != 0) {
    memcpy(input, pcm.mh.getPa(), pcm.size);
  } else {
    Serial.println("WARNING: Captured size is zero! @frontend_signal_input");
  }  
}

void signal_process(int16_t* stereo_input, int16_t* stereo_output, const uint32_t sample_size) {

  //***
  static bool pwm0 = true;
  static int16_t mono_A_input[SAMPLE_SIZE];
  static int16_t mono_B_input[SAMPLE_SIZE];
  int16_t *mono_input;

  int i = 0;
  for (int n = 0; n < SAMPLE_SIZE; ++n) {
    mono_A_input[n] = stereo_input[i++];
    mono_B_input[n] = stereo_input[i++];
  } 

  /* calculating a distance for the past frame */
  if (pwm0) mono_input = mono_B_input;
  else      mono_input = mono_A_input;
 
  uint32_t distance_time = 0; 
  uint32_t trigger_time  = 0;
  uint32_t alpha0 = 1280; // adjustment parameter for each device (microseconds)
  uint32_t alpha1 = 1060; // adjustment parameter for each device (microseconds)
  uint32_t offset_time = fire_time - last_time;
  if (pwm0) offset_time += alpha0;
  else      offset_time += alpha1; 

  int16_t max = 0;
  uint32_t t  = 0;
  for (uint32_t n = 0; n < sample_size; ++n) {
    int16_t amp = abs(mono_input[n]);
    if (amp > max) {
      max = amp; 
      t = n;
    }
  }

  trigger_time = (sample_size - t)*1000000/(192000); // microseconds
  distance_time = trigger_time - offset_time ; 

  // Serial.println("trigger_time: " + String(trigger_time) + " offset time: " + String(offset_time));

  const int ave_count = 1000;
  const float distance = 0.10;
  static uint32_t ave_distance_time0 = 0;
  static uint32_t ave_distance_time1 = 0;
  static uint32_t mem_distance_time0 = 0;
  static int counter0 = 0;
  static int counter1 = 0;
  if (distance_time > 0 && distance_time < 5333 /* about 2m */) {
    // Serial.println(String(offset_time) + "," + String(distance_cm));
    if (pwm0) {
      ave_distance_time0 += distance_time;
      if (counter0++ == ave_count) {
        ave_distance_time0 /= ave_count;
        mem_distance_time0 = ave_distance_time0;
        counter0 = 0;
        ave_distance_time0 = 0;
      }
      pwm0 = false;
    } else {
      ave_distance_time1 += distance_time;
      if (counter1++ == ave_count) {
        ave_distance_time1 /= ave_count;

        //**** Calucurate Wind Speed
        float wind_velocity  = (distance / 2)*(1000000./float(mem_distance_time0) - 1000000./float(ave_distance_time1));
        float sonic_velocity = (distance / 2)*(1000000./float(mem_distance_time0) + 1000000./float(ave_distance_time1));
        float calc_distance = sonic_velocity * float(mem_distance_time0 + ave_distance_time1) * 0.000001 / 2;
        Serial.println("time left to right: " + String(mem_distance_time0) + " usec");      
        Serial.println("time rihgt to left: " + String(ave_distance_time1) + " usec");      
        Serial.println("wind_velocity : " + String(wind_velocity) + " m/s");
        Serial.println("sonic_velocity: " + String(sonic_velocity) + " m/s");
        Serial.println("calc distance : " + String(calc_distance) + " m");
        //****

        counter1 = 0;
        ave_distance_time1 = 0;
      }
      pwm0 = true;
    }       
  }

  /* firing ultra sonic */
  //delayMicroseconds(1000);
  fire_time = micros();
  digitalWrite(4, HIGH);
  if (pwm0) {
    ioctl(fd0, PWMIOC_START, 0);
    delayMicroseconds(1000);
    ioctl(fd0, PWMIOC_STOP, 0);     
  } else {
    ioctl(fd1, PWMIOC_START, 0);
    delayMicroseconds(1000);
    ioctl(fd1, PWMIOC_STOP, 0);     
  }

  return;
}


void setup() {
  Serial.begin(115200);

  //***
  pinMode(4, OUTPUT);
  /* Setting PWM 40kHz for PWM0 (D06)*/
  Serial.println("setup /dev/pwm0");
  fd0 = open("/dev/pwm0", O_RDONLY);
  info0.frequency = 40000; // 40kHz
  info0.duty      = 0x7fff;
  ioctl(fd0, PWMIOC_SETCHARACTERISTICS, (unsigned long)((uintptr_t)&info0));
  ioctl(fd0, PWMIOC_STOP, 0);  

  /* Setting PWM 40kHz for PWM1 (D05)*/
  Serial.println("setup /dev/pwm1");
  fd1 = open("/dev/pwm1", O_RDONLY);
  info1.frequency = 40000; // 40kHz
  info1.duty      = 0x7fff;
  ioctl(fd1, PWMIOC_SETCHARACTERISTICS, (unsigned long)((uintptr_t)&info1));
  ioctl(fd1, PWMIOC_STOP, 0);  
  //***  


  /* Initialize memory pools and message libs */
  initMemoryPools();
  createStaticPools(MEM_LAYOUT_RECORDINGPLAYER);

  /* setup FrontEnd and Mixer */
  theFrontEnd = FrontEnd::getInstance();
  
  /* set clock mode */
  theFrontEnd->setCapturingClkMode(FRONTEND_CAPCLK_HIRESO);

  /* begin FrontEnd and OuputMixer */
  theFrontEnd->begin(frontend_attention_cb);
  Serial.println("Setup: FrontEnd began");

  /* activate FrontEnd and Mixer */
  theFrontEnd->setMicGain(0);
  theFrontEnd->activate(frontend_done_cb);
  delay(100); /* waiting for Mic startup */
  Serial.println("Setup: FrontEnd activated");

  /* Initialize FrontEnd */
  AsDataDest dst;
  dst.cb = frontend_pcm_cb;
  theFrontEnd->init(channel_num, bit_length, sample_size, AsDataPathCallback, dst);
  Serial.println("Setup: FrontEnd initialized");

  /* Unmute */
  board_external_amp_mute_control(false);
  theFrontEnd->start();
  Serial.println("Setup: FrontEnd started");
}

void loop() { 
  if (isErr == true) {
    board_external_amp_mute_control(true); 
    theFrontEnd->stop();
    theFrontEnd->deactivate();
    theFrontEnd->end();
    Serial.println("Capturing Process Terminated");
    while(1) {};
  }
}

```

パラメータの alpha0 と apha1 はそれぞれのセンサー毎に設定する必要があります。実際にオシロスコープで測定してみると感度が異なり振幅と信号の広がりが違うのが分かります。

**RX0の受信信号**
![NewFile1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/598183/90912a34-c179-979a-72b0-042e4c84391a.png)

**RX1の受信信号**
![NewFile2.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/598183/fff936c1-a66c-ddd5-d5f9-8a1a1b46c481.png)

調整は無風条件下で、”time left to right” と"time right to left"の値を一致させるのと室温に応じた音速にできるだけ近づけることで行います。音速は次の式で表せます。

C（音速m/s）= 331.5 + 0.6 x T (室温℃)


# 実際に測定をしてみる

測定をする際は、風が吹き抜けるようにセンサーを浮かせておく必要があります。机の上だと風が反射し乱流となり、思うように測定ができません。最初それに気が付かずにはまりました。

![DSC02560.JPG](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/598183/e97c92f4-c8eb-1488-36f5-00fd0cfc18f2.jpeg)

<br/>
測定は厳密に行わず、ファンの風を斜め上からあてるようにして測定をしました。こんな感じです。

![DSC02562.JPG](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/598183/6f3385b9-7295-ed12-65af-bcd7d63e7b25.jpeg)

測定結果は下記の通りです。

**無風時**
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/598183/1da1f69a-b371-b41a-4706-4914a3d7c647.png)


**順風時**
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/598183/1f032ed3-6c32-3cef-414c-25d7bdddcb75.png)


**逆風時**
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/598183/bfb154fe-adf3-724c-f76c-a5ea033a7b91.png)

<br/>

これらのグラフを見ればわかる通り、まったく使い物になりません。0.1 (m) で風速 1 (m/s) の違いを測定するには、

0.1 (m) ÷ 340 (m/s) = 0.000294117647 (s)　...①
0.1 (m) ÷ 341 (m/s) = 0.000293255132 (s)　...②

①-② = 0.000000862515 = 0.86 μ秒

なので、192kHz ではまったく精度が足りません。少なくとも1.2MHz以上のA/Dコンバータが必要ということです。う~ん、残念。とはいうものの何かしら手はないか考えてみたいところです。


# 作ってみた感想
超音波風速計は、超音波センサの軸合わせがかなりシビアでちょっとズレただけでも測定時間が大きくずれてしまいます。またセンサ毎の特性も異なるので、作るとなるとメカの精度や調整がかなり大変だなということは感じました。世の中の超音波風速計がかなりごつい作りになっているのも納得です。今回作ったものは、調整値の設定が大変だったのですが、どうもプログラムにも何か問題がありそうな気がしています。もう少し実験を繰り返して、傾向くらいは安定して出せるものに仕上げたいと思います。