Mercurial > dive4elements > river

comparison artifacts/src/main/java/org/dive4elements/river/artifacts/model/sq/Fitting.java @ 6777:48f6780c372d

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
S/Q relation: More Excel compat.
author Sascha L. Teichmann <teichmann@intevation.de>
date Wed, 07 Aug 2013 19:36:05 +0200
parents 9479cb7c8cd5
children b8f94e865875
comparison
equal deleted inserted replaced
6776:c146fd412f9d 6777:48f6780c372d
7 */ 7 */
8 8
9 package org.dive4elements.river.artifacts.model.sq; 9 package org.dive4elements.river.artifacts.model.sq;
10 10
11 import org.dive4elements.river.artifacts.math.fitting.Function; 11 import org.dive4elements.river.artifacts.math.fitting.Function;
12 import org.dive4elements.river.artifacts.math.fitting.Linear;
13 12
14 import java.util.ArrayList; 13 import java.util.ArrayList;
15 import java.util.List; 14 import java.util.List;
16 15
17 import org.apache.commons.math.MathException; 16 import org.apache.commons.math.MathException;
18 17
19 import org.apache.commons.math.optimization.fitting.CurveFitter; 18 import org.apache.commons.math.optimization.fitting.CurveFitter;
20 19
21 import org.apache.commons.math.optimization.general.AbstractLeastSquaresOptimizer;
22 import org.apache.commons.math.optimization.general.GaussNewtonOptimizer;
23 import org.apache.commons.math.optimization.general.LevenbergMarquardtOptimizer; 20 import org.apache.commons.math.optimization.general.LevenbergMarquardtOptimizer;
21 import org.apache.commons.math.stat.regression.SimpleRegression;
24 22
25 import org.apache.log4j.Logger; 23 import org.apache.log4j.Logger;
26 24
27 public class Fitting 25 public class Fitting
28 implements Outlier.Callback 26 implements Outlier.Callback
80 } 78 }
81 79
82 @Override 80 @Override
83 public void initialize(List<SQ> sqs) throws MathException { 81 public void initialize(List<SQ> sqs) throws MathException {
84 82
85 AbstractLeastSquaresOptimizer optimizer = getOptimizer(); 83 if (USE_NON_LINEAR_FITTING
84 || function.getInitialGuess().length != 2) {
85 nonLinearFitting(sqs);
86 }
87 else {
88 linearFitting(sqs);
89 }
90 }
91
92 protected void linearFitting(List<SQ> sqs) {
93
94 coeffs = linearRegression(sqs);
95
96 instance = function.instantiate(coeffs);
97 }
98
99 protected double [] linearRegression(List<SQ> sqs) {
100
101 SimpleRegression reg = new SimpleRegression();
102
103 int invalidPoints = 0;
104 for (SQ sq: sqs) {
105 double s = sq.getS();
106 double q = sq.getQ();
107 if (s <= 0d || q <= 0d) {
108 ++invalidPoints;
109 continue;
110 }
111 reg.addData(Math.log(q), Math.log(s));
112 }
113
114 if (sqs.size() - invalidPoints < 2) {
115 log.debug("not enough points");
116 return new double [] { 0, 0 };
117 }
118
119 double a = Math.exp(reg.getIntercept());
120 double b = reg.getSlope();
121
122 if (log.isDebugEnabled()) {
123 log.debug("invalid points: " +
124 invalidPoints + " (" + sqs.size() + ")");
125 log.debug("a: " + a + " (" + Math.log(a) + ")");
126 log.debug("b: " + b);
127 }
128
129 return new double [] { a, b };
130 }
131
132
133 protected void nonLinearFitting(List<SQ> sqs) throws MathException {
134
135 LevenbergMarquardtOptimizer optimizer =
136 new LevenbergMarquardtOptimizer();
86 137
87 CurveFitter cf = new CurveFitter(optimizer); 138 CurveFitter cf = new CurveFitter(optimizer);
88 double [] values = new double[2]; 139
89 for (SQ sq: sqs) { 140 for (SQ sq: sqs) {
90 values[0] = sq.getQ(); 141 cf.addObservedPoint(sq.getS(), sq.getQ());
91 values[1] = sq.getS();
92 transformInputValues(values);
93 cf.addObservedPoint(values[0], values[1]);
94 } 142 }
95 143
96 coeffs = cf.fit( 144 coeffs = cf.fit(
97 function, function.getInitialGuess()); 145 function, function.getInitialGuess());
98 146
99 transformCoeffsBack(coeffs);
100
101 instance = function.instantiate(coeffs); 147 instance = function.instantiate(coeffs);
102 148
103 chiSqr = optimizer.getChiSquare(); 149 chiSqr = optimizer.getChiSquare();
104 }
105
106 protected Function getFunction(Function function) {
107 return USE_NON_LINEAR_FITTING
108 ? function
109 : Linear.INSTANCE;
110 }
111
112 protected void transformInputValues(double [] values) {
113 if (!USE_NON_LINEAR_FITTING) {
114 for (int i = 0; i < values.length; ++i) {
115 values[i] = Math.log(values[i]);
116 }
117 }
118 }
119
120 protected AbstractLeastSquaresOptimizer getOptimizer() {
121 return USE_NON_LINEAR_FITTING
122 ? new LevenbergMarquardtOptimizer()
123 : new GaussNewtonOptimizer(false);
124 }
125
126 protected void transformCoeffsBack(double [] coeffs) {
127 if (!USE_NON_LINEAR_FITTING && coeffs.length > 0) {
128 coeffs[0] = Math.exp(coeffs[0]);
129 }
130 } 150 }
131 151
132 @Override 152 @Override
133 public double eval(SQ sq) { 153 public double eval(SQ sq) {
134 double s = instance.value(sq.q); 154 double s = instance.value(sq.q);
http://dive4elements.wald.intevation.org