Search in sources :

Example 11 with MismatchedDimensionException

use of org.opengis.geometry.MismatchedDimensionException in project jena by apache.

the class SpatialIndex method getGeometryLiteralIndexItems.

/**
 * @param model
 * @param srsURI
 * @return GeometryLiteral items prepared for adding to SpatialIndex.
 * @throws SpatialIndexException
 */
private static Collection<SpatialIndexItem> getGeometryLiteralIndexItems(Model model, String srsURI) throws SpatialIndexException {
    List<SpatialIndexItem> items = new ArrayList<>();
    StmtIterator stmtIt = model.listStatements(null, Geo.HAS_GEOMETRY_PROP, (Resource) null);
    while (stmtIt.hasNext()) {
        Statement stmt = stmtIt.nextStatement();
        Resource feature = stmt.getSubject();
        Resource geometry = stmt.getResource();
        ExtendedIterator<RDFNode> nodeIter = model.listObjectsOfProperty(geometry, Geo.HAS_SERIALIZATION_PROP);
        if (!nodeIter.hasNext()) {
            NodeIterator wktNodeIter = model.listObjectsOfProperty(geometry, Geo.AS_WKT_PROP);
            NodeIterator gmlNodeIter = model.listObjectsOfProperty(geometry, Geo.AS_GML_PROP);
            nodeIter = wktNodeIter.andThen(gmlNodeIter);
        }
        while (nodeIter.hasNext()) {
            Literal geometryLiteral = nodeIter.next().asLiteral();
            GeometryWrapper geometryWrapper = GeometryWrapper.extract(geometryLiteral);
            try {
                // Ensure all entries in the target SRS URI.
                GeometryWrapper transformedGeometryWrapper = geometryWrapper.convertSRS(srsURI);
                Envelope envelope = transformedGeometryWrapper.getEnvelope();
                SpatialIndexItem item = new SpatialIndexItem(envelope, feature);
                items.add(item);
            } catch (FactoryException | MismatchedDimensionException | TransformException ex) {
                throw new SpatialIndexException("Transformation Exception: " + geometryLiteral + ". " + ex.getMessage());
            }
        }
    }
    return items;
}
Also used : NodeIterator(org.apache.jena.rdf.model.NodeIterator) StmtIterator(org.apache.jena.rdf.model.StmtIterator) FactoryException(org.opengis.util.FactoryException) Statement(org.apache.jena.rdf.model.Statement) ArrayList(java.util.ArrayList) Resource(org.apache.jena.rdf.model.Resource) GeometryWrapper(org.apache.jena.geosparql.implementation.GeometryWrapper) TransformException(org.opengis.referencing.operation.TransformException) Envelope(org.locationtech.jts.geom.Envelope) MismatchedDimensionException(org.opengis.geometry.MismatchedDimensionException) Literal(org.apache.jena.rdf.model.Literal) RDFNode(org.apache.jena.rdf.model.RDFNode)

Example 12 with MismatchedDimensionException

use of org.opengis.geometry.MismatchedDimensionException in project jena by apache.

the class SearchEnvelope method build.

public static SearchEnvelope build(GeometryWrapper geometryWrapper, SRSInfo srsInfo, double radius, String unitsURI) {
    try {
        // Get the envelope of the target GeometryWrapper and convert to SRS URI, in case it is a complex polygon.
        GeometryWrapper envelopeGeometryWrapper = geometryWrapper.envelope();
        // Convert to SRS URI.
        GeometryWrapper srsGeometryWrapper = envelopeGeometryWrapper.convertSRS(srsInfo.getSrsURI());
        Envelope envelope = srsGeometryWrapper.getEnvelope();
        // Expand the envelope by the radius distance in all directions,
        // i.e. a bigger box rather than circle. More precise checks made later.
        SearchEnvelope searchEnvelope;
        if (srsInfo.isGeographic()) {
            double radiusMetres = UnitsOfMeasure.convertToMetres(radius, unitsURI, srsGeometryWrapper.getLatitude());
            searchEnvelope = expandGeographicEnvelope(envelope, radiusMetres, srsInfo);
        } else {
            String targetUnitsURI = srsInfo.getUnitsOfMeasure().getUnitURI();
            double targetRadius = UnitsOfMeasure.conversion(radius, unitsURI, targetUnitsURI);
            searchEnvelope = expandEnvelope(envelope, targetRadius, srsInfo);
        }
        return searchEnvelope;
    } catch (FactoryException | MismatchedDimensionException | TransformException ex) {
        throw new ExprEvalException(ex.getMessage() + ": " + geometryWrapper.asLiteral(), ex);
    }
}
Also used : FactoryException(org.opengis.util.FactoryException) GeometryWrapper(org.apache.jena.geosparql.implementation.GeometryWrapper) TransformException(org.opengis.referencing.operation.TransformException) Envelope(org.locationtech.jts.geom.Envelope) MismatchedDimensionException(org.opengis.geometry.MismatchedDimensionException) ExprEvalException(org.apache.jena.sparql.expr.ExprEvalException)

