refactoring ERA5 data mode - #396
uturuncoglu wants to merge 9 commits into
Conversation
|
@billsacks I am still working on it. Once tests are passed I'll let you know. |
|
Hi @uturuncoglu @billsacks . May I ask if this is still planned? I'm interested so I ran some tests in UFS. There is also interest in using ERA5 for datm-ocn-ice work, like for RTOFS development. fyi, AI came back with some review notes (NickSzapiro-NOAA#11 (review)). |
|
@uturuncoglu I got a link to this (draft) PR from @NickSzapiro-NOAA in the context of UFS DATM with the "data" atmosphere coming from ERA-5. 2 questions for either/both of you:
|
|
@sanAkel @NickSzapiro-NOAA The only thing that I need to do is to add new formulas to valuate saturation mixing ratio. I am plaining to add new option for it. We could also introduce namelist option to get those fractions from datm_in. I am also taking off starting from next week and back to US end of the month. So, I am not sure if I could do everything before love since I am also trying to finalize things for other projects but if I could do, I'll update you through the PR and and mark tis as ready to review. |
Thanks @uturuncoglu for ⬆️ update. If you have it ready before your leave, great! If not, we will wait or make do whatever we got! |
@uturuncoglu are you back? Hope you had a nice break! Checking if you have time to finish this work? 🙏 |
|
@sanAkel Yes, I am back. I am trying to catch other things but I could work on this next week and we could finalize it. |
|
@sanAkel BTW, i am assuming this is tested without any issue in your end without the remaining additional things. |
@uturuncoglu @raphaeldussin is going to set up a DATM run with the ERA5 forcing. If you happen to have a set of sample forcing files (from the ERA5 reanalysis), that you can point us to, we can test it at our end. Thanks again! |
|
@sanAkel i am trying the finalize the implementation but I think new algorithm to calculate qsat is wrong and give negative and very large values. I tested following code under Python like following, and it is giving |
|
AI suggests it is a units issue ... ie |
|
I used the Bolton approximation (which is fairly accurate at the observed range of temperatures), you can compare with my code here https://github.com/raphaeldussin/newERAtools/blob/e29d89707058cca52abb5d8015a85894756f9c37/compute_q2_ERA5.f90#L154 |
@uturuncoglu Please see ⬆️ response from @raphaeldussin. If your error persists, please check with @raphaeldussin. |
|
@sanAkel @raphaeldussin @NickSzapiro-NOAA I fixed the unit issue and pushed it to the fork. Now you could select type of algorithm in
Anyway, I tested this with UFS WM |
|
I also sync the branch with main. |
|
That's awesome @uturuncoglu 👏 ! @NickSzapiro-NOAA Can we run the same RT test @uturuncoglu ran on wcoss-2? To rephrase my question, do we have the same datasets for that test available? |
|
Thanks @uturuncoglu . I will start some tests tonight The same RT input data should be available on all platforms @sanAkel |
Thanks @NickSzapiro-NOAA
What do we need to do make into an "ocean" RT? Please list steps.. 🙏 |
Model configuration and input data with forcing, ICs, and fix files. I think there is quite some overlap with @uturuncoglu in datm-ocean-ice efforts, but I'm not sure this ESCOMP/CDEPS thread is right place to continue this topic. Happy to continue in ufs-community/ufs-weather-model#3294 , NOAA-EMC/RTOFS_GLO#147 , or wherever else as well |
|
To re-reiterate the ERA5 precipitation issue I mentioned at the UFS-coastal technical tag up today: This may be rather particular to our current DATA-CICE and CICE-SCHISM setups but, If ERA5 provides precipitation in units of meters (cumulative meters of rain per hour), rather than kg/m^2/s, it may be misinterpreted as a rate leading to 3.6X more precip than is actually occurring. (Code/line numbers are from datm/datm_datamode_era5_mod.F90 as of head e301493) But magnitude alone can't separate the two conventions the L682 test is choosing between. For true precip R mm/hr the file value is R/1000 as metres-of-1-hour-accumulation vs R/3600 as kg/m²/s — they differ by only 3.6×, and both sit in the same 0…~10⁻² range for ordinary precip. So the 0.01 cutoff effectively asks only "is there a ≥10 mm/hr cell somewhere in the domain at the first timestep?" Model reports ""already in kg/m^2/s unit!" to output log on L683, even though CDEPS doesn't actually check the units in the ERA5 input file. In our SCHISM CICE test case (prior to the implementation of Faxa_snow) this led to far too much snow falling on the sea ice, resulting in too thick ice that persisted too long into the spring. NOTE ALSO currently Faxa_snow is left out of the above conversion (L685-689) I think ideally CDEPS should read in units from the ERA5 file explicitly and only if units are not found should it make a guess. Barring that the guess should be reported in the log as an assumption that has some ambiguity associated with it. |
| ! Initialize value of export state | ||
| if (associated(Faxa_rain)) Faxa_rain(:) = 0.0_r8 | ||
| if (associated(Faxa_rainc)) Faxa_rainc(:) = 0.0_r8 | ||
| if (associated(Faxa_rainl)) Faxa_rainl(:) = 0.0_r8 | ||
| if (associated(Faxa_snow)) Faxa_snow(:) = 0.0_r8 | ||
| if (associated(Faxa_snowc)) Faxa_snowc(:) = 0.0_r8 | ||
| if (associated(Faxa_snowl)) Faxa_snowl(:) = 0.0_r8 |
There was a problem hiding this comment.
How do you decide which export state variables are initialized (to zero)?
Say, what about taux, tauy?
There was a problem hiding this comment.
@NickSzapiro-NOAA As i see the fields are created in dshr_fldlist_realize() call in dshr/dshr_fldlist_mod.F90 and it seems that there is no any initial value assigned to them. Maybe ESMF is initializing them internally but I am not %100 sure about it. We could use ESMF_FieldFill to set zero initially in that level to be in the safe side. @billsacks let me know what do you think?
There was a problem hiding this comment.
I did a very quick look through other CDEPS components and it looks like this kind of zeroing generally isn't done. I haven't looked at this ERA5 datamode enough to know if it might be needed for this one specially. My initial thought is that these variables should be set in the call to datm_datamode_era5_advance, so there's probably no need to initially zero them, but you shouldn't take my word for it because I'm not very familiar with this.
@NickSzapiro-NOAA Can you please help? |
@uturuncoglu It is not possible for me test these changes at this moment; see here for details. I will work with @NickSzapiro-NOAA on a way forward. Meanwhile you + others can proceed with your PR (test, merge, etc). |
|
@sanAkel Okay. Let me check on Hercules and try to run couple of UFS RTs. |
|
Hi @uturuncoglu . It seems #421 can be excepted compiler warning. Happy to test on UFS side or if anything else that would be helpful here |
|
@NickSzapiro-NOAA I think we could add fix for #421 to this PR. Does it work for you? |
|
It seems the suggested fixes are to either shorten the name or update compilers. Fixing would be nice, but if this is what's blocking this PR, I wouldn't worry so much about it because
|
|
@NickSzapiro-NOAA It would be nice to put the fix here. I am not sure about the timeline in here. People in here trying to finalize the CESM release. So, @billsacks is the right person to get some answer about it. This is just adding couple of namelist option and updates ERA5 data mode. It would be fine for CESM. BTW, is this tested with UFS WM? |
|
Thanks! I'll test again over the next few days ... I think it's just datm_cdeps_lnd_era5{_rst} as there is no more HAFS-related ERA5 RT |
|
@uturuncoglu and @NickSzapiro-NOAA - I haven't been following this discussion carefully. Can you clarify what the question is for me? |
|
@billsacks This PR basically restructure the ERA5 data mode to be used under CESM and also UFS WM. Since everybody is working on CESM release at this point, do we have bandwidth to test and merge this PR? I think @NickSzapiro-NOAA could test it on UFS WM side (let me know if I am wrong) but this needs to be tested with CESM too. So, at this point I am not sure about willingness to bring this to ESCOMP based on the work going on CESM side. |
billsacks
left a comment
There was a problem hiding this comment.
I have a bunch of specific comments below, but also have some general questions / comments about this.
My over-arching concern is that, by trying to provide so much flexibility, I'm concerned about both:
(1) The maintainability of this: it will be harder to make changes because you need to handle the possibility that any field could be present or absent
(2) The usability of this: if all checks are bypassed, the code will happily proceed if you have accidentally omitted some fields - filling in those fields with guesses that may not be what you want. For example, some of the field filling seems very arbitrary - like the swdn components, the assumed albedo, and the u/v wind components. As a user of this data mode, I could imagine preferring to be told that some required fields are missing that I need to add to my dataset, so I can then think about how to do it, rather than having fields silently filled in with what may be an oversimplified guess.
I am especially concerned because it will be (nearly) impossible to test all of the possible permutations here. My suggestion would be to only allow the bare minimum flexibility that is truly needed - so most fields are still required, and you perhaps just have flexibility in a few select fields. Would that be possible?
If you truly need to keep all of this flexibility, have you done a careful review – and ideally some testing – to confirm that any absent stream field is handled correctly?
I would also like to have reviews from:
- @mvertens - I'm hesitant to move ahead with this without at least a quick look from you, since you've been doing a lot of work with the data models lately
- Someone on the CESM side who intends to use this data mode. Would that be the ocean group? Ideally there would be reviews both from the science and software sides. @uturuncoglu do you know who would be appropriate here? - have you been working with anyone on the CESM side on this?
| if (associated(Sa_z)) then | ||
| lsize = size(Sa_z) | ||
| else | ||
| if (mainproc) write(logunit,*) subname,' Sa_z is not given. Try with Sa_pslv to get stream size' | ||
| end if | ||
| if (associated(Sa_pslv)) then | ||
| lsize = size(Sa_pslv) | ||
| else | ||
| if (mainproc) write(logunit,*) subname,' Sa_pslv is also not given. Exiting!' | ||
| call shr_log_error(subname//'ERROR: Sa_z and/or Sa_pslv must be associated for era5 datamode', rc=rc) | ||
| return | ||
| end if |
There was a problem hiding this comment.
Can you double-check this logic? It looks to me like the error return happens whenever Sa_pslv is unassociated, without considering whether Sa_z was associated. Should the associated(Sa_pslv) conditional be inside the else block for Sa_z?
| flds_preso3, & | ||
| export_all | ||
| export_all, & | ||
| skip_field_check, & |
There was a problem hiding this comment.
It seems potentially confusing that skip_field_check is only implemented for ERA5 data mode. If a user sets this for a different data mode, it could be confusing that it doesn't do anything. What would you think about reading this in all of the datm datamodes, but then in the non-era5 modes, aborting if this is set to true with a message like "skip_field_check not implemented for XXX datm mode"
| flds_preso3, & | ||
| export_all | ||
| export_all, & | ||
| skip_field_check, & |
There was a problem hiding this comment.
If skip_field_check will be needed in CESM, it will need to be added to the namelist definition xml file. It's not clear to me if this will be needed in CESM or not.
| end if | ||
|
|
||
| call shr_cal_date2julian(target_ymd, target_tod, rday, model_calendar) | ||
| rday = mod((rday - 1.0_R8),365.0_R8) |
There was a problem hiding this comment.
Should the 365 here and on the next line be made more general to get the current number of days per year (using something from the calendar) - to handle leap years and non-standard calendars?
|
|
||
| call shr_cal_date2julian(target_ymd, target_tod, rday, model_calendar) | ||
| rday = mod((rday - 1.0_R8),365.0_R8) | ||
| cosfactor = cos((2.0_R8*SHR_CONST_PI*rday)/365 - phs_c0) |
There was a problem hiding this comment.
Should this 365 be a real value rather than integer? See also comment above about generalizing the 365.
| if (associated(strm_Sa_pslv)) then | ||
| pbot = strm_Sa_pslv(n) | ||
| else if (associated(strm_Sa_pbot)) then | ||
| pbot = strm_Sa_pbot(n) | ||
| end if |
There was a problem hiding this comment.
Should there be an else here that aborts if neither is set - because then we won't have a pbot?
| if (associated(strm_Sa_t2m)) then | ||
| tbot = strm_Sa_t2m(n) | ||
| else if (associated(strm_Sa_tbot)) then | ||
| tbot = strm_Sa_tbot(n) | ||
| end if | ||
|
|
||
| if (associated(strm_Sa_pslv)) then | ||
| pbot = strm_Sa_pslv(n) | ||
| else if (associated(strm_Sa_pbot)) then | ||
| pbot = strm_Sa_pbot(n) | ||
| end if |
There was a problem hiding this comment.
The retrieval of tbot and pbot could be moved into if (.not. associated(strm_Sa_shum)), because I think they're only needed in that block.
| if (associated(strm_Sa_pbot) .and. associated(strm_Sa_tbot) .and. associated(strm_Sa_shum)) then | ||
| if (associated(Sa_dens)) then | ||
| Sa_dens(n) = strm_Sa_pbot(n)/(rdair*strm_Sa_tbot(n)*(1 + 0.608*strm_Sa_shum(n))) | ||
| end if | ||
| end if |
There was a problem hiding this comment.
It seems like this logic should be flipped around: first check if Sa_dens is associated; then check if the three needed stream variables are associated; if they aren't, then it's an error.
| if (associated(strm_Faxa_swdn)) then | ||
| if (associated(Faxa_swvdr)) Faxa_swvdr(n) = strm_Faxa_swdn(n)*(0.28_R8) | ||
| if (associated(Faxa_swndr)) Faxa_swndr(n) = strm_Faxa_swdn(n)*(0.31_R8) | ||
| if (associated(Faxa_swvdf)) Faxa_swvdf(n) = strm_Faxa_swdn(n)*(0.24_R8) | ||
| if (associated(Faxa_swndf)) Faxa_swndf(n) = strm_Faxa_swdn(n)*(0.17_R8) | ||
| end if |
There was a problem hiding this comment.
Again, I think this logic should be flipped so you can abort if the various Faxa things are associated but strm_Faxa_swdn is not associated.
| if (associated(strm_Faxa_swdn)) then | ||
| if (.not. associated(strm_Faxa_swnet) .and. associated(Faxa_swnet)) then | ||
| avg_alb = ( 0.069 - 0.011*cos(2.0_R8*yc(n)*degtorad ) ) | ||
| Faxa_swnet(n) = strm_Faxa_swdn(n)*(1.0_R8 - avg_alb) | ||
| Faxa_swnet(n) = Faxa_swnet(n)*3600.0_r8 ! to J/m^2 | ||
| end if | ||
| end if |
There was a problem hiding this comment.
Again, I think this logic should be flipped so that you can abort if Faxa_swnet is associated but strm_Faxa_swdn isn't available, so there's no way to set Faxa_swnet.
I'm fine merging this once it has undergone sufficient review. I did an initial - fairly quick, and not comprehensive - review, but I think it needs review from one or two people on the CESM side who have more knowledge of this data mode, if there is anyone. @uturuncoglu can you identify anyone on the CESM side who has this familiarity? If you're not sure, we could send an email to CSEG to see who could do this. |
I know both @alperaltuntas and @gustavo-marques use CESM in data mode. Maybe one of them could kindly review? |
|
@billsacks I create this PR since CESM folks are trying to replace JRA55 with ERA5. So, they asked me to make this modification. But, during this work restructuring PR came and break things also in UFS WM side. Adding controls for the export fields are great and prevents issues but under UFS WM we are mostly providing subset of variables for different data mode. It was working fine for us. So, those cases are not working at this point, since user needs to provide all the fields now. Then, I also put fix to this PR about it and we could keep continue to support regression tests and configurations under UFS WM. I think just removing control for some set of variables will not work since everybody has different needs. It would be hard to maintain. Anyway, at this point we have couple of options for now: (1) we could go with the implementation proposed in here, (2) we will force the existing users to change all the data coupled configurations and their input files, and (3) we could restructure CDEPS and make it more flexible (such as rather then hard coding those variables names to the data mode, we could ingest through the YAML file etc., so, every body will have their own configuration file and define subset of fields). All, works for me but I am not funded to do any of them. I could only make small fixes in here if it is requested. Anyway, I am happy to discuss more if you need. |
|
I understand your points, @uturuncoglu , and at the same time, I'm afraid that trying to add as much generality as is added here will be hard to maintain; and, worse, it may sometimes appear that things are working when in fact they're not working right, because there are missing fields that may or may not get filled correctly. It would help me to understand the actual fields present on the forcing datasets for CESM and for UFS for this data mode. If they are similar but with just some small differences, I think the best approach is to just give the necessary small amount of flexibility. If they are significantly different, then it feels like it would be cleanest and safest to have two separate data modes. If it's easier to move this to an in-person chat or video call, let's do that. @mvertens I'd also like to know if you'd like to be part of this discussion. |
|
@billsacks - I'd really like to be part of that discussion. I think there is interest in NorESM to run ERA5 forcing. @mvdebolskiy - do you want to join as well? I recall some change that is still on a CDEPS branch from a while ago to get the ERA5 working. |
|
@mvertens sure, let me know. Also, maybe jan ivar? |
|
I would like to be included as well, please. Maybe ERA5 datamode is particular since several usages want to share it |
|
I've run MOM6 with ERA5 with the FMS/coupler and I'm happy to share my 2 cents |
|
Please include me in the call as well. |
|
@uturuncoglu - I talked about this with @mvertens . She shares some of my concerns about bypassing the field checks and trying to fill in missing fields. When I started reviewing this, I was hoping this would be quick. Now I see that it is probably going to take significant discussion. With the upcoming CESM3 release, I don't have a lot of time to devote to this, so I'd prefer if we can wait to figure this out until after the current push towards the CESM3 release. But how time-critical is this for different applications? Do some people need this to move forward in some form? If this does need to move forward soon, then @mvertens and I think that the best approach could be to make a new data mode that works without needing to do something like skip_field_check. This datamode would handle exactly the set of fields that are provided in the dataset you're handling. There could still be a bit of flexibility where needed, but not so much flexibility that you need to bypass all of the field checks - because that leads to a situation where it can look like the model is working when in fact some fields have garbage, which is something we definitely want to avoid. This could mean that we need different datamodes for different flavors of ERA5; this feels better to us than trying to make a single generic datamode that tries to handle multiple very different situations. |
|
The OMWG would like this to be available in CESM 3. We are okay with introducing a new data mode if that means a safer solution. (In fact, that was our initial suggestion). |
|
@billsacks Again, any approach work for me but as I know the current version ESCOMP/CDEPS can not be synced with UFS WM easily and requires lots of work in their end. So, probably, they will maintain their on CDEPS fork and they will diverge from ESCOMP. If this is what NCAR wants, that is fine for me. Anyway, I do not have resource to work on this more. |
|
UFS would rather add more datamodes than further diverge from ESCOMP but recognize there is a cost to duplicating much of the same code and risk of different "ERA5"s diverging. Thanks for your effort and persistence with this @uturuncoglu in managing to unify different usages. On flexible ways to fill fields, maybe approaches like CDEPS/datm/datm_datamode_jra_mod.F90 Lines 261 to 265 in b539dff CDEPS/datm/datm_datamode_jra_mod.F90 Lines 297 to 308 in b539dff for rain/snow, partitioning shortwave, ... work in unified clean and narrow ways I would lean towards not advertising variables not in the input data rather than give fill values to a consumer |
|
I think that introducing a new data modes would help solve this problem. I know there will be NorESM specific changes that we want to bring in. Just as a thought - something like ERA5_UFS and ERA5_CESM and ERA5_NorESM. |
|
Thanks @NickSzapiro-NOAA and @raphaeldussin for participating in the discussion next week. My preference:
|
I hope we can find a solution that doesn't require a complete fork / divergence but still allows this to work for everyone. I'm thinking that separate ERA5 forcing modes for the different use cases might be most appropriate. My feelings agree closely with @NickSzapiro-NOAA 's recent comment. I'd be happy to participate in a meeting to discuss this further. Or if others are happy with introducing specific ERA5 forcing modes for now, that could be the simplest - @mvertens and I are in agreement that it will be easier to accept (for example) a reworked ERA5 mode for UFS if it is done in a way that it won't risk interfering with the NorESM / CESM needs. I'm also thinking: If we later find that there is significant duplication between the different ERA5 modes, it might be easier to consolidate them at that point, once we can all see more concretely where they are similar and where they differ. Is there someone with time to move this forward? If so, let us know how you'd like to move forward. |
|
I can work on a new datamode for UFS |

Description of changes
This PR aims to restructure ERA5 data mode to support:
(1) Using ERA5 as a replacement of JRA55 under CESM framework. This work does not include any necessary change needs to be done in CIME interface.
(2) Adding new data component configuration option for data atmosphere (DATM) as
skip_field_check. The default value is.false.and does not change the exiting behavior. By setting it.true.the field check part is skipped and user is able to provide just set of fields supported by the data mode.(3) The changes related to support UFS WM RTs with help of (2) and also some other changes.
Specific notes
Contributors other than yourself, if any:
CDEPS Issues Fixed (include github issue #):
Are there dependencies on other component PRs (if so list):
No
Are changes expected to change answers (bfb, different to roundoff, more substantial):
No
Any User Interface Changes (namelist or namelist defaults changes):
New optional
skip_field_checkoption for datmTesting performed (e.g. aux_cdeps, CESM prealpha, etc):
UFS WM RTs:
UFS Coastal:
Hashes used for testing: