[MATH-851] Added method MathArrays.convolve, thanks to Clemens Novak for the patch.

git-svn-id: https://svn.apache.org/repos/asf/commons/proper/math/trunk@1486982 13f79535-47bb-0310-9956-ffa450edef68
This commit is contained in:
Thomas Neidhart 2013-05-28 16:04:32 +00:00
parent 07f6ec88e9
commit 4fd920c291
3 changed files with 113 additions and 0 deletions

View File

@ -51,6 +51,10 @@ If the output is not quite correct, check for invisible trailing spaces!
</properties>
<body>
<release version="x.y" date="TBD" description="TBD">
<action dev="tn" type="add" issue="MATH-851" due-to="Clemens Novak">
Added method "MathArrays#convolve(double[], double[])" to compute the
discrete, linear convolution of two sequences.
</action>
<action dev="tn" type="add" issue="MATH-977">
Added low-discrepancy random generator "HaltonSequenceGenerator".
</action>

View File

@ -29,6 +29,7 @@ import org.apache.commons.math3.exception.DimensionMismatchException;
import org.apache.commons.math3.exception.MathArithmeticException;
import org.apache.commons.math3.exception.MathIllegalArgumentException;
import org.apache.commons.math3.exception.MathInternalError;
import org.apache.commons.math3.exception.NoDataException;
import org.apache.commons.math3.exception.NonMonotonicSequenceException;
import org.apache.commons.math3.exception.NotPositiveException;
import org.apache.commons.math3.exception.NotStrictlyPositiveException;
@ -1352,4 +1353,51 @@ public class MathArrays {
return array;
}
/**
* Calculates the convolution between two sequences.
* <p>
* The solution is obtained via straightforward computation of the convolution sum (and not via FFT; for longer sequences,
* the performance of this method might be inferior to an FFT-based implementation).
*
* @param x the first sequence (double array of length {@code N}); the sequence is assumed to be zero elsewhere
* (i.e. {x[i]}=0 for i<0 and i>={@code N}). Typically, this sequence will represent an input signal to a system.
* @param h the second sequence (double array of length {@code M}); the sequence is assumed to be zero elsewhere
* (i.e. {h[i]}=0 for i<0 and i>={@code M}). Typically, this sequence will represent the impulse response of the system.
* @return the convolution of {@code x} and {@code h} (double array of length {@code N} + {@code M} -1)
* @throws NullArgumentException if either {@code x} or {@code h} is null
* @throws NoDataException if either {@code x} or {@code h} is empty
*
* @see <a href="http://en.wikipedia.org/wiki/Convolution">Convolution (Wikipedia)</a>
* @since 4.0
*/
public static double[] convolve(double[] x, double[] h) throws NullArgumentException, NoDataException {
MathUtils.checkNotNull(x);
MathUtils.checkNotNull(h);
final int N = x.length;
final int M = h.length;
if (N == 0 || M == 0) {
throw new NoDataException();
}
// initialize the output array
final int totalLength = N + M - 1;
final double[] y = new double[totalLength];
// straightforward implementation of the convolution sum
for (int n = 0; n < totalLength; n++) {
double yn = 0;
for (int k = 0; k < M; k++) {
final int j = n - k;
if ((j > -1) && (j < N) ) {
yn = yn + x[j] * h[k];
}
}
y[n] = yn;
}
return y;
}
}

View File

@ -13,12 +13,15 @@
*/
package org.apache.commons.math3.util;
import static org.junit.Assert.fail;
import java.util.Arrays;
import org.apache.commons.math3.TestUtils;
import org.apache.commons.math3.exception.DimensionMismatchException;
import org.apache.commons.math3.exception.MathArithmeticException;
import org.apache.commons.math3.exception.MathIllegalArgumentException;
import org.apache.commons.math3.exception.NoDataException;
import org.apache.commons.math3.exception.NonMonotonicSequenceException;
import org.apache.commons.math3.exception.NotPositiveException;
import org.apache.commons.math3.exception.NotStrictlyPositiveException;
@ -833,4 +836,62 @@ public class MathArraysTest {
Assert.fail("expecting MathIllegalArgumentException");
} catch (MathIllegalArgumentException ex) {}
}
@Test
public void testConvolve() {
/* Test Case (obtained via SciPy)
* x=[1.2,-1.8,1.4]
* h=[1,0.8,0.5,0.3]
* convolve(x,h) -> array([ 1.2 , -0.84, 0.56, 0.58, 0.16, 0.42])
*/
double[] x1 = { 1.2, -1.8, 1.4 };
double[] h1 = { 1, 0.8, 0.5, 0.3 };
double[] y1 = { 1.2, -0.84, 0.56, 0.58, 0.16, 0.42 };
double tolerance = 1e-13;
double[] yActual = MathArrays.convolve(x1, h1);
Assert.assertArrayEquals(y1, yActual, tolerance);
double[] x2 = { 1, 2, 3 };
double[] h2 = { 0, 1, 0.5 };
double[] y2 = { 0, 1, 2.5, 4, 1.5 };
yActual = MathArrays.convolve(x2, h2);
Assert.assertArrayEquals(y2, yActual, tolerance);
try {
MathArrays.convolve(new double[]{1, 2}, null);
fail("an exception should have been thrown");
} catch (NullArgumentException e) {
// expected behavior
}
try {
MathArrays.convolve(null, new double[]{1, 2});
fail("an exception should have been thrown");
} catch (NullArgumentException e) {
// expected behavior
}
try {
MathArrays.convolve(new double[]{1, 2}, new double[]{});
fail("an exception should have been thrown");
} catch (NoDataException e) {
// expected behavior
}
try {
MathArrays.convolve(new double[]{}, new double[]{1, 2});
fail("an exception should have been thrown");
} catch (NoDataException e) {
// expected behavior
}
try {
MathArrays.convolve(new double[]{}, new double[]{});
fail("an exception should have been thrown");
} catch (NoDataException e) {
// expected behavior
}
}
}