"""
-----------------------------------------------------------------model_ini.py
ANUGA MODELLING TEMPLATE (Model Initialisation)
This initialisation script sets model run time flags and model-scenario-event

specific parameters and file references for use in the subsequent model_xx.py
scripts

To run in parallel you must use the folowing command line
mpirun -np n python model_run.py
where n is the total number of processors being used for the parallel run.

This script is imported by the model_data/build/checks/run/results scripts

<<<< THIS SCRIPT CONTAINS MOST OF THE CODE THAT WILL NORMALLY NEED TO  >>>>
<<<<       BE EDITED TO CREATE A PARTICULAR SIMULATION                 >>>>

Version         : 2.00 June 2019 flood template modified for tsunami modelling
Original Author : E Rigby (0437 250 500) ted.rigby@rienco.com.au
Modfied by      : C Kain, MRT, claire.kain@stategrowth.tas.gov.au

--------------------------------------------------------------------------------
"""
#################################################################################
#                                                                               #
#                     IMPORT THE REQUIRED PYTHON MODULES                        #
#                                                                               #
#################################################################################

# Standard python modules
import os
import sys
from os import sep, mkdir
from os.path import join, exists
from time import time,localtime, strftime


# Related ANUGA modules
import anuga
# Import the tools to manage the parallel interface
from anuga import distribute, myid, numprocs, finalize, barrier



#################################################################################
#                                                                               #
#           OPTIONALLY IDENTIFY THE USER & ORGANISATIONAL DETAILS               #
#       An optional block included in case more corporate ID desired            #
#################################################################################

org_name        = 'Mineral Resources Tasmania'
org_add1        = '30 Gordons Hill Road'
org_add2        = 'Rosny Park'
org_add3        = 'Tas 7018'
org_add4        = 'Australia'
org_phone       = '03 6165 4800'
org_web         = 'www.mrt.tas.gov.au'

user_name       = 'Ted'
user_email      = 'ted.rigby@rienco.com.au'

#################################################################################
#                                                                               #
#               SET THE CORE SIMULATION DESCRIPTORS FOR THIS RUN                #
#                                                                               #
#################################################################################

model        = 'NETAS'               # Set model name
scenario     = 'MRT2019_Scenario1'   # Set scenario
event_number = 'Puysegur3679'        # the event number - See suffix of nc scenario file from PTHA
runtype      = 'FINAL'               # (TRIAL, BASIC, STABLE, FINAL)
prefix       =  model+'_'+scenario+'_'+event_number+'_'+runtype   # run ID

anuga_verbose = True  # set to True to enable additional anuga internal output
model_verbose = True  # set to True to enable additional model specific output to screen
model_checks  = True  # set to True to enable additional model checks & plots to the checks dir
screen_log    = True # creates a logfile  of screen IO if true

# Note the anuga code uses final_time in evolve as duration rather than final_time
# This will eventually be corrected and will require mods to the following once implemented

starttime     = 0.0                     # start of simulation = sts start time for this event (approx 0 hrs))
duration      = 46800.0                 # simulation duration of 14400 secs (4 hrs) (was 14hrs for MRT/GA2009)
# note end simulation after a duration of 14400 secs (4 hrs after start at 0.0 hrs)
# note yieldstep is set later with  the run resolution reflecting coarse or fine simulations

#################################################################################
#                                                                               #
#             OPTIONALLY COPY SCREEN I-O TO LOGFILE & IDENTIFY BUILD            #
#                                                                               #
#################################################################################

# old screen log now depreciated - need to redevelop code for logging ???????????

# create the build identification of this run

user           = anuga.get_user_name()  # system names used in file stamping
host           = anuga.get_host_name()  # ditto
t0             = time()                 # record when the run started in secs
build_time     = strftime('%Y%m%d_%H%M%S',localtime())
build_ID       = build_time+'_build'+'_'+str(user)+'_on_'+str(host)
anuga_revision = anuga.get_revision_number()
anuga_ID       = str(anuga_revision)

if myid == 0:
    print ' '
    print prefix+ ' >>>> #########################################################################'
    print prefix+ ' >>>>'
    print prefix+ ' >>>> COMMENCING THIS MODEL SIMULATION USING ANUGA REVISION NUMBER %s  ' %(anuga_ID)
    print prefix+ ' >>>>'
    print prefix+ ' >>>>                           for the simulation named        '
    print prefix+ ' >>>>                          ', prefix
    print prefix+ ' >>>>                           and model Build ID  of '
    print prefix+ ' >>>>                          ', build_ID
    print prefix+ ' >>>> #########################################################################\n'
    print ' '

#################################################################################
#                                                                               #
#      SET THE REQUIRED ANUGA SPECIFIC RUN TIME FLAGS & VARIABLES               #
#                                                                               #
#################################################################################

minimum_allowed_height  = 0.001 # sets min depth (m) recognised in calcs
minimum_storable_height = 0.010 # writes d=0 (m)to sww if d < 0.010 - improves graphics
maximum_allowed_speed   = 5.00  # limits speed to 5m/sec max if d < minimum_allowed_height
flow_algorithm          = 'DE0' # set flow algorithm - refer anuga..config.py for options
dune_erosion            = True # simulates dune erosion in specified poly areas if true

# NOTES RE Flow_algorithm selection
# DE0 used Euler timestepping (and so is first order accurate in time).
# DE1 uses the average of two Euler steps to get a second order time stepping method.
# DE2 uses the average of three Euler steps to get a third order time stepping method.

# DE0 is very diffusive, so if you have a periodic forcing, the height of the wave will
# decrease  quite markedly. So not so good for complicated tsunami forcing problems.
# DE0 works quite well for flood modelling.
# Testing at GA suugests DE2 currently provides the best replication of tsunami behaviour

# DE2 will typically take 3 times as long as DE0 and DE1 takes 2 times as long as DE0.
# Recommend use of DE0 generally for speeed and switch to DE2 for final production runs
# when modelling a tsunami.

#################################################################################
#                                                                               #
#            SET THE SIMULATION SPECIFIC RUN TIME FLAGS & VARIABLES             #
#                                                                               #
#################################################################################

# assign names for the mesh and domain for this simulation

mesh_filename = prefix+'.tsh'   # msh for binary out and tsh for ascii output
domain_name   = prefix

##### Assign the different internal mesh resolutions/yieldsteps  to be used in this simulation

# Set the various model resolutions to suit the run type
if runtype  == 'TRIAL'           :  # only run with a extremely coarse mesh
    mesh_factor = 64.0
    yieldstep   = 600               # 10min
elif runtype ==  'BASIC'         :  # run with a relatively coarse mesh
    mesh_factor = 16.0
    yieldstep   = 240               # 6min
elif runtype == 'STABLE'         :  # run with reasonably fine mesh
    mesh_factor = 4.0
    yieldstep   = 120               # 2min
else                             :  # default is FINAL with unfactored res/yieldstep as below
    mesh_factor = 1.0
    yieldstep   = 60                # 1min

##### Assign the resolution class to factored resolution values for use in model_data

# As used for priority areas in GA2009
# vfine_res   = mesh_factor *    500.0  # m2 vrfine limit  32x32m tri  (typ high detail built up areas)
# As used for aos1 and aos2 areas in GA2009
# fine_res    = mesh_factor *   1500.0  # m2 fine limit    55x55m tri  (typ areas of interest generally)
# As used for sw area in GA2009
# medium_res  = mesh_factor *  30000.0  # m2 medium limit 245x245m tri (typ coastline generally and large estuaries)
# As used for default bounding poly res in GA2009
# coarse_res  = mesh_factor * 100000.0  # m2 coarse limit 450x450m tri (typ open deepwater areas)

# As used for resolution areas in MRT2016
xfine_res   = 50 * mesh_factor          # m2 10x10m triangles in high res plot area
vfine_res   = 200 * mesh_factor         # m2 20x20m triangles in std res plot areas of interest
fine_res    = 1250 * mesh_factor        # m2 50x50m triangles in coastal strip -10 to 10m AHD
medium_res  = 20000 * mesh_factor       # m2 200x200m triangles in moderate depth ocean -10 to -30m AMD
coarse_res  = 80000 * mesh_factor       # m2 400x400m triangles in deeper water -30m and deeper


default_res =  coarse_res               # applies to all bounding_poly pre internal res being added



##### Assign the initial default n to apply everywhere prior to applying spatial ns

default_n  = 0.010    # reflects estuary bed n - changes for land etc will be stamped over

##### Assign initial water level for this simulation (check??????????)

HAT = 0.8  #  NETAS
MSL = 0.0  #  NETAS
LAT = -0.8 #?  NETAS

# tide = xxx #  could be a timeseries
initial_water_level = HAT

################################################################################
#                                                                              #
#      ASSIGN AND CHECK THAT THE REQUIRED SIMULATION  DIRECTORIES EXIST        #
#                                                                              #
################################################################################

#-------------------------------------------------------------------------------
# Checks the current model directories  required to run the simulation are present
# This assumes a directory structure in which the  model run scripts,  model data files,
# model check files and model result files are in a directory called event - as outlined below.
#
# model\scenario\event\data    ( eg elevation, roughness, polygons, ICs,)
#                              \  Event \ boundary timeseries
#                     \checks  ( eg screen logs, check files and plots)
#                     \scripts ( py ini, data, build, check, run, results scripts)
#                     \results (sww result file)
#                             \ timeseries   (timeseries  files - (csv and plot files of eg stage H -v- time)
#                             \ spatial      (spatial eg asc max surface grid files)
#
#
# eg
# \SETAS2015\GA2009\58260\data     \ 58260 \ various timeseries files
#                                  \ various data files
#                        \ checks  \ various check files
#                        \ scripts \ python files
#                        \ results \ output.sww file
#                                  \ timeseries\ various timeseries files
#                                  \ spatial   \ various spatial files
#
#NOTE : \ scripts is the current directory when model_run.py is run
#         all directories are assigned relative to the scripts directory
#-------------------------------------------------------------------------------

##### Assign the required input directories

# ANUGA data data (input) directory. All input data is read from this directory
data_dir       = join('..' ,'data')

resolution_dir = join(data_dir, 'ResolutionPolygons2016_5revised')
gauges_dir     = data_dir+sep+'gauges'

# ANUGA results (output sww) directory and parent of timeseries ansd spatial output dirs
results_dir    = join('..','results')

# ANUGA results/timeseries (output) directory
timeseries_dir = results_dir+sep+'timeseries'

# ANUGA output spatial output files directory
spatial_dir    = results_dir+sep+'spatial'

# ANUGA output check  files directory
checks_dir   = join('..','checks')

##### Check the required output directories  are present and create if not


if not exists(results_dir)       :   # the required  results directory does not exist
    mkdir(results_dir)
    print prefix+ ' >>>> NOTE - Results directory %s was missing and has been created  ' %(results_dir)

if not exists(checks_dir)        :   # the required checks directory does not exist
    mkdir(checks_dir)
    print prefix+ ' >>>> NOTE - checks directory %s was missing and has been created  ' %(checks_dir)


if not exists(timeseries_dir)   :   # the required time series results directory does not exist
    mkdir(timeseries_dir)
    print prefix+ ' >>>> NOTE - Timeseries directory %s was missing and has been created  ' %(timeseries_dir)

if not exists(spatial_dir)      :   # the required spatial results directory does not exist
    mkdir(spatial_dir)
    print prefix+ ' >>>> NOTE - Spatial directory %s was missing and has been created  ' %(spatial_dir)

################################################################################
#                                                                              #
#             ASSIGN THE VARIOUS INPUT FILES FOR THIS SIMULATION               #
#                                                                              #
################################################################################

##### Assign the backdrop coastline polygon filename for model check plots
# Note - can not use shp for bounding poly as order in sts bdry points is
#        counterclockwise contray to GIS rules.

coastal_poly_filename  = 'netas_coastal_polygon.shp' # used only in the check plots

##### Assign the  bounding  polygon  for this simulation
# Note - can not use shp for bounding poly as order in sts bdry points is
#        counterclockwise contrary to GIS rules. Need to re-order sts if want to
#        use shp to input bounding poly in future!
bounding_poly_filename = 'bounding_poly_counterclockwise.csv'

##### Assign the internal polygon filenames for the MRT 2016 simulation

# Internal resolution polygons for MRT2016  - coarse is default set in bounding poly
int1_res_poly_filename = 'OceanApproachesNested.shp'      # medium res - ocean  ring 1 30m deep
int2_res_poly_filename = 'CoastalZoneNested.shp'          # fine res - coastal ring 2 10m deep - nested in ocean zone
# the CoastalZoneAbove10m is added in at this point to allow larger triangles in the higher ground of islands that
# are otherwise all at coastal res from Anuga's  need to have only simple polys (no holes)
int3_res_poly_filename = 'CoastalZoneAbove10mNested.shp'  # coarse res - land ring 3a > 10m AHD - nested in coastal zone
int4_res_poly_filename = 'PlotAreas.shp'                  # vfine res - plot areas merged ting 3b - nested in coastal zone
int5_res_poly_filename = 'StudyArea.shp'                  # High res - study area ring 4 - nested in plot areas

# same order as above ie int1_res_poly_filename polys will be modelled  at int1_res
int1_res = medium_res
int2_res = fine_res
int3_res = coarse_res
int4_res = vfine_res
int5_res = xfine_res

##### Assign the ANUGA input combined elevation spatial dataset filename for this simulation   (xyz,shp or asc)
elevation_points_filename  = 'elevation_netas.asc'      # assign the combined land/sea elevation points  txt filname

##### Assign the ANUGA input surface roughness (n) spatial dataset filename for this simulation
mannings_grid_filename     = 'mannings_n.asc'           # assign the combined land/sea roughness asc grid filename

#####  Assign the ANUGA input boundary conditions timeseries for this simulation
event_sts_filename        = 'puysegur_3679.sts'          # a preprepared netcdf of stage timeseries at various locations

