Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Height factor in configs and reading tiff. #55

Merged
merged 1 commit into from
Feb 15, 2024
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
Original file line number Diff line number Diff line change
Expand Up @@ -4,6 +4,9 @@ VAR_CODE = pd
CELL_METHODS = time: mean
UNIT = hPa
UNIT_CONVERSION = 1.0
HEIGHT_CORRECTION_FACTOR = 0.00025
# After the height correction allows removing negative values by setting them to zero
TRUNCATE_NEGATIVE_VALUES = True
VALUE_MIN = 0
VALUE_MAX = 50
VALUE_SCALE = 0.1
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -13,7 +13,7 @@ STANDARD_NAME = precipitation_amount
LONG_NAME = Daily Accumulated Precipitation
# 1303 - ERAinterim
# 1329 - ERA5-land
KIWIS_FILTER_PLUGIN_CLASSES = {'ObservationsKiwisFilter': {'1303': 100.0}}
KIWIS_FILTER_PLUGIN_CLASSES = {'ObservationsKiwisFilter': {'1303': 100.0, '1329': 100.0}}

[VAR_TIME]

Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -16,7 +16,7 @@ LONG_NAME = 6 Hourly Accumulated Precipitation
# 1295 - MARS
# 1303 - ERAinterim
# 1329 - ERA5-land
KIWIS_FILTER_PLUGIN_CLASSES = {'DowgradedObservationsKiwisFilter': {'1304': 1.0, '1302': 1.0, '1295': 1.0}, 'ObservationsKiwisFilter': {'1303': 100.0}}
KIWIS_FILTER_PLUGIN_CLASSES = {'DowgradedObservationsKiwisFilter': {'1304': 1.0, '1302': 1.0, '1295': 1.0}, 'ObservationsKiwisFilter': {'1303': 100.0, '1329': 100.0}}


[VAR_TIME]
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -4,6 +4,7 @@ VAR_CODE = ta
CELL_METHODS = time: mean
UNIT = degree_Celsius
UNIT_CONVERSION = 1.0
HEIGHT_CORRECTION_FACTOR = 0.006
VALUE_MIN = -50
VALUE_MAX = 50
VALUE_SCALE = 0.1
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -4,6 +4,7 @@ VAR_CODE = ta6
CELL_METHODS = time: mean
UNIT = degree_Celsius
UNIT_CONVERSION = 1.0
HEIGHT_CORRECTION_FACTOR = 0.006
VALUE_MIN = -52
VALUE_MAX = 50
VALUE_SCALE = 0.1
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -4,6 +4,7 @@ VAR_CODE = tn
CELL_METHODS = time: minimum
UNIT = degree_Celsius
UNIT_CONVERSION = 1.0
HEIGHT_CORRECTION_FACTOR = 0.006
VALUE_MIN = -55
VALUE_MAX = 50
VALUE_SCALE = 0.1
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -4,6 +4,7 @@ VAR_CODE = tx
CELL_METHODS = time: maximum
UNIT = degree_Celsius
UNIT_CONVERSION = 1.0
HEIGHT_CORRECTION_FACTOR = 0.006
VALUE_MIN = -55
VALUE_MAX = 55
VALUE_SCALE = 0.1
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -21,6 +21,10 @@ VAR_CODE = DUMMY_VAR
CELL_METHODS = time: mean
UNIT = DUMMY_UNIT
UNIT_CONVERSION = 1.0
# A value of 0.0 means no height correction will be done
HEIGHT_CORRECTION_FACTOR = 0.0
# After the height correction allows removing negative values by setting them to zero
TRUNCATE_NEGATIVE_VALUES = False
VALUE_MIN = 0
VALUE_MAX = 35
VALUE_SCALE = 0.1
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -4,6 +4,9 @@ VAR_CODE = pd
CELL_METHODS = time: mean
UNIT = hPa
UNIT_CONVERSION = 1.0
HEIGHT_CORRECTION_FACTOR = 0.00025
# After the height correction allows removing negative values by setting them to zero
TRUNCATE_NEGATIVE_VALUES = True
VALUE_MIN = 0
VALUE_MAX = 50
VALUE_SCALE = 0.1
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -4,6 +4,7 @@ VAR_CODE = ta
CELL_METHODS = time: mean
UNIT = degree_Celsius
UNIT_CONVERSION = 1.0
HEIGHT_CORRECTION_FACTOR = 0.006
VALUE_MIN = -50
VALUE_MAX = 50
VALUE_SCALE = 0.1
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -4,6 +4,7 @@ VAR_CODE = ta6
CELL_METHODS = time: mean
UNIT = degree_Celsius
UNIT_CONVERSION = 1.0
HEIGHT_CORRECTION_FACTOR = 0.006
VALUE_MIN = -52
VALUE_MAX = 50
VALUE_SCALE = 0.1
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -4,6 +4,7 @@ VAR_CODE = tn
CELL_METHODS = time: minimum
UNIT = degree_Celsius
UNIT_CONVERSION = 1.0
HEIGHT_CORRECTION_FACTOR = 0.006
VALUE_MIN = -55
VALUE_MAX = 50
VALUE_SCALE = 0.1
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -4,6 +4,7 @@ VAR_CODE = tx
CELL_METHODS = time: maximum
UNIT = degree_Celsius
UNIT_CONVERSION = 1.0
HEIGHT_CORRECTION_FACTOR = 0.006
VALUE_MIN = -55
VALUE_MAX = 55
VALUE_SCALE = 0.1
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -22,6 +22,10 @@ VAR_CODE = DUMMY_VAR
CELL_METHODS = time: mean
UNIT = DUMMY_UNIT
UNIT_CONVERSION = 1.0
# A value of 0.0 means no height correction will be done
HEIGHT_CORRECTION_FACTOR = 0.0
# After the height correction allows removing negative values by setting them to zero
TRUNCATE_NEGATIVE_VALUES = False
VALUE_MIN = 0
VALUE_MAX = 35
VALUE_SCALE = 0.1
Expand Down
37 changes: 27 additions & 10 deletions src/lisfloodutilities/gridding/generate_grids.py
Original file line number Diff line number Diff line change
Expand Up @@ -19,6 +19,7 @@
import os
import csv
from pathlib import Path
import numpy as np
from argparse import ArgumentParser, ArgumentTypeError
from datetime import datetime, timedelta
from lisfloodutilities.gridding.lib.utils import Printable, Dem, Config, FileUtils, GriddingUtils, KiwisLoader
Expand All @@ -45,8 +46,23 @@ def memory_save_mode_type(mode: str) -> str:
raise ArgumentTypeError(f'You must select a RAM save mode out of {list(Config.MEMORY_SAVE_MODES.keys())}.')
return mode

