Snap a Point to a Line¶
This example illustrates a common spatial operation: moving or snapping a point to lie on a nearby line. For instance, you might use this approach to align point locations from a mobile device with a map of roads. The JTS library is used for this sort of task day in and day out.
What you will learn:
Use of a spatial index to cache features in memory and search efficiently.
Going beyond the familiar JTS Geometry class methods and making direct use of other classes.
Related material:
To begin, we prompt the user for a shapefile containing line features and connect to it:
import java.io.File;
import java.util.List;
import java.util.Random;
import org.geotools.api.data.FeatureSource;
import org.geotools.api.data.FileDataStore;
import org.geotools.api.data.FileDataStoreFinder;
import org.geotools.api.feature.Feature;
import org.geotools.api.feature.FeatureVisitor;
import org.geotools.api.feature.simple.SimpleFeature;
import org.geotools.data.util.NullProgressListener;
import org.geotools.feature.FeatureCollection;
import org.geotools.geometry.jts.ReferencedEnvelope;
import org.geotools.swing.data.JFileDataStoreChooser;
import org.geotools.util.SuppressFBWarnings;
import org.locationtech.jts.geom.Coordinate;
import org.locationtech.jts.geom.Envelope;
import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.geom.LineString;
import org.locationtech.jts.geom.MultiLineString;
import org.locationtech.jts.index.SpatialIndex;
import org.locationtech.jts.index.strtree.STRtree;
import org.locationtech.jts.linearref.LinearLocation;
import org.locationtech.jts.linearref.LocationIndexedLine;
@SuppressFBWarnings("DMI_RANDOM_USED_ONLY_ONCE")
public class SnapToLine {
public static void main(String[] args) throws Exception {
/*
* Open a shapefile. You should choose one with line features
* (LineString or MultiLineString geometry)
*
*/
File file = JFileDataStoreChooser.showOpenFile("shp", null);
if (file == null) {
return;
}
FileDataStore store = FileDataStoreFinder.getDataStore(file);
FeatureSource source = store.getFeatureSource();
// Check that we have line features
Class<?> geomBinding = source.getSchema().getGeometryDescriptor().getType().getBinding();
boolean isLine =
geomBinding != null
&& (LineString.class.isAssignableFrom(geomBinding)
|| MultiLineString.class.isAssignableFrom(geomBinding));
if (!isLine) {
System.out.println("This example needs a shapefile with line features");
return;
}
You might be used to working with shapefiles as a streaming data source, i.e. reading features from disk as required. Here we optimize things by extracting the line geometries from the features and
caching them in a JTS SpatialIndex object. This gives us speed in two ways: we have the lines in
memory and can search for them efficiently by location:
final SpatialIndex index = new STRtree();
FeatureCollection features = source.getFeatures();
System.out.println("Slurping in features ...");
features.accepts(
new FeatureVisitor() {
@Override
public void visit(Feature feature) {
SimpleFeature simpleFeature = (SimpleFeature) feature;
Geometry geom = (MultiLineString) simpleFeature.getDefaultGeometry();
// Just in case: check for null or empty geometry
if (geom != null) {
Envelope env = geom.getEnvelopeInternal();
if (!env.isNull()) {
index.insert(env, new LocationIndexedLine(geom));
}
}
}
},
new NullProgressListener());
Notice that we wrapped each feature’s line geometry in a JTS LocationIndexedLine object which we will use to find the point on a line closest to a reference point. We could have loaded the lines directly into the spatial index, but this way we will avoid wrapping each line every time it is tested against a point.
Now let’s make some pretend point data:
/*
* For test data, we generate a large number of points placed randomly
* within the bounding rectangle of the features.
*/
final int NUM_POINTS = 10000;
ReferencedEnvelope bounds = features.getBounds();
Coordinate[] points = new Coordinate[NUM_POINTS];
Random rand = new Random(file.hashCode());
for (int i = 0; i < NUM_POINTS; i++) {
points[i] =
new Coordinate(
bounds.getMinX() + rand.nextDouble() * bounds.getWidth(),
bounds.getMinY() + rand.nextDouble() * bounds.getHeight());
}
At last we are ready to snap points. We create a search envelope of fixed size around each point and use this to query the lines in the spatial index.
In case your shapefile is large, we’ll set a time limit on how long snapping continues:
/*
* We defined the maximum distance that a line can be from a point
* to be a candidate for snapping (1% of the width of the feature
* bounds for this example).
*/
final double MAX_SEARCH_DISTANCE = bounds.getSpan(0) / 100.0;
// Maximum time to spend running the snapping process (milliseconds)
final long DURATION = 5000;
int pointsProcessed = 0;
int pointsSnapped = 0;
long elapsedTime = 0;
long startTime = System.currentTimeMillis();
while (pointsProcessed < NUM_POINTS
&& (elapsedTime = System.currentTimeMillis() - startTime) < DURATION) {
// Get point and create search envelope
Coordinate pt = points[pointsProcessed++];
Envelope search = new Envelope(pt);
search.expandBy(MAX_SEARCH_DISTANCE);
/*
* Query the spatial index for objects within the search envelope.
* Note that this just compares the point envelope to the line envelopes
* so it is possible that the point is actually more distant than
* MAX_SEARCH_DISTANCE from a line.
*/
@SuppressWarnings("unchecked")
List<LocationIndexedLine> lines = index.query(search);
// Initialize the minimum distance found to our maximum acceptable
// distance plus a little bit
double minDist = MAX_SEARCH_DISTANCE + 1.0e-6;
Coordinate minDistPoint = null;
for (LocationIndexedLine line : lines) {
LinearLocation here = line.project(pt);
Coordinate point = line.extractPoint(here);
double dist = point.distance(pt);
if (dist < minDist) {
minDist = dist;
minDistPoint = point;
}
}
if (minDistPoint == null) {
// No line close enough to snap the point to
System.out.println(pt + "- X");
} else {
System.out.printf("%s - snapped by moving %.4f\n", pt.toString(), minDist);
pointsSnapped++;
}
}
System.out.printf(
"Processed %d points (%.2f points per second). \n" + "Snapped %d points.\n\n",
pointsProcessed, 1000.0 * pointsProcessed / elapsedTime, pointsSnapped);
}
}
You can experiment with this code:
try some actual roads - using real data makes a difference
try using a
QuadTree- it is often much slower (but you can add and remove things from aQuadTreeat runtime)are you getting wacky results? Check if your
geometry.isValid()prior to using ittry simplifying the lines prior to creating the
LocationIndexedLine- it should be much fasterIf you are unsure how to do activities listed above consult the Secrets of JTS link provided at the top of the page.