Elevation Data as array

[FabianN]

Hi,

Is there a possibility to extract elevation data (DEM) within a given polygon (shape) - for example as a numpy array?

Thanks

[aldanstar]

Code: [Select]

import Metashape as ps
import os
import numpy as np
from PIL import Image
from zipfile import ZipFile
import xml.etree.ElementTree as ET

def read_elev():
    chunk = ps.app.document.chunk

    path = os.path.dirname(ps.app.document.path)
    name = os.path.basename(ps.app.document.path).split('.')[0]
    recurces_path = os.path.normpath(os.path.join(path, '{}.files'.format(name)))
    current_chunk_recurces_path = os.path.normpath(os.path.join(recurces_path, '{}'.format(chunk.key)))
    current_elevation_recurces_path = os.path.normpath(
        os.path.join(current_chunk_recurces_path, '0', 'elevation')) if chunk.elevation.key == 0 else os.path.normpath(
        os.path.join(current_chunk_recurces_path, '0', 'elevation.{}'.format(chunk.elevation.key)))

    if not os.path.exists(current_elevation_recurces_path): return

    elevation_zip_path = os.path.normpath(os.path.join(current_elevation_recurces_path, 'elevation.zip'))

    if not os.path.exists(elevation_zip_path): return

    with ZipFile(elevation_zip_path, 'r') as zipObj:
        listOfFileNames = zipObj.namelist()
        metaXml = zipObj.open(listOfFileNames[0])
        metaXmlData = metaXml.read()

    root = ET.fromstring(metaXmlData)

    tileWidth = int(root.findall("./params/dimensions/tileWidth")[0].text)
    tileHeight = int(root.findall("./params/dimensions/tileHeight")[0].text)

    width = int(root.findall("./params/dimensions/width")[0].text)
    height = int(root.findall("./params/dimensions/height")[0].text)

    tiles = []

    for child in root.findall("./tiles/tile"):
        try:
            child.attrib['level']
        except KeyError:
            tile = {}
            tile['x'] = int(child.attrib['x']) * tileWidth
            tile['y'] = int(child.attrib['y']) * tileHeight
            tile['path'] = os.path.normpath(
                os.path.join(current_elevation_recurces_path, child.findall("./path")[0].text))
            tiles.append(tile)

    arr = np.zeros((height, width))

    for tile in tiles:
        im = Image.open(tile['path'])
        imarray = np.array(im)
        x = tile['x']
        y = tile['y']
        arr[y:y + tileHeight, x:x + tileWidth] = imarray

    result = Image.fromarray(arr)

    result.show()