chore(sf|h3): reimplement polyfill h3 functions

README.md

@@ -62,6 +62,43 @@ Right now the only way to get access the Analytics toolbox is by installing it d
 | Postgres | [README.md](./clouds/postgres/README.md) |
 | Databricks | [README.md](./clouds/databricks/README.md) |
 
+### Useful make commands
+
+To run tests, switch to a specific cloud directory. For example, Showflake: `cd clouds/snowflake`.  
+
+```
+# All tests
+make test
+
+# Specific module(s)
+make test modules=h3
+make test modules=h3,transformations
+
+# Specific function(s)
+make test functions=H3_POLYFILL
+make test functions=H3_POLYFILL,ST_BUFFER
+```
+
+### Rebuild h3-js 3.7.2 dependency
+First, ensure you have yarn and docker installed.  
+
+```
+wget https://github.com/uber/h3-js/releases/tag/v3.7.2
+unzip h3-js-3.7.2.zip
+cd h3-js-3.7.2
+yarn docker-boot && yarn build-emscripten
+```
+
+Remove all the unneeded bindings from the `lib/bindings.js`  
+
+Then run:  
+```
+yarn docker-emscripten-run
+```
+
+Your new library file is available at `out/a.out.js`. Copy it to the correct location with the new filename. For example: `cp out/a.out.js  ~/development/analytics-toolbox/core/clouds/snowflake/libraries/javascript/src/h3/h3_polyfill/libh3_custom.js`. Ensure it is named `libh3_custom.js`.
+
+
 ## Contribute
 
 This project is public. We are more than happy of receiving feedback and contributions. Feel free to open a ticket with a bug, a doubt or a discussion, or open a pull request with a fix or a new feature.

clouds/snowflake/libraries/javascript/libs/h3_kring_distances.js

@@ -2,4 +2,4 @@ import { h3Distance } from '../src/h3/h3_kring_distances/h3core_custom';
 
 export default {
     h3Distance
-};
+};

clouds/snowflake/libraries/javascript/libs/h3_polyfill.js

@@ -1,7 +1,11 @@
-import { polyfill } from '../src/h3/h3_polyfill/h3core_custom';
-import { bboxClip } from '@turf/turf';
+import { booleanContains, booleanIntersects, intersect, polygon, multiPolygon } from '@turf/turf';
+import { h3ToGeoBoundary } from '../src/h3/h3_polyfill/h3core_custom';
 
 export default {
-    bboxClip,
-    polyfill
+    booleanContains,
+    booleanIntersects,
+    intersect,
+    polygon,
+    multiPolygon,
+    h3ToGeoBoundary
 };

clouds/snowflake/libraries/javascript/test/h3.test.js

