seareport_meshes's Issues
Generate station.in per mesh
For each mesh we should also create a corresponding station.in
.
The following snippet is merging STOFS-2D and IOC stations in a single file. The results is similar to this:
$ head station.in
1 0 0 0 0 0 0 0 0 ! https://schism-dev.github.io/schism/master/input-output/optional-inputs.html#stationin-bp-format
3323 ! number of stations
1 144.2806072235 44.0247266972 0 ! abas 2269385 144.2900000000 44.0200000000 13.0000000000 916.6594922518
2 -2.0464618098 57.1249802456 0 ! abed 962967 -2.0800000000 57.1400000000 24.0000000000 2623.9968649481
3 131.4139218534 31.5887292407 0 ! abur 1854768 131.4100000000 31.5800000000 26.0000000000 1039.3086788890
4 -89.8377180390 13.5701359897 0 ! acaj 2746534 -89.8381280000 13.5737920000 9.0000000000 408.9377541903
5 -99.9228904646 16.8436445908 0 ! acap 2681197 -99.9166000000 16.8333000000 15.0000000000 1330.9068274144
6 -99.9228904646 16.8436445908 0 ! acap2 2681197 -99.9030000000 16.8379333000 15.0000000000 2210.0747947232
7 -74.4102481571 39.3729424866 0 ! acnj 1727855 -74.4183000000 39.3550000000 -0.0000000000 2111.7822884706
8 -74.4102481571 39.3729424866 0 ! acnj2 1727855 -74.4183000000 39.3550000000 -0.0000000000 2111.7822884706
$ tail station.in
3314 -133.2311735220 54.4803069628 0 ! "UJ845 SOUS00 SA845" 2077033 -133.2200000000 54.4600000000 345.0000000000 2370.6561082907
3315 -56.6849517955 54.4792339399 0 ! "UJ846 SOUS00 SA846" 2014221 -56.6800000000 54.4600000000 103.0000000000 2162.5208329908
3316 -36.8209103440 54.4617331441 0 ! "UJ847 SOUS00 SA847" 2527088 -36.8400000000 54.4600000000 2839.0000000000 1248.7824430343
3317 -17.1246075148 54.5407370460 0 ! "UJ848 SOUS00 SA848" 492125 -17.0000000000 54.4600000000 2442.0000000000 12055.4642648129
3318 147.4342946221 58.9668723570 0 ! "UJ849 SOUS00 SA849" 2908985 147.4100000000 58.9600000000 111.0000000000 1588.6814030428
3319 -178.5729671215 58.9347728367 0 ! "UJ850 SOUS00 SA850" 1232699 -178.5800000000 58.9600000000 2763.0000000000 2833.9881290880
3320 -144.4913277716 58.9552577946 0 ! "UJ851 SOUS00 SA851" 474508 -144.5700000000 58.9600000000 3772.0000000000 4541.7870224183
3321 -68.0504057985 58.9667668588 0 ! "UJ852 SOUS00 SA852" 1587142 -68.0200000000 58.9600000000 48.0000000000 1898.6426578148
3322 -50.7457083278 58.9051883153 0 ! "UJ853 SOUS00 SA853" 2699858 -51.0100000000 58.9600000000 3433.0000000000 16344.3736645492
3323 -25.5662962363 58.9674260336 0 ! "UJ854 SOUS00 SA854" 2101694 -25.5000000000 58.9600000000 2517.0000000000 3889.4252078064
import os
import geopandas as gpd
import holoviews as hv
import hvplot.pandas
import numpy as np
import pandas as pd
import panel as pn
import searvey
import sklearn.neighbors
from pyposeidon.utils.cfl import parse_hgrid_nodes
from pyposeidon.utils.cpoint import find_nearest_nodes
from pyposeidon.utils.obs import serialize_stations
EARTH_RADIUS = 6371000 # in meters
@pn.cache
def get_stofs2d_meta():
stofs2d = pd.read_csv(
"https://polar.ncep.noaa.gov/stofs/data/stofs_2d_glo_elev_stat_v2_1_0",
names=["coords", "name"],
sep="!",
header=None,
skiprows=1
)
stofs2d = stofs2d.assign(
lon=stofs2d.coords.str.split("\t", n=1).str[0].astype(float),
lat=stofs2d.coords.str.strip().str.rsplit("\t", n=1).str[1].astype(float),
stofs2d_name=stofs2d.name.str.strip(),
).drop(columns=["coords", "name"])
return stofs2d
@pn.cache
def get_ioc_meta() -> gpd.GeoDataFrame:
meta_web = searvey.get_ioc_stations().drop(columns=["lon", "lat"])
meta_api = (
pd.read_json(
"http://www.ioc-sealevelmonitoring.org/service.php?query=stationlist&showall=all"
)
.drop_duplicates()
.drop(columns=["lon", "lat"])
.rename(columns={"Code": "ioc_code", "Lon": "lon", "Lat": "lat"})
)
merged = pd.merge(
meta_web,
meta_api[["ioc_code", "lon", "lat"]].drop_duplicates(),
on=["ioc_code"],
)
return merged
def find_closest_nodes(mesh_nodes: pd.DataFrame, points: pd.DataFrame):
# The only requirement is that both `mesh_nodes and `points` have `lon/lat` columns
tree = sklearn.neighbors.BallTree(
np.radians(mesh_nodes[["lat", "lon"]]).values,
metric="haversine",
)
distances, indices = tree.query(np.radians(points[["lat", "lon"]].values))
closest_nodes = (
mesh_nodes
.rename(columns={"lon": "mesh_lon", "lat": "mesh_lat"})
.iloc[indices.flatten()]
.assign(distance=(distances.flatten() * EARTH_RADIUS))
.reset_index(names=["mesh_index"])
)
return pd.concat((points, closest_nodes), axis="columns")
def merge_ioc_and_stofs(ioc: pd.DataFrame, stofs2d: pd.DataFrame) -> pd.DataFrame:
stations = pd.concat((ioc, stofs2d), ignore_index=True)
stations = stations.assign(unique_id=stations.ioc_code.combine_first(stations.stofs2d_name))
return stations
mesh_nodes = parse_hgrid_nodes("/home/panos/Prog/poseidon/seareport_meshes/v2.2/schism/hgrid.gr3")
ioc = get_ioc_meta()
stofs2d = get_stofs2d_meta()
stations = merge_ioc_and_stofs(ioc=ioc, stofs2d=stofs2d)
closest_nodes = find_nearest_nodes(mesh_nodes=mesh_nodes, points=stations)
serialize_stations(closest_nodes, path="./station.in")
Some of the functions live in the dev
branch of pyposeidon but we could easily inline them.
Remove hardcoded paths
There are hardcoded paths:
$ ag 'home\/tomsail'
mesh_v0.1.py
55: wind = glob.glob("/home/tomsail/Documents/work/python/pyPoseidon/Tutorial/data/era5_202307_uvp.nc")
62: fdem = "/home/tomsail/Documents/work/python/pyPoseidon/Tutorial/data/ETOPO_0.03.nc"
mesh_v0.py
56: wind = glob.glob("/home/tomsail/Documents/work/python/pyPoseidon/Tutorial/data/era5_202307_uvp.nc")
63: fdem = "/home/tomsail/Documents/work/python/pyPoseidon/Tutorial/data/ETOPO_0.03.nc"
Generate different mesh files for both SCHISM and Telemac
Telemac requires to remove the elements that have a node on the north pole. To simplify usage, we should generate two different mesh files, one for each solver.
generate to2d.npy, too
I am not 100% sure, but I think that the to2d.npy should probably be generated alongside the hgrid.gr3 file. Or at the very least there should be an option to generate it. @brey thoughts?
Create a CI that runs every week / month
This will test the capacity of the package to produce reproducible models, including the following steps:
- meshing
- modeling
- post-processing
For that a link with the data is necessary:
- for the coastlines (I am using the parquet version of GSHHS
f
resolution by default) - for the bathy (for now I am using ETOPO)
- for the wind forcing (ERA5 as we cannot distribute ECMWF operational data)
Bathymetry in `gr3` files needs to be positive
The bathymetry in the gr3
files need to be positive numbers. This means that:
- the
v0.1
needs to be multiplied by -1. - the check in
v0.2
needs to be updated toz[z < 20] = 20
instead of:
https://github.com/seareport/seareport_meshes/blob/b5708949fad3d8f9e7e25c01f3c03f1a9a6e05c6/scripts/mesh_v0.2.py#L48C1-L48C21
v0
does not have this issue.
Add a lock file for conda
The yaml file is not enough to reliably reproduce an environment. We should use either conda-lock
or at least the output of conda env export --explicit
Add a poetry.lock/requirements.txt file for each mesh
As it is there are no lock files whatsoever. Furthermore with the current structure whenever we update one of the dependencies, the scripts might generate a different mesh. This is definitely not what we want. Each mesh version should have its own lock file.
`scripts/mesh_v2.2.py` is missing.
Rename `v0` to `v0.0`
v2.1 is missing
`gshhg_{resolution}_nosea.parquet`
There is no reference to this "nosea" coastlines file. Should we add it in https://github.com/seareport/coastlines ?
Recommend Projects
-
React
A declarative, efficient, and flexible JavaScript library for building user interfaces.
-
Vue.js
๐ Vue.js is a progressive, incrementally-adoptable JavaScript framework for building UI on the web.
-
Typescript
TypeScript is a superset of JavaScript that compiles to clean JavaScript output.
-
TensorFlow
An Open Source Machine Learning Framework for Everyone
-
Django
The Web framework for perfectionists with deadlines.
-
Laravel
A PHP framework for web artisans
-
D3
Bring data to life with SVG, Canvas and HTML. ๐๐๐
-
Recommend Topics
-
javascript
JavaScript (JS) is a lightweight interpreted programming language with first-class functions.
-
web
Some thing interesting about web. New door for the world.
-
server
A server is a program made to process requests and deliver data to clients.
-
Machine learning
Machine learning is a way of modeling and interpreting data that allows a piece of software to respond intelligently.
-
Visualization
Some thing interesting about visualization, use data art
-
Game
Some thing interesting about game, make everyone happy.
Recommend Org
-
Facebook
We are working to build community through open source technology. NB: members must have two-factor auth.
-
Microsoft
Open source projects and samples from Microsoft.
-
Google
Google โค๏ธ Open Source for everyone.
-
Alibaba
Alibaba Open Source for everyone
-
D3
Data-Driven Documents codes.
-
Tencent
China tencent open source team.