The Negligible Lab

Astray in the forest of electrical and electronic circuits. Adrift in the gap between time and frequency domains. 独立独歩を是とする。

M5Stackで遊ぶ ～浮遊静電容量で系統周波数を測る～

電子工作 M5Stack Googleスプレッドシート

はじめに

概論

電力系統の周波数は，言わずもがな，東日本では50 Hz，西日本では60 Hzです。しかし，いつもピッタリ50 Hzや60 Hzとなっているわけではなく，それぞれ50 Hz，60 Hzを中心としながらも，常に揺らいでいます。発電が消費よりも大きければ周波数は上昇し，逆に消費が発電よりも大きければ周波数は低下します。現状，水車や蒸気タービンで回転している同期発電機が大勢を占めていることから，負荷が重ければ減速し，軽ければ加速する──とも言い換えられますね*1。

もちろん，50 Hz，60 Hzから離れ過ぎると問題が生じるので，地域によって異なりますが，±0.2 ~ ±0.3 Hzに収まるよう（±0.1 Hzの滞在率が概ね99%以上*2），水車に流す水量や，ボイラに投入する燃料を調整したり，周波数を観測して特定の発電所を制御したり，時刻，曜日，季節に応じて運転・停止する発電所を選択したりといった措置が取られています。

実際どうなの？

諸外国には，電力系統の周波数をリアルタイムにインターネット上に公開しているTSO (transmission system operator)もあります。例えば，スウェーデンのSvenska Kraftnätなどがそうです*3。

しかし，日本ではこれは実施されていません。電力広域的運営推進機関(OCCTO)ホームページでも，地域間連系線（周波数変換所や直流送電を含む）の融通量を見ることはできますが，周波数は見ることができません。実際，どういう風に周波数が変動しているのか？

もちろん，インターネット上には先行研究があります。例えば，とものさんのブログでは同じく50 Hz系統の周波数を測定した例が記載されています。

tomono.tokyo

tomono.tokyo

とは言え，やっぱり自分で作った装置で見てみたいではありませんか…？

本記事では，M5Stackにて系統周波数を測る試みについて記載する他，M5StackからGoogleスプレッドシートにデータを送信・蓄積する方法について，最小構成のプログラムを示しながら概説します。完成イメージは図1のようになります。「浮遊静電容量式M5Stack電力系統周波数監視装置」，通称「FreqFloat」と呼ぶことにしました。

f:id:s-inoue2010:20210907182350j:plain:w600
図1: M5Stackで電力系統周波数をフローティングで測る

www.youtube.com
動画1: 浮遊静電容量式M5Stack電力系統周波数監視装置FreqFloatが動作している様子

やっぱり，アナログメータも欲しいですよね。

目次

アイディア

100 Vのコンセントから変圧器やフォトカプラなどを通してマイコンで電圧を検知し，周期を計測したりPLL (phase-locked loop)したりすることで周波数を測ることができそうです。しかし…100 Vなどという超高電圧に関わる工作など，危ないではないか*4…！ 何とか，系統に繋がないで周波数を測ることができないかと考えていたところ，こんな現象を思い出しました。

f:id:s-inoue2010:20210829013938j:plain:w600
図2: 浮いているプローブを触った場合のオシロスコープの波形

図2は，オシロスコープのプローブが何処にも触っていない（浮いている）時に，プローブ先端を指で触った場合の波形です。オシロスコープを触った経験のある方であれば覚えがあるのではと思います。そう，系統電圧（と同じ周波数と思われる波形）が見えています…！ オシロスコープを見ていると，振幅は指に加える力にも依存しますが，周波数は変わりません。これなら少なくとも周波数は検知できる…！ 恐らく，人体が持つ（対地）静電容量が関係しているのではないかと推測しますが，ある程度の大きさの金属塊など，同様の静電容量を用意すれば，人間が触っている必要はないと思われます。

また，オーディオ機器などで音声入力の端子を手で触った際にブ～ンといういわゆるハムノイズが聞こえることがありますね。これも同じ現象かと考えます。

マイコンのアナログ入力（に設定した）ピンも，オシロスコープと同様に高い入力抵抗と小さな入力容量を持つはずですので，マイコンで電力系統の周波数を常時監視できるのではないか…？ 考えるより先に作ってみよう…！

M5Stackとの出会い

