https://github.com/minorua/Qgis2threejs
Tip revision: 744bd03195167190b0f1fc5229ad210aea432fc6 authored by Minoru Akagi on 31 October 2018, 07:34:05 UTC
version 2.2
version 2.2
Tip revision: 744bd03
builddem.py
# -*- coding: utf-8 -*-
"""
/***************************************************************************
Qgis2threejs
A QGIS plugin
export terrain data, map canvas image and vector data to web browser
-------------------
begin : 2014-01-16
copyright : (C) 2014 Minoru Akagi
email : akaginch@gmail.com
***************************************************************************/
/***************************************************************************
* *
* This program is free software; you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by *
* the Free Software Foundation; either version 2 of the License, or *
* (at your option) any later version. *
* *
***************************************************************************/
"""
import struct
from PyQt5.QtCore import QByteArray, QSize
from qgis.core import QgsPoint, QgsProject
from .conf import DEBUG_MODE
from .datamanager import MaterialManager
from .buildlayer import LayerBuilder
from .geometry import PolygonGeometry, TriangleMesh, IndexedTriangles2D, dissolvePolygonsOnCanvas
from .propertyreader import DEMPropertyReader
from .rotatedrect import RotatedRect
class DEMLayerBuilder(LayerBuilder):
def __init__(self, settings, imageManager, layer, pathRoot=None, urlRoot=None, progress=None):
"""if both pathRoot and urlRoot are None, object is built in all_in_dict mode."""
LayerBuilder.__init__(self, settings, imageManager, layer, pathRoot, urlRoot, progress)
self.provider = settings.demProviderByLayerId(layer.layerId)
self.prop = DEMPropertyReader(layer.layerId, layer.properties)
def build(self, build_blocks=False):
if self.provider is None:
return None
d = {
"type": "layer",
"id": self.layer.jsLayerId,
"properties": self.layerProperties()
}
if DEBUG_MODE:
d["PROPERTIES"] = self.properties
# DEM block
if build_blocks:
d["data"] = [block.build() for block in self.blocks()]
else:
d["data"] = []
return d
def layerProperties(self):
p = LayerBuilder.layerProperties(self)
p["type"] = "dem"
p["shading"] = self.properties.get("checkBox_Shading", True)
return p
def blocks(self):
mapTo3d = self.settings.mapTo3d()
baseExtent = self.settings.baseExtent
center = baseExtent.center()
rotation = baseExtent.rotation()
base_grid_size = self.prop.demSize(self.settings.mapSettings.outputSize())
# clipping
clip_geometry = None
clip_option = self.properties.get("checkBox_Clip", False)
if clip_option:
clip_layerId = self.properties.get("comboBox_ClipLayer")
clip_layer = QgsProject.instance().mapLayer(clip_layerId) if clip_layerId else None
if clip_layer:
clip_geometry = dissolvePolygonsOnCanvas(self.settings, clip_layer)
# surroundings
surroundings = self.properties.get("checkBox_Surroundings", False) #TODO: [GSIElevProvider] if prop.layerId else False
roughening = self.properties["spinBox_Roughening"] if surroundings else 1
size = self.properties["spinBox_Size"] if surroundings else 1
size2 = size * size
blks = []
for i in range(size2):
sx = i % size - (size - 1) // 2
sy = i // size - (size - 1) // 2
dist2 = sx * sx + sy * sy
blks.append([dist2, i, sx, sy])
for dist2, blockIndex, sx, sy in sorted(blks):
#self.progress(20 * i / size2 + 10)
is_center = (sx == 0 and sy == 0)
if is_center:
extent = baseExtent
grid_size = base_grid_size
else:
block_center = QgsPoint(center.x() + sx * baseExtent.width(), center.y() + sy * baseExtent.height())
extent = RotatedRect(block_center, baseExtent.width(), baseExtent.height()).rotate(rotation, center)
grid_size = QSize(max(2, (base_grid_size.width() - 1) // roughening + 1),
max(2, (base_grid_size.height() - 1) // roughening + 1))
block = DEMBlockBuilder(self.settings,
self.imageManager,
self.layer,
blockIndex,
self.provider,
grid_size,
extent,
mapTo3d.planeWidth,
mapTo3d.planeHeight,
offsetX=mapTo3d.planeWidth * sx,
offsetY=mapTo3d.planeHeight * sy,
edgeRougheness=roughening if is_center else 1,
clip_geometry=clip_geometry if is_center else None,
pathRoot=self.pathRoot,
urlRoot=self.urlRoot)
yield block
class DEMBlockBuilder:
def __init__(self, settings, imageManager, layer, blockIndex, provider, grid_size, extent, planeWidth, planeHeight, offsetX=0, offsetY=0, edgeRougheness=1, clip_geometry=None, pathRoot=None, urlRoot=None):
self.settings = settings
self.imageManager = imageManager
self.materialManager = MaterialManager(settings.materialType())
self.layer = layer
self.properties = layer.properties
self.blockIndex = blockIndex
self.provider = provider
self.grid_size = grid_size
self.extent = extent
self.planeWidth = planeWidth
self.planeHeight = planeHeight
self.offsetX = offsetX
self.offsetY = offsetY
self.edgeRougheness = edgeRougheness
self.clip_geometry = clip_geometry
self.pathRoot = pathRoot
self.urlRoot = urlRoot
def build(self):
if self.edgeRougheness == 1:
ba = self.provider.read(self.grid_size.width(), self.grid_size.height(), self.extent)
else:
grid_values = list(self.provider.readValues(self.grid_size.width(), self.grid_size.height(), self.extent))
self.processEdges(grid_values, self.edgeRougheness)
ba = struct.pack("{0}f".format(self.grid_size.width() * self.grid_size.height()), *grid_values)
# write grid values to an external binary file
if self.pathRoot is not None:
with open(self.pathRoot + "{0}.bin".format(self.blockIndex), "wb") as f:
f.write(ba)
# block data
g = {"width": self.grid_size.width(),
"height": self.grid_size.height()}
if self.urlRoot is None:
g["binary"] = QByteArray(ba)
# g["array"] = grid_values
else:
g["url"] = self.urlRoot + "{0}.bin".format(self.blockIndex)
# material
material = self.material()
mapTo3d = self.settings.mapTo3d()
b = {"type": "block",
"layer": self.layer.jsLayerId,
"block": self.blockIndex,
"grid": g,
"width": self.planeWidth,
"height": self.planeHeight,
"translate": [self.offsetX, self.offsetY, mapTo3d.verticalShift * mapTo3d.multiplierZ],
"zShift": mapTo3d.verticalShift,
"zScale": mapTo3d.multiplierZ,
"material": material}
# clipped with polygon layer
if self.clip_geometry:
b["clip"] = self.clipped()
# sides and bottom
if self.properties.get("checkBox_Sides", False):
b["sides"] = True
# frame
if self.properties.get("checkBox_Frame", False) and not self.properties.get("checkBox_Clip", False):
b["frame"] = True
return b
def material(self):
# properties
texture_scale = self.properties.get("comboBox_TextureSize", 100) // 100
opacity = self.properties.get("spinBox_Opacity", 100) / 100
transp_background = self.properties.get("checkBox_TransparentBackground", False)
# display type
canvas_size = self.settings.mapSettings.outputSize()
if self.properties.get("radioButton_MapCanvas", False):
#if texture_scale == 1:
# mi = self.materialManager.getCanvasImageIndex(opacity, transp_background)
#else:
mi = self.materialManager.getMapImageIndex(canvas_size.width() * texture_scale, canvas_size.height() * texture_scale, self.extent, opacity, transp_background)
elif self.properties.get("radioButton_LayerImage", False):
layerids = self.properties.get("layerImageIds", [])
mi = self.materialManager.getLayerImageIndex(layerids, canvas_size.width() * texture_scale, canvas_size.height() * texture_scale, self.extent, opacity, transp_background)
elif self.properties.get("radioButton_ImageFile", False):
filepath = self.properties.get("lineEdit_ImageFile", "")
mi = self.materialManager.getImageFileIndex(filepath, opacity, transp_background, True)
else: #.get("radioButton_SolidColor", False)
mi = self.materialManager.getMeshMaterialIndex(self.properties.get("colorButton_Color", ""), opacity, True)
#elif self.properties.get("radioButton_Wireframe", False):
# mi = self.materialManager.getWireframeIndex(self.properties["lineEdit_Color"], opacity)
# build material
filepath = None if self.pathRoot is None else "{0}{1}.png".format(self.pathRoot, self.blockIndex)
url = None if self.urlRoot is None else "{0}{1}.png".format(self.urlRoot, self.blockIndex)
return self.materialManager.build(mi, self.imageManager, filepath, url, self.settings.base64)
def clipped(self):
mapTo3d = self.settings.mapTo3d()
z_func = lambda x, y: 0
transform_func = lambda x, y, z: mapTo3d.transform(x, y, z)
# create triangle mesh
hw = 0.5 * mapTo3d.planeWidth
hh = 0.5 * mapTo3d.planeHeight
tmesh = TriangleMesh(-hw, -hh,
hw, hh,
self.grid_size.width() - 1, self.grid_size.height() - 1)
# split polygons with triangle mesh
geom = PolygonGeometry.fromQgsGeometry(self.clip_geometry, z_func, transform_func)
geom.splitPolygon(tmesh, z_func)
triangles = IndexedTriangles2D()
split_polygons = []
for polygon in geom.split_polygons:
boundary = polygon[0]
if len(polygon) == 1 and len(boundary) == 4:
triangles.addTriangle(boundary[0], boundary[2], boundary[1]) # vertex order should be counter-clockwise
else:
bnds = [[[pt.x, pt.y] for pt in bnd] for bnd in polygon]
split_polygons.append(bnds)
return {"polygons": geom.asList2(),
"triangles": {"v": [[pt.x, pt.y] for pt in triangles.vertices],
"f": triangles.faces},
"split_polygons": split_polygons}
def processEdges(self, grid_values, roughness):
grid_width = self.grid_size.width()
grid_height = self.grid_size.height()
rg_grid_width = (grid_width - 1) // roughness + 1
rg_grid_height = (grid_height - 1) // roughness + 1
ii = range(roughness)[1:]
for x0 in range(rg_grid_width - 1):
# top edge
ix0 = x0 * roughness
z0 = grid_values[ix0]
z1 = grid_values[ix0 + roughness]
for i in ii:
grid_values[ix0 + i] = (z1 - z0) * i / roughness + z0
# bottom edge
iy0 = grid_width * (grid_height - 1)
z0 = grid_values[iy0 + ix0]
z1 = grid_values[iy0 + ix0 + roughness]
for i in ii:
grid_values[iy0 + ix0 + i] = (z1 - z0) * i / roughness + z0
rw = roughness * grid_width
for y0 in range(rg_grid_height - 1):
# left edge
iy0 = y0 * grid_width
z0 = grid_values[iy0]
z1 = grid_values[iy0 + rw]
for i in ii:
grid_values[iy0 + i * grid_width] = (z1 - z0) * i / roughness + z0
# right edge
iy0 += grid_width - 1
z0 = grid_values[iy0]
z1 = grid_values[iy0 + rw]
for i in ii:
grid_values[iy0 + i * grid_width] = (z1 - z0) * i / roughness + z0
def getValue(self, x, y):
def _getValue(gx, gy):
return self.grid_values[gx + self.grid_width * gy]
if 0 <= x and x <= self.grid_width - 1 and 0 <= y and y <= self.grid_height - 1:
ix, iy = int(x), int(y)
sx, sy = x - ix, y - iy
z11 = _getValue(ix, iy)
z21 = 0 if x == self.grid_width - 1 else _getValue(ix + 1, iy)
z12 = 0 if y == self.grid_height - 1 else _getValue(ix, iy + 1)
z22 = 0 if x == self.grid_width - 1 or y == self.grid_height - 1 else _getValue(ix + 1, iy + 1)
return (1 - sx) * ((1 - sy) * z11 + sy * z12) + sx * ((1 - sy) * z21 + sy * z22) # bilinear interpolation
return 0 # as safe null value
def gridPointToPoint(self, x, y):
x = self.rect.xMinimum() + self.rect.width() / (self.grid_width - 1) * x
y = self.rect.yMaximum() - self.rect.height() / (self.grid_height - 1) * y
return x, y
def pointToGridPoint(self, x, y):
x = (x - self.rect.xMinimum()) / self.rect.width() * (self.grid_width - 1)
y = (self.rect.yMaximum() - y) / self.rect.height() * (self.grid_height - 1)
return x, y
class DEMBlocks:
def __init__(self):
self.blocks = []
def appendBlock(self, block):
self.blocks.append(block)
def appendBlocks(self, blocks):
self.blocks += blocks
def processEdges(self):
"""for now, this function is designed for simple resampling mode with surroundings"""
count = len(self.blocks)
if count < 9:
return
ci = (count - 1) // 2
size = int(count ** 0.5)
center = self.blocks[0]
blocks = self.blocks[1:ci + 1] + [center] + self.blocks[ci + 1:]
grid_width, grid_height, grid_values = center.grid_width, center.grid_height, center.grid_values
for istop, neighbor in enumerate([blocks[ci - size], blocks[ci + size]]):
if grid_width == neighbor.grid_width:
continue
y = grid_height - 1 if not istop else 0
for x in range(grid_width):
gx, gy = center.gridPointToPoint(x, y)
gx, gy = neighbor.pointToGridPoint(gx, gy)
grid_values[x + grid_width * y] = neighbor.getValue(gx, gy)
for isright, neighbor in enumerate([blocks[ci - 1], blocks[ci + 1]]):
if grid_height == neighbor.grid_height:
continue
x = grid_width - 1 if isright else 0
for y in range(grid_height):
gx, gy = center.gridPointToPoint(x, y)
gx, gy = neighbor.pointToGridPoint(gx, gy)
grid_values[x + grid_width * y] = neighbor.getValue(gx, gy)
def stats(self):
if len(self.blocks) == 0:
return {"max": 0, "min": 0}
block = self.blocks[0]
stats = {"max": block.orig_stats["max"], "min": block.orig_stats["min"]}
for block in self.blocks[1:]:
stats["max"] = max(block.orig_stats["max"], stats["max"])
stats["min"] = min(block.orig_stats["min"], stats["min"])
return stats
def dummyProgress(progress=None, statusMsg=None):
pass