Repacking

[wachsylon]

Hi,

I just learned "repacking" should be integrated into CMOR which would be great in a sense, since, imo, noone would understand why we need two tools for data standardization and publication.

Is there a plan when that can be realized? We may wait until then with CMIP7 publication.

Best,
Fabi

[matthew-mizielinski]

For context here, there will be some communication from the WIP about tools to optimise packing of data in netCDF files to optimise speed of extraction.

If I recall correctly then the changes made to support this are;

• single chunk for time points (CMOR has a default chunksize of 512)
• single chunk for time bounds (I think this needs looking at in CMOR)
• chunking data into 4 MB blocks.

I'm a little hazy about the details as to how the repacking tool works, so @davidhassell would you be able to discuss with @mauzey1 and @durack1 the details that would allow CMOR implement this change or get as close to it as possible?

[davidhassell]

Yes! Just a bit a brain-dump for now:

Doing things exactly like cmip7repack is necessary - a repacked file just has to pass the tests that ESGF will apply to it, namely:

1. the /time variable has one chunk
2. the time bounds variable has one chunk
3. the data variable has chunks of "at least" 4 MiB (the slight vagary here is that a chunk size under 4 MiB is OK if increasing the chunks by one element along its leading dimension would put it over 4 MiB)
4. For each variable, its first data chunk must appear later in the file that that variables last B-tree (internal metadata) fragment.

cmip7repack achieves number 4 by moving all of the B-tree fragments (for every variable) to the top of the file, which is optimal, but clumping the B-tree fragments together on a per-variable basis, wherever in the file, is still fine. This is a typical h5repack command built and run by cmip7repack:

h5repack --metadata_block_size=236570
    -f /time:SHUF -f /time:GZIP=4 -f /time:FLET -l /time:CHUNK=1800
    -f /time_bnds:SHUF -f /time_bnds:GZIP=4 -f /time_bnds:FLET -l /time_bnds:CHUNK=1800x2
    -f /pr:SHUF -f /pr:GZIP=4 -f /pr:FLET -l /pr:CHUNK=37x144x192
    file.nc file_repacked.nc

(It is possible for number 4 to be true and for the file to still be badly packed, but that to make a file like this would obtuse and hard to do we can discount it as a problem!)

[wachsylon]

@davidhassell @matthew-mizielinski

Is there any information why these 4 points are really necessary? For me, 1. and 2. sound reasonable, but also not at all like a reason to reject data in CMIP.

Isnt the optimal size of a chunk depending on the compressing algorithm? Note that CMOR is also allowing and supporting zstandard which is to be much faster than zlib. Looking at the h5repack tool, I am not sure if that is supported at all. Please do not force data producers to use an old compressor.

Also, the question that always come regarding rechunking: how can you guarantee that the repacked chunk shape is optimal for the most important use case of the data? Maybe for a specific experiment, one would more often analyse the full time series for a given coordinate point, while for others one would analyse a full record for a given time step. If model producers knew the most often use case and did create chunk shapes on purpose (unlikely, but who knows) which in that scenario would not allow to create chunks >4MiB, why would you take the risk forcing a potential non-optimal rechunking?

Not only, but especially the B-tree thingy sounds like sth very hdf5 or NC internal. Why is this B-tree fragmentation not the default in these formats?

If 4. is a necessary condition to have performant range request at all, it would be way easier to force ESGF nodes to create virtual datasets accessible through Zarr instead of forcing data producers into a legacy format. I bet that this will be done by many ESGF nodes anyway, or at least cloud providers or Pangeo supporters. From my experience, chunked index coordinates are not an issue for kerchunked based data access from ESGF.

[davidhassell]

Hello @wachsylon,

Thanks for your questions and feedback. The purpose of all this is to allow CMIP7 to be accessible via direct remote access. For this, parsing the netCDF file requires a range get for every internal file metadata fragment (sometimes called B-tree fragments) , and every accessing every data chunk also requires a range get. There are CMIP datasets with 1 million chunks, which makes them entirely unusable for direct remote access. The repacking process reduces the number of chunks, and hence B-trees, and reduces the number of range gets needed to access the data. The internal B-tree structure of a file is related directly to the chunk layout, and the order in which chunks were written to the file, so is more a function of use, rather than a library deficiency (although there are techniques, such as "paging", that can be used at original file creation time to improve things).

I agree that there is no one-size-fits-all optimal chunking strategy for the data variable, but there are most definitely poor chunking strategies. Much of the CMIP6 archive was chunked into one chunk per time step, which is bad every time. The strategy of chunking along the time axis to create chunks of at least 4 MiB (or larger) chunks is our compromise. The ESGF test will be that data variable chunks are no smaller than 4MiB, but they may be larger if the user wants. If the user has already chunked their data into > 4 MiB chunks, the cmip7repack tools will not change it.

The only tool we are aware of for easily changing the internal structure of the dataset is h5repack, and our output files have to be netCDF4, so we can only use compression filters provided by base libhdf5 that are supported by netCDF. So unfortunately zstandard is off limits. Any suggestions welcome, on this front.

