SRTMデータの可視化

SRTM(Shuttle Radar Topography Mission)は、シャトルに搭載される合成開口レーダーを用いたリモートセンシング技術により取得された地表のレーダー画像です。

可視化ソフトウェアParaViewを使用して、このSRTMデータを可視化してみましょう。

(注:簡単な変換で可視化するため、変換後の座標などは不正確です)

(注:データとプログラム(Windows用)を最下段のhgtToCSV.zipに置きました。run.batを修正して、実行すると動作します)

SRTMデータの取得

SRTMデータは、Shuttle Radar Toplography Missionからダウンロードできます。

今回は、この中の https://dds.cr.usgs.gov/srtm/version2_1/SRTM3/Africa/ から、適当に選んだN32E020.hgt.zipを使用します。

このファイルをダウンロード、解凍すると

N32E020.hgt

という、HGT形式のファイルが格納されています。

SRTM(HGT形式)からCSVファイルへの変換

SRTMのHGT形式は、ビッグエンディアン(Big endian)で記入された1201x1201個の座標点からなるバイナリファイルです。

以下の簡単なPythonプログラムでCSVファイルに変換します。

hgtToCSV.py

Everything is expanded.Everything is shortened.
  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
import sys
 
nrows = ncols = 1201
xdim = ydim = 0.000833
 
args = sys.argv
filename = args[1]
name, ext = os.path.splitext(filename)
outname = name + '.csv'
 
longitude = float(name[1:3])
latitude = float(name[4:7])
 
print('longitude: ' + name[0] + name[1:3])
print('latitude: ' + name[3] + name[4:7])
 
ulxmap = longitude + (xdim/2.0)
ulymap = (latitude+1.0) - (ydim/2.0)
 
a = array.array('h')
f = open(filename, 'rb')
 
a.fromfile(f, nrows*ncols)
f.close()
 
a.byteswap()
alist = a.tolist()
 
out = open(outname, 'w')
 
count = 0
for row in range(nrows):
  for col in range(ncols):
    x = ulxmap + xdim * col
    y = ulymap - ydim * row
    z = max(0.0, alist[count])
    out.write(str(x) + ',' + str(y) + ',' + str(z/10000.0) + '\n')
    count += 1
 
    if count%100000 == 0:
        print(str(int(100.0 * count/(nrows*ncols))) + '%')
 
print('100%')
print('node number: ' + str(count))
 
out.close()

以下のコマンド

 > python hgtToCSV.py N32E020.hgt

を実行すると、N32E020.csvに変換されます。

注意:座標は不正確です。実際に使用するときは、正しい座標に変換してください。

CSVファイルからVTSファイルへの変換

次に、CSVファイルをParaViewで読み込めるVTS形式(VTKの構造格子用データ形式)に変換します。

こちらも簡単なPythonスクリプトを使用します。なお、標高データは変数zに格納しています。

csvToVts.py

Everything is expanded.Everything is shortened.
  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 
 
-
!
-
!
 
 
 
 
 
 
 
-
!
 
 
 
 
 
 
-
!
 
 
 
 
 
-
!
 
 
 
 
 
 
 
 
 
-
!
 
 
 
 
 
import os.path
import sys
#### import the simple module from the paraview
from paraview.simple import *
#### disable automatic camera reset on 'Show'
paraview.simple._DisableFirstRenderCameraReset()
 
args = sys.argv
filename = args[1]
 
workingDir = './'
fileName = args[1]
 
# create a new 'CSV Reader'
csvReader = CSVReader(FileName=[workingDir + fileName + '.csv'])
csvReader.DetectNumericColumns = 1
csvReader.UseStringDelimiter = 1
csvReader.HaveHeaders = 0
csvReader.FieldDelimiterCharacters = ','
csvReader.MergeConsecutiveDelimiters = 0
 
# create a new 'Table To Structured Grid'
tableToStructuredGrid = TableToStructuredGrid(Input=csvReader)
tableToStructuredGrid.WholeExtent = [0, 1200, 0, 1200, 0, 0]
tableToStructuredGrid.XColumn = 'Field 0'
tableToStructuredGrid.YColumn = 'Field 1'
tableToStructuredGrid.ZColumn = 'Field 2'
 
# create a new 'Calculator'
calculator = Calculator(Input=tableToStructuredGrid)
calculator.AttributeMode = 'Point Data'
calculator.CoordinateResults = 0
calculator.ResultNormals = 0
calculator.ResultTCoords = 0
calculator.ResultArrayName = 'z'
calculator.Function = 'coordsZ'
calculator.ReplaceInvalidResults = 1
calculator.ReplacementValue = 0.0
 
# save data
SaveData(workingDir + fileName + '.vts', proxy=calculator, 
           Writealltimestepsasfileseries=0, 
           DataMode='Appended', 
           EncodeAppendedData=1, 
           CompressorType='ZLib')

実行は、ParaViewのインストールフォルダ

  • C:/Program Files/ParaView 5.2.0-Qt4-OpenGL2-Windows-64bit/bin/pvpython.exe

にあるpvpythonを使用して実行します。

 > pvpython csvToVts.py N32E020

変換が終了すると、N32E020.vtsが出力されます。

ParaViewによる可視化

上記で作成されたN32E020.vtsファイルをParaViewに読み込むと、以下のように可視化できます。

N32E020

ダウンロード

hgtToCSVプログラム一式:


添付ファイル: filehgtToCSV.zip 218件 [詳細] fileN32E020.png 258件 [詳細]

トップ   編集 凍結 差分 バックアップ 添付 複製 名前変更 リロード   新規 一覧 単語検索 最終更新   ヘルプ   最終更新のRSS
Last-modified: 2017-02-17 (金) 08:15:48 (461d)