Revision 6a557d2bd737748d622d1ae12c6537f5a2139543 authored by Minoru Akagi on 30 June 2016, 00:50:24 UTC, committed by Minoru Akagi on 30 June 2016, 00:50:24 UTC
1 parent 2c38ad6
demblock.py
# -*- coding: utf-8 -*-
"""
/***************************************************************************
DEMBlock
-------------------
begin : 2015-05-23
copyright : (C) 2015 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. *
* *
***************************************************************************/
"""
from qgis.core import QgsRectangle
import gdal2threejs
from geometry import PolygonGeometry, Triangles
from qgis2threejstools import pyobj2js
class DEMBlock:
def __init__(self, dem_width, dem_height, dem_values, plane_width, plane_height, offsetX, offsetY):
self.dem_width = dem_width
self.dem_height = dem_height
self.dem_values = dem_values
self.plane_width = plane_width
self.plane_height = plane_height
self.offsetX = offsetX
self.offsetY = offsetY
self.orig_stats = {"max": max(dem_values), "min": min(dem_values)}
self.rect = QgsRectangle(offsetX - plane_width * 0.5, offsetY - plane_height * 0.5,
offsetX + plane_width * 0.5, offsetY + plane_height * 0.5)
self.properties = {"width": dem_width, "height": dem_height}
self.properties["plane"] = {"width": plane_width, "height": plane_height,
"offsetX": offsetX, "offsetY": offsetY}
self.clip_geometry = None
def set(self, key, value):
"""set property"""
self.properties[key] = value
def setClipGeometry(self, geometry):
self.clip_geometry = geometry
def zShift(self, shift):
if shift != 0:
self.dem_values = map(lambda x: x + shift, self.dem_values)
def zScale(self, scale):
if scale != 1:
self.dem_values = map(lambda x: x * scale, self.dem_values)
def write(self, writer):
mapTo3d = writer.settings.mapTo3d()
writer.write("bl = lyr.addBlock({0}, {1});\n".format(pyobj2js(self.properties), pyobj2js(bool(self.clip_geometry))))
writer.write("bl.data = [{0}];\n".format(",".join(map(gdal2threejs.formatValue, self.dem_values))))
# clipped with polygon layer
if self.clip_geometry:
z_func = lambda x, y: 0
transform_func = lambda x, y, z: mapTo3d.transform(x, y, z)
geom = PolygonGeometry.fromQgsGeometry(self.clip_geometry, z_func, transform_func)
geom.splitPolygon(writer.triangleMesh(self.dem_width, self.dem_height))
polygons = []
for polygon in geom.polygons:
bnds = []
for boundary in polygon:
bnds.append(map(lambda pt: [pt.x, pt.y], boundary))
polygons.append(bnds)
writer.write("bl.clip = {};\n")
writer.write("bl.clip.polygons = {0};\n".format(pyobj2js(polygons)))
triangles = Triangles()
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 = [map(lambda pt: [pt.x, pt.y], bnd) for bnd in polygon]
polygons.append(bnds)
vf = {"v": map(lambda pt: [pt.x, pt.y], triangles.vertices), "f": triangles.faces}
writer.write("bl.clip.triangles = {0};\n".format(pyobj2js(vf)))
writer.write("bl.clip.split_polygons = {0};\n".format(pyobj2js(polygons)))
def getValue(self, x, y):
def _getValue(gx, gy):
return self.dem_values[gx + self.dem_width * gy]
if 0 <= x and x <= self.dem_width - 1 and 0 <= y and y <= self.dem_height - 1:
ix, iy = int(x), int(y)
sx, sy = x - ix, y - iy
z11 = _getValue(ix, iy)
z21 = 0 if x == self.dem_width - 1 else _getValue(ix + 1, iy)
z12 = 0 if y == self.dem_height - 1 else _getValue(ix, iy + 1)
z22 = 0 if x == self.dem_width - 1 or y == self.dem_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.dem_width - 1) * x
y = self.rect.yMaximum() - self.rect.height() / (self.dem_height - 1) * y
return x, y
def pointToGridPoint(self, x, y):
x = (x - self.rect.xMinimum()) / self.rect.width() * (self.dem_width - 1)
y = (self.rect.yMaximum() - y) / self.rect.height() * (self.dem_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:]
dem_width, dem_height, dem_values = center.dem_width, center.dem_height, center.dem_values
for istop, neighbor in enumerate([blocks[ci - size], blocks[ci + size]]):
if dem_width == neighbor.dem_width:
continue
y = dem_height - 1 if not istop else 0
for x in range(dem_width):
gx, gy = center.gridPointToPoint(x, y)
gx, gy = neighbor.pointToGridPoint(gx, gy)
dem_values[x + dem_width * y] = neighbor.getValue(gx, gy)
for isright, neighbor in enumerate([blocks[ci - 1], blocks[ci + 1]]):
if dem_height == neighbor.dem_height:
continue
x = dem_width - 1 if isright else 0
for y in range(dem_height):
gx, gy = center.gridPointToPoint(x, y)
gx, gy = neighbor.pointToGridPoint(gx, gy)
dem_values[x + dem_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 write(self, writer, separated=False):
for block in self.blocks:
if separated:
writer.nextFile(True)
block.write(writer)
![swh spinner](/static/img/swh-spinner.gif)
Computing file changes ...