Skip to content

Commit 6edbe09

Browse files
authored
Merge pull request #65 from ec-jrc/lisvap1h
Issue #64: Sub-daily Run
2 parents 5d92fac + fe65fa7 commit 6edbe09

Some content is hidden

Large Commits have some content hidden by default. Use the searchbox below for content that may be hidden.

53 files changed

+303
-74
lines changed

.gitignore

Lines changed: 3 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -18,7 +18,10 @@ lisflood_lisvap.egg-info/
1818
!tests/data/tests_cordex.xml
1919
!/tests/data/tests_glofas.xml
2020
!tests/data/tests_efas_1arcmin.xml
21+
!tests/data/tests_efas_1arcmin_6hourly.xml
22+
!tests/data/tests_efas_1arcmin_hourly.xml
2123
!tests/data/tests_efas_1arcmin_yearly.xml
24+
!tests/data/tests_efas_1arcmin_yearly_output.xml
2225
!tests/data/tests_efas_1arcmin_360days_calendar.xml
2326

2427
# unignoring reference and input data for tests

Dockerfile

Lines changed: 4 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -1,5 +1,5 @@
1-
FROM python:3.7-stretch
2-
MAINTAINER Domenico Nappo <[email protected]>
1+
FROM python:3.7-buster
2+
MAINTAINER Goncalo Gomes <[email protected]>
33

44
ENV no_proxy=jrc.it,localhost,ies.jrc.it,127.0.0.1,jrc.ec.europa.eu
55
ENV ftp_proxy=http://10.168.209.72:8012
@@ -24,12 +24,13 @@ RUN wget -q http://pcraster.geo.uu.nl/pcraster/4.2.1/pcraster-4.2.1.tar.bz2 && t
2424

2525
COPY requirements.txt /
2626
RUN /usr/local/bin/pip3.7 install -U pip && /usr/local/bin/pip3.7 install -r /requirements.txt \
27-
&& cd /usr/lib/x86_64-linux-gnu/ && ln -s libboost_python-py35.so libboost_python3.so
27+
&& cd /usr/lib/x86_64-linux-gnu/ && [ ! -f libboost_python3.so ] && ln -s libboost_python-py35.so libboost_python3.so || true
2828

2929
WORKDIR /opt/pcraster-4.2.1/build
3030
RUN cmake -DFERN_BUILD_ALGORITHM:BOOL=TRUE -DCMAKE_INSTALL_PREFIX:PATH=/opt/pcraster -DPYTHON_EXECUTABLE:FILEPATH=/usr/local/bin/python3.7 ../ \
3131
&& cmake --build ./ && make install
3232

33+
WORKDIR /
3334
COPY tests/. /tests/
3435
COPY basemaps/. /basemaps/
3536
COPY src/. /