@@ -13,7 +13,7 @@ test('h3 library defined', () => {
     expect(hexRingLib.hexRing).toBeDefined();
     expect(hexRingLib.h3IsValid).toBeDefined();
     expect(isPentagonLib.h3IsPentagon).toBeDefined();
-    expect(h3PolyfillLib.polyfill).toBeDefined();
+    expect(h3PolyfillLib.h3ToGeoBoundary).toBeDefined();
     expect(boundaryLib.h3ToGeoBoundary).toBeDefined();
     expect(boundaryLib.h3IsValid).toBeDefined();
     expect(kringDistancesLib.h3Distance).toBeDefined();

clouds/snowflake/modules/doc/h3/H3_POLYFILL.md

@@ -1,15 +1,21 @@
 ## H3_POLYFILL
 
 ```sql:signature
-H3_POLYFILL(geography, resolution)
+H3_POLYFILL(geography, resolution)  
+H3_POLYFILL(geography, resolution, mode)
 ```
 
 **Description**
 
-Returns an array with all the H3 cell indexes **with centers** contained in a given polygon. It will return `null` on error (invalid geography type or resolution out of bounds).
+Returns an array with all H3 cell indicies contained in the given polygon. There are three modes which decide if a H3 cell is contained in the polygon:  
 
-* `geography`: `GEOGRAPHY` **polygon** or **multipolygon** representing the shape to cover.
-* `resolution`: `INT` number between 0 and 15 with the [H3 resolution](https://h3geo.org/docs/core-library/restable).
+- **center** (Default) - The center point of the H3 cell must be within the polygon
+- **contains** - the H3 cell must be fully contained within the polygon (least inclusive)
+- **intersects** - The H3 cell intersects in any way with the polygon (most inclusive)
+
+- `geography`: `GEOGRAPHY` **polygon** or **multipolygon** representing the shape to cover. **GeometryCollections** are also allowed but they should contain **polygon** or **multipolygon** geographies. Non-Polygon types will not raise an error but will be ignored instead.
+- `resolution`: `INT` number between 0 and 15 with the [H3 resolution](https://h3geo.org/docs/core-library/restable).
+- `mode`: `STRING` `<center|contains|intersects>`. Optional. Defaults to 'center' mode.
 
 **Return type**
 
@@ -19,7 +25,7 @@ Returns an array with all the H3 cell indexes **with centers** contained in a gi
 
 ```sql
 SELECT carto.H3_POLYFILL(
-    TO_GEOGRAPHY('POLYGON ((30 10, 40 40, 20 40, 10 20, 30 10))'), 4);
+    TO_GEOGRAPHY('POLYGON ((30 10, 40 40, 20 40, 10 20, 30 10))'), 4, 'center');
 -- 842da29ffffffff
 -- 843f725ffffffff
 -- 843eac1ffffffff

clouds/snowflake/modules/sql/h3/H3_POLYFILL.sql

@@ -1,83 +1,177 @@
-----------------------------
--- Copyright (C) 2021 CARTO
-----------------------------
+--------------------------------
+-- Copyright (C) 2021-2024 CARTO
+--------------------------------
 
-CREATE OR REPLACE FUNCTION @@SF_SCHEMA@@._H3_POLYFILL
-(geojson STRING, input_resolution DOUBLE)
-RETURNS ARRAY
+CREATE OR REPLACE FUNCTION @@SF_SCHEMA@@._HAS_POLYGON_JS
+(geojson STRING)
+RETURNS BOOLEAN
 LANGUAGE JAVASCRIPT
 IMMUTABLE
 AS $$
-    if (!GEOJSON || INPUT_RESOLUTION == null) {
-        return [];
+    let inputGeoJSON = JSON.parse(GEOJSON);
+    let geometries = inputGeoJSON.geometries ? inputGeoJSON.geometries : [inputGeoJSON] // geometrycollection or regular feature geometry
+    for (let g of geometries) {
+        if (g.type === 'Polygon' || g.type === 'MultiPolygon') {
+	    return true
+        }	
     }
+    return false
+$$;
+
+CREATE OR REPLACE FUNCTION @@SF_SCHEMA@@._HAS_POLYGON
+(geog GEOGRAPHY)
+RETURNS BOOLEAN
+LANGUAGE SQL
+IMMUTABLE
+AS $$
+	@@SF_SCHEMA@@._HAS_POLYGON_JS(CAST(ST_ASGEOJSON(GEOG) AS STRING))
+$$;
+
+CREATE OR REPLACE FUNCTION @@SF_SCHEMA@@._FILTER_GEOG_JS
+(geojson STRING)
+RETURNS STRING
+LANGUAGE JAVASCRIPT
+IMMUTABLE
+AS $$
+    // remove non-polygons and split polygons >= 180 degrees
+    // output is a always MULTIPOLYGON
+
+    let inputGeoJSON = JSON.parse(GEOJSON);
 
     @@SF_LIBRARY_H3_POLYFILL@@
 
-    const resolution = Number(INPUT_RESOLUTION);
-    if (resolution < 0 || resolution > 15) {
-        return [];
-    }
+    const westernHemisphere = h3PolyfillLib.polygon([[ [-180, 90], [0, 90], [0, -90], [-180, -90], [-180, 90]]]);
+    const easternHemisphere = h3PolyfillLib.polygon([[ [0, 90], [180, 90], [180, -90], [0, -90], [0, 90] ]]);
 
-    const bboxA = [-180, -90, 0, 90]
-    const bboxB = [0, -90, 180, 90]
-    const featureGeometry = JSON.parse(GEOJSON)
-    let polygonCoordinatesA = [];
-    let polygonCoordinatesB = [];
-    switch(featureGeometry.type) {
-        case 'GeometryCollection':
-            featureGeometry.geometries.forEach(function (geom) {
-                if (geom.type === 'MultiPolygon') {
-                    var clippedGeometryA = h3PolyfillLib.bboxClip(geom, bboxA).geometry;
-                    polygonCoordinatesA = polygonCoordinatesA.concat(clippedGeometryA.coordinates);
-                    var clippedGeometryB = h3PolyfillLib.bboxClip(geom, bboxB).geometry;
-                    polygonCoordinatesB = polygonCoordinatesB.concat(clippedGeometryB.coordinates);
-                } else if (geom.type === 'Polygon') {
-                    var clippedGeometryA = h3PolyfillLib.bboxClip(geom, bboxA).geometry;
-                    polygonCoordinatesA = polygonCoordinatesA.concat([clippedGeometryA.coordinates]);
-                    var clippedGeometryB = h3PolyfillLib.bboxClip(geom, bboxB).geometry;
-                    polygonCoordinatesB = polygonCoordinatesB.concat([clippedGeometryB.coordinates]);
-                }
-            });
-        break;
-        case 'MultiPolygon':
-            var clippedGeometryA = h3PolyfillLib.bboxClip(featureGeometry, bboxA).geometry;
-            polygonCoordinatesA = clippedGeometryA.coordinates;
-            var clippedGeometryB = h3PolyfillLib.bboxClip(featureGeometry, bboxB).geometry;
-            polygonCoordinatesB = clippedGeometryB.coordinates;
-        break;
-        case 'Polygon':
-            var clippedGeometryA = h3PolyfillLib.bboxClip(featureGeometry, bboxA).geometry;
-            polygonCoordinatesA = [clippedGeometryA.coordinates];
-            var clippedGeometryB = h3PolyfillLib.bboxClip(featureGeometry, bboxB).geometry;
-            polygonCoordinatesB = [clippedGeometryB.coordinates];
-        break;
-        default:
-            return [];
-    }
+    let polygons = [];
+    let geometries = inputGeoJSON.geometries ? inputGeoJSON.geometries : [inputGeoJSON]
 
-    if (polygonCoordinatesA.length + polygonCoordinatesB.length === 0) {
-        return [];
-    }
+    geometries.forEach(g => {
+        if (g.type === 'Polygon') {
+	    polygons.push({type: 'Feature', geometry: g})	
+        }	
+        else if (g.type === 'MultiPolygon') {
+	    g.coordinates.forEach(ring => polygons.push({type: 'Feature', geometry: {type: 'Polygon', coordinates: ring}}))
+        }
+    });
+
+
+    let intersections = [];
+
+    let intersectAndPush = (hemisphere, poly) => {
+        const intersection = h3PolyfillLib.intersect(poly, hemisphere);
+	if (intersection) {
+            if (intersection.geometry.type === 'Polygon') {
+	        intersections.push(intersection);
+            }	
+            else if (intersection.geometry.type === 'MultiPolygon') {
+	        intersection.geometry.coordinates.forEach(ring => intersections.push({type: 'Feature', geometry: {type: 'Polygon', coordinates: ring}}))
+            }
+	}
+    };
+
+    polygons.forEach(p => {
+        intersectAndPush(westernHemisphere, p);
+        intersectAndPush(easternHemisphere, p);
+    })
+
+    return JSON.stringify(h3PolyfillLib.multiPolygon(intersections.map(i => i.geometry.coordinates)).geometry)
+$$;
+
+CREATE OR REPLACE FUNCTION @@SF_SCHEMA@@._FILTER_GEOG
+(geog GEOGRAPHY)
+RETURNS GEOGRAPHY
+LANGUAGE SQL
+IMMUTABLE
+AS $$
+	TO_GEOGRAPHY(@@SF_SCHEMA@@._FILTER_GEOG_JS(CAST(ST_ASGEOJSON(GEOG) AS STRING)))
+$$;
+
+CREATE OR REPLACE FUNCTION @@SF_SCHEMA@@._H3_POLYFILL_CONTAINS
+(geojson STRING, indicies ARRAY)
+RETURNS ARRAY
+LANGUAGE JAVASCRIPT
+IMMUTABLE
+AS $$
+    let results = []
+    let inputGeoJSON = JSON.parse(GEOJSON);
+
+    @@SF_LIBRARY_H3_POLYFILL@@
+
+    // @@SF_SCHEMA@@.ST_BUFFER demotes MULTIPOLYGONs to POLYGON if it only has one ring. So we check again here.
+    let polygons = inputGeoJSON.type === 'MultiPolygon' ? inputGeoJSON.coordinates.map(ring => h3PolyfillLib.polygon(ring)) : [h3PolyfillLib.polygon(inputGeoJSON.coordinates)]
+
+
+    INDICIES.forEach(h3Index => {
+	polygons.some(p => {
+	    if (h3PolyfillLib.booleanContains(p, h3PolyfillLib.polygon([h3PolyfillLib.h3ToGeoBoundary(h3Index, true)]))) {
+	        results.push(h3Index)	
+	    }
+	})	
+    })
+    return results
+$$;
+
+CREATE OR REPLACE SECURE FUNCTION @@SF_SCHEMA@@._H3_POLYFILL_INTERSECTS_FILTER
+(h3Indicies ARRAY, geojson STRING)
+RETURNS ARRAY
+LANGUAGE JAVASCRIPT
+IMMUTABLE
+AS $$
+
+    let results = []
+    let inputGeoJSON = JSON.parse(GEOJSON);
+
+    @@SF_LIBRARY_H3_POLYFILL@@
+
+    // @@SF_SCHEMA@@.ST_BUFFER demotes MULTIPOLYGONs to POLYGON if it only has one ring. So we check again here.
+    let polygons = inputGeoJSON.type === 'MultiPolygon' ? inputGeoJSON.coordinates.map(ring => h3PolyfillLib.polygon(ring)) : [h3PolyfillLib.polygon(inputGeoJSON.coordinates)]
+
+    H3INDICIES.forEach(h3Index => {
+	if (polygons.some(p => h3PolyfillLib.booleanIntersects(p, h3PolyfillLib.polygon([h3PolyfillLib.h3ToGeoBoundary(h3Index, true)])))) {
+	    results.push(h3Index)
+	}
+    })
+    return [...new Set(results)]
+$$;
 
-    let hexesA = polygonCoordinatesA.reduce(
-        (acc, coordinates) => acc.concat(h3PolyfillLib.polyfill(coordinates, resolution, true)),
-        []
-    ).filter(h => h != null);
-    let hexesB = polygonCoordinatesB.reduce(
-        (acc, coordinates) => acc.concat(h3PolyfillLib.polyfill(coordinates, resolution, true)),
-        []
-    ).filter(h => h != null);
-    hexes = [...hexesA, ...hexesB];
-    hexes = [...new Set(hexes)];
-
-    return hexes;
+CREATE OR REPLACE FUNCTION @@SF_SCHEMA@@._CHECK_TOO_WIDE(geo GEOGRAPHY)
+RETURNS BOOLEAN
+AS
+$$
+    CASE
+        WHEN ST_XMax(geo) < ST_XMin(geo) THEN
+            -- Adjusts for crossing the antimeridian
+            360 + ST_XMax(geo) - ST_XMin(geo) >= 180
+        ELSE
+            ST_XMax(geo) - ST_XMin(geo) >= 180
+    END
 $$;
 
 CREATE OR REPLACE SECURE FUNCTION @@SF_SCHEMA@@.H3_POLYFILL
 (geog GEOGRAPHY, resolution INT)
 RETURNS ARRAY
 IMMUTABLE
 AS $$
-    @@SF_SCHEMA@@._H3_POLYFILL(CAST(ST_ASGEOJSON(GEOG) AS STRING), CAST(RESOLUTION AS DOUBLE))
+    IFF(
+        GEOG IS NOT NULL AND RESOLUTION >= 0 AND RESOLUTION <= 15 AND @@SF_SCHEMA@@._HAS_POLYGON(GEOG),
+	COALESCE(H3_POLYGON_TO_CELLS_STRINGS(@@SF_SCHEMA@@._FILTER_GEOG(GEOG), RESOLUTION), []),
+	[]
+    ) 
+$$;
+
+CREATE OR REPLACE SECURE FUNCTION @@SF_SCHEMA@@.H3_POLYFILL
+(geog GEOGRAPHY, resolution INT, mode STRING)
+RETURNS ARRAY
+IMMUTABLE
+AS $$
+    CASE WHEN GEOG IS NULL OR RESOLUTION < 0 OR RESOLUTION > 15 OR NOT @@SF_SCHEMA@@._HAS_POLYGON(GEOG) THEN []
+	WHEN MODE = 'center' THEN @@SF_SCHEMA@@.H3_POLYFILL(GEOG, RESOLUTION)
+	WHEN MODE = 'intersects' THEN
+	    CASE WHEN @@SF_SCHEMA@@._CHECK_TOO_WIDE(GEOG) THEN @@SF_SCHEMA@@._H3_POLYFILL_INTERSECTS_FILTER(H3_COVERAGE_STRINGS(@@SF_SCHEMA@@._FILTER_GEOG(GEOG), RESOLUTION), CAST(ST_ASGEOJSON(GEOG) AS STRING))
+	        ELSE H3_COVERAGE_STRINGS(@@SF_SCHEMA@@.ST_BUFFER(@@SF_SCHEMA@@._FILTER_GEOG(GEOG), CAST(0.00000001 AS DOUBLE)), RESOLUTION)
+	    END
+	WHEN MODE = 'contains' THEN @@SF_SCHEMA@@._H3_POLYFILL_CONTAINS(CAST(ST_ASGEOJSON(@@SF_SCHEMA@@.ST_BUFFER(@@SF_SCHEMA@@._FILTER_GEOG(GEOG), CAST(0.00000001 AS DOUBLE)) ) AS STRING), @@SF_SCHEMA@@.H3_POLYFILL(GEOG, RESOLUTION))
+	ELSE []
+    END
 $$;