# HG changeset patch # User Sascha L. Teichmann # Date 1337706567 0 # Node ID a441be7f15892d521997f4dc1b631e8d466749b1 # Parent 306b9d0f0fb3149e8cbd20560dfef3ba4661c263 Added Fix calculation. flys-artifacts/trunk@4462 c6561f87-3c4e-4783-a992-168aeb5c3f6f diff -r 306b9d0f0fb3 -r a441be7f1589 flys-artifacts/src/main/java/de/intevation/flys/artifacts/model/CalculationResult.java --- a/flys-artifacts/src/main/java/de/intevation/flys/artifacts/model/CalculationResult.java Tue May 22 15:29:54 2012 +0000 +++ b/flys-artifacts/src/main/java/de/intevation/flys/artifacts/model/CalculationResult.java Tue May 22 17:09:27 2012 +0000 @@ -14,6 +14,10 @@ public CalculationResult() { } + public CalculationResult(Calculation report) { + this(null, report); + } + /** * @param report report (e.g. error messages). */ diff -r 306b9d0f0fb3 -r a441be7f1589 flys-artifacts/src/main/java/de/intevation/flys/artifacts/model/FixingsColumn.java --- a/flys-artifacts/src/main/java/de/intevation/flys/artifacts/model/FixingsColumn.java Tue May 22 15:29:54 2012 +0000 +++ b/flys-artifacts/src/main/java/de/intevation/flys/artifacts/model/FixingsColumn.java Tue May 22 17:09:27 2012 +0000 @@ -28,22 +28,26 @@ } public boolean getW(double km, double [] w) { + return getW(km, w, 0); + } + + public boolean getW(double km, double [] w, int index) { if (kms.length == 0 || km < kms[0] || km > kms[kms.length-1]) { - w[0] = Double.NaN; + w[index] = Double.NaN; return true; } int idx = Arrays.binarySearch(kms, km); if (idx >= 0) { - w[0] = ws[idx]; + w[index] = ws[idx]; return true; } idx = -idx - 1; - w[0] = Linear.linear(km, kms[idx], kms[idx+1], ws[idx], ws[idx+1]); + w[index] = Linear.linear(km, kms[idx], kms[idx+1], ws[idx], ws[idx+1]); return false; } diff -r 306b9d0f0fb3 -r a441be7f1589 flys-artifacts/src/main/java/de/intevation/flys/artifacts/model/FixingsColumnFactory.java --- a/flys-artifacts/src/main/java/de/intevation/flys/artifacts/model/FixingsColumnFactory.java Tue May 22 15:29:54 2012 +0000 +++ b/flys-artifacts/src/main/java/de/intevation/flys/artifacts/model/FixingsColumnFactory.java Tue May 22 17:09:27 2012 +0000 @@ -43,6 +43,10 @@ private FixingsColumnFactory() { } + public static FixingsColumnFactory getInstance() { + return INSTANCE; + } + public FixingsColumn getColumnData(Fixing.Column column) { boolean debug = log.isDebugEnabled(); diff -r 306b9d0f0fb3 -r a441be7f1589 flys-artifacts/src/main/java/de/intevation/flys/artifacts/model/Parameters.java --- a/flys-artifacts/src/main/java/de/intevation/flys/artifacts/model/Parameters.java Tue May 22 15:29:54 2012 +0000 +++ b/flys-artifacts/src/main/java/de/intevation/flys/artifacts/model/Parameters.java Tue May 22 17:09:27 2012 +0000 @@ -2,43 +2,11 @@ import gnu.trove.TDoubleArrayList; -import java.util.Map; - import java.io.Serializable; public class Parameters implements Serializable { - public static class Parameter - implements Serializable - { - public String name; - public double value; - - public Parameter() { - } - - public Parameter(String name, double value) { - this.name = name; - this.value = value; - } - public String getName() { - return name; - } - - public void setName(String name) { - this.name = name; - } - - public double getValue() { - return value; - } - - public void setValue(double value) { - this.value = value; - } - } // Parameter - protected String [] columnNames; protected TDoubleArrayList [] columns; @@ -65,23 +33,7 @@ return -1; } - public void add(Parameter [] parameters) { - - int N = columns[0].size(); - - for (int i = 0; i < columns.length; ++i) { - columns[i].add(Double.NaN); - } - - for (Parameter parameter: parameters) { - int index = columnIndex(parameter.getName()); - if (index >= 0) { - columns[index].setQuick(N, parameter.getValue()); - } - } - } - - public void add(Map parameters) { + public int newRow() { int N = columns[0].size(); @@ -89,31 +41,7 @@ columns[i].add(Double.NaN); } - for (Map.Entry entry: parameters.entrySet()) { - int index = columnIndex(entry.getKey()); - Double v = entry.getValue(); - if (index >= 0 && v != null) { - columns[index].setQuick(N, v); - } - } - } - - public Parameter [] get(int i) { - Parameter [] parameters = new Parameter[columns.length]; - for (int j = 0; i < parameters.length; ++j) { - parameters[j] = new Parameter(columnNames[j], Double.NaN); - } - return get(i, parameters); - } - - public Parameter [] get(int i, Parameter [] parameters) { - for (Parameter parameter: parameters) { - int index = columnIndex(parameter.getName()); - if (index >= 0) { - parameter.setValue(columns[index].getQuick(i)); - } - } - return parameters; + return N; } public double get(int i, int index) { @@ -127,6 +55,13 @@ : Double.NaN; } + public void set(int i, String columnName, double value) { + int idx = columnIndex(columnName); + if (idx >= 0) { + columns[i].setQuick(idx, value); + } + } + public int size() { return columns[0].size(); } diff -r 306b9d0f0fb3 -r a441be7f1589 flys-artifacts/src/main/java/de/intevation/flys/artifacts/model/fixings/FixCalculation.java --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/flys-artifacts/src/main/java/de/intevation/flys/artifacts/model/fixings/FixCalculation.java Tue May 22 17:09:27 2012 +0000 @@ -0,0 +1,252 @@ +package de.intevation.flys.artifacts.model.fixings; + +import de.intevation.flys.artifacts.FixationArtifactAccess; + +import de.intevation.flys.artifacts.math.fitting.Function; +import de.intevation.flys.artifacts.math.fitting.FunctionFactory; + +import de.intevation.flys.artifacts.model.Calculation; +import de.intevation.flys.artifacts.model.CalculationResult; +import de.intevation.flys.artifacts.model.FixingsColumn; +import de.intevation.flys.artifacts.model.FixingsColumnFactory; +import de.intevation.flys.artifacts.model.FixingsOverview; +import de.intevation.flys.artifacts.model.FixingsOverviewFactory; +import de.intevation.flys.artifacts.model.Parameters; + +import de.intevation.flys.artifacts.model.FixingsOverview.Fixing; +import de.intevation.flys.artifacts.model.FixingsOverview.IdFilter; + +import de.intevation.flys.utils.DoubleUtil; + +import java.util.ArrayList; +import java.util.List; + +import org.apache.commons.math.MathException; + +import org.apache.commons.math.optimization.fitting.CurveFitter; + +import org.apache.commons.math.optimization.general.LevenbergMarquardtOptimizer; + +import org.apache.log4j.Logger; + +public class FixCalculation +extends Calculation +{ + private static Logger logger = Logger.getLogger(FixCalculation.class); + + protected String river; + protected double from; + protected double to; + protected double step; + protected boolean preprocessing; + protected String function; + protected int [] events; + + public FixCalculation() { + } + + public FixCalculation(FixationArtifactAccess access) { + + String river = access.getRiver(); + Double from = access.getFrom(); + Double to = access.getTo(); + Double step = access.getStep(); + String function = access.getFunction(); + int [] events = access.getEvents(); + + if (river == null) { + // TODO: i18n + addProblem("fix.missing.river"); + } + + if (from == null) { + // TODO: i18n + addProblem("fix.missing.from"); + } + + if (to == null) { + // TODO: i18n + addProblem("fix.missing.to"); + } + + if (step == null) { + // TODO: i18n + addProblem("fix.missing.step"); + } + + if (function == null) { + // TODO: i18n + addProblem("fix.missing.function"); + } + + if (events == null || events.length < 1) { + // TODO: i18n + addProblem("fix.missing.events"); + } + + if (!hasProblems()) { + this.river = river; + this.from = from; + this.to = to; + this.step = step; + this.function = function; + this.events = events; + } + } + + public CalculationResult calculate() { + + FixingsOverview overview = + FixingsOverviewFactory.getOverview(river); + + if (overview == null) { + addProblem("fix.no.overview.available"); + } + + Function func = FunctionFactory.getInstance() + .getFunction(function); + + if (func == null) { + // TODO: i18n + addProblem("fix.invalid.function.name"); + } + + if (hasProblems()) { + return new CalculationResult(this); + } + + FixingsColumnFactory fcf = FixingsColumnFactory.getInstance(); + + List dataColumns = + new ArrayList(events.length); + + for (int eventId: events) { + IdFilter idFilter = new IdFilter(eventId); + + List columns = overview.filter(null, idFilter); + if (columns.isEmpty()) { + // TODO: i18n + addProblem("fix.missing.column", eventId); + continue; + } + FixingsColumn dataColumn = fcf.getColumnData(columns.get(0)); + if (dataColumn == null) { + // TODO: i18n + addProblem("fix.cannot.load.data", eventId); + continue; + } + dataColumns.add(dataColumn); + } + + if (dataColumns.size() < 2) { + // TODO: i18n + addProblem("fix.too.less.data.columns"); + return new CalculationResult(this); + } + + double [] kms = DoubleUtil.explode(from, to, step); + + double [] ws = new double[dataColumns.size()]; + double [] qs = new double[ws.length]; + + String [] parameterNames = func.getParameterNames(); + + Parameters results = + new Parameters(createColumnNames(parameterNames)); + + boolean invalid = false; + + for (int i = 0; i < kms.length; ++i) { + double km = kms[i]; + + for (int j = 0; j < ws.length; ++j) { + FixingsColumn column = dataColumns.get(j); + qs[j] = column.getQ(km); + boolean interpolated = column.getW(km, ws, j); + // TODO: mark as interpolated. + } + + // TODO: Do preprocessing here! + double [] parameters = fit(func, km, ws, qs); + if (parameters == null) { // Problems are reported already. + continue; + } + + int row = results.newRow(); + + results.set(row, "km", km); + for (int j = 0; j < parameters.length; ++j) { + if (Double.isNaN(parameters[j])) { + invalid = true; + } + else { + results.set(row, parameterNames[j], parameters[j]); + } + } + // TODO: Calculate statistics, too! + } + + if (invalid) { + // TODO: i18n + addProblem("fix.invalid.values"); + results.removeNaNs(); + } + + // TODO: Calculate Delta W/t, too. + return new CalculationResult(results, this); + } + + protected static String [] createColumnNames(String [] parameters) { + String [] result = new String[parameters.length + 1]; + result[0] = "km"; + // TODO: Add statistic columns, too. + System.arraycopy(parameters, 0, result, 1, parameters.length); + return result; + } + + protected double [] fit( + Function function, + double km, + double [] ws, + double [] qs + ) { + LevenbergMarquardtOptimizer lmo = new LevenbergMarquardtOptimizer(); + CurveFitter cf = new CurveFitter(lmo); + + boolean missingWs = false; + boolean missingQs = false; + + for (int i = 0; i < ws.length; ++i) { + boolean ignore = false; + if (Double.isNaN(ws[i])) { + ignore = true; + if (!missingWs) { + missingWs = true; + // TODO: i18n + addProblem(km, "fix.missing.w"); + } + } + if (Double.isNaN(qs[i])) { + ignore = true; + if (!missingQs) { + missingQs = true; + // TODO: i18n + addProblem(km, "fix.missing.q"); + } + } + if (!ignore) { + cf.addObservedPoint(ws[i], qs[i]); + } + } + + try { + return cf.fit(function, function.getInitialGuess()); + } + catch (MathException me) { + addProblem(km, "fix.fitting.failed"); + } + + return null; + } +} +// vim:set ts=4 sw=4 si et sta sts=4 fenc=utf8 :