在 python3.x 中计算 ASPECT、SLOPE(matlab gradientm 函数)

calculating ASPECT, SLOPE in python3.x (matlab gradientm function)

gradientm是一个计算梯度的matlab函数,slopeaspect数据网格

语法

[ASPECT, SLOPE, gradN, gradE] = gradientm(Z, R)

描述

[ASPECT, SLOPE, gradN, gradE] = gradientm(Z, R) 计算常规数据网格 Z 相对于参考 R 的梯度的斜率、纵横比以及北和东分量。如果网格包含以米为单位的高程,由此产生的坡向和坡度以从北顺时针方向和从水平方向向上的度数为单位。北部和东部梯度分量是地图变量在北部和东部方向每米距离的变化。计算对默认地球椭球体上的地图变量使用有限差分。

问题

我想要一些东西equivalent in python。 我已经创建了这个 python 库 PyDEM 但它是用于 python 2.x 但我正在使用 python 3.x

这是使用 DEMProcessorDEM 文件计算 slopeaspect 的方法 PyDEM:

# needs to match above command
filename = 'Shasta-30m-DEM.tif'

# instantiate a processor object
processor = DEMProcessor(filename)

# get magnitude of slope and aspect
mag, aspect = processor.calc_slopes_directions()

第一步:安装

第二步:安装rasterio

第 3 步:像 py 一样简单:

from osgeo import gdal
import numpy as np
import rasterio


def calculate_slope(DEM):
    gdal.DEMProcessing('slope.tif', DEM, 'slope')
    with rasterio.open('slope.tif') as dataset:
        slope=dataset.read(1)
    return slope

def calculate_aspect(DEM):
    gdal.DEMProcessing('aspect.tif', DEM, 'aspect')
    with rasterio.open('aspect.tif') as dataset:
        aspect=dataset.read(1)
    return aspect

slope=calculate_slope('DEM.tif')
aspect=calculate_aspect('DEM.tif')

print(type(slope))
print(slope.dtype)
print(slope.shape)