Search in sources :

Example 6 with GridGeometry

use of org.apache.sis.coverage.grid.GridGeometry in project sis by apache.

the class CoverageReadConsistency method readAndCompareRandomRegions.

/**
 * Implementation of methods testing reading in random sub-regions with random sub-samplings.
 *
 * @param  label  a label for the test being run.
 * @throws DataStoreException if an error occurred while using the resource.
 */
private void readAndCompareRandomRegions(final String label) throws DataStoreException {
    randomConfigureResource();
    final GridGeometry gg = resource.getGridGeometry();
    final int dimension = gg.getDimension();
    final long[] low = new long[dimension];
    final long[] high = new long[dimension];
    final int[] subsampling = new int[dimension];
    final int[] subOffsets = new int[dimension];
    final int numBands = resource.getSampleDimensions().size();
    /*
         * We will collect statistics on execution time only if the
         * test is executed in a more verbose mode than the default.
         */
    final Statistics durations = (VERBOSE || !failOnMismatch) ? new Statistics(label) : null;
    int failuresCount = 0;
    for (int it = 0; it < numIterations; it++) {
        final GridGeometry domain = randomDomain(gg, low, high, subsampling);
        final int[] selectedBands = randomRange(numBands);
        /*
             * Read a coverage containing the requested sub-domain. Note that the reader is free to read
             * more data than requested. The extent actually read is `actualReadExtent`. It shall contain
             * fully the requested `domain`.
             */
        final long startTime = System.nanoTime();
        final GridCoverage subset = resource.read(domain, selectedBands);
        final GridExtent actualReadExtent = subset.getGridGeometry().getExtent();
        if (failOnMismatch) {
            assertEquals("Unexpected number of dimensions.", dimension, actualReadExtent.getDimension());
            for (int d = 0; d < dimension; d++) {
                if (subsampling[d] == 1) {
                    assertTrue("Actual extent is too small.", actualReadExtent.getSize(d) > high[d] - low[d]);
                    assertTrue("Actual extent is too small.", actualReadExtent.getLow(d) <= low[d]);
                    assertTrue("Actual extent is too small.", actualReadExtent.getHigh(d) >= high[d]);
                }
            }
        }
        /*
             * If subsampling was enabled, the factors selected by the reader may be different than
             * the subsampling factors that we specified. The following block updates those values.
             */
        if (allowSubsampling && full != null) {
            final GridDerivation change = full.getGridGeometry().derive().subgrid(subset.getGridGeometry());
            System.arraycopy(change.getSubsampling(), 0, subsampling, 0, dimension);
            System.arraycopy(change.getSubsamplingOffsets(), 0, subOffsets, 0, dimension);
        }
        /*
             * Iterate over all dimensions greater than 2. In the common case where we are reading a
             * two-dimensional image, the following loop will be executed only once. If reading a 3D
             * or 4D image, the loop is executed for all possible two-dimensional slices in the cube.
             */
        final int sd = actualReadExtent.getDimension();
        final long[] sliceMin = new long[sd];
        final long[] sliceMax = new long[sd];
        for (int i = 0; i < sd; i++) {
            sliceMin[i] = actualReadExtent.getLow(i);
            sliceMax[i] = actualReadExtent.getHigh(i);
        }
        nextSlice: for (; ; ) {
            System.arraycopy(sliceMin, BIDIMENSIONAL, sliceMax, BIDIMENSIONAL, dimension - BIDIMENSIONAL);
            final PixelIterator itr = iterator(full, sliceMin, sliceMax, subsampling, subOffsets, allowSubsampling);
            final PixelIterator itc = iterator(subset, sliceMin, sliceMax, subsampling, subOffsets, false);
            if (itr != null) {
                assertEquals(itr.getDomain().getSize(), itc.getDomain().getSize());
                final double[] expected = new double[selectedBands.length];
                double[] reference = null, actual = null;
                while (itr.next()) {
                    assertTrue(itc.next());
                    reference = itr.getPixel(reference);
                    actual = itc.getPixel(actual);
                    for (int i = 0; i < selectedBands.length; i++) {
                        expected[i] = reference[selectedBands[i]];
                    }
                    if (!Arrays.equals(expected, actual)) {
                        failuresCount++;
                        if (!failOnMismatch)
                            break;
                        final Point pr = itr.getPosition();
                        final Point pc = itc.getPosition();
                        final StringBuilder message = new StringBuilder(100).append("Mismatch at position (").append(pr.x).append(", ").append(pr.y).append(") in full image and (").append(pc.x).append(", ").append(pc.y).append(") in tested sub-image");
                        findMatchPosition(itr, pr, selectedBands, actual, message);
                        assertArrayEquals(message.toString(), expected, actual, STRICT);
                    /*
                             * POSSIBLE CAUSES FOR TEST FAILURE (known issues):
                             *
                             *   - If the `GridGeometry` has no `gridToCRS` transform, then `GridDerivation` manages
                             *     to save the scales (subsampling factors) anyway but the translations (subsampling
                             *     offsets) are lost. It causes an image shift if the offsets were not zero. Because
                             *     `gridToCRS` should never be null with spatial data, we do not complexify the code
                             *     for what may be a non-issue.
                             */
                    }
                }
                assertFalse(itc.next());
            } else {
                // Unable to create a reference image. Just check that no exception is thrown.
                double[] actual = null;
                while (itc.next()) {
                    actual = itc.getPixel(actual);
                }
            }
            /*
                 * Move to the next two-dimensional slice and read again.
                 * We stop the loop after we have read all 2D slices.
                 */
            for (int d = dimension; --d >= BIDIMENSIONAL; ) {
                if (sliceMin[d]++ <= actualReadExtent.getHigh(d))
                    continue nextSlice;
                sliceMin[d] = actualReadExtent.getLow(d);
            }
            break;
        }
        if (durations != null) {
            durations.accept((System.nanoTime() - startTime) / (double) StandardDateFormat.NANOS_PER_MILLISECOND);
        }
    }
    /*
         * Show statistics only if the test are executed with the `VERBOSE` flag set,
         * or if this `CoverageReadConsistency` is used for benchmark.
         */
    if (durations != null) {
        if (statistics == null) {
            statistics = new ArrayList<>();
        }
        statistics.add(durations);
        final int totalCount = durations.count();
        out.println("Number of failures: " + failuresCount + " / " + totalCount + " (" + (failuresCount / (totalCount / 100f)) + "%)");
    }
}
Also used : GridGeometry(org.apache.sis.coverage.grid.GridGeometry) GridExtent(org.apache.sis.coverage.grid.GridExtent) PixelIterator(org.apache.sis.image.PixelIterator) Point(java.awt.Point) Statistics(org.apache.sis.math.Statistics) Point(java.awt.Point) GridDerivation(org.apache.sis.coverage.grid.GridDerivation) GridCoverage(org.apache.sis.coverage.grid.GridCoverage)

