variable frequency in atmospheric forcing fields

General scientific issues regarding ROMS

Moderators: arango, robertson

Post Reply
Message
Author
stef
Posts: 175
Joined: Tue Mar 13, 2007 6:38 pm
Location: Independent researcher
Contact:

variable frequency in atmospheric forcing fields

#1 Unread post by stef »

Hi,

I am comparing a simulations forced with ERA5, GDAS/FNL and GFS fields using bulk fluxes (I know it can be problematic [1]).

ERA5 has radiation fluxes at the same time steps as wind, albeit temporally staggered because it's an averaged quantity, but effectively all fields have the same resolution.

GDAS/FNL and GFS have e.g. wind every 3 hours if one counts the +00h forecast/analysis and the +03h forecasts. However, the radiation fluxes are effectively only available in 6-hourly resolution (?), because they are averaged quantities and

* the +00h fields do not contain fluxes at all
* the +03h steps are centered at +1.5, and are effectively not usable because there are no averaged fluxed between +03h and +06h
* the +06h steps are averaged over 6 hours and hence centered at +03h.

So the only way to get conservative fluxes is to use the last option +06h centered at +03h, of which there are 4 available per day, so it's effectively a 6-hourly field.

Is this correct? And if so, would you use bulk fluxes where some fields are 3-hourly and some are 6-hourly? Or do you think this is a problem? Should I use consistently 6-hours instead of mixed 3 and 6 hours?

Thanks for your help!!

[1] viewtopic.php?t=3003

User avatar
wilkin
Posts: 875
Joined: Mon Apr 28, 2003 5:44 pm
Location: Rutgers University
Contact:

Re: variable frequency in atmospheric forcing fields

#2 Unread post by wilkin »

The MATLAB codes I wrote (roms_get_era5_NCARds633_bulkflux.m and roms_write_era5_NCARds633_frcfile.m) at https://github.com/johnwilkin/roms_wilkin extract hourly data for all fields from the ERA5 dataset in the Research Data Archive at NCAR (rda.ucar.edu).
John Wilkin: DMCS Rutgers University
71 Dudley Rd, New Brunswick, NJ 08901-8521, USA. ph: 609-630-0559 jwilkin@rutgers.edu

stef
Posts: 175
Joined: Tue Mar 13, 2007 6:38 pm
Location: Independent researcher
Contact:

Re: variable frequency in atmospheric forcing fields

#3 Unread post by stef »

Thanks for your reply!

I have seen the ERA5 scripts, thanks! However, my question was specific to GDAS/FNL and GFS. In that repository you linked to, do you have scripts that extract hourly data from those as well?

The question was whether
1) the radiative fluxes for GDAS/FNL and GFS really are (as I believe) available only at 6-hourly resolution, i.e. the +06 forecasts which contain averaged fluxes for the period +00 to +06 (six hours) for each synoptic time
2) if 1) is true, would you use those 6-hourly fluxes in bulk formulations together with 3-hourly fields for other variables (wind, sea level pressure etc.)

The reason I like to experiment with GDAS/FNL is because GFS makes forecasts publicly available for free via the NCEP NOMADS server, in contrast to ECMWF which provides forecasts (for all required fields) only for payment.

stef
Posts: 175
Joined: Tue Mar 13, 2007 6:38 pm
Location: Independent researcher
Contact:

Re: variable frequency in atmospheric forcing fields

#4 Unread post by stef »

Hmm, I guess one could easily compute the 3-hourly resolution from the available data:

F_a: average flux from +00 to +03 (available from the resources)
F_b: average flux from +03 to +06 (missing in the resources)
F_c: average flux from +00 to +06 (available from the resources)

Can't one compute F_b as

F_b = 2*F_c-F_a

because
E_c=E_a+E_b, where E_c (units Joule per m^2) is F_c (Watts per m^2) integrated over 6 hours, i.e. E_c = 6*F_c
and so
6F_c = 3F_a + 3F_b

Is it really as simple as that? I should test if it makes a difference to using the 6-hourly fluxes.

User avatar
wilkin
Posts: 875
Joined: Mon Apr 28, 2003 5:44 pm
Location: Rutgers University
Contact:

Re: variable frequency in atmospheric forcing fields

#5 Unread post by wilkin »