netCDF is not a legacy format (a vast number of petabytes are about to be created, afterall :)). A Zarr copy of a netCDF dataset often appears to perform better, but that is usually because it is chunked differently (i.e. better) than the original netCDF. Well-chunked (with nice B-trees) netCDF performs as well as Zarr. This repacking procedure makes netCDF "cloud optimised", and it will perform as well as Zarr. A virtual Zarr view of badly a formed netCDF file is itself poorly performing (a kerchunk file with 1 million entries is not fun!), which is why actual copies were made in (e.g. pangeo) - so they could be rechunked. Given that CMIP7 datasets will be already cloud optimised, because of this repacking, I wonder if there will be the appetite for creating Zarr view or copies? We shall see :) (Note that thread-safe parallel reads to HDF5 files in python are now possible with pyfive (which is also very soon to be a new backend option to h5netcdf).

All the best,
David

[wachsylon]

Thanks David for the quick and helpful response. I understand the benefits of good chunking and would like to see the 4 points built in in CMOR.

If however there is no way around a second tool, I suggest to move at least 2 things into that new package which was always an issue between data producers and data publication in CMIP6:

• PID creation: CMOR creates, with the help of the prefix for each project, a PID-like tracking_id. The issue is: this is only a valid PID if it is registered. CMOR on the other hand is just a tool producing files with such ids that may or may not be registered. PID manager do not like the idea of files laying around with PIDs which are not PIDs. @TobiasWeigel

• version directory:
  The version directory is created by the data producer. The rule is:

  approximate date when files were created

  which lead to problems in CMIP6 as there were multiple directories for the same dataset. Aside from this, it would make sense to rather have a directory date that represents the date of publication. Between creation and publication, there can be months or even years of validation. So an old creation version directory looks like an outdated dataset however it can brand new. That could be corrected with a tool that is actually applied by the data publishers.

[valeriupredoi]

A virtual Zarr view of badly a formed netCDF file is itself poorly performing (a kerchunk file with 1 million entries is not fun!), which is why actual copies were made in (e.g. pangeo) - so they could be rechunked.

This is an excellent point @davidhassell ! Also, hi Matt @matthew-mizielinski 🍺

[valeriupredoi]

just a quick note, we now moved the repacking tool to its permanent location https://github.com/NCAS-CMS/cmip7repack

[valeriupredoi]

Hi @wachsylon too - sorry, I've not remembered your last name from EGU 😃 We can dissect Zarr over a pint or two at EGU BTW 🍺

[wachsylon]

The purpose of all this is to allow CMIP7 to be accessible via direct remote access

As explained, there are multiple ways to achieve that. One could also have considered allowing Zarr as a format, which is compatible through NCZarr with NetCDF afterall. You want to make it accessible with netCDF

parsing the netCDF file requires a range get for every internal file metadata fragment (sometimes called B-tree fragments)

With a virtual dataset approach, one would do that once and store the information outside of the file. And one will do it before publication. How would netcdf realize that? Would people again and again query the metadata information?

There are CMIP datasets with 1 million chunks, which makes them entirely unusable for direct remote access

Not true. Have a look on lazyReferenceMapper. The kerchunk-parquet format allows users to only open the chunk tables they need. If parquet tables are small, indexing chunks is not an issue.

Then, it is just a question whether the infrastructure holds the chunk requests. From my experience, large amounts of small requests is not an issue but rather any request to large chunks because these will be buffered or stored in memory on the server.

While you can ignore all my arguments up to now, you may want to set a maximum chunk size as well: 49.9MB. That is because browsers only cache chunk sizes <50MB per default.

That being said, even 10MB chunks would mean for users that use a default python http client, they would request 100*10MB chunks at once - which may result in large requests queues and up to 1GB memory usage on the server.

which is bad every time

No. What about creating a plot for the entire globe?

is our compromise

A compromise that ... which community agreed on?

So unfortunately zstandard is off limits

Imo, we would need proof that improved range request is more important than flexiblity in the compression algorithm before we unable CMOR again to use zstandard and force CMIP into repacking. For me, it is hard to believe these exist. I guess most CMIP users work on HPC for which data centers will repliate data once and who need high-performant data access. This is not possible with zlib anymore.

netCDF is not a legacy format

sorry, that was an unnecessary shot fired :) we use it a lot too, of course.

(a kerchunk file with 1 million entries is not fun!)

No, it is indeed fun if you use the right output format. Either parquet, or just recently Icechunk.

I wonder if there will be the appetite for creating Zarr view or copies?

The real magic of virtual datasets is not really the cloud-optimized storage. It is the high aggregation. Imagine users who analyse an ensemble of the same model given by multiple members of an experiment. They would need to aggregate these ensemble members on a dimension to work nicely with the data. This exact preparation step can be made with virtualization and it turns out to be usefull for everyone because cutting is easier than glueing. Finding these ensemble members, validate that they can be merged on an ensemble dimension - these tasks are coastful and would have to be done for each ensemble analysist. In contrast, a single sel on a high aggregation cost nothing.

And that were just the advantages of users that you also put into foreground. For data providers, we get another inherent quality assurance because virtual datasets can only be created if aggregation is possible which is only possible if the data is consistent across ensemble members. We can also edit the json metadata file of the virtual dataset and enrich it - e.g. with STAC metadata.

