EFH2HydrologyModel.java [src/java/m/efh2] Revision: ed263028ceff14df667d73424780ec801704c6f3  Date: Tue Sep 09 13:00:29 MDT 2014
package m.efh2;

import m.efh2.EFH2RainfallCoefficients.RfC;
import java.util.List;

/**
 * Storm runoff model based on conventions in Engineering Field Handbook. This
 * uses a lookup method to compute total runoff, based on a number of measures
 * of the drained field.
 */
public class EFH2HydrologyModel {

    /**
     * storm type parameter, from values in EFH2RainfallCoefficients
     */
    private String stormType = "II";
    /**
     * total precipitation parameter, inches
     */
    private double precip = 1.0;
    /**
     * watershed length parameter, feet
     */
    private double watershedLength = 200.0;
    /**
     * watershed slope parameter, percent
     */
    private double watershedSlope = 0.1;
    /**
     * runoff curve number parameter, for lookups
     */
    private int runoffCurveNumber = 40;
    // outputs
    /**
     * runoff depth output, in
     */
    private double runoffQ = Double.NaN;
    /**
     * time of concentration, hours
     */
    private double Tc;
    /**
     * unit peak discharge output, cfs/ac-in
     */
    private double qu;

    public void simulate() {
        // lookup Q based on P and CN, using Table 2.2
        runoffQ = EFH2RunoffTable.lookupQ(precip, runoffCurveNumber);
        computeTc();
        computeQu();
    }

    private void computeTc() {
        // compute Tc using Eq. 2-5
        Tc = Math.pow(watershedLength, 0.8) * Math.pow((1000.0 / runoffCurveNumber) - 9.0, 0.7)
                / (1140.0 * Math.sqrt(watershedSlope));
        if (Tc < 0.1) {
            Tc = 0.1;
        } else if (Tc > 10) {
            throw new IllegalArgumentException(
                    "Computed Tc > 10 hours; EFH2 model is not valid for this range.\n" + "Please estimate peak discharge using NEH-4 procedure, TR-20 program");
        }
    }

    private void computeQu() {
        // lookup rainfall coefficients for qu computation
        double Ia = EFH2IaTable.lookupIa(runoffCurveNumber);
        double SAPG = Ia / precip;    // Ia/P
        List<RfC> rfcs = EFH2RainfallCoefficients.getRainfallCoeff(stormType, SAPG);
        if (rfcs.size() < 1) {
            //Logger logger = Logger.getLogger(this.getClass());
            System.err.println("Failed to lookup rainfall coefficients for storm " + stormType + " and Ia/P = " + SAPG);
            qu = Double.NaN;
            return;
        }

        // compute qu, per algorithm in Basic source
        double logTc = Math.log10(Tc);
        RfC rfc = rfcs.get(0);
        double logLower = rfc.getRfC2() + (rfc.getRfC3() + rfc.getRfC4() * logTc) * logTc;
        double lowerQug = Math.pow(10.0, logLower);
        if (rfcs.size() == 1) {

            // Ia/P matched one of the standard values, use just one set of coeffs
            qu = lowerQug / 640.0;

        } else {
            // Ia/P fell between two standard values, read second RfC coeffs
            double lowerSAPG = rfc.getSAPG();
            rfc = rfcs.get(1);
            double upperSAPG = rfc.getSAPG();

            // compute second Qug result and interpolate
            double logUpper = rfc.getRfC2() + (rfc.getRfC3() + rfc.getRfC4() * logTc) * logTc;
            double upperQug = Math.pow(10.0, logUpper);
            double Qug = lowerQug + (upperQug - lowerQug) * (SAPG - lowerSAPG) / (upperSAPG - lowerSAPG);
            qu = Qug / 640.0;
        }
    }

    /**
     * get storm type parameter, see enumerated values
     */
    public String getStormType() {
        return stormType;
    }

    /**
     * set storm type parameter
     *
     * @param stormType
     */
    public void setStormType(String stormType) {
        this.stormType = stormType;
    }

    public List<String> getStormTypes() {
        return EFH2RainfallCoefficients.getStormTypes();
    }

    /**
     * get total precipitation parameter
     *
     * @return precip in inches
     */
    public double getPrecip() {
        return precip;
    }

    /**
     * set total precipitation parameter
     *
     * @param precip total storm precip in inches
     */
    public void setPrecip(double precip) {
        if (precip < 1.0 || precip > 15) {
            throw new IllegalArgumentException("Precip out of range: (1.0 .. 15.0) : " + precip);
        }
        this.precip = Math.abs(precip);
    }

    /**
     * get watershed length parameter, feet
     */
    public double getWatershedLength() {
        return watershedLength;
    }

    /**
     * set watershed length parameter
     *
     * @param watershedLength feet
     */
    public void setWatershedLength(double watershedLength) {
        if (this.watershedLength < 200 || this.watershedLength > 26000) {
            throw new IllegalArgumentException(
                    "Watershed flow length " + watershedLength
                    + " out of bounds; please use TR-55 program to compute Tc and qu");
        }
        this.watershedLength = Math.abs(watershedLength);
    }

    /**
     * get watershed slope parameter
     *
     * @return slope, ft/ft
     */
    public double getWatershedSlope() {
        return watershedSlope;
    }

    /**
     * set watershed slope parameter
     *
     * @param watershedSlope ft/ft
     */
    public void setWatershedSlope(double watershedSlope) {
        if (Math.abs(watershedSlope) > 64.0 || Math.abs(watershedSlope) < 0.1) {
            throw new IllegalArgumentException(
                    "Watershed slope must be in percent: 0.1% .. 64%");
        }
        this.watershedSlope = Math.abs(watershedSlope);
    }

    /**
     * get runoff curve number parameter
     *
     * @return curve number
     */
    public int getRunoffCurveNumber() {
        return runoffCurveNumber;
    }

    /**
     * set runoff curve number parameter
     *
     * @param runoffCurveNumber curve number
     */
    public void setRunoffCurveNumber(int runoffCurveNumber) {
        if (runoffCurveNumber < 40 || runoffCurveNumber > 95) {
            throw new IllegalArgumentException("Runoff Curve Number: 40 .. 95");
        }
        this.runoffCurveNumber = runoffCurveNumber;
    }

    /**
     * get computed runoff, in inches
     */
    public double getRunoffQ() {
//        simulate();
        return runoffQ;
    }

    /**
     * get computed Tc time of concentration, hours
     */
    public double getTimeOfConcentration() {
//        simulate();
        return Tc;
    }

    /**
     * get computed peack discharge rate, cfs/ac-in??
     */
    public double getUnitPeakDischarge() {
//        simulate();
        return qu;
    }
}