I have added codes roms_{get,write}_gfs_NCARds083_bulkflux.m for accessing the NCAR archive of GFS data to my tools https://github.com/johnwilkin/roms_wilkin. They run similarly to the ERA5 downloads, though there is presently an issue with updating the GFS THREDDS aggregation that the NCAR Helpdesk are looking into.

The codes return everything at 3-hour intervals because NCAR have unpacked those 6-hour accumulated variables before adding to the aggregation.
John Wilkin: DMCS Rutgers University
71 Dudley Rd, New Brunswick, NJ 08901-8521, USA. ph: 609-630-0559 jwilkin@rutgers.edu

stef
Posts: 175
Joined: Tue Mar 13, 2007 6:38 pm
Location: Independent researcher
Contact:

Re: variable frequency in atmospheric forcing fields

#6 Unread post by stef »

Wow, great!
The codes return everything at 3-hour intervals because NCAR have unpacked those 6-hour accumulated variables before adding to the aggregation.
That's interesting, may I ask how you found this? I didn't see it in the manuals. Do you ask at their helpdesk for this type of information?

Hmm, but they didn't update the time variables and variable names accordingly? Not sure if I understand correctly, but e.g. the variable 'Downward_Short-Wave_Radiation_Flux_surface_Mixed_intervals_Average' in [1] has dimensions of

Code: Select all

reftime12 time12 lat lon
reftime12 is (see [2]):
0.0, 0.0, 6.0, 6.0, 12.0, 12.0, 18.0, 18.0, ...
time12 is (see [3]):
3.0, 6.0, 9.0, 12.0, 15.0, 18.0, 21.0, 24.0, ...

Both are hours, and so the difference valid-time minus reftime in hours is
3.0, 6.0, 3.0, 6.0, 3.0, 6.0, ....

The lags of 3 and 6 hours, together with the fact that the variable is explicitly called "Mixed intervals average" would suggest that they didn't unpack it? I tried to compare it with recent GFS data from NOMADS, but the last timestep in "time12" is 52593 hours since 2015-07-08 UTC, which is 2021-07-07 09:00 I think. Can't get it from NOMADS, so I compared the following:

* Downward_Short-Wave_Radiation_Flux_surface_6_Hour_Average for 2021-01-02 +06 in [4]

* Downward_Short-Wave_Radiation_Flux_surface_Mixed_intervals_Average for time12=2021-01-02 06:00:00+00:00 in [5]


The lon/lat values are roughly GMT+6. The figures look very similar.

[1] https://rda.ucar.edu/thredds/dodsC/aggr ... /Best.html
[2] https://rda.ucar.edu/thredds/dodsC/aggr ... ?reftime12
[3] https://rda.ucar.edu/thredds/dodsC/aggr ... cii?time12

[4] rda.ucar.edu/thredds/dodsC/files/g/ds084.1/2021/20210102/gfs.0p25.2021010200.f006.grib2?time[0:0],lat[371:461],lon[299:401],Downward_Short-Wave_Radiation_Flux_surface_6_Hour_Average[0:0][371:461][299:401]

[5] rda.ucar.edu/thredds/dodsC/aggregations/g/ds083.3/1/Best?time12[16041:16041],lat[371:461],lon[299:401],Downward_Short-Wave_Radiation_Flux_surface_Mixed_intervals_Average[16041:16041][371:461][299:401]

User avatar
wilkin
Posts: 875
Joined: Mon Apr 28, 2003 5:44 pm
Location: Rutgers University
Contact:

Re: variable frequency in atmospheric forcing fields

#7 Unread post by wilkin »

As I understand it, the reftime tells when that forecast cycle started. There is a new forecast every 6 hours, and data from that cycle is saved at 3-hour intervals. Hence, reftime repeats in pairs 0, 0, 6, 6, whereas the time coordinate advances by 3 hours for each record.

This is a Forecast Model Run Collection (FMRC). If the coordinates attribute for a variable is coordinates: reftime12 time12 lat lon then reading right to left this tells that the time coordinate variable is time12.

The units for time12 are clear units: Hour since 2015-02-22T00:00:00Z so it's straightforward to match up time with the data.