So for me, yes, there will be virtual datasets. No matter whether they reference data on local storage or cloud storage or somewhere else. Just have to be created once - afterwards, it is a string replace.

[florianziemen]

Hmm, I don't really see the conflict between using zstd and netCDF, but maybe that's because I consider netcdf as something netcdf-c can open. In principle, netcdf-c can simply use blosc or a plain zstd (I'd go with blosc), as long as you provide the codecs as plugins. netcdf-c can also handle plain zarr2, so I'll say hdf5 for netcdf stored as hdf5 in the following.

Regarding the chunk sizes:
Isn't the good thing about putting the data in hdf5 that you can put many chunks into a file, and thus can use small chunks? We have built massive kerchunk tables with parquet without problems, the main problem we have is handling datasets that have big chunks that don't perfectly align with the users' interests, and thus require reading way too much bycatch (cheers to GRIB).

What one could do to defuse the zarr / hdf5 debate (and I think this could be useful, independent of the debate), is to make the datasets conform to both standards:
hdf5 allows you to attach random data at the end of your files - what's not listed in the headers / trees is not seen. Zarr3 sharding usually has the shard info at the back of a file, so what you could do is to attach zarr shard tables at the ends of the hdf5 files (and provide a zarr.json file and stores for the dimension variables, and matching directory hierarchy - if need be via links). This would allow friends of zarr to use the data as zarr, and friends of hdf5 to use it as hdf5. The alias/link game does not really work with object stores, but one could write a pretty minimal web-server that provides a zarr view with redirects to the object store.

libnetcdf currenly is a bit slow in keeping up with the zarr developments (a bit more staff would probably help there), so pure zarr3 with sharding via https not yet an option there.

[florianziemen]

Ah, sorry, the double selling probably won't work as it (currently) requires the chunks / shards to be equal shape, and that conflicts with files based on months / years with changing lengths. Then you'd have to do hdf5 files for things like 30 / 360 / ... days, and that's inconvenient for users that use file names and copy files around.

Rectangular chunk grids are progressing though, so the kerchunk variant could fly in the mid-term future.

[durack1]

@cofinoa just tagging you on this issue - considering previous tweaks captured in #601 - would be great to get your eyes on the incoming PR when it lands

[mauzey1]

Update on changes for repacking

I have made a change to the time and time bounds chunking where if the time axis length is defined before writing to the file, the chunk size of time will be the length of the axis with the bounds having the chunk dimensions axis_length x 2. If the time axis is not defined before writing data, such as when you are writing one time step at a time, then a default chunk size of 512 will be used for time with chunk dimensions of 512 x 2 for the bounds. Relevant code section below.

cmor/Src/cmor.c

Lines 5559 to 5599 in 52c936c

	// time and time_bnds chunking
	for (i = 0; i < pVar->ndims; i++) {
	if (cmor_axes[pVar->axes_ids[i]].axis == 'T') {
	// Set chunking size of the time axis to the size of the
	// time axis if it is greater than zero. Otherwise, set it
	// to the default time chunking size.
	size_t time_chunk_size = CMOR_DEFAULT_TIME_CHUNK_SIZE;
	if (cmor_axes[pVar->axes_ids[i]].length > 0) {
	time_chunk_size = cmor_axes[pVar->axes_ids[i]].length;
	}
	ierr = nc_def_var_chunking(ncid, nc_vars[i], NC_CHUNKED, &time_chunk_size);
	if (ierr != NC_NOERR) {
	cmor_handle_error_var_variadic(
	"NetCDF Error (%i: %s) defining time coordinate chunking "
	"for variable '%s'",
	CMOR_NORMAL, var_id,
	ierr, nc_strerror(ierr), pVar->id);
	}
	// Set chunking size of the time bounds axis with the first
	// dimension set to time's chunking size and the second dimension
	// set to 2.
	if (nc_bnds_vars[i] != -1) {
	size_t time_bnds_chunk_sizes[2] = {time_chunk_size, 2};
	ierr = nc_def_var_chunking(ncafid, nc_bnds_vars[i],
	NC_CHUNKED, time_bnds_chunk_sizes);
	if (ierr != NC_NOERR) {
	cmor_handle_error_var_variadic(
	"NetCDF Error (%i: %s) defining time_bnds chunking "
	"for variable '%s'",
	CMOR_NORMAL, var_id,
	ierr, nc_strerror(ierr), pVar->id);
	}
	}
	break;
	}
	}

The current chunking strategy for the data variable is to use a chunk dimension of 1 time step and the lengths of the other axes. I plan to make a chunking strategy where chunks of ~4 MB will be used if the total size of the data variable is >4 MB and the size of each time step of the data is <4 MB.

Heads up: I won't be available next week but I will return to my normal schedule on December 1.

[mauzey1]

Should we also chunk the other axis variables (latitude, longitude, etc.) rather than have them be contiguous?

Currently, CMOR will make these axes contiguous but will make their bounds, like lat_bnds and lon_bnds, chunked by the size of the axis x 2.