Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

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

gradientm is a matlab function that Calculate gradient, slope, and aspect of data grid

Syntax

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

Description

[ASPECT, SLOPE, gradN, gradE] = gradientm(Z, R) computes the slope, aspect, and north and east components of the gradient for a regular data grid Z with respect to reference R. If the grid contains elevations in meters, the resulting aspect and slope are in units of degrees clockwise from north and up from the horizontal. The north and east gradient components are the change in the map variable per meter of distance in the north and east directions. The computation uses finite differences for the map variable on the default earth ellipsoid.

Problem

I want something equivalent in python. I have founded this python library PyDEM but it is for python 2.x but I am working with python 3.x

this is The way of Use the DEMProcessor to compute slope and aspect for the DEM files by 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()
like image 957
pd shah Avatar asked Dec 06 '22 12:12

pd shah


1 Answers

step1: install GDAL

step2: install rasterio

step3: easy as 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)
like image 136
pd shah Avatar answered May 16 '23 07:05

pd shah