docker-entrypoint.sh

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -6,7 +6,7 @@ if [[ "$1" = 'usecases' ]]; then
66
mkdir -p /input/efas
77
mkdir -p /input/basemaps
88
echo "Copying basemaps to /input/..."
9-
cp /basemaps/* /input/basemaps/
9+
cp -r /basemaps/* /input/basemaps/
1010
echo "copy input files to /input/cordex......"
1111
cp /tests/data/input/cordex/*.nc /input/cordex/
1212
echo "copy cordex settings to /input/..."

settings_tpl.xml

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -35,7 +35,7 @@ You can use $(SettingsDir) or $(SettingsPath) to refer directory containing the
3535

3636
<textvar name="DtSec" value="86400">
3737
<comment>
38-
time step [seconds] ALWAYS USE 86400!!
38+
time step [seconds] USE 86400 for daily, 43200 for 12hourly, 21600 for 6hourly and 3600 for hourly!!
3939
</comment>
4040
</textvar>
4141

src/lisvap/__init__.py

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -1,7 +1,7 @@
1-
version = version = (1, 2, 9)
1+
version = version = (1, 3, 0)
22
__authors__ = "Peter Burek, Johan van der Knijff, Ad de Roo"
33
__version__ = '.'.join(list(map(str, version)))
4-
__date__ = "22/01/2024"
4+
__date__ = "07/03/2024"
55
__copyright__ = "Copyright 2020, Lisflood Open Source"
66
__maintainers__ = "Goncalo Gomes, Domenico Nappo, Valerio Lorini, Lorenzo Mentaschi, Lorenzo Alfieri"
77
__status__ = "Operation"

src/lisvap/lisvapdynamic.py

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -135,7 +135,7 @@ def dynamic(self):
135135
# CM: get time for operation "Start dynamic"
136136
tp.timemeasure('Start dynamic')
137137
# CM: date corresponding to the model time step (yyyy-mm-dd hh:mm:ss)
138-
self.calendar_date = self.calendar_day_start + datetime.timedelta(days=(self.currentTimeStep()) * self.DtDay)
138+
self.calendar_date = self.calendar_day_start + datetime.timedelta(seconds=(self.currentTimeStep()) * self.DtSec)
139139
# CM: day of the year corresponding to the model time step
140140
self.calendar_day = int(self.calendar_date.strftime("%j"))
141141

@@ -146,7 +146,7 @@ def dynamic(self):
146146
self.time_since_start = self.currentTimeStep() - self.firstTimeStep() + 1
147147

148148
if settings.flags['loud']:
149-
print("%-6i %10s" % (self.currentTimeStep(), self.calendar_date.strftime("%d/%m/%Y")))
149+
print("%-6i %10s" % (self.currentTimeStep(), self.calendar_date.strftime("%d/%m/%Y %H:%M")))
150150
elif not settings.flags['checkfiles']:
151151
if settings.flags['quiet'] and not settings.flags['veryquiet']:
152152
sys.stdout.write(".")

src/lisvap/utils/readers.py

Lines changed: 10 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -275,9 +275,17 @@ def netcdf_step(averageyearflag, nf1, timestampflag, timestep, splitIO=False):
275275
# Time step, expressed as fraction of day (same as self.var.DtSec and self.var.DtDay)
276276
# get date of current simulation step
277277
current_date = calendar(timestep)
278+
278279
if not isinstance(current_date, datetime.datetime):
279-
current_date_number = date2num(begin, units=t_unit, calendar=t_cal) + ((current_date - 1) * DtDay)
280-
current_date = num2date(current_date_number, t_unit, t_cal)
280+
current_date_number = current_date * int(settings.binding['internal.time.frequency'])
281+
init_t_unit = settings.binding['internal.time.unit']
282+
init_t_cal = settings.binding['internal.time.calendar']
283+
begin_date_number = date2num(begin, units=init_t_unit, calendar=init_t_cal)
284+
cur_date = num2date(current_date_number, init_t_unit, init_t_cal)
285+
next_date = cur_date - datetime.timedelta(seconds=DtSec)
286+
cur_step = date2num(next_date, units=init_t_unit, calendar=init_t_cal)
287+
current_date_number = begin_date_number + cur_step
288+
current_date = num2date(current_date_number, init_t_unit, init_t_cal)
281289
# if reading from an average year NetCDF stack, ignore the year in current simulation date and change it to the netCDF time unit year
282290
if averageyearflag:
283291
# CM: get year from time unit in case average year is used

src/lisvap/utils/tools.py

Lines changed: 8 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -29,6 +29,7 @@
2929
from pcraster.operations import scalar, defined, maptotal, ifthenelse, mapminimum, mapmaximum
3030

3131
from . import LisfloodError
32+
from . import LisSettings
3233
from .decorators import counted
3334

3435

@@ -97,10 +98,16 @@ def get_calendar_configuration(netcdf_file_obj, settings=None):
9798
except AttributeError: # Attribute does not exist
9899
t_unit = u'hours since 1990-01-01 06:00:00'
99100
t_cal = u'gregorian' # Use standard calendar
101+
try:
102+
t_frequency = int(netcdf_file_obj.variables['time'].frequency) # get frequency
103+
except AttributeError: # Attribute does not exist
104+
t_frequency = 1
100105
if settings is not None and not ('internal.time.unit' in settings.binding and
101-
'internal.time.calendar' in settings.binding):
106+
'internal.time.calendar' in settings.binding and
107+
'internal.time.frequency' in settings.binding):
102108
settings.binding['internal.time.unit'] = t_unit
103109
settings.binding['internal.time.calendar'] = t_cal
110+
settings.binding['internal.time.frequency'] = t_frequency
104111
return t_unit, t_cal
105112

106113

@@ -142,7 +149,6 @@ def date_to_int(date_in, both=False):
142149
the number of steps as integer is returned
143150
:return: number of steps as integer and input date as string
144151
"""
145-
from lisvap.utils import LisSettings
146152
settings = LisSettings.instance()
147153
binding = settings.binding
148154
# CM: get reference date to be used with step numbers from 'CalendarDayStart' in Settings.xml file
@@ -175,7 +181,6 @@ def checkdate(start, end):
175181
:param start: start date for model run (# or date as string)
176182
:param end: end date for model run (# or date as string)
177183
"""
178-
from . import LisSettings
179184
settings = LisSettings.instance()
180185
binding = settings.binding
181186
# CM: calendar date start (CalendarDayStart)

src/lisvap/utils/writers.py

Lines changed: 60 additions & 50 deletions
Original file line numberDiff line numberDiff line change
@@ -63,55 +63,53 @@ def coordinates_range(start=0, nelems=1, step=1):
6363
return elem_array
6464

6565

66-
def get_output_parameters_monthly(start_date, timestep, current_output_index):
66+
def get_output_parameters_monthly(start_date, timestep, time_frequency, timestep_stride, current_output_index):
6767
output_index = current_output_index
68-
start_yearmonth = start_date.strftime('%Y%m')
69-
current_date = start_date + datetime.timedelta(days=timestep-1)
70-
current_yearmonth = current_date.strftime('%Y%m')
71-
filename_suffix = current_yearmonth
72-
if start_yearmonth != current_yearmonth:
73-
first_day_current_month = current_date.replace(day=1, hour=0, minute=0, second=0, microsecond=0)
74-
current_day = current_date.replace(hour=0, minute=0, second=0, microsecond=0)
75-
days_between = (current_day - first_day_current_month).days
76-
if days_between == 0:
77-
# Get last month because the 1st of the month belongs in the last month file
78-
last_day_last_month = first_day_current_month - datetime.timedelta(days=1)
79-
filename_suffix = (last_day_last_month).strftime('%Y%m')
80-
num_days_last_month = last_day_last_month.day
81-
# If last month was complete the idx needs to be set to the last idx of the previous file
82-
if output_index > num_days_last_month:
83-
output_index = num_days_last_month - 1
84-
else:
85-
# Since the 1st day belongs to the last year file, the remaining days need to start on index 0...
86-
output_index = days_between - 1
68+
current_date = start_date + datetime.timedelta(seconds=((timestep - 1) * timestep_stride))
69+
filename_suffix = current_date.strftime('%Y%m')
70+
71+
first_date_current_month = start_date.replace(year=current_date.year, month=current_date.month)
72+
seconds_between = (current_date - first_date_current_month).total_seconds()
73+
num_steps_done_in_current_month = int(seconds_between / timestep_stride) + 1
74+
last_date_last_month = first_date_current_month - datetime.timedelta(seconds=timestep_stride)
75+
day_inside_last_month = last_date_last_month - datetime.timedelta(seconds=timestep_stride)
76+
77+
if current_date == first_date_current_month:
78+
output_index = 0
79+
elif current_date == last_date_last_month:
80+
filename_suffix = day_inside_last_month.strftime('%Y%m')
81+
first_date_last_month = day_inside_last_month.replace(day=1) + datetime.timedelta(seconds=(2 * timestep_stride))
82+
num_steps_done_in_last_month = int((last_date_last_month - first_date_last_month).total_seconds() / timestep_stride) + 1
83+
output_index = num_steps_done_in_last_month
84+
elif current_output_index >= num_steps_done_in_current_month:
85+
output_index = num_steps_done_in_current_month - 1
8786
return filename_suffix, output_index
8887

8988

90-
def get_output_parameters_yearly(start_date, timestep, current_output_index):
89+
def get_output_parameters_yearly(start_date, timestep, time_frequency, timestep_stride, current_output_index):
9190
output_index = current_output_index
92-
start_year = start_date.strftime('%Y')
93-
current_date = start_date + datetime.timedelta(days=timestep-1)
94-
current_year = current_date.strftime('%Y')
95-
filename_suffix = current_year
96-
if start_year != current_year:
97-
first_day_current_year = current_date.replace(month=1, day=1, hour=0, minute=0, second=0, microsecond=0)
98-
current_day = current_date.replace(hour=0, minute=0, second=0, microsecond=0)
99-
days_between = (current_day - first_day_current_year).days
100-
if days_between == 0:
101-
# Get last year because the 1st of January belongs in the last year file
102-
last_day_last_year = first_day_current_year - datetime.timedelta(days=1)
103-
filename_suffix = (last_day_last_year).strftime('%Y')
104-
num_days_last_year = last_day_last_year.timetuple().tm_yday
105-
# If last year was complete the idx needs to be set to the last idx of the previous file
106-
if output_index > num_days_last_year:
107-
output_index = num_days_last_year - 1
108-
else:
109-
# Since the 1st day belongs to the last year file, the remaining days need to start on index 0...
110-
output_index = days_between - 1
91+
current_date = start_date + datetime.timedelta(seconds=((timestep - 1) * timestep_stride))
92+
filename_suffix = current_date.strftime('%Y')
93+
94+
first_date_current_year = start_date.replace(year=current_date.year)
95+
seconds_between = (current_date - first_date_current_year).total_seconds()
96+
num_steps_done_in_current_year = int(seconds_between / timestep_stride) + 1
97+
last_date_last_year = first_date_current_year - datetime.timedelta(seconds=timestep_stride)
98+
day_inside_last_year = last_date_last_year - datetime.timedelta(seconds=timestep_stride)
99+
100+
if current_date == first_date_current_year:
101+
output_index = 0
102+
elif current_date == last_date_last_year:
103+
filename_suffix = day_inside_last_year.strftime('%Y')
104+
first_date_last_year = first_date_current_year.replace(year=current_date.year - 1)
105+
num_steps_done_in_last_year = int((last_date_last_year - first_date_last_year).total_seconds() / timestep_stride) + 1
106+
output_index = num_steps_done_in_last_year - 1
107+
elif current_output_index >= num_steps_done_in_current_year:
108+
output_index = num_steps_done_in_current_year - 1
111109
return filename_suffix, output_index
112110

113111

114-
def get_output_parameters(settings, netcdf_output_file, start_date, timestep, current_output_index):
112+
def get_output_parameters(settings, netcdf_output_file, start_date, timestep, time_frequency, timestep_stride, current_output_index):
115113
output_index = current_output_index
116114
p = Path(netcdf_output_file)
117115
netfile = Path(p.parent) / Path('{}.nc'.format(p.name) if not p.name.endswith('.nc') else p.name)
@@ -120,9 +118,9 @@ def get_output_parameters(settings, netcdf_output_file, start_date, timestep, cu
120118
if splitIO:
121119
monthlyIO = settings.get_option('monthlyOutput')
122120
if monthlyIO:
123-
filename_suffix, output_index = get_output_parameters_monthly(start_date, timestep, current_output_index)
121+
filename_suffix, output_index = get_output_parameters_monthly(start_date, timestep, time_frequency, timestep_stride, current_output_index)
124122
else:
125-
filename_suffix, output_index = get_output_parameters_yearly(start_date, timestep, current_output_index)
123+
filename_suffix, output_index = get_output_parameters_yearly(start_date, timestep, time_frequency, timestep_stride, current_output_index)
126124
netfile = Path(p.parent) / Path('{}_{}.nc'.format(p.name, filename_suffix) if not p.name.endswith('.nc') else p.name)
127125
return prefix, netfile, output_index
128126

@@ -140,6 +138,16 @@ def set_time_dimension(settings, netcdf_obj, time_variable_name, start_date, var
140138
start_date_6hourly = start_date - datetime.timedelta(hours=18)
141139
time.units = 'hours since %s' % start_date_6hourly.strftime('%Y-%m-%d %H:%M:%S.0')
142140
time.frequency = 6
141+
elif 'internal.time.unit' in settings.binding:
142+
internal_time_units = settings.binding['internal.time.unit']
143+
# Separate the different parts of time units, ex:
144+
# "seconds since 1970-01-01 00:00:00"
145+
# "days since 1970-01-01 00:00:00"
146+
# "hours since 1970-01-01 00:00:00"
147+
# "minutes since 1970-01-01 00:00:00"
148+
time_units_parts = internal_time_units.split(' ')
149+
time.units = '%s since %s' % (time_units_parts[0], start_date.strftime('%Y-%m-%d %H:%M:%S.0'))
150+
time.frequency = int(settings.binding['internal.time.frequency'])
143151
else:
144152
time.units = 'days since %s' % start_date.strftime('%Y-%m-%d %H:%M:%S.0')
145153
time.frequency = 1
@@ -276,19 +284,21 @@ def writenet(current_output_index, inputmap, netcdf_output_file, current_timeste
276284
netfile: output netcdf filename
277285
timestep:
278286
"""
279-
start_date = calendar_day_start + datetime.timedelta(days=1)
287+
settings = LisSettings.instance()
288+
289+
timestep_stride = int(settings.binding['DtSec'])
290+
time_frequency = int(settings.binding['internal.time.frequency'])
291+
start_date = calendar_day_start + datetime.timedelta(seconds=timestep_stride)
280292
timestep = current_timestep
281293

282294
cutmap = CutMap.instance()
283295
nrows = np.abs(cutmap.cuts[3] - cutmap.cuts[2])
284296
ncols = np.abs(cutmap.cuts[1] - cutmap.cuts[0])
285297

286-
settings = LisSettings.instance()
287-
288298
time_variable = 'time'
289299
output6hourly = settings.get_option('output6hourly')
290-
prefix, netfile, output_index = get_output_parameters(settings, netcdf_output_file, start_date,
291-
timestep, current_output_index)
300+
prefix, netfile, output_index = get_output_parameters(settings, netcdf_output_file, start_date, timestep,
301+
time_frequency, timestep_stride, current_output_index)
292302

293303
# Create and setup the netcdf file when the first map needs to be stored
294304
if output_index == 0 or not netfile.exists():
@@ -304,15 +314,15 @@ def writenet(current_output_index, inputmap, netcdf_output_file, current_timeste
304314
mapnp[np.isnan(mapnp)] = (nan_value - add_offset) * scale_factor
305315
if flag_time:
306316
# In case output6hourly==True, replicate four daily maps to get the 6 hourly output (EFCC-2316)
307-
# The timestep need to increase by 4
317+
# The timestep needs to increase by 4
308318
if output6hourly:
309319
time_frequency = 6
310320
for i in range(4):
311321
map_idx = output_index * 4 + i
312322
nf1.variables[time_variable][map_idx] = (timestep * 4 - 4 + i) * time_frequency
313323
nf1.variables[prefix][map_idx, :, :] = mapnp
314324
else: # Generate daily output
315-
nf1.variables[time_variable][output_index] = timestep - 1
325+
nf1.variables[time_variable][output_index] = (timestep - 1) * time_frequency # timestep - time_frequency
316326
nf1.variables[prefix][output_index, :, :] = mapnp
317327
else:
318328
nf1.variables[prefix][:, :] = mapnp

src/lisvap1.py

Lines changed: 3 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -50,10 +50,11 @@ def lisvapexe(settings):
5050
tp = TimeProfiler()
5151
step_start = settings.binding['StepStart']
5252
step_end = settings.binding['StepEnd']
53+
timestep_stride = int(settings.binding['DtSec'])
5354
start_date, end_date = datetime.datetime.strptime(step_start, '%d/%m/%Y %H:%M'), datetime.datetime.strptime(step_end, '%d/%m/%Y %H:%M')
5455
start_date_simulation = datetime.datetime.strptime(settings.binding['CalendarDayStart'], '%d/%m/%Y %H:%M')
55-
timestep_start = (start_date - start_date_simulation).days + 1
56-
timestep_end = (end_date - start_date_simulation).days + 1
56+
timestep_start = int((start_date - start_date_simulation).total_seconds() / timestep_stride) + 1
57+
timestep_end = int((end_date - start_date_simulation).total_seconds() / timestep_stride) + 1
5758
checkdate('StepStart', 'StepEnd')
5859
print('Start date: {} ({}) - End date: {} ({})'.format(step_start, timestep_start, step_end, timestep_end))
5960

0 commit comments

Comments
 (0)