How to use GeometricShapeFactory in geoTools to create a Circle on map

Question

I am currently using the below code to create a GeoJson Polygon. this gives me a bad circle which is not valid...

in this case RADIUS = 1609.34 , which is 1 mile in meters.

        public  GeoJsonPolygon createRadiusPolygon(  Point point,double RADIUS) {       

              GeometricShapeFactory shapeFactory = new GeometricShapeFactory();
              shapeFactory.setNumPoints(32);
              shapeFactory.setCentre(new com.vividsolutions.jts.geom.Coordinate(point.getX(), point.getY()));
              shapeFactory.setSize(RADIUS * 2);
              com.vividsolutions.jts.geom.Geometry circle = shapeFactory.createCircle();
              List<Point> points = new ArrayList<Point>();
                for (com.vividsolutions.jts.geom.Coordinate coordinate : circle.getCoordinates()) {
                    Point lngLatAtl = new Point(coordinate.x, coordinate.y);
                    points.add(lngLatAtl);
                }
                Collections.reverse(points);
                return new GeoJsonPolygon(points);
            }

referenced: http://docs.geotools.org/stable/userguide/library/jts/geometry.html

currently if i use Point(-73.87,40.84) RADIUS = 1609.34, i get the below link. https://gist.githubusercontent.com/VanitySoft/56c4ce0f5c1c7e7fe0461ed46fd5ed11/raw/94544750a140d81780ebe9206395a21ab88bb1f7/circle

===SOLVED== from @Ian answer: Using method in his answer. RADIUS is in miles, to get the Circle used to create the GeoJson.

...
  com.vividsolutions.jts.geom.Point jtsPoint =  new GeometryFactory().createPoint(new com.vividsolutions.jts.geom.Coordinate(point.getY(), point.getX()));
                  javax.measure.Measure measure = Measure.valueOf(RADIUS, NonSI.MILE);
                  com.vividsolutions.jts.geom.Geometry circle = createCircleRadis(measure,CRS.decode("epsg:4326"),jtsPoint );
...

...

Answer 1

If you just want a pretty circle on a flat map, this worked for me:

import java.io.ByteArrayOutputStream;
import java.io.IOException;
import java.util.ArrayList;
import java.util.List;

import org.geotools.geojson.geom.GeometryJSON;
import org.locationtech.jts.geom.Coordinate;
import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.geom.GeometryFactory;
import org.locationtech.jts.geom.Point;
import org.locationtech.jts.geom.PrecisionModel;


public class GeoUtils {

    public class GeoCoordDTO {
        float lat;

        public float getLat() {
            return lat;
        }

        float lon;

        public float getLon() {
            return lon;
        }
    }

    public static String createGeoJsonCircle(
        GeoCoordDTO g1,
        double radiusInKm
    ) throws IOException {
        Geometry shape = new GeometryFactory().createMultiPoint(createCircle(64, g1, radiusInKm).toArray(new Point[] {}));
        GeometryJSON g = new GeometryJSON(17);
        ByteArrayOutputStream baos = new ByteArrayOutputStream();
        g.write(shape, baos);
        return baos.toString();
    }

    private static double toRadians(float angleInDegrees) {
        return (angleInDegrees * Math.PI) / 180;
    }

    private static double toDegrees(double angleInRadians) {
        return (angleInRadians * 180) / Math.PI;
    }


    private static Point createPoint(
        GeoCoordDTO c1,
        double distance,
        double bearing
    ) {
        int earthRadius = 6378137;
        double lat1 = toRadians(c1.getLat());
        double lon1 = toRadians(c1.getLon());
        double dByR = distance / earthRadius;
        double lat = Math.asin(Math.sin(lat1) * Math.cos(dByR) + Math.cos(lat1) * Math.sin(dByR) * Math.cos(bearing));
        double lon = lon1 + Math.atan2(Math.sin(bearing) * Math.sin(dByR) * Math.cos(lat1), Math.cos(dByR) - Math.sin(lat1) * Math.sin(lat));
        return new Point(new Coordinate(toDegrees(lon), toDegrees(lat)), new PrecisionModel(PrecisionModel.FLOATING), 0);
    }

