sascha@3220: package de.intevation.flys.exports.fixings; sascha@3220: sascha@3220: import de.intevation.artifacts.CallMeta; sascha@3220: sascha@3220: import de.intevation.flys.artifacts.math.fitting.Function; sascha@3220: sascha@3220: import de.intevation.flys.artifacts.model.Parameters; sascha@3220: sascha@3220: import de.intevation.flys.artifacts.resources.Resources; sascha@3220: sascha@3220: import de.intevation.flys.exports.ATWriter; sascha@3220: sascha@3220: import java.io.IOException; sascha@3220: import java.io.PrintWriter; sascha@3220: import java.io.Writer; sascha@3220: sascha@3220: import java.util.Locale; sascha@3220: sascha@3220: import org.apache.log4j.Logger; sascha@3220: felix@3905: /** Export Fixation Analysis Results to AT. */ sascha@3220: public class FixATWriter sascha@3220: { felix@3905: /** Private logger. */ sascha@3220: private static Logger log = Logger.getLogger(FixATWriter.class); sascha@3220: sascha@3220: public static final String I18N_HEADER_KEY = sascha@3220: "fix.export.at.header"; sascha@3220: sascha@3220: public static final String I18N_HEADER_DEFAULT = sascha@3220: "Exported fixings discharge curve for {0} {0}-km: {1}"; sascha@3220: sascha@3220: public static final String [] Q_MAX_COLUMN = new String [] { "max_q" }; sascha@3220: sascha@3641: private static final int MAX_ITERATIONS = 10000; sascha@3641: private static final double EPSILON = 1e-8; sascha@3643: private static final double MIN_Q = 1e-4; sascha@3220: sascha@3220: protected Function function; sascha@3220: protected Parameters parameters; sascha@3220: sascha@3220: public FixATWriter() { sascha@3220: } sascha@3220: sascha@3220: public FixATWriter(Function function, Parameters parameters) { sascha@3220: this.function = function; sascha@3220: this.parameters = parameters; sascha@3220: } sascha@3220: sascha@3220: public void write( sascha@3220: Writer writer, sascha@3220: CallMeta meta, sascha@3220: String river, sascha@3220: double km sascha@3220: ) sascha@3220: throws IOException { sascha@3220: PrintWriter out = new PrintWriter(writer); sascha@3220: printHeader(out, meta, river, km); sascha@3220: sascha@3220: double [] coeffs = parameters.interpolate( sascha@3220: "km", km, function.getParameterNames()); sascha@3220: sascha@3220: double [] qMax = parameters.interpolate( sascha@3220: "km", km, Q_MAX_COLUMN); sascha@3220: sascha@3220: if (coeffs == null || qMax == null) { sascha@3220: log.debug("No data found at km " + km + "."); sascha@3220: return; sascha@3220: } sascha@3220: sascha@3220: de.intevation.flys.artifacts.math.Function funcInst = sascha@3220: function.instantiate(coeffs); sascha@3220: sascha@3642: // Increase Q max about 5%. sascha@3642: qMax[0] += Math.abs(qMax[0])*0.05; sascha@3642: sascha@3220: double wMax = funcInst.value(qMax[0]); sascha@3220: sascha@3220: if (Double.isNaN(wMax) || wMax < 0d) { sascha@3220: log.debug("function '" + function.getName() + sascha@3220: "' eval failed at " + wMax); sascha@3220: return; sascha@3220: } sascha@3220: sascha@3220: Function inverse = function.getInverse(); sascha@3220: sascha@3220: de.intevation.flys.artifacts.math.Function invInst = sascha@3220: inverse.instantiate(coeffs); sascha@3220: sascha@3641: double wMin = minW(invInst, wMax, qMax[0]); sascha@3220: sascha@3641: double wMinCM = wMin * 100d; sascha@3641: double wMaxCM = wMax * 100d; sascha@3220: sascha@3641: int wRow = ((int)wMinCM / 10) * 10; sascha@3220: sascha@3641: if ((wMinCM - (int)wMinCM) > 0d) { sascha@3641: wMinCM = (int)wMinCM + 1d; sascha@3641: } sascha@3641: sascha@3641: double w = wMinCM / 100.0; sascha@3641: sascha@3641: int wcm = ((int)wMinCM) % 10; sascha@3641: sascha@3641: if (log.isDebugEnabled()) { sascha@3641: log.debug("wMinCM: " + wMinCM); sascha@3641: log.debug("wMaxCM: " + wMaxCM); sascha@3641: log.debug("wcm: " + wcm); sascha@3641: } sascha@3641: sascha@3641: out.printf(Locale.US, "%8d", wRow); sascha@3641: sascha@3641: if (wcm > 0) { sascha@3641: int rest = 10 - wcm; sascha@3641: while (rest-- > 0) { sascha@3641: out.print(ATWriter.EMPTY); sascha@3220: } sascha@3220: } sascha@3220: sascha@3641: for (;;) { sascha@3641: while (wcm++ < 10) { sascha@3641: if (w > wMax) { sascha@3641: break; sascha@3641: } sascha@3641: double q = invInst.value(w); sascha@3641: if (Double.isNaN(w)) { sascha@3641: out.print(ATWriter.EMPTY); sascha@3641: } sascha@3641: else { sascha@3641: ATWriter.printQ(out, q); sascha@3641: } sascha@3641: w += 0.01d; sascha@3641: } sascha@3641: out.println(); sascha@3641: if (w > wMax) { sascha@3641: break; sascha@3641: } sascha@3641: out.printf(Locale.US, "%8d", wRow += 10); sascha@3641: wcm = 0; sascha@3220: } sascha@3220: sascha@3220: out.flush(); sascha@3220: } sascha@3220: sascha@3220: protected void printHeader( sascha@3220: PrintWriter out, sascha@3220: CallMeta meta, sascha@3220: String river, sascha@3220: double km sascha@3220: ) { sascha@3220: out.println(Resources.format( sascha@3220: meta, sascha@3220: I18N_HEADER_KEY, sascha@3220: I18N_HEADER_DEFAULT, sascha@3220: river, km)); sascha@3220: } sascha@3641: sascha@3641: private static double minW( sascha@3641: de.intevation.flys.artifacts.math.Function function, sascha@3641: double maxW, sascha@3641: double maxQ sascha@3641: ) { sascha@3641: double stepWidth = 10d; sascha@3641: sascha@3641: double lastW = maxW; sascha@3641: double lastQ = maxQ; sascha@3641: sascha@3641: for (int i = 0; i < MAX_ITERATIONS; ++i) { sascha@3641: double w = lastW - stepWidth; sascha@3641: double q = function.value(w); sascha@3641: sascha@3643: if (Double.isNaN(q) || q > lastQ || q < MIN_Q) { sascha@3641: if (stepWidth < EPSILON) { sascha@3641: break; sascha@3641: } sascha@3641: stepWidth *= 0.5d; sascha@3641: continue; sascha@3641: } sascha@3641: sascha@3641: lastW = w; sascha@3641: lastQ = q; sascha@3641: } sascha@3641: sascha@3641: return lastW; sascha@3641: } sascha@3220: } sascha@3220: // vim:set ts=4 sw=4 si et sta sts=4 fenc=utf8 :