1 /* 2 * Licensed to the Apache Software Foundation (ASF) under one or more 3 * contributor license agreements. See the NOTICE file distributed with 4 * this work for additional information regarding copyright ownership. 5 * The ASF licenses this file to You under the Apache License, Version 2.0 6 * (the "License"); you may not use this file except in compliance with 7 * the License. You may obtain a copy of the License at 8 * 9 * http://www.apache.org/licenses/LICENSE-2.0 10 * 11 * Unless required by applicable law or agreed to in writing, software 12 * distributed under the License is distributed on an "AS IS" BASIS, 13 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. 14 * See the License for the specific language governing permissions and 15 * limitations under the License. 16 */ 17 package org.apache.commons.math3.fitting; 18 19 import java.util.ArrayList; 20 import java.util.List; 21 import java.util.Random; 22 23 import org.apache.commons.math3.analysis.function.HarmonicOscillator; 24 import org.apache.commons.math3.exception.MathIllegalStateException; 25 import org.apache.commons.math3.exception.NumberIsTooSmallException; 26 import org.apache.commons.math3.util.FastMath; 27 import org.apache.commons.math3.util.MathUtils; 28 import org.junit.Assert; 29 import org.junit.Test; 30 31 public class HarmonicCurveFitterTest { 32 /** 33 * Zero points is not enough observed points. 34 */ 35 @Test(expected=NumberIsTooSmallException.class) testPreconditions1()36 public void testPreconditions1() { 37 HarmonicCurveFitter.create().fit(new WeightedObservedPoints().toList()); 38 } 39 40 @Test testNoError()41 public void testNoError() { 42 final double a = 0.2; 43 final double w = 3.4; 44 final double p = 4.1; 45 final HarmonicOscillator f = new HarmonicOscillator(a, w, p); 46 47 final WeightedObservedPoints points = new WeightedObservedPoints(); 48 for (double x = 0.0; x < 1.3; x += 0.01) { 49 points.add(1, x, f.value(x)); 50 } 51 52 final HarmonicCurveFitter fitter = HarmonicCurveFitter.create(); 53 final double[] fitted = fitter.fit(points.toList()); 54 Assert.assertEquals(a, fitted[0], 1.0e-13); 55 Assert.assertEquals(w, fitted[1], 1.0e-13); 56 Assert.assertEquals(p, MathUtils.normalizeAngle(fitted[2], p), 1e-13); 57 58 final HarmonicOscillator ff = new HarmonicOscillator(fitted[0], fitted[1], fitted[2]); 59 for (double x = -1.0; x < 1.0; x += 0.01) { 60 Assert.assertTrue(FastMath.abs(f.value(x) - ff.value(x)) < 1e-13); 61 } 62 } 63 64 @Test 65 public void test1PercentError() { 66 final Random randomizer = new Random(64925784252L); 67 final double a = 0.2; 68 final double w = 3.4; 69 final double p = 4.1; 70 final HarmonicOscillator f = new HarmonicOscillator(a, w, p); 71 72 final WeightedObservedPoints points = new WeightedObservedPoints(); 73 for (double x = 0.0; x < 10.0; x += 0.1) { 74 points.add(1, x, f.value(x) + 0.01 * randomizer.nextGaussian()); 75 } 76 77 final HarmonicCurveFitter fitter = HarmonicCurveFitter.create(); 78 final double[] fitted = fitter.fit(points.toList()); 79 Assert.assertEquals(a, fitted[0], 7.6e-4); 80 Assert.assertEquals(w, fitted[1], 2.7e-3); 81 Assert.assertEquals(p, MathUtils.normalizeAngle(fitted[2], p), 1.3e-2); 82 } 83 84 @Test 85 public void testTinyVariationsData() { 86 final Random randomizer = new Random(64925784252L); 87 88 final WeightedObservedPoints points = new WeightedObservedPoints(); 89 for (double x = 0.0; x < 10.0; x += 0.1) { 90 points.add(1, x, 1e-7 * randomizer.nextGaussian()); 91 } 92 93 final HarmonicCurveFitter fitter = HarmonicCurveFitter.create(); 94 fitter.fit(points.toList()); 95 96 // This test serves to cover the part of the code of "guessAOmega" 97 // when the algorithm using integrals fails. 98 } 99 100 @Test 101 public void testInitialGuess() { 102 final Random randomizer = new Random(45314242L); 103 final double a = 0.2; 104 final double w = 3.4; 105 final double p = 4.1; 106 final HarmonicOscillator f = new HarmonicOscillator(a, w, p); 107 108 final WeightedObservedPoints points = new WeightedObservedPoints(); 109 for (double x = 0.0; x < 10.0; x += 0.1) { 110 points.add(1, x, f.value(x) + 0.01 * randomizer.nextGaussian()); 111 } 112 113 final HarmonicCurveFitter fitter = HarmonicCurveFitter.create() 114 .withStartPoint(new double[] { 0.15, 3.6, 4.5 }); 115 final double[] fitted = fitter.fit(points.toList()); 116 Assert.assertEquals(a, fitted[0], 1.2e-3); 117 Assert.assertEquals(w, fitted[1], 3.3e-3); 118 Assert.assertEquals(p, MathUtils.normalizeAngle(fitted[2], p), 1.7e-2); 119 } 120 121 @Test 122 public void testUnsorted() { 123 Random randomizer = new Random(64925784252L); 124 final double a = 0.2; 125 final double w = 3.4; 126 final double p = 4.1; 127 final HarmonicOscillator f = new HarmonicOscillator(a, w, p); 128 129 // Build a regularly spaced array of measurements. 130 final int size = 100; 131 final double[] xTab = new double[size]; 132 final double[] yTab = new double[size]; 133 for (int i = 0; i < size; i++) { 134 xTab[i] = 0.1 * i; 135 yTab[i] = f.value(xTab[i]) + 0.01 * randomizer.nextGaussian(); 136 } 137 138 // shake it 139 for (int i = 0; i < size; i++) { 140 int i1 = randomizer.nextInt(size); 141 int i2 = randomizer.nextInt(size); 142 double xTmp = xTab[i1]; 143 double yTmp = yTab[i1]; 144 xTab[i1] = xTab[i2]; 145 yTab[i1] = yTab[i2]; 146 xTab[i2] = xTmp; 147 yTab[i2] = yTmp; 148 } 149 150 // Pass it to the fitter. 151 final WeightedObservedPoints points = new WeightedObservedPoints(); 152 for (int i = 0; i < size; ++i) { 153 points.add(1, xTab[i], yTab[i]); 154 } 155 156 final HarmonicCurveFitter fitter = HarmonicCurveFitter.create(); 157 final double[] fitted = fitter.fit(points.toList()); 158 Assert.assertEquals(a, fitted[0], 7.6e-4); 159 Assert.assertEquals(w, fitted[1], 3.5e-3); 160 Assert.assertEquals(p, MathUtils.normalizeAngle(fitted[2], p), 1.5e-2); 161 } 162 163 @Test(expected=MathIllegalStateException.class) 164 public void testMath844() { 165 final double[] y = { 0, 1, 2, 3, 2, 1, 166 0, -1, -2, -3, -2, -1, 167 0, 1, 2, 3, 2, 1, 168 0, -1, -2, -3, -2, -1, 169 0, 1, 2, 3, 2, 1, 0 }; 170 final List<WeightedObservedPoint> points = new ArrayList<WeightedObservedPoint>(); 171 for (int i = 0; i < y.length; i++) { 172 points.add(new WeightedObservedPoint(1, i, y[i])); 173 } 174 175 // The guesser fails because the function is far from an harmonic 176 // function: It is a triangular periodic function with amplitude 3 177 // and period 12, and all sample points are taken at integer abscissae 178 // so function values all belong to the integer subset {-3, -2, -1, 0, 179 // 1, 2, 3}. 180 new HarmonicCurveFitter.ParameterGuesser(points); 181 } 182 } 183