気象庁XMLの台風解析・予報情報電文を使用して台風情報を描画する備忘録 その1

この記事は、[防災アプリ Advent Calendar 2022](https://adventar.org/calendars/7488) 22日目の記事です。


## はじめに
気象庁が発表する気象情報や地震情報をプログラムで処理したくなったときに便利なのが、気象庁防災情報 XML（以下、気象庁 XML）です。
ありがたいことに [Atom フィード](http://xml.kishou.go.jp/xmlpull.html)[^1] が提供されており手軽に利用させていただくことができます。

気象庁 XML ではどのような情報が提供されているのか、XML のフォーマットはどうなっているのかといった気象庁 XML を取り扱う際に必要な情報は、気象庁の [気象庁防災情報XMLフォーマット 技術資料](http://xml.kishou.go.jp/tec_material.html) のページにあります。（ここまでテンプレ）

台風解析・予報情報電文[^2]を使用してブラウザ上で台風の現在位置や進路予報図を描画したので処理方法をメモしておこうと思います。
本当にメモ程度なので説明不足な点が大いにあることを予めお断りしておきます。
また、サンプルコードのつくりも適当なので各自で改良してお使いください。

今回は、台風の現在位置を表示するところまでを説明します。


## 現在の強風域、暴風域、中心位置を描画
`Body/MeteorologicalInfos[0]/MeteorologicalInfo[0]/Item[0]/Kind` のうち `Property/Type` が「中心」の要素「風」の要素を使用します。
それぞれの `Property[0]` 要素を `center` と `wind` としておきます。

---

はじめに強風域を描画します。
`wind/WarningAreaPart[@type='強風域']/Circle/Axes/Axis` を `axis` と置いておきます。
強風域が存在しないとき[^3]は `axis[0]/Radius` の `@condition` と `@description` に「なし」が設定されるため、この場合は以下の処理を行わないように注意が必要です。
ここまでで強風域の中心位置と大きさが分かったため描画できるかと思いきや、**強風域が非対称な形をしている場合がある**ため状況に応じて処理を分ける必要があります。

まず、強風域が台風の中心位置から等距離で表現されている場合は `axis[0]/Direction` は空タグになり、`@condition` と `@description` に「全域」が入ります。
`axis` の要素数も1つだけになります。
このときは、`center/CenterPart/Coordinate[@type='中心位置（度）']` を中心位置として `axis[0]/Radius[@unit='km']` km の円を描画すれば OK です。
中心位置の経緯度は ISO 6709 形式の文字列になっているため緯度と経度に分解します。
以下は JavaScript のサンプルコードです（厳密に ISO 6709 形式に従っているわけではなく、本電文に出現する経緯度の形式に対応する処理です）。
```js
function splitCoordinate(coordinateString) {
  const match = coordinateString.match(/^(?<latitude>[+-]\d{2}(\.\d+)?)(?<longitude>[+-]\d{3}(\.\d+)?)\/$/);
  if (!match) return null;
  return [parseFloat(match.groups.latitude), parseFloat(match.groups.longitude)];
}
```

次に、強風域が台風の中心距離から等距離で表現ではない場合は、楕円を正円に近似して地図上に表示する必要があります。
この場合 `axis` の要素数が2つになり、長径、短径の順に要素が出現します。

以下は JavaScript のサンプルコードです。
```js
/*
 * coreLatitude: center/CenterPart/Coordinate[@type='中心位置（度）'] で取得した緯度
 * coreLongitude: center/CenterPart/Coordinate[@type='中心位置（度）'] で取得した経度
 * majorAxisRadius: axis[0]/Radius[@unit='km'] の値
 * minorAxisRadius: axis[1]/Radius[@unit='km'] の値
 * majorDirection: axis[0]/Direction の値
 * minorDirection: axis[1]/Direction の値
 * isComposite: 関数内でのみ使用（false）
 */
function calcCenterCoordinate(coreLatitude, coreLongitude, majorAxisRadius, minorAxisRadius, majorDirection, minorDirection, isComposite = false) {
  const earthRadius = 6371; // 地球の半径 (km)
  const circumference = 2 * Math.PI * earthRadius;
  const latitudePerKilometer = 360 / circumference;
  const longitudePerKilometer = 360 / (circumference * Math.cos(coreLatitude * Math.PI / 180));
  const averageRadius = (majorAxisRadius + minorAxisRadius) / 2;
  let offsetRadius = majorAxisRadius - averageRadius;
  if (isComposite) offsetRadius /= Math.sqrt(2);

  if (majorDirection.length > 1) {
    // 「北東」のように単純な4方位ではない場合は「北」と「東」のようにそれぞれの要素に分解して再帰的に処理
    const result1 = calcCenterCoordinate(coreLatitude, coreLongitude, majorAxisRadius, minorAxisRadius, majorDirection[0], minorDirection[0], true);
    const result2 = calcCenterCoordinate(coreLatitude, coreLongitude, majorAxisRadius, minorAxisRadius, majorDirection[1], minorDirection[1], true);
    return { latitude: result1.latitude, longitude: result2.longitude, radius: result2.radius };
  } else {
    const vector = direction => {
      switch (direction) {
        case '北': return { x: 0, y: latitudePerKilometer * offsetRadius }
        case '東': return { x: longitudePerKilometer * offsetRadius, y: 0 }
        case '南': return { x: 0, y: -latitudePerKilometer * offsetRadius }
        case '西': return { x: -longitudePerKilometer * offsetRadius, y: 0 }
      }
    };
    return { latitude: coreLatitude + vector(majorDirection[0]).y, longitude: coreLongitude + vector(majorDirection[0]).x, radius: averageRadius };
  }
}
```
返ってきた値を使って `latitude` と `longitude` を中心位置とする `radius` km の円を描画します。

---

暴風域は上記の処理を適宜読み替えて同じ処理を行います。

最後に、台風の中心を描画する場合は `center/CenterPart/Coordinate[@type='中心位置（度）']` を使用すれば OK です。


## おわりに
本当は進路予報図を描画するところまでを書きたかったのですが、想像以上にまとめるのが大変だったため別の記事で書きます。[（書きました）](https://qiita.com/oruponu/items/817d949f944e090b45cb)
地図投影法の話[^4]なども盛り込みたかったですが、こちらも機会があれば書きます（気が向いたら追加するかもしれません）。

最後に、今回作成したプログラムをもとに令和3年台風第14号の経路（速報解析による速報値）をアニメーション表示してみたので紹介します。
<blockquote class="twitter-tweet"><p lang="ja" dir="ltr">台風14号の経路（速報値） <a href="https://t.co/6uDlOKrV6w">pic.twitter.com/6uDlOKrV6w</a></p>&mdash; のた (@oruponu) <a href="https://twitter.com/oruponu/status/1439205028665847810?ref_src=twsrc%5Etfw">September 18, 2021</a></blockquote> <script async src="https://platform.twitter.com/widgets.js" charset="utf-8"></script>


[^1]: 情報発表から更新されるまでに数分かかりますが、リクエストの頻度に気を付けたり `If-Modified-Since` を使ったりするなど、お行儀よく取得させていただきましょう
[^2]: より正確には「台風解析・予報情報（５日進路・強度予報）電文」（2022年12月現在）
[^3]: 強風域が存在しないのは気象擾乱が台風の条件を満たしていないとき
[^4]: 気象庁ホームページの台風情報は、リニューアル前はポーラーステレオ図法で表示されていましたが、リニューアル後はメルカトル図法になりました