Processingでお手軽3次元データビジュアライズ（フライトデータ編）

本記事は、[Processing Advent Calendar 2020](https://adventar.org/calendars/4953)の13日目の記事です。

[昨年に続き](https://qiita.com/masakiyamabe/items/5c5ebc1878865b524d4b)、Processingでお手軽に3次元データビジュアライズするシリーズです。

# はじめに
新型コロナウイルスの世界的なパンデミックによって、飛行機に乗る機会が全くなくなってしまいました。飛行機＆空港好きの自分としては寂しい限りです。

再び気軽に飛行機に乗れるようになることを願いながら、フライトデータの3次元的な可視化に挑戦してみてはいかがでしょうか？

本記事では、[flightradar24](https://www.flightradar24.com/data/flights)でダウンロードしたフライトデータを元に、Processingで可視化します。(注：ご自身で任意のフライトデータをダウンロードするには有料プランが必要です。フリートライアル有）

フライトデータを可視化するアートとしては、アーロン・コブリンの「[Flight Patterns](http://www.aaronkoblin.com/work/flightpatterns/)」が有名です。

今回は、飛行機１機分の飛行経路を可視化するのでここまでハイクオリティではないですが、サクッとやってみたい方には参考になるかと思います。また、他の地理空間移動情報の可視化にも参考になるかと思います。

この記事では紹介しきれませんがサンプルコードの中では、時刻（Timestamp）も使って、飛行をアニメーションで再現したり、スライダーを使って時間を操作するサンプルもありますので、興味のある方はご覧ください（詳細は記事末尾に）。ちなみ使用したデータは、[羽田発福岡行きJAL321便](https://www.flightradar24.com/data/flights/jl321)です。

サンプルコードは[こちらからダウンロード](https://files.masakiyamabe.com/qiita/visualization_3D_time.zip)できます。
文中のサンプルコードのファイル名と対応しているので適宜参照しながらお読みください。

__完成図__

sketch_flight.pde
![sketch_flight.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/80141/c0fe7dc5-d0bf-67e3-4275-044c45e7520f.png)



# 前提条件
Processingのバージョンと使用するライブラリは次の通りです。マウスで視点操作を行うので、ライブラリはPeasyCamを使います。

* Processing3.5.4
* PeasyCam 302

# 素材の準備

## データ

![flightdata.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/80141/248e06d2-64b8-783c-30df-1ae08c9a8968.png)

Flightradar24からダウンロードしたデータは、時刻（Timestamp、UTC）、飛行位置（緯度、経度）、コールサイン（Callsign）、高度（Altitude）、速度（Speed）、方向（Direction）を持った時系列のcsvファイルです。この記事では、緯度、経度、高度を元に飛行位置をXYZ座標でプロットして線でつなぎ飛行経路を可視化します。単位はあまり気にしなくても大丈夫ですが、気になる方は以下を参考にしてください。
高度の単位は「フィート」（1フィート＝0.3048m）
速度の単位は「ノット」（1ノット＝時速1.852km）
方向の単位は「°」（北＝0°東＝90°南＝180°西＝270°）　

### 緯度と経度

まず緯度と経度のおさらいです。

__緯度：latitude__

赤道上は緯度0°、北に向かって緯度は増加、北極点は緯度（北緯）90°

__経度：longitude__

イギリスのグリニッジは経度0°、東に向かって経度は増加、グリニッジの真裏の太平洋上で経度180°東京の緯度経度は、緯度：35.69°　経度：139.69°
[さらに詳しい説明](https://www.gsi.go.jp/KIDS/PAMPHLET/p9.htm)

## 地図画像

使用する地図の画像はProcessingで読み込めればどのようなものでも大丈夫ですが、緯度、経度、それぞれどんな範囲の画像なのかあらかじめ知っておくと、緯度経度から、x座標、y座標へ変換が簡単になります。

例えば、今回使用している画像のサイズと緯度経度の範囲の関係は次のようになっています。

![1000-3.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/80141/6c395e95-a71e-e568-b539-0dc38d8f7ae2.png)

この関係から、例えば、東京の緯度・経度（35.69,139.69）はmap関数を使って次のように地図上のx座標、y座標に変換できます。この変換式は後で使います。

```java

//東京の経度から地図上のx座標への変換
float x = map(139.69, 122.0, 148.0, 0, 822);
//東京の緯度から地図上のy座標への変換
float y = map(35.69, 24.0, 46.0, 850, 0);
```

この変換方法はあくまで近似的なもので、厳密に正確な位置にマッピングしたい場合や、広範囲をカバーしたい場合には向きません。特にメルカトル図法では緯度とY座標は非線形なので高緯度の場所ではズレが大きくなります。今回のような日本全体で大まかに見るような場合は十分かなと思います。

# 実装

## 準備:座標軸と図の関係

具体的な実装の紹介の前に、座標軸と図の関係を整理しておきます。下図のように、x軸（西→東）、y軸（北→南）、z軸（地表→上空）となります。

![1000-4.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/80141/10c40b1e-c1e4-f1d5-10be-b1436f53f20f.png)

## 地図をテクスチャとして張り込む

データや画像などの素材を読み込んでおき、地図をテクスチャとした面の描画部分は次の通りです。
サンプルコード：sketch_flight_step1

```java:sketch_flight_step1.pde
import peasy.*; //PeasyCamライブラリをインストールしてインポート
PeasyCam cam;
Table table;
PImage japanmap;
PFont font;

//経度の範囲
float lonFrom = 122.0;
float lonTo = 148.0;

//緯度の範囲
float latFrom = 24.0;
float latTo = 46.0;

void setup() {
  size(1024, 768, P3D);

  //CSVファイルをTableに読み込む
  table = loadTable("JL321_229092b2.csv", "header");

  //地図画像を読み込む
  japanmap = loadImage("map-japan.png");

  cam = new PeasyCam(this, 700);

  colorMode(HSB);
  font = createFont("Osaka", 12);
  textFont(font);
}


void draw() {
  background(0);
  translate(-japanmap.width/2, -height/2);

  //地図の描画
  noStroke();
  beginShape();
  texture(japanmap);//テクスチャに読み込んだ画像を指定
  vertex(0, 0, 0, 0, 0);//vertex(x,y,z,textureのu,textureのv)
  vertex(japanmap.width, 0, 0, japanmap.width, 0);
  vertex(japanmap.width, japanmap.height, 0, japanmap.width, japanmap.height);
  vertex(0, japanmap.height, 0, 0, japanmap.height);
  endShape();
```
テクスチャを指定した後のvertex関数は引数が5つありますが、前の3つはxyz座標で後の２つは面の頂点とテクスチャ画像の対応する位置を関連づけるために必要です。図にすると次のようになります。

![1000-7.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/80141/ad42afb7-3401-6736-cfd4-5c1a9ede420f.png)

## 緯度・経度をXY座標に変換する方法

例えば、東京の緯度経度をxy座標に変換し、地図上にboxを表示するような場合、以下のように変換を行います。

```java:
 //東京の緯度は
 float latTokyo = 35.69;
 //東京の経度は
 float lonTokyo = 139.69;
 
 //東京のx座標、y座標は
 float x = map(lonTokyo, lonFrom, lonTo, 0, japanmap.width);
 float y = map(latTokyo, latFrom, latTo, japanmap.height,0);
 
 fill(255,0,0);
 
 //東京の場所にboxを表示
 pushMatrix();
 translate(x, y, 0);
 box(10,10,10);
 popMatrix(); 
```
実行結果
![1000-8.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/80141/90681298-4c1e-fa1d-2ba2-1cefc6b06603.png)

## csvからデータを取り出して経路を表示

フライトデータのcsvから取り出した緯度、経度、高度を対応するXYZ座標に変換して、経路全体を描きます。通過した点を白いドットで表示し、ドット間は線で繋ぎます。その際、線の色は速度の値を元に色情報に変換します。

サンプルコード「sketch_flight_step1.pde」

```java:sketch_flight_step1.pde
//データから全ての行を取り出す

  noFill();
  strokeWeight(1);
  beginShape();

  //データのある位置を線でつなぐ。線の色は高度から色相に変換。
  for (int i=0; i<table.getRowCount(); i++) {
    float lat = table.getRow(i).getFloat("Lat");
    float lon = table.getRow(i).getFloat("Lon");
    float alt = table.getRow(i).getFloat("Altitude");
    float speed = table.getRow(i).getFloat("Speed");

    float x = map(lon, lonFrom, lonTo, 0, japanmap.width);
    float y = map(lat, latFrom, latTo, japanmap.height, 0);
    float z = map(alt, 0, 40000, 0, 60); //高度をzに変換
    float hue = map(speed, 0, 500, 192, 0); //速度を色相に変換
    stroke(hue, 255, 255);
    vertex(x, y, z);
  }

  endShape();
  
  fill(0,0,255);
  noStroke();
  //データのある位置をプロットする
  for (int i=0; i<table.getRowCount(); i++) {
    float lat = table.getRow(i).getFloat("Lat");
    float lon = table.getRow(i).getFloat("Lon");
    float alt = table.getRow(i).getFloat("Altitude");

    float x = map(lon, lonFrom, lonTo, 0, japanmap.width);
    float y = map(lat, latFrom, latTo, japanmap.height, 0);
    float z = map(alt, 0, 40000, 0, 60);
    pushMatrix();
    translate(x,y,z);
    sphere(0.5);
    popMatrix();
  }
  
```
実行結果
![sketch_flight.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/80141/c0fe7dc5-d0bf-67e3-4275-044c45e7520f.png)


## 完成

完成したら、グリグリ動かしてさまざまな角度からお楽しみください。

![sketch_flight_2.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/80141/c4b261dd-b139-48f8-d417-a05efa708ef8.png)

![sketch_flight_3.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/80141/82e0f480-e32a-6fc8-f2ea-e19e653c54cb.png)


# 発展形

よりフライト気分を味わうためにアニメーションやスライダーを使うサンプルも入っています。動画は最終形のサンプルです。（sketch_flight_final.pde）

![5f5d6266268863d589527903928e0e6a.gif](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/80141/339e7e55-5a9f-09b5-659a-2cfdd1c67f7e.gif)

サンプルコードには以下のような、さまざまな発展形を入れてあります。

* step1: データのある位置をプロットする（上で解説済み）
* step2: 通過済みのポイントは赤くする。
    * Timestampを利用する。Unixtimeとは？[Flightradarの参考記事](https://www.flightradar24.com/blog/using-the-new-flightradar24-kml-and-csv-export-tools/)
* step3: ポイント10個毎に通過時刻を表示
    * Unixtimeと時刻の変換
* step4: 飛んでいる位置に合わせて点を移動させる
    * 出発してからの時間t（秒）を増やしつつ、tから現在地を特定する
    * 特定した位置を青いプロットで表現
* step5: スライダーで操作で任意の時間の飛行位置を表示する
    * UIライブラリ「controlP5」を使ってスライダーを作る
* step6: オートプレイとスライダー操作の切り分け方
    * スペースキーでオートプレイとポーズを切り替える
* step7:現在地に飛行機のアイコン画像を表示
    * 飛行機のアイコン画像の読み込みと表示
    * 方向（Direction）を使って進行方向とアイコンを一致させる
* step8:滑らかな移動の表現
    * 位置の線形補間（lerp関数）を使う
* final: いろいろ装飾してみる。（上記動画）

いろいろと発展させて、快適な空の旅をお楽しみください。



# 終わりに

一刻も早く、自由な移動や旅行が楽しめるようになるといいですね。