/******************************************************************************* * Copyright (c) 2015 Voyager Search and MITRE * All rights reserved. This program and the accompanying materials * are made available under the terms of the Apache License, Version 2.0 which * accompanies this distribution and is available at * http://www.apache.org/licenses/LICENSE-2.0.txt ******************************************************************************/ package org.locationtech.spatial4j.shape.impl; import org.locationtech.spatial4j.context.SpatialContext; import org.locationtech.spatial4j.distance.DistanceUtils; import org.locationtech.spatial4j.shape.Point; import org.locationtech.spatial4j.shape.Rectangle; import org.locationtech.spatial4j.shape.SpatialRelation; import java.util.Formatter; import java.util.Locale; /** * A circle as it exists on the surface of a sphere. */ public class GeoCircle extends CircleImpl { private GeoCircle inverseCircle;//when distance reaches > 1/2 way around the world, cache the inverse. private double horizAxisY;//see getYAxis public GeoCircle(Point p, double radiusDEG, SpatialContext ctx) { super(p, radiusDEG, ctx); assert ctx.isGeo(); init(); } @Override public void reset(double x, double y, double radiusDEG) { super.reset(x, y, radiusDEG); init(); } private void init() { if (radiusDEG > 90) { //--spans more than half the globe assert enclosingBox.getWidth() == 360; double backDistDEG = 180 - radiusDEG; if (backDistDEG > 0) { double backRadius = 180 - radiusDEG; double backX = DistanceUtils.normLonDEG(getCenter().getX() + 180); double backY = DistanceUtils.normLatDEG(getCenter().getY() + 180); //Shrink inverseCircle as small as possible to avoid accidental overlap. // Note that this is tricky business to come up with a value small enough // but not too small or else numerical conditioning issues become a problem. backRadius -= Math.max(Math.ulp(Math.abs(backY)+backRadius), Math.ulp(Math.abs(backX)+backRadius)); if (inverseCircle != null) { inverseCircle.reset(backX, backY, backRadius); } else { inverseCircle = new GeoCircle(ctx.makePoint(backX, backY), backRadius, ctx); } } else { inverseCircle = null;//whole globe } horizAxisY = getCenter().getY();//although probably not used } else { inverseCircle = null; double _horizAxisY = ctx.getDistCalc().calcBoxByDistFromPt_yHorizAxisDEG(getCenter(), radiusDEG, ctx); //some rare numeric conditioning cases can cause this to be barely beyond the box if (_horizAxisY > enclosingBox.getMaxY()) { horizAxisY = enclosingBox.getMaxY(); } else if (_horizAxisY < enclosingBox.getMinY()) { horizAxisY = enclosingBox.getMinY(); } else { horizAxisY = _horizAxisY; } //assert enclosingBox.relate_yRange(horizAxis,horizAxis,ctx).intersects(); } } @Override protected double getYAxis() { return horizAxisY; } /** * Called after bounding box is intersected. * @param bboxSect INTERSECTS or CONTAINS from enclosingBox's intersection * @return DISJOINT, CONTAINS, or INTERSECTS (not WITHIN) */ @Override protected SpatialRelation relateRectanglePhase2(Rectangle r, SpatialRelation bboxSect) { if (inverseCircle != null) { return inverseCircle.relate(r).inverse(); } //if a pole is wrapped, we have a separate algorithm if (enclosingBox.getWidth() == 360) { return relateRectangleCircleWrapsPole(r, ctx); } //This is an optimization path for when there are no dateline or pole issues. if (!enclosingBox.getCrossesDateLine() && !r.getCrossesDateLine()) { return super.relateRectanglePhase2(r, bboxSect); } //Rectangle wraps around the world longitudinally creating a solid band; there are no corners to test intersection if (r.getWidth() == 360) { return SpatialRelation.INTERSECTS; } //do quick check to see if all corners are within this circle for CONTAINS int cornersIntersect = numCornersIntersect(r); if (cornersIntersect == 4) { //ensure r's x axis is within c's. If it doesn't, r sneaks around the globe to touch the other side (intersect). SpatialRelation xIntersect = r.relateXRange(enclosingBox.getMinX(), enclosingBox.getMaxX()); if (xIntersect == SpatialRelation.WITHIN) return SpatialRelation.CONTAINS; return SpatialRelation.INTERSECTS; } //INTERSECT or DISJOINT ? if (cornersIntersect > 0) return SpatialRelation.INTERSECTS; //Now we check if one of the axis of the circle intersect with r. If so we have // intersection. /* x axis intersects */ if ( r.relateYRange(getYAxis(), getYAxis()).intersects() // at y vertical && r.relateXRange(enclosingBox.getMinX(), enclosingBox.getMaxX()).intersects() ) return SpatialRelation.INTERSECTS; /* y axis intersects */ if (r.relateXRange(getXAxis(), getXAxis()).intersects()) { // at x horizontal double yTop = getCenter().getY()+ radiusDEG; assert yTop <= 90; double yBot = getCenter().getY()- radiusDEG; assert yBot >= -90; if (r.relateYRange(yBot, yTop).intersects())//back bottom return SpatialRelation.INTERSECTS; } return SpatialRelation.DISJOINT; } private SpatialRelation relateRectangleCircleWrapsPole(Rectangle r, SpatialContext ctx) { //This method handles the case where the circle wraps ONE pole, but not both. For both, // there is the inverseCircle case handled before now. The only exception is for the case where // the circle covers the entire globe, and we'll check that first. if (radiusDEG == 180)//whole globe return SpatialRelation.CONTAINS; //Check if r is within the pole wrap region: double yTop = getCenter().getY() + radiusDEG; if (yTop > 90) { double yTopOverlap = yTop - 90; assert yTopOverlap <= 90; if (r.getMinY() >= 90 - yTopOverlap) return SpatialRelation.CONTAINS; } else { double yBot = point.getY() - radiusDEG; if (yBot < -90) { double yBotOverlap = -90 - yBot; assert yBotOverlap <= 90; if (r.getMaxY() <= -90 + yBotOverlap) return SpatialRelation.CONTAINS; } else { //This point is probably not reachable ?? assert yTop == 90 || yBot == -90;//we simply touch a pole //continue } } //If there are no corners to check intersection because r wraps completely... if (r.getWidth() == 360) return SpatialRelation.INTERSECTS; //Check corners: int cornersIntersect = numCornersIntersect(r); // (It might be possible to reduce contains() calls within nCI() to exactly two, but this intersection // code is complicated enough as it is.) double frontX = getCenter().getX(); if (cornersIntersect == 4) {//all double backX = frontX <= 0 ? frontX + 180 : frontX - 180; if (r.relateXRange(backX, backX).intersects()) return SpatialRelation.INTERSECTS; else return SpatialRelation.CONTAINS; } else if (cornersIntersect == 0) {//none if (r.relateXRange(frontX, frontX).intersects()) return SpatialRelation.INTERSECTS; else return SpatialRelation.DISJOINT; } else//partial return SpatialRelation.INTERSECTS; } /** Returns either 0 for none, 1 for some, or 4 for all. */ private int numCornersIntersect(Rectangle r) { //We play some logic games to avoid calling contains() which can be expensive. boolean bool;//if true then all corners intersect, if false then no corners intersect // for partial, we exit early with 1 and ignore bool. bool = (contains(r.getMinX(),r.getMinY())); if (contains(r.getMinX(),r.getMaxY())) { if (!bool) return 1;//partial } else { if (bool) return 1;//partial } if (contains(r.getMaxX(),r.getMinY())) { if (!bool) return 1;//partial } else { if (bool) return 1;//partial } if (contains(r.getMaxX(),r.getMaxY())) { if (!bool) return 1;//partial } else { if (bool) return 1;//partial } return bool?4:0; } @Override public String toString() { //Add distance in km, which may be easier to recognize. double distKm = DistanceUtils.degrees2Dist(radiusDEG, DistanceUtils.EARTH_MEAN_RADIUS_KM); //instead of String.format() so that we get consistent output no matter the locale String dStr = new Formatter(Locale.ROOT).format("%.1f\u00B0 %.2fkm", radiusDEG, distKm).toString(); return "Circle(" + point + ", d=" + dStr + ')'; } }