H3 cells in Cesium view

h3_cells_visible.js

import * as Cesium from "cesium";
import * as h3 from "https://cdnjs.cloudflare.com/ajax/libs/h3-js/4.3.0/h3-js.es.min.js";

/**************************************************************************************************
 * "Patch" for Cesium camera computeViewRectangle.
 * 
 * This path provides a modified implementation of computeViewRectangle that
 * returns the correct view rectangle when the camera is oriented towards the south
 * and the horizon is visible. The default implementation truncates the view rectangle
 * at a lower corner of the view which can significantly reduce the reported view
 * rectangle. 
 * 
 * Load this in sandcastle.
 */

const VIEWER = new Cesium.Viewer(
    'cesiumContainer', {
        timeline: false,
        animation: false,
        baseLayerPicker: false,
        terrain: Cesium.Terrain.fromWorldTerrain(),
    });

const TOOLBAR = document.getElementById('toolbar');

const someColors = [
    Cesium.Color.GOLD,
    Cesium.Color.AQUA,
    Cesium.Color.BLUEVIOLET,
    Cesium.Color.CHARTREUSE,
    Cesium.Color.CORAL,
    Cesium.Color.DEEPPINK
];

// The global bounding box latitude and longitude values.
const GLOBAL_RECT = new Cesium.Rectangle(-Cesium.Math.PI, -Cesium.Math.PI_OVER_TWO, Cesium.Math.PI, Cesium.Math.PI_OVER_TWO);

let bbEntity = null;


/**
 * Scans vertically at x pixels min_y to max_y to find a pixel that intersect 
 * the ellipsoid. 
 *
 * In screen space canvas, 0,0 is the top left corner.
 * Returns the cartographic coordinates of the picked point or null if 
 * there is no intersection.
 */
function computeHorizonPointY(camera, ellips, x, min_y, max_y) {
    for (let j=min_y; j < max_y; j += 5) {
        const cp = camera.pickEllipsoid(new Cesium.Cartesian2(x, j), ellips);
        if (cp) {
            return ellips.cartesianToCartographic(cp);
        }
    }
    return null;
}


/**
 * Compute if point is visible in current view
 * p: Cartographic3
 */
function isPositionVisible(camera, ellips, cartesian) {
    const frustum = camera.frustum;
    const cullingVolume = frustum.computeCullingVolume(
        camera.position,
        camera.direction,
        camera.up
    );
    const intersection = cullingVolume.computeVisibility(new Cesium.BoundingSphere(cartesian, 0.0));
    if (intersection === Cesium.Intersect.INSIDE) {
        const globeBoundingSphere = new Cesium.BoundingSphere(
            Cesium.Cartesian3.ZERO,
            ellips.minimumRadius
        );
        const occluder = new Cesium.Occluder(
            globeBoundingSphere,
            camera.position
        )
        return occluder.isPointVisible(cartesian);
    }
    return false;
}


/**
 * Return 0 if neither visible, -1 for south, 1 for north
 */
 function isPoleVisible(camera, ellips) {
  const NORTH_POLE = Cesium.Cartographic.toCartesian(
      Cesium.Cartographic.fromDegrees(0.0, 90.0, 1.0, new Cesium.Cartographic()), 
      ellips,
      new Cesium.Cartesian3()
  );
  if (isPositionVisible(camera, ellips, NORTH_POLE)) {
      return 1;
  }
  const SOUTH_POLE = Cesium.Cartographic.toCartesian(
    Cesium.Cartographic.fromDegrees(0.0, -90.0, 1.0, new Cesium.Cartographic()), 
    ellips,
    new Cesium.Cartesian3()
);
if (isPositionVisible(camera, ellips, SOUTH_POLE)) {
      return -1;
  }
  return 0;
}


/**
 * Computes five cartographic points corresponding with the canvas top left,
 * top middle, top right, lower right, and lower left.
 * 
 * The top three points are computed at the intersection of that column of pixels
 * with the ellpsoid.
 * 
 * null is returned for any point that does not intersect with the ellpsoid.
 */
function computeViewHorizonPoints(canvas, camera, ellips) {
    const xs = [0, Math.floor(canvas.width/2), canvas.width];
    const points = [null, null, null, null, null];
    for (let i=0; i<3; i++) {
        points[i] = computeHorizonPointY(camera, ellips, xs[i], 0, canvas.height);        
    }
    points[3] = computeHorizonPointY(camera, ellips, xs[0], canvas.height-1, canvas.height);
    points[4] = computeHorizonPointY(camera, ellips, xs[2], canvas.height-1, canvas.height);
    return points;
}


/**
 * Given a list of cartographic points, compute
 * the bounding rectangle for the points.
 */
