Introduce "gadf" format for RegionNDMap

[adonath]

Description

This pull request introduces a proof of concept for a "gadf" format for RegionNDMap. It serialises as following:

• Main HDU as an ImageHDU with the data, but no WCS info (to be defined / documented)
• Region HDU as an BinTableHDU in the ogip like format (see https://fits.gsfc.nasa.gov/registry/region.html)
• A bands HDU equivalent to what we have for normal maps (see https://gamma-astro-data-formats.readthedocs.io/en/latest/skymaps/index.html#bands-hdu)

The main question here is whether this is needed at all for now. The main use case I hand in mind was to be able to serialise a RegionNDMap with multiple axes. However this currently only applies to the EDispKernelMap and PSFMap, for both we can fall back to different solutions: they can be serialised as a WCS map with one spatial bin.

Beside the extra axes the serialisation will also keep all the region information on the maps, which might be interesting for provenance and it provides a simple solution to serialise both SpectrumDataset as well as SpectrumDatasetOnOff in single FITS file, which is probably a convenient alternative to the current ogip format.

The goal here is not to define any official DL4 data format, but just provide a simple and convenient way to serialise a RegionNDMap of any kind. The format will not be optimised by any means, but contain all the info needed to allow for round-tripping of a RegionNDMap, without loss of information.

If we decide this is a useful addition to our serialisation functionality, the only remaining thing is to discuss the structure, that is used to store the data. Beside the ImageHDU, with no associated WCS, it could be an BinTableHDU as well, similar to the IRF DL3 formats.

[codecov]

Codecov Report

Merging #3231 (376eaeb) into master (ded60ed) will decrease coverage by 0.00%.
The diff coverage is 94.62%.

@@            Coverage Diff             @@
##           master    #3231      +/-   ##
==========================================
- Coverage   93.91%   93.90%   -0.01%     
==========================================
  Files         144      144              
  Lines       17450    17478      +28     
==========================================
+ Hits        16388    16413      +25     
- Misses       1062     1065       +3     

Impacted Files	Coverage Δ
gammapy/maps/core.py	86.92% <50.00%> (-0.17%)	⬇️
gammapy/maps/regionnd.py	92.22% <91.89%> (+<0.01%)	⬆️
gammapy/maps/geom.py	91.05% <96.77%> (-0.07%)	⬇️
gammapy/datasets/io.py	95.10% <100.00%> (ø)
gammapy/irf/core.py	91.10% <100.00%> (ø)
gammapy/maps/region.py	90.57% <100.00%> (+0.28%)	⬆️

---

Continue to review full report at Codecov.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update ded60ed...376eaeb. Read the comment docs.

[registerrier]

Thanks @adonath for starting this.

This pull request introduces a proof of concept for a "gadf" format for RegionNDMap. It serialises as following:

Technically, we need a specific format for map serialization but it should decoupled from data level issues. I am not sure it is relevant to have Map serialization as part of higher level data. After all, this is really basic FITS convention on how to store multidimensional quantities that is not fully supported by WCS.

• Main HDU as an ImageHDU with the data, but no WCS info (to be defined / documented)

How about a simple table HDU?

• Region HDU as an BinTableHDU in the ogip like format (see https://fits.gsfc.nasa.gov/registry/region.html)

Fine. As a side note, the region convention is a regular FITS convention, it is not connected to ogip, no?

• A bands HDU equivalent to what we have for normal maps (see https://gamma-astro-data-formats.readthedocs.io/en/latest/skymaps/index.html#bands-hdu)

OK

The main question here is whether this is needed at all for now. The main use case I hand in mind was to be able to serialise a RegionNDMap with multiple axes. However this currently only applies to the EDispKernelMap and PSFMap, for both we can fall back to different solutions: they can be serialised as a WCS map with one spatial bin.

I think we need it. While OGIP format can support the 1D spectral analysis. It is not fully compliant with all use cases (e.g. what if I want to serialize a SpectrumDataset with a FoVBackgroundModel).

Beside the extra axes the serialisation will also keep all the region information on the maps, which might be interesting for provenance and it provides a simple solution to serialise both SpectrumDataset as well as SpectrumDatasetOnOff in single FITS file, which is probably a convenient alternative to the current ogip format.

+1

The goal here is not to define any official DL4 data format, but just provide a simple and convenient way to serialise a RegionNDMap of any kind. The format will not be optimised by any means, but contain all the info needed to allow for round-tripping of a RegionNDMap, without loss of information.

Indeed. And again something like WcsNDMap or ``RegionNDMap` does not really belong to DL4. It is not describing a data format but rather a possible FITS convention.

If we decide this is a useful addition to our serialisation functionality, the only remaining thing is to discuss the structure, that is used to store the data. Beside the ImageHDU, with no associated WCS, it could be an BinTableHDU as well, similar to the IRF DL3 formats.

[adonath]

Thanks for the feedback @registerrier, I'll explore the use of a BinTableHDU for the data instead...

[adonath]

@registerrier Following your suggestion I changed to use a BinTableHDU and decided to merge the bands HDU with the data hdu (similar to the ogip specresp). The table has now the following structure:

I think together with the region HDU this is a very clean structure now. I'll start adding more tests and documentation now.

[registerrier]

Thanks @adonath . This looks good!
The proposed format is nice. I just wonder whether we should call it gadf for now.
I have left a few comments inline.

[registerrier]

Should we call it gadf, if it is still internal?

[adonath]

See my comment below...

[adonath]

Choosing a different name for this format is possible but introduces the problem of defining multiple formats in MapDataset.write(). Consider the case where dataset.psf is a RegionGeom (assuming we want to support this) and one uses MapDataset.write(format="gadf"). I guess we don't really want to introduce something like ,write(format={"counts": "gadf", "background": "fgst-ccube", "psf": "gp-region"}), this would make the number of possible single fits file definitions for a MapDataset explode...

[registerrier]

Technically you should also have the ogip-rmf for the Region EdispKernel no?

[adonath]

Yes, this is true. My proposal would be that I enable the format for IRF maps such as PSFMap and EDispKernelMap in a follow up PR.

[registerrier]

Shouldn't this be tested?

[adonath]

Indeed, I'll add a test...

[adonath]

Thanks a lot @registerrier, I have made a request here: open-gamma-ray-astro/gamma-astro-data-formats#170, if it is accepted, I think we can keep the "gadf" name, otherwise I'll choose a different one.

[adonath]

Following a suggestion by @registerrier in the developer call, I now changed to introducing a dedicated bin table HDU for the data as well. This allows to share a common bands HDU if multiple "spectra" are stored within the same FITS file (see e.g. serialisation of #3075)

[adonath]

@registerrier In order to make progress, I'll take the risk and merge this PR as is. The PR in the gammapy astro data format repo is open and any user can go there to look at the definition. I assume it will be merged in a some version close to the proposed one. If we decide the definition of the FITS conventions is not part of GADF, I guess we have to adapt the tag in Gammapy consistently anyway.