6月のある日，髪を切りに1,000円カットのお店に出掛けました。私は3ヶ月に1回，1,000円カットに通っているだけです。野球少年かというほどバサッと切り，そしてボサボサになるまで伸ばし続けます。このように，毎年ほぼ4回しか理髪店に行かずに済んでいるのですが*5，その結果，毎月5,000円の美容室に通っている人と比べると，年間数万円もの差額が生じるではないか！ と言うことで，その1,000円カットの帰りに秋月電子通商まで足を延ばしてみると，いつの間にかM5Stack Basicを購入していました。

M5Stackは中国深圳に本拠を置くM5Stack社が開発した，言わば液晶ディスプレイ & ケース付きマイコンボードです。CPU (SoC)には同じく中国のEspressif Systems社によるESP32-D0WDQ6が用いられています。ESP32は240 MHzのデュアルコアで，しかもWiFiとBluetoothに対応しています。Raspberry Piとは異なり，Linuxなどの一般的なOSは載っておらず，ベアメタルでの開発に主眼を置いているようです*6。プログラムの開発にはArduino IDEやUIFlow*7が使えます。

まずは，下島健彦先生の「みんなのM5Stack入門」で勉強させて頂きました。

books.rakuten.co.jp

そうこうしているうちに，何故かM5Stack Grayが増殖してしまいました😅

f:id:s-inoue2010:20210830183033j:plain:w600
図3: M5Stack GrayとM5Stack Basic

個人的には，「あ，こんなの試してみたい…」と（妙な）アイディアが出てきた際に，サクっと試すことができる“実験場”として，とても面白いのではと思っています。

予備実験

アナログ入力ピンから見える波形を確かめる

さて，図2と同じようなことをM5Stackで試してみました。320 × 240の液晶ディスプレイで2周期分(40 ms)を表示するように，40 ms / 320 = 125 μs周期でサンプリングします*8。ここでは，アナログ入力端子としてGPIO35を選びました。結果，図4のように，半波整流のような波形を観測できました。う～ん，当たり前ですが，電源が片極性なのでA-D変換器も正の電圧しか見えません。正負対称の波形は見えず，半波整流のような波形しか取得できないようです。しかしこれでも周波数を測ること自体は可能ですね。また，指を離しても，振幅は下がるものの，場所によっては十分な大きさの波形を観測できます。

f:id:s-inoue2010:20210830183303j:plain:w600
図4: 浮いているアナログ入力端子を指で触った場合

f:id:s-inoue2010:20210901184730j:plain:w600
図5: M5Stackの向かい合ったピンは接続されている

図5に示すように，M5Stackの裏面にはピンの説明が書かれています（ゴム足を換えています）。左右側面の向かい合う端子は内部で接続されており，また，上下面の向かい合う端子も内部で接続されています。しかも，左と下はピンヘッダ，右と上はピンソケットになっており，外部回路に接続する際に，ジャンパワイヤのオス・メスのいずれでも接続できるようになっています*9。例えば，I²Cについては下面にSDA，SCLと書かれており，ピンヘッダが出ておりますが，上面のGPIO21，22のピンソケットとも繋がっており，いずれに接続しても機能します。

周期から周波数を測ってみる

とりあえず，0 Vからの立ち上がりを検知することで周期の始まりを捉え，これを50回繰り返すことで50周期分の時間を測ります。その逆数をとって50を掛けれれば，周波数を計算できるはずです。Arduino IDEではmicros()という関数が用意されています。これを使うと，いわゆるtick*10をμs単位で取得できます。50周期の始まりと終わりでtickを取得し，その差分から50周期分の時間をカウントします。図6は，実際にやってみた結果です*11。それらしい値を表示していることが判ります。

f:id:s-inoue2010:20210830185346j:plain:w600
図6: 浮いていても周波数を測れる

置く場所を工夫すれば，人間が触っていなくとも周波数を検知し続けることが可能でした。この辺り，もうすこし研究が必要そうです。

アナログメータを作る

M5Stackには7セグメントLED/LCDを模したフォントがデフォルトで組み込まれており，M5.Lcd.setCursor()関数の3番目の引数を7にすることで使うことができます。図6がまさに使用例です。しかし，やっぱり何か電気的な量を測るのであれば，アナログメータが欲しくなりませんか？ M5Stackには液晶ディスプレイに図形を描画する関数（メソッド）が豊富に用意されています*12。これを使って， 車輪の再発明になるかもしれませんが自分でアナログメータを実装することにしました。図7にアナログメータを表示しているM5Stackを示します。

f:id:s-inoue2010:20210907183246j:plain:w300
図7: アナログメータ

さて，プログラムではひとつひとつの図形を地道に描いていくしかありません。そこで，図8のように，まずは紙の上で考えてみました。

