search

NCAR / DART

Data Assimilation Research Testbed
https://dart.ucar.edu/
Apache License 2.0
195 stars 145 forks source link

Converting sea ice freeboard observations #32

Closed timhoar closed 3 years ago

timhoar commented 4 years ago

Ed note: This issue came as an email to the DART team.

I’m a graduate student at UW working with CC Bitz (if she has already reached out to you about this, I apologize for doubling down). I am interested using DART to assimilate sea ice freeboard measurements from ICESat-2. I have the data (currently as .h5 files) and from what I understand, I will need to create a new observation TYPE and/or QUANTITY and then write/adapt an observation converter and forward operator. The DART website kindly encouraged me to reach out for pointers and advice regarding how best to go about these steps. I’ve found the website and the tutorials extremely helpful, but I am new to DART and would like to make sure that I set up the converter and forward operator correctly.

timhoar commented 4 years ago

Thanks for going through the tutorials - you are starting off on the right foot.

The observation QUANTITYs are defined in _assimilation_code/modules/observations/DEFAULT_obs_kindmod.F90 and you're in luck ... there is already a QTY_SEAICE_AGREG_FREEBOARD

The observation TYPEs are defined in the forward operator code. Since there are many logical collections of forward operators, and there is no need to support ocean observations if you are doing Lorenz_63, for example, there are many Fortran modules that provide support for the myriad observations. You're in luck again ... _observations/forward_operators/obs_def_cicemod.f90 has the infrastructure to support SAT_SEAICE_AGREG_FREEBOARD observations by calling the _get_expected_agregfreeboard() routine.

You should look through the _obs_def_cice.f90:get_expected_agregfreeboard() routine to see if it employs the algorithm you want to use to convert the model state to the freeboard quantity (pay attention to the units of the observation and the units of the model state - you have not said what model you are using ... CESM CICE?)

Your luck runs out (or the fun begins, depending on your personality) with converting the HDF5-format observations to a DART observation sequence file. There is an _observations/obsconverters/cice directory that has some hints about converting observations, but not in HDF5 format. Actually, there is a comment in _modis_ist_to_obsnetcdf.f90 that indicates Yongfei simply converted the HDF5 to netCDF using an NCL program, and then the converter reads the netCDF to create a DART observation sequence file. I'm not a fan of the two-step process, but ...

The _observations/obsconverters/AIRS directory reads JPLs HDF5 files for AIRS observations and requires the HDF-EOS libraries to read the data files.

So ... my question is ... how much HDF5 do you know? Where are you getting the observations? (from a DAAC? do they have an example read program?) You may be in luck if they have not been too tricky with the HDF5 ... since netCDF uses the HDF libraries 'under the hood', some (many, actually) Fortran netCDF read routines actually work on HDF5 files.

12:22 PM (1 hour ago)

Thanks so much! This was all great news!

I have a couple of questions (and clarifications):

I am planning to use CICE within CESM. I’m taking a look at the freeboard algorithm now, and it seems good to go from my perspective for now (though I am new to Fortran, so I will also check with CC).

I did see Yongfei’s comment about converting first to netCDF using the ncl converter, which I am happy to do, though it would help my peace of mind to avoid jumping between dataset formats, for sure. I don’t know very much HDF5, though I understand the structure and idea behind the data format. The ICESat-2 data was downloaded from the NSIDC DAAC, which does provide netCDF 3 and 4 reformatting capabilities and some tutorials on reading HDF5 files, but primarily (only!) in Python, as far as I am aware. From what I can tell, I don’t think the ICESat-2 team was tricky with the HDF5; most of the emphasis I’ve gotten from them is that the format is very similar to netCDF.

I have no objections to the GitHub Issue! I imagine it might be quite helpful, and I would appreciate it.

Just for my own understanding: my conversations with CC have generally implied that we would need to write a new forward operator for real freeboard observation assimilations; is this not the case? Using the cice forward operator with the get_expected_agreg_freeboard() subroutine would do the job? Is it deriving some freeboard from the model state and then updating that in the model with the assimilated observations? If so, would it be possible to adapt that code to assimilate the freeboard and update both sea ice thickness and snow depth directly? Please let me know if I’ve misunderstood!

Thanks so much for your welcome and your help!

timhoar commented 4 years ago

I believe one of us may have misunderstood.

The existing _get_expected_agregfreeboard() routine will convert the model state at any location to the estimate of the freeboard. The only part of the observation that is of any interest to the forward operator is the observation location (and time, naturally).

The reason to have a TYPE of observation is that there are many platforms or instruments that measure something. The model states have (fundamental?) QUANTITIES. The forward operators encounter an observation of a particular TYPE and ask the model for the associated QUANTITY. The association happens in the _obs_def_cicemod.f90 (in this case).

Any observation (real or synthetic) with a TYPE of SAT_SEAICE_AGREG_FREEBOARD will invoke the one and only forward operator routine _get_expected_agregfreeboard() . Now ... if there is more than one satellite that measures freeboard, you could add a new TYPE (say ... ICESAT2_SEAICE_FREEBOARD) and then you would have to add that TYPE to the case statements so that type invokes the _get_expected_agregfreeboard() routine - and several other appropriate routines. You could also have that TYPE invoke some other forward operator code that you would have to supply. That's up to you.

The reason to discriminate among satellites (or any other observing platform) is that you can evaluate one set of observations and assimilate (or ignore) the other observations.

mollymwieringa commented 4 years ago

I think the misunderstanding was mine regarding the role of the forward operator. This makes more sense now, and the distinction between TYPE and QUANTITY is helpful.

Sounds like we're good to go on the QUANTITY front, and can just move forward adding a new ICESAT2 freeboard TYPE.

timhoar commented 4 years ago

Let me know when you have it working and want to issue a pull request.

Tim Hoar Data Assimilation Research Section National Center for Atmospheric Research thoar@ucar.edu 303.497.1708

On Sun, Jun 28, 2020 at 10:42 AM mollymwieringa notifications@github.com wrote:

I think the misunderstanding was mine regarding the role of the forward operator. This makes more sense now, and the distinction between TYPE and QUANTITY is helpful.

Sounds like we're good to go on the QUANTITY front, and can just move forward adding a new ICESAT2 freeboard TYPE.

— You are receiving this because you authored the thread. Reply to this email directly, view it on GitHub https://github.com/NCAR/DART/issues/32#issuecomment-650791590, or unsubscribe https://github.com/notifications/unsubscribe-auth/ADTRYOWN6SQANEBOKERQQCDRY5XIFANCNFSM4OF7AGDA .