Examples with ReducedGridCoverage org.geotoolkit.coverage.ReducedGridCoverage used on opensource projects

Example 1 with ReducedGridCoverage

use of org.geotoolkit.coverage.ReducedGridCoverage in project geotoolkit by Geomatys.

the class AbstractCoverageSymbolizerRenderer method getObjectiveCoverage.

/**
 * Returns expected {@linkplain GridCoverage elevation coverage} or {@linkplain GridCoverage coverage}
 * from given {@link ProjectedCoverage}.
 *
 * TODO: add a margin or interpolation parameter. To properly interpolate border on output canvas, we need extra
 *  lines/columns on source image. Their number depends on applied interpolation (bilinear, bicubic, etc.).
 *
 * @param ref coverage resource
 * @param canvasGrid Rendering canvas grid geometry
 * @param isElevation {@code true} if we want elevation coverage, else ({@code false}) for features coverage.
 * @param sourceBands coverage source bands to features
 * @return expected {@linkplain GridCoverage elevation coverage} or {@linkplain GridCoverage coverage}
 * @throws org.geotoolkit.coverage.io.CoverageStoreException if problem during coverage reading.
 * @throws org.opengis.referencing.operation.TransformException if problem during {@link Envelope} transformation.
 * @throws org.opengis.util.FactoryException if problem during {@link Envelope} study.
 * @throws org.geotoolkit.process.ProcessException if problem during resampling processing.
 * @see ProjectedCoverage#getElevationCoverage(org.geotoolkit.coverage.io.GridCoverageReadParam)
 * @see ProjectedCoverage#getCoverage(org.geotoolkit.coverage.io.GridCoverageReadParam)
 */
protected GridCoverage getObjectiveCoverage(final GridCoverageResource ref, GridGeometry canvasGrid, final boolean isElevation, int[] sourceBands) throws DataStoreException, TransformException, FactoryException, ProcessException {
    ArgumentChecks.ensureNonNull("projectedCoverage", ref);
    final InterpolationCase interpolation = InterpolationCase.BILINEAR;
    final GridGeometry refGG = ref.getGridGeometry();
    // fast envelope intersection in 2D
    if (refGG.isDefined(GridGeometry.ENVELOPE)) {
        Envelope bbox = renderingContext.getCanvasObjectiveBounds2D();
        Envelope refEnv = Envelopes.transform(refGG.getEnvelope(), bbox.getCoordinateReferenceSystem());
        if (!AbstractEnvelope.castOrCopy(bbox).intersects(refEnv, true)) {
            throw new DisjointExtentException("Coverage resource envelope do not intersect canvas");
        }
    }
    final GridGeometry baseGG = trySubGrid(refGG, canvasGrid);
    final GridGeometry slice = extractSlice(baseGG, canvasGrid, computeMargin2D(interpolation), true);
    if (sourceBands != null && sourceBands.length < 1)
        sourceBands = null;
    GridCoverage coverage = ref.read(slice, sourceBands);
    if (coverage instanceof GridCoverageStack) {
        Logger.getLogger("org.geotoolkit.display2d.primitive").log(Level.WARNING, "Coverage reader return more than one slice.");
    }
    while (coverage instanceof GridCoverageStack) {
        // pick the first slice
        coverage = ((GridCoverageStack) coverage).coverageAtIndex(0);
    }
    // we remove all other dimension to simplify any following operation
    if (coverage.getCoordinateReferenceSystem().getCoordinateSystem().getDimension() > 2) {
        coverage = new ReducedGridCoverage(coverage, 0, 1);
    }
    final CoordinateReferenceSystem crs2d = CRS.getHorizontalComponent(canvasGrid.getCoordinateReferenceSystem());
    if (Utilities.equalsIgnoreMetadata(crs2d, coverage.getCoordinateReferenceSystem())) {
        return coverage;
    } else {
        coverage = prepareCoverageToResampling(coverage, symbol);
        // resample
        final double[] fill = new double[coverage.getSampleDimensions().size()];
        Arrays.fill(fill, Double.NaN);
        // ///// HACK FOR 0/360 /////////////////////////////////////////
        GeneralEnvelope ge = new GeneralEnvelope(coverage.getGridGeometry().getEnvelope());
        try {
            GeneralEnvelope cdt = GeneralEnvelope.castOrCopy(Envelopes.transform(coverage.getGridGeometry().getEnvelope(), CommonCRS.WGS84.normalizedGeographic()));
            cdt.normalize();
            if (!cdt.isEmpty()) {
                ge = cdt;
            }
        } catch (ProjectionException ex) {
            LOGGER.log(Level.INFO, ex.getMessage(), ex);
        }
        GridGeometry resampleGrid = canvasGrid;
        try {
            resampleGrid = resampleGrid.derive().rounding(GridRoundingMode.ENCLOSING).subgrid(ge).build().reduce(0, 1);
        } catch (DisjointExtentException ex) {
        // don't log, still continue
        } catch (IllegalGridGeometryException ex) {
            LOGGER.log(Level.INFO, ex.getMessage(), ex);
        }
        resampleGrid = CoverageUtilities.forceLowerToZero(resampleGrid);
        // ///// HACK FOR 0/360 /////////////////////////////////////////
        final MathTransform gridToCRS = coverage.getGridGeometry().getGridToCRS(PixelInCell.CELL_CENTER);
        if (isNonLinear(gridToCRS)) {
            final GridGeometry slice2 = extractSlice(refGG, canvasGrid, computeMargin2D(interpolation), false);
            coverage = ref.read(slice2, sourceBands);
            // we remove all other dimension to simplify any following operation
            if (coverage.getCoordinateReferenceSystem().getCoordinateSystem().getDimension() > 2) {
                coverage = new ReducedGridCoverage(coverage, 0, 1);
            }
            return forwardResample(coverage, resampleGrid);
        } else {
            ResampleProcess process = new ResampleProcess(coverage, crs2d, resampleGrid, interpolation, fill);
            // do not extrapolate values, can cause large areas of incorrect values
            process.getInput().parameter(ResampleDescriptor.IN_BORDER_COMPORTEMENT_TYPE.getName().getCode()).setValue(ResampleBorderComportement.FILL_VALUE);
            return process.executeNow();
        }
    }
}