In the preamble comments to roms_get_gfs_NCARds083_bulkflux.m I've done my best to explain what I consider the best match-up between the GFS outputs and the ROMS variables for bulk fluxes. I compared all outputs downloaded to maps of the same fields from ERA5 (which are less ambiguous) and I'm comfortable that this passed that basic sanity check.

It's my understanding that variables denoted "average" are just that ... they are averages over the 3-hour interval as opposed to snapshots at a time instant. This is because are derived from accumulated variables that integrate over time. So W m-2 become J m-2. Much older MATLAB code I wrote to process files directly downloaded from ECMWF has to unpack those J m-2 etc. accumulations, which took some work because the data were integral 0 to 3 hours and integral 0 to 6 hours. But it's easy to see if you have it right in shortwave radiation because if you get it wrong you don't have a "sinusoidal" diurnal cycle but a saw-tooth. I know, because I got it wrong at first.

In the NCAR RDA archive, all heat flux components have units W m-2, not J m-2, so it's clear they have removed the time integration.

If you think my code is matching up the wrong variables let me know. But I'm optimistic it's right.

As I noted in my earlier post, and you discovered yourself, while there are files for 2023 for GFS at rda.ucar.edu they are not all added to the THREDDS aggregations. I have lodged a ticket with the Helpdesk asking them to update the aggregation since my code works on those single end-points, not the individual files.
John Wilkin: DMCS Rutgers University
71 Dudley Rd, New Brunswick, NJ 08901-8521, USA. ph: 609-630-0559 jwilkin@rutgers.edu

stef
Posts: 175
Joined: Tue Mar 13, 2007 6:38 pm
Location: Independent researcher
Contact:

Re: variable frequency in atmospheric forcing fields

#8 Unread post by stef »

It's my understanding that variables denoted "average" are just that ... they are averages over the 3-hour interval as opposed to snapshots at a time instant.
So the question is whether they are 3-hour averages, or mixed 3-hour averages and 6-hour averages, right? Earlier you wrote
The codes return everything at 3-hour intervals because NCAR have unpacked those 6-hour accumulated variables before adding to the aggregation.
When you say "unpack", do you mean they turned accumulated variables into averages (energy into radiant fluxes in the case of radiation), or do you mean they went a step further and calculated 3-hour averages from 6-hour avarages?
If you think my code is matching up the wrong variables let me know. But I'm optimistic it's right.
Did I understand your code correctly in that it assumes rda.ucar staff have computed 3-hour averages from 6-hour averages? Or did I miss something and your code *does* make this conversion?
As I noted in my earlier post, and you discovered yourself, while there are files for 2023 for GFS at rda.ucar.edu they are not all added to the THREDDS aggregations. I have lodged a ticket with the Helpdesk asking them to update the aggregation since my code works on those single end-points, not the individual files.
Understood. Thanks, I might benefit from this too!

I made some plots comparing downward shortwave datiation datasets:

1) GFS from NOMADS server
2) rda.ucar ds084.1 (NCEP GFS 0.25 Degree Global Forecast Grids Historical Archive)
3) rda.ucar ds083.3 aggregations Best (NCEP GDAS/FNL 0.25 Degree Global Tropospheric Analyses and Forecast Grids)

All comparisons are for 18N,90W (Bay of Bengal, roughly GMT+6). Dates are UTC.

First is a comparison for recents dates, just to see if rda.ucar does any modifications to the 6-hour averages available from NOMADS.

The second is a comparison for older dates which are available in ds083.3 aggregations Best. The line for 2) is the same as before (but different dates), and the markers 'o' show reftime12+03 and reftime12+06.

So it looks like reftime12+06 shows 6-hour averages, right? So we have to account for this and turn them into 3-hour averages, right?

rad_gfs_vs_gdas_6days_1d_1.png
Attachments
rad_gfs_vs_gdas_6days_1d.png

User avatar
wilkin
Posts: 875
Joined: Mon Apr 28, 2003 5:44 pm
Location: Rutgers University
Contact:

Re: variable frequency in atmospheric forcing fields

#9 Unread post by wilkin »

Looks to me like your magenta dots are simply half way between the orange dots, and if you plot them together you get a nice diurnal cycle resolved at 3-hour intervals? I think the files you are reading have correctly parsed data from the 3 and 6 hour fields in each cycle. No further tinkering required.

