use of org.opengis.coverage.CannotEvaluateException in project geotoolkit by Geomatys.
the class BandedCoverageResource method sample.
/**
* @param resource
* @param domain desired grid extent and resolution, not null.
* @param range 0-based indices of sample dimensions to read, or {@code null} or an empty sequence for reading them all.
* @return
* @throws CannotEvaluateException
* @throws DataStoreException
*/
public static GridCoverage sample(BandedCoverageResource resource, final GridGeometry domain, int... range) throws CannotEvaluateException, DataStoreException {
ArgumentChecks.ensureNonNull("resource", resource);
ArgumentChecks.ensureNonNull("domain", domain);
try {
final BandedCoverage coverage = resource.read(domain, range);
if (coverage.getSampleDimensions().size() != 1) {
throw new CannotEvaluateException("Only single band sampling supported in current implementation.");
}
final GridExtent extent = domain.getExtent();
final int width = Math.toIntExact(extent.getSize(0));
final int height = Math.toIntExact(extent.getSize(1));
final long lowX = extent.getLow(0);
final long lowY = extent.getLow(1);
final MathTransform gridToCRS = domain.getGridToCRS(PixelInCell.CELL_CENTER);
// Verify no overflow is possible before allowing any array
final int nbPts = Math.multiplyExact(width, height);
final int xyLength = Math.multiplyExact(nbPts, 2);
final double[] xyGrid = new double[xyLength];
final double[] z = new double[nbPts];
for (int y = 0; y < height; y++) {
for (int x = 0; x < width; x++) {
int idx = (y * width + x) * 2;
xyGrid[idx] = lowX + x;
xyGrid[idx + 1] = lowY + y;
}
}
// convert to resulting coverage crs
final double[] xyResult;
final CoordinateReferenceSystem crs2d = CRS.getHorizontalComponent(resource.getEnvelope().get().getCoordinateReferenceSystem());
final CoordinateReferenceSystem gridCrs2d = CRS.getHorizontalComponent(domain.getCoordinateReferenceSystem());
if (!Utilities.equalsIgnoreMetadata(gridCrs2d, crs2d)) {
MathTransform trs = CRS.findOperation(gridCrs2d, crs2d, null).getMathTransform();
trs = MathTransforms.concatenate(gridToCRS, trs);
xyResult = xyGrid;
trs.transform(xyResult, 0, xyResult, 0, xyResult.length / 2);
} else {
gridToCRS.transform(xyGrid, 0, xyGrid, 0, xyGrid.length / 2);
xyResult = xyGrid;
}
final ThreadLocal<Evaluator> tl = new ThreadLocal<>();
IntStream.range(0, xyResult.length / 2).parallel().forEach((int i) -> {
Evaluator evaluator = tl.get();
if (evaluator == null) {
evaluator = coverage.evaluator();
evaluator.setNullIfOutside(true);
tl.set(evaluator);
}
final Vector3d dp = new Vector3d();
dp.x = xyResult[i * 2];
dp.y = xyResult[i * 2 + 1];
final double[] sample = evaluator.apply(dp);
if (sample != null) {
z[i] = sample[0];
} else {
z[i] = Double.NaN;
}
});
return new GridCoverageBuilder().setDomain(domain).setRanges(coverage.getSampleDimensions()).setValues(new DataBufferDouble(z, z.length), new Dimension(width, height)).build();
} catch (TransformException | FactoryException ex) {
throw new CannotEvaluateException(ex.getMessage(), ex);
}
}
use of org.opengis.coverage.CannotEvaluateException in project geotoolkit by Geomatys.
the class GridCoverageStack method seek.
/**
* Loads coverages required for a linear interpolation at the specified <var>z</var> value.
* The loaded coverages will be stored in {@link #lower} and {@link #upper} fields. It is
* possible that the same coverage is given to those two fields, if this method determine
* that no interpolation is necessary.
*
* @param z The z value.
* @return {@code true} if data were found.
* @throws PointOutsideCoverageException if the <var>z</var> value is outside the allowed range.
* @throws CannotEvaluateException if the operation failed for some other reason.
*/
private boolean seek(final double z) throws CannotEvaluateException {
assert Thread.holdsLock(this);
/*
* Check if currently loaded coverages
* are valid for the requested z value.
*/
if ((z >= lowerZ && z <= upperZ) || (isNaN(z) && isNaN(lowerZ) && isNaN(upperZ))) {
return true;
}
/*
* Currently loaded coverages are not valid for the requested z value.
* Search for the coverage to use as upper bounds ({@link #upper}).
*/
final Number Z = Double.valueOf(z);
int index;
try {
index = Arrays.binarySearch(elements, Z, COMPARATOR);
} catch (BackingStoreException exception) {
// TODO: localize
throw new CannotEvaluateException("Can't fetch coverage properties.", exception.unwrapOrRethrow(IOException.class));
}
try {
if (index >= 0) {
/*
* An exact match has been found.
* Load only this coverage and exit.
*/
load(index);
return true;
}
// Insertion point (note: ~ is NOT the minus sign).
index = ~index;
if (index == elements.length) {
if (--index >= 0) {
/*
* The requested z is after the last coverage's central z.
* Maybe it is not after the last coverage's upper z. Check...
*/
if (elements[index].getZRange().containsAny(Z)) {
load(index);
return true;
}
}
// fall through the exception at this method's end.
} else if (index == 0) {
/*
* The requested z is before the first coverage's central z.
* Maybe it is not before the first coverage's lower z. Check...
*/
if (elements[index].getZRange().containsAny(Z)) {
load(index);
return true;
}
// fall through the exception at this method's end.
} else {
/*
* An interpolation between two coverages seems possible.
* Checks if there is not a z lag between both.
*/
final Element lowerElement = elements[index - 1];
final Element upperElement = elements[index];
final NumberRange<?> lowerRange = lowerElement.getZRange();
final NumberRange<?> upperRange = upperElement.getZRange();
final double lowerEnd = lowerRange.getMaxDouble();
final double upperStart = upperRange.getMinDouble();
if (lowerEnd + lagTolerance >= upperStart) {
if (interpolationEnabled) {
load(lowerElement, upperElement);
} else {
if (Math.abs(getZ(upperElement) - z) > Math.abs(z - getZ(lowerElement))) {
index--;
}
load(index);
}
return true;
}
if (lowerRange.containsAny(Z)) {
load(index - 1);
return true;
}
if (upperRange.containsAny(Z)) {
load(index);
return true;
}
// Missing data.
return false;
}
} catch (IOException exception) {
String message = exception.getLocalizedMessage();
if (message == null) {
message = Classes.getShortClassName(exception);
}
throw new CannotEvaluateException(message, exception);
}
final Object Zp;
if (zCRS instanceof TemporalCRS) {
Zp = DefaultTemporalCRS.castOrCopy((TemporalCRS) zCRS).toDate(z);
} else {
Zp = Z;
}
throw new PointOutsideCoverageException(Errors.format(Errors.Keys.ZvalueOutsideCoverage_2, "stack", Zp));
}
use of org.opengis.coverage.CannotEvaluateException in project geotoolkit by Geomatys.
the class GridCoverageStack method snap.
/**
* Snaps the specified coordinate point to the coordinate of the nearest voxel available in
* this coverage. Invoking any {@code evaluate(...)} method with snapped coordinates will
* return non-interpolated values.
* <p>
* <ul>
* <li>First, this method locates the {@linkplain Element coverage element} at or near
* the last ordinate value (the <var>z</var> value). If no coverage is available at
* the specified <var>z</var> value, then the nearest one is selected.</li>
* <li>Next, this method locates the pixel under the {@code point} coordinate in the
* coverage element.</li>
* <li>The {@code point} is then set to the pixel center coordinate and to the
* <var>z</var> value of the selected coverage element.</li>
* </ul>
*
* @param point The point to snap.
* @throws IOException if an I/O operation was required but failed.
*/
public void snap(final DirectPosition point) throws IOException {
// No synchronization needed.
double z = point.getOrdinate(zDimension);
int index;
try {
index = Arrays.binarySearch(elements, Double.valueOf(z), COMPARATOR);
} catch (BackingStoreException exception) {
throw exception.unwrapOrRethrow(IOException.class);
}
if (index < 0) {
/*
* There is no exact match for the z value.
* Snap it to the closest coverage element.
*/
index = ~index;
if (index == elements.length) {
if (index == 0) {
// No elements in this coverage
return;
}
z = getZ(elements[--index]);
} else if (index == 0) {
z = getZ(elements[index]);
} else {
final double lowerZ = getZ(elements[index - 1]);
final double upperZ = getZ(elements[index]);
// Use !(...) in order to accept NaN values.
assert !(z <= lowerZ || z >= upperZ) : z;
if (isNaN(upperZ) || z - lowerZ < upperZ - z) {
index--;
z = lowerZ;
} else {
z = upperZ;
}
}
point.setOrdinate(zDimension, z);
}
/*
* Now that we know the coverage element,
* snap the spatial coordinate point.
*/
final Element element = elements[index];
final GridGeometry geometry = element.getGridGeometry();
if (geometry != null) {
final GridExtent range = geometry.getExtent();
final MathTransform transform = geometry.getGridToCRS(PixelInCell.CELL_CENTER);
final int dimension = transform.getSourceDimensions();
DirectPosition position = new GeneralDirectPosition(dimension);
for (int i = dimension; --i >= 0; ) {
// Copy only the first dimensions (may not be up to crs.dimension)
position.setOrdinate(i, point.getOrdinate(i));
}
try {
position = transform.inverse().transform(position, position);
for (int i = dimension; --i >= 0; ) {
position.setOrdinate(i, Math.max(range.getLow(i), Math.min(range.getHigh(i), (int) Math.rint(position.getOrdinate(i)))));
}
position = transform.transform(position, position);
for (int i = Math.min(dimension, zDimension); --i >= 0; ) {
// Do not touch the z-value, copy the other coordinates.
point.setOrdinate(i, position.getOrdinate(i));
}
} catch (TransformException exception) {
throw new CannotEvaluateException(cannotEvaluate(point), exception);
}
}
}
Aggregations