Commit abab4a51 authored by luroth's avatar luroth
Browse files

Strickhof processing

parent c2821068
from Common import GISfunctions
import os
import shutil
import numpy as np
from scipy import ndimage
import gdal
import cv2
from pathlib import Path
def process_df(path_dtm, path_dsm, path_output):
GISfunctions.create_elevation_diff_GeoTiff(path_dtm, path_dsm, path_output)
if __name__ == "__main__":
path_dtm = r'O:\UAV\_Processed_\DSP_strickhof_DSWW003\20181109\28m_M600P\orthoJPEG_DEM.tif'
path_dtm_smoothed = r'O:\UAV\_Processed_\DSP_strickhof_DSWW005\20181109\28m_M600P\orthoJPEG_DEM_mv5median.tif'
if not Path(path_dtm_smoothed).exists():
# Read DTM
print("Reading DTM")
dtm_ds = gdal.Open(path_dtm)
dtm_band1 = dtm_ds.GetRasterBand(1)
dtm_arr = dtm_band1.ReadAsArray()
dtm_transform = dtm_ds.GetGeoTransform()
print("Smoothing DTM")
dtm_smoothed = cv2.medianBlur(dtm_arr, 5)
print("Saving smoothed DTM")
driver = gdal.GetDriverByName('GTiff')
ds_smoothed = driver.Create(path_dtm_smoothed, xsize=dtm_arr.shape[1], ysize=dtm_arr.shape[0], bands=1,
eType=gdal.GDT_Float32)
band1 = ds_smoothed.GetRasterBand(1)
band1.WriteArray(dtm_smoothed)
ds_smoothed.SetGeoTransform(dtm_transform)
ds_smoothed.SetProjection(dtm_ds.GetProjection())
ds_smoothed = None
path_dtm = path_dtm_smoothed
paths = [
r'O:\UAV\_Processed_\DSP_strickhof_DSWW003\20190215\28m_M600P',
r'O:\UAV\_Processed_\DSP_strickhof_DSWW003\20190222\28m_M600P',
r'O:\UAV\_Processed_\DSP_strickhof_DSWW003\20190322\28m_M600P',
r'O:\UAV\_Processed_\DSP_strickhof_DSWW003\20190411\28m_M600P',
r'O:\UAV\_Processed_\DSP_strickhof_DSWW003\20190502\28m_M600P',
r'O:\UAV\_Processed_\DSP_strickhof_DSWW003\20190509\28m_M600P',
r'O:\UAV\_Processed_\DSP_strickhof_DSWW003\20190517\28m_M600P',
r'O:\UAV\_Processed_\DSP_strickhof_DSWW003\20190523\28m_M600P',
r'O:\UAV\_Processed_\DSP_strickhof_DSWW003\20190527\28m_M600P',
r'O:\UAV\_Processed_\DSP_strickhof_DSWW003\20190530\28m_M600P',
r'O:\UAV\_Processed_\DSP_strickhof_DSWW003\20190601\28m_M600P',
r'O:\UAV\_Processed_\DSP_strickhof_DSWW003\20190605\28m_M600P',
r'O:\UAV\_Processed_\DSP_strickhof_DSWW003\20190607\28m_M600P',
r'O:\UAV\_Processed_\DSP_strickhof_DSWW003\20190612\28m_M600P',
r'O:\UAV\_Processed_\DSP_strickhof_DSWW003\20190614\28m_M600P',
r'O:\UAV\_Processed_\DSP_strickhof_DSWW003\20190617\28m_M600P',
r'O:\UAV\_Processed_\DSP_strickhof_DSWW003\20190624\28m_M600P',
r'O:\UAV\_Processed_\DSP_strickhof_DSWW003\20190628\28m_M600P',
r'O:\UAV\_Processed_\DSP_strickhof_DSWW003\20190703\28m_M600P',
r'O:\UAV\_Processed_\DSP_strickhof_DSWW003\20190705\28m_M600P',
r'O:\UAV\_Processed_\DSP_strickhof_DSWW003\20190708\28m_M600P',
]
for path in paths:
print("Process", path)
path_dsm = os.path.join(path, 'orthoJPEG_DEM.tif')
path_output = os.path.join(path, 'PH_DEM.tif')
process_df(path_dtm, path_dsm, path_output)
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment