Skip to content
Snippets Groups Projects
Commit 2010fd02 authored by Lene Wasskog's avatar Lene Wasskog
Browse files

Merge branch 'main' into release

parents 84f5996d b3d6aded
No related branches found
No related tags found
No related merge requests found
...@@ -16,7 +16,8 @@ ...@@ -16,7 +16,8 @@
along with this program. If not, see <https://www.gnu.org/licenses/>. along with this program. If not, see <https://www.gnu.org/licenses/>.
""" """
import subprocess, os import subprocess
import os
from dotenv import load_dotenv from dotenv import load_dotenv
from datetime import datetime from datetime import datetime
import pytz import pytz
...@@ -30,6 +31,7 @@ load_dotenv() ...@@ -30,6 +31,7 @@ load_dotenv()
config = configparser.ConfigParser() config = configparser.ConfigParser()
config.read("PSILARTEMP.cfg") config.read("PSILARTEMP.cfg")
year = os.getenv("YEAR")
weatherdata_path = os.getenv("WEATHER_DATA_DIR") weatherdata_path = os.getenv("WEATHER_DATA_DIR")
tmp_path = "tmp/" tmp_path = "tmp/"
out_path = os.getenv("DATA_DIR") out_path = os.getenv("DATA_DIR")
...@@ -37,19 +39,32 @@ out_path = os.getenv("DATA_DIR") ...@@ -37,19 +39,32 @@ out_path = os.getenv("DATA_DIR")
# TODO: Put timezone out in .env # TODO: Put timezone out in .env
local_timezone = pytz.timezone("Europe/Oslo") local_timezone = pytz.timezone("Europe/Oslo")
#""" # """
# Calculate cumulated degree days above 5 degrees after 1st of April # Calculate cumulated degree days above 5 degrees after 1st of April
# Remove all values before April 1st # Remove all values before April 1st
subprocess.run(f"cdo -selname,TM -seldate,2023-04-01T00:00:00,2023-12-31T00:00:00 {weatherdata_path}met_1_0km_nordic-2023.nc {tmp_path}TM_from_april.nc", shell=True) subprocess.run(
f"cdo -selname,TM -seldate,{year}-04-01T00:00:00,{year}-12-31T00:00:00 {weatherdata_path}met_1_0km_nordic-{year}.nc {tmp_path}TM_from_april.nc",
shell=True,
)
# Subtracting 5 deg C from all cells # Subtracting 5 deg C from all cells
subprocess.run(f"cdo -subc,5 {tmp_path}TM_from_april.nc {tmp_path}TM_minus_5.nc", shell=True) subprocess.run(
f"cdo -subc,5 {tmp_path}TM_from_april.nc {tmp_path}TM_minus_5.nc", shell=True
)
# Create an .nc file with all cells containing 1 if cell in tgminus5.nc is >=0, 0 otherwise # Create an .nc file with all cells containing 1 if cell in tgminus5.nc is >=0, 0 otherwise
subprocess.run(f"cdo -gtc,0 {tmp_path}TM_minus_5.nc {tmp_path}TM_minus_5_gtc.nc", shell=True) subprocess.run(
f"cdo -gtc,0 {tmp_path}TM_minus_5.nc {tmp_path}TM_minus_5_gtc.nc", shell=True
)
# Multiplying tgminus5.nc with tgminus5gtc.nc, so that all sub zero values are set to 0 # Multiplying tgminus5.nc with tgminus5gtc.nc, so that all sub zero values are set to 0
subprocess.run(f"cdo -mul {tmp_path}TM_minus_5.nc {tmp_path}TM_minus_5_gtc.nc {tmp_path}TM_minus_5_nozero.nc", shell=True) subprocess.run(
f"cdo -mul {tmp_path}TM_minus_5.nc {tmp_path}TM_minus_5_gtc.nc {tmp_path}TM_minus_5_nozero.nc",
shell=True,
)
# Accumulate day degrees with the corrected values # Accumulate day degrees with the corrected values
subprocess.run(f"cdo -timcumsum {tmp_path}TM_minus_5_nozero.nc {tmp_path}DD_unmasked.nc", shell=True) subprocess.run(
f"cdo -timcumsum {tmp_path}TM_minus_5_nozero.nc {tmp_path}DD_unmasked.nc",
shell=True,
)
# Mask results using a CSV file with polygons # Mask results using a CSV file with polygons
...@@ -57,7 +72,10 @@ subprocess.run(f"cdo -timcumsum {tmp_path}TM_minus_5_nozero.nc {tmp_path}DD_unma ...@@ -57,7 +72,10 @@ subprocess.run(f"cdo -timcumsum {tmp_path}TM_minus_5_nozero.nc {tmp_path}DD_unma
if os.getenv("MASK_FILE") is not None: if os.getenv("MASK_FILE") is not None:
mask_file = os.getenv("MASK_FILE") mask_file = os.getenv("MASK_FILE")
print(f"Applying mask file {mask_file} to result.nc") 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) subprocess.run(
f"cdo -maskregion,{mask_file} {tmp_path}DD_unmasked.nc {tmp_path}DD.nc",
shell=True,
)
else: else:
os.rename(f"{tmp_path}DD_unmasked.nc", f"{tmp_path}DD.nc") os.rename(f"{tmp_path}DD_unmasked.nc", f"{tmp_path}DD.nc")
...@@ -68,66 +86,87 @@ else: ...@@ -68,66 +86,87 @@ else:
# 4: 360 < DD <= 560 # Peak flight period # 4: 360 < DD <= 560 # Peak flight period
# 5: 560 < DD # 1st generation flight period ended # 5: 560 < DD # 1st generation flight period ended
# (A>0)*(A<=260)*2 + (A>260)*(A<=360)*3 + (A>360)*(A<=560)*4 + (A>560)*5 # (A>0)*(A<=260)*2 + (A>260)*(A<=360)*3 + (A>360)*(A<=560)*4 + (A>560)*5
subprocess.run(f'cdo -aexpr,"WARNING_STATUS=(TM>0)*(TM<=260)*2 + (TM>260)*(TM<=360)*3 + (TM>360)*(TM<=560)*4 + (TM>560)*5" {tmp_path}DD.nc {tmp_path}result.nc', shell=True) subprocess.run(
f'cdo -aexpr,"WARNING_STATUS=(TM>0)*(TM<=260)*2 + (TM>260)*(TM<=360)*3 + (TM>360)*(TM<=560)*4 + (TM>560)*5" {tmp_path}DD.nc {tmp_path}result.nc',
shell=True,
)
#""" # """
# Split the combined file into daily .nc files again, with YYYY-MM-DD in the filename. Convert to corresponding GeoTIFF files # Split the combined file into daily .nc files again, with YYYY-MM-DD in the filename. Convert to corresponding GeoTIFF files
# Variables that needs discrete classification, must be integers in order for mapserver to work properly (Don't ask!) # Variables that needs discrete classification, must be integers in order for mapserver to work properly (Don't ask!)
# Since we need WARNING_STATUS to be discretely classified, we need to create a separate GeoTIFF file for it # Since we need WARNING_STATUS to be discretely classified, we need to create a separate GeoTIFF file for it
result = nc.Dataset(f'{tmp_path}result.nc', 'r') result = nc.Dataset(f"{tmp_path}result.nc", "r")
timesteps = result.variables["time"][:] timesteps = result.variables["time"][:]
#print(timesteps) # print(timesteps)
timestep_index = 1 timestep_index = 1
timestep_dates = [] timestep_dates = []
for timestep in timesteps: for timestep in timesteps:
timestep_date = datetime.fromtimestamp(timestep).astimezone(local_timezone) timestep_date = datetime.fromtimestamp(timestep).astimezone(local_timezone)
file_date = timestep_date.astimezone(local_timezone).strftime("%Y-%m-%d") 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) timestep_dates.append(file_date)
# We must delete the GeoTIFF file before merging # 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 {tmp_path}result_*{file_date}_lcc.tif", shell=True)
subprocess.run(f'rm {out_path}result_*{file_date}.tif', shell=True) subprocess.run(f"rm {out_path}result_*{file_date}.tif", shell=True)
with open("/dev/null", "w") as devnull: 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(
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) 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 # 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(
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) 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 timestep_index = timestep_index + 1
if len(timestep_dates) != len(set(timestep_dates)): if len(timestep_dates) != len(set(timestep_dates)):
print("ERROR: Something is wrong with your timesteps") print("ERROR: Something is wrong with your timesteps")
#print(timestep_dates) # print(timestep_dates)
#print(len(timestep_dates)) # print(len(timestep_dates))
# Generate mapfile # Generate mapfile
# Building data sets for language specific legends # Building data sets for language specific legends
languages = [] languages = []
language_codes = config["i18n"]["languages"].split(","); language_codes = config["i18n"]["languages"].split(",")
for language_code in language_codes: for language_code in language_codes:
language = {"language_code": language_code} language = {"language_code": language_code}
if ("i18n.%s" % language_code) in config: if ("i18n.%s" % language_code) in config:
for keyword in config["i18n.%s" % language_code]: for keyword in config["i18n.%s" % language_code]:
language[keyword] = config["i18n.%s" % language_code][keyword] language[keyword] = config["i18n.%s" % language_code][keyword]
languages.append(language) languages.append(language)
# The paths should be set in a .env file # The paths should be set in a .env file
env = Environment(loader=FileSystemLoader('.')) env = Environment(loader=FileSystemLoader("."))
template = env.get_template("mapfile/template.j2") template = env.get_template("mapfile/template.j2")
output = template.render({ output = template.render(
"timestep_dates": timestep_dates, {
"mapserver_data_dir": os.getenv("MAPSERVER_DATA_DIR"), "timestep_dates": timestep_dates,
"mapserver_mapfile_dir": os.getenv("MAPSERVER_MAPFILE_DIR"), "mapserver_data_dir": os.getenv("MAPSERVER_DATA_DIR"),
"mapserver_log_file": os.getenv("MAPSERVER_LOG_FILE"), "mapserver_mapfile_dir": os.getenv("MAPSERVER_MAPFILE_DIR"),
"mapserver_image_path": os.getenv("MAPSERVER_IMAGE_PATH"), "mapserver_log_file": os.getenv("MAPSERVER_LOG_FILE"),
"mapserver_extent": os.getenv("MAPSERVER_EXTENT"), "mapserver_image_path": os.getenv("MAPSERVER_IMAGE_PATH"),
"languages": languages, "mapserver_extent": os.getenv("MAPSERVER_EXTENT"),
"language_codes": language_codes "languages": languages,
}) "language_codes": language_codes,
}
)
mapfile_outdir = os.getenv("MAPFILE_DIR") mapfile_outdir = os.getenv("MAPFILE_DIR")
with open(f"{mapfile_outdir}/PSILARTEMP.map", 'w') as f: with open(f"{mapfile_outdir}/PSILARTEMP.map", "w") as f:
f.write(output) f.write(output)
\ No newline at end of file
# Use this as an example to create your own .env file # Use this as an example to create your own .env file
# Year to run for
YEAR=2024
# Where your application resides # Where your application resides
HOME_DIR=/home/foo/ HOME_DIR=/home/foo/
# Path to the weather data # Path to the weather data
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment