Apache
/
commons-math
mirror of https://github.com/apache/commons-math.git
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

MATH-1014 

Created "HarmonicCurveFitter" as a replacement for "HarmonicFitter".


git-svn-id: https://svn.apache.org/repos/asf/commons/proper/math/trunk@1520807 13f79535-47bb-0310-9956-ffa450edef68
This commit is contained in:
Gilles Sadowski 2013-09-07 20:37:49 +00:00
parent 6429fc8d89
commit f2b1cf5e5d
3 changed files with 625 additions and 0 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

439
src/main/java/org/apache/commons/math3/fitting/HarmonicCurveFitter.java Normal file
View File

@ -0,0 +1,439 @@
/*
* 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.commons.math3.fitting;
import java.util.Collection;
import java.util.List;
import java.util.ArrayList;
import org.apache.commons.math3.analysis.function.HarmonicOscillator;
import org.apache.commons.math3.exception.ZeroException;
import org.apache.commons.math3.exception.NumberIsTooSmallException;
import org.apache.commons.math3.exception.MathIllegalStateException;
import org.apache.commons.math3.exception.util.LocalizedFormats;
import org.apache.commons.math3.fitting.leastsquares.LevenbergMarquardtOptimizer;
import org.apache.commons.math3.fitting.leastsquares.WithStartPoint;
import org.apache.commons.math3.fitting.leastsquares.WithMaxIterations;
import org.apache.commons.math3.linear.DiagonalMatrix;
import org.apache.commons.math3.util.FastMath;
/**
* Fits points to a {@link
* org.apache.commons.math3.analysis.function.HarmonicOscillator.Parametric harmonic oscillator}
* function.
* <br/>
* The {@link #withStartPoint(double[]) initial guess values} must be passed
* in the following order:
* <ul>
* <li>Amplitude</li>
* <li>Angular frequency</li>
* <li>phase</li>
* </ul>
* The optimal values will be returned in the same order.
*
* @version $Id$
* @since 3.3
*/
public class HarmonicCurveFitter extends AbstractCurveFitter<LevenbergMarquardtOptimizer>
implements WithStartPoint<HarmonicCurveFitter>,
WithMaxIterations<HarmonicCurveFitter> {
/** Parametric function to be fitted. */
private static final HarmonicOscillator.Parametric FUNCTION = new HarmonicOscillator.Parametric();
/** Initial guess. */
private final double[] initialGuess;
/** Maximum number of iterations of the optimization algorithm. */
private final int maxIter;
/**
* Contructor used by the factory methods.
*
* @param initialGuess Initial guess. If set to {@code null}, the initial guess
* will be estimated using the {@link ParameterGuesser}.
* @param maxIter Maximum number of iterations of the optimization algorithm.
*/
private HarmonicCurveFitter(double[] initialGuess,
int maxIter) {
this.initialGuess = initialGuess;
this.maxIter = maxIter;
}
/**
* Creates a default curve fitter.
* The initial guess for the parameters will be {@link ParameterGuesser}
* computed automatically, and the maximum number of iterations of the
* optimization algorithm is set to {@link Integer#MAX_VALUE}.
*
* @return a curve fitter.
*
* @see #withStartPoint(double[])
* @see #withMaxIterations(int)
*/
public static HarmonicCurveFitter create() {
return new HarmonicCurveFitter(null, Integer.MAX_VALUE);
}
/** {@inheritDoc} */
public HarmonicCurveFitter withStartPoint(double[] start) {
return new HarmonicCurveFitter(start.clone(),
maxIter);
}
/** {@inheritDoc} */
public HarmonicCurveFitter withMaxIterations(int max) {
return new HarmonicCurveFitter(initialGuess,
max);
}
/** {@inheritDoc} */
@Override
protected LevenbergMarquardtOptimizer getOptimizer(Collection<WeightedObservedPoint> observations) {
// Prepare least-squares problem.
final int len = observations.size();
final double[] target = new double[len];
final double[] weights = new double[len];
int i = 0;
for (WeightedObservedPoint obs : observations) {
target[i] = obs.getY();
weights[i] = obs.getWeight();
++i;
}
final AbstractCurveFitter.TheoreticalValuesFunction model
= new AbstractCurveFitter.TheoreticalValuesFunction(FUNCTION,
observations);
final double[] startPoint = initialGuess != null ?
initialGuess :
// Compute estimation.
new ParameterGuesser(observations).guess();
// Return a new optimizer set up to fit a Gaussian curve to the
// observed points.
return LevenbergMarquardtOptimizer.create()
.withMaxEvaluations(Integer.MAX_VALUE)
.withMaxIterations(maxIter)
.withStartPoint(startPoint)
.withTarget(target)
.withWeight(new DiagonalMatrix(weights))
.withModelAndJacobian(model.getModelFunction(),
model.getModelFunctionJacobian());
}
/**
* This class guesses harmonic coefficients from a sample.
* <p>The algorithm used to guess the coefficients is as follows:</p>
*
* <p>We know f (t) at some sampling points t<sub>i</sub> and want to find a,
* &omega; and &phi; such that f (t) = a cos (&omega; t + &phi;).
* </p>
*
* <p>From the analytical expression, we can compute two primitives :
* <pre>
* If2 (t) = &int; f<sup>2</sup> = a<sup>2</sup> &times; [t + S (t)] / 2
* If'2 (t) = &int; f'<sup>2</sup> = a<sup>2</sup> &omega;<sup>2</sup> &times; [t - S (t)] / 2
* where S (t) = sin (2 (&omega; t + &phi;)) / (2 &omega;)
* </pre>
* </p>
*
* <p>We can remove S between these expressions :
* <pre>
* If'2 (t) = a<sup>2</sup> &omega;<sup>2</sup> t - &omega;<sup>2</sup> If2 (t)
* </pre>
* </p>
*
* <p>The preceding expression shows that If'2 (t) is a linear
* combination of both t and If2 (t): If'2 (t) = A &times; t + B &times; If2 (t)
* </p>
*
* <p>From the primitive, we can deduce the same form for definite
* integrals between t<sub>1</sub> and t<sub>i</sub> for each t<sub>i</sub> :
* <pre>
* If2 (t<sub>i</sub>) - If2 (t<sub>1</sub>) = A &times; (t<sub>i</sub> - t<sub>1</sub>) + B &times; (If2 (t<sub>i</sub>) - If2 (t<sub>1</sub>))
* </pre>
* </p>
*
* <p>We can find the coefficients A and B that best fit the sample
* to this linear expression by computing the definite integrals for
* each sample points.
* </p>
*
* <p>For a bilinear expression z (x<sub>i</sub>, y<sub>i</sub>) = A &times; x<sub>i</sub> + B &times; y<sub>i</sub>, the
* coefficients A and B that minimize a least square criterion
* &sum; (z<sub>i</sub> - z (x<sub>i</sub>, y<sub>i</sub>))<sup>2</sup> are given by these expressions:</p>
* <pre>
*
* &sum;y<sub>i</sub>y<sub>i</sub> &sum;x<sub>i</sub>z<sub>i</sub> - &sum;x<sub>i</sub>y<sub>i</sub> &sum;y<sub>i</sub>z<sub>i</sub>
* A = ------------------------
* &sum;x<sub>i</sub>x<sub>i</sub> &sum;y<sub>i</sub>y<sub>i</sub> - &sum;x<sub>i</sub>y<sub>i</sub> &sum;x<sub>i</sub>y<sub>i</sub>
*
* &sum;x<sub>i</sub>x<sub>i</sub> &sum;y<sub>i</sub>z<sub>i</sub> - &sum;x<sub>i</sub>y<sub>i</sub> &sum;x<sub>i</sub>z<sub>i</sub>
* B = ------------------------
* &sum;x<sub>i</sub>x<sub>i</sub> &sum;y<sub>i</sub>y<sub>i</sub> - &sum;x<sub>i</sub>y<sub>i</sub> &sum;x<sub>i</sub>y<sub>i</sub>
* </pre>
* </p>
*
*
* <p>In fact, we can assume both a and &omega; are positive and
* compute them directly, knowing that A = a<sup>2</sup> &omega;<sup>2</sup> and that
* B = - &omega;<sup>2</sup>. The complete algorithm is therefore:</p>
* <pre>
*
* for each t<sub>i</sub> from t<sub>1</sub> to t<sub>n-1</sub>, compute:
* f (t<sub>i</sub>)
* f' (t<sub>i</sub>) = (f (t<sub>i+1</sub>) - f(t<sub>i-1</sub>)) / (t<sub>i+1</sub> - t<sub>i-1</sub>)
* x<sub>i</sub> = t<sub>i</sub> - t<sub>1</sub>
* y<sub>i</sub> = &int; f<sup>2</sup> from t<sub>1</sub> to t<sub>i</sub>
* z<sub>i</sub> = &int; f'<sup>2</sup> from t<sub>1</sub> to t<sub>i</sub>
* update the sums &sum;x<sub>i</sub>x<sub>i</sub>, &sum;y<sub>i</sub>y<sub>i</sub>, &sum;x<sub>i</sub>y<sub>i</sub>, &sum;x<sub>i</sub>z<sub>i</sub> and &sum;y<sub>i</sub>z<sub>i</sub>
* end for
*
* |--------------------------
* \ | &sum;y<sub>i</sub>y<sub>i</sub> &sum;x<sub>i</sub>z<sub>i</sub> - &sum;x<sub>i</sub>y<sub>i</sub> &sum;y<sub>i</sub>z<sub>i</sub>
* a = \ | ------------------------
* \| &sum;x<sub>i</sub>y<sub>i</sub> &sum;x<sub>i</sub>z<sub>i</sub> - &sum;x<sub>i</sub>x<sub>i</sub> &sum;y<sub>i</sub>z<sub>i</sub>
*
*
* |--------------------------
* \ | &sum;x<sub>i</sub>y<sub>i</sub> &sum;x<sub>i</sub>z<sub>i</sub> - &sum;x<sub>i</sub>x<sub>i</sub> &sum;y<sub>i</sub>z<sub>i</sub>
* &omega; = \ | ------------------------
* \| &sum;x<sub>i</sub>x<sub>i</sub> &sum;y<sub>i</sub>y<sub>i</sub> - &sum;x<sub>i</sub>y<sub>i</sub> &sum;x<sub>i</sub>y<sub>i</sub>
*
* </pre>
* </p>
*
* <p>Once we know &omega;, we can compute:
* <pre>
* fc = &omega; f (t) cos (&omega; t) - f' (t) sin (&omega; t)
* fs = &omega; f (t) sin (&omega; t) + f' (t) cos (&omega; t)
* </pre>
* </p>
*
* <p>It appears that <code>fc = a &omega; cos (&phi;)</code> and
* <code>fs = -a &omega; sin (&phi;)</code>, so we can use these
* expressions to compute &phi;. The best estimate over the sample is
* given by averaging these expressions.
* </p>
*
* <p>Since integrals and means are involved in the preceding
* estimations, these operations run in O(n) time, where n is the
* number of measurements.</p>
*/
public static class ParameterGuesser {
/** Amplitude. */
private final double a;
/** Angular frequency. */
private final double omega;
/** Phase. */
private final double phi;
/**
* Simple constructor.
*
* @param observations Sampled observations.
* @throws NumberIsTooSmallException if the sample is too short.
* @throws ZeroException if the abscissa range is zero.
* @throws MathIllegalStateException when the guessing procedure cannot
* produce sensible results.
*/
public ParameterGuesser(Collection<WeightedObservedPoint> observations) {
if (observations.size() < 4) {
throw new NumberIsTooSmallException(LocalizedFormats.INSUFFICIENT_OBSERVED_POINTS_IN_SAMPLE,
observations.size(), 4, true);
}
final WeightedObservedPoint[] sorted
= sortObservations(observations).toArray(new WeightedObservedPoint[0]);
final double aOmega[] = guessAOmega(sorted);
a = aOmega[0];
omega = aOmega[1];
phi = guessPhi(sorted);
}
/**
* Gets an estimation of the parameters.
*
* @return the guessed parameters, in the following order:
* <ul>
* <li>Amplitude</li>
* <li>Angular frequency</li>
* <li>Phase</li>
* </ul>
*/
public double[] guess() {
return new double[] { a, omega, phi };
}
/**
* Sort the observations with respect to the abscissa.
*
* @param unsorted Input observations.
* @return the input observations, sorted.
*/
private List<WeightedObservedPoint> sortObservations(Collection<WeightedObservedPoint> unsorted) {
final List<WeightedObservedPoint> observations = new ArrayList<WeightedObservedPoint>(unsorted);
// Since the samples are almost always already sorted, this
// method is implemented as an insertion sort that reorders the
// elements in place. Insertion sort is very efficient in this case.
WeightedObservedPoint curr = observations.get(0);
final int len = observations.size();
for (int j = 1; j < len; j++) {
WeightedObservedPoint prec = curr;
curr = observations.get(j);
if (curr.getX() < prec.getX()) {
// the current element should be inserted closer to the beginning
int i = j - 1;
WeightedObservedPoint mI = observations.get(i);
while ((i >= 0) && (curr.getX() < mI.getX())) {
observations.set(i + 1, mI);
if (i-- != 0) {
mI = observations.get(i);
}
}
observations.set(i + 1, curr);
curr = observations.get(j);
}
}
return observations;
}
/**
* Estimate a first guess of the amplitude and angular frequency.
*
* @param observations Observations, sorted w.r.t. abscissa.
* @throws ZeroException if the abscissa range is zero.
* @throws MathIllegalStateException when the guessing procedure cannot
* produce sensible results.
* @return the guessed amplitude (at index 0) and circular frequency
* (at index 1).
*/
private double[] guessAOmega(WeightedObservedPoint[] observations) {
final double[] aOmega = new double[2];
// initialize the sums for the linear model between the two integrals
double sx2 = 0;
double sy2 = 0;
double sxy = 0;
double sxz = 0;
double syz = 0;
double currentX = observations[0].getX();
double currentY = observations[0].getY();
double f2Integral = 0;
double fPrime2Integral = 0;
final double startX = currentX;
for (int i = 1; i < observations.length; ++i) {
// one step forward
final double previousX = currentX;
final double previousY = currentY;
currentX = observations[i].getX();
currentY = observations[i].getY();
// update the integrals of f<sup>2</sup> and f'<sup>2</sup>
// considering a linear model for f (and therefore constant f')
final double dx = currentX - previousX;
final double dy = currentY - previousY;
final double f2StepIntegral =
dx * (previousY * previousY + previousY * currentY + currentY * currentY) / 3;
final double fPrime2StepIntegral = dy * dy / dx;
final double x = currentX - startX;
f2Integral += f2StepIntegral;
fPrime2Integral += fPrime2StepIntegral;
sx2 += x * x;
sy2 += f2Integral * f2Integral;
sxy += x * f2Integral;
sxz += x * fPrime2Integral;
syz += f2Integral * fPrime2Integral;
}
// compute the amplitude and pulsation coefficients
double c1 = sy2 * sxz - sxy * syz;
double c2 = sxy * sxz - sx2 * syz;
double c3 = sx2 * sy2 - sxy * sxy;
if ((c1 / c2 < 0) || (c2 / c3 < 0)) {
final int last = observations.length - 1;
// Range of the observations, assuming that the
// observations are sorted.
final double xRange = observations[last].getX() - observations[0].getX();
if (xRange == 0) {
throw new ZeroException();
}
aOmega[1] = 2 * Math.PI / xRange;
double yMin = Double.POSITIVE_INFINITY;
double yMax = Double.NEGATIVE_INFINITY;
for (int i = 1; i < observations.length; ++i) {
final double y = observations[i].getY();
if (y < yMin) {
yMin = y;
}
if (y > yMax) {
yMax = y;
}
}
aOmega[0] = 0.5 * (yMax - yMin);
} else {
if (c2 == 0) {
// In some ill-conditioned cases (cf. MATH-844), the guesser
// procedure cannot produce sensible results.
throw new MathIllegalStateException(LocalizedFormats.ZERO_DENOMINATOR);
}
aOmega[0] = FastMath.sqrt(c1 / c2);
aOmega[1] = FastMath.sqrt(c2 / c3);
}
return aOmega;
}
/**
* Estimate a first guess of the phase.
*
* @param observations Observations, sorted w.r.t. abscissa.
* @return the guessed phase.
*/
private double guessPhi(WeightedObservedPoint[] observations) {
// initialize the means
double fcMean = 0;
double fsMean = 0;
double currentX = observations[0].getX();
double currentY = observations[0].getY();
for (int i = 1; i < observations.length; ++i) {
// one step forward
final double previousX = currentX;
final double previousY = currentY;
currentX = observations[i].getX();
currentY = observations[i].getY();
final double currentYPrime = (currentY - previousY) / (currentX - previousX);
double omegaX = omega * currentX;
double cosine = FastMath.cos(omegaX);
double sine = FastMath.sin(omegaX);
fcMean += omega * currentY * cosine - currentYPrime * sine;
fsMean += omega * currentY * sine + currentYPrime * cosine;
}
return FastMath.atan2(-fsMean, fcMean);
}
}
}

