How to see rHEALPix with QGIS?

Question

Using the spatial reference +proj=rhealpix and comparing with "pure WGS84", in separate tables to avoid confusion in the QGIS renderization. But no visualization... is a QGIS bug? For example:

CREATE TABLE test1_b AS    -- Pure WGS84
  SELECT 1 AS gid,
         ST_SetSRID(ST_GeomFromGeoHash('d2g6b'), 4326) AS geom
;
CREATE TABLE test1_c AS   -- rHEALPix to WGS84 (double transform)
  SELECT 1 AS gid, 
      ST_Transform(
        ST_Transform(
           ST_SetSRID(ST_GeomFromGeoHash('d2g6c'),4326),
           '+proj=rhealpix +ellps=WGS84 +south_square=0 +north_square=2'
        ),
        4326
      ) AS geom
; -- BUG!?  ERROR:  ST_Transform: Input geometry has unknown (0) SRID
CREATE TABLE test1_f AS   -- rHEALPix
  SELECT 1 AS gid, 
        ST_Transform(
           ST_SetSRID(ST_GeomFromGeoHash('d2g6f'),4326),
           '+proj=rhealpix +ellps=WGS84 +south_square=0 +north_square=2'
        ) AS geom
;

QGIS v3 only shows test1_b:

Notes

Geohashes d2g6b, d2g6c and d2g6f are small boxes and horizontal neighbors.
In a related question, "How to see ISEA with QGIS?", no clues about visualization, except that rHEALPix is fine because has inverse — as Proj guide say about available forms, "Forward and inverse, spherical and ellipsoidal".
Double transform also not works.
Dump to text seems good: select st_astext(geom) from test1_f results the WKT below

 POLYGON((
  -8238641.204453768 612551.682534819,
  -8238641.204453768 618263.7731066671,
  -8233754.703146263 618263.7731066671,
  -8233754.703146263 612551.682534819,
  -8238641.204453768 612551.682534819
))

Using WGS84 as option after alert,

Is it an "old-PROJ problem"? – Peter Krauss Feb 17 '22 at 11:21 — Peter Krauss, Feb 17 '22 at 11:21

Peter Krauss · Answer 1 · 2022-02-17T22:48:07.393

The PostGIS inverse only works for registered transform

The PostGIS error, "ST_Transform: Input geometry has unknown (0) SRID" is sayng that we must to INSERT the projection in the spatial_ref_sys table. A new and non-standard projection must use the "free range" of the srid column, you can use any arbitrary integer after SELECT MAX(srid) FROM spatial_ref_sys of the PostGIS standard installation, is 104992. So we can use 955001.

INSERT INTO spatial_ref_sys (srid, auth_name, auth_srid, proj4text) VALUES
(
  955001, -- arbitrary number after 105000
  'DGGS:rHEALPix', -- arbitrary label
  9000,   -- arbitrary number
  '+proj=rhealpix +ellps=WGS84 +south_square=0 +north_square=2'
);

Repeat with SRID instead string

And not need one table per projection, QGIS manages the projection diversity. Here adding a name for each.

CREATE TABLE test2 AS
    SELECT 1 AS gid, 'WGS84' as name,
           ST_SetSRID(ST_GeomFromGeoHash('d2g6b'), 4326) AS geom
UNION ALL
  SELECT 2, 'back',
      ST_Transform(    -- PostGIS inversion
        ST_Transform(
           ST_SetSRID(ST_GeomFromGeoHash('d2g6c'),4326),
           955001 -- rhealpix
        ),
        4326      -- back to "pure WGS84" (PostGIS invertion) 
      )
UNION ALL
SELECT 3 AS gid, 'rHEALPix', 
        ST_Transform(
           ST_SetSRID(ST_GeomFromGeoHash('d2g6f'),4326),
           955001 -- rhealpix
        )
;

Result on QGIS

Discussion and conclusion

QGIS say nothing, but the PostGIS error was a good clue. QGIS need to use the inverse transform to show the geometry, and it only works for registered SRIDs.