f:id:s-inoue2010:20210908135819j:plain:w600
図8: 紙の上でアナログメータの描き方を考える

M5Stackの液晶ディスプレイ(LCD)は320 × 240（懐かしのQVGA）ですが，目盛り部分は画面外の点である(159, 319)を中心とする中心角56°の円弧とすることにしました。目盛りの中央をもちろん50 Hzとします。針は中ほど少し下から出て，その時の値を指します。また，単位「Hz」を中央に表示します。

これをどのように実装するかですが，C++らしくクラスとして作ることにしました。また，この周波数監視装置以外にも使えるように目盛りの数や単位を変えられるようにすることにしました。まず，表示する数値をfloat型で格納した配列と，単位として表示する文字列を引数として，AnalogMeterというクラスからインスタンスを生成します。ここでは，range = {49.0, 49.5, 50.0, 50.5, 51.0}という配列とHzという単位を指定しました。

生成したインスタンスに対してdrawGauge()というメソッドを実行すると，指定した49 Hzから51 Hzに至る円弧をcos関数とsin関数で描きながら*13，一定間隔毎に目盛りとなる線を引き，その上に49, 49.5, 50, 50.5, 51の数値を書き入れます。ここで，小数点以下が0であった場合，int型に変換して小数点以下を取り除き，さらにフォントを大きくしてから書き入れます。これで準備は完了です。

表示したい値がある場合，生成したインスタンスに対してupdate()というメソッドを適用します（引数はfloat型）。まず，前回描いた針を背景色（デフォルトでは黒）で上書きした後，単位(Hz)を書き入れ，最後に指定した値を指し示す針を（デフォルトでは）赤色で描きます。現在の作りでは，上記のdrawGauge()を実行せずともupdate()を実行できてしまうので，そうなると針だけが表示されます💦

privateメソッドとして，drawLinePolar()というメソッドを作りました。これは，指定した点(x, y)を原点として，極座標で角度，描き始める半径，線の長さを指定します。上述の「目盛り」や「針」を描く際に，それぞれのメソッドから呼ばれます。これは「みんなのM5Stack入門」のスケッチ5.6に記載のdrawlinebyangle()関数を参考にさせて頂きました。

作り始めると，見た目をより綺麗にしようと欲が出てしまいました。文字や針の色などの設定にはセッターやゲッターを使わず，public変数を直に操作するようになっています。改善が必要かもしれませんね。

あ，度からラジアンに変換するradians()という関数が用意されていたのか…💦 使えばよかった…😅

測定結果

付録にて後述するように，測定した周波数のデータをGoogleスプレッドシート*14に蓄積することにしました。図9，10に，Googleスプレッドシートで記録した周波数の測定例を示します。青色は10秒毎のデータ，橙色は3分間の移動平均です。想像していた以上に，数秒 ~ 十数秒単位での周波数変動が起きており，±0.1 Hzの範囲の中では大きく動いていることが判りました。こんな風になっているとは思っておらず，当初，ノイズやプログラムの誤りかと疑いましたが，2台のM5Stack（と別途購入したESP32-DevKitC）で同じ測定結果が得られることなどから，そう大きく誤ってはいないと推測しています。

f:id:s-inoue2010:20210901174529p:plain
図9: 60分間の周波数の測定例

f:id:s-inoue2010:20210901174606p:plain
図10: 24時間の周波数の測定例

毎日観測を続けていると，同じ時刻に同じような挙動を示すことが判ってきます。例えば，正午になると周波数が50.1 Hzを超えるか超えないかまで急上昇することが多いです。これは，いろいろな産業がお昼休みに入り，負荷が軽くなるためと思われます。また，何故だかわかりませんが，15時の直前に周波数が高くなり，その後，だらだらと（3分移動平均が）50 Hz未満の期間が続いたりします。また，当然これは季節によっても変わってくるのではと予想します。

まとめ

系統に繋がないで系統の周波数を測るため，M5Stackのアナログ入力ピンに人体や金属塊などの（対地）浮遊静電容量を接続し，得られた半波整流波形から周波数を測定できることを確認しました。また，得られたデータをGoogleスプレッドシートに格納し，グラフとして観測できるようにしました。

今後も観測を続け*15，例えば冬の重負荷時などに備えたい（？）と考えます。今のところ，地震や大停電に伴う周波数変動は観測されていません。何かが観測出来たらまたご報告致します。

追記: 地震発生に伴う周波数変動について