Example 13 with MismatchedDimensionException

use of org.opengis.geometry.MismatchedDimensionException in project ddf by codice.

the class GeospatialUtil method createBufferedCircleFromPoint.

private static Geometry createBufferedCircleFromPoint(Measure distance, CoordinateReferenceSystem origCRS, Geometry point) {
    Geometry pointGeo = point;
    Unit<?> unit = distance.getUnit();
    UnitConverter unitConverter = null;
    if (!(origCRS instanceof ProjectedCRS)) {
        double x = point.getCoordinate().x;
        double y = point.getCoordinate().y;
        // CRS code for UTM
        String crsCode = "AUTO:42001," + x + "," + y;
        try {
            CoordinateReferenceSystem utmCrs = CRS.decode(crsCode);
            MathTransform toTransform = CRS.findMathTransform(DefaultGeographicCRS.WGS84, utmCrs);
            MathTransform fromTransform = CRS.findMathTransform(utmCrs, DefaultGeographicCRS.WGS84);
            pointGeo = JTS.transform(point, toTransform);
            return JTS.transform(pointGeo.buffer(distance.doubleValue()), fromTransform);
        } catch (MismatchedDimensionException | TransformException | FactoryException e) {
            LOGGER.debug("Unable to create buffered circle from point.", e);
        }
    } else {
        try {
            unitConverter = unit.getConverterToAny(origCRS.getCoordinateSystem().getAxis(0).getUnit());
        } catch (IncommensurableException e) {
            LOGGER.debug("Unable to create unit converter.", e);
        }
    }
    if (unitConverter != null) {
        return pointGeo.buffer(unitConverter.convert(distance.doubleValue()));
    }
    return pointGeo.buffer(distance.doubleValue());
}
Also used : MathTransform(org.opengis.referencing.operation.MathTransform) FactoryException(org.opengis.referencing.FactoryException) IncommensurableException(javax.measure.IncommensurableException) TransformException(org.opengis.referencing.operation.TransformException) MismatchedDimensionException(org.opengis.geometry.MismatchedDimensionException) Geometry(org.locationtech.jts.geom.Geometry) ProjectedCRS(org.opengis.referencing.crs.ProjectedCRS) UnitConverter(javax.measure.UnitConverter) CoordinateReferenceSystem(org.opengis.referencing.crs.CoordinateReferenceSystem)

Example 14 with MismatchedDimensionException

use of org.opengis.geometry.MismatchedDimensionException in project hale by halestudio.

the class StreamGmlWriter method partitionByExtent.

private void partitionByExtent(ProgressIndicator progress, IOReporter reporter) throws IOException {
    int maxNodes = getParameter(PARAM_PARTITION_BY_EXTENT_MAX_NODES).as(Integer.class, 1000);
    String mode = getParameter(PARAM_PARTITION_BY_EXTENT_MODE).as(String.class, PARTITION_BY_EXTENT_MODE_DATASET);
    final SubtaskProgressIndicator qtProgress = new SubtaskProgressIndicator(progress) {

        @Override
        protected String getCombinedTaskName(String taskName, String subtaskName) {
            return taskName + " (" + subtaskName + ")";
        }
    };
    // Map for instances that either contain no or multiple geometries
    Map<String, InstanceReference> unhandledInstances = new HashMap<>();
    QuadtreeBuilder<Point, InstanceReference> builder = new QuadtreeBuilder<>();
    try (ResourceIterator<Instance> it = getInstances().iterator()) {
        qtProgress.begin("Collecting geometries", getInstances().size());
        final XMLInspector gadget = new XMLInspector();
        int i = 0;
        while (it.hasNext()) {
            Instance inst = it.next();
            InstanceReference instRef = getInstances().getReference(inst);
            InstanceTraverser traverser = new DepthFirstInstanceTraverser();
            GeometryFinder finder = new GeometryFinder(getTargetCRS());
            traverser.traverse(inst, finder);
            List<GeometryProperty<?>> geoms = finder.getGeometries();
            if (geoms.isEmpty() || geoms.size() > 1) {
                unhandledInstances.put(gadget.getIdentity(inst), instRef);
            } else {
                GeometryProperty<?> geomProperty = geoms.get(0);
                Geometry geom = geomProperty.getGeometry();
                Point centroid;
                switch(mode) {
                    case PARTITION_BY_EXTENT_MODE_WORLD:
                        CoordinateReferenceSystem sourceCrs = geomProperty.getCRSDefinition().getCRS();
                        CodeDefinition wgs84 = new CodeDefinition("EPSG:4326");
                        try {
                            MathTransform toWgs84 = CRS.findMathTransform(sourceCrs, wgs84.getCRS());
                            Geometry geomWgs84 = JTS.transform(geom, toWgs84);
                            centroid = geomWgs84.getCentroid();
                        } catch (FactoryException | MismatchedDimensionException | TransformException e) {
                            log.error("Unable to transform geometry to WGS 84", e);
                            throw new IllegalStateException(e.getMessage(), e);
                        }
                        break;
                    case PARTITION_BY_EXTENT_MODE_DATASET:
                    // fall through to default
                    default:
                        centroid = geom.getCentroid();
                }
                builder.add(centroid, new IdentifiableInstanceReference(instRef, gadget.getIdentity(inst)));
            }
            qtProgress.advance(1);
            if (++i % 100 == 0) {
                qtProgress.setCurrentTask(MessageFormat.format("{0} instances processed", i));
            }
        }
        qtProgress.setCurrentTask("Building quadtree");
        FixedBoundaryQuadtree<InstanceReference> qt;
        switch(mode) {
            case PARTITION_BY_EXTENT_MODE_DATASET:
                qt = builder.build(maxNodes);
                break;
            case PARTITION_BY_EXTENT_MODE_WORLD:
                Envelope world = new Envelope(-180, 180, -90, 90);
                qt = builder.build(maxNodes, world);
                break;
            default:
                log.error(MessageFormat.format("Unrecognized extent partitioning mode \"{0}\", using dataset boundaries", mode));
                qt = builder.build(maxNodes);
        }
        qtProgress.setCurrentTask("Performing spatial partitioning");
        final Map<String, String> idToKeyMapping = new HashMap<>();
        final Map<String, Collection<InstanceReference>> keyToRefsMapping = new HashMap<>();
        // Instances without geometry or with multiple geometries
        keyToRefsMapping.put(ExtentPartsHandler.KEY_NO_GEOMETRY, unhandledInstances.values());
        unhandledInstances.keySet().stream().forEach(id -> idToKeyMapping.put(id, ExtentPartsHandler.KEY_NO_GEOMETRY));
        buildMappings(qt, idToKeyMapping, keyToRefsMapping);
        // Partition source instances based on quadtree tiles
        Iterator<InstanceCollection> collIt = new Iterator<InstanceCollection>() {

            private final Queue<String> keySet = new LinkedList<>(keyToRefsMapping.keySet());

            @Override
            public boolean hasNext() {
                return !keySet.isEmpty();
            }

            @Override
            public InstanceCollection next() {
                String key = keySet.poll();
                Collection<InstanceReference> refs = keyToRefsMapping.get(key);
                InstanceCollection instColl = new DefaultInstanceCollection(refs.stream().map(ref -> getInstances().getInstance(IdentifiableInstanceReference.getRootReference(ref))).collect(Collectors.toList()));
                return new ExtentPartsHandler.TreeKeyDecorator(instColl, key);
            }
        };
        final Map<String, URI> keyToTargetMapping = new HashMap<>();
        keyToRefsMapping.keySet().stream().forEach(k -> keyToTargetMapping.put(k, new File(ExtentPartsHandler.getTargetFilename(k, getTarget().getLocation())).toURI()));
        final ExtentPartsHandler handler = new ExtentPartsHandler(keyToTargetMapping, idToKeyMapping);
        qtProgress.end();
        try {
            writeParts(collIt, handler, progress, reporter);
        } catch (XMLStreamException e) {
            throw new IOException(e.getMessage(), e);
        }
    }
}
Also used : MathTransform(org.opengis.referencing.operation.MathTransform) HashMap(java.util.HashMap) Instance(eu.esdihumboldt.hale.common.instance.model.Instance) GeometryFinder(eu.esdihumboldt.hale.common.instance.geometry.GeometryFinder) FactoryException(org.opengis.referencing.FactoryException) IdentifiableInstanceReference(eu.esdihumboldt.hale.common.instance.model.IdentifiableInstanceReference) Envelope(org.locationtech.jts.geom.Envelope) DefaultInstanceCollection(eu.esdihumboldt.hale.common.instance.model.impl.DefaultInstanceCollection) MismatchedDimensionException(org.opengis.geometry.MismatchedDimensionException) URI(java.net.URI) CodeDefinition(eu.esdihumboldt.hale.common.instance.geometry.impl.CodeDefinition) ResourceIterator(eu.esdihumboldt.hale.common.instance.model.ResourceIterator) Iterator(java.util.Iterator) CoordinateReferenceSystem(org.opengis.referencing.crs.CoordinateReferenceSystem) Queue(java.util.Queue) DepthFirstInstanceTraverser(eu.esdihumboldt.hale.common.instance.helper.DepthFirstInstanceTraverser) InstanceTraverser(eu.esdihumboldt.hale.common.instance.helper.InstanceTraverser) GeometryProperty(eu.esdihumboldt.hale.common.schema.geometry.GeometryProperty) PerTypeInstanceCollection(eu.esdihumboldt.hale.common.instance.model.ext.impl.PerTypeInstanceCollection) DefaultInstanceCollection(eu.esdihumboldt.hale.common.instance.model.impl.DefaultInstanceCollection) InstanceCollection(eu.esdihumboldt.hale.common.instance.model.InstanceCollection) TransformException(org.opengis.referencing.operation.TransformException) SubtaskProgressIndicator(eu.esdihumboldt.hale.common.core.io.impl.SubtaskProgressIndicator) QuadtreeBuilder(eu.esdihumboldt.util.geometry.quadtree.QuadtreeBuilder) Point(org.locationtech.jts.geom.Point) IOException(java.io.IOException) Point(org.locationtech.jts.geom.Point) DepthFirstInstanceTraverser(eu.esdihumboldt.hale.common.instance.helper.DepthFirstInstanceTraverser) Geometry(org.locationtech.jts.geom.Geometry) XMLStreamException(javax.xml.stream.XMLStreamException) InstanceReference(eu.esdihumboldt.hale.common.instance.model.InstanceReference) IdentifiableInstanceReference(eu.esdihumboldt.hale.common.instance.model.IdentifiableInstanceReference) XMLInspector(eu.esdihumboldt.hale.common.instance.graph.reference.impl.XMLInspector) PerTypeInstanceCollection(eu.esdihumboldt.hale.common.instance.model.ext.impl.PerTypeInstanceCollection) DefaultInstanceCollection(eu.esdihumboldt.hale.common.instance.model.impl.DefaultInstanceCollection) Collection(java.util.Collection) InstanceCollection(eu.esdihumboldt.hale.common.instance.model.InstanceCollection) File(java.io.File)

Example 15 with MismatchedDimensionException

use of org.opengis.geometry.MismatchedDimensionException in project sis by apache.

the class Line method fit.

/**
 * Given a sequence of points, fits them to a straight line
 * <var>y</var> = <var>slope</var>⋅<var>x</var> + <var>y₀</var> in a least-squares senses.
 * Points shall be two dimensional with ordinate values in the (<var>x</var>,<var>y</var>) order.
 * This method assumes that the <var>x</var> values are precise and all uncertainty is in <var>y</var>.
 * {@link Double#NaN} ordinate values are ignored.
 *
 * @param  points  the two-dimensional points.
 * @return estimation of the correlation coefficient. The closer this coefficient is to +1 or -1, the better the fit.
 * @throws MismatchedDimensionException if a point is not two-dimensional.
 */
public double fit(final Iterable<? extends DirectPosition> points) {
    ArgumentChecks.ensureNonNull("points", points);
    int i = 0, n = 0;
    final DoubleDouble mean_x = new DoubleDouble();
    final DoubleDouble mean_y = new DoubleDouble();
    for (final DirectPosition p : points) {
        final int dimension = p.getDimension();
        if (dimension != DIMENSION) {
            throw new MismatchedDimensionException(Errors.format(Errors.Keys.MismatchedDimension_3, "points[" + i + ']', DIMENSION, dimension));
        }
        i++;
        final double x, y;
        if (// Test first the dimension which is most likely to contain NaN.
        !isNaN(y = p.getOrdinate(1)) && !isNaN(x = p.getOrdinate(0))) {
            mean_x.add(x);
            mean_y.add(y);
            n++;
        }
    }
    mean_x.divide(n, 0);
    mean_y.divide(n, 0);
    /*
         * We have to solve two equations with two unknowns:
         *
         *   1)    mean(y)   = m⋅mean(x) + y₀
         *   2)    mean(x⋅y) = m⋅mean(x²) + y₀⋅mean(x)
         *
         * Those formulas lead to a quadratic equation. However, the formulas become very simples
         * if we set 'mean(x) = 0'. We can achieve this result by computing instead of (2):
         *
         *   2b)   mean(Δx⋅y) = m⋅mean(Δx²)
         *
         * where dx = x-mean(x). In this case mean(Δx) == 0.
         */
    final DoubleDouble a = new DoubleDouble();
    final DoubleDouble b = new DoubleDouble();
    final DoubleDouble mean_x2 = new DoubleDouble();
    final DoubleDouble mean_y2 = new DoubleDouble();
    final DoubleDouble mean_xy = new DoubleDouble();
    for (final DirectPosition p : points) {
        final double y;
        if (// Test first the dimension which is most likely to contain NaN.
        !isNaN(y = p.getOrdinate(1)) && !isNaN(a.value = p.getOrdinate(0))) {
            a.error = 0;
            // Δx = x - mean_x
            a.subtract(mean_x);
            b.setFrom(a);
            // (Δx)² = Δx * Δx
            b.multiply(a);
            // mean_x² += (Δx)²
            mean_x2.add(b);
            // xy = Δx * y
            a.multiply(y);
            // mean_xy += xy
            mean_xy.add(a);
            // y² = y * y
            a.setToProduct(y, y);
            // mean_y² += y²
            mean_y2.add(a);
        }
    }
    mean_x2.divide(n, 0);
    mean_y2.divide(n, 0);
    mean_xy.divide(n, 0);
    /*
         * Assuming that 'mean(x) == 0', then the correlation
         * coefficient can be approximate by:
         *
         * R = mean(xy) / sqrt( mean(x²) * (mean(y²) - mean(y)²) )
         */
    a.setFrom(mean_xy);
    // slope = mean_xy / mean_x²
    a.divide(mean_x2);
    b.setFrom(mean_x);
    b.multiply(a);
    b.negate();
    // y₀ = mean_y - mean_x * slope
    b.add(mean_y);
    setEquation(a, b);
    /*
         * Compute the correlation coefficient:
         * mean_xy / sqrt(mean_x2 * (mean_y2 - mean_y * mean_y))
         */
    a.setFrom(mean_y);
    a.multiply(mean_y);
    a.negate();
    a.add(mean_y2);
    a.multiply(mean_x2);
    a.sqrt();
    a.inverseDivide(mean_xy);
    return a.value;
}
Also used : DirectPosition(org.opengis.geometry.DirectPosition) MismatchedDimensionException(org.opengis.geometry.MismatchedDimensionException) DoubleDouble(org.apache.sis.internal.util.DoubleDouble)

Aggregations

MismatchedDimensionException (org.opengis.geometry.MismatchedDimensionException)29 TransformException (org.opengis.referencing.operation.TransformException)26 GeometryWrapper (org.apache.jena.geosparql.implementation.GeometryWrapper)19 FactoryException (org.opengis.util.FactoryException)19 ExprEvalException (org.apache.jena.sparql.expr.ExprEvalException)15 DatatypeFormatException (org.apache.jena.datatypes.DatatypeFormatException)10 Envelope (org.locationtech.jts.geom.Envelope)8 MathTransform (org.opengis.referencing.operation.MathTransform)6 FactoryException (org.opengis.referencing.FactoryException)4 ArrayList (java.util.ArrayList)3 Literal (org.apache.jena.rdf.model.Literal)3 Resource (org.apache.jena.rdf.model.Resource)3 Geometry (org.locationtech.jts.geom.Geometry)3 DirectPosition (org.opengis.geometry.DirectPosition)3 CoordinateReferenceSystem (org.opengis.referencing.crs.CoordinateReferenceSystem)3 GeometryFinder (eu.esdihumboldt.hale.common.instance.geometry.GeometryFinder)2 CodeDefinition (eu.esdihumboldt.hale.common.instance.geometry.impl.CodeDefinition)2 DepthFirstInstanceTraverser (eu.esdihumboldt.hale.common.instance.helper.DepthFirstInstanceTraverser)2 InstanceTraverser (eu.esdihumboldt.hale.common.instance.helper.InstanceTraverser)2 GeometryProperty (eu.esdihumboldt.hale.common.schema.geometry.GeometryProperty)2