sis-commits mailing list archives

Site index · List index
Message view « Date » · « Thread »
Top « Date » · « Thread »
From desruisse...@apache.org
Subject [sis] 01/01: First attempt to implement bicubic interpolation. Current code produces wrong results; more work is needed before potential merge.
Date Mon, 06 Apr 2020 17:25:59 GMT
This is an automated email from the ASF dual-hosted git repository.

desruisseaux pushed a commit to branch feat/bicubic
in repository https://gitbox.apache.org/repos/asf/sis.git

commit c18d424196f43f30a72d4298e3ce932bfb1041f1
Author: Martin Desruisseaux <martin.desruisseaux@geomatys.com>
AuthorDate: Mon Apr 6 19:24:32 2020 +0200

    First attempt to implement bicubic interpolation. Current code produces wrong results;
more work is needed before potential merge.
    
    https://issues.apache.org/jira/browse/SIS-494
---
 .../org/apache/sis/image/BicubicInterpolation.java |  96 +++++++++++++++++
 .../java/org/apache/sis/image/Interpolation.java   |  21 ++++
 .../sis/image/BicubicReferenceInterpolation.java   | 119 +++++++++++++++++++++
 .../org/apache/sis/image/ResampledImageTest.java   |  44 ++++++++
 4 files changed, 280 insertions(+)

