Sentinel-1 VV VH backscatter DA

[KUL-RSDA]

Description

At KU Leuven (Michel Bechtold, Alex Gruber, Gabrielle De Lannoy) and CNR IRPI (Sara Mondanesi, Christian Massari) we developed plugins for a Sentinel-1 backscatter data assimilation with the water cloud model as RTM.
There are two DA cases included for NoahMP 3.6:

• VV backscatter for soil moisture updating, and
• VV&VH backscatter for soil moisture and LAI updating

The PR includes bug fixes related to

• the use of multiple observation types which we needed for the joint VV and VH assimilation issue (Wrong use of DA instance index k in write out of spread and incr. #1012)
  We found this bug when assimilating VV and VH in two different instances (the approach of using two separate instances we eventually replaced by using the multiple observation type option).

The PR also includes a bug fix we encountered in the map_utils that occurs due to longitude 0/180 for the latlong projection. We adopted the correct approach that was already present in llij_gauss.

First paper making use of the Sentinel-1 data assimilation:
https://hess.copernicus.org/preprints/hess-2022-61/

Resolves #1192
Resolves #1012

Testcase

LISF/lis/testcases/dataassim/S1_backscatter

[emkemp]

We will begin working on this pull request. Two issues that need particular attention before merging:

• Need to confirm the bug fixes are acceptable
• Need to better handle netCDF file in test case (we prefer not to include binary files in GitHub)

[cmclaug2]

Here is an update on this pull request:

• We have extracted the testcase.
• We have run the testcase and confirmed we can produce output.
• We are still reviewing the code to make sure all changes are acceptable, specifically the compilation and interpolation routines.
• We need to receive example output from you to compare our results.

Do you have a preferred way to share data with us? Perhaps you could push a temporary repository to github that we can pull from? Or, if you have access to a shared data portal on discover, you could save the output there.

Thanks

[cmclaug2]

Hi developers,

Please see the message above. We still need to compare our output to yours. Thank you!

[Rajainoubli]

Hello everyone,
I used the wcm to obtain simulated backscatter and I need to invert the model to obtain soil moisture using the ANN. Can you please help me to know the required steps and inputs to achieve my goal. Thanks in advance for your help.

[cmclaug2]

@KUL-RSDA Just a nudge! We have made progress reviewing this pull request. We will need example output for your testcase.

@Rajainoubli Thank you for your message. If this question is best answered by the developers of this code, I suggest reaching out to them outside of this thread.

[mbechtold]

@KUL-RSDA Just a nudge! We have made progress reviewing this pull request. We will need example output for your testcase.

@Rajainoubli Thank you for your message. If this question is best answered by the developers of this code, I suggest reaching out to them outside of this thread.

My apologies for the delay in our response! Unfortunately, we did not receive a notification via our KUL-RSDA github account. We will reply asap after next week's EGU.

[cmclaug2]

Hi @mbechtold,

We'd love to work towards merging your contribution. Once we have example output from your testcase we can take the next step.

[mbechtold]

Hi @cmclaug2
so sorry for not getting back to you earlier! Part of our HPC cluster was down and I couldn't actually access all files.
I created a temporary github repo with the reference output for the testcase:
https://github.com/mbechtold/S1_LIS_PR_output_temporary
Let me know if you need more from me. Directly address a comment to @mbechtold and not to @KUL-RSDA since we have problems in getting forwarded github messages with our group account email address. Thx!

[cmclaug2]

@mbechtold

Thank you for creating that temporary directory! I ran your testcases again and compared our output. I am noticing some differences, and I was wondering if you could double check something for me -

In your lis.config file on lines 152-155, you point to the nccs data portal noahmp_params folder. I notice from a note in noah_README.txt that this folder has since been deleted, but it pointed to noahmp36_params so it was no issue to just edit the config file.

However, line 152 of your config file points to a nonexistent noahmp36 land use parameter table, VEGPARM_UMD.TBL. In its place, I ran your testcase with the existing parameter table VEGPARM.TBL, but I am noticing differences. I do not have enough background on these tables to know what happened to VEGPARM_UMD.TBL or if it even had differences compared to VEGPARM.TBL

In short, can you please rerun your testcase but edit line 152 of your config file to the new VEGPARM table? Or can you help point me to the VEGPARM table you are using?

Thanks so much for helping us review this PR.

[mbechtold]

@cmclaug2
I am not exactly sure what you mean by 'new VEGPARM.TBL'.
The VEGPARM_UMD.TBL we are using is in line with what I can find here: https://portal.nccs.nasa.gov/lisdata_pub/data/PARAMETERS/noahmp36_parms/noahmp_params/
I also putshed it on https://github.com/mbechtold/S1_LIS_PR_output_temporary
If you let me know which VEGPARAM.TBL I should use, I'm also fine to run the test case again!

[cmclaug2]

@mbechtold,

Thanks for the response! The issue was on my end. I was having trouble connecting to the data portal. All is well and I ran your testcase so we can move on to the next step.

I am noticing differences in your output. We are going to troubleshoot that on our end and see what we can find. In the meantime, it would not hurt if you rerun your testcase. Can you also push /lis/make/configure.lis? This will help make sure we are running your testcase using the same configuration settings as you.

[mbechtold]

@cmclaug2
I added configure.lis
and also the other input TBL files I used.
For finding the difference it might be also useful to directly run with my
lis_input.d01.nc which I now also put. Then we can exclude that it is related to ldt.
Reg. the rerun, what do you want me to change? So far I didn't rerun anything and only provided you all files I used for the run that I put on github.

[cmclaug2]

Hi @mbechtold,

Sorry for the delay. We did some investigation and cannot seem to pinpoint where the differences are coming from. In all variables (expect TotalPrecip_tavg and Water_tableD_tave) we are noticing differences across both time steps. Most of these differences are statistically insignificant and may be attributed to using different compilers. These are not a concern.

However, a few of the differences caught out attention. LAI_inst, SoilMoist_tavg, and SoilTemp_tavg all seem to have differences that are approaching statistical significance. For context, here are two plots of SoilTemp_tavg at the second time step. Plot 1 shows your result. Plot 2 shows the difference between your results and our results.

Plot 1) Your Data