There were bugs in the code that I hastily added to roms_wilkin to extract GFS from rda.ucar.edu - very embarrassing. But I got it going and what I show in the plot below is swrad during January 2020 in the middle of the Bay of Bengal. ERA5 and GFS. ERA5 is hourly and GFS 3-hourly.
swradBoB.png
My code reported reading these variables:

Code: Select all

 Reading Downward_Short-Wave_Radiation_Flux_surface_Mixed_intervals_Average
   Reading Upward_Short-Wave_Radiation_Flux_surface_Mixed_intervals_Average
    got units W.m-2
and took their difference to get swrad. This is reading from the THREDDS aggregation https://rda.ucar.edu/thredds/dodsC/aggr ... /Best.html
John Wilkin: DMCS Rutgers University
71 Dudley Rd, New Brunswick, NJ 08901-8521, USA. ph: 609-630-0559 jwilkin@rutgers.edu

stef
Posts: 175
Joined: Tue Mar 13, 2007 6:38 pm
Location: Independent researcher
Contact:

Re: variable frequency in atmospheric forcing fields

#10 Unread post by stef »

and if you plot them together you get a nice diurnal cycle resolved at 3-hour intervals?
So I'm comparing this with ERA5 data, and also what I would get from applying a (hypothetically necessary) conversion from the 6-hour averages to 3-hour averages of the RDA aggregation/Best fluxes.

Below figure shows

* black diamonds show modified ERA5 'Mean surface downward short-wave radiation flux' (msdwswrf, W/m^2). The raw data are hourly averages mapped to the end of the avaraging interval, and my modification is simply to take 3-hour avarages, mapped to the end of the 3-hour interval.

* green squares show ds083.3 aggregation/Best mixed averages, except that the values as +06 are modified by taking twice their value minus the preceeding +03 value (as per hypothesis in the post above)

The green squares look more similar to the black diamonds than the yellow circles, so my conclusion is that the conversin from 6-hour averages to 3-hour averages has to be made by the user, and is not done by rda.ucar.

The non-zero yellow circles (non-zero 6-hour averages) occur at UTC 06:00 and UTC 12:00, which at 18N,90W (Bay of Bengal, roughly GMT+6) is roughly 12:00 and 18:00 local time. The fact that they are roughly equally large agrees with the notion that they represent fluxes during the entire noon/afternoon. The green squares are nonzero at UTC 03:00, 06:00, 09:00, 12:00 corresponding to local times 9:00, 12:00, 15:00, 18:00. The points at 12:00 and 15:00 local time represent averages 9-12 and 12-15, and intuitively they should be roughly equally large, which they are (the caveat is I have no clue about the diurnal cycle of deep convection at that point, and how that would affect short-wave radiation.... ? but that aside, it agrees with intuition).


Note that all quantities are mapped to the end of the averaging interval, but for actual forcing fields, they should be mapped to the center of the interval. I didn't do this here because it involves yet another manipulation of the raw data.
rad_gfs_vs_gdas_6days_1d.png

User avatar
wilkin
Posts: 875
Joined: Mon Apr 28, 2003 5:44 pm
Location: Rutgers University
Contact:

Re: variable frequency in atmospheric forcing fields

#11 Unread post by wilkin »

I think we might be speaking at cross purposes here.

You are plotting values extracted from individual files, which have 3- and 6-hour average values in each cycle. I was plotting "3-hour" data from the THREDDS FMRC aggregations. I may be wrong, but I assumed whoever configured the FMRC correctly accounted for the different averages. You have the workflow set-up to test that, so I invite you to check that out.

Or, since the folks at NCAR are usually super helpful, perhaps log a ticket with their helpdesk asking for clarification.

On the time coordinate ...
Note that all quantities are mapped to the end of the averaging interval, but for actual forcing fields, they should be mapped to the center of the interval. I didn't do this here because it involves yet another manipulation of the raw data.
... all you need to do is define a new time coordinate for your adjusted variables that reports 1.5 hour offset middle of the average interval.

ROMS is going to linearly interpolate based on each variables time coordinate, so let it do the work instead of yourself.

Keep the original time coordinate for the other variables (winds etc.) because they are snapshots of the state, not averages.
John Wilkin: DMCS Rutgers University
71 Dudley Rd, New Brunswick, NJ 08901-8521, USA. ph: 609-630-0559 jwilkin@rutgers.edu

