RoutingPreProc.py 3.5 KB
Newer Older
1 2 3 4 5 6 7 8 9 10 11 12
#
#
#
import numpy as np
import os, sys
from mpi4py import MPI
import gc
#
import matplotlib as mpl
mpl.use('Agg')
from matplotlib.backends.backend_pdf import PdfPages
import matplotlib.pyplot as plt
13
import time
14 15 16 17
#
# Gert the information from the configuration file.
#
import configparser
18
config=configparser.ConfigParser()
19
config.read("run.def")
POLCHER Jan's avatar
POLCHER Jan committed
20
nbasmax=config.getint("OverAll", "nbasmax")
21
numop=config.getint("OverAll", "numop", fallback=100)
22
OutGraphFile=config.get("OverAll","GraphFile")
23
DumpHydroSuper=config.getboolean("OverAll","DumpHydroSuper",fallback=False)
24
wfile=config.get("OverAll","WeightFile",fallback=None)
25 26 27 28 29 30
#
import ModelGrid as MG
import Partition as PA
import RPPtools as RPP
import HydroGrid as HG
import Interface as IF
31
from Truncate import trunc as TR
32 33 34 35 36 37 38 39 40 41 42 43 44 45
#
# Logging in MPI : https://groups.google.com/forum/#!topic/mpi4py/SaNzc8bdj6U
#
import getargs
log_master, log_world = getargs.getLogger()
INFO, DEBUG, ERROR = log_master.info, log_master.debug, log_world.error
INFO_ALL, DEBUG_ALL = log_world.info, log_world.debug
#
# Read full grid and partition the domain.
#
gg = MG.GlobalGrid()
#
comm = MPI.COMM_WORLD
szhalo=1
46
nbcore=comm.Get_size()
47 48 49 50
rank=comm.Get_rank()
#
# Region of grid to be treated
#
51
part = PA.partition(gg.ni, gg.nj, gg.land, comm, nbcore, szhalo, rank)
52 53
#
modelgrid=MG.ModelGrid(part.ihstart+gg.igstart,part.nih,part.jhstart+gg.jgstart,part.njh)
54
#
55 56
INFO("Longitude interval on proc "+str(rank)+" = "+str(modelgrid.box_land[0][0])+" : "+str(modelgrid.box_land[0][1]))
INFO("Latitude interval on proc "+str(rank)+" = "+str(modelgrid.box_land[1][0])+" : "+str(modelgrid.box_land[1][1]))
57
#
58
hydrogrid=HG.HydroGrid(modelgrid.box_land, wfile)
59
#
60
# Computes weights of overlap of modelgrid and hydrogrid
61
w = RPP.compweights(wfile, part, modelgrid, hydrogrid)
62 63 64
#
#
#
65 66
nbpt = len(w.index)
nbvmax = part.domainmax(max(w.hpts))
67 68 69
INFO("nbpt : {0}".format(nbpt))
INFO("nbvmax : {0}".format(nbvmax))
INFO("nbasmax : {0}".format(nbasmax))
70 71 72
#
# Extract hydo data from file
#
73
INFO("=================== HYDRODATA ====================")
74
hydrodata = HG.HydroData(hydrogrid.ncfile, hydrogrid.box, w.index)
75

76
INFO("=================== INITATMGRID ====================")
77
IF.initatmgrid(rank, nbcore, nbpt, modelgrid)
78

79
INFO("=================== HOVERLAP ====================")
80
hoverlap = IF.HydroOverlap(nbpt, nbvmax, w.hpts, w.index, w.area, w.lon, w.lat, part, modelgrid, hydrodata)
81 82 83
#
# Do some memory management and synchronize procs.
#
84

85
del w
86
gc.collect()
87
comm.Barrier()
88

89
hsuper = IF.HydroSuper(nbvmax, hydrodata, hoverlap, nbasmax, part)
90
#
91
INFO("=================== LINKUP ====================")
92
hsuper.linkup(hydrodata)
93

94 95 96
#
# Do some memory management and synchronize procs.
#
97

98
del hoverlap
99
gc.collect()
100
comm.Barrier()
101 102 103
#
hsuper.correct_outflows(part)
#
104
INFO("=================== Compute fetch ====================")
105
t = time.time()
POLCHER Jan's avatar
POLCHER Jan committed
106
hsuper.fetch(part)
107 108
comm.Barrier()
t1 = time.time()
109
print("Time for fetch: {:0.2f} s.".format(t1-t))
110 111
comm.Barrier()

112
if DumpHydroSuper :
113
    INFO("Dumping HydroSuper")
114
    hsuper.dumpnetcdf(OutGraphFile.replace(".nc","_HydroSuper.nc"), gg, modelgrid, part)
115

Anthony's avatar
Anthony committed
116
INFO("=================== Truncate ====================")
117
print("Rank : {0} - Basin_count Before Truncate : ".format(part.rank), hsuper.basin_count)
118
hs = TR(hsuper, part, comm, modelgrid, numop = numop)
119 120
print("Rank : {0} - Basin_count After Truncate : ".format(part.rank), hsuper.basin_count)

Anthony's avatar
Anthony committed
121
INFO("=================== HydroGraph ====================")
122 123
hgraph = IF.HydroGraph(nbasmax, hsuper, part, modelgrid)
hgraph.dumpnetcdf(OutGraphFile, gg, modelgrid, part)
124 125

IF.closeinterface(comm)