    private static List<Point> createCircle(
        int noOfPoints,
        GeoCoordDTO g1,
        double radiusInKm
    ) {
        List<Point> coordinates = new ArrayList<>();
        for (int i = 0; i < noOfPoints; ++i) {
            coordinates.add(createPoint(g1, radiusInKm * 1000, (2 * Math.PI * -i) / noOfPoints));
        }
        coordinates.add(coordinates.get(0));
        return coordinates;
    }
}

Answer 2

    double latitude = 40.689234d;
    double longitude = -74.044598d;
    double diameterInMeters = 2000d; //2km

    GeometricShapeFactory shapeFactory = new GeometricShapeFactory();
    shapeFactory.setNumPoints(64); // adjustable
    shapeFactory.setCentre(new Coordinate(latitude, longitude));
    // Length in meters of 1° of latitude = always 111.32 km
    shapeFactory.setWidth(diameterInMeters/111320d);
    // Length in meters of 1° of longitude = 40075 km * cos( latitude ) / 360
    shapeFactory.setHeight(diameterInMeters / (40075000 * Math.cos(Math.toRadians(latitude)) / 360));

    Polygon circle = shapeFactory.createEllipse();

Answer 3

Your output circle is valid, it just happens to exceed the diameter of the Earth's surface so your GIS may have issues drawing it! The problem is that you are mixing degrees and meters indiscriminately and GeoTools has no clue what you want it to do.

You need to add some information about the coordinate reference system of the point to the program, and if that projection is geographic (i.e. in degrees), transform the problem to a projection that is in meters.

public Geometry bufferPoint(Measure<Double, Length> distance, CoordinateReferenceSystem origCRS, Geometry geom) {
    Geometry pGeom = geom;
    MathTransform toTransform, fromTransform = null;
    // reproject the geometry to a local projection
    Unit<Length> unit = distance.getUnit();
    if (!(origCRS instanceof ProjectedCRS)) {

      double x = geom.getCoordinate().x;
      double y = geom.getCoordinate().y;

      String code = "AUTO:42001," + x + "," + y;
      // System.out.println(code);
      CoordinateReferenceSystem auto;
      try {
        auto = CRS.decode(code);
        toTransform = CRS.findMathTransform(DefaultGeographicCRS.WGS84, auto);
        fromTransform = CRS.findMathTransform(auto, DefaultGeographicCRS.WGS84);
        pGeom = JTS.transform(geom, toTransform);
        unit = SI.METER;
      } catch (MismatchedDimensionException | TransformException | FactoryException e) {
        // TODO Auto-generated catch block
        e.printStackTrace();
      }

    } else {
      unit = (Unit<Length>) origCRS.getCoordinateSystem().getAxis(0).getUnit();

    }


    // buffer
    Geometry out = pGeom.buffer(distance.doubleValue(unit));
    Geometry retGeom = out;
    // reproject the geometry to the original projection
    if (!(origCRS instanceof ProjectedCRS)) {
      try {
        retGeom = JTS.transform(out, fromTransform);

      } catch (MismatchedDimensionException | TransformException e) {
        // TODO Auto-generated catch block
        e.printStackTrace();
      }
    }
    return retGeom;
  }

AUTO:42001,x,y is a special projection centred on the point x,y in meters that allows us to use the JTS buffer method which is easier than the circle operation you are using.

For your inputs this gives me an ellipse over New York, note this is expected and is due to the distorting effects of using unprojected Lat/Lon coordinates on a curved Earth.

You can call it using:

//Measure<Double, Length> dist = Measure.valueOf(50.0, SI.KILOMETER);
Measure<Double, Length> dist = Measure.valueOf(1.0, NonSI.MILE);
GeometryFactory gf = new GeometryFactory();
Point p = gf.createPoint(new Coordinate(-73.87,40.84));
buf.bufferPoint(dist, DefaultGeographicCRS.WGS84, p);