stef
Posts: 175
Joined: Tue Mar 13, 2007 6:38 pm
Location: Independent researcher
Contact:

Re: variable frequency in atmospheric forcing fields

#12 Unread post by stef »

You are plotting values extracted from individual files
Well I'm plotting the files that are accessed by the scripts you mention above, where you wrote
roms_{get,write}_gfs_NCARds083_bulkflux.m for accessing the NCAR archive of GFS data to my tools https://github.com/johnwilkin/roms_wilkin.
Specifically in 'roms_get_gfs_NCARds083_bulkflux.m' the dataset is currently

Code: Select all

dataurl = 'https://rda.ucar.edu/thredds/dodsC/aggregations/g/ds083.3/1/Best'
and the downward shortwave variable name is

Code: Select all

Downward_Short-Wave_Radiation_Flux_surface_Mixed_intervals_Average
This is one of the datasets/variables I plotted above, and it seems the converion from 6-hourly averages to 3-hourly averages is up to the user (not done by rda.ucar).

I was plotting "3-hour" data from the THREDDS FMRC aggregations.
Hmm, not sure which dataset you mean. Is this hosted by rda.ucar or a different organization?

https://rda.ucar.edu/datasets/ds083.3/? ... on#!access

(login necessary) I can't see FMRC.

stef
Posts: 175
Joined: Tue Mar 13, 2007 6:38 pm
Location: Independent researcher
Contact:

Re: variable frequency in atmospheric forcing fields

#13 Unread post by stef »

Or, since the folks at NCAR are usually super helpful, perhaps log a ticket with their helpdesk asking for clarification.
Ok, I'll log a ticket! I'll update this post once I get a response.

User avatar
wilkin
Posts: 875
Joined: Mon Apr 28, 2003 5:44 pm
Location: Rutgers University
Contact:

Re: variable frequency in atmospheric forcing fields

#14 Unread post by wilkin »

I see your point.

The individual files, e.g., https://rda.ucar.edu/thredds/dodsC/file ... grib2.html (one every 3 hours)
have variable Downward_Short-Wave_Radiation_Flux_surface_3_Hour_Average

Perhaps that is what should have been aggregated in THREDDS and was not.

If you can reconcile these 3_Hour_Averages with the mixed interval averages, I can add code to unpack them as you propose. But I'll wait until you can confirm with certainty that you have sorted out what is going on.

I wonder if the rain is also based on an average and needs an adjusted time coordinate.
John Wilkin: DMCS Rutgers University
71 Dudley Rd, New Brunswick, NJ 08901-8521, USA. ph: 609-630-0559 jwilkin@rutgers.edu

stef
Posts: 175
Joined: Tue Mar 13, 2007 6:38 pm
Location: Independent researcher
Contact:

Re: variable frequency in atmospheric forcing fields

#15 Unread post by stef »

The individual files, e.g., https://rda.ucar.edu/thredds/dodsC/file ... grib2.html (one every 3 hours)
have variable Downward_Short-Wave_Radiation_Flux_surface_3_Hour_Average
Well, I guess not all of them. The +09 file you linked above

https://rda.ucar.edu/thredds/dodsC/file ... grib2.html

does indeed have a Downward_Short-Wave_Radiation_Flux_surface_3_Hour_Average, but the +06 file:

https://rda.ucar.edu/thredds/dodsC/file ... grib2.html

only has a Downward_Short-Wave_Radiation_Flux_surface_6_Hour_Average. So apparently the +03 files have 3-hour averages, the +06 files have 6-hour averages, and the +09 files have 3-hour averages again.

That's interesting, I have never looked at the +09 files.
I wonder if the rain is also based on an average and needs an adjusted time coordinate.
I see, so we have to look at all accumulations/averages.

stef
Posts: 175
Joined: Tue Mar 13, 2007 6:38 pm
Location: Independent researcher
Contact:

Re: variable frequency in atmospheric forcing fields

#16 Unread post by stef »

It's been about a week since I asked at the rda.ucar helpdesk, no answer yet. They seem to re-organize their website, maybe the ticket got lost. It may also be possible that the helpdesk tickets I lodge have lower priority than others (due to my non-institutional email address, no affiliation with large organizations etc.). Perhaps somebody else want's to try their luck...?

Post Reply