diff --git a/core/sis-feature/src/main/java/org/apache/sis/image/BicubicInterpolation.java
b/core/sis-feature/src/main/java/org/apache/sis/image/BicubicInterpolation.java
new file mode 100644
index 0000000..5ad27d9
--- /dev/null
+++ b/core/sis-feature/src/main/java/org/apache/sis/image/BicubicInterpolation.java
@@ -0,0 +1,96 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one or more
+ * contributor license agreements.  See the NOTICE file distributed with
+ * this work for additional information regarding copyright ownership.
+ * The ASF licenses this file to You under the Apache License, Version 2.0
+ * (the "License"); you may not use this file except in compliance with
+ * the License.  You may obtain a copy of the License at
+ *
+ *     http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+package org.apache.sis.image;
+
+import java.awt.Dimension;
+import java.nio.DoubleBuffer;
+
+
+/**
+ * Bicubic interpolation.
+ *
+ * @author  Rémi Marechal (Geomatys)
+ * @author  Martin Desruisseaux (Geomatys)
+ * @version 1.1
+ *
+ * @see <a href="https://en.wikipedia.org/wiki/Bicubic_interpolation">Bicubic interpolation
on Wikipedia</a>
+ *
+ * @since 1.1
+ * @module
+ */
+final class BicubicInterpolation implements Interpolation {
+    /**
+     * The number of pixel required horizontally and vertically for performing this interpolation.
+     */
+    static final int SUPPORT = 4;
+
+    /**
+     * Creates a new interpolation.
+     */
+    BicubicInterpolation() {
+    }
+
+    /**
+     * Interpolation name for debugging purpose.
+     */
+    @Override
+    public String toString() {
+        return "BICUBIC";
+    }
+
+    /**
+     * Size of the area over which to provide values.
+     */
+    @Override
+    public Dimension getSupportSize() {
+        return new Dimension(SUPPORT, SUPPORT);
+    }
+
+    /**
+     * Applies bicubic interpolation on a 4×4 window.
+     */
+    @Override
+    public boolean interpolate(final DoubleBuffer source, final int numBands,
+            double xfrac, double yfrac, final double[] writeTo, int writeToOffset)
+    {
+        xfrac++;    // TODO: update coefficients in `interpolate(…)` method instead.
+        yfrac++;
+        final double[] y = new double[SUPPORT];
+        for (int b=0; b<numBands; b++) {
+            int p = source.position() + b;
+            for (int j=0; j<SUPPORT; j++) {
+                y[j] = interpolate(xfrac, source.get(p            ),
+                                          source.get(p += numBands),
+                                          source.get(p += numBands),
+                                          source.get(p += numBands));
+                p += numBands;
+            }
+            writeTo[writeToOffset++] = interpolate(yfrac, y[0], y[1], y[2], y[3]);
+        }
+        return true;
+    }
+
+    /**
+     * Applies a bicubic interpolation on a row or a column.
+     */
+    private static double interpolate(final double t, final double s0, final double s1, final
double s2, final double s3) {
+        final double a1 = ( 1d/3)*s3 + (-3d/2)*s2 + ( 3d  )*s1 + (-11d/6)*s0;
+        final double a2 = (-1d/2)*s3 + ( 2d  )*s2 + (-5d/2)*s1 +          s0;
+        final double a3 = ( 1d/6)*s3 + (-1d/2)*s2 + ( 1d/2)*s1 + ( -1d/6)*s0;
+        return s0 + (a1 + (a2 + a3*t)*t)*t;
+    }
+}
diff --git a/core/sis-feature/src/main/java/org/apache/sis/image/Interpolation.java b/core/sis-feature/src/main/java/org/apache/sis/image/Interpolation.java
index 6cf39e8..e60f326 100644
--- a/core/sis-feature/src/main/java/org/apache/sis/image/Interpolation.java
+++ b/core/sis-feature/src/main/java/org/apache/sis/image/Interpolation.java
@@ -164,6 +164,27 @@ public interface Interpolation {
     };
 
     /**
+     * Bicubic interpolation.
+     * The bicubic convolution algorithm uses a kernel with coefficients computed as below:
+     *
+     * {@preformat text
+     *   W(x)  =  (a+2)|x|³ − (a+3)|x|²         + 1       for 0 ≤ |x| ≤ 1
+     *   W(x)  =      a|x|³ −    5a|x|² + 8a|x| − 4a      for 1 < |x| < 2
+     *   W(x)  =  0                                       otherwise
+     * }
+     *
+     * The <var>a</var> value is typically −0.5, −0.75 or −1.
+     * This interpolation instance uses <var>a</var> = −0.5.
+     *
+     * <div class="note"><b>Reference:</b>
+     * Digital Image Warping, George Wolberg, 1990, pp 129-131, IEEE Computer Society Press,
ISBN 0-8186-8944-7.
+     * </div>
+     *
+     * @see <a href="https://en.wikipedia.org/wiki/Bicubic_interpolation">Bicubic interpolation
on Wikipedia</a>
+     */
+    Interpolation BICUBIC = new BicubicInterpolation();
+
+    /**
      * Lanczos interpolation for photographic images.
      * This interpolation is not recommended for images that may contain NaN values.
      * The Lanczos reconstruction kernel is:
diff --git a/core/sis-feature/src/test/java/org/apache/sis/image/BicubicReferenceInterpolation.java
b/core/sis-feature/src/test/java/org/apache/sis/image/BicubicReferenceInterpolation.java
new file mode 100644
index 0000000..3b08c8f
--- /dev/null
+++ b/core/sis-feature/src/test/java/org/apache/sis/image/BicubicReferenceInterpolation.java
@@ -0,0 +1,119 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one or more
+ * contributor license agreements.  See the NOTICE file distributed with
+ * this work for additional information regarding copyright ownership.
+ * The ASF licenses this file to You under the Apache License, Version 2.0
+ * (the "License"); you may not use this file except in compliance with
+ * the License.  You may obtain a copy of the License at
+ *
+ *     http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+package org.apache.sis.image;
+
+import java.awt.Dimension;
+import java.nio.DoubleBuffer;
+
+import static org.apache.sis.image.BicubicInterpolation.SUPPORT;
+
+
+/**
+ * A {@link BicubicInterpolation} to be used as a reference implementation.
+ * This class implements the bicubic convolution algorithm as documented on
+ * <a href="https://en.wikipedia.org/wiki/Bicubic_interpolation#Bicubic_convolution_algorithm">Wikipedia</a>,
+ * also found in Java Advanced Imaging (JAI) as below:
+ *
+ * {@preformat text
+ *   W(x)  =  (a+2)|x|³ − (a+3)|x|²         + 1       for 0 ≤ |x| ≤ 1
+ *   W(x)  =      a|x|³ −    5a|x|² + 8a|x| − 4a      for 1 < |x| < 2
+ *   W(x)  =  0                                       otherwise
+ * }
+ *
+ * The <var>a</var> value is typically −0.5 or −1, as recommended by Rifman
and Keys respectively (Reference:
+ * Digital Image Warping, George Wolberg, 1990, pp 129-131, IEEE Computer Society Press,
ISBN 0-8186-8944-7).
+ * The −1 value often (not always) produces sharper results than the −0.5 value.
+ *
+ * <p>The {@link BicubicInterpolation} class develops coefficients for the <var>a</var>
= −0.5 case.
+ * This class can be used for comparing {@link BicubicInterpolation} results with results
calculated
+ * from original formulas.</p>
+ *
+ * @author  Rémi Marechal (Geomatys)
+ * @author  Martin Desruisseaux (Geomatys)
+ * @version 1.1
+ * @since   1.1
+ * @module
+ */
+public final strictfp class BicubicReferenceInterpolation implements Interpolation {
+    /**
+     * The <var>a</var> coefficient used in bicubic convolution algorithm.
+     * Typical values are −0.5, −0.75 or −1.
+     */
+    private final double a;
+
+    /**
+     * Creates a new interpolation.
+     */
+    BicubicReferenceInterpolation(final double a) {
+        this.a = a;
+    }
+
+    /**
+     * Interpolation name for debugging purpose.
+     */
+    @Override
+    public String toString() {
+        return "BicubicReference[a=" + a + ']';
+    }
+
+    /**
+     * Size of the area over which to provide values.
+     */
+    @Override
+    public Dimension getSupportSize() {
+        return new Dimension(SUPPORT, SUPPORT);
+    }
+
+    /**
+     * Applies bicubic interpolation on a 4×4 window.
+     */
+    @Override
+    public boolean interpolate(final DoubleBuffer source, final int numBands,
+            final double xfrac, final double yfrac, final double[] writeTo, int writeToOffset)
+    {
+        for (int b=0; b<numBands; b++) {
+            double sum = 0;
+            int p = source.position() + b;
+            for (int j=0; j<SUPPORT; j++) {
+                double s = 0;
+                for (int i=0; i<SUPPORT; i++) {
+                    s += convolutionValue(xfrac - i) * source.get(p);
+                    p += numBands;
+                }
+                sum += convolutionValue(yfrac - j) * s;
+            }
+            writeTo[writeToOffset++] = sum;
+        }
+        return true;
+    }
+
+    /**
+     * Compute a value of kernel filter to apply on current pixel value.
+     *
+     * @param  t  difference between interpolation position and pixel position.
+     */
+    private double convolutionValue(double t) {
+        t = StrictMath.abs(t);
+        if (t <= 1) {
+            return ((a+2)*t - (a+3))*(t*t) + 1;     // (a+2)|x|³ − (a+3)|x|² + 1
+        } else if (t < 2) {
+            return (((t - 5)*t + 8)*t - 4)*a;       // a|x|³ − 5a|x|² + 8a|x| − 4a
+        } else {
+            return 0;
+        }
+    }
+}
diff --git a/core/sis-feature/src/test/java/org/apache/sis/image/ResampledImageTest.java b/core/sis-feature/src/test/java/org/apache/sis/image/ResampledImageTest.java
index 546c652..3065ef1 100644
--- a/core/sis-feature/src/test/java/org/apache/sis/image/ResampledImageTest.java
+++ b/core/sis-feature/src/test/java/org/apache/sis/image/ResampledImageTest.java
@@ -235,6 +235,28 @@ public final strictfp class ResampledImageTest extends TestCase {
     }
 
     /**
+     * Tests {@link Interpolation#BICUBIC} on floating point values.
+     */
+    @Test
+    public void testBicubicOnFloats() {
+        source = createRandomImage(DataBuffer.TYPE_FLOAT);
+        interpolation = Interpolation.BICUBIC;
+        createScaledByTwo(-15, -12);
+        verifyAtIntegerPositions();
+    }
+
+    /**
+     * Tests {@link Interpolation#BICUBIC} on integer values.
+     */
+    @Test
+    public void testBicubicOnIntegers() {
+        source = createRandomImage(DataBuffer.TYPE_SHORT);
+        interpolation = Interpolation.BICUBIC;
+        createScaledByTwo(12, 7);
+        verifyAtIntegerPositions();
+    }
+
+    /**
      * Resamples a single-tiled and single-banded image with values 1 everywhere except in
center.
      * The {@linkplain #source} is a 3×3 or 4×4 image with the following values:
      *
@@ -316,4 +338,26 @@ public final strictfp class ResampledImageTest extends TestCase {
             1.111111111f, 1f, 0.888888888f, 0.777777777f, 0.666666666f, 0.777777777f, 0.888888888f,
1f, 1.111111111f
         }, resampleSimpleImage(3), (float) STRICT);
     }
+
+    /**
+     * Checks all values of a bicubic interpolation on a simple image.
+     */
+    @Test
+    public void verifyBicubicResults() {
+        interpolation = Interpolation.BICUBIC;
+        assertArrayEquals(new float[] {
+            0f, 0f, 0f,                0f,                0f,                0f,        
       0f,                0f,                0f,                0f,                  0f, 0f,
+            0f, 1f, 1f,                1f,                1f,                1f,        
       1f,                1f,                1f,                1f,                  1f, 0f,
+            0f, 1f, 1.19753086419753f, 1.34567901234568f, 1.44444444444444f, 1.49382716049383f,
1.49382716049383f, 1.44444444444444f, 1.34567901234568f, 1.19753086419753f,   1f, 0f,
+            0f, 1f, 1.34567901234568f, 1.60493827160494f, 1.77777777777777f, 1.86419753086420f,
1.86419753086420f, 1.77777777777777f, 1.60493827160494f, 1.34567901234568f,   1f, 0f,
+            0f, 1f, 1.44444444444444f, 1.77777777777777f, 2f,                2.11111111111111f,
2.11111111111111f, 2f,                1.77777777777777f, 1.44444444444444f,   1f, 0f,
+            0f, 1f, 1.49382716049383f, 1.86419753086420f, 2.11111111111111f, 2.23456790123457f,
2.23456790123457f, 2.11111111111111f, 1.86419753086420f, 1.49382716049382f,   1f, 0f,
+            0f, 1f, 1.49382716049383f, 1.86419753086420f, 2.11111111111111f, 2.23456790123457f,
2.23456790123457f, 2.11111111111111f, 1.86419753086420f, 1.49382716049382f,   1f, 0f,
+            0f, 1f, 1.44444444444444f, 1.77777777777777f, 2f,                2.11111111111111f,
2.11111111111111f, 2f,                1.77777777777777f, 1.44444444444444f,   1f, 0f,
+            0f, 1f, 1.34567901234568f, 1.60493827160494f, 1.77777777777777f, 1.86419753086420f,
1.86419753086420f, 1.77777777777777f, 1.60493827160494f, 1.34567901234568f,   1f, 0f,
+            0f, 1f, 1.19753086419753f, 1.34567901234568f, 1.44444444444444f, 1.49382716049383f,
1.49382716049383f, 1.44444444444444f, 1.34567901234568f, 1.19753086419753f,   1f, 0f,
+            0f, 1f, 1f,                1f,                1f,                1f,        
       1f,                1f,                1f,                1f,                  1f, 0f,
+            0f, 0f, 0f,                0f,                0f,                0f,        
       0f,                0f,                0f,                0f,                  0f, 0f
+        }, resampleSimpleImage(4), (float) STRICT);
+    }
 }


Mime
View raw message