2021年10月7日22時41分頃に千葉県北西部を震源とするマグニチュード5.9の地震があり，埼玉と東京の一部で震度5強を観測しました。我が家でも震度5弱の揺れでした。この地震が発生した際に，50 Hz系統の周波数に擾乱が生じ，本記事のFreqFloatでも観測できました（今ではM5StackではなくESP32 DevKitCで動いています）。図11に観測結果を示します。

f:id:s-inoue2010:20211015153902j:plain
図11: 2021年10月7日22時41分頃の地震に伴う電力系統周波数の変動

気象庁の報道発表によれば，地震波検知は22時41分34.8秒でした*16。図13を見ると，22時41分40秒から（いや，変化が現れるのはその次の22時41分50秒からですね）周波数が上昇し始め，22時42分10秒に最大値50.21 Hzに達しました。大都市圏で地震が発生したため，電源（発電所）というよりも負荷（需要家）が解列したことが原因で周波数が上昇したものと推測します。その後，数分を掛けて擾乱は落ち着いております。もし，電源地域で大きな地震が起きていたら，発電所の解列に伴って周波数は下がっていくのではないかと考えます。

なお，ESP32 DevKit CはNTPプロトコルでNICTのNTPサーバと同期しているので時刻は正確です*17。

付録1: Googleスプレッドシートにデータを蓄積する

Googleスプレッドシートとスクリプトの概要

M5Stackの液晶ディスプレイに表示するだけでは，その瞬間の系統周波数しか分かりません。M5Stack単体でグラフを表示しても良いのですが，パソコンやスマートフォンなどからもデータにアクセスできた方が便利です。そんな時に使えるのが，Googleスプレッドシートです。

GoogleスプレッドシートにはExcel VBAに相当する「App Script」と呼ばれる機能があります。App Scriptで作成したスクリプトでセルの値やグラフを操作することができるのですが，Google Drive上にスプレッドシートを置いた場合，Web APIを作ることができます。要するに，あるURLをHTTPでGETしたりPOSTしたりすると，スクリプトにデータを渡すことができます（データを受け取ることも可能と思われますが，試したことがありません😅）。私は，M5Stackの中で時刻と周波数のデータをJSONとして生成し，これをWeb APIにPOSTすることによって，データをGoogleスプレッドシートに渡すことにしました。

スプレッドシートのスクリプトでは，渡されたデータから時刻と周波数を抽出し，これをデータ格納用のシートの最終行に追加します。また，グラフの範囲を更新し，過去60分間と過去24時間の周波数変動を示す2つのグラフを表示するようにしました。

最小構成のプログラム例

これだけでは分かり難いので，Googleスプレッドシートにデータを渡すための最小構成のM5Stack側，Googleスプレッドシート側のプログラムを示します。

この例では，loop()関数内でint型の変数iを0から2秒毎にカウントアップしていき，これをJSONにしてGoogleスプレッドシートのWebアプリにPOSTします。送っている変数iの値は，M5Stackのディスプレイに表示します。Googleスプレッドシート側では，POSTされてきたJSONを解釈し，変数iに相当する数値をアクティブなシートの最終行に追加します。int型の変数iの箇所をいろいろと変えてあげれば，つまり，JSONを必要に応じて作り変えれば，様々なデータをGoogleスプレッドシートに蓄積できます。

M5Stack (Arduino IDE)側のスケッチ

M5Stack側(Arduino IDE)は下記のようになります。プログラム中の「XXXXXXXXXX」「YYYYYYYYYY」「ZZZZZZZZZZ」には，皆さまの環境に合う文字列をはめ込んで下さい。

#include <M5Stack.h>
#include <WiFi.h>
#include <HTTPClient.h>
#include <ArduinoJson.h>

// Declare WiFi constants
const char* ssid = "XXXXXXXXXX";
const char* password = "YYYYYYYYYY";

// Declare a function to report data to Google Spreadsheet
void reportGSS(int dataToReport);
const char* apiURL = "https://script.google.com/macros/s/ZZZZZZZZZZ/exec";
int httpCode;

// Obtain HTTP client handler
HTTPClient http;

/*
   The setup function
*/
void setup() {
  // Initialize M5Stack
  M5.begin();
  M5.Power.begin();

  // Connect to WiFi
  WiFi.begin(ssid, password);
  while (WiFi.status() != WL_CONNECTED);

  // Setup display
  M5.Lcd.setTextSize(0);
}

/*
   The loop function
*/
void loop() {
  static int i = 0;

  // Show what number is to be reported
  M5.Lcd.setCursor(32, 32, 7);
  M5.Lcd.print(i);
  
  // Store data from 0 to 99 to Google Spreadsheet every 30 sec
  reportGSS(i);
  i = (i + 1) % 100;

  delay(2000);
}

