schism-dev / pyschism

Python interface for handling the SCHISM model.
https://schism-dev.github.io/schism/master/getting-started/pre-processing-with-pyschism/overview.html
Apache License 2.0
24 stars 20 forks source link

Issue Using ERA5 Forcing #30

Open SorooshMani-NOAA opened 2 years ago

SorooshMani-NOAA commented 2 years ago

When I try to use ERA5 as sflux1 I get an error in the nws.nws2.nws2.NWS2.write method. Looking at the code it seems that ERA5 is not following the same API as other sflux subclasses for write method. For example GFS accepts a level argument while ERA5 does not. This results in ERA5's write to get 1 as positional argument for start_date and causes issue due to passing start_date as keyword arg.

I'm using https://github.com/schism-dev/pyschism/commit/ec00dffb251e49759a1162eff4e911c82a7d0cf2, but it seems that tag v0.1.0 (https://github.com/schism-dev/pyschism/commit/e04ee75b1b0bdb9a895f02dcd9aadf9a4722f8c7) has the same issue.

cuill commented 2 years ago

ERA5 was not added to NWS2 class, which writes data twice, first writing to tmpdir then rewriting as sflux_*.nc. Instead, ERA5 only writes once. The following script will generate sflux_1 from ERA5.

from datetime import datetime from time import time import pathlib

from pyschism.forcing.nws.nws2.era5 import ERA5

from pyschism.mesh.hgrid import Hgrid

startdate=datetime(2017, 4, 26) rnday=180

t0=time() hgrid=Hgrid.open('./hgrid.gr3',crs='EPSG:4326') bbox = hgrid.get_bbox('EPSG:4326', output_type='bbox')

er=ERA5() outdir = pathlib.Path('./') er.write(outdir=outdir, start_date=startdate, rnday=rnday, air=True, rad=True, prc=True, bbox=bbox, overwrite=True)print(f'It took {(time()-t0)/60} minutes to generate {rnday} days')

SorooshMani-NOAA commented 2 years ago

I see. I can update my client script to use ERA5 directly, but is there any plan for updating ERA5 to work like GFS and HRRR?

SorooshMani-NOAA commented 2 years ago

@cuill, writing ERA5 alone doesn't solve the issue completely. I also need to manually modify model_driver_obj.param.opt.wtiminc and model_driver_obj.param.opt.nws since it won't be automatically updated based on wind forcing in ModelDriver.__init__ at https://github.com/schism-dev/pyschism/blob/ec00dffb251e49759a1162eff4e911c82a7d0cf2/pyschism/driver.py#L348-L355

So I think it'd be good idea to make the API similar to what NWS2 superclass expects. @jreniel can you please comment on this?

jreniel commented 2 years ago

This is what happens when we prioritize "make it fast" over "make it right". I've been aware of this issue for a long time and will eventually fix it. If you want to provide a fix, you are welcome to open a PR. Thanks.

SorooshMani-NOAA commented 2 years ago

I need to first familiarize myself with the nws module in pyschism. If I get to it I'll open a PR.

SorooshMani-NOAA commented 2 years ago

This might be obvious, but I'm just writing it for the sake of documentation:

I noticed that since ERA5 is bypassing NWS2, other than updating relevant parameters, the user needs to write out sflux_inputs.txt and windrot_geo2proj.gr3 manually as well. In my case I'm just writing an empty namelist file for sflux_inputs.txt and also creating and writing Windrot to the SCHISM input dir (following the logic from NWS2.write).

huyquangtranaus commented 2 years ago

Hi all, have this been solved yet? I think the issue related topyschism.forcing.nws.nws2.sflux because when I made a modified version of your ERA5 code without callingSfluxDataset, and everything seems to work perfectly (for my purpose). What I did is to callinventory = ERA5DataInventory(self.start_date, self.rnday, bbox)directly then using put_sflux_fields to write sflux files.

There is a minor concern: While https://cds.climate.copernicus.eu/api/v2/resources/reanalysis-era5-single-levels) provides hourly data, I noticed that the "sfluxair*.nc files are 3-hourly data.. It's a bit surprising to me, then I have a closer look at put_sflux_fields (/forcing/nws/nws2/era5.py) https://github.com/schism-dev/pyschism/blob/d2b474e8a32e22eea28f7da462feae8424ab0ca4/pyschism/forcing/nws/nws2/era5.py#L95 and it seems to be correct. Can you please confirm this?

Thank you :)

cuill commented 2 years ago

@huyquangtranaus Yes, it generates 3-hourly netCDF files. Originally, it was designed to satisfy our group's needs. I think it's better to add the output interval as a parameter to let the user decide. I will consider adding this feature. Thanks.

huyquangtranaus commented 2 years ago

@cuill Yeah... that is useful too. Thank you.

Also, there is likely an issue in the subset forcing area is in the southern hemisphere (SH) because you use ::-1 to reverse the list, for [example] (https://github.com/schism-dev/pyschism/blob/d2b474e8a32e22eea28f7da462feae8424ab0ca4/pyschism/forcing/nws/nws2/era5.py#L131)

So the 2d matrix is flipped like this:
Screenshot from 2022-05-27 20-19-52

Have you seen this?

cuill commented 2 years ago

@huyquangtranaus The old version (before Jan 31, 2022) of schism only allowed the counter-clockwise direction in the sflux horizontal grid. So the reversed list flips the variable in the latitude direction: dev/pyschism/blob/df803edb53184625b12399f38a8bd26a022abbc1/pyschism/forcing/nws/nws2/era5.py#L78

However, since this restriction has been removed in the latest schism (https://github.com/schism-dev/schism/commit/8ecf95efc028136f1358351ce223a32e589c944f), we may consider removing the reversing.