Many resources are needed to download a project. Please understand that we have to compensate our server costs. Thank you in advance.
Project price only 1 $
You can buy this project and download/modify it how often you want.
eu.europa.ec.eurostat.jgiscotools.algo.base.NodingUtil Maven / Gradle / Ivy
/**
*
*/
package eu.europa.ec.eurostat.jgiscotools.algo.base;
import java.util.ArrayList;
import java.util.Collection;
import java.util.HashSet;
import java.util.List;
import org.apache.logging.log4j.LogManager;
import org.apache.logging.log4j.Logger;
import org.locationtech.jts.geom.Coordinate;
import org.locationtech.jts.geom.Envelope;
import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.geom.LineSegment;
import org.locationtech.jts.geom.LineString;
import org.locationtech.jts.geom.LinearRing;
import org.locationtech.jts.geom.MultiLineString;
import org.locationtech.jts.geom.MultiPoint;
import org.locationtech.jts.geom.MultiPolygon;
import org.locationtech.jts.geom.Point;
import org.locationtech.jts.geom.Polygon;
import org.locationtech.jts.geom.TopologyException;
import org.locationtech.jts.index.SpatialIndex;
import org.locationtech.jts.index.quadtree.Quadtree;
import org.locationtech.jts.index.strtree.STRtree;
import eu.europa.ec.eurostat.jgiscotools.feature.Feature;
import eu.europa.ec.eurostat.jgiscotools.feature.FeatureUtil;
import eu.europa.ec.eurostat.jgiscotools.feature.JTSGeomUtil;
/**
* @author julien Gaffuri
*
*/
public class NodingUtil {
public final static Logger LOGGER = LogManager.getLogger(NodingUtil.class.getName());
public enum NodingIssueType { PointPoint, LinePoint }
public static class NodingIssue{
public NodingIssueType type;
public Coordinate c;
public double distance;
public NodingIssue(NodingIssueType type, Coordinate c, double distance) { this.type=type; this.c=c; this.distance=distance; }
public String toString() { return "NodingIssue "+type+" c="+c+" d="+distance; }
}
private static boolean checkLPNodingIssue(Coordinate c, Coordinate c1, Coordinate c2, double nodingResolution) {
return
c.distance(c1) > nodingResolution
&& c.distance(c2) > nodingResolution
&& new LineSegment(c1,c2).distance(c) <= nodingResolution
;
}
private static boolean checkPPNodingIssue(Coordinate c, Coordinate c_, double nodingResolution) {
double d = c.distance(c_);
return(d!=0 && d <= nodingResolution);
}
//get noding issues for multi-polygonal features
public static Collection getNodingIssues(NodingIssueType type, Collection mpfs, double nodingResolution) {
STRtree index = type==NodingIssueType.LinePoint? FeatureUtil.getIndexSTRtreeCoordinates(mpfs) : getSTRtreeCoordinatesForPP(mpfs, nodingResolution);
Collection nis = new HashSet();
for(Feature mpf : mpfs)
nis.addAll(getNodingIssues(type, mpf, index, nodingResolution));
return nis;
}
public static Collection getNodingIssues(NodingIssueType type, Feature mpf, SpatialIndex index, double nodingResolution) {
return getNodingIssues(type, (MultiPolygon)mpf.getGeometry(), index, nodingResolution);
}
public static Collection getNodingIssues(NodingIssueType type, MultiPolygon mp, SpatialIndex index, double nodingResolution) {
Collection out = new HashSet();
for(int i=0; i getNodingIssues(NodingIssueType type, Polygon p, SpatialIndex index, double nodingResolution) {
Collection out = new HashSet();
for(LineString lr : JTSGeomUtil.getRings(p))
out.addAll( getNodingIssues(type,lr,index,nodingResolution) );
return out;
}
private static Collection getNodingIssues(NodingIssueType type, LineString ls, SpatialIndex index, double nodingResolution) {
Collection out = new HashSet();
if(type == NodingIssueType.LinePoint ) {
//go through segments of l1
Coordinate[] c1s = ls.getCoordinates();
Coordinate c1 = c1s[0], c2;
for(int i=1; i)index.query(new Envelope(c1,c2))) {
NodingIssue ni = getLinePointNodingIssues(c,c1,c2,nodingResolution);
if(ni != null) out.add(ni);
}
c1 = c2;
}
} else if(type == NodingIssueType.PointPoint ) {
//go through coordinates of l1
Coordinate[] c1s = ls.getCoordinates();
for(int i=0; i)index.query(env)) {
NodingIssue ni = getPointPointNodingIssues(c,c_,nodingResolution);
if(ni != null) out.add(ni);
}
}
}
return out;
}
public static NodingIssue getLinePointNodingIssues(Coordinate c, Coordinate c1, Coordinate c2, double nodingResolution) {
if( checkLPNodingIssue(c,c1,c2,nodingResolution) )
return new NodingIssue(NodingIssueType.LinePoint,c,new LineSegment(c1,c2).distance(c));
else
return null;
}
public static NodingIssue getPointPointNodingIssues(Coordinate c, Coordinate c_, double nodingResolution) {
if( checkPPNodingIssue(c,c_,nodingResolution) )
return new NodingIssue(NodingIssueType.PointPoint,c,c.distance(c_));
else
return null;
}
public static void fixNoding(NodingIssueType type, Collection mpfs, double nodingResolution) {
STRtree index = type==NodingIssueType.LinePoint? FeatureUtil.getIndexSTRtreeCoordinates(mpfs) : getSTRtreeCoordinatesForPP(mpfs, nodingResolution);
for(Feature mpf : mpfs)
fixNoding(type, mpf, index, nodingResolution);
}
private static void fixNoding(NodingIssueType type, Feature mpf, SpatialIndex index, double nodingResolution) {
MultiPolygon mp = fixNoding(type, (MultiPolygon) mpf.getGeometry(), index, nodingResolution);
mpf.setGeometry(mp);
}
public static MultiPolygon fixNoding(NodingIssueType type, MultiPolygon mp, SpatialIndex index, double nodingResolution) {
Polygon[] ps = new Polygon[mp.getNumGeometries()];
for(int i=0; i nis = getNodingIssues(type, ls, index, nodingResolution);
while(nis.size() != 0) {
NodingIssue ni = nis.iterator().next();
out = fixNoding(ni.type, out, ni.c, nodingResolution);
nis = getNodingIssues(type, out, index, nodingResolution);
}
return out;
}
public static LineString fixNoding(NodingIssueType type, LineString ls, Coordinate c, double nodingResolution) {
LineString out = null;
if(type == NodingIssueType.PointPoint)
out = fixPPNoding(ls, c, nodingResolution);
if(type == NodingIssueType.LinePoint)
out = fixLPNoding(ls, c, nodingResolution);
return out;
}
//fix a noding issue by moving a coordinate (or several for closed lines) to a target position
public static LineString fixPPNoding(LineString ls, Coordinate c, double nodingResolution) {
Coordinate[] cs = ls.getCoordinates();
Coordinate[] csOut = new Coordinate[cs.length];
boolean found = false;
for(int i=0; i fs, double nodingResolution) {
Collection geoms = new HashSet();
for(Feature f : fs) geoms.add(f.getGeometry());
return getSTRtreeCoordinatesForPPG(geoms, nodingResolution);
}
/*private static SpatialIndex getSTRtreeCoordinatesForPP(double nodingResolution, Geometry... geoms) {
Collection gs = new HashSet();
for(Geometry g : geoms) gs.add(g);
return getSTRtreeCoordinatesForPPG(gs, nodingResolution);
}*/
private static STRtree getSTRtreeCoordinatesForPPG(Collection gs, double nodingResolution) {
//build index of all coordinates, ensuring newly added coordinates are not within a radius of nodingResolution of other ones.
Quadtree index = new Quadtree();
boolean found;
for(Geometry g : gs) {
for(Coordinate c : g.getCoordinates()) {
//try to find a coordinate of the index within nodingResolution radius of c
found = false;
Envelope env = new Envelope(c); env.expandBy(nodingResolution*1.01);
for(Coordinate c2 : (List)index.query(env )) {
if(c.distance(c2) > nodingResolution) continue;
found = true;
break;
}
if(found) continue;
index.insert(new Envelope(c), c);
}
}
STRtree index_ = new STRtree();
for(Coordinate c : (List)index.queryAll()) index_.insert(new Envelope(c), c);
return index_;
}
//public static void main(String[] args) {
//LOGGER.info("Start");
/*
double nodingResolution = 1e-3;
Polygon p1 = JTSGeomUtil.createPolygon(0,0, 1,0, 0,1, 0,0);
Polygon p2 = JTSGeomUtil.createPolygon(1,0, 0.5,0.5, 1,1, 1,0);
//SpatialIndex index = getCoordinatesSpatialIndex(p1, p2);
SpatialIndex index = getSTRtreeCoordinatesForPP(nodingResolution, p1, p2);
System.out.println(p1);
System.out.println(p2);
for(NodingIssue ni : getNodingIssues(NodingIssueType.LinePoint, p1, index, nodingResolution)) System.out.println(ni);
p1 = fixNoding(NodingIssueType.LinePoint, p1, index, nodingResolution);
p2 = fixNoding(NodingIssueType.LinePoint, p2, index, nodingResolution);
System.out.println(p1);
System.out.println(p2);
for(NodingIssue ni : getNodingIssues(NodingIssueType.LinePoint, p1, index, nodingResolution)) System.out.println(ni);
*/
/*
Polygon p1 = JTSGeomUtil.createPolygon(0,1, 0,0, 1.00001,0, 0,1);
Polygon p2 = JTSGeomUtil.createPolygon(1,0, 0,1, 1,1, 1,0);
SpatialIndex index = getSTRtreeCoordinatesForPP(nodingResolution, p1, p2);
System.out.println(p1);
System.out.println(p2);
for(NodingIssue ni : getNodingIssues(NodingIssueType.PointPoint, p1, index, nodingResolution)) System.out.println(ni);
p1 = fixNoding(NodingIssueType.PointPoint, p1, index, nodingResolution);
p2 = fixNoding(NodingIssueType.PointPoint, p2, index, nodingResolution);
System.out.println(p1);
System.out.println(p2);
for(NodingIssue ni : getNodingIssues(NodingIssueType.PointPoint, p1, index, nodingResolution)) System.out.println(ni);
*/
//LOGGER.info("End");
//}
//node features with linear geoemtries intersecting
public static void fixLineStringsIntersectionNoding(Collection fs) {
//make spatial index
Quadtree si = FeatureUtil.getIndexQuadtree(fs);
boolean b;
//go through pairs of features to check their intersection
for(Feature f1 : fs) {
Geometry g1 = f1.getGeometry();
for(Object f2_ : si.query(g1.getEnvelopeInternal())) {
if(f1==f2_) continue;
Feature f2 = (Feature) f2_;
if(f1.getID().compareTo(f2.getID()) < 0) continue;
Geometry g2 = f2.getGeometry();
if(! g1.getEnvelopeInternal().intersects(g2.getEnvelopeInternal())) continue;
Geometry inter = null;
try {
inter = g1.intersection(g2);
} catch (Exception e) {
try {
inter = g2.intersection(g1);
} catch (TopologyException e1) {
LOGGER.error("Could not compute intersection");
LOGGER.error(e1.getMessage());
continue;
}
}
if(inter.isEmpty()) continue;
//case when intersection has line component
if(inter.getLength() != 0) {
g1 = g1.difference(g2);
inter = g1.intersection(g2);
}
if(inter.isEmpty()) continue;
//intersection can only be with points
Coordinate[] cs = inter.getCoordinates();
for(Coordinate c : cs) {
g1 = insertCoordinate(g1, c);
g2 = insertCoordinate(g2, c);
}
//update both features
b = si.remove(f1.getGeometry().getEnvelopeInternal(), f1);
if(!b) LOGGER.warn("Error when removing feature in spatial index - fixLineStringsIntersectionNoding (1)");
b = si.remove(f2.getGeometry().getEnvelopeInternal(), f2);
if(!b) LOGGER.warn("Error when removing feature in spatial index - fixLineStringsIntersectionNoding (2)");
f1.setGeometry(g1);
f2.setGeometry(g2);
si.insert(f1.getGeometry().getEnvelopeInternal(), f1);
si.insert(f2.getGeometry().getEnvelopeInternal(), f2);
}
}
}
public static Geometry insertCoordinate(Geometry g, Coordinate c) {
if(g instanceof Point)
return insertCoordinate((Point)g, c);
if(g instanceof MultiPoint)
return insertCoordinate((MultiPoint)g, c);
if(g instanceof LineString)
return insertCoordinate((LineString)g, c);
if(g instanceof MultiLineString)
return insertCoordinate((MultiLineString)g, c);
LOGGER.warn("Method insertCoordinate not supported for geometry type "+g.getClass().getSimpleName());
return g;
}
public static Geometry insertCoordinate(Point p, Coordinate c) {
return p.union(p.getFactory().createPoint(c));
}
public static Geometry insertCoordinate(MultiPoint mp, Coordinate c) {
return mp.union(mp.getFactory().createPoint(c));
}
public static LineString insertCoordinate(LineString ls, Coordinate c) {
double d = ls.distance(ls.getFactory().createPoint(c));
if(d == 0) return ls;
//System.out.println(d + " " + c);
return fixLPNoding(ls, c, d * 1.1);
}
public static MultiLineString insertCoordinate(MultiLineString mls, Coordinate c) {
if(mls == null || mls.isEmpty()) return mls;
//get component closest to p
ArrayList lss = new ArrayList<>( JTSGeomUtil.getGeometries(mls) );
Geometry lsMin = null;
double dMin = Double.MAX_VALUE;
Point pt = mls.getFactory().createPoint(c);
for(Geometry ls : lss) {
double d = ls.distance(pt);
if(d
