Extraction of ERA 5 data for custom grid

[neelbasak08]

Hello people,

I wanted to use ERA-5 dataset to create an extraction data that can be used with ROMS . But I was unable to find any script that could extract data from ERA-5 and generate an output netcdf file either in MATLAB format . Scripts were available for ERA interim but I guess obsolete now.

I tried my luck searching the internet but couldn't find any.

Is there any MATLAB script available that could do the job ?

Any help in any form is highly appreciated.

Kind regards
-Neel

[wilkin]

I have MATLAB code roms_get_era5_NCARds633_bulkflux.m that reads ECMWF ERA-5 meteorological reanalysis from the NCAR Research Data Archive (RDA) dataset ds633.0 https://rda.ucar.edu/datasets/ds633.0. Access to this dataset at RDA requires authentication with a login, so you need to be a registered user.

The NCAR dataset is accessed via OPeNDAP, so my function is designed to extract a regional subset instead of having to download entire global files and then subset them. The output from the script is passed to another function roms_write_era5_NCARds633_frcfile.m to write the forcing files.

I have only set this up to make BULK_FLUXES forcing files. It would take some work to extend to stresses and net heat fluxes. It's on my long to-do list.

Both these functions are in my roms_wilkin repo on github https://github.com/johnwilkin/roms_wilkin

[neelbasak08]

Thanks a lot Mr. Wilkin.

Although I faced an error(I am using MATLAB 2019b):

roms_get_era5_NCARds633_bulkflux
Using username:password for rda.ucar.edu from $HOME/.netrc

bbox =

-1.92 -2.24 49.75 49.65


Processing variable mean sea-level pressure
Error using ncread
Expected count to be positive.

