0

I wanted to build a simple grid over the earth's surface. I need this functionality for a computer program where I quickly need to check wheter a given coordinate in (latitude, longitude) is inside a cluster.
These cluster are therefore represented as index-pairs (latitude-index, longitude-index). Which I was computing quite naively, independent in each dimension. Here is some Java code on how I compute these indices:

public GeoIndex getIndex(double latitude, double longitude, double cellSize) {
    if (!isPointInBoundaries(latitude, longitude)) {
        return null;
    }
    long latitudeIndex = (long) Math.floor(haversine(boundaries.latBottomLeft, boundaries.longBottomLeft,
                                                    latitude, boundaries.longBottomLeft) / cellSize);
    long longitudeIndex = (long) Math.floor(haversine(boundaries.latBottomLeft, boundaries.longBottomLeft,
                                                    boundaries.latBottomLeft, longitude) / cellSize);

return new GeoIndex(latitudeIndex, longitudeIndex);
}

Here the cellSize is supposed to be the grid length for each grid-cell. And the boundaries object gives a point that is the 0-point in my coordinate system, so to say. The lowest, left-most point from where the coordinate system will be originated. I want to cover regions like, Europe, North America, China...

First question, is this a solid approach and will points always map to the same grid cells? Or is this too simplistic.

My real issue started when I was trying to compute the boundary points of each grid-cell. I though I could just recompute the boundaries. But it seems to be an issue that the earth radius shrinks the farther north you go. Is there a way to recompute the boundaries of each cell? I was trying the following approach:

static toCoordinates(geoIndex) {
    const coordinates = [];
    const cells = geoCluster.getAiActivatedMapCells();
    const distanceInKm = cellSize / 1000.0;
    const latitudePoint = GeoUtils.calculateEndPointWithDistance(
        boundaries.latBottomLeft,
        boundaries.longBottomLeft,
        distanceInKm * geoIndex.latitude,
        0.0
      );
      const longitudePoint = GeoUtils.calculateEndPointWithDistance(
        boundaries.latBottomLeft,
        boundaries.longBottomLeft,
        distanceInKm * geoIndex.longitude,
        90.0
      );
      const bottomLeftPoint = {
        lat: latitudePoint.lat,
        lng: longitudePoint.lng,
      };

      const topLeftPoint = GeoUtils.calculateEndPointWithDistance(
        latitudePoint.lat,
        longitudePoint.lng,
        distanceInKm,
        0.0
      );
      const bottomRightPoint = GeoUtils.calculateEndPointWithDistance(
        latitudePoint.lat,
        longitudePoint.lng,
        distanceInKm,
        90.0
      );
      const topRightPoint = {
        lat: topLeftPoint.lat,
        lng: bottomRightPoint.lng,
      };

      return {
        topLeft: topLeftPoint,
        topRight: topRightPoint,
        bottomRight: bottomRightPoint,
        bottomLeft: bottomLeftPoint,
      };
}

static calculateEndPointWithDistance(
    startLatitude,
    startLongitude,
    distanceInKm,
    bearing
  ) {
    const distanceRatio = distanceInKm / GeoUtils.EARTH_RADIUS;
    const startLat = GeoUtils.deg2rad(startLatitude);
    const startLon = GeoUtils.deg2rad(startLongitude);

    const endLat = Math.asin(
      Math.sin(startLat) * Math.cos(distanceRatio) +
        Math.cos(startLat) *
          Math.sin(distanceRatio) *
          Math.cos(GeoUtils.deg2rad(bearing))
    );
    const endLon =
      startLon +
      Math.atan2(
        Math.sin(GeoUtils.deg2rad(bearing)) *
          Math.sin(distanceRatio) *
          Math.cos(startLat),
        Math.cos(distanceRatio) - Math.sin(startLat) * Math.sin(endLat)
      );

    return { lat: GeoUtils.rad2deg(endLat), lng: GeoUtils.rad2deg(endLon) };
}

Sorry that this is JavaScript. But basically I try to recompute the latitude and longitude of the lower left point of the grid cell independently again. With the function calculateEndPointWithDistance which is an implementation of the the function I found here: https://www.movable-type.co.uk/scripts/latlong.html (Destination point given distance and bearing from start point). And the compute the grid boarders.

However, the longitude is always way off. I don't really understand way. I understand why the grid cells are to big and overlap, but I don't understand why they are way off to west at all. Can you explain me why?

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.