Precip units do not match.

[jbusecke]

Raised by Sid on Slack:

The virtual store 'TotalPrecip_Percentiles' has the wrong units (should be kg m-2, but is m3/s).

Open Questions

• Are just the units wrong, or is this pointing to the wrong files?

Log

• I confirmed that all the files have the expected units
• The individual virtual datasets have the correct units (debugging the script)
• So the issue must be in the serialization or happened after the creation of the virtual store.

[jsignell]

Could the issue also be in the concatenation?

[jbusecke]

Ok so I produced a debug version (only one variable, modified script here) that can be accessed like this:

import icechunk
import xarray as xr

storage = icechunk.s3_storage(
    bucket="nasa-veda-scratch",
    prefix=f"jbusecke/RASI/debug/HISTORICAL/", #replace HISTORICAL with SSP245/SSP585 for future scenarios
    from_env=True
)

chunk_url = "s3://nasa-waterinsight/RASI/"
virtual_credentials = icechunk.containers_credentials({
    chunk_url: icechunk.s3_anonymous_credentials()
})

repo = icechunk.Repository.open(
    storage=storage,
    authorize_virtual_chunk_access=virtual_credentials,
)

session = repo.readonly_session('main')
ds = xr.open_zarr(session.store, consolidated=False, zarr_version=3, chunks={})
ds

This seems to have the appropriate units:

Next ill rewrite the full store as a test again.

[jbusecke]

Could the issue also be in the concatenation?

The test above would suggest otherwise, but let me run the full thing and confirm.

[jbusecke]

Ok I have confirmed that this is an issue with the merging logic. Let me see if I can find a quick fix.

[jbusecke]

I wasn't able to get this working today, but will continue tomorrow!

[jbusecke]

Ok I think I got this working. The solution for now was to build a dataset per variable, then manually constructing a new dataset and writing that to icechunk.
My prelim conclusion is that there is something going wrong with the combine logic in xarray, which seems weird to me. Ill try to pin this down a bit more, but will not have too much time for this today.

[jbusecke]

Gotta switch gears for a bit, but Id like to come back to this and produce a small reproducer:

• ideally taking out the virtualization entirely
• Only a few files

Wondering if this is actually a bug that I should raise over at xarray.

[jbusecke]

Ok I noticed another issue with the dataset. The variable 'Streamflow_percentiles' is missing a year (first year in HISTORICAL, and last year in SSP*).

Example file listings for HISTORICAL:

HISTORICAL RZSM_percentiles
120
['nasa-waterinsight/RASI/LSM/HISTORICAL/PERCENTILES/RZSM_percentiles_195001.nc', 'nasa-waterinsight/RASI/LSM/HISTORICAL/PERCENTILES/RZSM_percentiles_195002.nc', 'nasa-waterinsight/RASI/LSM/HISTORICAL/PERCENTILES/RZSM_percentiles_195003.nc']
['nasa-waterinsight/RASI/LSM/HISTORICAL/PERCENTILES/RZSM_percentiles_195910.nc', 'nasa-waterinsight/RASI/LSM/HISTORICAL/PERCENTILES/RZSM_percentiles_195911.nc', 'nasa-waterinsight/RASI/LSM/HISTORICAL/PERCENTILES/RZSM_percentiles_195912.nc']
================================================================================
================================================================================
HISTORICAL Streamflow_percentiles
108
['nasa-waterinsight/RASI/ROUTING/HISTORICAL/PERCENTILES/Streamflow_percentiles_195101.nc', 'nasa-waterinsight/RASI/ROUTING/HISTORICAL/PERCENTILES/Streamflow_percentiles_195102.nc', 'nasa-waterinsight/RASI/ROUTING/HISTORICAL/PERCENTILES/Streamflow_percentiles_195103.nc']
['nasa-waterinsight/RASI/ROUTING/HISTORICAL/PERCENTILES/Streamflow_percentiles_195910.nc', 'nasa-waterinsight/RASI/ROUTING/HISTORICAL/PERCENTILES/Streamflow_percentiles_195911.nc', 'nasa-waterinsight/RASI/ROUTING/HISTORICAL/PERCENTILES/Streamflow_percentiles_195912.nc']

Ill factor this out into a separate issue and work on a fix for now.

[jbusecke]

Ok I just naively tried to do an 'outer' join (filling streamflow data with nans) and this does not work on Manifest arrays (yet, due to limitations on the indexing).

I think the most straightforward fix for now is to handle this at the file search level, omitting the extra years from all variables but 'Streamflow_...', this will be equivalent to doing an 'inner' join.

[jbusecke]

Oh boy there is all kinds of things wrong with the streamflow data.

The netcdf data variable is actually called TotalPrecip_Percentiles!

Virtualizing var='Streamflow_percentiles' # this is parsed from the filename
Single variable virtual Dataset: <xarray.Dataset> Size: 327MB
Dimensions:                  (time: 816, percentile: 5, lat: 152, lon: 132)
Coordinates:
  * percentile               (percentile) int32 20B 10 25 50 75 90
  * lat                      (lat) float32 608B 23.95 24.05 ... 38.95 39.05
  * lon                      (lon) float32 528B -109.5 -109.5 ... -96.55 -96.45
  * time                     (time) datetime64[ns] 7kB 2015-01-01 ... 2082-12-01
Data variables:
    TotalPrecip_Percentiles  (time, percentile, lat, lon) float32 327MB Manif...

This now totally explains why the units were mangled! I think we should for now not virtualize Streamflow_percentiles until we have clarity in #7? What do you all think?

EDIT: I have confirmed these are not matching both for HISTORICAL and SSP245