Plot 2) Difference

Here are a few ideas to see what is going on:

1. Can you rerun your test double check that all of the updates in the PR were applied when you test was run? Also, I think it could be a good idea to lengthen your testcase to a week. This way we can see if the differences grow over time.

2. If possible, can you push all of the input data you used for this testcase to your https://github.com/mbechtold/S1_LIS_PR_output_temporary github directory? I just want to make sure changes were not made to the data on our data portal that could have occurred after you grabbed it.

3. I am not an expert on the topics in this PR. I am wondering what you think about these differences. Do they appear significant to you? We generally want to make sure we can fully replicate results before we add outside contributions. I think we will know more after we test over a longer time period to see if the differences grow.

Thank you in advance.

[mbechtold]

Hi @cmclaug2
I finally started with your idea 1. Could you compare your output against the new reference output that I pushed to
https://github.com/mbechtold/S1_LIS_PR_output_temporary
Note: If it's now consistent with your results then sth went wrong on my end. I feel already very sorry for it. It seems my testcase run linked not to the Sentinel-1 reference data of the PR but to another Sentinel-1 dataset with a different processing. I now took the reduced testcase S1 dataset of the PR and the output is slightly different in the direction of what you observed. Let me know if it's true. Otherwise, we should continue with your idea 2. On idea 3., I think the differences are too large and not acceptable, so we need to figure out the reason.
Thanks, Michel

[cmclaug2]

@mbechtold,

Thank you for rerunning your testcase. The results in the SURFACEMODEL/ folder now differ from mine on an order of magnitude consistent with roundoff errors. That is good!

However, I am noticing differences on the order of 5%-10% in the output in the RTM/ folder. I am going to discuss this with @jvgeiger next week and let you know what we think about this.

Also, jotting this down for a reminder to ourselves - before merging, we will need to do the following on our end:

• Extract the testcase
• Add WCM documentation to user_cfg_table.adoc and build.adoc.

[mbechtold]