Error in internal.matlab.imagesci.nc/read (line 562)
validateattributes(count,{'numeric'},...

Error in ncread (line 66)
vardata = ncObj.read(varName, varargin{:});

Error in roms_get_era5_NCARds633_bulkflux (line 295)
lon = ncread(url,'longitude',Is,Ilen)-360;

I read through you code and found that you had mentioned the following as well:

Code: Select all

 % NOTE on Matlab versions:
%   This code runs successfully with Matlab Release 2020b (9.9.0.1467703).
%   Some Matlab versions throw an error in ncread when reading flux 
%   variables (longwave radiation etc.). Circumstantial evidence links this
%   to ncread built with NetCDF 4.6.1. You can check the version of NetCDF 
%   that your Matlab is using with: >> libvers = netcdf.inqLibVers 

When I changed the bbox to the default West Atlantic model one bbox = [-110 -30 0 55] i got no error


Is there a way to solve this error and make it work in my MATLAB 2019b or I need to download MATLAB 2020b version for this?
Is it a problem with my bbox variable which I assigned bbox=[-2.246 -1.92 49.65 49.75] ([long_left long_right lat_bottom lat_top])?

Also, when I ran the script on MATLAb online 2020b the default West atlantic model it ran successfully but when I ran the roms_write_era5_NCARds633_frcfile.m it showed the following error:

Unable to resolve the name E.time.data.

Error in roms_write_era5_NCARds633_frcfile (line 30)
time = E.time.data - Time0;

I had the following in roms_write_era5_NCARds633_frcfile:

Code: Select all

My Time0 = datenum(2011,1,1);

in my roms_get_era5_ncar_ds633.m I had the following :

Code: Select all

  epoch = datenum(1900,1,1);
yyyy = 2018;
mm = 05;
bbox = [-110 -30 0 55];
YYYYMM = [sprintf('%d',yyyy) sprintf('%02d',mm);];
YYYY = sprintf('%d',yyyy); 

So if I am not wrong, after extraction file is created , I need to define bulk_fluxes header in my roms application header file to get it running with my grid in ROMS model. Am I right ?


Extreme regards
-Neel

[neelbasak08]

I have MATLAB code roms_get_era5_NCARds633_bulkflux.m that reads ECMWF ERA-5 meteorological reanalysis from the NCAR Research Data Archive (RDA) dataset ds633.0 https://rda.ucar.edu/datasets/ds633.0. Access to this dataset at RDA requires authentication with a login, so you need to be a registered user.

The NCAR dataset is accessed via OPeNDAP, so my function is designed to extract a regional subset instead of having to download entire global files and then subset them. The output from the script is passed to another function roms_write_era5_NCARds633_frcfile.m to write the forcing files.

I have only set this up to make BULK_FLUXES forcing files. It would take some work to extend to stresses and net heat fluxes. It's on my long to-do list.

Both these functions are in my roms_wilkin repo on github https://github.com/johnwilkin/roms_wilkin

Dear Dr. Wilkin, Thanks a lot for the script you graciously provided me with. I tried to run the roms_get_era5_NCARds633_bulkflux.m script for default West Atlantic Model case as in the script in MATLAB 2020b as you had stated previously and it ran error free. But as soon as I try to execute roms_write_era5_NCARds633_frcfile.m I still the following error:

Unable to resolve the name E.time.data.

Error in roms_write_era5_NCARds633_frcfile (line 30)
time = E.time.data - Time0;

I had read your script cautiously and it stated the following :

Code: Select all

 NOTE: This is an example script. The user needs to configure the ROMS
% time coordinate basedate (Time0) and the output filename (ncname) and
% title and any other pertinent metadata.
%
% This script uses routines from the myroms.org Matlab tools
% (roms_metadata, nc_create, nc_write, etc.). 

I use MATLAB 2020b and it still shows this error. I am unable to understand what is meant by settimg Time0 time cordinate basedate. I want my simulation from August 1, 2018 for 30 days. I had set the following parameter in my roms_write_era5_NCARds633_frcfile.m:

Code: Select all

 Time0 = datenum(2018,5,1);
time = E.time.data - Time0;
Tname = 'time'; % need to reset some parameters from roms_metadata
 

And the following in my roms_get_era5_NCARds633_bulkflux.m:

Code: Select all

 epoch = datenum(1900,1,1);
yyyy = 2018;
mm = 05;
bbox = [-110 -30 0 55];
YYYYMM = [sprintf('%d',yyyy) sprintf('%02d',mm);];
YYYY = sprintf('%d',yyyy); 

Am I doing something incorrect or am I at fault at understanding something ? Why is E.time.data showing error when roms_get_era5_NCARds633_bulkflux.m ran successfully? I am not able to figure it out.

Any help in any kind is extremely appreciated by me.

Kind regards
-Neel

[wilkin]

In the preamble comments for these two scripts it says:

In roms_write_era5_NCARds633_frcfile.m:

% First run e.g. E = roms_get_era5_ncar_ds633(year,month,bounding_box...)
% and then run this script to create the forcing netcdf file.

In roms_get_era5_NCARds633_bulkflux.m:

% Use this function to extract the data for a chosen month, then use
% roms_write_era5_NCARds633_frcfile to create and write the ROMS format
% surface forcing netcdf file.

So, first run E = roms_get_era5_NCARds633_bulkflux ... for the year, month and bounding box you want, and then run roms_write_era5_NCARds633_frcfile, which will process the data in structure E.

I did this as shown below:

Code: Select all

>> bbox = [-80 -67 33 42]
>> E = roms_get_era5_NCARds633_bulkflux(2018,5,bbox)
Using username:password for rda.ucar.edu from $HOME/.netrc
 
Processing variable mean sea-level pressure
 Reading 1 file ...
  in 2 16-day chunks:   Elapsed time is 6.399692 seconds.
 
Processing variable 2 meter temperature
 Reading 1 file ...
  in 2 16-day chunks:   Elapsed time is 7.523408 seconds.
 
Processing variable 2 meter dewpoint temperature
 Reading 1 file ...
  in 2 16-day chunks:   Elapsed time is 6.080859 seconds.
 
Processing variable 10-m wind (u)
 Reading 1 file ...

and then ...

Code: Select all

>> roms_write_era5_NCARds633_frcfile

Creating ROMS NetCDF forcing file: frc_MY_APPLICATION_ERA5_bulkflux_201805.nc

Wrote lon                 Min=-8.00000e+01 Max=-6.72500e+01
Wrote lat                 Min= 3.32500e+01 Max= 4.20000e+01
Wrote Uwind               into record: , Min=-1.22101e+01 Max= 1.47395e+01
Wrote Vwind               into record: , Min=-1.17272e+01 Max= 2.09668e+01
Wrote Pair                into record: , Min= 1.00461e+03 Max= 1.03120e+03

[neelbasak08]

In the preamble comments for these two scripts it says:

In roms_write_era5_NCARds633_frcfile.m:

% First run e.g. E = roms_get_era5_ncar_ds633(year,month,bounding_box...)
% and then run this script to create the forcing netcdf file.

In roms_get_era5_NCARds633_bulkflux.m:

% Use this function to extract the data for a chosen month, then use
% roms_write_era5_NCARds633_frcfile to create and write the ROMS format
% surface forcing netcdf file.

So, first run E = roms_get_era5_NCARds633_bulkflux ... for the year, month and bounding box you want, and then run roms_write_era5_NCARds633_frcfile, which will process the data in structure E.

I did this as shown below:

Code: Select all

>> bbox = [-80 -67 33 42]
>> E = roms_get_era5_NCARds633_bulkflux(2018,5,bbox)
Using username:password for rda.ucar.edu from $HOME/.netrc
 
Processing variable mean sea-level pressure
 Reading 1 file ...
  in 2 16-day chunks:   Elapsed time is 6.399692 seconds.
 
Processing variable 2 meter temperature
 Reading 1 file ...
  in 2 16-day chunks:   Elapsed time is 7.523408 seconds.
 
Processing variable 2 meter dewpoint temperature
 Reading 1 file ...
  in 2 16-day chunks:   Elapsed time is 6.080859 seconds.
 
Processing variable 10-m wind (u)
 Reading 1 file ...

and then ...

Code: Select all

>> roms_write_era5_NCARds633_frcfile

Creating ROMS NetCDF forcing file: frc_MY_APPLICATION_ERA5_bulkflux_201805.nc

Wrote lon                 Min=-8.00000e+01 Max=-6.72500e+01
Wrote lat                 Min= 3.32500e+01 Max= 4.20000e+01
Wrote Uwind               into record: , Min=-1.22101e+01 Max= 1.47395e+01
Wrote Vwind               into record: , Min=-1.17272e+01 Max= 2.09668e+01
Wrote Pair                into record: , Min= 1.00461e+03 Max= 1.03120e+03

Thanks Dr. Wilkin for your kindest reply. I ran the the scripts in a similar manner to yours but faced an error for some error:

Code: Select all

   E = roms_get_era5_NCARds633_bulkflux(2018,5,bbox)
Using username:password for rda.ucar.edu from $HOME/.netrc
 
Processing variable mean sea-level pressure
 Reading 1 file ...
  in 2 16-day chunks:   Elapsed time is 18.282090 seconds.
 
Processing variable 2 meter temperature
 Reading 1 file ...
  in 2 16-day chunks:   Elapsed time is 19.090334 seconds.
 
Processing variable 2 meter dewpoint temperature
 Reading 1 file ...
  in 2 16-day chunks:   Elapsed time is 14.885951 seconds.
 
Processing variable 10-m wind (u)
 Reading 1 file ...
  in 2 16-day chunks:   Elapsed time is 14.918444 seconds.
 
Processing variable 10-m wind (v)
 Reading 1 file ...
  in 2 16-day chunks:   Elapsed time is 23.771228 seconds.
 
Processing variable mean surface downward long-wave radiation flux
 Reading file 1 of 3 ...  Elapsed time is 4.213215 seconds.
  Reading file 2 of 3 ...  Elapsed time is 22.879712 seconds.
   Reading file 3 of 3 ...  Elapsed time is 22.151422 seconds.
 
Processing variable mean surface net long-wave radiation flux
 Reading file 1 of 3 ...  Elapsed time is 4.787892 seconds.
  Reading file 2 of 3 ...  Elapsed time is 33.216153 seconds.
   Reading file 3 of 3 ...  Elapsed time is 28.094873 seconds.
 
Processing variable mean surface net short-wave radiation flux
 Reading file 1 of 3 ...  Elapsed time is 4.530642 seconds.
  Reading file 2 of 3 ...  Elapsed time is 47.138120 seconds.
   Reading file 3 of 3 ...  Elapsed time is 17.842173 seconds.
 
Processing variable mean total precipitation rate
 Reading file 1 of 3 ...  Elapsed time is 4.505570 seconds.
  Reading file 2 of 3 ...  Elapsed time is 20.921956 seconds.
   Reading file 3 of 3 ...  Elapsed time is 20.904085 seconds.

E = 

  struct with fields:

            u10: [1×1 struct]
            v10: [1×1 struct]
             d2: [1×1 struct]
             t2: [1×1 struct]
            msl: [1×1 struct]
          msshf: [1×1 struct]
          mslhf: [1×1 struct]
       msdwlwrf: [1×1 struct]
        msnlwrf: [1×1 struct]
        msnswrf: [1×1 struct]
           mtpr: [1×1 struct]
           time: [1×1 struct]
            lon: [1×1 struct]
            lat: [1×1 struct]
           yyyy: 2018
             mm: 5
    description: 'https://rda.ucar.edu/datasets/ds633.0'
       citation: 'European Centre for Medium-Range Weather Forecasts, 2019,updated monthly. ERA5 Reanalysis (0.25 Degree Latitude-Longitude Grid). Research Data Archive at the National Center for Atmospheric Research, Computational and Information Systems Laboratory. https://doi.org/10.5065/BH6N-5N20. Accessed 29-Jan-2021 20:40:44'
 

But here lies the error:

Code: Select all

 >> roms_write_era5_NCARds633_frcfile

Creating ROMS NetCDF forcing file: frc_watl_ERA5_bulkflux_201805.nc

Unrecognized function or variable 'nc_constant'.

Error in roms_write_era5_NCARds633_frcfile (line 122)
S.Dimensions(3).Length = nc_constant('nc_unlimited');
 
>>  

I checked there is nothing wrong with my netcdf version and I also have snctools and mexcdf added to my matlab paths. What might be causing this issue in your opinion ? Is there a remedy to this ?

Any help is deeply appreciated by me.

Kind regards
-Neel

[ymamoutos]

Greetings

nc_constant.m is part of ROMS matlab toolbox.
You can download it with svn using your ROMS credentials
from here:

https://www.myroms.org/svn/src/matlab/

Giannis

[wilkin]

Script roms_write_era5_NCARds633_frcfile.m notes that ...

Code: Select all

% This script uses routines from the myroms.org Matlab tools
% (roms_metadata, nc_create, nc_write, etc.).

nc_constant is part of that suite of MATLAB codes available from myroms.org.

BTW, when you post here there is no need to quote a previous post in its entirety ... we can see it in the thread.

[neelbasak08]

Thanks a lot Dr. Wilkin and Mr. Ymamoutos for your generous help.

The script is working fine for me. Sorry from my end for being not so good while posting in ROMS discussion. I will keep that advice.

Another query I had , as per the script instruction I shall need to define only BULKFLUX in my header file of application right ? Do I need to parallelly activate any other cppdefs option in order to use the output file ?

Again , lots and lots of thanks and gratitude for your generous help and time.

Kind regards
-Neel

[wilkin]

#define BULK_FLUXES is basically what you need.

You have other choices though.

#define LONGWAVE_OUT if you wish ROMS to compute upward longwave radiation from the model SST. ROMS will read lwrad_down from the ERA-5 file. The default if this is not defined is to use lwrad (net) from the ERA-5 file.

#define EMINUSP if you wish ROMS to compute evaporation from latent heat and combine with precipitation from ERA-5 to get the freshwater flux.

Do NOT #define DIRUNAL_SRFLUX because these ERA-5 data are hourly and resolve the diurnal cycle of shortwave radiation (swrad).

[neelbasak08]

Sorry to bother you once again Mr. Wilkin but I had a small issue which I noticed while extraction of ERA5 data using the script which you generously provided.

I am using MATLAB 2020b software. When I am passing on the argument yyyy as 2019 in function, the output file generated still has the yyyy as 2018. I haven't changed anything in the script.
Is there any specific reason to that or it is an error ? I'm not sure. Also it's taking a lot of time to process.

I am trying to extract ERA5 data for entire Philippines.

bbox = [114 129 4 19]

bbox =

114 129 4 19

>> E = roms_get_era5_NCARds633_bulkflux(2019,5,bbox)
Using username:password for rda.ucar.edu from $HOME/.netrc

Processing variable mean sea-level pressure
Reading 1 file ...
in 2 16-day chunks: Elapsed time is 96.768261 seconds.

Processing variable 2 meter temperature
Reading 1 file ...
in 2 16-day chunks: Elapsed time is 315.813836 seconds.

Processing variable 2 meter dewpoint temperature
Reading 1 file ...
in 2 16-day chunks: Elapsed time is 238.946152 seconds.

Processing variable 10-m wind (u)
Reading 1 file ...
in 2 16-day chunks: Elapsed time is 458.363336 seconds.

Processing variable 10-m wind (v)
Reading 1 file ...
in 2 16-day chunks: Elapsed time is 170.625708 seconds.

Processing variable mean surface downward long-wave radiation flux
Reading file 1 of 3 ... Elapsed time is 7.369953 seconds.
Reading file 2 of 3 ... Elapsed time is 57.551023 seconds.
Reading file 3 of 3 ... Elapsed time is 70.367142 seconds.

Processing variable mean surface net long-wave radiation flux
Reading file 1 of 3 ... Elapsed time is 8.145141 seconds.
Reading file 2 of 3 ... Elapsed time is 74.051057 seconds.
Reading file 3 of 3 ... Elapsed time is 78.556211 seconds.

Processing variable mean surface net short-wave radiation flux
Reading file 1 of 3 ... Elapsed time is 7.436637 seconds.
Reading file 2 of 3 ... Elapsed time is 50.002450 seconds.
Reading file 3 of 3 ... Elapsed time is 57.693464 seconds.

Processing variable mean total precipitation rate
Reading file 1 of 3 ... Elapsed time is 7.705115 seconds.
Reading file 2 of 3 ... Elapsed time is 127.911613 seconds.
Reading file 3 of 3 ... Elapsed time is 76.806878 seconds.

E =

struct with fields:

u10: [1×1 struct]
v10: [1×1 struct]
d2: [1×1 struct]
t2: [1×1 struct]
msl: [1×1 struct]
msshf: [1×1 struct]
mslhf: [1×1 struct]
msdwlwrf: [1×1 struct]
msnlwrf: [1×1 struct]
msnswrf: [1×1 struct]
mtpr: [1×1 struct]
time: [1×1 struct]
lon: [1×1 struct]
lat: [1×1 struct]
yyyy: 2018
mm: 5
description: 'https://rda.ucar.edu/datasets/ds633.0'
citation: 'European Centre for Medium-Range Weather Forecasts, 2019,updated monthly. ERA5 Reanalysis (0.25 Degree Latitude-Longitude Grid). Research Data Archive at the National Center for Atmospheric Research, Computational and Information Systems Laboratory. https://doi.org/10.5065/BH6N-5N20. Accessed 06-Feb-2021 13:04:55'

>>

Is there any solution to this issue ?

Kind regards
-Neel

[neelbasak08]

Also when I tried to run the script on MATLAB 2020b online , I was thrown the following error :

Cordinates I want to use: Left longitude = 114 Right Longitude= 129 Top Latitude= 19 Bottom Latitude = 4. Date of simulation May, 2018 Duration 30 days.

bbox = [114 129 4 19]

bbox =

114 129 4 19

>> E = roms_get_era5_NCARds633_bulkflux(2018,5,bbox)
Using username:password for rda.ucar.edu from $HOME/.netrc

Processing variable mean sea-level pressure
Error using ncread
Expected start to be nonempty.

Error in internal.matlab.imagesci.nc/read (line 557)
validateattributes(start,{'numeric'},...

Error in ncread (line 66)
vardata = ncObj.read(varName, varargin{:});

Error in roms_get_era5_NCARds633_bulkflux (line 292)
lon = ncread(url,'longitude',Is,Ilen)-360;

I am not sure why it's showing this error as my coordinates are right and in the format as well if I am not wrong.

I tried it using MATLAB 2020b software as well but it's showing up a different issue which I have written just 1 post above the current one.

Is there any way to solve this issue or is it something I am doing wrong ?

Any help is highly appreciated.

Kind regards
-Neel

[wilkin]

Your two questions:

(1) Wrong yyyy in filename:

I have no idea how you got this:

>> E = roms_get_era5_NCARds633_bulkflux(2019,5,bbox)
...
E =
struct with fields:
...
yyyy: 2018
mm: 5

Did you hack something inside the code while you were testing? Use the MATLAB debugger to step through roms_get_era5_NCARds633_bulkflux and see where yyyy get altered. This doesn't happen in my present code.

(2) ncread error

Find the line of code with

Code: Select all

lon = ncread(url,'longitude')-360;

and take out the minus 360.

I work in the West Atlantic region and typically use negative west longitudes in my coordinates, so I wrote the code accordingly. Maybe one day I will update to trap and remedy this, but for now just modify your code to work for you.

[wilkin]

Your two questions:

(1) Wrong yyyy in filename:

I have no idea how you got this:

>> E = roms_get_era5_NCARds633_bulkflux(2019,5,bbox)
...
E =
struct with fields:
...
yyyy: 2018
mm: 5

Did you hack something inside the code while you were testing? Use the MATLAB debugger to step through roms_get_era5_NCARds633_bulkflux and see where yyyy get altered. This doesn't happen in my present code.

(2) ncread error

Find the line of code with

Code: Select all

lon = ncread(url,'longitude')-360;

and take out the minus 360.

I work in the West Atlantic region and typically use negative west longitudes in my coordinates, so I wrote the code accordingly. Maybe one day I will update to trap and remedy this, but for now just modify your code to work for you.

[wilkin]

I updated my code on github to automatically detect if the user requested longitudes as negative or positive in BBOX and adjust accordingly.

[neelbasak08]

Dear Dr. Wilkin

The debug worked fine. I have some issue with running of model based on this era5 data. I was thrown the following error while running ROMS model:

NETCDF_INQ_VAR - requested variable: srf_time
not found in input file: frc_watl_ERA5_bulkflux_alderney201905_new.nc
call from: ROMS/Utility/get_cycle.F

Found Error: 02 Line: 140 Source: ROMS/Utility/get_cycle.F
Found Error: 02 Line: 321 Source: ROMS/Utility/inquiry.F
Found Error: 02 Line: 133 Source: ROMS/Utility/get_2dfld.F
Found Error: 02 Line: 341 Source: ROMS/Nonlinear/get_data.F
Found Error: 02 Line: 659 Source: ROMS/Nonlinear/initial.F
Found Error: 02 Line: 191 Source: ROMS/Drivers/nl_ocean.h, ROMS_initialize

I read on the forum that swrad:time='srf_time' needs to be equated to Unwind:time='time'. But I'm not sure how to do this in your current set of MATLAB scripts. Is there a way to make these changes ?

Kind regards
-Neel

[donoso]

Hello Mr. John,

I hava a conceptual question about roms_get_era5_NCARds633_bulkflux.m script. In that file, why you use forecast fields to 'msdwlwrf','msnlwrf','msnswrf','mtpr' and analysis fields to 'msl','t2','d2','u10','v10'? Wouldn't it be better to just use analysis fields?

Thanks in advance.

Cheers,
David.

I have MATLAB code roms_get_era5_NCARds633_bulkflux.m that reads ECMWF ERA-5 meteorological reanalysis from the NCAR Research Data Archive (RDA) dataset ds633.0 https://rda.ucar.edu/datasets/ds633.0. Access to this dataset at RDA requires authentication with a login, so you need to be a registered user.

The NCAR dataset is accessed via OPeNDAP, so my function is designed to extract a regional subset instead of having to download entire global files and then subset them. The output from the script is passed to another function roms_write_era5_NCARds633_frcfile.m to write the forcing files.

I have only set this up to make BULK_FLUXES forcing files. It would take some work to extend to stresses and net heat fluxes. It's on my long to-do list.

Both these functions are in my roms_wilkin repo on github https://github.com/johnwilkin/roms_wilkin

[wilkin]

Responses to two posts in this thread:

(1) From user neelbasak08 regarding the error:

Code: Select all

NETCDF_INQ_VAR - requested variable: srf_time
not found in input file: frc_watl_ERA5_bulkflux_alderney201905_new.nc
call from: ROMS/Utility/get_cycle.F

Sorry, that was a bug. I grab the basic metadata for each forcing variable using roms_metadata, but that uses

Code: Select all

time:srf_time

etc. which I have to change to be

Code: Select all

time:time

to use the common time coordinate in the file.

I had assumed that S.Variables(*).Attributes(3) was always "time", but for the radiation variables it is S.Variables(*).Attributes(5). I've updated the code to use my function findstrinstruct to detect the correct attribute to change.

(2) From user donoso

I had the same concern about mixing analysis and forecast fields, so I asked the UCAR helpdesk (who are super helpful and responsive, by the way). I posed it this way:

Do you see any severe logistical disconnect in me using some variables from analysis and some from forecast? For my purposes I doubt it.

David Stepaniak who curates this data set for NCAR replied:

Yes, it's OK to mix analysis and forecast parameters — in particular, radiation parameters are derived at least in part from the radiation transfer scheme (RTTOV611) that can only be run as part of the forecast in physical space.

[donoso]

Thank you John! I'm not sure to understand that answer pretty well. But, in my interpretation, the background reason to request different data sources (analysis vs forecast) would be related to store physical limitation at the NCAR. That is correct or have a sense for you?

Cheers,
David.

Responses to two posts in this thread:

(1) From user neelbasak08 regarding the error:

Code: Select all

NETCDF_INQ_VAR - requested variable: srf_time
not found in input file: frc_watl_ERA5_bulkflux_alderney201905_new.nc
call from: ROMS/Utility/get_cycle.F

Sorry, that was a bug. I grab the basic metadata for each forcing variable using roms_metadata, but that uses

Code: Select all

time:srf_time

etc. which I have to change to be

Code: Select all

time:time

to use the common time coordinate in the file.

I had assumed that S.Variables(*).Attributes(3) was always "time", but for the radiation variables it is S.Variables(*).Attributes(5). I've updated the code to use my function findstrinstruct to detect the correct attribute to change.

(2) From user donoso

I had the same concern about mixing analysis and forecast fields, so I asked the UCAR helpdesk (who are super helpful and responsive, by the way). I posed it this way:

Do you see any severe logistical disconnect in me using some variables from analysis and some from forecast? For my purposes I doubt it.

David Stepaniak who curates this data set for NCAR replied:

Yes, it's OK to mix analysis and forecast parameters — in particular, radiation parameters are derived at least in part from the radiation transfer scheme (RTTOV611) that can only be run as part of the forecast in physical space.

[wilkin]

would be related to store physical limitation at the NCAR. That is correct or have a sense for you?

No. My understanding is that the radiation variables are not part of the state space of the analysis. They are computed diagnostically by the radiation transfer scheme (RTTOV611) while the forecast is running dependent on the model cloud solution.

[donoso]

Oh, that sound completely new to me. If you access ERA5 datasets at the Climate Data Store of Copernicus, you will not see the concept of 'analysis' and 'forecast'. You need to request all variables, including radiation variables (fluxes), like 'analysis' fields. However, that doesn't mean all variables available there are truly type ''analysis'. I guess that is the origin of my confusion. Anyway, I need to read more about how ERA5 reanalysis works.

Thank you John!

Cheers,
David.

would be related to store physical limitation at the NCAR. That is correct or have a sense for you?

No. My understanding is that the radiation variables are not part of the state space of the analysis. They are computed diagnostically by the radiation transfer scheme (RTTOV611) while the forecast is running dependent on the model cloud solution.

[stef]

I'm trying to learn how to use atmospheric forcing. Thanks for the matlab scripts on github! They are great.

I saw in [1] that e.g. the mean evaporation rate in the era5 reanalysis is

This parameter is a mean over a particular time period (the processing period) which depends on the data extracted. For the reanalysis, the processing period is over the 1 hour ending at the validity date and time.

Do I understand correctly that validity times are always full hours? Does this mean that we should use 00:30, 01:30, etc. for the time axis of the mean rates and fluxes? When I use grib_ls e.g. on mean evap. rate it says:

Code: Select all

dataDate     stepRange    dataType
20151130     5-6              fc

while for the analyzed fields, e.g. mean sea level pressure

Code: Select all

dataDate     stepRange    dataType 
20151201     0                 an

The date and time specification are documented in Table 0 of [5]. In [4], era5 time terminology is further documented, but I'm not sure what 'stepRange' means. Does 5-6 mean that it starts at forecast hour 5 and ends at forecast hour 6? So it should really be centered in between the forecast hours, so at XX:30?


Another question, which is less important as I think I could use any of the existing formulas, but I'm still curious:
I have not found a reference for the formula of vapor pressure you are using (I browsed through [2] and [3]):

VP = 6.11 .* 10.0 .^ (7.5 .* tdew ./ (237.7 + tdew))

I guess this is in hectopascal? I compared it to the August-Roche-Magnus formula on Wikipedia and it is very similar. Do you happen to have a reference for this formula? Thanks so much for your help!


[1] https://cds.climate.copernicus.eu/cdsap ... b=overview
[2] https://en.wikipedia.org/wiki/Vapour_pr ... n_formulas
[3] http://climate.envsci.rutgers.edu/pdf/L ... ntBAMS.pdf
[4] https://confluence.ecmwf.int/pages/view ... d=85402030
[5] https://confluence.ecmwf.int/display/CK ... cification

[stef]

I just saw in [1] it says

In HRES atmospheric (hourly steps), mean rate precipitation at time=06:00, step=3 represents the average precipitation rate in the 1-hour period from 08:00 to 09:00,

so that's clear now. I'm downloading from "cds.climate.copernicus.eu". Not sure which time axis "rda.ucar" uses.

[1] https://confluence.ecmwf.int/pages/view ... parameters

[wilkin]

I would look at the documentation at NCAR for this specific version of the dataset, not other versions hosted by Copernicus.

NCAR provide some variables as mean fluxes per unit time and some as instantaneous time analysis.

The structure E in my code shows which variables are loaded from the meanflux catalog and which are instantaneous from the surface analysis catalog ... e.g.

Code: Select all

E.mer.long = 'mean evaporation rate'; % opposite sign convention to ROMS
E.mer.name = 'MER';
E.mer.code = '235_043';
E.mer.units = 'kg m-2 s-1';
E.mer.set = 'e5.oper.fc.sfc.meanflux';  % 3 files
E.mer.v = 'mer';

The NCAR mean fluxes is a different format from the time accumulated values that were so cumbersome with ERA Interim at ECMWF (and for which there is out-of-date emcwf2roms Matlab code).

At rad.ucar.edu I went to https://rda.ucar.edu/datasets/ds633.0/# ... html?_do=y and clicked on the entry for mean evaporation rate kg m-2 s-1 which pops up a box that indicates the data are 1-hour Average (initial+0 to initial+1), etc.

The code I wrote computes the time as forecast_initial_time plus forecast_hour, so if the data are indeed the 1-hour average then -- to answer your question -- one should probably add (or is it subtract?) 30 minutes to the time for any of the variables downloaded from the meanflux catalog. I'd need to dig a bit deeper to absolutely convince myself whether it is add or subtract 30 minutes. For 0.25 degree data, do you think 30 mins can possibly matter?

Regarding your question about the calculation of relative humidity, the equation is derived from the Clausius-Clapeyron equation. I forget the precise provenance of this formulation in terms of powers of 10. Most versions of this are written in terms of exponents of e (natural logarithm). For example, see ANNEX 4.B. FORMULAE FOR THE COMPUTATION OF MEASURES OF HUMIDITY in WMO publication Guide to Instruments and Methods of Observation, Volume I: Measurement of Meteorological Variables, WMO No. 8, 2018 https://library.wmo.int/index.php?lvl=n ... iGSwf5KiUk. There you see equations for saturation water vapor in the form ew(t) = 6.11 exp (...). Using this formulation, relative humidity is 100*ew(t_dew point)/ew(t_2m). You can verify that my code in the format 10^() and WMO equations in the format exp() give the same result.

[stef]

Thanks for the WMO paper, I should read it! These are great.

Yes we should look at the rda.ucar documentation. Unfortunately I did not realize that rda.ucar offers the ERA5 data set until it was too late. I like to download the netcdfs via thredds/opendap on rda.ucar. But I had already invested some time in the copernicus version before I saw your scripts.

Regarding the question whether 30 mins matter for a 0.25 degree data set: I don't know. But it's an hourly dataset, so I thought I might as well get a 30min phase shift right.

[wilkin]

I have updated my routines roms_get_era5_NCARds633_bulkflux.m and roms_write_era5_NCARds633_frcfile.m at
https://github.com/johnwilkin/roms_wilkin

Relative humidity is now computed using the formula documented by WMO-No.8 report (2018). In my test comparison the values differed at most by a factor of 1.008.

I have also updated the code (I did this a while ago) to allow the option of building ROMS forcing files for either #define BULK_FLUXES or #undef BULK_FLUXES. In the latter case, ROMS will look for stresses (surface momentum flux), net heat flux, net shortwave (if #define SOLAR_SOURCE) and net freshwater flux "swflux"(evaporation minus precipitation).

There is the option 'allfluxes' to build forcing files that are complete for both options, so that the user can switch between them without changing the FRCFILE list in roms.in.

However, something we have noticed in tests is that near the coast, in regions where the ERA5 land/sea mask is a fractional value and not exactly 0 (which is strictly water), that the stresses can be very large. Much larger than typical values of sea water (like 0.3 Pa in Long Island Sound). I'm not meteorologist, but I assume this indicates that over land the drag coefficient is much larger than over the ocean (seems plausible). This difference is sufficient to cause dramatic excessive coastal upwelling. So, I would recommend users stick with BULK_FLUXES in conjunction with ERA5 if they have a coastal ocean in the domain.

[stef]

Fantastic! Thanks! I'm sure a lot of us appreciate the work you put into this. It seems to be a very widely used data set, and having the options for direct momenum/heat/salt fluxes and BULK_FLUXES is great.

However, something we have noticed in tests is that near the coast, in regions where the ERA5 land/sea mask is a fractional value and not exactly 0 (which is strictly water), that the stresses can be very large. Much larger than typical values of sea water (like 0.3 Pa in Long Island Sound). I'm not meteorologist, but I assume this indicates that over land the drag coefficient is much larger than over the ocean (seems plausible). This difference is sufficient to cause dramatic excessive coastal upwelling. So, I would recommend users stick with BULK_FLUXES in conjunction with ERA5 if they have a coastal ocean in the domain.

Wow, thanks for the warning. So ERA5's "Mean east/westward turbulent surface stress" are problematic at the coast, and it's less of a problem if we compute the stresses via BULK_FLUXES? Thanks for the tip.

[Dan_chan]

Hi all,

I met a new wrong when using roms_get_era5_NCARds633_bulkflux function. Actually, it worked well yesterday and I got some forcing nc. But, today the error shows like:

Code: Select all

Processing variable mean sea-level pressure
Error using matlab.internal.imagesci.netcdflib
The NetCDF library encountered an error during execution of 'open' function - 'I/O failure (NC_EIO)'. If the data
source was an  OPeNDAP URL, see the OPeNDAP Troubleshooting section in the Users Guide.

Error in netcdf.open (line 67)
        [varargout{:}] = matlab.internal.imagesci.netcdflib ( 'open', filename, varargin{1} );

Error in internal.matlab.imagesci.nc/openToRead (line 1280)
            this.ncRootid = netcdf.open(this.Filename,'NOWRITE');

Error in internal.matlab.imagesci.nc (line 124)
                    this.openToRead();

Error in ncread (line 61)
ncObj   = internal.matlab.imagesci.nc(ncFile);

Error in roms_get_era5_NCARds633_bulkflux (line 293)
  lon = ncread(url,'longitude',Is,Ilen)-360;

Error in test1_forcing (line 12)
E = roms_get_era5_NCARds633_bulkflux(yyyy,mm,bbox,user);

I don't know what's wrong with it?

Kind Regards,
Dan

[wilkin]

Is the error repeatable?

It looks like a failure on the server side which might be a temporary problem at NCAR.

[Dan_chan]

Hi Wilkin,

You are right, the error may be fixed automatically. Now it works well.

Thank you

[mariefrahman]

Hi Wilkin,

do you know why we can used the roms_get_era5_NCARds633_bulkflux.m to get the era5 dataset? I tried the script and the error come like this:

Error in roms_get_era5_NCARds633_bulkflux (line 410)
lon = ncread(url,'longitude');

and i tried to look into the rda.ucar.edu/thredds and look inside ds633.0, and its empty. do you have idea about this?

thank you.
arief

[wilkin]

It appears NCAR changed the hostname for the THREDDS server. It is now https://thredds.rda.ucar.edu not just https://rda.ucar.edu.

If you update the line in the code that sets urlbase to read ...

Code: Select all

urlbase = 'thredds.rda.ucar.edu/thredds/dodsC/files/g/ds633.0';

... it should work.

It also appears that username:password authentication is no longer required for the ERA5 data, so all the code about parsing the user's credentials for rda.ucar.edu from .netrc is unnecessary. I'll clean and update the code on GitHub in a while, but for now that change to urlbase should do it.

[mariefrahman]

Hi Wilkin, thank you

I changed the urlbased into 'thredds.rda.ucar.edu/thredds/dodsC/files/g/ds633.0' as you mention before, but I get another error:

Error using matlab.internal.imagesci.netcdflib
The NetCDF library encountered an error during execution of 'open' function - 'file not found (-90)'.

Error in netcdf.open (line 77)
[varargout{:}] = matlab.internal.imagesci.netcdflib('open', ...

Error in internal.matlab.imagesci.nc/openToRead (line 1290)
this.ncRootid = netcdf.open(this.Filename,'NOWRITE');

Error in internal.matlab.imagesci.nc (line 124)
this.openToRead();

Error in ncread (line 69)
ncObj = internal.matlab.imagesci.nc(ncFile);

Error in roms_get_era5_NCARds633_bulkflux_update (line 478)
itime = double(ncread(url,'forecast_initial_time'))/24 + epoch;

I think that the script cannot read the nc from the url, am i right?

[wilkin]

If you got to that place in the code ...

Error in roms_get_era5_NCARds633_bulkflux_update (line 478)
itime = double(ncread(url,'forecast_initial_time'))/24 + epoch;

... then it must have read the first set of data files, right? So the urlbase is not wrong. The File Not Found error in this case possibly means you lost connection to the server, or your query is too large and causing it to time out.

That said, I get errors when the function gets to downward longwave radiation. This is read from a different catalog (meanflux) than the previous variables (sfc analysis). There might be something wrong with the catalog. In a few days time I can take a look at this more closely.

[wilkin]

Please try running this ... modified ... code again. The NCAR helpdesk claim the issue is resolved.