Improve the logic for determine grid/image gtype based on upstream GM…

…T codes

pygmt/datatypes/header.py

@@ -250,6 +250,16 @@ def gtype(self) -> int:
         "lon"/"lat" or have units "degrees_east"/"degrees_north", then the grid is
         assumed to be geographic.
         """
+        gtype = 0  # Cartesian by default
+
         dims = self.dims
-        gtype = 1 if dims[0] == "lat" and dims[1] == "lon" else 0
+        if dims[0] == "lat" and dims[1] == "lon":
+            # Check dimensions for grids that following CF-conventions
+            gtype = 1
+        elif self.ProjRefPROJ4 is not None:
+            # Check ProjRefPROJ4 for images imported via GDAL.
+            # The logic comes from GMT's `gmtlib_read_image_info` function.
+            projref = self.ProjRefPROJ4.decode()
+            if "longlat" in projref or "latlong" in projref:
+                gtype = 1
         return gtype

pygmt/datatypes/image.py

@@ -166,7 +166,7 @@ def to_dataarray(self) -> xr.DataArray:
             axis:          Y
             actual_range:  [-90.  90.]
         >>> da.gmt.registration, da.gmt.gtype
-        (1, 0)
+        (1, 1)
         """
         # The image header
         header = self.header.contents