diff --git a/src/main/java/org/orekit/rugged/intersection/duvenhage/MinMaxTreeTile.java b/src/main/java/org/orekit/rugged/intersection/duvenhage/MinMaxTreeTile.java
index 1a6924d1191a32afc08eeeedfcad933377ac915f..7681e2c76693fd6ed9e8cf823ac50a879ffcbbc6 100644
--- a/src/main/java/org/orekit/rugged/intersection/duvenhage/MinMaxTreeTile.java
+++ b/src/main/java/org/orekit/rugged/intersection/duvenhage/MinMaxTreeTile.java
@@ -20,6 +20,7 @@ import org.apache.commons.math3.analysis.BivariateFunction;
import org.apache.commons.math3.analysis.function.Max;
import org.apache.commons.math3.analysis.function.Min;
import org.apache.commons.math3.util.FastMath;
+import org.orekit.rugged.errors.DumpManager;
import org.orekit.rugged.raster.SimpleTile;
/** Simple implementation of a {@link org.orekit.rugged.raster.Tile}
@@ -147,6 +148,13 @@ public class MinMaxTreeTile extends SimpleTile {
final int levelJ = j >> colShift;
final int levelC = 1 + ((getLongitudeColumns() - 1) >> colShift);
+ if (DumpManager.isActive()) {
+ // the following call is performed only for its side-effect
+ // of dumping the ancestor elevation, during the final
+ // call to getElevationAtIndices
+ findAncestor(i, j, level, true);
+ }
+
return minTree[start[level] + levelI * levelC + levelJ];
}
@@ -170,10 +178,112 @@ public class MinMaxTreeTile extends SimpleTile {
final int levelJ = j >> colShift;
final int levelC = 1 + ((getLongitudeColumns() - 1) >> colShift);
+ if (DumpManager.isActive()) {
+ // the following call is performed only for its side-effect
+ // of dumping the ancestor elevation, during the final
+ // call to getElevationAtIndices
+ findAncestor(i, j, level, false);
+ }
+
return maxTree[start[level] + levelI * levelC + levelJ];
}
+ /** Find a min or max ancestor cell.
+ * <p>
+ * The min/max ancestor is the cell that reaches the min (resp. max)
+ * elevation for the sub-tile containing cell (i, j) at specified level.
+ * </p>
+ * @param i row index of the cell
+ * @param j column index of the cell
+ * @param level tree level of the sub-tile considered
+ * @param isMin if true, the min ancestor cell is desired, otherwise the max ancestor
+ * cell is desired
+ * @return row/column indices of the min/max ancestor cell
+ */
+ public int[] findAncestor(final int i, final int j, final int level, final boolean isMin) {
+
+ final int k = start.length - level;
+ int rowShift = k / 2;
+ int colShift = (k + 1) / 2;
+ int levelI = i >> rowShift;
+ int levelJ = j >> colShift;
+ int levelR = 1 + ((getLatitudeRows() - 1) >> rowShift);
+ int levelC = 1 + ((getLongitudeColumns() - 1) >> colShift);
+
+ // track the cell ancestors from merged tree at specified level up to tree at level 1
+ for (int l = level; l + 1 < start.length; ++l) {
+
+ if (isColumnMerging(l + 1)) {
+
+ --colShift;
+ levelC = 1 + ((getLongitudeColumns() - 1) >> colShift);
+ levelJ = levelJ << 1;
+
+ if (levelJ + 1 < levelC) {
+ // the cell results from a regular merging of two columns
+ if (isMin) {
+ if (minTree[start[l] + levelI * levelC + levelJ] <=
+ minTree[start[l] + levelI * levelC + levelJ + 1]) {
+ levelJ++;
+ }
+ } else {
+ if (maxTree[start[l] + levelI * levelC + levelJ] >=
+ maxTree[start[l] + levelI * levelC + levelJ + 1]) {
+ levelJ++;
+ }
+ }
+ }
+
+ } else {
+
+ --rowShift;
+ levelR = 1 + ((getLatitudeRows() - 1) >> rowShift);
+ levelI = levelI << 1;
+
+ if (levelI + 1 < levelR) {
+ // the cell results from a regular merging of two rows
+ if (isMin) {
+ if (minTree[start[l] + levelI * levelC + levelJ] <=
+ minTree[start[l] + (levelI + 1) * levelC + levelJ]) {
+ levelI++;
+ }
+ } else {
+ if (maxTree[start[l] + levelI * levelC + levelJ] >=
+ maxTree[start[l] + (levelI + 1) * levelC + levelJ]) {
+ levelI++;
+ }
+ }
+ }
+
+ }
+
+ }
+
+ // we are now at first merge level, which always results from a column merge of raw data
+ levelJ = levelJ << 1;
+ levelC = getLongitudeColumns();
+ if (levelJ + 1 < levelC) {
+ // the cell results from a regular merging of two columns
+ if (isMin) {
+ if (getElevationAtIndices(levelI, levelJ) <=
+ getElevationAtIndices(levelI, levelJ + 1)) {
+ levelJ++;
+ }
+ } else {
+ if (getElevationAtIndices(levelI, levelJ) >=
+ getElevationAtIndices(levelI, levelJ + 1)) {
+ levelJ++;
+ }
+ }
+ }
+
+ return new int[] {
+ levelI, levelJ
+ };
+
+ }
+
/** Get the deepest level at which two cells are merged in the same min/max sub-tile.
* @param i1 row index of first cell
* @param j1 column index of first cell
diff --git a/src/test/java/org/orekit/rugged/intersection/duvenhage/MinMaxTreeTileTest.java b/src/test/java/org/orekit/rugged/intersection/duvenhage/MinMaxTreeTileTest.java
index e41dcb1ee833a4e55a14405a8fd87f35ac9bc6b2..84549d1b3dda444c74949a73ecbf8b96743df142 100644
--- a/src/test/java/org/orekit/rugged/intersection/duvenhage/MinMaxTreeTileTest.java
+++ b/src/test/java/org/orekit/rugged/intersection/duvenhage/MinMaxTreeTileTest.java
@@ -18,6 +18,8 @@ package org.orekit.rugged.intersection.duvenhage;
import java.lang.reflect.Field;
+import org.apache.commons.math3.random.RandomGenerator;
+import org.apache.commons.math3.random.Well1024a;
import org.apache.commons.math3.util.FastMath;
import org.junit.Assert;
import org.junit.Test;
@@ -130,6 +132,39 @@ public class MinMaxTreeTileTest {
}
}
+ @Test
+ public void testAncestor() throws RuggedException {
+ RandomGenerator random = new Well1024a(0xca9883209c6e740cl);
+ for (int nbRows = 1; nbRows < 25; nbRows++) {
+ for (int nbColumns = 1; nbColumns < 25; nbColumns++) {
+
+ MinMaxTreeTile tile = new MinMaxTreeTileFactory().createTile();
+ tile.setGeometry(1.0, 2.0, 0.1, 0.2, nbRows, nbColumns);
+ for (int i = 0; i < nbRows; ++i) {
+ for (int j = 0; j < nbColumns; ++j) {
+ tile.setElevation(i, j, 1000.0 * random.nextDouble());
+ }
+ }
+ tile.tileUpdateCompleted();
+
+ for (int i = 0; i < tile.getLatitudeRows(); ++i) {
+ for (int j = 0; j < tile.getLongitudeColumns(); ++j) {
+ for (int l = 0; l < tile.getLevels(); ++l) {
+ int[] minAncestor = tile.findAncestor(i, j, l, true);
+ Assert.assertEquals(tile.getMinElevation(i, j, l),
+ tile.getElevationAtIndices(minAncestor[0], minAncestor[1]),
+ 1.0e-10);
+ int[] maxAncestor = tile.findAncestor(i, j, l, false);
+ Assert.assertEquals(tile.getMaxElevation(i, j, l),
+ tile.getElevationAtIndices(maxAncestor[0], maxAncestor[1]),
+ 1.0e-10);
+ }
+ }
+ }
+ }
+ }
+ }
+
@Test
public void testMergeLarge() throws RuggedException {
MinMaxTreeTile tile = createTile(1201, 1201);