Hi @cmclaug2
thanks, good to hear about the better match ... except for the RTM. This is difficult to explain for me since the rest seems fine. I just checked but I am sure I pushed the right files. So you don't see any difference for DAOBS and also not for the innov files, but the RTM output is inconsitent? This shouldn't happen since one can actually calculate the right RTM output from DAOBS and innov.
I am on holidays for two weeks from tomorrow night onward. Will take care of the documentation afterwards.

[cmclaug2]

@mbechtold,

The DAOBS files are identical, but I do see small differences in the innov files. We will investigate more next week and hope to have a plan when you return. As of now, we will take care of updating the documentation.

If you have time to push all your input data, that would be very helpful for us.

Have a great holiday!

[mbechtold]

Hi @cmclaug2
I put everything I used for the LIS run (including the lis_input.d01.nc) on
https://github.com/mbechtold/S1_LIS_PR_output_temporary
Talk to you soon! Hope the remaining difference can be sorted out soon. Thx!

[cmclaug2]

Update - I am able to reproduce your sample output for your testcase.

[cmclaug2]

Hi @mbechtold,

Summary of where we are:

• I reproduced your output within an acceptable tolerance.
• We are working on a rebased branch of your PR where we have extracted the testcase and performed clean up in the files lis/arch/Config.pl and lis/make/default.cfg.
• We appreciate the documentation of your bug fix Wrong use of DA instance index k in write out of spread and incr. #1012. Thank you! We are performing a final check that the bug fix is acceptable with a final review of lis/dataassim/algorithm/enkf/enkf_Mod.F90.
• We appreciate your implementation of more robust code for the latlon projection in lis/interp/map_utils.F90.

Questions:

1. In lis/core/LIS_historyMod.F90, two checks for NaNs were added on lines 8723 and 8857. Is there a reason this check for NaNs was added? I have run your code without the NaN checks and found identical results.
2. For the writing of the innovation files in lis/core/LIS_DAobservationsMod.F90, indexing was changed on lines 2047-2053. Why was this necessary for multiple observation types?

[KUL-RSDA]

Hi @cmclaug2,

thanks, great to hear about your updates! I am glad that you understood and appreciated the two general bug fixes.

On 1.: This check for NaN seems to be obsolete indeed. We were not immediately able to reconstruct why this was once added. It goes back to the initial phase of developing the WCM plugin. KUL-RSDA@47e5278 We will run one more test at our end to be sure that it can be removed.

On 2.: This indexing change was for sure necessary. The use of multiple observation types was not functional in its previous implementation. I assume this affected all DA cases in LIS with multiple observation types. How much was this feature used before? I will look asap into the debugger when running the 'VVVH' (number of obs types = 2) case to give you a more detailed reply on how the indexing was exactly wrong originally. I don't recall right now.

One update from our side. We were running for the first time the Sentinel-1 DA over Australia and encountered a problem with the Sentinel-1 observation reader due to the file reading on a day-by-day basis. Some change in handling the local time is still needed.

I have time to work on all three issues in about a week from now and let you know about the progress.

[cmclaug2]

@KUL-RSDA @mbechtold,

At long last we are ready to merge this PR. Since we needed to rebase this PR to extract the testcase, here is what will happen.

1. I will open a new PR (Sentinel-1 VV VH backscatter DA (REWORKED) #1415) that is identical to this one except:

   • Does not contain the testcase.
   • Removes the NaN check from LIS_historyMod.F90.
   • Does not contain your updates to Config.pl. (lis.config.adoc will specifiy the need to configure with CRTM and CMEM. WCM will no longer be a configuration setting but will remain ON/TRUE in default.cfg)
   • Adds documentation updates to user_cfg_tables.adoc, default.cfg, and lis.config.adoc.

2. I will follow up on issue Wrong use of DA instance index k in write out of spread and incr. #1012. We appreciate your bug fix to enkf_Mod.F90. We will eventually need to propagate that bug fix through the other DA algorithms.

3. When you fix the newly identified bug related to the Sentinel-1 observation reader handles local time, submit a new pull request. In the meantime, please open an issue so we can track the problem.

4. When the new PR is opened and merged in, this PR will be closed.