function pointsToBoundingRectangle(camera, ellips, points, result) {
    for (let i=0; i < points.length; i++) {
        if (!points[i]) {
            return GLOBAL_RECT;
        }
    }    
    result = Cesium.Rectangle.fromCartographicArray(points, result);
    const pvisible = isPoleVisible(camera, ellips);
    if (pvisible === 1) {
        result.west = -Cesium.Math.PI;
        result.east = Cesium.Math.PI;
        result.north = Cesium.Math.PI_OVER_TWO;
    } else if (pvisible === -1) {
        result.west = -Cesium.Math.PI;
        result.east = Cesium.Math.PI;
        result.south = -Cesium.Math.PI_OVER_TWO;
    }
    return result;
}


/**
 * Compute the (approximate) bounding rectangle of the camera view.
 * 
 * returns Cesium.Rectangle
 */
Cesium.Camera.prototype.computeViewRectangle2 = function(ellips, result) {
  const points = computeViewHorizonPoints(this._scene.canvas, this, ellips);
  return pointsToBoundingRectangle(this, ellips, points, result);
}

function getH3IndexesInBondingBox(minLat, minLng, maxLat, maxLng, resolution) {
  // Define the bounding box as a GeoJSON-like polygon:	
  // The structure is [outer ring, ...holes], and each ring is an array of [lat, lng]	
  const polygon = [	
    [	
      [maxLat, maxLng], // Top Right	
      [minLat, maxLng], // Bottom Right	
      [minLat, minLng], // Bottom Left	
      [maxLat, minLng], // Top Left	
      [maxLat, maxLng]  // Closing the loop (same as Top Right)	
    ]	
  ];	

  // The h3.polygonToCells function (or h3.polyfill in older v3 library)	
  // returns an array of H3 indexes whose centroids fall within the polygon.	
  const h3Indexes = h3.polygonToCells(polygon, resolution);	
  return h3Indexes;	  
}

function distanceToLevel(d) {
  // really rough estimate of distance to center of view vs zoom level.
  const ranges = [
    5655046,
    2534653,
    1123998,
    326074,
    78393,
    30721,
    18239,
    4639,
    1102,
    500,
    250
  ];
  if (d < 0) {
    return 0;
  }
  for (let i=0; i < ranges.length; i++) {
    if (d > ranges[i]) {
      return i+1;
    }
  }
  return 12;
}

function createH3GroundPrimitives(h3Indexes, alpha=0.3) {
  return h3Indexes.map(index => {
    const color = someColors[Math.floor(Math.random()*someColors.length)].withAlpha(alpha);
    const boundary = h3.cellToBoundary(index); 
    const positions = boundary.map(coords => {
      return Cesium.Cartesian3.fromDegrees(coords[1], coords[0]);
    });
    const instance = new Cesium.GeometryInstance({
      geometry: new Cesium.PolygonGeometry({
        polygonHierarchy: new Cesium.PolygonHierarchy(positions),
        height: 0 // Keep at ground level
      }),
      attributes: {
        color: Cesium.ColorGeometryInstanceAttribute.fromColor(color)
      }
    });
    return new Cesium.GroundPrimitive({
      geometryInstances: instance,
      appearance: new Cesium.PerInstanceColorAppearance({
        flat: true // Flat shading often looks better for data overlays
      })
    });
  });
}


const handler = new Cesium.ScreenSpaceEventHandler(VIEWER.canvas);
handler.setInputAction((click) => {
    const bb = VIEWER.camera.computeViewRectangle2(VIEWER.scene.globe.ellipsoid, new Cesium.Rectangle());

    const minLng = Cesium.Math.toDegrees(bb.west);
    const minLat = Cesium.Math.toDegrees(bb.south);
    const maxLat = Cesium.Math.toDegrees(bb.north)
    const maxLng = Cesium.Math.toDegrees(bb.east);

    let h3cells = [];
    let distanceToTerrain = -1;
    if (bb === GLOBAL_RECT) {
      h3cells = h3.getRes0Cells();
    } else {
      const ray = VIEWER.camera.getPickRay(
        new Cesium.Cartesian2(
          VIEWER.canvas.width / 2, 
          VIEWER.canvas.height / 2
        )
      );
      const intersection = VIEWER.scene.globe.pick(ray, VIEWER.scene);
      if (intersection) {
          distanceToTerrain = Cesium.Cartesian3.distance(VIEWER.camera.position, intersection);
      }
      h3cells = getH3IndexesInBondingBox(minLat, minLng, maxLat, maxLng, distanceToLevel(distanceToTerrain));
    }

    TOOLBAR.innerHTML = `<pre>
        west: ${minLng.toFixed(2)}
        south:${minLat.toFixed(2)}
        east:${maxLng.toFixed(2)}
        north:${maxLat.toFixed(2)}
        distance:${distanceToTerrain.toFixed(0)}
        ncells:${h3cells.length}</pre>`;
    VIEWER.scene.primitives.removeAll();
    if (h3cells.length > 0) {
      const primitives = createH3GroundPrimitives(h3cells.slice(0,2000));
      primitives.forEach(p => VIEWER.scene.primitives.add(p));
    }
}, Cesium.ScreenSpaceEventType.LEFT_UP);

Example in sandcastle.