public static void iterateOverTriangles(final Polygon polygon, final Consumer<Geometry> action) { final double elevation = getContourCoordinates(polygon).averageZ(); final double sizeTol = FastMath.sqrt(polygon.getArea()) / 100.0; final DelaunayTriangulationBuilder dtb = new DelaunayTriangulationBuilder(); final PreparedGeometry buffered = PREPARED_GEOMETRY_FACTORY.create(polygon.buffer(sizeTol, 5, 0)); final Envelope3D env = Envelope3D.of(buffered.getGeometry()); try { dtb.setSites(polygon); dtb.setTolerance(sizeTol); applyToInnerGeometries(dtb.getTriangles(GEOMETRY_FACTORY), (gg) -> { final ICoordinates cc = getContourCoordinates(gg); if (cc.isCoveredBy(env) && buffered.covers(gg)) { cc.setAllZ(elevation); gg.geometryChanged(); action.accept(gg); } }); } catch (final LocateFailureException | ConstraintEnforcementException e) { final IScope scope = GAMA.getRuntimeScope(); GamaRuntimeException.warning("Impossible to triangulate: " + new WKTWriter().write(polygon), scope); iterateOverTriangles((Polygon) DouglasPeuckerSimplifier.simplify(polygon, 0.1), action); return; } }
private static IShape difference(final IScope scope, final IShape geom, final List<IShape> geoms, final PreparedGeometry ref) { if (geom == null) { return null; } IShape gR = new GamaShape(geom); for (final IShape g : geoms) { if (g != null && geom.intersects(g)) { gR = Spatial.Operators.minus(scope, gR, g); if (gR == null) return null; if (gR.getGeometries().size() > 1) { for (final IShape sh : gR.getGeometries()) { if (!ref.disjoint(sh.getInnerGeometry())) { gR = sh; break; } } } else if (ref.disjoint(gR.getInnerGeometry())) { return null; } } } return gR; }
public static void main(String[] args) throws Exception { Geometry circle = createCircle(); PreparedGeometry prepCircle = PreparedGeometryFactory.prepare(circle); int count = 0; int inCount = 0; for (int i = 0; i < MAX_ITER; i++) { count++; Point randPt = createRandomPoint(); if (prepCircle.intersects(randPt)) { inCount++; } //System.out.println("Approximation to PI: " + (4.0 * inCount / (double) count)); } double approxPi = 4.0 * inCount / (double) count; double approxDiffPct = 1.0 - approxPi / Math.PI; System.out.println("Approximation to PI: " + approxPi + " ( % difference from actual = " + 100 * approxDiffPct + " )" ); }
/** * Finds all {@link PreparedGeometry}s which intersect a given {@link Geometry} * * @param g the geometry to query by * @return a list of intersecting PreparedGeometrys */ public List intersects(Geometry g) { List result = new ArrayList(); List candidates = query(g); for (Iterator it = candidates.iterator(); it.hasNext(); ) { PreparedGeometry prepGeom = (PreparedGeometry) it.next(); if (prepGeom.intersects(g)) { result.add(prepGeom); } } return result; }
/** */ public List<VectorFeature> intersectingFeatures (Geometry geom) { final PreparedGeometry pGeom = PreparedGeometryFactory.prepare(geom); final List<VectorFeature> result = new ArrayList<VectorFeature>(); _spatialIndex.query(geom.getEnvelopeInternal(), new ItemVisitor() { public void visitItem (Object item) { VectorFeature feature = (VectorFeature)item; if (pGeom.intersects(feature.getGeometry())) { result.add(feature); } } }); return result; }
private static TimedPosition findRegionCrossingPoint(Shapefile region, Fix fix1, Fix fix2) { Coordinate[] coords = new Coordinate[] { new Coordinate(fix1.lon(), fix1.lat()), new Coordinate(fix2.lon(), fix2.lat()) }; LineString line = new GeometryFactory().createLineString(coords); for (PreparedGeometry g : region.geometries()) { if (g.crosses(line)) { Geometry intersection = g.getGeometry().intersection(line); // expecting just one point Coordinate coord = intersection.getCoordinate(); double longitude = coord.x; double latitude = coord.y; Position a = Position.create(fix1.lat(), fix1.lon()); Position b = Position.create(fix2.lat(), fix2.lon()); Position c = Position.create(latitude, longitude); double ac = a.getDistanceToKm(c); double bc = b.getDistanceToKm(c); if (ac == 0) { return new TimedPosition(fix1.lat(), fix1.lon(), fix1.time()); } else if (bc == 0) { return new TimedPosition(fix2.lat(), fix2.lon(), fix2.time()); } else { // predict the timestamp based on distance from a and b long diff = fix2.time() - fix1.time(); long t = Math.round(fix1.time() + ac * diff / (ac + bc)); return new TimedPosition(latitude, longitude, t); } } } throw new RuntimeException("crossing not found"); }
public static boolean contains(GeometryFactory gf, Collection<PreparedGeometry> geometries, double lat, double lon) { Point point = gf.createPoint(new Coordinate(lon, lat)); for (PreparedGeometry g : geometries) { if (g.contains(point)) return true; } return false; }
private Shapefile load() { if (state.compareAndSet(NOT_LOADED, LOADED)) { try { final List<PreparedGeometry> geometries = new ArrayList<>(); for (String typeName : datastore.getTypeNames()) { SimpleFeatureSource source = datastore.getFeatureSource(typeName); crs = source.getBounds().getCoordinateReferenceSystem(); final SimpleFeatureCollection features = source.getFeatures(); SimpleFeatureIterator it = features.features(); while (it.hasNext()) { SimpleFeature feature = it.next(); Geometry g = (Geometry) feature.getDefaultGeometry(); if (bufferDistance > 0) g = g.buffer(bufferDistance); geometries.add(PreparedGeometryFactory.prepare(g)); this.features.add(feature); } it.close(); } this.geometries = geometries; } catch (IOException e) { throw new RuntimeException(e); } } else if (state.get() == CLOSED) throw new RuntimeException("Shapefile is closed and can't be accessed"); return this; }
public Rect mbr() { // TODO assumes that shapefile is using WGS84? load(); Rect r = null; for (PreparedGeometry g : geometries) { Coordinate[] v = g.getGeometry().getEnvelope().getCoordinates(); System.out.println(Arrays.toString(v)); Rect rect = new Rect(v[0].y, v[0].x, v[2].y, v[2].x); if (r == null) r = rect; else r = r.add(rect); } return r; }
private static boolean contains(GeometryFactory gf, Collection<PreparedGeometry> geometries, double lat, double lon) { for (PreparedGeometry g : geometries) { if (g.contains(gf.createPoint(new Coordinate(lon, lat)))) return true; } return false; }
private void linkBasinWithNetwork() throws Exception { FeatureExtender fExt = new FeatureExtender(inNetwork.getSchema(), new String[]{NetworkChannel.NETNUMNAME}, new Class[]{Integer.class}); DefaultFeatureCollection newCollection = new DefaultFeatureCollection(); SimpleFeatureIterator hillslopeFeatures = inHillslope.features(); while( hillslopeFeatures.hasNext() ) { SimpleFeature hFeature = hillslopeFeatures.next(); Object netNum = hFeature.getAttribute(NetworkChannel.NETNUMNAME); Geometry hGeometry = (Geometry) hFeature.getDefaultGeometry(); PreparedGeometry preparedHGeometry = PreparedGeometryFactory.prepare(hGeometry); SimpleFeatureIterator netFeatures = inNetwork.features(); while( netFeatures.hasNext() ) { SimpleFeature nFeature = netFeatures.next(); Geometry geometry = (Geometry) nFeature.getDefaultGeometry(); if (geometry.getNumGeometries() != 1) { throw new ModelsRuntimeException("The network geometries have to be single lines.", this); } LineString nLine = (LineString) geometry.getGeometryN(0); Point startPoint = nLine.getStartPoint(); if (preparedHGeometry.contains(startPoint)) { SimpleFeature extendFeature = fExt.extendFeature(nFeature, new Object[]{netNum}); newCollection.add(extendFeature); break; } } } inNetwork = newCollection; }
@Override public boolean compare( final Geometry dataGeometry, final PreparedGeometry constraintGeometry ) { // This method is same as Geometry.equalsTopo which is // computationally expensive. // See equalsExact for quick structural equality return constraintGeometry.getGeometry().equals( dataGeometry); }
private PreparedFeature( final Feature feature, final PreparedGeometry preparedGeometry) { this.feature = feature; this.preparedGeometry = preparedGeometry; }
@Override public void setup( final Configuration configuration, final List<ColumnInterface> columnList) throws IOException { m_spatialIndex = new STRtree(); m_buffer = configuration.getFloat(GeoEnrichmentJob.KEY_BUFFER, 0.000001F); final URL url = getUrl(configuration); final ShapefileDataStoreFactory factory = new ShapefileDataStoreFactory(); final Map params = new HashMap(); params.put(ShapefileDataStoreFactory.URLP.key, url); final DataStore dataStore = factory.createDataStore(params); try { final String[] typeNames = dataStore.getTypeNames(); final SimpleFeatureSource featureSource = dataStore.getFeatureSource(typeNames[0]); m_geometryName = featureSource.getSchema().getGeometryDescriptor().getLocalName(); final SimpleFeatureCollection featureCollection = featureSource.getFeatures(); featureCollection.accepts(new FeatureVisitor() { public void visit(final Feature feature) { final Geometry geometry = (Geometry) feature.getProperty(m_geometryName).getValue(); final PreparedGeometry preparedGeometry = PreparedGeometryFactory.prepare(geometry); m_spatialIndex.insert(geometry.getEnvelopeInternal(), new PreparedFeature(feature, preparedGeometry)); } }, new NullProgressListener()); } finally { dataStore.dispose(); } }
private static Geometry intersectionSim( PreparedGeometry pg, Geometry g2) { PreparedPolygon ppoly = (PreparedPolygon) pg; FastSegmentSetIntersectionFinder intf = ppoly.getIntersectionFinder(); Coordinate[] pts = g2.getCoordinates(); List segStrings = SegmentStringUtil.extractSegmentStrings(g2); intf.intersects(segStrings ); return g2; }
private static PreparedGeometry cacheFind(Geometry g) { PreparedGeometry pg = cache.get(g); if (pg == null) { pg = PreparedGeometryFactory.prepare(g); cache.put(g, pg); //System.out.println("cache size = " + cache.size()); } return pg; }
/** * Creates a prepared geometry. <p> Prepared geometries make operations like intersection must faster. </p> */ public static PreparedGeometry prepare(Geometry g) { return PreparedGeometryFactory.prepare(g); }
@operator ( value = "masked_by", category = { IOperatorCategory.SPATIAL }, concept = { IConcept.GEOMETRY, IConcept.SPATIAL_COMPUTATION, IConcept.OBSTACLE }) @doc ( examples = { @example ( value = "perception_geom masked_by obstacle_list", equals = "the geometry representing the part of perception_geom visible from the agent position considering the list of obstacles obstacle_list.", isExecutable = false) }) public static IShape masked_by(final IScope scope, final IShape source, final IContainer<?, IShape> obstacles, Integer precision) { final IAgent a = scope.getAgent(); final List<IShape> obst = obstacles == null ? new ArrayList<IShape>() : obstacles.listValue(scope, Types.GEOMETRY, false); final ILocation location = a != null ? a.getLocation() : new GamaPoint(0, 0); if (precision == null) { precision = 120; } final Geometry visiblePercept = GeometryUtils.GEOMETRY_FACTORY.createGeometry(source.getInnerGeometry()); final boolean isPoint = source.isPoint(); if (obstacles != null && !obstacles.isEmpty(scope)) { final Geometry pt = GeometryUtils.GEOMETRY_FACTORY.createPoint(GeometryUtils.toCoordinate(location)); final Geometry locG = pt.buffer(0.01).getEnvelope(); double percep_dist = 0; for (final ILocation p : source.getPoints()) { final double dist = location.euclidianDistanceTo(p); if (dist > percep_dist) { percep_dist = dist; } } final Geometry gbuff = pt.buffer(percep_dist, precision / 4); final List<IShape> geoms = new ArrayList<IShape>(); for (int k = 1; k < gbuff.getNumPoints(); k++) { final IList coordinates = GamaListFactory.create(Types.POINT, 4); coordinates.add(location); coordinates.add(new GamaPoint(gbuff.getCoordinates()[k - 1])); coordinates.add(new GamaPoint(gbuff.getCoordinates()[k])); coordinates.add(location); final IShape gg = Spatial.Operators.inter(scope, source, Spatial.Creation.polygon(scope, coordinates)); if (gg != null && (isPoint || !gg.isPoint())) { geoms.add(new GamaShape(gg)); } } final GamaList geomsVisible = (GamaList) GamaListFactory.create(); final PreparedGeometry ref = PreparedGeometryFactory.prepare(locG); for (final IShape geom : geoms) { if (!intersection(geom, obst)) { geomsVisible.addValue(scope, geom); } else { final IShape perceptReal = difference(scope, geom, obst, ref); if (perceptReal != null && (isPoint || !perceptReal.isPoint())) { geomsVisible.addValue(scope, perceptReal); } } } if (!geomsVisible.isEmpty(scope)) { IShape result = Cast.asGeometry(scope, geomsVisible, false); if (result.getInnerGeometry() instanceof GeometryCollection) { result = Spatial.Transformations.enlarged_by(scope, result, 0.1); } return result; } return null; } return new GamaShape(visiblePercept); }
public List<PreparedGeometry> geometries() { load(); return geometries; }
@Execute public void process() throws Exception { checkNull(inRaster); ISurfaceInterpolator interpolator; if (pMode.equals(IDW)) { interpolator = new IDWInterpolator(pBuffer); } else { interpolator = new TPSInterpolator(pBuffer); } RegionMap regionMap = CoverageUtilities.getRegionParamsFromGridCoverage(inRaster); int rows = regionMap.getRows(); int cols = regionMap.getCols(); WritableRaster outWR = CoverageUtilities.renderedImage2WritableRaster(inRaster.getRenderedImage(), false); WritableRandomIter outIter = CoverageUtilities.getWritableRandomIterator(outWR); GridGeometry2D gridGeometry = inRaster.getGridGeometry(); PreparedGeometry preparedRoi = null; if (inROI != null) { List<Geometry> roiList = FeatureUtilities.featureCollectionToGeometriesList(inROI, false, null); GeometryCollection gc = new GeometryCollection(roiList.toArray(GeometryUtilities.TYPE_GEOMETRY), gf); preparedRoi = PreparedGeometryFactory.prepare(gc); } pm.beginTask("Filling holes...", cols - 2); for( int c = 1; c < cols - 1; c++ ) { for( int r = 1; r < rows - 1; r++ ) { if (pm.isCanceled()) { return; } double value = outIter.getSampleDouble(c, r, 0); if (isNovalue(value)) { DirectPosition worldPosition = gridGeometry.gridToWorld(new GridCoordinates2D(c, r)); double[] coordinate = worldPosition.getCoordinate(); Coordinate pointCoordinate = new Coordinate(coordinate[0], coordinate[1]); Point point = gf.createPoint(pointCoordinate); if (preparedRoi == null || preparedRoi.intersects(point)) { // TODO this could be done considering more points and more far away points. // For now, this works. List<Coordinate> surroundingValids = getValidSurroundingPoints(outIter, gridGeometry, c, r); if (surroundingValids.size() > 3) { double newValue = interpolator.getValue(surroundingValids.toArray(new Coordinate[0]), pointCoordinate); outIter.setSample(c, r, 0, newValue); } } } } pm.worked(1); } pm.done(); outIter.done(); outRaster = CoverageUtilities.buildCoverage("nulled", outWR, regionMap, inRaster.getCoordinateReferenceSystem()); }
@SuppressWarnings("unchecked") @Execute public void process() throws Exception { checkNull(inVector, pMaxOverlap, inRaster); RandomIter rasterIter = CoverageUtilities.getRandomIterator(inRaster); GridGeometry2D gridGeometry = inRaster.getGridGeometry(); double[] tm_utm_tac = new double[3]; STRtree circlesTree = FeatureUtilities.featureCollectionToSTRtree(inVector); List<SimpleFeature> circlesList = FeatureUtilities.featureCollectionToList(inVector); DefaultFeatureCollection outFC = new DefaultFeatureCollection(); for( SimpleFeature circleFeature : circlesList ) { Geometry geometry = (Geometry) circleFeature.getDefaultGeometry(); Polygon circle = (Polygon) geometry.getGeometryN(0); PreparedGeometry preparedCircle = PreparedGeometryFactory.prepare(circle); List<SimpleFeature> circlesAround = circlesTree.query(circle.getEnvelopeInternal()); List<Geometry> intersectedCircles = new ArrayList<Geometry>(); for( SimpleFeature circleAround : circlesAround ) { if (circleAround.equals(circleFeature)) { continue; } Geometry circleAroundGeometry = (Geometry) circleAround.getDefaultGeometry(); if (preparedCircle.intersects(circleAroundGeometry)) { intersectedCircles.add(circleAroundGeometry); } } Point centroid = circle.getCentroid(); int intersectionsCount = intersectedCircles.size(); if (intersectionsCount != 0) { // check how many circles overlapped if (intersectionsCount > pMaxOverlapCount) { continue; } // check if the circles overlap too much, i.e. cover their baricenter boolean intersected = false; for( Geometry intersectedCircle : intersectedCircles ) { if (intersectedCircle.intersects(centroid)) { intersected = true; break; } } if (intersected) { continue; } } // check if the center has a raster value, i.e. is not empty double value = CoverageUtilities.getValue(inRaster, centroid.getCoordinate()); if (!HMConstants.isNovalue(value)) { continue; } // check if the inner part of the circle is indeed rather empty // min, max, mean, var, sdev, activeCellCount, passiveCellCount double[] stats = OmsZonalStats.polygonStats(circle, gridGeometry, rasterIter, false, tm_utm_tac, 0, pm); // if we have many more active cells than passive cells, that is not a circle double activeCells = stats[5]; double novalues = stats[6]; if (activeCells * 1.5 > novalues) { continue; } // take it as valid circle outFC.add(circleFeature); } outCircles = outFC; rasterIter.done(); }
@Execute public void process() throws Exception { checkNull(inMap1, inMap2); outMap = new DefaultFeatureCollection(); if (!doKeepFirstAttributes) { SimpleFeatureCollection inMapTmp = inMap1; inMap1 = inMap2; inMap2 = inMapTmp; } List<Geometry> geometries = FeatureUtilities.featureCollectionToGeometriesList(inMap2, false, null); GeometryCollection geometryCollection = new GeometryCollection(geometries.toArray(new Geometry[geometries.size()]), gf); Geometry intersectionGeometry = geometryCollection.buffer(0); PreparedGeometry preparedIntersectionGeometry = PreparedGeometryFactory.prepare(intersectionGeometry); List<SimpleFeature> mainFeatures = FeatureUtilities.featureCollectionToList(inMap1); if (mainFeatures.size() == 0) { throw new ModelsIllegalargumentException("No features found in the layer.", this); } EGeometryType geometryType = EGeometryType.forGeometry((Geometry) mainFeatures.get(0).getDefaultGeometry()); Class< ? > multiClazz = geometryType.getMultiClazz(); EGeometryType newGeometryType = EGeometryType.forClass(multiClazz); FeatureGeometrySubstitutor sub = new FeatureGeometrySubstitutor(inMap1.getSchema(), multiClazz); pm.beginTask("Performing intersection...", mainFeatures.size()); for( SimpleFeature feature : mainFeatures ) { Geometry geometry = (Geometry) feature.getDefaultGeometry(); if (preparedIntersectionGeometry.intersects(geometry)) { Geometry intersection = geometry.intersection(intersectionGeometry); EGeometryType intersectionGeometryType = EGeometryType.forGeometry(intersection); if (intersectionGeometryType.isCompatibleWith(newGeometryType)) { SimpleFeature newFeature = sub.substituteGeometry(feature, intersection); ((DefaultFeatureCollection) outMap).add(newFeature); } else { pm.errorMessage("Could not add intersection result geometry to layer due to incompatibility: " + intersection); } } pm.worked(1); } pm.done(); }
/** * Retrieve all the trees envelopes that intersect the geometry. * * @param checkGeom the {@link com.vividsolutions.jts.geom.Geometry} to use to check. * @param doOnlyEnvelope check for the geom envelope instead of a intersection with it. * @param minMaxZ an array to be filled with the min and max z to be used as style. * @return the list of envelopes contained in the supplied geometry. * @throws Exception */ @Override public synchronized List<Geometry> getEnvelopesInGeometry( Geometry checkGeom, boolean doOnlyEnvelope, double[] minMaxZ ) throws Exception { checkOpen(); ArrayList<Geometry> envelopeListForTile = new ArrayList<Geometry>(); Envelope env = checkGeom.getEnvelopeInternal(); PreparedGeometry preparedGeometry = null; if (!doOnlyEnvelope) { preparedGeometry = PreparedGeometryFactory.prepare(checkGeom); } double min = Double.POSITIVE_INFINITY; double max = Double.NEGATIVE_INFINITY; List< ? > filesList = mainLasFolderIndex.query(env); for( Object fileName : filesList ) { if (fileName instanceof String) { String name = (String) fileName; File lasFile = new File(lasFolder, name); File lasIndexFile = FileUtilities.substituteExtention(lasFile, "lasfix"); String absolutePath = lasIndexFile.getAbsolutePath(); STRtreeJGT lasIndex = fileName2IndexMap.get(absolutePath); if (lasIndex == null) { lasIndex = OmsLasIndexReader.readIndex(absolutePath); fileName2IndexMap.put(absolutePath, lasIndex); } List< ? > queryBoundables = lasIndex.queryBoundables(env); for( Object object : queryBoundables ) { if (object instanceof ItemBoundable) { ItemBoundable itemBoundable = (ItemBoundable) object; double[] item = (double[]) itemBoundable.getItem(); if (item.length > 0) { Envelope bounds = (Envelope) itemBoundable.getBounds(); Polygon envelopePolygon = LasIndexer.envelopeToPolygon(bounds); envelopePolygon.setUserData(new double[]{item[2], item[3]}); if (minMaxZ != null) { min = Math.min(min, item[2]); max = Math.max(max, item[2]); } if (doOnlyEnvelope) { envelopeListForTile.add(envelopePolygon); } else { if (preparedGeometry.intersects(envelopePolygon)) { envelopeListForTile.add(envelopePolygon); } } } } } } } if (minMaxZ != null) { minMaxZ[0] = min; minMaxZ[1] = max; } return envelopeListForTile; }
private PreparedGeometry getFloodingArea( SimpleFeatureCollection inFloodingAreas ) throws Exception { List<Geometry> geometriesList = FeatureUtilities.featureCollectionToGeometriesList(inFloodingAreas, true, null); Geometry polygonUnion = CascadedPolygonUnion.union(geometriesList); PreparedGeometry preparedGeometry = PreparedGeometryFactory.prepare(polygonUnion); return preparedGeometry; }
public boolean compare( final Geometry dataGeometry, final PreparedGeometry constraintGeometry );
@Override public boolean compare( final Geometry dataGeometry, final PreparedGeometry constraintGeometry ) { return constraintGeometry.contains(dataGeometry); }
@Override public boolean compare( final Geometry dataGeometry, final PreparedGeometry constraintGeometry ) { return constraintGeometry.overlaps(dataGeometry); }
@Override public boolean compare( final Geometry dataGeometry, final PreparedGeometry constraintGeometry ) { return constraintGeometry.intersects(dataGeometry); }
@Override public boolean compare( final Geometry dataGeometry, final PreparedGeometry constraintGeometry ) { return constraintGeometry.touches(dataGeometry); }
@Override public boolean compare( final Geometry dataGeometry, final PreparedGeometry constraintGeometry ) { return constraintGeometry.within(dataGeometry); }
@Override public boolean compare( final Geometry dataGeometry, final PreparedGeometry constraintGeometry ) { return constraintGeometry.disjoint(dataGeometry); }
@Override public boolean compare( final Geometry dataGeometry, final PreparedGeometry constraintGeometry ) { return constraintGeometry.crosses(dataGeometry); }
public GeometryImage( final PreparedGeometry preparedGeometry ) { super(); this.preparedGeometry = preparedGeometry; geometryBinary = GeometryUtils.geometryToBinary(preparedGeometry.getGeometry()); }
public PreparedGeometry getGeometry() { return preparedGeometry; }
