Wrap GMT's standard data type GMT_IMAGE for images

[seisman]

Description of proposed changes

This PR is a subset of PR #3128. In this PR, we focus on wrapping the GMT_IMAGE structure and will implement the to_dataarray method in the future.

Reminders

• Run make format and make check to make sure the code follows the style guide.
• Add tests for new features or tests that would have caught the bug that you're fixing.
• Add new public functions/methods/classes to doc/api/index.rst.
• Write detailed docstrings for all functions/methods.
• If wrapping a new module, open a 'Wrap new GMT module' issue and submit reasonably-sized PRs.
• If adding new functionality, add an example to docstrings or tutorials.
• Use underscores (not hyphens) in names of Python files and directories.

Slash Commands

You can write slash commands (/command) in the first line of a comment to perform
specific operations. Supported slash command is:

• /format: automatically format and lint the code

[weiji14]

Getting this error at https://github.com/GenericMappingTools/pygmt/actions/runs/10017032324/job/27690729251?pr=3338#step:8:1556:

ERROR 4: `/home/runner/.gmt/server/earth/earth_relief/earth_relief_01d_p.grd' not recognized as being in a supported file format. It could have been recognized by driver HDF5, but plugin gdal_HDF5.so is not available in your installation. You may install it with 'conda install -c conda-forge libgdal-hdf5'
Error: ession [ERROR]: Unable to open earth_relief_01d_p.grd.
Error: ession [ERROR]: ERROR reading image with gdalread.
pygmt-session (gmtapi_import_image): Could not read from file [earth_relief_01d_p.grd]
[Session pygmt-session (198)]: Error returned from GMT API: GMT_IMAGE_READ_ERROR (22)
[Session pygmt-session (198)]: Error returned from GMT API: GMT_IMAGE_READ_ERROR (22)

So do we need to add libgdal-hdf5 to gmt-feedstock, similar to conda-forge/gmt-feedstock#291?

[seisman]

I think yes

[seisman]

Working in conda-forge/gmt-feedstock#293.

Edit:

Tried to run gmt grdcut @earth_relief_01d -R-10/-9/3/5 -Greliefcut.nc in conda-forge/gmt-feedstock@ff8d995 and it worked. It means reading the netCDF file usually doesn't need the HDF5 library, but here we're actually reading the netCDF file as an image (kind="image" which reads the file into a GMT_IMAGE data container), which calls the GDAL library to read rather than the netCDF library.

[seisman]

This test still fails even if libgdal-hdf5 is installed, see https://github.com/GenericMappingTools/pygmt/actions/runs/10018688146/job/27694518606?pr=3338.

The expected header.wesn is [-180, 180, -90, 90], but the actual results are platform-dependent. For Python 3.12 + Ubuntu, the result is [-179.5, 179.5, -89.5, 89.5] but for Python 3.10 + Ubuntu, the result is [0.5, 359.5, 0.5, 179.5].

    def test_clib_read_data_image_actual_grid():
        """
        Test the Session.read_data method for image, but actually the file is a grid.
        """
        with Session() as lib:
            data_ptr = lib.read_data(
                "@earth_relief_01d_p", kind="image", mode="GMT_CONTAINER_ONLY"
            )
            image = data_ptr.contents
            header = image.header.contents
            assert header.n_rows == 180
            assert header.n_columns == 360
>           assert header.wesn[:] == [-179.5, 179.5, -89.5, 89.5]
E           assert [0.5, 359.5, 0.5, 179.5] == [-179.5, 179.5, -89.5, 89.5]

The test added in conda-forge/gmt-feedstock#293 also fails (see https://github.com/conda-forge/gmt-feedstock/pull/293/checks?check_run_id=27696860313, https://dev.azure.com/conda-forge/feedstock-builds/_build/results?buildId=983672&view=logs&jobId=656edd35-690f-5c53-9ba3-09c10d0bea97&j=656edd35-690f-5c53-9ba3-09c10d0bea97&t=986b1512-c876-5f92-0d81-ba851554a0a3).

Not interested in digging down the rabbit hole, so I prefer to removing this test.

[weiji14]

Hmm, or we can keep the test and mark it with xfail? Might have something to do with different GDAL versions or something. I almost think we should include GDAL in pygmt.show_versions() to more easily check things.

[seisman]

I've added the test back but the macOS + Python 3.12 CI job keeps crashing (https://github.com/GenericMappingTools/pygmt/actions/runs/10024733960/job/27707337991?pr=3338).

[seisman]

header.z_min and header.z_max are not checked here, because they have invalid values like 1.79769313486231570000e+308 (i.e., DBL_MAX).

$ gmt grdinfo @earth_day_01d_p                               
/Users/seisman/.gmt/server/earth/earth_day/earth_day_01d_p.tif: Title: Grid imported via GDAL
/Users/seisman/.gmt/server/earth/earth_day/earth_day_01d_p.tif: Command:
/Users/seisman/.gmt/server/earth/earth_day/earth_day_01d_p.tif: Remark:
/Users/seisman/.gmt/server/earth/earth_day/earth_day_01d_p.tif: Pixel node registration used [Geographic grid]
/Users/seisman/.gmt/server/earth/earth_day/earth_day_01d_p.tif: Grid file format: gd = Import/export through GDAL
/Users/seisman/.gmt/server/earth/earth_day/earth_day_01d_p.tif: x_min: -180 x_max: 180 x_inc: 1 name: x n_columns: 360
/Users/seisman/.gmt/server/earth/earth_day/earth_day_01d_p.tif: y_min: -90 y_max: 90 y_inc: 1 name: y n_rows: 180
/Users/seisman/.gmt/server/earth/earth_day/earth_day_01d_p.tif: v_min: 1.79769313486e+308 v_max: -1.79769313486e+308 name: z
/Users/seisman/.gmt/server/earth/earth_day/earth_day_01d_p.tif: scale_factor: 1 add_offset: 0
/Users/seisman/.gmt/server/earth/earth_day/earth_day_01d_p.tif: Default CPT:
+proj=longlat +R=6378137 +no_defs

[seisman]

Will merge this PR in 24 hours.

[weiji14]

Ok, we can continue more in #3128