Example 7 with GridGeometry

use of org.apache.sis.coverage.grid.GridGeometry in project sis by apache.

the class CoverageQueryTest method testWithSubgrid.

/**
 * Tests with a sub-grid geometry.
 *
 * @throws DataStoreException if query execution failed.
 */
@Test
public void testWithSubgrid() throws DataStoreException {
    final GridGeometry subGrid = createSubGrid(0);
    final CoverageQuery query = new CoverageQuery();
    query.setSelection(subGrid);
    final GridCoverageResource subset = resource.subset(query);
    assertEquals(subGrid, subset.getGridGeometry());
    verifyRead(subset, 0);
}
Also used : GridGeometry(org.apache.sis.coverage.grid.GridGeometry) Test(org.junit.Test)

Example 8 with GridGeometry

use of org.apache.sis.coverage.grid.GridGeometry in project sis by apache.

the class CoverageQueryTest method verifyRead.

/**
 * Verifies that the read operation adds the expected margins.
 */
private void verifyRead(final GridCoverageResource subset, final int expansion) throws DataStoreException {
    /*
         * Test reading the whole image. The grid geometry of returned coverage
         * must be the same than the grid geometry of the GridCoverageResource,
         * which has been verified by the caller to contain the expansion.
         */
    assertEquals(subset.getGridGeometry(), subset.read(null).getGridGeometry());
    /*
         * Request for a smaller area and verify that the request has the expected size,
         * including expansion.
         */
    GridGeometry request = createSubGrid(-4);
    final GridCoverage coverage = subset.read(request);
    if (expansion != 0) {
        request = createSubGrid(-4 + expansion);
    }
    assertEquals(request, coverage.getGridGeometry());
}
Also used : GridGeometry(org.apache.sis.coverage.grid.GridGeometry) GridCoverage(org.apache.sis.coverage.grid.GridCoverage)

