Search code examples
javaalgorithmgeotools

Building GeoTools Geometry "segments" from route coordinates


From a set of coordinates that define a route, I want to draw a Geometry that mimicks a theoretical highway of that track, given an arbitrary number of meters wide (e.g. 20)

I don't know if GeoTools provides tools for constructing a Geometry with such inputs, so my initial idea is to split the track coordinates (several thousands) in pairs (coord0, coord1), (coord1, coord2), ...., (coordN, coordN-1) and, with each pair, build a rectangle assuming that the two points are the midpoints of a 20m wide segment (as in Knowing two points of a rectangle, how can I figure out the other two?), and joining all the resulting geometries.

Maybe it's overkill but I haven't found a cheaper way to do this

Any ideas would be greatly appreciated!


Solution

  • The easy way to do this is to use a 20m buffer around a line created from your points. So some code like this to create a line from the points (:

    String[] wkt = {
        "Point (-0.13666168754467312 50.81919869153657743)",
        "Point (-0.13622277073931291 50.82205165077141373)",
        "Point (-0.13545466632993253 50.82512406840893959)",
        "Point (-0.13457683271921211 50.82687973563037787)",
        "Point (-0.13413791591385191 50.82907431965718104)",
        "Point (-0.13951464677951447 50.8294035072611976)",
        "Point (-0.14346489802775639 50.83082998687861931)",
        "Point (-0.14697623247063807 50.83072025767727808)",
        "Point (-0.15004865010815954 50.83390240451614517)",
        "Point (-0.15740050659794308 50.8349996965295432)",
        "Point (-0.16486209228906662 50.83741373895902171)",
        "Point (-0.17276259478555042 50.83894994777778464)",
        "Point (-0.18549118214099652 50.8387304893751022)"
        };
    
        //build line
        WKTReader2 reader = new WKTReader2();
        GeometryFactory gf = new GeometryFactory();
        Coordinate[] points = new Coordinate[wkt.length];
        int i=0;
        for(String w:wkt) {
          Point p;
      try {
        p = (Point) reader.read(w);
        points[i++]=p.getCoordinate();
      } catch (ParseException e) {
        // TODO Auto-generated catch block
        e.printStackTrace();
      }
    
        }
        LineString line = gf.createLineString(points);
        SimpleFeatureBuilder builder = new SimpleFeatureBuilder(schema);
        builder.set("locations", line);
        SimpleFeature feature = builder.buildFeature("1");
    

    And then a BufferLine method like:

    public SimpleFeature bufferFeature(SimpleFeature feature, Measure<Double, Length> distance) {
        // extract the geometry
        GeometryAttribute gProp = feature.getDefaultGeometryProperty();
        CoordinateReferenceSystem origCRS = gProp.getDescriptor().getCoordinateReferenceSystem();
    
        Geometry geom = (Geometry) feature.getDefaultGeometry();
        Geometry pGeom = geom;
        MathTransform toTransform, fromTransform = null;
        // reproject the geometry to a local projection
        if (!(origCRS instanceof ProjectedCRS)) {
    
            Point c = geom.getCentroid();
            double x = c.getCoordinate().x;
            double y = c.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);
            } catch (MismatchedDimensionException | TransformException | FactoryException e) {
                // TODO Auto-generated catch block
                e.printStackTrace();
            }
    
        }
    
        // buffer
        Geometry out = buffer(pGeom, distance.doubleValue(SI.METER));
        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 a new feature containing the geom
        SimpleFeatureType schema = feature.getFeatureType();
        SimpleFeatureTypeBuilder ftBuilder = new SimpleFeatureTypeBuilder();
        ftBuilder.setCRS(origCRS);
        // ftBuilder.setDefaultGeometry("buffer");
        ftBuilder.addAll(schema.getAttributeDescriptors());
        ftBuilder.setName(schema.getName());
    
        SimpleFeatureType nSchema = ftBuilder.buildFeatureType();
        SimpleFeatureBuilder builder = new SimpleFeatureBuilder(nSchema);
        List<Object> atts = feature.getAttributes();
        for (int i = 0; i < atts.size(); i++) {
            if (atts.get(i) instanceof Geometry) {
                atts.set(i, retGeom);
            }
        }
        SimpleFeature nFeature = builder.buildFeature(null, atts.toArray());
        return nFeature;
    }
    
    /**
     * create a buffer around the geometry, assumes the geometry is in the same
     * units as the distance variable.
     * 
     * @param geom
     *          a projected geometry.
     * @param dist
     *          a distance for the buffer in the same units as the projection.
     * @return
     */
    private Geometry buffer(Geometry geom, double dist) {
    
        Geometry buffer = geom.buffer(dist);
    
        return buffer;
    
    }
    

    The tricky part is reprojecting into a locally flat CRS so that you can use metres for the buffer size. If you know of a locally good projection you could just use that (in this case we could have used OSGB (EPSG:27700) for better results).

    This gives the following map:

    enter image description here