Skip to content

Commit 13eab8b

Browse files
Merge pull request #36 from MET-OM/norkyst800
Norkyst800 url varies by year
2 parents afa7715 + 528d000 commit 13eab8b

File tree

2 files changed

+126
-63
lines changed

2 files changed

+126
-63
lines changed

dnora/ocr/read_metno.py

+110-47
Original file line numberDiff line numberDiff line change
@@ -1,10 +1,11 @@
1-
from abc import ABC, abstractmethod
1+
from abc import ABC, abstractmethod
22
from copy import copy
33
import numpy as np
44
import xarray as xr
55
from subprocess import call
66
import os, glob
77
import time
8+
89
# Import objects
910
from ..grd.grd_mod import Grid
1011

@@ -13,7 +14,13 @@
1314

1415
# Import aux_funcsiliry functions
1516
from .. import msg
16-
from ..aux_funcs import create_time_stamps, u_v_from_dir, expand_area, lon_in_km, pyfimex
17+
from ..aux_funcs import (
18+
create_time_stamps,
19+
u_v_from_dir,
20+
expand_area,
21+
lon_in_km,
22+
pyfimex,
23+
)
1724

1825

1926
class NorKyst800(OceanCurrentReader):
@@ -27,7 +34,14 @@ class NorKyst800(OceanCurrentReader):
2734
No. 1 : User Manual and technical descriptions.
2835
"""
2936

30-
def __init__(self, stride: int=24, hours_per_file: int=24, last_file: str='', lead_time: int=0, program: str='pyfimex'):
37+
def __init__(
38+
self,
39+
stride: int = 24,
40+
hours_per_file: int = 24,
41+
last_file: str = "",
42+
lead_time: int = 0,
43+
program: str = "pyfimex",
44+
):
3145
"""The data is currently in daily files. Do not change the default
3246
setting unless you have a good reason to do so.
3347
"""
@@ -39,18 +53,27 @@ def __init__(self, stride: int=24, hours_per_file: int=24, last_file: str='', le
3953
self.program = program
4054
return
4155

42-
def __call__(self, grid: Grid, start_time: str, end_time: str, expansion_factor: float):
56+
def __call__(
57+
self, grid: Grid, start_time: str, end_time: str, expansion_factor: float
58+
):
4359
"""Reads in all grid points between the given times and at for the given indeces"""
4460
self.start_time = start_time
4561
self.end_time = end_time
4662

4763
start_times, end_times, file_times = create_time_stamps(
48-
start_time, end_time, self.stride, self.hours_per_file, self.last_file, self.lead_time)
64+
start_time,
65+
end_time,
66+
self.stride,
67+
self.hours_per_file,
68+
self.last_file,
69+
self.lead_time,
70+
)
4971

5072
msg.info(
51-
f"Getting ocean_current forcing from Norkyst800 from {self.start_time} to {self.end_time}")
73+
f"Getting ocean_current forcing from Norkyst800 from {self.start_time} to {self.end_time}"
74+
)
5275

53-
temp_folder = 'dnora_ocr_temp'
76+
temp_folder = "dnora_ocr_temp"
5477
if not os.path.isdir(temp_folder):
5578
os.mkdir(temp_folder)
5679
print("Creating folder %s..." % temp_folder)
@@ -59,70 +82,110 @@ def __call__(self, grid: Grid, start_time: str, end_time: str, expansion_factor:
5982
for f in glob.glob("dnora_ocr_temp/*MetNo_Norkyst800.nc"):
6083
os.remove(f)
6184

62-
6385
# Define area to search in
64-
lon_min, lon_max, lat_min, lat_max = expand_area(min(grid.lon()), max(grid.lon()), min(grid.lat()), max(grid.lat()), expansion_factor)
86+
lon_min, lon_max, lat_min, lat_max = expand_area(
87+
min(grid.lon()),
88+
max(grid.lon()),
89+
min(grid.lat()),
90+
max(grid.lat()),
91+
expansion_factor,
92+
)
6593

6694
# Setting resolution to roughly 0.8 km
67-
dlat = 0.8/111
68-
mean_lon_in_km = (lon_in_km(grid.lat()[0])+lon_in_km(grid.lat()[-1]))*0.5
69-
dlon = 0.8/mean_lon_in_km
95+
dlat = 0.8 / 111
96+
mean_lon_in_km = (lon_in_km(grid.lat()[0]) + lon_in_km(grid.lat()[-1])) * 0.5
97+
dlon = 0.8 / mean_lon_in_km
7098

7199
ocr_list = []
72-
print('Apply >>> '+ self.program)
100+
print("Apply >>> " + self.program)
73101
for n in range(len(file_times)):
74102
url = self.get_url(file_times[n])
75103

76104
msg.from_file(url)
77105
msg.plain(
78-
f"Reading ocean_current forcing data: {start_times[n]}-{end_times[n]}")
106+
f"Reading ocean_current forcing data: {start_times[n]}-{end_times[n]}"
107+
)
79108

80-
nc_fimex = f'dnora_ocr_temp/ocean_current_{n:04.0f}_MetNo_Norkyst800.nc'
109+
nc_fimex = f"dnora_ocr_temp/ocean_current_{n:04.0f}_MetNo_Norkyst800.nc"
81110
# Apply pyfimex or fimex
82-
if self.program == 'pyfimex':
83-
pyfimex(input_file=url,output_file=nc_fimex,
111+
if self.program == "pyfimex":
112+
pyfimex(
113+
input_file=url,
114+
output_file=nc_fimex,
84115
projString="+proj=latlong +ellps=sphere +a=6371000 +e=0",
85-
xAxisValues=np.arange(lon_min,lon_max+dlon,dlon),
86-
yAxisValues=np.arange(lat_min,lat_max+dlat,dlat),
87-
selectVariables=['u', 'v'],
88-
reduceTime_start = start_times[n].strftime('%Y-%m-%dT%H:%M:%S'),
89-
reduceTime_end = end_times[n].strftime('%Y-%m-%dT%H:%M:%S'))
90-
elif self.program == 'fimex':
91-
fimex_command = ['fimex', '--input.file='+url,
92-
'--interpolate.method=bilinear',
93-
'--interpolate.projString=+proj=latlong +ellps=sphere +a=6371000 +e=0',
94-
'--interpolate.xAxisValues='+ str(lon_min)+','+str(lon_min+dlon)+ ',...,'+str(lon_max)+'',
95-
'--interpolate.yAxisValues=' + str(lat_min)+','+str(lat_min+dlat)+ ',...,'+str(lat_max)+'',
96-
'--interpolate.xAxisUnit=degree', '--interpolate.yAxisUnit=degree',
97-
'--process.rotateVector.all',
98-
'--extract.selectVariables=u', '--extract.selectVariables=v',
99-
'--extract.reduceTime.start=' + \
100-
start_times[n].strftime('%Y-%m-%dT%H:%M:%S'),
101-
'--extract.reduceTime.end=' + \
102-
end_times[n].strftime('%Y-%m-%dT%H:%M:%S'),
103-
'--process.rotateVector.direction=latlon',
104-
#'--extract.reduceDimension.name=depth',
105-
#'--extract.reduceDimension.start=0',
106-
#'--extract.reduceDimension.end=0',
107-
'--output.file='+nc_fimex]
116+
xAxisValues=np.arange(lon_min, lon_max + dlon, dlon),
117+
yAxisValues=np.arange(lat_min, lat_max + dlat, dlat),
118+
selectVariables=["u", "v"],
119+
reduceTime_start=start_times[n].strftime("%Y-%m-%dT%H:%M:%S"),
120+
reduceTime_end=end_times[n].strftime("%Y-%m-%dT%H:%M:%S"),
121+
)
122+
elif self.program == "fimex":
123+
fimex_command = [
124+
"fimex",
125+
"--input.file=" + url,
126+
"--interpolate.method=bilinear",
127+
"--interpolate.projString=+proj=latlong +ellps=sphere +a=6371000 +e=0",
128+
"--interpolate.xAxisValues="
129+
+ str(lon_min)
130+
+ ","
131+
+ str(lon_min + dlon)
132+
+ ",...,"
133+
+ str(lon_max)
134+
+ "",
135+
"--interpolate.yAxisValues="
136+
+ str(lat_min)
137+
+ ","
138+
+ str(lat_min + dlat)
139+
+ ",...,"
140+
+ str(lat_max)
141+
+ "",
142+
"--interpolate.xAxisUnit=degree",
143+
"--interpolate.yAxisUnit=degree",
144+
"--process.rotateVector.all",
145+
"--extract.selectVariables=u",
146+
"--extract.selectVariables=v",
147+
"--extract.reduceTime.start="
148+
+ start_times[n].strftime("%Y-%m-%dT%H:%M:%S"),
149+
"--extract.reduceTime.end="
150+
+ end_times[n].strftime("%Y-%m-%dT%H:%M:%S"),
151+
"--process.rotateVector.direction=latlon",
152+
#'--extract.reduceDimension.name=depth',
153+
#'--extract.reduceDimension.start=0',
154+
#'--extract.reduceDimension.end=0',
155+
"--output.file=" + nc_fimex,
156+
]
108157
call(fimex_command)
109158
ocr_list.append(xr.open_dataset(nc_fimex).squeeze())
110159

111160
oceancurrent_forcing = xr.concat(ocr_list, dim="time")
112161
# Rename X/Y to lon/lat
113-
oceancurrent_forcing = oceancurrent_forcing.rename_dims({'Y': 'lat', 'X': 'lon'})
114-
oceancurrent_forcing = oceancurrent_forcing.rename_vars({'Y': 'lat', 'X': 'lon'})
162+
oceancurrent_forcing = oceancurrent_forcing.rename_dims(
163+
{"Y": "lat", "X": "lon"}
164+
)
165+
oceancurrent_forcing = oceancurrent_forcing.rename_vars(
166+
{"Y": "lat", "X": "lon"}
167+
)
115168
# Select depth = 0 m
116169
oceancurrent_forcing = oceancurrent_forcing.sel(depth=0)
117170

118-
oceancurrent_forcing['u'] = oceancurrent_forcing['u'].fillna(0)
119-
oceancurrent_forcing['v'] = oceancurrent_forcing['v'].fillna(0)
171+
oceancurrent_forcing["u"] = oceancurrent_forcing["u"].fillna(0)
172+
oceancurrent_forcing["v"] = oceancurrent_forcing["v"].fillna(0)
120173

121-
oceancurrent_forcing = oceancurrent_forcing.transpose('time', 'lat', 'lon')
174+
oceancurrent_forcing = oceancurrent_forcing.transpose("time", "lat", "lon")
122175

123176
return oceancurrent_forcing
124177

125178
def get_url(self, time_stamp):
126-
filename = 'NorKyst-800m_ZDEPTHS_his.an.'+time_stamp.strftime('%Y%m%d')+'00.nc'
127-
url = 'https://thredds.met.no/thredds/dodsC/sea/norkyst800mv0_1h/'+ filename
179+
filename = (
180+
"NorKyst-800m_ZDEPTHS_his.an." + time_stamp.strftime("%Y%m%d") + "00.nc"
181+
)
182+
if time_stamp.year < 2018:
183+
184+
url = (
185+
"https://thredds.met.no/thredds/dodsC/sea/norkyst800mv0_1h/" + filename
186+
)
187+
else:
188+
url = (
189+
"https://thredds.met.no/thredds/dodsC/fou-hi/norkyst800m-1h/" + filename
190+
)
128191
return url

examples/example_mdl_swan.py

+16-16
Original file line numberDiff line numberDiff line change
@@ -2,53 +2,52 @@
22
# IMPORT dnora
33
# =============================================================================
44
from dnora import grd, mdl, wnd, bnd, inp, wlv, ocr, run
5+
56
# =============================================================================
67
# DEFINE GRID OBJECT
78
# =============================================================================
89
# Set grid definitions
9-
grid = grd.Grid(lon=(5.35, 5.6), lat=(59.00, 59.17), name='Boknafjorden')
10+
grid = grd.Grid(lon=(5.35, 5.6), lat=(59.00, 59.17), name="Boknafjorden")
1011

1112
# Set spacing and boundary points
1213
grid.set_spacing(dm=250)
1314

1415
# Import topography and mesh it down to the grid definitions
1516
# Options for grd.readers: EMODNET2020, KartverketNo50m, GEBCO2021
16-
topo_reader = grd.read.EMODNET2020(
17-
tile='*5', folder='/home/konstac/bathy')
17+
topo_reader = grd.read.EMODNET2020(tile="*5", folder="/home/konstac/bathy")
1818
grid.import_topo(topo_reader=topo_reader)
1919

2020
# This can be used to get an empty topography for testing
21-
#grid.import_topo(topo_reader=grd.read.EmptyTopo(grid=grid))
21+
# grid.import_topo(topo_reader=grd.read.EmptyTopo(grid=grid))
2222
#
2323
grid.mesh_grid()
2424
#
2525
# Set the boundaries
26-
bnd_set = grd.boundary.EdgesAsBoundary(edges=[ 'W', 'S'], step = 10)
26+
bnd_set = grd.boundary.EdgesAsBoundary(edges=["W", "S"], step=10)
2727
grid.set_boundary(boundary_setter=bnd_set)
2828
#
2929
#
3030
# =============================================================================
3131
# DEFINE MODEL OBJECT (mdl.)
3232
# =============================================================================
3333
# Options for mdl: SWAN_NORA3, SWAN_ERA5, SWAN_WW3_4km, SWAN_WAM4km
34-
model = mdl.SWAN_NORA3(grid, start_time='2020-12-17T12:00',
35-
end_time='2020-12-17T18:00')
34+
model = mdl.SWAN_NORA3(grid, start_time="2020-12-17T12:00", end_time="2020-12-17T18:00")
3635
# =============================================================================
3736
# IMPORT BOUNDARIES AND FORCING
3837
# =============================================================================
3938
# Boundary Spectra
4039
model.import_boundary(write_cache=True, read_cache=False)
41-
#model.import_boundary(bnd.read_metno.NORA3(source='lustre'),write_cache=True, read_cache=False) # for internal
40+
# model.import_boundary(bnd.read_metno.NORA3(source='lustre'),write_cache=True, read_cache=False) # for internal
4241

4342
# Wind Forcing
4443
model.import_forcing(write_cache=True, read_cache=False)
45-
#model.import_forcing(wnd.read_metno.NORA3(source='lustre'),write_cache=True, read_cache=False) # for internal
44+
# model.import_forcing(wnd.read_metno.NORA3(source='lustre'),write_cache=True, read_cache=False) # for internal
4645

4746
# Ocean Current
48-
#model.import_oceancurrent(ocr.read_metno.NorKyst800())
47+
# model.import_oceancurrent(ocr.read_metno.NorKyst800())
4948

5049
# Water Level
51-
#model.import_waterlevel(wlv.read_ec.GTSM_ERA5())
50+
# model.import_waterlevel(wlv.read_ec.GTSM_ERA5())
5251
# =============================================================================
5352
# PLOT GRID, FORCING AND BOUNDARIES
5453
# =============================================================================
@@ -59,11 +58,12 @@
5958
model.export_grid()
6059
model.export_boundary()
6160
model.export_forcing()
62-
#model.export_oceancurrent()
63-
#model.export_waterlevel()
64-
model.write_input_file(input_file_writer=inp.SWAN(
65-
spec_points=[(5.50, 59.16), (5.55, 59.15)]))
61+
# model.export_oceancurrent()
62+
# model.export_waterlevel()
63+
model.write_input_file(
64+
input_file_writer=inp.SWAN(spec_points=[(5.50, 59.16), (5.55, 59.15)])
65+
)
6666
# =============================================================================
6767
# SWAN RUN
6868
# =============================================================================
69-
model.run_model(model_executer = run.SWAN(nproc=15))
69+
model.run_model(model_executer=run.SWAN(nproc=15))

0 commit comments

Comments
 (0)