From 2ca256db9e6b93f373bc08ee4b59deda51fd22e4 Mon Sep 17 00:00:00 2001
From: Luc Maisonobe <luc@orekit.org>
Date: Tue, 3 Feb 2015 22:08:02 +0100
Subject: [PATCH] Prepared a way to dump min/max ancestors.

The new test for the ancestor method fails, but the failure is probably
related to issue#189 (in fact the issue was identified while testing
this). So we will first solve issue#189 in another branch and merge it
once solved, before continuing working on this branch.
---
 .../duvenhage/MinMaxTreeTile.java             | 110 ++++++++++++++++++
 .../duvenhage/MinMaxTreeTileTest.java         |  35 ++++++
 2 files changed, 145 insertions(+)

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 1a6924d1..7681e2c7 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 e41dcb1e..84549d1b 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);
-- 
GitLab