Examples with Matrix org.opengis.referencing.operation.Matrix used on opensource projects

Example 21 with Matrix

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

the class Shapes2D method transform.

/**
 * Implementation of {@link #transform(MathTransform2D, Rectangle2D, Rectangle2D)} with the
 * opportunity to save the projected center coordinate. This method sets {@code point} to
 * the center of the source envelope projected to the target CRS.
 */
@SuppressWarnings("fallthrough")
private static Rectangle2D transform(final MathTransform2D transform, final Rectangle2D envelope, Rectangle2D destination, final double[] point) throws TransformException {
    if (envelope == null) {
        return null;
    }
    double xmin = Double.POSITIVE_INFINITY;
    double ymin = Double.POSITIVE_INFINITY;
    double xmax = Double.NEGATIVE_INFINITY;
    double ymax = Double.NEGATIVE_INFINITY;
    /*
         * Notation (as if we were applying a map projection, but this is not necessarily the case):
         *   - (λ,φ) are ordinate values before projection.
         *   - (x,y) are ordinate values after projection.
         *   - D[00|01|10|11] are the ∂x/∂λ, ∂x/∂φ, ∂y/∂λ and ∂y/∂φ derivatives respectively.
         *   - Variables with indice 0 are for the very first point in iteration order.
         *   - Variables with indice 1 are for the values of the previous iteration.
         *   - Variables with indice 2 are for the current values in the iteration.
         *   - P1-P2 form a line segment to be checked for curvature.
         */
    double x0 = 0, y0 = 0, λ0 = 0, φ0 = 0;
    double x1 = 0, y1 = 0, λ1 = 0, φ1 = 0;
    Matrix D0 = null, D1 = null, D2 = null;
    // x2 and y2 defined inside the loop.
    boolean isDerivativeSupported = true;
    final CurveExtremum extremum = new CurveExtremum();
    for (int i = 0; i <= 8; i++) {
        /*
             * Iteration order (center must be last):
             *
             *   (6)────(5)────(4)
             *    |             |
             *   (7)    (8)    (3)
             *    |             |
             *   (0)────(1)────(2)
             */
        double λ2, φ2;
        switch(i) {
            case 0:
            case 6:
            case 7:
                λ2 = envelope.getMinX();
                break;
            case 1:
            case 5:
            case 8:
                λ2 = envelope.getCenterX();
                break;
            case 2:
            case 3:
            case 4:
                λ2 = envelope.getMaxX();
                break;
            default:
                throw new AssertionError(i);
        }
        switch(i) {
            case 0:
            case 1:
            case 2:
                φ2 = envelope.getMinY();
                break;
            case 3:
            case 7:
            case 8:
                φ2 = envelope.getCenterY();
                break;
            case 4:
            case 5:
            case 6:
                φ2 = envelope.getMaxY();
                break;
            default:
                throw new AssertionError(i);
        }
        point[0] = λ2;
        point[1] = φ2;
        try {
            D1 = D2;
            D2 = Envelopes.derivativeAndTransform(transform, point, point, 0, isDerivativeSupported && i != 8);
        } catch (TransformException e) {
            if (!isDerivativeSupported) {
                // Derivative were already disabled, so something went wrong.
                throw e;
            }
            isDerivativeSupported = false;
            D2 = null;
            point[0] = λ2;
            point[1] = φ2;
            transform.transform(point, 0, point, 0, 1);
            // Log only if the above call was successful.
            Envelopes.recoverableException(Shapes2D.class, e);
        }
        double x2 = point[0];
        double y2 = point[1];
        if (x2 < xmin)
            xmin = x2;
        if (x2 > xmax)
            xmax = x2;
        if (y2 < ymin)
            ymin = y2;
        if (y2 > ymax)
            ymax = y2;
        switch(i) {
            case 0:
                {
                    // Remember the first point.
                    λ0 = λ2;
                    x0 = x2;
                    φ0 = φ2;
                    y0 = y2;
                    D0 = D2;
                    break;
                }
            case 8:
                {
                    // Close the iteration with the first point.
                    // Discard P2 because it is the rectangle center.
                    λ2 = λ0;
                    // Discard P2 because it is the rectangle center.
                    x2 = x0;
                    φ2 = φ0;
                    y2 = y0;
                    D2 = D0;
                    break;
                }
        }
        /*
             * At this point, we expanded the rectangle using the projected points. Now try
             * to use the information provided by derivatives at those points, if available.
             * For the following block, notation is:
             *
             *   - s  are ordinate values in the source space (λ or φ)
             *   - t  are ordinate values in the target space (x or y)
             *
             * They are not necessarily in the same dimension. For example would could have
             * s=λ while t=y. This is typically the case when inspecting the top or bottom
             * line segment of the rectangle.
             *
             * The same technic is also applied in the transform(MathTransform, Envelope) method.
             * The general method is more "elegant", at the cost of more storage requirement.
             */
        if (D1 != null && D2 != null) {
            final int srcDim;
            // Ordinate values in source space (before projection)
            final double s1, s2;
            switch(i) {
                // Horizontal segment
                case 1:
                // Horizontal segment
                case 2:
                // Horizontal segment
                case 5:
                // Horizontal segment
                case 6:
                    {
                        assert φ2 == φ1;
                        srcDim = 0;
                        s1 = λ1;
                        s2 = λ2;
                        break;
                    }
                // Vertical segment
                case 3:
                // Vertical segment
                case 4:
                // Vertical segment
                case 7:
                // Vertical segment
                case 8:
                    {
                        assert λ2 == λ1;
                        srcDim = 1;
                        s1 = φ1;
                        s2 = φ2;
                        break;
                    }
                default:
                    throw new AssertionError(i);
            }
            final double min, max;
            if (s1 < s2) {
                min = s1;
                max = s2;
            } else {
                min = s2;
                max = s1;
            }
            int tgtDim = 0;
            do {
                // Executed exactly twice, for dimensions 0 and 1 in the projected space.
                extremum.resolve(s1, (tgtDim == 0) ? x1 : y1, D1.getElement(tgtDim, srcDim), s2, (tgtDim == 0) ? x2 : y2, D2.getElement(tgtDim, srcDim));
                /*
                     * At this point we found the extremum of the projected line segment
                     * using a cubic curve t = A + Bs + Cs² + Ds³ approximation.  Before
                     * to add those extremum into the projected bounding box, we need to
                     * ensure that the source ordinate is inside the the original
                     * (unprojected) bounding box.
                     */
                boolean isP2 = false;
                do {
                    // Executed exactly twice, one for each point.
                    final double se = isP2 ? extremum.ex2 : extremum.ex1;
                    if (se > min && se < max) {
                        final double te = isP2 ? extremum.ey2 : extremum.ey1;
                        if ((tgtDim == 0) ? (te < xmin || te > xmax) : (te < ymin || te > ymax)) {
                            /*
                                 * At this point, we have determined that adding the extremum point
                                 * to the rectangle would have expanded it. However we will not add
                                 * that point directly, because maybe its position is not quite right
                                 * (since we used a cubic curve approximation). Instead, we project
                                 * the point on the rectangle border which is located vis-à-vis the
                                 * extremum. Our tests show that the correction can be as much as 50
                                 * metres.
                                 */
                            final double oldX = point[0];
                            final double oldY = point[1];
                            if (srcDim == 0) {
                                point[0] = se;
                                // == φ2 since we have an horizontal segment.
                                point[1] = φ1;
                            } else {
                                // == λ2 since we have a vertical segment.
                                point[0] = λ1;
                                point[1] = se;
                            }
                            transform.transform(point, 0, point, 0, 1);
                            final double x = point[0];
                            final double y = point[1];
                            if (x < xmin)
                                xmin = x;
                            if (x > xmax)
                                xmax = x;
                            if (y < ymin)
                                ymin = y;
                            if (y > ymax)
                                ymax = y;
                            point[0] = oldX;
                            point[1] = oldY;
                        }
                    }
                } while ((isP2 = !isP2) == true);
            } while (++tgtDim == 1);
        }
        λ1 = λ2;
        x1 = x2;
        φ1 = φ2;
        y1 = y2;
        D1 = D2;
    }
    if (destination != null) {
        destination.setRect(xmin, ymin, xmax - xmin, ymax - ymin);
    } else {
        destination = new IntervalRectangle(xmin, ymin, xmax, ymax);
    }
    /*
         * Note: a previous version had an "assert" statement here comparing our calculation
         * with the calculation performed by the more general method working on Envelope. We
         * verified that the same values (coordinate points and derivatives) were ultimately
         * passed to the CurveExtremum.resolve(…) method, so we would expect the same result.
         * However the iteration order is different. The result seems insensitive to iteration
         * order most of the time, but not always. However, it seems that the cases were the
         * results are different are the cases where the methods working with CoordinateOperation
         * object wipe out that difference anyway.
         */
    return destination;
}