3
src/main/java/org/apache/commons/math3/fitting/HarmonicFitter.java
View File

@ -35,7 +35,10 @@ import org.apache.commons.math3.util.FastMath;
*
* @version $Id: HarmonicFitter.java 1416643 2012-12-03 19:37:14Z tn $
* @since 2.0
* @deprecated As of 3.3. Please use {@link HarmonicCurveFitter} and
* {@link WeightedObservedPoints} instead.
*/
@Deprecated
public class HarmonicFitter extends CurveFitter<HarmonicOscillator.Parametric> {
/**
* Simple constructor.

183
src/test/java/org/apache/commons/math3/fitting/HarmonicCurveFitterTest.java Normal file
View File

@ -0,0 +1,183 @@
/*
* 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.commons.math3.fitting;
import java.util.Random;
import java.util.List;
import java.util.ArrayList;
import org.apache.commons.math3.optim.nonlinear.vector.jacobian.LevenbergMarquardtOptimizer;
import org.apache.commons.math3.analysis.function.HarmonicOscillator;
import org.apache.commons.math3.exception.NumberIsTooSmallException;
import org.apache.commons.math3.exception.MathIllegalStateException;
import org.apache.commons.math3.util.FastMath;
import org.apache.commons.math3.util.MathUtils;
import org.junit.Test;
import org.junit.Assert;
public class HarmonicCurveFitterTest {
/**
* Zero points is not enough observed points.
*/
@Test(expected=NumberIsTooSmallException.class)
public void testPreconditions1() {
HarmonicCurveFitter.create().fit(new WeightedObservedPoints().toList());
}
@Test
public void testNoError() {
final double a = 0.2;
final double w = 3.4;
final double p = 4.1;
final HarmonicOscillator f = new HarmonicOscillator(a, w, p);
final WeightedObservedPoints points = new WeightedObservedPoints();
for (double x = 0.0; x < 1.3; x += 0.01) {
points.add(1, x, f.value(x));
}
final HarmonicCurveFitter fitter = HarmonicCurveFitter.create();
final double[] fitted = fitter.fit(points.toList());
Assert.assertEquals(a, fitted[0], 1.0e-13);
Assert.assertEquals(w, fitted[1], 1.0e-13);
Assert.assertEquals(p, MathUtils.normalizeAngle(fitted[2], p), 1e-13);
final HarmonicOscillator ff = new HarmonicOscillator(fitted[0], fitted[1], fitted[2]);
for (double x = -1.0; x < 1.0; x += 0.01) {
Assert.assertTrue(FastMath.abs(f.value(x) - ff.value(x)) < 1e-13);
}
}
@Test
public void test1PercentError() {
final Random randomizer = new Random(64925784252L);
final double a = 0.2;
final double w = 3.4;
final double p = 4.1;
final HarmonicOscillator f = new HarmonicOscillator(a, w, p);
final WeightedObservedPoints points = new WeightedObservedPoints();
for (double x = 0.0; x < 10.0; x += 0.1) {
points.add(1, x, f.value(x) + 0.01 * randomizer.nextGaussian());
}
final HarmonicCurveFitter fitter = HarmonicCurveFitter.create();
final double[] fitted = fitter.fit(points.toList());
Assert.assertEquals(a, fitted[0], 7.6e-4);
Assert.assertEquals(w, fitted[1], 2.7e-3);
Assert.assertEquals(p, MathUtils.normalizeAngle(fitted[2], p), 1.3e-2);
}
@Test
public void testTinyVariationsData() {
final Random randomizer = new Random(64925784252L);
final WeightedObservedPoints points = new WeightedObservedPoints();
for (double x = 0.0; x < 10.0; x += 0.1) {
points.add(1, x, 1e-7 * randomizer.nextGaussian());
}
final HarmonicCurveFitter fitter = HarmonicCurveFitter.create();
final double[] fitted = fitter.fit(points.toList());
// This test serves to cover the part of the code of "guessAOmega"
// when the algorithm using integrals fails.
}
@Test
public void testInitialGuess() {
final Random randomizer = new Random(45314242L);
final double a = 0.2;
final double w = 3.4;
final double p = 4.1;
final HarmonicOscillator f = new HarmonicOscillator(a, w, p);
final WeightedObservedPoints points = new WeightedObservedPoints();
for (double x = 0.0; x < 10.0; x += 0.1) {
points.add(1, x, f.value(x) + 0.01 * randomizer.nextGaussian());
}
final HarmonicCurveFitter fitter = HarmonicCurveFitter.create()
.withStartPoint(new double[] { 0.15, 3.6, 4.5 });
final double[] fitted = fitter.fit(points.toList());
Assert.assertEquals(a, fitted[0], 1.2e-3);
Assert.assertEquals(w, fitted[1], 3.3e-3);
Assert.assertEquals(p, MathUtils.normalizeAngle(fitted[2], p), 1.7e-2);
}
@Test
public void testUnsorted() {
Random randomizer = new Random(64925784252L);
final double a = 0.2;
final double w = 3.4;
final double p = 4.1;
final HarmonicOscillator f = new HarmonicOscillator(a, w, p);
// Build a regularly spaced array of measurements.
final int size = 100;
final double[] xTab = new double[size];
final double[] yTab = new double[size];
for (int i = 0; i < size; i++) {
xTab[i] = 0.1 * i;
yTab[i] = f.value(xTab[i]) + 0.01 * randomizer.nextGaussian();
}
// shake it
for (int i = 0; i < size; i++) {
int i1 = randomizer.nextInt(size);
int i2 = randomizer.nextInt(size);
double xTmp = xTab[i1];
double yTmp = yTab[i1];
xTab[i1] = xTab[i2];
yTab[i1] = yTab[i2];
xTab[i2] = xTmp;
yTab[i2] = yTmp;
}
// Pass it to the fitter.
final WeightedObservedPoints points = new WeightedObservedPoints();
for (int i = 0; i < size; ++i) {
points.add(1, xTab[i], yTab[i]);
}
final HarmonicCurveFitter fitter = HarmonicCurveFitter.create();
final double[] fitted = fitter.fit(points.toList());
Assert.assertEquals(a, fitted[0], 7.6e-4);
Assert.assertEquals(w, fitted[1], 3.5e-3);
Assert.assertEquals(p, MathUtils.normalizeAngle(fitted[2], p), 1.5e-2);
}
@Test(expected=MathIllegalStateException.class)
public void testMath844() {
final double[] y = { 0, 1, 2, 3, 2, 1,
0, -1, -2, -3, -2, -1,
0, 1, 2, 3, 2, 1,
0, -1, -2, -3, -2, -1,
0, 1, 2, 3, 2, 1, 0 };
final List<WeightedObservedPoint> points = new ArrayList<WeightedObservedPoint>();
for (int i = 0; i < y.length; i++) {
points.add(new WeightedObservedPoint(1, i, y[i]));
}
// The guesser fails because the function is far from an harmonic
// function: It is a triangular periodic function with amplitude 3
// and period 12, and all sample points are taken at integer abscissae
// so function values all belong to the integer subset {-3, -2, -1, 0,
// 1, 2, 3}.
new HarmonicCurveFitter.ParameterGuesser(points);
}
}