Example 9 with GridGeometry

use of org.apache.sis.coverage.grid.GridGeometry in project sis by apache.

the class GridMapping method adaptGridCRS.

/**
 * Creates the grid geometry from the {@link #crs} and {@link #gridToCRS} field,
 * completing missing information with the given template.
 *
 * @param  variable  the variable for which to create a grid geometry.
 * @param  template  template to use for completing missing information.
 * @param  anchor    whether we computed "grid to CRS" transform relative to pixel center or pixel corner.
 * @return the grid geometry with modified CRS and "grid to CRS" transform, or {@code null} in case of failure.
 */
GridGeometry adaptGridCRS(final Variable variable, final GridGeometry template, final PixelInCell anchor) {
    CoordinateReferenceSystem givenCRS = crs;
    int firstAffectedCoordinate = 0;
    boolean isSameGrid = true;
    if (template.isDefined(GridGeometry.CRS)) {
        final CoordinateReferenceSystem templateCRS = template.getCoordinateReferenceSystem();
        if (givenCRS == null) {
            givenCRS = templateCRS;
        } else {
            /*
                 * The CRS built by Grid may have a different axis order than the CRS specified by grid mapping attributes.
                 * Check which axis order seems to fit, then replace grid CRS by given CRS (potentially with swapped axes).
                 * This is where the potential difference between EPSG axis order and grid axis order is handled. If we can
                 * not find where to substitute the CRS, assume that the given CRS describes the first dimensions. We have
                 * no guarantees that this later assumption is right, but it seems to match common practice.
                 */
            final CoordinateSystem cs = templateCRS.getCoordinateSystem();
            CoordinateSystem subCS = givenCRS.getCoordinateSystem();
            firstAffectedCoordinate = AxisDirections.indexOfColinear(cs, subCS);
            if (firstAffectedCoordinate < 0) {
                givenCRS = AbstractCRS.castOrCopy(givenCRS).forConvention(AxesConvention.RIGHT_HANDED);
                subCS = givenCRS.getCoordinateSystem();
                firstAffectedCoordinate = AxisDirections.indexOfColinear(cs, subCS);
                if (firstAffectedCoordinate < 0) {
                    firstAffectedCoordinate = 0;
                    if (!isEPSG) {
                        // If specified by WKT, use the given CRS verbatim.
                        givenCRS = crs;
                        subCS = givenCRS.getCoordinateSystem();
                    }
                }
            }
            /*
                 * Replace the grid CRS (or a component of it) by the CRS parsed from WKT or EPSG code with same (if possible)
                 * axis order. If the grid CRS contains more axes (for example elevation or time axis), we try to keep them.
                 */
            try {
                givenCRS = new GeodeticObjectBuilder(variable.decoder, variable.decoder.listeners.getLocale()).replaceComponent(templateCRS, firstAffectedCoordinate, givenCRS);
            } catch (FactoryException e) {
                canNotCreate(variable, Resources.Keys.CanNotCreateCRS_3, e);
                return null;
            }
            isSameGrid = templateCRS.equals(givenCRS);
            if (isSameGrid) {
                // Keep existing instance if appropriate.
                givenCRS = templateCRS;
            }
        }
    }
    /*
         * Perform the same substitution than above, but in the "grid to CRS" transform. Note that the "grid to CRS"
         * is usually not specified, so the block performing substitution will rarely be executed. If executed, then
         * then we need to perform selection in target dimensions (not source dimensions) because the first affected
         * coordinate computed above is in CRS dimension, which is the target of "grid to CRS" transform.
         */
    MathTransform givenG2C = gridToCRS;
    if (template.isDefined(GridGeometry.GRID_TO_CRS)) {
        final MathTransform templateG2C = template.getGridToCRS(anchor);
        if (givenG2C == null) {
            givenG2C = templateG2C;
        } else
            try {
                int count = 0;
                MathTransform[] components = new MathTransform[3];
                final TransformSeparator sep = new TransformSeparator(templateG2C, variable.decoder.getMathTransformFactory());
                if (firstAffectedCoordinate != 0) {
                    sep.addTargetDimensionRange(0, firstAffectedCoordinate);
                    components[count++] = sep.separate();
                    sep.clear();
                }
                components[count++] = givenG2C;
                final int next = firstAffectedCoordinate + givenG2C.getTargetDimensions();
                final int upper = templateG2C.getTargetDimensions();
                if (next != upper) {
                    sep.addTargetDimensionRange(next, upper);
                    components[count++] = sep.separate();
                }
                components = ArraysExt.resize(components, count);
                givenG2C = MathTransforms.compound(components);
                if (templateG2C.equals(givenG2C)) {
                    // Keep using existing instance if appropriate.
                    givenG2C = templateG2C;
                } else {
                    isSameGrid = false;
                }
            } catch (FactoryException e) {
                canNotCreate(variable, Resources.Keys.CanNotCreateGridGeometry_3, e);
                return null;
            }
    }
    /*
         * At this point we finished to compute the grid geometry components.
         * If any of them have changed, create the new grid geometry.
         */
    if (isSameGrid) {
        return template;
    } else {
        return new GridGeometry(template.getExtent(), anchor, givenG2C, givenCRS);
    }
}
Also used : GridGeometry(org.apache.sis.coverage.grid.GridGeometry) MathTransform(org.opengis.referencing.operation.MathTransform) FactoryException(org.opengis.util.FactoryException) CoordinateSystem(org.opengis.referencing.cs.CoordinateSystem) GeodeticObjectBuilder(org.apache.sis.internal.referencing.GeodeticObjectBuilder) CoordinateReferenceSystem(org.opengis.referencing.crs.CoordinateReferenceSystem) TransformSeparator(org.apache.sis.referencing.operation.transform.TransformSeparator)

