Mercurial > dive4elements > river
changeset 6787:51eb6491c537
S/Q: Excel compat completed: Now the data is linearized before fitting. This can be prevented by setting the system property S/Q: Excel compat completed: Now the data is linearized before fitting. This can be prevented by setting the system property "minfo.sq.calcution.non.linear.fitting" to true.
| author | Sascha L. Teichmann <teichmann@intevation.de> |
|---|---|
| date | Thu, 08 Aug 2013 18:14:38 +0200 |
| parents | 70b440dc5317 |
| children | 6587058498f1 |
| files | artifacts/src/main/java/org/dive4elements/river/artifacts/model/sq/Fitting.java artifacts/src/main/java/org/dive4elements/river/artifacts/model/sq/SQRelationCalculation.java |
| diffstat | 2 files changed, 175 insertions(+), 55 deletions(-) [+] |
line wrap: on
line diff
--- a/artifacts/src/main/java/org/dive4elements/river/artifacts/model/sq/Fitting.java Thu Aug 08 16:19:59 2013 +0200 +++ b/artifacts/src/main/java/org/dive4elements/river/artifacts/model/sq/Fitting.java Thu Aug 08 18:14:38 2013 +0200 @@ -8,6 +8,7 @@ package org.dive4elements.river.artifacts.model.sq; +import org.dive4elements.artifacts.common.utils.StringUtils; import org.dive4elements.river.artifacts.math.fitting.Function; import java.util.ArrayList; @@ -83,7 +84,7 @@ public void initialize(List<SQ> sqs) throws MathException { if (USE_NON_LINEAR_FITTING - || function.getInitialGuess().length != 2) { + || function.getParameterNames().length != 2) { nonLinearFitting(sqs); } else { @@ -92,43 +93,48 @@ } protected void linearFitting(List<SQ> sqs) { - - coeffs = linearRegression(sqs); - + coeffs = linearRegression(sqs); instance = function.instantiate(coeffs); } protected double [] linearRegression(List<SQ> sqs) { + String [] pns = function.getParameterNames(); + double [] result = new double[pns.length]; + + if (sqs.size() < 2) { + log.debug("not enough points"); + return result; + } + SimpleRegression reg = new SimpleRegression(); - int invalidPoints = 0; for (SQ sq: sqs) { - double s = sq.getS(); - double q = sq.getQ(); - if (s <= 0d || q <= 0d) { - ++invalidPoints; - continue; - } - reg.addData(Math.log(q), Math.log(s)); + double s = sqView.getS(sq); + double q = sqView.getQ(sq); + reg.addData(q, s); } - if (sqs.size() - invalidPoints < 2) { - log.debug("not enough points"); - return new double [] { 0, 0 }; - } - - double a = Math.exp(reg.getIntercept()); + double m = reg.getIntercept(); double b = reg.getSlope(); if (log.isDebugEnabled()) { - log.debug("invalid points: " + - invalidPoints + " (" + sqs.size() + ")"); - log.debug("a: " + a + " (" + Math.log(a) + ")"); + log.debug("m: " + m); log.debug("b: " + b); } - return new double [] { a, b }; + int mIdx = StringUtils.indexOf("m", pns); + int bIdx = StringUtils.indexOf("b", pns); + + if (mIdx == -1 || bIdx == -1) { + log.error("index not found: " + mIdx + " " + bIdx); + return result; + } + + result[bIdx] = m; + result[mIdx] = b; + + return result; } @@ -140,7 +146,7 @@ CurveFitter cf = new CurveFitter(optimizer); for (SQ sq: sqs) { - cf.addObservedPoint(sq.getS(), sq.getQ()); + cf.addObservedPoint(sqView.getQ(sq), sqView.getS(sq)); } coeffs = cf.fit( @@ -179,7 +185,7 @@ chiSqr); } - public boolean fit(List<SQ> sqs, String method, Callback callback) { + public boolean fit(List<SQ> sqs, String method, Callback callback) { if (sqs.size() < 2) { log.warn("Too less points for fitting.");
--- a/artifacts/src/main/java/org/dive4elements/river/artifacts/model/sq/SQRelationCalculation.java Thu Aug 08 16:19:59 2013 +0200 +++ b/artifacts/src/main/java/org/dive4elements/river/artifacts/model/sq/SQRelationCalculation.java Thu Aug 08 18:14:38 2013 +0200 @@ -8,6 +8,7 @@ package org.dive4elements.river.artifacts.model.sq; +import org.dive4elements.artifacts.common.utils.StringUtils; import org.dive4elements.river.artifacts.access.SQRelationAccess; import org.dive4elements.river.artifacts.math.fitting.Function; @@ -21,6 +22,7 @@ import org.dive4elements.river.backend.SedDBSessionHolder; import java.util.ArrayList; +import java.util.Arrays; import java.util.List; import org.apache.log4j.Logger; @@ -30,7 +32,11 @@ private static final Logger log = Logger.getLogger(SQRelationCalculation.class); - public static final String SQ_FUNCTION_NAME = "sq-pow"; + public static final boolean NON_LINEAR_FITTING = + Boolean.getBoolean("minfo.sq.calcution.non.linear.fitting"); + + public static final String SQ_POW_FUNCTION_NAME = "sq-pow"; + public static final String SQ_LIN_FUNCTION_NAME = "linear"; protected String river; protected double location; @@ -49,8 +55,6 @@ Double outliers = access.getOutliers(); String method = access.getOutlierMethod(); - //river = "Rhein"; - if (river == null) { // TODO: i18n addProblem("sq.missing.river"); @@ -102,20 +106,76 @@ } } + public interface TransformCoeffs { + double [] transform(double [] coeffs); + } + + public static final TransformCoeffs IDENTITY_TRANS = + new TransformCoeffs() { + @Override + public double [] transform(double [] coeffs) { + return coeffs; + } + }; + + public static final TransformCoeffs LINEAR_TRANS = + new TransformCoeffs() { + @Override + public double [] transform(double [] coeffs) { + log.debug("before transform: " + Arrays.toString(coeffs)); + if (coeffs.length == 2) { + coeffs = new double [] { Math.exp(coeffs[1]), coeffs[0] }; + } + log.debug("after transform: " + Arrays.toString(coeffs)); + return coeffs; + } + }; + protected CalculationResult internalCalculate() { - Function function = FunctionFactory + Function powFunction = FunctionFactory .getInstance() - .getFunction(SQ_FUNCTION_NAME); + .getFunction(SQ_POW_FUNCTION_NAME); - if (function == null) { - log.error("No '" + SQ_FUNCTION_NAME + "' function found."); + if (powFunction == null) { + log.error("No '" + SQ_POW_FUNCTION_NAME + "' function found."); // TODO: i18n addProblem("sq.missing.sq.function"); + return new CalculationResult(new SQResult[0], this); } - SQ.View sqView = SQ.SQ_VIEW; - SQ.Factory sqFactory = SQ.SQ_FACTORY; + Function function; + SQ.View sqView; + SQ.Factory sqFactory; + ParameterCreator pc; + + + if (NON_LINEAR_FITTING) { + log.debug("Use non linear fitting."); + sqView = SQ.SQ_VIEW; + sqFactory = SQ.SQ_FACTORY; + function = powFunction; + pc = new ParameterCreator( + powFunction.getParameterNames(), + powFunction.getParameterNames()); + } + else { + log.debug("Use linear fitting."); + sqView = LogSQ.LOG_SQ_VIEW; + sqFactory = LogSQ.LOG_SQ_FACTORY; + function = FunctionFactory + .getInstance() + .getFunction(SQ_LIN_FUNCTION_NAME); + if (function == null) { + log.error("No '" + SQ_LIN_FUNCTION_NAME + "' function found."); + // TODO: i18n + addProblem("sq.missing.sq.function"); + return new CalculationResult(new SQResult[0], this); + } + pc = new LinearParameterCreator( + powFunction.getParameterNames(), + function.getParameterNames()); + } Measurements measurements = MeasurementFactory.getMeasurements( @@ -124,13 +184,14 @@ SQFractionResult [] fractionResults = new SQFractionResult[SQResult.NUMBER_FRACTIONS]; + for (int i = 0; i < fractionResults.length; ++i) { List<SQ> sqs = measurements.getSQs(i); SQFractionResult fractionResult; List<SQFractionResult.Iteration> iterations = - doFitting(function, sqs, sqView); + doFitting(function, sqs, sqView, pc); if (iterations == null) { // TODO: i18n @@ -152,9 +213,10 @@ } protected List<SQFractionResult.Iteration> doFitting( - final Function function, - List<SQ> sqs, - SQ.View sqView + final Function function, + List<SQ> sqs, + SQ.View sqView, + final ParameterCreator pc ) { final List<SQFractionResult.Iteration> iterations = new ArrayList<SQFractionResult.Iteration>(); @@ -171,8 +233,7 @@ double standardDeviation, double chiSqr ) { - Parameters parameters = createParameters( - function.getParameterNames(), + Parameters parameters = pc.createParameters( coeffs, standardDeviation, chiSqr); @@ -186,22 +247,75 @@ return success ? iterations : null; } - public static final Parameters createParameters( - String [] names, - double [] values, - double standardDeviation, - double chiSqr - ) { - String [] columns = new String[names.length + 2]; - columns[0] = "chi_sqr"; - columns[1] = "std_dev"; - System.arraycopy(names, 0, columns, 2, names.length); - Parameters parameters = new Parameters(columns); - int row = parameters.newRow(); - parameters.set(row, names, values); - parameters.set(row, "chi_sqr", chiSqr); - parameters.set(row, "std_dev", standardDeviation); - return parameters; + public static class ParameterCreator { + + protected String [] origNames; + protected String [] proxyNames; + + public ParameterCreator(String [] origNames, String [] proxyNames) { + this.origNames = origNames; + this.proxyNames = proxyNames; + } + + protected double [] transformCoeffs(double [] coeffs) { + return coeffs; + } + + public Parameters createParameters( + double [] coeffs, + double standardDeviation, + double chiSqr + ) { + String [] columns = new String[origNames.length + 2]; + columns[0] = "chi_sqr"; + columns[1] = "std_dev"; + System.arraycopy(origNames, 0, columns, 2, origNames.length); + Parameters parameters = new Parameters(columns); + int row = parameters.newRow(); + parameters.set(row, origNames, transformCoeffs(coeffs)); + parameters.set(row, "chi_sqr", chiSqr); + parameters.set(row, "std_dev", standardDeviation); + return parameters; + } + } + + /** We need to transform the coeffs back to the original function. */ + public static class LinearParameterCreator extends ParameterCreator { + + public LinearParameterCreator( + String [] origNames, + String [] proxyNames + ) { + super(origNames, proxyNames); + } + + @Override + protected double [] transformCoeffs(double [] coeffs) { + + int bP = StringUtils.indexOf("m", proxyNames); + int mP = StringUtils.indexOf("b", proxyNames); + + int aO = StringUtils.indexOf("a", origNames); + int bO = StringUtils.indexOf("b", origNames); + + if (bP == -1 || mP == -1 || aO == -1 || bO == -1) { + log.error("index not found: " + + bP + " " + mP + " " + + aO + " " + bO); + return coeffs; + } + + double [] ncoeffs = (double [])coeffs.clone(); + ncoeffs[aO] = Math.exp(coeffs[mP]); + ncoeffs[bO] = coeffs[bP]; + + if (log.isDebugEnabled()) { + log.debug("before transform: " + Arrays.toString(coeffs)); + log.debug("after transform: " + Arrays.toString(ncoeffs)); + } + + return ncoeffs; + } } } // vim:set ts=4 sw=4 si et sta sts=4 fenc=utf-8 :