台つき もう少し広範囲の３次元メッシュ（xml）

　絵が見づらいです。5mの三宅島で、こちらも台つきバージョンです。

　Pythonのスクリプトは次のようになります。ここに生pyファイルがあります。

import bpy import math from xml.etree.ElementTree import ElementTree import glob # # 台つきもう少し広範囲の３次元メッシュ（xml） # demmulti.py # 2016/11/17 Shizuka-Toikanbetsu # XmlFile = "FG-GML-*.xml" tblheight = 2 # 台の部分の高さ # # 既存カーブの削除 for item in bpy.context.scene.objects: if item.type == 'MESH': bpy.context.scene.objects.unlink(item) for item in bpy.data.objects: if item.type == 'MESH': bpy.data.objects.remove(item) for item in bpy.data.meshes: bpy.data.meshes.remove(item) # # namespaceを定義 ns = { 'xml': 'http://fgd.gsi.go.jp/spec/2008/FGD_GMLSchema', 'gml': 'http://www.opengis.net/gml/3.2' } # fact = 36000.0 # # 世界座標の基準緯度と世界座標の目標座標の距離をメートル単位で返す def GetLatOffset(pBaseLat, pLat): return ( pLat - pBaseLat ) * 40007000.0 / 360.0 # # 世界座標の基準経度と世界座標の目標座標の距離をメートル単位で返す def GetLonOffset(pBaseLon, pLat, pLon): return ( pLon - pBaseLon ) * 40075000.0 * math.cos( math.radians( pLat ) ) / 360.0 # class CMinMax: def __init__(self): self.Xmin = 100000000.0 self.Ymin = self.Xmin self.Xmax = - self.Xmin self.Ymax = self.Xmax self.px = 0 self.py = 0 self.Xspan = 0 self.Yspan = 0 # # ２次メッシュ、３次メッシュ専用緯度経度補正関数 def Adjust2(Item): return round( Item * fact ) # def GetMapRange(pXml): print(pXml) NowRange = CMinMax() tree = ElementTree() tree.parse( pXml ) xpath = ".//gml:lowerCorner" lowerCorner = tree.find( xpath, ns ) tempstr = lowerCorner.text item = tempstr.split( ' ' ) NowRange.Ymin = Adjust2( float( item[ 0 ] ) ) NowRange.Xmin = Adjust2( float( item[ 1 ] ) ) xpath = ".//gml:upperCorner" upperCorner = tree.find( xpath, ns ) tempstr = upperCorner.text item = tempstr.split( ' ' ) NowRange.Ymax = Adjust2( float( item[ 0 ] ) ) NowRange.Xmax = Adjust2( float( item[ 1 ] ) ) xpath = ".//gml:high" high = tree.find( xpath, ns ) tempstr = high.text item = tempstr.split( ' ' ) NowRange.px = int( item[ 0 ] ) + 1 NowRange.py = int( item[ 1 ] ) + 1 NowRange.Xspan = ( NowRange.Xmax - NowRange.Xmin ) / NowRange.px NowRange.Yspan = ( NowRange.Ymax - NowRange.Ymin ) / NowRange.py return NowRange # print( "\nStart Script" ) # destminmax = CMinMax() # get xml files files = glob.glob( XmlFile ) if len( files ): print( 'Now get range from all xml.' ) for file in files: NowMapRange = GetMapRange( file ) destminmax.Xmin = min( NowMapRange.Xmin, destminmax.Xmin ) destminmax.Xmax = max( NowMapRange.Xmax, destminmax.Xmax ) destminmax.Ymin = min( NowMapRange.Ymin, destminmax.Ymin ) destminmax.Ymax = max( NowMapRange.Ymax, destminmax.Ymax ) destminmax.px = round( ( destminmax.Xmax - destminmax.Xmin ) / NowMapRange.Xspan ) destminmax.py = round( ( destminmax.Ymax - destminmax.Ymin ) / NowMapRange.Yspan ) destminmax.Xspan = int( NowMapRange.Xspan ) destminmax.Yspan = int( NowMapRange.Yspan ) wx = ( destminmax.Xmin + destminmax.Xmax ) / 2.0 / fact wy = ( destminmax.Ymin + destminmax.Ymax ) / 2.0 / fact BlockWidth = ( NowMapRange.Xmax - NowMapRange.Xmin ) BlockHeight = ( NowMapRange.Ymax - NowMapRange.Ymin ) BlockXspan = int( BlockWidth / NowMapRange.px ) BlockYspan = int( BlockHeight / NowMapRange.py ) # print( "add grid." ) bpy.ops.mesh.primitive_grid_add( x_subdivisions = destminmax.px, y_subdivisions = destminmax.py, location = ( 0, 0, 0 ) ) Grid = bpy.data.objects[ 'Grid' ] Grid.scale[0] = 0.001 Grid.scale[1] = 0.001 Grid.scale[2] = 0.001 bpy.context.scene.objects.active = Grid # print( "realighn grid" ) bgnY = int( destminmax.Ymin + destminmax.Yspan / 2 ) btmY = destminmax.Ymax ptr = 0 bgn = int( destminmax.Ymin + destminmax.Yspan / 2 ) for ly in range( bgnY, btmY, destminmax.Yspan ): bgnX = int( destminmax.Xmin + destminmax.Xspan / 2 ) btmX = destminmax.Xmax posY = ly / fact for lx in range( bgnX, btmX, destminmax.Xspan ): posX = lx / fact dx = GetLonOffset( wx, posY, posX ) dy = GetLatOffset( wy, posY ) bpy.context.object.data.vertices[ ptr ].co[ 0 ] = dx bpy.context.object.data.vertices[ ptr ].co[ 1 ] = dy ptr = ptr + 1 # bpy.ops.object.modifier_add(type='SOLIDIFY') bpy.context.object.modifiers["Solidify"].thickness = tblheight * 1000 bpy.context.object.modifiers["Solidify"].thickness_clamp = 0 bpy.ops.object.modifier_apply(apply_as='DATA', modifier="Solidify") # print( "\n\nDecording XML." ) for file in files: print( file ) tree = ElementTree() tree.parse( file ) # NowRange = CMinMax() xpath = ".//gml:lowerCorner" lowerCorner = tree.find( xpath, ns ) tempstr = lowerCorner.text item = tempstr.split( ' ' ) NowRange.Ymin = Adjust2( float( item[ 0 ] ) ) NowRange.Xmin = Adjust2( float( item[ 1 ] ) ) xpath = ".//gml:upperCorner" upperCorner = tree.find( xpath, ns ) tempstr = upperCorner.text item = tempstr.split( ' ' ) NowRange.Ymax = Adjust2( float( item[ 0 ] ) ) NowRange.Xmax = Adjust2( float( item[ 1 ] ) ) xpath = ".//gml:tupleList" temp = tree.find( xpath, ns ) tupleListStr = temp.text xpath = ".//gml:startPoint" startPoint = tree.find( xpath, ns ) tempstr = startPoint.text item = tempstr.split( ' ' ) sx = int( item[ 0 ] ) sy = int( item[ 1 ] ) # # 標高データの文字列を分解し、実数型の準備データを作成 AltSrc = [] AltList = tupleListStr.split( '\n' ) for AltLine in AltList: if AltLine: # 空行でなければ以下を実行 LineData = AltLine.split( ',' ) if LineData: AltSrc.append( float( LineData[ 1 ] ) ) # # 準備データを本番データに転記 Alt = [] for destptr in range( 0, NowMapRange.px * NowMapRange.py ): Alt.append( -9999.0 ) srcptr = 0 destptr = sy * NowMapRange.px + sx for ptr in range( 0, len( AltSrc ) ): Alt[ destptr ] = AltSrc[ srcptr ] destptr = destptr + 1 srcptr = srcptr + 1 Xblock = int( ( NowRange.Xmin - destminmax.Xmin ) / BlockWidth ) Yblock = int( ( NowRange.Ymin - destminmax.Ymin ) / BlockHeight ) AltPtr = 0 for Yptr in range( 0, NowMapRange.py ): Ybase = ( NowMapRange.py - 1 ) - Yptr Ybase = Yblock * NowMapRange.py + Ybase for Xptr in range( 0, NowMapRange.px ): Xbase = Xblock * NowMapRange.px + Xptr destindex = Ybase * destminmax.px + Xbase vert = bpy.context.object.data.vertices[ destindex ] temp = Alt[ AltPtr ] if temp == -9999.0: temp = 0 vert.co[ 2 ] = temp temp = destindex + destminmax.px * destminmax.py AltPtr = AltPtr + 1 Grid.location[2] = tblheight # print( 'Exit' )

<トップページに戻る>