オープンソースの点群ライブラリ、PDALの利用方法について紹介します。PDALは点群データを扱うのに様々な機能を提供しており、複雑なプログラミングを組まなくても利用が可能です。また、Python、JavaやScalaなどから呼び出し、点群データを簡単に扱うことができます。
点群データについて
点群（ポイントクラウド）データは、X、Y、Zの基本的な位置情報や色、分類ラベル（植生や建築物など）などの情報を持つデータで、1行に各列にそれらの情報がそれぞれ取りまとめられています。データフォーマットはテキスト（txtやCSV）だったり、LAS、LAZなどがあります。
点群データの主な取得方法は、LiDARによる計測が代表的です。LiDARのセンサーを地上、ドローン、航空機、衛星などに搭載し、地表物を３次元に計測します。他にも写真測量による３次元モデルなど点群データのソースは様々です。
点群データは、ドローンや衛星から取得する画像よりも、情報量や容量が大きくなる傾向にあります。自前プログラミングで処理ソフトを作ることも可能ですが、高度な数学・プログラミング知識が求められます。そのせいか、点群データを処理するソフトは高価になる場合があります。PDALは、それらの課題を解決する素晴らしいライブラリです。
PDALのインストール
PDALをインストールするには、Anacondaを利用します。Windows OSを利用する場合、本家からインストーラをダウンロードし、インストールしてください。LinuxやMacなどUNIX OSを利用の場合、以下のコマンドでAnacondaをインストールしてください。
# pyenvのインストールと設定
git clone git://github.com/yyuu/pyenv.git $HOME/.pyenv
git clone https://github.com/yyuu/pyenv-virtualenv.git $HOME/.pyenv/plugins/pyenv-virtualenv
export PYENV_ROOT=$HOME/.pyenv
export PATH=$PYENV_ROOT/shims:$PYENV_ROOT/bin:$PATH
# Bashを利用の場合、~./bashrc
echo 'eval "$(pyenv init -)"' >> ~/.zshrc
# Anacondaをインストール。バージョンは任意に変更
pyenv install anaconda3-2019.10
conda install -c conda-forge pdal python-pdal gdalWindowsの場合、Anacondaを起動し、pdal、python-pdalとgdalをインストールするようにしてください。
コマンドによるPDALの利用
PDALを利用す流にはいくつかの方法があります。一つ目がコマンドによるものです。PDALには様々な機能がすでに実装されており、自ら処理内容をプログラミングする必要がありません。例えば、フィルタリングによるノイズ除去、点群の分類、ファイルサイズの圧縮、分割、統合、投影変換などです。機能の一覧はコチラで確認できます。
点群データを処理する際、複数のステップを行う必要があります。PDALは、自ら定義した一連の処理を順番通りに行うpipelineという機能があります。処理の流れはJSONフォーマットで定義します。
今回は以下のようなJSONファイルを用意しました。
{
"pipeline":[
"input.las",
{
"type": "filters.outlier",
"method": "statistical",
"multiplier": 2.2,
"mean_k": 20
},
{
"type": "filters.range",
"limits": "Classification![7:7],Z[-100:3000]"
},
{
"type": "writers.las",
"compression": "true",
"dataformat_id": "0",
"filename":"output.las"
}
]
}
上記の内容のJSONファイルを用意したら、以下のコマンドで実行することができます。処理時間はファイルサイズやマシンのスペックに依存します。
pdal pipeline process.jsonまず、”filters.outlier”では、点群データに含まれる異常値を除去します。次に”filters.range”では、z軸の値の範囲を設定し、それ以外を除去します。最後の”writers.las”は、処理結果を出力します。”compression”がtrueとなっており、圧縮され、ファイルサイズが小さくなります。
以下がPDALの処理前になります。地表面より上の上空にいくつかのノイズのような点があることがわかります。
以下が処理後のデータです。上空にあったノイズがなくなっています。
PythonによるPDALの利用
２つ目のPDALの利用方法はプログラミングで呼び出し、実行することです。PDALは、PythonやScalaから呼び出し、利用することが可能です。プログラミング言語で利用する場合、主に２つの利用方法が用意されています。
１つ目は、上記のコマンドラインの様に、パイプライン処理を実行する方法です。２つ目は、点群データを配列データのようにメモリ上で点群処理をする方法です。何かしらの高度な処理を行いたい場合、２つ目の方法で利用するのが良いでしょう。
上記のパイプラインと同じように実行する場合、以下のようなPythonコードになります。
import pdal
from glob import glob
from os.path import join as path_join
from os.path import basename
pipeline_json = """{{
"pipeline": [
{{
"filename":"{input_las}",
"type":"readers.las",
}},
{{
"type":"filters.outliers",
"method":"statistical",
"multiplier":2.2,
"mean_k": 20
}},
{{
"type": "filters.range",
"limits": "Classification![7:7],Z[-100:3000]"
}}
{{
"type":"writers.las",
"compression":"true",
"filename":"{output_las}"
}}
]
}}"""
def pdal_filter(input_filepath: str, output_folder: str):
output_filepath = input_filepath.replace('.las', '_reprojected.las')
output_filepath = path_join(output_folder, basename(output_filepath))
# noinspection PyStringFormat
modified_json = pipeline_json.format(input_las=input_filepath, output_las=output_filepath)
pipeline = pdal.Pipeline(modified_json)
count = pipeline.execute()
return output_filepath
また、点群データをPythonで配列処理したい場合、以下のコードでデータを取り込むことが可能です。
json = """
[
"filename.las",
{
"type": "filters.sort",
"dimension": "X"
}
]
"""
import pdal
pipeline = pdal.Pipeline(json)
count = pipeline.execute()
arrays = pipeline.arrays
pipeline.execute()は配列読み込み前に実行が必要です。配列は以下のようなlistにnumpy配列が格納されています。
[array([(364431.73710739, 3749236.48392952, 971.46795866, 230, 1, 2, 1, 0, 2, 0., 0, 0, 447956.75597013),
(364431.73710739, 3749236.48392952, 971.46795866, 230, 2, 2, 1, 0, 2, 0., 0, 0, 447956.75597013),
(364431.76510739, 3749236.47092952, 971.46395866, 230, 1, 2, 1, 0, 2, 0., 0, 0, 447955.15434912),
...,
(364480.05810739, 3749195.91192952, 1007.22895866, 172, 1, 1, 1, 0, 2, 0., 0, 0, 447960.52295206),
(364480.05910739, 3749195.88392952, 1007.20495866, 172, 1, 1, 1, 0, 2, 0., 0, 0, 447950.61333706),
(364480.05910739, 3749195.88592952, 1007.20495866, 172, 1, 1, 1, 0, 2, 0., 0, 0, 447951.11374606)],
dtype=[('X', '<f8'), ('Y', '<f8'), ('Z', '<f8'), ('Intensity', '<u2'), ('ReturnNumber', 'u1'), ('NumberOfReturns', 'u1'), ('ScanDirectionFlag', 'u1'), ('EdgeOfFlightLine', 'u1'), ('Classification', 'u1'), ('ScanAngleRank', '<f4'), ('UserData', 'u1'), ('PointSourceId', '<u2'), ('GpsTime', '<f8')])]
ScalaによるPDALの利用
まず、ScalaでPDALを利用する場合、build.sbtに以下を記述してPDALのパッケージをインストールします。
libraryDependencies ++= Seq(
"com.typesafe.akka" %% "akka-actor" % "2.6.14",
// "org.scala-lang.modules" %% "scala-swing" % "3.0.0",
"org.scalafx" %% "scalafx-extras" % "0.3.6",
// "org.openjfx" % "javafx" % "12.0.2" pomOnly(),
"io.pdal" %% "pdal" % "2.2.0",
"io.pdal" %% "pdal-scala" % "2.2.0",
"io.pdal" % "pdal-native" % "2.2.0"
)ScalaでPythonのようにJSON形式でpipelineを組んで処理を実行することができます。しかし、Scalaでは、もっとプログラミングっぽく（？）記述して、実行することができます。その例が以下のコードです。
def process(input_filepath: String, output_filepath: String): Unit = {
val exp =
ReadLas(input_filepath, spatialreference = Some(inSrs)) ~
FilterOutlier(method = Some("statistical"), multiplier = Some(2.2), meanK = Some(20)) ~
FilterRange(limits = Some("Classification![7:7],Z[-100:3000]")) ~
WriteLas(filename = output_filepath, compression = Some("true"))
val pipe = exp.toPipeline
pipe.validate()
pipe.setLogLevel(8)
pipe.execute()
}
Scalaでも点群データの配列を取り込むことができます。Scalaで点群データを処理する理由は、Sparkなどの分散処理だと思います。点群データは、データポイント数が膨大のため、高速に処理する機能が重要になります。分散処理を適用することで、巨大な点群データでも短時間で処理を完了させることが期待できます。
まとめ
以上がPDALになります。点群データを扱うのに、高価なソフトウェアが必要かと思います。しかし、PDALを使えば、多くの処理を自動で行うことができます。昨日の詳細はPDALのウェブサイトに紹介されています。
Leave a Reply