Example 22 with Matrix

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

the class TransformResultComparator method derivative.

/**
 * Delegates to the tested implementation and verifies that the value is equals
 * to the one provided by the reference implementation.
 */
@Override
public Matrix derivative(DirectPosition point) throws MismatchedDimensionException, TransformException {
    final Matrix value = tested.derivative(point);
    assertMatrixEquals("derivative", reference.derivative(point), value, tolerance);
    return value;
}

Example 23 with Matrix

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

the class PassThroughTransformTest method testIdentity.

/**
 * Tests the pass through transform using an identity transform.
 * The "pass-through" of such transform shall be itself the identity transform.
 *
 * @throws TransformException should never happen.
 */
@Test
public void testIdentity() throws TransformException {
    final Matrix matrix = new Matrix3();
    runTest(MathTransforms.linear(matrix), IdentityTransform.class);
}

Example 24 with Matrix

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

the class ProjectiveTransformTest method ensureImplementRightInterface.

/**
 * Executed after every test in order to ensure that the {@linkplain #transform transform}
 * implements the {@link MathTransform1D} or {@link MathTransform2D} interface as needed.
 * In addition, all Apache SIS classes for linear transforms shall implement
 * {@link LinearTransform} and {@link Parameterized} interfaces.
 */
@After
public final void ensureImplementRightInterface() {
    if (transform instanceof TransformResultComparator) {
        transform = ((TransformResultComparator) transform).tested;
    }
    /*
         * Below is a copy of MathTransformTestCase.validate(), with minor modifications
         * due to the fact that this class does not extend MathTransformTestCase.
         */
    assertNotNull("The 'transform' field shall be assigned a value.", transform);
    Validators.validate(transform);
    final int dimension = transform.getSourceDimensions();
    if (transform.getTargetDimensions() == dimension && !skipInterfaceCheckForDimension(dimension)) {
        assertEquals("MathTransform1D", dimension == 1, (transform instanceof MathTransform1D));
        assertEquals("MathTransform2D", dimension == 2, (transform instanceof MathTransform2D));
    } else {
        assertFalse("MathTransform1D", transform instanceof MathTransform1D);
        assertFalse("MathTransform2D", transform instanceof MathTransform2D);
    }
    assertInstanceOf("Parameterized", Parameterized.class, transform);
    /*
         * End of MathTransformTestCase.validate(). Remaining is specific to LinearTransform implementations.
         */
    assertInstanceOf("Not a LinearTransform.", LinearTransform.class, transform);
    final Matrix tm = ((LinearTransform) transform).getMatrix();
    assertTrue("The matrix declared by the MathTransform is not equal to the one given at creation time.", Matrices.equals(matrix, tm, tolerance, false));
    assertSame("ParameterDescriptor", Affine.getProvider(transform.getSourceDimensions(), transform.getTargetDimensions(), true).getParameters(), ((Parameterized) transform).getParameterDescriptors());
}