##### Assign the ANUGA  gauge locations dataset filename at which  time-series output from this simulation is required
gauges_filename           = 'time_series_points.csv'     # assign the gauges (points) filename

if dune_erosion :
    ##### Assign the ANUGA erosion zone polygons filename containing polygons within which dune erosion will be simulated
    erosion_zones_filename     = 'NETAS_dune_erosion_zones.shp'              # assign the erosion zones (polygons) shp filename


#################################################################################
#                                                                               #
#       CREATE AND CHECK THE REQUIRED PATHNAMES ARE PRESENT                     #
#                                                                               #
#################################################################################


# create the full pathname for the output msh file
mesh_pathname = join(checks_dir,mesh_filename)                     # dump the mesh into the checks directory

# create the full pathname for coastal_filename and check exists
coastal_poly_pathname = join(data_dir,coastal_poly_filename)
if not exists(coastal_poly_pathname) :                 # check and stop run if absent
    print prefix+ ' >>>> ERROR - specified coastal poly  file %s missing ' %(coastal_poly_pathname)
    sys.exit(0)                                        # terminate the program

# create the full pathname for bounding_poly_filename and check exists
bounding_poly_pathname = join(data_dir,bounding_poly_filename)
if not exists(bounding_poly_pathname) :                # check and stop run if absent
    print prefix+ ' >>>> ERROR - specified bounding poly  file %s missing ' %(bounding_poly_pathname)
    sys.exit(0) 	                                   # terminate the program

# create the full pathname for the 5int res polys and check exists
int1_res_poly_pathname = join(data_dir,int1_res_poly_filename)  # Coarsest of the interior polys first
int2_res_poly_pathname = join(data_dir,int2_res_poly_filename)
int3_res_poly_pathname = join(data_dir,int3_res_poly_filename)
int4_res_poly_pathname = join(data_dir,int4_res_poly_filename)
int5_res_poly_pathname = join(data_dir,int5_res_poly_filename)


if not exists(int1_res_poly_pathname) :                # check and stop run if absent
    print prefix+ ' >>>> ERROR - specified int res poly  file %s missing ' %(int1_res_poly_pathname)
    sys.exit(0)
if not exists(int2_res_poly_pathname) :                # check and stop run if absent
    print prefix+ ' >>>> ERROR - specified int res poly  file %s missing ' %(int2_res_poly_pathname)
    sys.exit(0)
if not exists(int3_res_poly_pathname) :                # check and stop run if absent
    print prefix+ ' >>>> ERROR - specified int res poly  file %s missing ' %(int3_res_poly_pathname)
    sys.exit(0)
if not exists(int4_res_poly_pathname) :                # check and stop run if absent
    print prefix+ ' >>>> ERROR - specified int res poly  file %s missing ' %(int4_res_poly_pathname)
    sys.exit(0)
if not exists(int5_res_poly_pathname) :                # check and stop run if absent
    print prefix+ ' >>>> ERROR - specified int res poly  file %s missing ' %(int5_res_poly_pathname)
    sys.exit(0)

# Create the full pathname for the elevation points dataset and check exists
elevation_points_pathname = join (data_dir, elevation_points_filename)
if not exists(elevation_points_pathname) :              # check and stop run if absent
    print prefix+ ' >>>> ERROR - specified elevation  file %s missing ' %(elevation_points_pathname)
    sys.exit(0)

# Create the full pathname for the mannings roughness polys  and check exists
mannings_grid_pathname =  join(data_dir,mannings_grid_filename)
if not exists(mannings_grid_pathname) :                 # check and stop run if absent
    print prefix+ ' >>>> ERROR - specified mannings roughness file %s missing ' %(mannings_grid_pathname)
    sys.exit(0)

# Create the full pathname for the netcdf sts boundary conditions file  and check exists
event_sts_pathname =  join(data_dir,event_sts_filename)
if not exists(event_sts_pathname)    :                  # check and stop run if absent
    print prefix+ ' >>>> ERROR - specified  boundary conditions file %s missing ' %(event_sts_pathname)
    sys.exit(0)

# Create the full pathname for the file defining the output gauge lo0cations and check exists
gauges_pathname =  join(gauges_dir,gauges_filename)
if not exists(gauges_pathname)       :                  # check and stop run if absent
    print prefix+ ' >>>> ERROR - specified Gauges file %s missing ' %(gauges_pathname)
    sys.exit(0)

# Conditionally create the full pathname for the file defining the dune erosion zones and check exists
if dune_erosion :
    erosion_zones_pathname =  join(data_dir,erosion_zones_filename)
    if not exists(erosion_zones_pathname)       :                  # check and stop run if absent
        print prefix+ ' >>>> ERROR - specified erosion zone file %s missing ' %(erosion_zones_pathname)
        sys.exit(0)


