Extract statistics from .res11 or .res1d-files

Here is a code example on how to extract statistics from a res11 or res1d-files. Copy the code and the replace the file names and paths.

import os
from hydroinform.ts_tools import HydroTs
from hydroinform import DFS, Mike_river
import pandas as pd
import numpy as np


def get_single_point_ts(filename, rivername, topo_id, chainage):
    """
    Gets a single timeseries from a single point in a river
    Note, this method reads the entire output file, so use with care
    """
    with Mike_river.River_setup(filename) as res:
        reach = next((r for r in res.reaches if r.name == 'rivername' and r.topo_id==topo_id), None)
        if reach != None:
            gp = next((p for p in reach.grid_points if p.chainage == chainage), None)
            ts = res.get_timeseries(gp)
            return ts

def statistics(filename):
    """
    Returns a dataframe with statistics on the water flow from either a .res11- or res1d-file
    """
    res11=Mike_river.River_setup.from_file(filename)

    columns=['Branch','Chainage','X','Y', 'Min', 'Max', 'Average', 'DFFV_EQR', 'DVFI_EQR', 'DVPI_EQR', 'Q95', 'Q75', 'Q50', 'Q25', 'MedianMin']    
    nbqpoints = np.sum([len(r.grid_points_q) for r in res11.reaches])
    data = pd.DataFrame(columns=columns,index = range(nbqpoints))
   
    index =0
    for reach in res11.reaches: #Loop all reaches
        for gp in reach.grid_points_q: #Loop all q-points
            ts = res11.get_timeseries(gp) # get the discharge time series
            hts = HydroTs(ts) # Create a hydro time series
            #set all values
            data.at[index, 'Branch']=reach.name
            data.at[index, 'Chainage']=gp.chainage
            data.at[index, 'X']=gp.x
            data.at[index, 'Y']=gp.y
            data.at[index, 'Min']=ts.min()
            data.at[index, 'Max']=ts.max()
            data.at[index, 'Average']=np.average(ts)
            data.at[index, 'DFFV_EQR']=hts.DFFV_EQR()
            data.at[index, 'DVFI_EQR']=hts.DVFIEQR()
            data.at[index, 'DVPI_EQR']=hts.DVPIEQR()
            data.at[index, 'Q95']=hts.Q95
            data.at[index, 'Q75']=hts.Q75
            data.at[index, 'Q50']=hts.Q50
            data.at[index, 'Q25']=hts.Q25
            data.at[index, 'MedianMin']=hts.medianmin()
            index+=1
    res11.dispose()
    return data



if __name__ == "__main__":

    #get_single_point_ts(r'QPoints2\R201803_vandweb\m11_res\abs\DK1_Abs.res11','MERNAA', 'novana_model', 8251).to_csv(r'c:\temp\testData.csv')

    import shapefile
    df =pd.DataFrame()

    MikeHydroFiles=[r'C:\projects\GEUS\DKM2019\Model\Result\DK1\DK1_mh_2018.res1d'
           r'C:\projects\GEUS\DKM2019\Model\Result\DK2\DK2_mh_2018.res1d',
           r'C:\projects\GEUS\DKM2019\Model\Result\DK3\DK3_mh_2018.res1d',
           r'C:\projects\GEUS\DKM2019\Model\Result\DK4\DK4_mh_2018.res1d',
           r'C:\projects\GEUS\DKM2019\Model\Result\DK5\DK5_mh_2018.res1d',
           r'C:\projects\GEUS\DKM2019\Model\Result\DK6\DK6_mh_2018.res1d'
           ]

    #Loop the files
    for file in MikeHydroFiles:
        print('Reading: ' + file)
        #Join the dataframes
        df=df.append(statistics(file), ignore_index=True)
    
    #Write to shapefile
    with shapefile.Writer(r'c:\temp\output.shp', shapeType=shapefile.POINT) as w:
        w.autoBalance=True
        #First column is text (Branch name)
        w.field(df.columns[0],'C')
        #Remaining columns are numbers
        for c in df.columns[1:]:
            w.field(c, 'F', decimal=10)

        for i, row in df.iterrows():
            w.record(*row)
            w.point(row['X'],row['Y'])