Time MFD aspect calculations

time-mfd-aspect.py

# encoding=UTF-8
"""time-mfd-aspect.py

Time MFD mean aspect ratio generation routine over 100 runs.

Requires:
 - pygeoprocessing==1.9.2
 - git+https://github.com/phargogh/invest.git@8ac2575#egg=natcap.invest
"""

import os
import time

import pygeoprocessing
import pygeoprocessing.routing
from natcap.invest.sdr import sdr_core

WORKSPACE = 'test-mfd-aspect'
DEM = 'natcap/invest-sample-data/Base_Data/Freshwater/dem/hdr.adf'

if not os.path.exists(WORKSPACE):
    os.mkdir(WORKSPACE)

print('filling dem')
FILLED_DEM = os.path.join(WORKSPACE, 'filled_dem.tif')
pygeoprocessing.routing.fill_pits(
    (DEM, 1), FILLED_DEM)

print('calculating flow dir')
FLOW_DIR = os.path.join(WORKSPACE, 'flow_dir_mfd.tif')
pygeoprocessing.routing.flow_dir_mfd(
    (FILLED_DEM, 1), FLOW_DIR)

AVG_ASPECT = os.path.join(WORKSPACE, 'avg_aspect_%s.tif')

print('calculating avg aspect')
start_time = time.time()
N_REPS = 100
for i in range(N_REPS):
    if i % 10 == 0:
        print(i)
    sdr_core.calculate_average_aspect(FLOW_DIR, AVG_ASPECT % str(i))
print(round((time.time() - start_time) / N_REPS, 5))