Example 25 with Matrix

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

the class ServicesForMetadata method toFormattableObject.

/**
 * Converts the given object in a {@code FormattableObject} instance. Callers should verify that the given
 * object is not already an instance of {@code FormattableObject} before to invoke this method. This method
 * returns {@code null} if it can not convert the object.
 *
 * @param  object    the object to wrap.
 * @param  internal  {@code true} if the formatting convention is {@code Convention.INTERNAL}.
 * @return the given object converted to a {@code FormattableObject} instance, or {@code null}.
 *
 * @since 0.6
 */
@Override
public FormattableObject toFormattableObject(final MathTransform object, boolean internal) {
    Matrix matrix;
    final ParameterValueGroup parameters;
    if (internal && (matrix = MathTransforms.getMatrix(object)) != null) {
        parameters = Affine.parameters(matrix);
    } else if (object instanceof Parameterized) {
        parameters = ((Parameterized) object).getParameterValues();
    } else {
        matrix = MathTransforms.getMatrix(object);
        if (matrix == null) {
            return null;
        }
        parameters = Affine.parameters(matrix);
    }
    return new FormattableObject() {

        @Override
        protected String formatTo(final Formatter formatter) {
            WKTUtilities.appendParamMT(parameters, formatter);
            return WKTKeywords.Param_MT;
        }
    };
}