Example 10 with GridGeometry

use of org.apache.sis.coverage.grid.GridGeometry in project sis by apache.

the class GridMapping method createGridCRS.

/**
 * Creates a new grid geometry for the given extent.
 * This method should be invoked only when no existing {@link GridGeometry} can be used as template.
 */
GridGeometry createGridCRS(final Variable variable) {
    final List<Dimension> dimensions = variable.getGridDimensions();
    final long[] upper = new long[dimensions.size()];
    for (int i = 0; i < upper.length; i++) {
        // Convert CRS dimension to netCDF dimension.
        final int d = (upper.length - 1) - i;
        upper[i] = dimensions.get(d).length();
    }
    return new GridGeometry(new GridExtent(null, null, upper, false), PixelInCell.CELL_CENTER, gridToCRS, crs);
}
Also used : GridGeometry(org.apache.sis.coverage.grid.GridGeometry) GridExtent(org.apache.sis.coverage.grid.GridExtent)

Aggregations

GridGeometry (org.apache.sis.coverage.grid.GridGeometry)25 GridExtent (org.apache.sis.coverage.grid.GridExtent)10 SampleDimension (org.apache.sis.coverage.SampleDimension)6 CoordinateReferenceSystem (org.opengis.referencing.crs.CoordinateReferenceSystem)6 GridCoverage (org.apache.sis.coverage.grid.GridCoverage)5 MathTransform (org.opengis.referencing.operation.MathTransform)5 Test (org.junit.Test)4 TransformException (org.opengis.referencing.operation.TransformException)4 GridDerivation (org.apache.sis.coverage.grid.GridDerivation)3 Point (java.awt.Point)2 RenderedImage (java.awt.image.RenderedImage)2 GridCoverage2D (org.apache.sis.coverage.grid.GridCoverage2D)2 LinearTransform (org.apache.sis.referencing.operation.transform.LinearTransform)2 DataStoreContentException (org.apache.sis.storage.DataStoreContentException)2 FactoryException (org.opengis.util.FactoryException)2 AffineTransform (java.awt.geom.AffineTransform)1 BufferedImage (java.awt.image.BufferedImage)1 DataBuffer (java.awt.image.DataBuffer)1 RasterFormatException (java.awt.image.RasterFormatException)1 IOException (java.io.IOException)1