Skip to content
Snippets Groups Projects
Commit ca1984c1 authored by Tor-Einar Skog's avatar Tor-Einar Skog
Browse files

Pipe gdal stdout to /dev/null. Make mask optional

parent 1bd23968
No related branches found
No related tags found
No related merge requests found
......@@ -48,10 +48,14 @@ subprocess.run(f"cdo -mul {tmp_path}TM_minus_5.nc {tmp_path}TM_minus_5_gtc.nc {t
subprocess.run(f"cdo -timcumsum {tmp_path}TM_minus_5_nozero.nc {tmp_path}DD_unmasked.nc", shell=True)
# Mask the output with Norway land borders
# TODO make optional
subprocess.run(f"cdo -maskregion,Norge_landomrader.csv {tmp_path}DD_unmasked.nc {tmp_path}DD.nc", shell=True)
# Mask results using a CSV file with polygons
# Env variable MASK_FILE must be set
if os.getenv("MASK_FILE") is not None:
mask_file = os.getenv("MASK_FILE")
print(f"Applying mask file {mask_file} to result.nc")
subprocess.run(f'cdo -maskregion,{mask_file} {tmp_path}DD_unmasked.nc {tmp_path}DD.nc', shell=True)
else:
os.rename(f"{tmp_path}DD_unmasked.nc", f"{tmp_path}DD.nc")
# Add the DD threshold classification => warning status
# (A>0)*(A<=260)*2 + (A>260)*(A<=360)*3 + (A>360)*(A<=560)*4 + (A>560)*0
......@@ -68,23 +72,24 @@ timestep_dates = []
for timestep in timesteps:
timestep_date = datetime.fromtimestamp(timestep).astimezone(local_timezone)
file_date = timestep_date.astimezone(local_timezone).strftime("%Y-%m-%d")
print(f"{timestep_index}: {file_date}")
#print(f"{timestep_index}: {file_date}")
timestep_dates.append(file_date)
# We must delete the GeoTIFF file before merging
subprocess.run(f'rm {tmp_path}result_*{file_date}_lcc.tif', shell=True)
subprocess.run(f'rm {out_path}result_*{file_date}.tif', shell=True)
subprocess.run(f'gdal_translate -ot Int16 -of GTiff -b {timestep_index} NETCDF:"{tmp_path}result.nc":WARNING_STATUS {tmp_path}result_WARNING_STATUS_{file_date}_lcc.tif', shell=True)
subprocess.run(f'gdal_translate -ot Float32 -of GTiff -b {timestep_index} NETCDF:"{tmp_path}result.nc":TM {tmp_path}result_{file_date}_lcc.tif', shell=True)
# Need to reproject the files, to ensure we have the projection given in the generted mapfile. We always use EPSG:4326 for this
subprocess.run(f'gdalwarp -t_srs EPSG:4326 {tmp_path}result_WARNING_STATUS_{file_date}_lcc.tif {out_path}result_WARNING_STATUS_{file_date}.tif', shell=True)
subprocess.run(f'gdalwarp -t_srs EPSG:4326 {tmp_path}result_{file_date}_lcc.tif {out_path}result_{file_date}.tif', shell=True)
with open("/dev/null", "w") as devnull:
subprocess.run(f'gdal_translate -ot Int16 -of GTiff -b {timestep_index} NETCDF:"{tmp_path}result.nc":WARNING_STATUS {tmp_path}result_WARNING_STATUS_{file_date}_lcc.tif', shell=True, stdout=devnull)
subprocess.run(f'gdal_translate -ot Float32 -of GTiff -b {timestep_index} NETCDF:"{tmp_path}result.nc":TM {tmp_path}result_{file_date}_lcc.tif', shell=True, stdout=devnull)
# Need to reproject the files, to ensure we have the projection given in the generted mapfile. We always use EPSG:4326 for this
subprocess.run(f'gdalwarp -t_srs EPSG:4326 {tmp_path}result_WARNING_STATUS_{file_date}_lcc.tif {out_path}result_WARNING_STATUS_{file_date}.tif', shell=True, stdout=devnull)
subprocess.run(f'gdalwarp -t_srs EPSG:4326 {tmp_path}result_{file_date}_lcc.tif {out_path}result_{file_date}.tif', shell=True, stdout=devnull)
timestep_index = timestep_index + 1
if len(timestep_dates) != len(set(timestep_dates)):
print("ERROR: Something is wrong with your timesteps")
#print(timestep_dates)
print(len(timestep_dates))
#print(len(timestep_dates))
# Generate mapfile
# The paths should be set in a .env file
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment