Search in sources :

Example 1 with NoninvertibleTransformException

use of org.opengis.referencing.operation.NoninvertibleTransformException in project sis by apache.

the class Shapes2D method transform.

/**
 * Transforms a rectangular envelope using the given coordinate operation.
 * The transformation is only approximative: the returned envelope may be bigger
 * than the smallest possible bounding box, but should not be smaller in most cases.
 *
 * <p>This method can handle the case where the rectangle contains the North or South pole,
 * or when it cross the ±180° longitude.</p>
 *
 * @param  operation    the operation to use. Source and target dimension must be 2.
 * @param  envelope     the rectangle to transform (may be {@code null}).
 * @param  destination  the destination rectangle (may be {@code envelope}).
 *         If {@code null}, a new rectangle will be created and returned.
 * @return {@code destination}, or a new rectangle if {@code destination} was non-null and {@code envelope} was null.
 * @throws TransformException if a transform failed.
 *
 * @see #transform(MathTransform2D, Rectangle2D, Rectangle2D)
 * @see Envelopes#transform(CoordinateOperation, Envelope)
 */
@SuppressWarnings("null")
public static Rectangle2D transform(final CoordinateOperation operation, final Rectangle2D envelope, Rectangle2D destination) throws TransformException {
    ArgumentChecks.ensureNonNull("operation", operation);
    if (envelope == null) {
        return null;
    }
    final MathTransform transform = operation.getMathTransform();
    if (!(transform instanceof MathTransform2D)) {
        throw new MismatchedDimensionException(Errors.format(Errors.Keys.IllegalPropertyValueClass_3, "transform", MathTransform2D.class, MathTransform.class));
    }
    MathTransform2D mt = (MathTransform2D) transform;
    final double[] center = new double[2];
    destination = transform(mt, envelope, destination, center);
    /*
         * If the source envelope crosses the expected range of valid coordinates, also projects
         * the range bounds as a safety. See the comments in transform(Envelope, ...).
         */
    final CoordinateReferenceSystem sourceCRS = operation.getSourceCRS();
    if (sourceCRS != null) {
        final CoordinateSystem cs = sourceCRS.getCoordinateSystem();
        if (cs != null && cs.getDimension() == 2) {
            // Paranoiac check.
            CoordinateSystemAxis axis = cs.getAxis(0);
            double min = envelope.getMinX();
            double max = envelope.getMaxX();
            Point2D.Double pt = null;
            for (int i = 0; i < 4; i++) {
                if (i == 2) {
                    axis = cs.getAxis(1);
                    min = envelope.getMinY();
                    max = envelope.getMaxY();
                }
                final double v = (i & 1) == 0 ? axis.getMinimumValue() : axis.getMaximumValue();
                if (!(v > min && v < max)) {
                    continue;
                }
                if (pt == null) {
                    pt = new Point2D.Double();
                }
                if ((i & 2) == 0) {
                    pt.x = v;
                    pt.y = envelope.getCenterY();
                } else {
                    pt.x = envelope.getCenterX();
                    pt.y = v;
                }
                destination.add(mt.transform(pt, pt));
            }
        }
    }
    /*
         * Now take the target CRS in account.
         */
    final CoordinateReferenceSystem targetCRS = operation.getTargetCRS();
    if (targetCRS == null) {
        return destination;
    }
    final CoordinateSystem targetCS = targetCRS.getCoordinateSystem();
    if (targetCS == null || targetCS.getDimension() != 2) {
        // It should be an error, but we keep this method tolerant.
        return destination;
    }
    /*
         * Checks for singularity points. See the Envelopes.transform(CoordinateOperation, Envelope)
         * method for comments about the algorithm. The code below is the same algorithm adapted for
         * the 2D case and the related objects (Point2D, Rectangle2D, etc.).
         *
         * The 'border' variable in the loop below is used in order to compress 2 dimensions
         * and 2 extremums in a single loop, in this order: (xmin, xmax, ymin, ymax).
         */
    TransformException warning = null;
    Point2D sourcePt = null;
    Point2D targetPt = null;
    // A bitmask for each (dimension, extremum) pairs.
    int includedBoundsValue = 0;
    for (int border = 0; border < 4; border++) {
        // 2 dimensions and 2 extremums compacted in a flag.
        // The dimension index being examined.
        final int dimension = border >>> 1;
        final CoordinateSystemAxis axis = targetCS.getAxis(dimension);
        if (axis == null) {
            // Should never be null, but check as a paranoiac safety.
            continue;
        }
        final double extremum = (border & 1) == 0 ? axis.getMinimumValue() : axis.getMaximumValue();
        if (Double.isInfinite(extremum) || Double.isNaN(extremum)) {
            continue;
        }
        if (targetPt == null) {
            try {
                mt = mt.inverse();
            } catch (NoninvertibleTransformException exception) {
                Envelopes.recoverableException(Shapes2D.class, exception);
                return destination;
            }
            targetPt = new Point2D.Double();
        }
        switch(dimension) {
            case 0:
                targetPt.setLocation(extremum, center[1]);
                break;
            case 1:
                targetPt.setLocation(center[0], extremum);
                break;
            default:
                throw new AssertionError(border);
        }
        try {
            sourcePt = mt.transform(targetPt, sourcePt);
        } catch (TransformException exception) {
            if (warning == null) {
                warning = exception;
            } else {
                warning.addSuppressed(exception);
            }
            continue;
        }
        if (envelope.contains(sourcePt)) {
            destination.add(targetPt);
            includedBoundsValue |= (1 << border);
        }
    }
    /*
         * Iterate over all dimensions of type "WRAPAROUND" for which minimal or maximal axis
         * values have not yet been included in the envelope. We could inline this check inside
         * the above loop, but we don't in order to have a chance to exclude the dimensions for
         * which the point have already been added.
         *
         * See transform(CoordinateOperation, Envelope) for more comments about the algorithm.
         */
    if (includedBoundsValue != 0) {
        /*
             * Bits mask transformation:
             *   1) Swaps the two dimensions               (YyXx  →  XxYy)
             *   2) Insert a space between each bits       (XxYy  →  X.x.Y.y.)
             *   3) Fill the space with duplicated values  (X.x.Y.y.  →  XXxxYYyy)
             *
             * In terms of bit positions 1,2,4,8 (not bit values), we have:
             *
             *   8421  →  22881144
             *   i.e. (ymax, ymin, xmax, xmin)  →  (xmax², ymax², xmin², ymin²)
             *
             * Now look at the last part: (xmin², ymin²). The next step is to perform a bitwise
             * AND operation in order to have only both of the following conditions:
             *
             *   Borders not yet added to the envelope: ~(ymax, ymin, xmax, xmin)
             *   Borders in which a singularity exists:  (xmin, xmin, ymin, ymin)
             *
             * The same operation is repeated on the next 4 bits for (xmax, xmax, ymax, ymax).
             */
        int toTest = ((includedBoundsValue & 1) << 3) | ((includedBoundsValue & 4) >>> 1) | ((includedBoundsValue & 2) << 6) | ((includedBoundsValue & 8) << 2);
        // Duplicate the bit values.
        toTest |= (toTest >>> 1);
        toTest &= ~(includedBoundsValue | (includedBoundsValue << 4));
        /*
             * Forget any axes that are not of kind "WRAPAROUND". Then get the final
             * bit pattern indicating which points to test. Iterate over that bits.
             */
        if ((toTest & 0x33333333) != 0 && !CoordinateOperations.isWrapAround(targetCS.getAxis(0)))
            toTest &= 0xCCCCCCCC;
        if ((toTest & 0xCCCCCCCC) != 0 && !CoordinateOperations.isWrapAround(targetCS.getAxis(1)))
            toTest &= 0x33333333;
        while (toTest != 0) {
            final int border = Integer.numberOfTrailingZeros(toTest);
            final int bitMask = 1 << border;
            // Clear now the bit, for the next iteration.
            toTest &= ~bitMask;
            final int dimensionToAdd = (border >>> 1) & 1;
            final CoordinateSystemAxis toAdd = targetCS.getAxis(dimensionToAdd);
            final CoordinateSystemAxis added = targetCS.getAxis(dimensionToAdd ^ 1);
            double x = (border & 1) == 0 ? toAdd.getMinimumValue() : toAdd.getMaximumValue();
            double y = (border & 4) == 0 ? added.getMinimumValue() : added.getMaximumValue();
            if (dimensionToAdd != 0) {
                final double t = x;
                x = y;
                y = t;
            }
            targetPt.setLocation(x, y);
            try {
                sourcePt = mt.transform(targetPt, sourcePt);
            } catch (TransformException exception) {
                if (warning == null) {
                    warning = exception;
                } else {
                    warning.addSuppressed(exception);
                }
                continue;
            }
            if (envelope.contains(sourcePt)) {
                destination.add(targetPt);
            }
        }
    }
    /*
         * At this point we finished envelope transformation. Verify if some ordinates need to be "wrapped around"
         * as a result of the coordinate operation.   This is usually the longitude axis where the source CRS uses
         * the [-180 … +180]° range and the target CRS uses the [0 … 360]° range, or the converse. In such case we
         * set the rectangle to the full range (we do not use the mechanism documented in Envelope2D) because most
         * Rectangle2D implementations do not support spanning the anti-meridian. This results in larger rectangle
         * than what would be possible with GeneralEnvelope or Envelope2D, but we try to limit the situation where
         * this expansion is applied.
         */
    final Set<Integer> wrapAroundChanges;
    if (operation instanceof AbstractCoordinateOperation) {
        wrapAroundChanges = ((AbstractCoordinateOperation) operation).getWrapAroundChanges();
    } else {
        wrapAroundChanges = CoordinateOperations.wrapAroundChanges(sourceCRS, targetCS);
    }
    for (int dim : wrapAroundChanges) {
        // Empty in the vast majority of cases.
        final CoordinateSystemAxis axis = targetCS.getAxis(dim);
        final double minimum = axis.getMinimumValue();
        final double maximum = axis.getMaximumValue();
        final double o1, o2;
        if (dim == 0) {
            o1 = destination.getMinX();
            o2 = destination.getMaxX();
        } else {
            o1 = destination.getMinY();
            o2 = destination.getMaxY();
        }
        if (o1 < minimum || o2 > maximum) {
            final double span = maximum - minimum;
            if (dim == 0) {
                destination.setRect(minimum, destination.getY(), span, destination.getHeight());
            } else {
                destination.setRect(destination.getX(), minimum, destination.getWidth(), span);
            }
        }
    }
    if (warning != null) {
        Envelopes.recoverableException(Shapes2D.class, warning);
    }
    return destination;
}
Also used : MathTransform(org.opengis.referencing.operation.MathTransform) CoordinateSystem(org.opengis.referencing.cs.CoordinateSystem) NoninvertibleTransformException(org.opengis.referencing.operation.NoninvertibleTransformException) TransformException(org.opengis.referencing.operation.TransformException) CoordinateSystemAxis(org.opengis.referencing.cs.CoordinateSystemAxis) MismatchedDimensionException(org.opengis.geometry.MismatchedDimensionException) NoninvertibleTransformException(org.opengis.referencing.operation.NoninvertibleTransformException) Point2D(java.awt.geom.Point2D) MathTransform2D(org.opengis.referencing.operation.MathTransform2D) CoordinateReferenceSystem(org.opengis.referencing.crs.CoordinateReferenceSystem) AbstractCoordinateOperation(org.apache.sis.referencing.operation.AbstractCoordinateOperation)

Example 2 with NoninvertibleTransformException

use of org.opengis.referencing.operation.NoninvertibleTransformException in project sis by apache.

the class FranceGeocentricInterpolation method createMathTransform.

/**
 * Creates a transform from the specified group of parameter values.
 * This method creates the transform from <em>target</em> to <em>source</em>
 * (which is the direction that use the interpolation grid directly without iteration),
 * then inverts the transform.
 *
 * @param  factory  the factory to use if this constructor needs to create other math transforms.
 * @param  values   the group of parameter values.
 * @return the created math transform.
 * @throws ParameterNotFoundException if a required parameter was not found.
 * @throws FactoryException if an error occurred while loading the grid.
 */
@Override
public MathTransform createMathTransform(final MathTransformFactory factory, final ParameterValueGroup values) throws ParameterNotFoundException, FactoryException {
    boolean withHeights = false;
    final Parameters pg = Parameters.castOrWrap(values);
    final Integer dim = pg.getValue(Molodensky.DIMENSION);
    if (dim != null)
        switch(dim) {
            case 2:
                break;
            case 3:
                withHeights = true;
                break;
            default:
                throw new InvalidParameterValueException(Errors.format(Errors.Keys.IllegalArgumentValue_2, "dim", dim), "dim", dim);
        }
    final Path file = pg.getMandatoryValue(FILE);
    final DatumShiftGridFile<Angle, Length> grid = getOrLoad(file, isRecognized(file) ? new double[] { TX, TY, TZ } : null, PRECISION);
    MathTransform tr = createGeodeticTransformation(factory, createEllipsoid(pg, Molodensky.TGT_SEMI_MAJOR, Molodensky.TGT_SEMI_MINOR, // GRS 1980 ellipsoid
    CommonCRS.ETRS89.ellipsoid()), createEllipsoid(pg, Molodensky.SRC_SEMI_MAJOR, Molodensky.SRC_SEMI_MINOR, // Clarke 1880 (IGN) ellipsoid
    null), withHeights, grid);
    try {
        tr = tr.inverse();
    } catch (NoninvertibleTransformException e) {
        // Should never happen.
        throw new FactoryException(e);
    }
    return tr;
}
Also used : Path(java.nio.file.Path) NoninvertibleTransformException(org.opengis.referencing.operation.NoninvertibleTransformException) Parameters(org.apache.sis.parameter.Parameters) Angle(javax.measure.quantity.Angle) MathTransform(org.opengis.referencing.operation.MathTransform) InvalidParameterValueException(org.opengis.parameter.InvalidParameterValueException) Length(javax.measure.quantity.Length) FactoryException(org.opengis.util.FactoryException)

Example 3 with NoninvertibleTransformException

use of org.opengis.referencing.operation.NoninvertibleTransformException in project sis by apache.

the class FranceGeocentricInterpolation method getOrLoad.

/**
 * Returns the grid of the given name. This method returns the cached instance if it still exists,
 * or load the grid otherwise.
 *
 * @param  file      name of the datum shift grid file to load.
 * @param  averages  an "average" value for the offset in each dimension, or {@code null} if unknown.
 * @param  scale     the factor by which to multiply each compressed value before to add to the average value.
 */
@SuppressWarnings("null")
static DatumShiftGridFile<Angle, Length> getOrLoad(final Path file, final double[] averages, final double scale) throws FactoryException {
    final Path resolved = DataDirectory.DATUM_CHANGES.resolve(file).toAbsolutePath();
    DatumShiftGridFile<?, ?> grid = DatumShiftGridFile.CACHE.peek(resolved);
    if (grid == null) {
        final Cache.Handler<DatumShiftGridFile<?, ?>> handler = DatumShiftGridFile.CACHE.lock(resolved);
        try {
            grid = handler.peek();
            if (grid == null) {
                try (BufferedReader in = Files.newBufferedReader(resolved)) {
                    DatumShiftGridLoader.log(FranceGeocentricInterpolation.class, file);
                    final DatumShiftGridFile.Float<Angle, Length> g = load(in, file);
                    grid = DatumShiftGridCompressed.compress(g, averages, scale);
                } catch (IOException | NoninvertibleTransformException | RuntimeException e) {
                    // NumberFormatException, ArithmeticException, NoSuchElementException, possibly other.
                    throw DatumShiftGridLoader.canNotLoad(HEADER, file, e);
                }
                grid = grid.useSharedData();
            }
        } finally {
            handler.putAndUnlock(grid);
        }
    }
    return grid.castTo(Angle.class, Length.class);
}
Also used : Path(java.nio.file.Path) IOException(java.io.IOException) NoninvertibleTransformException(org.opengis.referencing.operation.NoninvertibleTransformException) Angle(javax.measure.quantity.Angle) Length(javax.measure.quantity.Length) BufferedReader(java.io.BufferedReader) Cache(org.apache.sis.util.collection.Cache)

Example 4 with NoninvertibleTransformException

use of org.opengis.referencing.operation.NoninvertibleTransformException in project sis by apache.

the class GeocentricAffineBetweenGeographic method createMathTransform.

/**
 * Creates a math transform from the specified group of parameter values.
 * This method wraps the affine operation into Geographic/Geocentric conversions.
 *
 * @param  factory  the factory to use for creating concatenated transforms.
 * @param  values   the group of parameter values.
 * @return the created math transform.
 * @throws FactoryException if a transform can not be created.
 */
@Override
public MathTransform createMathTransform(final MathTransformFactory factory, final ParameterValueGroup values) throws FactoryException {
    final Parameters pv = Parameters.castOrWrap(values);
    final MathTransform affine = super.createMathTransform(factory, pv);
    /*
         * Create a "Geographic to Geocentric" conversion with ellipsoid axis length units converted to metres
         * (the unit implied by SRC_SEMI_MAJOR) because it is the unit of Bursa-Wolf parameters that we created above.
         */
    MathTransform toGeocentric = EllipsoidToCentricTransform.createGeodeticConversion(factory, pv.doubleValue(SRC_SEMI_MAJOR), pv.doubleValue(SRC_SEMI_MINOR), Units.METRE, getSourceDimensions() >= 3, EllipsoidToCentricTransform.TargetType.CARTESIAN);
    /*
         * Create a "Geocentric to Geographic" conversion with ellipsoid axis length units converted to metres
         * because this is the unit of the Geocentric CRS used above.
         */
    MathTransform toGeographic = EllipsoidToCentricTransform.createGeodeticConversion(factory, pv.doubleValue(TGT_SEMI_MAJOR), pv.doubleValue(TGT_SEMI_MINOR), Units.METRE, getTargetDimensions() >= 3, EllipsoidToCentricTransform.TargetType.CARTESIAN);
    try {
        toGeographic = toGeographic.inverse();
    } catch (NoninvertibleTransformException e) {
        // Should never happen with SIS implementation.
        throw new FactoryException(e);
    }
    /*
         * The  Geocentric → Affine → Geographic  chain.
         */
    return factory.createConcatenatedTransform(toGeocentric, factory.createConcatenatedTransform(affine, toGeographic));
}
Also used : NoninvertibleTransformException(org.opengis.referencing.operation.NoninvertibleTransformException) Parameters(org.apache.sis.parameter.Parameters) MathTransform(org.opengis.referencing.operation.MathTransform) FactoryException(org.opengis.util.FactoryException)

Example 5 with NoninvertibleTransformException

use of org.opengis.referencing.operation.NoninvertibleTransformException in project sis by apache.

the class NTv2 method getOrLoad.

/**
 * Returns the grid of the given name. This method returns the cached instance if it still exists,
 * or load the grid otherwise.
 *
 * @param  file  name of the datum shift grid file to load.
 */
@SuppressWarnings("null")
static DatumShiftGridFile<Angle, Angle> getOrLoad(final Path file) throws FactoryException {
    final Path resolved = DataDirectory.DATUM_CHANGES.resolve(file).toAbsolutePath();
    DatumShiftGridFile<?, ?> grid = DatumShiftGridFile.CACHE.peek(resolved);
    if (grid == null) {
        final Cache.Handler<DatumShiftGridFile<?, ?>> handler = DatumShiftGridFile.CACHE.lock(resolved);
        try {
            grid = handler.peek();
            if (grid == null) {
                try (ReadableByteChannel in = Files.newByteChannel(resolved)) {
                    DatumShiftGridLoader.log(NTv2.class, file);
                    final Loader loader = new Loader(in, file);
                    grid = loader.readGrid();
                    loader.reportWarnings();
                } catch (IOException | NoninvertibleTransformException | RuntimeException e) {
                    throw DatumShiftGridLoader.canNotLoad("NTv2", file, e);
                }
                grid = grid.useSharedData();
            }
        } finally {
            handler.putAndUnlock(grid);
        }
    }
    return grid.castTo(Angle.class, Angle.class);
}
Also used : Path(java.nio.file.Path) NoninvertibleTransformException(org.opengis.referencing.operation.NoninvertibleTransformException) ReadableByteChannel(java.nio.channels.ReadableByteChannel) IOException(java.io.IOException) Cache(org.apache.sis.util.collection.Cache)

Aggregations

NoninvertibleTransformException (org.opengis.referencing.operation.NoninvertibleTransformException)13 MathTransform (org.opengis.referencing.operation.MathTransform)7 Path (java.nio.file.Path)4 FactoryException (org.opengis.util.FactoryException)4 IOException (java.io.IOException)3 Parameters (org.apache.sis.parameter.Parameters)3 Cache (org.apache.sis.util.collection.Cache)3 CoordinateSystem (org.opengis.referencing.cs.CoordinateSystem)3 ReadableByteChannel (java.nio.channels.ReadableByteChannel)2 Angle (javax.measure.quantity.Angle)2 Length (javax.measure.quantity.Length)2 AbstractCoordinateOperation (org.apache.sis.referencing.operation.AbstractCoordinateOperation)2 CoordinateReferenceSystem (org.opengis.referencing.crs.CoordinateReferenceSystem)2 CoordinateSystemAxis (org.opengis.referencing.cs.CoordinateSystemAxis)2 MathTransform1D (org.opengis.referencing.operation.MathTransform1D)2 TransformException (org.opengis.referencing.operation.TransformException)2 Point2D (java.awt.geom.Point2D)1 BufferedReader (java.io.BufferedReader)1 ByteBuffer (java.nio.ByteBuffer)1 FloatBuffer (java.nio.FloatBuffer)1