/*
   Report to Google Spreadsheet
*/
void reportGSS(int dataToReport)
{
  char pubMessage[64];
  
  // Create JSON message
  StaticJsonDocument<500> doc;
  JsonObject object = doc.to<JsonObject>();
  object["data"] = dataToReport;
  serializeJson(doc, pubMessage);

  // Report to Google Spreadsheet API
  http.begin(apiURL);
  httpCode = http.POST(pubMessage);
}
Googleスプレッドシート側のスクリプト

Googleスプレッドシート側のスクリプトは下記のようになります。

function doPost(e) {
  // Extract POST'ed data
  var params = JSON.parse(e.postData.getDataAsString());
  var data = params.data;

  // Obtain spreadsheet and sheet handler
  var ss = SpreadsheetApp.getActiveSpreadsheet();
  var sheet = ss.getActiveSheet();

  // Append data to the last row
  sheet.appendRow([data]);
}

実行結果

図12にM5Stackの画面を，図13にGoogleスプレッドシートの画面を示します。

f:id:s-inoue2010:20210901181635j:plain:w600
図12: 最小構成のGoogleスプレッドシートのプログラム例を実行しているM5Stackの画面

f:id:s-inoue2010:20210901181733p:plain
図13: Googleスプレッドシートに書き込まれた数値

M5Stackに「12」と表示されているので，int型の変数i = 12であり，Googleスプレッドシートも12まで書き込まれていることが判ります。このまま放置しておくと，i = 99までカウントアップしてまたi = 0に戻ります。

付録2: 本記事のプログラム

本記事のプログラムをGitHubに置きました。拙いですが，M5Stack本体と僅かなジャンパワイヤがあれば動作を試すことができます。ご興味のある方はどうぞご覧下さい。

github.com

付録3: 同様の試みについて

2021年9月に入ってから，本記事と同じようにM5Stackを電力系統に接続することなく系統周波数を計測しようという試みを見つけました。

protopedia.net

のるさんによるこの作品では，同じようにGoogleスプレッドシートにデータを蓄積している他，M5Stack単体でもグラフを表示できるようになっています。2021年2月ごろに発生した停電をきっかけとして，かつ，ハムノイズを利用して，M5Stack単体で周波数を計測する作品を作られたようです。本記事よりかなり早い時期からの着眼❗ う～ん，やられました❗ しかし，同じことを考えている方がいらっしゃったというのは面白いです。

*1:太陽光や風力など，系統連系インバータを介した発電設備が増えてくるとこの挙動が変わってくるとして，その対策に関する研究が最近では盛んです。

*2:https://www.occto.or.jp/iinkai/chouseiryoku/jukyuchousei/2020/files/jukyu_shijyo_19_02_02.pdf

*3:The control room | Svenska kraftnät

*4:実際には，100 V/6.3 Vの変圧器を持っていたりするのですが，工作する時間が取れません…。

*5:ええ，貧乏性です…💦

*6:もちろんOSを載せることもできますが，Raspberry Piのようなデスクトップとしての使用は難しいでしょう。また，Arduino IDEで使う場合，FreeRTOSが無意識のうちに載っています。

*7:これも試してみたいです。

*8:定周期割り込みについては別途書きたいと思います。

*9:Arduinoはピンソケット派，Raspberry Piはピンヘッダ派ですよね。一方，M5Stackは液晶ディスプレイやmicroSDのスロットなど，もともとマイコン(ESP32)から見た外付けデバイスが多いため，ユーザーが自由に使えるピン数は先の2者に比べると少なめです。ピンヘッダとピンソケットの両方を使えるのはそれを逆手にとっていますね。なお，Nucleoでは上面にArduino互換の配置をもつピンソケットと，Nucleo独自の配置をもつピンヘッダが用意されており，やはり両方使えるようになっています。

*10:電源を投入してからの経過時間です。

*11:図5は後述のGoogleスプレッドシートへのデータ送信の機能も追加した後の画面です。したがって，時計やIPアドレスをも表示しています。

*12:らびやんさんのLovyanGFXを使わせて頂いた方が良かったのかもしれません。

*13:もちろん，math.hをincludeします。

*14:正確には「Google スプレッドシート」と間にスペースが入るようですが，日本語として不自然なので，スペースなしの表記で書かせて頂きます。

*15:M5Stackを使ってしまうのはもったいないので，ESP32-DevKitCにプログラムを移設しました

*16:気象庁｜報道発表資料

*17:プログラム開始時にsetup()関数内で一度だけだから少しずつズレていくかも…。