DFLOW_lowFlowStats.java [src/java/m/cfa/timeseries] Revision: e92755c243161b026a60fe7fb06b32b04402a8dc  Date: Thu Jan 25 10:37:07 MST 2018
package m.cfa.timeseries;

import csip.Executable;
import java.io.BufferedReader;
import java.io.File;
import java.io.FileNotFoundException;
import java.io.FileReader;
import java.io.FileWriter;
import m.cfa.DoubleArray;
import m.cfa.DoubleMath;
import java.io.IOException;
import java.io.PrintWriter;
import java.text.ParseException;
import java.util.ArrayList;

/**
* Last Updated: 5-June-2017
* @author Tyler Wible
* @since 17-May-2017
*/
public class DFLOW_lowFlowStats {
    /**
     * Calculate the DFLOW "extreme value" design flow (see DFLOW user manual) 
     * @param e the DFLOW executable
     * @param directory  the working directory for the executable
     * @param flowData  flow data, column1 = dates (format yyyy-mm-dd), column2 = flow values
     * @param M  the number of days to average for annual low flow analysis (m-day average)
     * @param R  the desired return period of the m-day low flow (in years)
     * @return
     * @throws java.io.IOException
     */
    public double DFLOW_ExtremeValue(Executable e,
                                     String directory,
                                     String[][] flowData,
                                     int M,
                                     double R) throws IOException{
        double designFlow = -1;
        if(flowData.length == 0){
            return designFlow;
        }
        
        //Check inputs
        if(M < 1 || M > 30){
            throw new IOException("Error encountered. Please see the following message for details: \n"
                    + "Input out of allowable range\n Current averaging length (M): " + M + "\nAllowable range: 1 to 30");
        }
        if(R < 1 || R > 500){
            throw new IOException("Error encountered. Please see the following message for details: \n"
                    + "Input out of allowable range\n Current return period (R): " + R + "\nAllowable range: 1 to 500");
        }
        
        //Write inputs
        writeDFLOWinp(directory, 2, 2, M, R, 120, 5);
        writeDFLOWflo(directory, flowData);
        
        //Call DFLOW
        e.redirectOutput("dflow.out");
        e.exec();
        
        //Expected Output: "dflow.out"
        if (!new File(directory, "dflow.out").exists()) {
            throw new FileNotFoundException("dflow.out");
        }
        
        //Parse DFLOW output
        FileReader file_to_read = new FileReader(directory + File.separator + "dflow.out");
        BufferedReader reader = new BufferedReader(file_to_read);
        String currentLine;
        while((currentLine = reader.readLine()) !=null){
            if(currentLine.contains("DESIGN FLOWS FOR USGS GAGE")){
                currentLine = reader.readLine();//skip header
                currentLine = reader.readLine();//skip period of record summary
                currentLine = reader.readLine();//design flow
                String flow_str = currentLine.substring(currentLine.indexOf(":") + 1, currentLine.indexOf("CFS")).trim();
                designFlow = Double.parseDouble(flow_str);
                break;
            }
        }
        reader.close();
        file_to_read.close();
        
        //Round design flows
        designFlow = DoubleMath.round(designFlow, 2);
        
        return designFlow;
    }
    /**
     * Calculate the CDPHE "biologically-based" design flow (see DFLOW user manual) 
     * which is an 'm'-day harmonic average low flow based on a certain excursion count (exceedance?), 
     * performs this calculation for the entire flow record as well as each month of the year
     * @param e the DFLOW executable
     * @param directory  the working directory for the executable
     * @param flowData  a String[][] of flow data, column1 = dates (format yyyy-mm-dd), column2 = flow values
     * @param M  the number of days to average for annual low flow analysis (m-day average)
     * @param R  the desired return period of the m-day low flow (in years)
     * @param clusterLength  the length of time for the definition of a 'cluster' of excursion 
     * periods (default should be 120). All excursion periods within this amount of 
     * time of each other will be counted, subject to the clusterCountMax below.
     * @param clusterCountMax  the upper count limit on how many excursion periods 
     * will be counted within an excursion cluster (default should be 5)
     * @return
     * @throws IOException
     * @throws ParseException 
     */
    public double[][] DFLOW_Biological(Executable e,
                                     String directory,
                                     String[][] flowData,
                                     int M,
                                     int R,
                                     int clusterLength,
                                     int clusterCountMax) throws IOException, ParseException{
        double[] designFlows = {-1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1};
        double[] designYears = {-1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1};
        if(flowData.length == 0){
            double[][] returnArray = {designFlows, designYears};
            return returnArray;
        }
        
        //Check inputs
        if(M < 1 || M > 30){
            throw new IOException("Error encountered. Please see the following message for details: \n"
                    + "Input out of allowable range\nCurrent averaging length (M): " + M + "\nAllowable range: 1 to 30");
        }
        if(R < 1 || R > 500){
            throw new IOException("Error encountered. Please see the following message for details: \n"
                    + "Input out of allowable range\nCurrent return period (R): " + R + "\nAllowable range: 1 to 500");
        }
        if(clusterLength < 1){
            throw new IOException("Error encountered. Please see the following message for details: \n"
                    + "Input out of allowable range\nLength of excursion clustering period (XMI): " + clusterLength + "\nAllowable range: greater than 1");
        }
        if(clusterCountMax < 1){
            throw new IOException("Error encountered. Please see the following message for details: \n"
                    + "Input out of allowable range\nMaximum number of excursions (XME): " + clusterCountMax + "\nAllowable range: greater than 1");
        }
        int biologicalType = 2;//i.e. M and R are custom
        if(M == 30 && R == 3){
            biologicalType = 3;//i.e. M and R are 30Q3
        }
        
        //Write inputs
        writeDFLOWinp(directory, 1, biologicalType, M, R, clusterLength, clusterCountMax);
        writeDFLOWflo(directory, flowData);
        
        //Call DFLOW
        e.redirectOutput("dflow.out");
        e.exec();
        
        //Expected Output: "dflow.out"
        if (!new File(directory, "dflow.out").exists()) {
            throw new FileNotFoundException("dflow.out");
        }
        
        //Parse DFLOW output
        FileReader file_to_read = new FileReader(directory + File.separator + "dflow.out");
        BufferedReader reader = new BufferedReader(file_to_read);
        String currentLine;
        while((currentLine = reader.readLine()) !=null){
            if(currentLine.contains("DESIGN FLOWS FOR USGS GAGE")){
                //Parse annual flow
                currentLine = reader.readLine();//skip header
                currentLine = reader.readLine();//skip period of record summary
                currentLine = reader.readLine();//skip number of excursions
                currentLine = reader.readLine();//skip length of flow averaging period (M)
                currentLine = reader.readLine();//skip average interval between excursions
                currentLine = reader.readLine();//design flow
                
                String flow_str = currentLine.substring(currentLine.indexOf(":") + 1, currentLine.indexOf("CFS")).trim();
                designFlows[0] = Double.parseDouble(flow_str);
                
                //Parse monthly flows
                currentLine = reader.readLine();//skip blank line
                for(int i=1; i<=12; i++){
                    currentLine = reader.readLine();
                    String flow_tmp = currentLine.substring(currentLine.indexOf(":") + 1, currentLine.indexOf("cfs")).trim();
                    designFlows[i] = Double.parseDouble(flow_tmp);
                    
                    String year = currentLine.substring(currentLine.indexOf(":")-4,currentLine.indexOf(":"));
                    designYears[i] = Double.parseDouble(year);
                }
                break;
            }
        }
        reader.close();
        file_to_read.close();
        
        //Round design flows
        designFlows = DoubleMath.roundColumn(designFlows, 2);
        double[][] returnArray = {designFlows, designYears};
        return returnArray;
    }
    /**
     * Calculate the CDPHE "human health" design flow (see DFLOW user manual) 
     * which is the harmonic mean of the flows
     * @param e the DFLOW executable
     * @param directory  the working directory for the executable
     * @param flowData  flow data, column1 = dates (format yyyy-mm-dd), column2 = flow values
     * @return
     * @throws IOException
     * @throws ParseException 
     */
    public double DFLOW_HumanHealth(Executable e, String directory, String[][] flowData) throws IOException, ParseException{
        double designFlow = -1;
        if(flowData.length == 0){
            return designFlow;
        }
        
        //Write inputs
        writeDFLOWinp(directory, 3, 1, 1, 1, 120, 5);
        writeDFLOWflo(directory, flowData);
        
        //Call DFLOW
        e.redirectOutput("dflow.out");
        e.exec();
        
        //Expected Output: "dflow.out"
        if (!new File(directory, "dflow.out").exists()) {
            throw new FileNotFoundException("dflow.out");
        }
        
        //Parse DFLOW output
        FileReader file_to_read = new FileReader(directory + File.separator + "dflow.out");
        BufferedReader reader = new BufferedReader(file_to_read);
        String currentLine;
        while((currentLine = reader.readLine()) !=null){
            if(currentLine.contains("DESIGN FLOWS FOR USGS GAGE")){
                currentLine = reader.readLine();//skip header
                currentLine = reader.readLine();//skip period of record summary
                currentLine = reader.readLine();//design flow
                String flow_str = currentLine.substring(currentLine.indexOf(":") + 1, currentLine.indexOf("CFS")).trim();
                designFlow = Double.parseDouble(flow_str);
                break;
            }
        }
        reader.close();
        file_to_read.close();
        
        //Round design flows
        designFlow = DoubleMath.round(designFlow, 2);
        
        return designFlow;
    }
    /**
     * @param e the DFLOW executable
     * @param directory  the working directory for the executable
     * @param flowData  a String[][] of flow data, column1 = dates (format yyyy-mm-dd), column2 = flow values
     * @return
     * @throws IOException
     * @throws ParseException 
     */
    public String DFLOW_REG31(Executable e, String directory, String[][] flowData) throws IOException, ParseException{
        ArrayList<Integer> monthlyRowIndex = new ArrayList<>();
        monthlyRowIndex.add(1);
        monthlyRowIndex.add(2);
        monthlyRowIndex.add(3);
        monthlyRowIndex.add(4);
        monthlyRowIndex.add(5);
        monthlyRowIndex.add(6);
        monthlyRowIndex.add(7);
        monthlyRowIndex.add(8);
        monthlyRowIndex.add(9);
        monthlyRowIndex.add(10);
        monthlyRowIndex.add(11);
        monthlyRowIndex.add(12);
        
        //Calculate acute biologically based low flows (1-day 3-year)
        System.out.println("Calcualting Reg. 31 1E3 Acute Monthly Low Flows...");
        double[][] monthlyBiological_1day = DFLOW_Biological(e, directory, flowData, 1, 3, 120, 5);
        double monlyBiological_1day_min = DoubleMath.min(DoubleArray.getRows(DoubleArray.getColumn(monthlyBiological_1day, 0), monthlyRowIndex));
        
        //Calculate chronic biologically based low flows (7-day 3-year)
        System.out.println("Calcualting Reg. 31 7E3 Chronic Monthly Low Flows...");
        double[][] monthlyBiological_7day = DFLOW_Biological(e, directory, flowData, 7, 3, 120, 5);
        double monlyBiological_7day_min = DoubleMath.min(DoubleArray.getRows(DoubleArray.getColumn(monthlyBiological_7day, 0), monthlyRowIndex));
        
        //Calculate chronic biologically based low flows (30-day 3-year)
        System.out.println("Calcualting Reg. 31 30E3 Chronic Monthly Low Flows...");
        double[][] monthlyBiological_30day = DFLOW_Biological(e, directory, flowData, 30, 3, 120, 5);
        double monlyBiological_30day_min = DoubleMath.min(DoubleArray.getRows(DoubleArray.getColumn(monthlyBiological_30day,0), monthlyRowIndex));
        
        //Calculate annual median of data with a 5-year return period
        System.out.println("Calcualting Reg. 31 1E5 Annual Median of Daily Average Flows...");
        CDPHE_lowFlowStats cdphe_lowFlowStats = new CDPHE_lowFlowStats();
        double median_R_yearFlow = cdphe_lowFlowStats.annualMedianReturnPeriod(flowData, 5);
        median_R_yearFlow = DoubleMath.round(median_R_yearFlow, 2);
        
        //Format the results into a summary table
        String summaryTable = "\t1E3 Acute Monthly Low Flows\t\t7E3 Chronic Monthly Low Flows\t\t30E3 Chronic Monthly Low Flows\t" + "\r\nMonth\tYear\tFlow [cfs]\tYear \tFlow [cfs]\tYear\tFlow [cfs]";
        String[] months = {"Entire Record", "January", "February", "March", "April", "May", "June", "July", "August", "September", "October", "November", "December", "Lowest Monthly Flow"};
        for(int i=0; i<13; i++){
            int year_1day = (int) monthlyBiological_1day[1][i];
            int year_7day = (int) monthlyBiological_7day[1][i];
            int year_30day = (int) monthlyBiological_30day[1][i];
            summaryTable = summaryTable + "\r\n" + months[i] + "\t" + year_1day + "\t" + monthlyBiological_1day[0][i] + "\t" + year_7day + "\t" + monthlyBiological_7day[0][i] + "\t" + year_30day + "\t" + monthlyBiological_30day[0][i];
        }
        summaryTable = summaryTable + "\r\n" + months[13] + "\t" + monlyBiological_1day_min + "\t" + monlyBiological_7day_min + "\t" + monlyBiological_30day_min;
        summaryTable = summaryTable + "\r\n" + median_R_yearFlow;
        
        return summaryTable;
    }
    /**
     * Creates the input file for DFLOW parameters
     * @param directory  the working directory of the DFLOW executable
     * @param designFlowMethod  the design flow method (1=biologically based, 2=extreme value, 3=human health)
     * @param biologicalType  the parameter set for biologically based design flows (1=default, 2=custom, 3=EPA 3Q30)
     * @param M  The desired return period (for biological and extreme value)
     * @param R  the desired average number of years between excursions (for biological and extreme value)
     * @param clusterLength  the length of excursion clustering period (for biological)
     * @param clusterCountMax the maximum number of excursions (for biological)
     * @throws IOException 
     */
    private void writeDFLOWinp(String directory, int designFlowMethod, int biologicalType, int M, double R, int clusterLength, int clusterCountMax) throws IOException{
        //Output data to text file
        String path = directory + File.separator + "dflow_input.inp";
        FileWriter writer =  new FileWriter(path, false);
        PrintWriter print_line = new PrintWriter(writer);
        
        //Add fake header id to the file
        print_line.printf("     1     |jtype |format of input '.flo' file (1=DFLOW, 2=USGS)\r\n");
        print_line.printf("%7d   |I     |design flow method: (1=biologially based, 2=extreme value, 3=human health harmonic mean)\r\n", designFlowMethod);
        print_line.printf("%7d   |I1    |biological: type (1=default, 2=custom, 3=EPA default)\r\n", biologicalType);
        print_line.printf("%7d   |M     |biological/extreme val: number of days in flow averaging period (Biological: 1=acute, 4=chronic, Extreme Value: 7=default, must be between 1 and 30)\r\n", M);
        print_line.printf("%7.1f   |XMR   |biological/extreme val: average number of years between excursions (Biological: 3=default, Extreme Value: 10=default, must be between 1 and 500)\r\n", R);
        print_line.printf("%7d   |XMI   |biological: Length of excursion clustering period (120=default, must be greater than 0)\r\n", clusterLength);
        print_line.printf("%7d   |XME   |biological: Maximum number of excursions (5=default, must be greater than 0)\r\n", clusterCountMax);
        
        print_line.close();
        writer.close();
        System.out.println("Text File located at:\t" + path);
    }
    private void writeDFLOWflo(String directory, String[][] flowData) throws IOException{
        //Output data to text file
        String path = directory + File.separator + "river.flo";
        FileWriter writer =  new FileWriter(path, false);
        PrintWriter print_line = new PrintWriter(writer);
        
        //Add fake header id to the file
        print_line.printf("%s" + "\r\n", " 1234567");
        
        //Write timeseries data
        for(int i=0; i<flowData.length; i++){
            String year = flowData[i][0].substring(0,4);
            String month = flowData[i][0].substring(5,7);
            String day = flowData[i][0].substring(8);
            double flow = Double.parseDouble(flowData[i][1]);
            print_line.printf("%s    %s      %s            %11.2f\r\n", year, month, day, flow);
        }
        
        print_line.close();
        writer.close();
        System.out.println("Text File located at:\t" + path);
    }
}