def get_grid(grid_utils: GriddingUtils, filename: Path, tiff_filepath: Path=None, get_existing_tiff: bool=False) -> np.ndarray:
grid_data = None
# read the existing TIFF instead of interpolating the grid
if get_existing_tiff:
try:
grid_data = grid_utils.read_tiff(tiff_filepath)
except Exception as e:
print_msg(f"Warning: The file {tiff_filepath} is probably missing or is corrupt. It will be regenerated by this tool.")
# Interpolate the grid out of the txt input file
grid_data = grid_utils.generate_grid(filename)
else:
# Interpolate the grid out of the txt input file
grid_data = grid_utils.generate_grid(filename)
return grid_data

def run(config_filename: str, infolder: str, output_file: str, processing_dates_file: str, file_utils: FileUtils,
output_tiff: bool, output_netcdf: bool, overwrite_output: bool, use_existing_file: bool,
output_tiff: bool, output_netcdf: bool, overwrite_output: bool, use_existing_file: bool, get_existing_tiff: bool,
start_date: datetime = None, end_date: datetime = None, interpolation_mode: str = 'adw',
use_broadcasting: bool = False, memory_save_mode: str = None):
"""
Expand Down Expand Up @@ -77,26 +93,23 @@ def run(config_filename: str, infolder: str, output_file: str, processing_dates_
str(conf.dem_max_y)))

print_msg('Start reading files')
# inwildcard = conf.var_code + FileUtils.FILES_WILDCARD

netcdf_offset_file_date = int(conf.get_config_field('VAR_TIME','OFFSET_FILE_DATE'))

outfile = output_file
if output_tiff:
output_writer_tiff = GDALWriter(conf, overwrite_output, quiet_mode)
if output_netcdf:
output_writer_netcdf = NetCDFWriter(conf, overwrite_output, quiet_mode)
output_writer_netcdf.open(Path(outfile))
output_writer_netcdf.open(Path(output_file))
file_loader = KiwisLoader(conf, Path(infolder), dates_to_process, overwrite_output, use_existing_file, quiet_mode)
for filename, kiwis_timestamp_str in file_loader:
# file_timestamp = file_utils.get_timestamp_from_filename(filename) + timedelta(days=netcdf_offset_file_date)
kiwis_timestamp = datetime.strptime(kiwis_timestamp_str, FileUtils.DATE_PATTERN_CONDENSED_SHORT)
file_timestamp = kiwis_timestamp + timedelta(days=netcdf_offset_file_date)
print_msg(f'Processing file: {filename}')
tiff_filepath = filename.with_suffix('.tiff')
if output_tiff:
outfilepath = filename.with_suffix('.tiff')
output_writer_tiff.open(outfilepath)
grid_data = grid_utils.generate_grid(filename)
output_writer_tiff.open(tiff_filepath)
grid_data = get_grid(grid_utils, filename, tiff_filepath, get_existing_tiff)
if output_netcdf:
output_writer_netcdf.write(grid_data, file_timestamp)
if output_tiff:
Expand Down Expand Up @@ -142,6 +155,7 @@ def main(argv):
out_netcdf=False,
overwrite_output=False,
use_existing_file=False,
get_existing_tiff=False,
start_date='',
end_date=END_DATE_DEFAULT,
interpolation_mode='adw',
Expand Down Expand Up @@ -181,6 +195,8 @@ def main(argv):
help="Force write to existing file. TIFF files will be overwritten and netCDF file will be appended. [default: %(default)s]")
parser.add_argument("-u", "--useexisting", dest="use_existing_file", action="store_true",
help="Force to use existing point/txt filenames, so these files will be used for gridding. [default: %(default)s]")
parser.add_argument("-g", "--getexistingtiff", dest="get_existing_tiff", action="store_true",
help="Force to use existing tiff files instead of interpolating the values again. [default: %(default)s]")
parser.add_argument("-m", "--mode", dest="interpolation_mode", required=False, type=interpolation_mode_type,
help="Set interpolation mode. [default: %(default)s]",
metavar=f"{list(Config.INTERPOLATION_MODES.keys())}")
Expand Down Expand Up @@ -240,15 +256,16 @@ def main(argv):
print_msg(f"Input Folder: {args.in_file_or_folder}")
print_msg(f"Overwrite Output: {args.overwrite_output}")
print_msg(f"Use Existing Point Files: {args.use_existing_file}")
print_msg(f"Use Existing TIFF Files: {args.get_existing_tiff}")
print_msg(f"Interpolation Mode: {args.interpolation_mode}")
print_msg(f"RAM Save Mode: {args.memory_save_mode}")
print_msg(f"Broadcasting: {args.use_broadcasting}")
print_msg(f"Processing Dates File: {args.processing_dates_file}")
print_msg(f"Config File: {config_filename}")

run(config_filename, args.in_file_or_folder, args.output_file, args.processing_dates_file, file_utils, args.out_tiff,
args.out_netcdf, args.overwrite_output, args.use_existing_file, start_date, end_date, args.interpolation_mode,
args.use_broadcasting, args.memory_save_mode)
args.out_netcdf, args.overwrite_output, args.use_existing_file, args.get_existing_tiff, start_date, end_date,
args.interpolation_mode, args.use_broadcasting, args.memory_save_mode)
except Exception as e:
indent = len(program_name) * " "
sys.stderr.write(program_name + ": " + repr(e) + "\n")
Expand Down
28 changes: 21 additions & 7 deletions src/lisfloodutilities/gridding/lib/utils.py
Original file line number Diff line number Diff line change
Expand Up @@ -21,6 +21,7 @@
import numpy as np
import numpy.ma as ma
import pandas as pd
import rasterio
from decimal import *
from datetime import datetime, timedelta
from collections import OrderedDict
Expand Down Expand Up @@ -199,7 +200,6 @@ def __init__(self, config_filename: str, start_date: datetime = None, end_date:
self.interpolation_mode = interpolation_mode
self.memory_save_mode = memory_save_mode
self.__setup_variable_config(config_filename)
self.__HEIGHT_CORRECTION_FACTOR = {"tx": 0.006, "tn": 0.006, "ta": 0.006, "ta6": 0.006, "pd": 0.00025}

self._dem_file = Path(self.config_path).joinpath('dem.nc')
self._dem_interpolation_file = Path(self.config_path).joinpath('dem.wgs')
Expand Down Expand Up @@ -352,11 +352,15 @@ def var_code(self) -> str:

@property
def do_height_correction(self) -> bool:
return self.var_code in self.__HEIGHT_CORRECTION_FACTOR
return self.height_correction_factor != 0.0

@property
def height_correction_factor(self) -> float:
return self.__HEIGHT_CORRECTION_FACTOR[self.var_code]
return float(self.get_config_field('PROPERTIES', 'HEIGHT_CORRECTION_FACTOR'))

@property
def truncate_negative_values(self) -> bool:
return str(self.get_config_field('PROPERTIES', 'TRUNCATE_NEGATIVE_VALUES')).lower() in ('true', '1')

@property
def neighbours_near(self) -> int:
Expand Down Expand Up @@ -419,11 +423,11 @@ def reset_height(self, result: ma.MaskedArray) -> np.ndarray:
mask=self.conf.dem_mask,
fill_value=self.conf.VALUE_NAN)
result -= grid_cells_height * self.conf.height_correction_factor
if self.conf.var_code == "pd":
# Eliminating the eventual negative values
# from Vapor Pressure resulting from height correction
result[np.where(result<0)] = 0
self.print_msg('Finish reseting height correction')
if self.conf.truncate_negative_values:
# Eliminating the eventual negative values like
# for Vapor Pressure resulting from height correction
result[np.where(result<0)] = 0
result = result.filled()
return result

Expand Down Expand Up @@ -453,6 +457,16 @@ def check_grid_nan(self, filename: Path, result: np.ndarray):
if count_dem_nan != count_grid_nan:
print(f"WARNING: The grid interpolated from file {filename.name} contains different NaN values ({count_grid_nan}) than the DEM ({count_dem_nan}). diff: {count_grid_nan - count_dem_nan}")

def read_tiff(self, tiff_filepath: Path) -> np.ndarray:
# This method could be implemented on a GDALReader class
src = rasterio.open(tiff_filepath)
data = src.read(1)
scale = src.scales[0]
offset = src.offsets[0]
data = data * scale + offset
src.close()
return self.prepare_grid(data, self.conf.dem.lons.shape)

def generate_grid(self, filename: Path) -> np.ndarray:
mv_target = self.conf.dem.mv
mv_source = None
Expand Down
9 changes: 6 additions & 3 deletions tests/test_gridding.py
Original file line number Diff line number Diff line change
Expand Up @@ -57,6 +57,7 @@ def test_generate_netcdf(self):
infolder = input_folder
overwrite_output = True
use_existing_file = True
get_existing_tiff = False
outfolder_or_file = output_netcdf
processing_dates_file = None
interpolation_mode = 'adw'
Expand All @@ -76,7 +77,7 @@ def test_generate_netcdf(self):
end_date = datetime.strptime(end_date_str, FileUtils.DATE_PATTERN_CONDENSED)

run(config_filename, infolder, outfolder_or_file, processing_dates_file,
file_utils, out_tiff, out_netcdf, overwrite_output, use_existing_file, start_date, end_date,
file_utils, out_tiff, out_netcdf, overwrite_output, use_existing_file, get_existing_tiff, start_date, end_date,
interpolation_mode=interpolation_mode, use_broadcasting=use_broadcasting, memory_save_mode=memory_save_mode)

reference = Dataset(reference_output)
Expand Down Expand Up @@ -108,6 +109,7 @@ def test_update_netcdf(self):
infolder = input_folder
overwrite_output = True
use_existing_file = True
get_existing_tiff = False
outfolder_or_file = output_netcdf
interpolation_mode = 'adw'
use_broadcasting = False
Expand All @@ -125,7 +127,7 @@ def test_update_netcdf(self):
end_date = datetime.strptime(end_date_str, FileUtils.DATE_PATTERN_CONDENSED)

run(config_filename, infolder, outfolder_or_file, processing_dates_file,
file_utils, out_tiff, out_netcdf, overwrite_output, use_existing_file, start_date, end_date,
file_utils, out_tiff, out_netcdf, overwrite_output, use_existing_file, get_existing_tiff, start_date, end_date,
interpolation_mode=interpolation_mode, use_broadcasting=use_broadcasting, memory_save_mode=memory_save_mode)

reference = Dataset(reference_output)
Expand Down Expand Up @@ -159,6 +161,7 @@ def test_generate_tiff(self):
infolder = input_folder
overwrite_output = True
use_existing_file = True
get_existing_tiff = False
outfolder_or_file = str(Path(output_folder, 'output.nc'))
processing_dates_file = None
interpolation_mode = 'adw'
Expand All @@ -177,7 +180,7 @@ def test_generate_tiff(self):
end_date = datetime.strptime(end_date_str, FileUtils.DATE_PATTERN_CONDENSED)

run(config_filename, infolder, outfolder_or_file, processing_dates_file,
file_utils, out_tiff, out_netcdf, overwrite_output, use_existing_file, start_date, end_date,
file_utils, out_tiff, out_netcdf, overwrite_output, use_existing_file, get_existing_tiff, start_date, end_date,
interpolation_mode=interpolation_mode, use_broadcasting=use_broadcasting, memory_save_mode=memory_save_mode)

# Remove netcdf output file since we are only interested in the output tiffs
Expand Down
Loading