diff --git a/DocumentationInterface b/DocumentationInterface
index 499d46d1f82f5e4167b0ee67ef12c4ac1fa56932..4c4d1c54d134f9fc7bc8514e73f1e031b498b8fa 100644
--- a/DocumentationInterface
+++ b/DocumentationInterface
@@ -186,7 +186,7 @@ inflow_number : rank-2 array('i') with bounds (nbpt,8 * nbxmax_in)
inflow_grid : rank-3 array('i') with bounds (nbpt,8 * nbxmax_in,8 * nbxmax_in)
inflow_basin : rank-3 array('i') with bounds (nbpt,8 * nbxmax_in,8 * nbxmax_in)
====================================================================
-outflow_grid_new = fetch(basin_count,basin_area,basin_id,basin_hierarchy,basin_fac,outflow_grid,outflow_basin,fetch_basin,[nbpt,nwbas])
+outflow_uparea = fetch(basin_count,basin_area,basin_id,basin_hierarchy,basin_fac,outflow_grid,outflow_basin,partial_sum,fetch_basin,[nbpt,nwbas])
Wrapper for ``fetch``.
@@ -199,6 +199,7 @@ basin_hierarchy : input rank-2 array('f') with bounds (nbpt,nwbas)
basin_fac : input rank-2 array('f') with bounds (nbpt,nwbas)
outflow_grid : in/output rank-2 array('i') with bounds (nbpt,nwbas)
outflow_basin : input rank-2 array('i') with bounds (nbpt,nwbas)
+partial_sum : input rank-2 array('f') with bounds (nbpt,nwbas)
fetch_basin : in/output rank-2 array('f') with bounds (nbpt,nwbas)
Other Parameters
@@ -210,15 +211,16 @@ nwbas : input int, optional
Returns
-------
-outflow_grid_new : rank-2 array('i') with bounds (nbpt,nwbas)
+outflow_uparea : rank-1 array('f') with bounds (nbpt*nwbas)
====================================================================
-routing_area,routing_cg,topo_resid,route_togrid,route_tobasin,route_nbintobas,global_basinid,route_outlet,route_type,origin_nbintobas,routing_fetch = truncate(nbasmax,basin_count,basin_notrun,basin_area,basin_cg,basin_topoind,fetch_basin,basin_id,basin_coor,basin_type,basin_flowdir,outflow_grid,outflow_basin,inflow_number,inflow_grid,inflow_basin,[nbpt,nbxmax_in,nwbas])
+routing_area,routing_cg,topo_resid,route_togrid,route_tobasin,route_nbintobas,global_basinid,route_outlet,route_type,origin_nbintobas,routing_fetch = truncate(nbasmax,num_largest,basin_count,basin_notrun,basin_area,basin_cg,basin_topoind,fetch_basin,basin_id,basin_coor,basin_type,basin_flowdir,outflow_grid,outflow_basin,inflow_number,inflow_grid,inflow_basin,[nbpt,nbxmax_in,nwbas])
Wrapper for ``truncate``.
Parameters
----------
nbasmax : input int
+num_largest : input int
basin_count : in/output rank-1 array('i') with bounds (nbpt)
basin_notrun : in/output rank-1 array('i') with bounds (nbpt)
basin_area : in/output rank-2 array('f') with bounds (nbpt,nwbas)
diff --git a/F90subroutines/routing_interface.f90 b/F90subroutines/routing_interface.f90
index 4d72ec328e9b1d9d0b7380cad3c2b93751f31583..67d52162c3098a9593e6a44870d8c9f3e96e4a99 100644
--- a/F90subroutines/routing_interface.f90
+++ b/F90subroutines/routing_interface.f90
@@ -318,8 +318,8 @@ SUBROUTINE areanorm(nbpt, nwbas, basin_count, basin_area, basin_area_norm)
!
END SUBROUTINE areanorm
-SUBROUTINE fetch(nbpt, nwbas, basin_count, basin_area, basin_id, basin_hierarchy, basin_fac,&
- & outflow_grid, outflow_basin, fetch_basin, outflow_grid_new)
+SUBROUTINE fetch(nbpt, nwbas, basin_count, basin_area, basin_id, basin_hierarchy, basin_fac, &
+ & outflow_grid, outflow_basin, partial_sum, fetch_basin, outflow_uparea)
!
USE ioipsl
USE grid
@@ -336,23 +336,56 @@ SUBROUTINE fetch(nbpt, nwbas, basin_count, basin_area, basin_id, basin_hierarchy
INTEGER(i_std), DIMENSION(nbpt,nwbas), INTENT(in) :: basin_id !!
REAL(r_std), DIMENSION(nbpt,nwbas), INTENT(in) :: basin_hierarchy
REAL(r_std), DIMENSION(nbpt,nwbas), INTENT(in) :: basin_fac
- INTEGER(i_std), DIMENSION(nbpt,nwbas), INTENT(inout) :: outflow_grid !! Type of outflow on the grid box (unitless)
INTEGER(i_std), DIMENSION(nbpt,nwbas), INTENT(in) :: outflow_basin !!
+ REAL(r_std), DIMENSION(nbpt,nwbas), INTENT(in) :: partial_sum
!
!! IN-OUTPUT VARIABLES
+ INTEGER(i_std), DIMENSION(nbpt,nwbas), INTENT(inout) :: outflow_grid !! Type of outflow on the grid box (unitless)
REAL(r_std), DIMENSION(nbpt,nwbas), INTENT(inout) :: fetch_basin !!
! Output
- INTEGER(i_std), DIMENSION(nbpt,nwbas), INTENT(out) :: outflow_grid_new !! Type of outflow on the grid box (unitless)
- !
- outflow_grid_new(:,:) = outflow_grid(:,:)
+ REAL(r_std), DIMENSION(nbpt*nwbas), INTENT(out) :: outflow_uparea
!
CALL routing_reg_fetch(nbpt, area, contfrac, nwbas, basin_count, basin_area, basin_id, basin_hierarchy, basin_fac,&
- & outflow_grid_new, outflow_basin, fetch_basin)
+ & outflow_grid, outflow_basin, partial_sum, fetch_basin, outflow_uparea)
!
END SUBROUTINE fetch
-SUBROUTINE truncate(nbpt, nbxmax_in, nbasmax, nwbas, basin_count, basin_notrun, basin_area, basin_cg, basin_topoind,&
+SUBROUTINE rivclassification(nbpt, nwbas, basin_count, outflow_grid, outflow_basin, fetch_basin, largest_rivarea, num_largest)
+ !
+ USE defprec
+ !
+ !! INPUT VARIABLES
+ INTEGER(i_std), INTENT (in) :: nbpt !! Domain size (unitless)
+ INTEGER(i_std), INTENT (in) :: nwbas !!
+ INTEGER(i_std), DIMENSION(nbpt), INTENT(in) :: basin_count !!
+ INTEGER(i_std), DIMENSION(nbpt,nwbas), INTENT(inout) :: outflow_grid
+ INTEGER(i_std), DIMENSION(nbpt,nwbas), INTENT(in) :: outflow_basin !!
+ REAL(r_std), DIMENSION(nbpt,nwbas), INTENT(in) :: fetch_basin !!
+ REAL(r_std), INTENT(in) :: largest_rivarea
+ INTEGER(i_std), INTENT (out) :: num_largest
+ !
+ !! LOCAL VARIABLES
+ INTEGER(i_std) :: ib, ij
+ !
+ !
+ num_largest = 0
+ !
+ DO ib=1,nbpt
+ !
+ DO ij=1,basin_count(ib)
+ !
+ IF (outflow_grid(ib,ij) .LT. 0 .AND. fetch_basin(ib,ij) .GE. largest_rivarea) THEN
+ num_largest = num_largest + 1
+ outflow_grid(ib,ij) = -1
+ ENDIF
+ !
+ ENDDO
+ ENDDO
+ !
+END SUBROUTINE rivclassification
+
+SUBROUTINE truncate(nbpt, nbxmax_in, nbasmax, nwbas, num_largest, basin_count, basin_notrun, basin_area, basin_cg, basin_topoind,&
& fetch_basin, basin_id, basin_coor, basin_type, basin_flowdir, outflow_grid, outflow_basin, inflow_number,&
& inflow_grid, inflow_basin, routing_area, routing_cg, topo_resid, route_togrid, route_tobasin, route_nbintobas, &
& global_basinid, route_outlet, route_type, origin_nbintobas, routing_fetch)
@@ -368,7 +401,8 @@ SUBROUTINE truncate(nbpt, nbxmax_in, nbasmax, nwbas, basin_count, basin_notrun,
INTEGER(i_std), INTENT (in) :: nbpt !! Domain size (unitless)
INTEGER(i_std), INTENT (in) :: nbxmax_in, nbasmax
INTEGER(i_std), INTENT (in) :: nwbas !!
-
+ INTEGER(i_std), INTENT (in) :: num_largest
+ !
INTEGER(i_std), DIMENSION(nbpt), INTENT(inout) :: basin_count !!
INTEGER(i_std), DIMENSION(nbpt), INTENT(inout) :: basin_notrun !!
INTEGER(i_std), DIMENSION(nbpt,nwbas), INTENT(inout) :: basin_id !!
@@ -407,7 +441,7 @@ SUBROUTINE truncate(nbpt, nbxmax_in, nbasmax, nwbas, basin_count, basin_notrun,
STOP
ENDIF
- CALL routing_reg_truncate(nbpt, nbasmax, area, contfrac, nwbas, basin_count, basin_notrun, basin_area, basin_cg, &
+ CALL routing_reg_truncate(nbpt, nbasmax, area, contfrac, nwbas, num_largest, basin_count, basin_notrun, basin_area, basin_cg, &
& basin_topoind, fetch_basin, basin_id, basin_coor, basin_type, basin_flowdir, outflow_grid, outflow_basin, inflow_number,&
& inflow_grid, inflow_basin)
diff --git a/F90subroutines/routing_reg.f90 b/F90subroutines/routing_reg.f90
index 6f4ad0125ced42b5807f6c439615c3513ceea0bd..6ed789ebc8fadeac5b4507d0dee950bb068d9315 100644
--- a/F90subroutines/routing_reg.f90
+++ b/F90subroutines/routing_reg.f90
@@ -10,7 +10,6 @@ MODULE routing_reg
INTEGER(i_std), SAVE :: nbvmax=440 !! The maximum number of basins we can handle at any time during the generation of the maps (unitless)
INTEGER(i_std), SAVE :: nbxmax=63 !! The maximum number of points in one direction (NS or WE) using in routing_reg_getgrid (unitless)
INTEGER(i_std), SAVE :: maxpercent=2 !! The maximum area percentage of a sub-basin in a grib-box (default = 2%)
- INTEGER(i_std), SAVE :: num_largest = 50
INTEGER(i_std), SAVE :: nbpt_save
INTEGER(i_std), SAVE :: nbasmax_save
@@ -2692,7 +2691,7 @@ SUBROUTINE routing_reg_linkup(nbpt, neighbours, nwbas, basin_count, basin_area,
!_ ================================================================================================================================
SUBROUTINE routing_reg_fetch(nbpt, gridarea, contfrac, nwbas, basin_count, basin_area, basin_id, basin_hierarchy, basin_fac,&
- & outflow_grid, outflow_basin, fetch_basin)
+ & outflow_grid, outflow_basin, partial_sum, fetch_basin, outflow_uparea)
!
IMPLICIT NONE
!
@@ -2709,10 +2708,12 @@ SUBROUTINE routing_reg_fetch(nbpt, gridarea, contfrac, nwbas, basin_count, basin
REAL(r_std), DIMENSION(nbpt,nwbas), INTENT(in) :: basin_hierarchy
REAL(r_std), DIMENSION(nbpt,nwbas), INTENT(in) :: basin_fac
INTEGER(i_std), DIMENSION(nbpt,nwbas), INTENT(in) :: outflow_basin !!
+ REAL(r_std), DIMENSION(nbpt,nwbas), INTENT(in) :: partial_sum
!
!! INPUT/OUTPUT VARIABLES
INTEGER(i_std), DIMENSION(nbpt,nwbas), INTENT(inout) :: outflow_grid !! Type of outflow on the grid box (unitless)
REAL(r_std), DIMENSION(nbpt,nwbas), INTENT(inout) :: fetch_basin !!
+ REAL(r_std), DIMENSION(nbpt*nwbas), INTENT(out) :: outflow_uparea
!
!! LOCAL VARIABLES
INTEGER(i_std) :: ib, ij, ff(1), it, itt, igrif, ibasf, nboutflow !! Indices (unitless)
@@ -2720,41 +2721,90 @@ SUBROUTINE routing_reg_fetch(nbpt, gridarea, contfrac, nwbas, basin_count, basin
INTEGER(i_std), DIMENSION(nbpt*nbvmax,2) :: tmpindex !!
!
INTEGER(i_std), PARAMETER :: maxiterations=10000
-
+ !
!_ ================================================================================================================================
!
+ ! Propagate first the partial sums from the halo into the domain
!
DO ib=1,nbpt
+ !
+ IF ( SUM(partial_sum(ib,:)) > 0 ) THEN
+ WRITE(numout,*) ib, "Have some partial system to add in point", lalo(ib,2), lalo(ib,1)
+ ENDIF
!
DO ij=1,basin_count(ib)
!
- fetch_basin(ib, ij) = fetch_basin(ib, ij) + basin_area(ib,ij)
- !
- igrif = outflow_grid(ib,ij)
- ibasf = outflow_basin(ib,ij)
- !
- itt = 0
+ IF (partial_sum(ib,ij) > EPSILON(partial_sum)) THEN
+ !
+ fetch_basin(ib, ij) = MAX(fetch_basin(ib, ij), partial_sum(ib,ij))
+ !
+ igrif = outflow_grid(ib,ij)
+ ibasf = outflow_basin(ib,ij)
+ !
+ itt = 0
+ !
+ DO WHILE (igrif .GT. 0)
+ !
+ fetch_basin(igrif,ibasf) = MAX(fetch_basin(igrif,ibasf), partial_sum(ib, ij))
+ !
+ it = outflow_grid(igrif, ibasf)
+ ibasf = outflow_basin(igrif, ibasf)
+ igrif = it
+ itt = itt + 1
+ IF ( itt .GT. maxiterations) THEN
+ IF ( itt .GT. maxiterations+20) THEN
+ CALL ipslerr_p(3,'routing_reg_fetch','Problem in propagating partial sum','','')
+ ENDIF
+ ENDIF
+ ENDDO
+ ENDIF
+ ENDDO
+ ENDDO
+ !
+ ! Add the areas contributed by the core region of the domain.
+ !
+ DO ib=1,nbpt
+ !
+ DO ij=1,basin_count(ib)
!
- DO WHILE (igrif .GT. 0)
+ IF (partial_sum(ib,ij) <= 0.01) THEN
+ !
+ fetch_basin(ib, ij) = fetch_basin(ib, ij) + basin_area(ib,ij)
!
- fetch_basin(igrif,ibasf) = fetch_basin(igrif,ibasf) + basin_area(ib, ij)
- it = outflow_grid(igrif, ibasf)
- ibasf = outflow_basin(igrif, ibasf)
- igrif = it
- itt = itt + 1
- IF ( itt .GT. maxiterations) THEN
- WRITE(numout,&
- "('Grid ',I5, ' and basin ',I5, ' did not converge after iteration ',I5)") ib, ij, itt
- WRITE(numout, "('==> We flow through grid ',I5,' and basin ',I5)") igrif, ibasf
- WRITE(numout,*) '==> Basin ID :', basin_id(igrif,ibasf), "Hierarchy :", basin_hierarchy(igrif,ibasf)
- WRITE(numout,*) "==> Flow accumulation :", basin_fac(igrif,ibasf)
- WRITE(numout,*) "==> Coordinates : ", lalo_g(igrif,2), lalo_g(igrif,1)
+ igrif = outflow_grid(ib,ij)
+ ibasf = outflow_basin(ib,ij)
+ !
+ itt = 0
+ !
+ DO WHILE (igrif .GT. 0 )
!
- IF ( itt .GT. maxiterations+20) THEN
- CALL ipslerr_p(3,'routing_reg_fetch','Problem...','','')
+ IF (partial_sum(igrif,ibasf) <= 0.01) THEN
+ ! In the core regions
+ fetch_basin(igrif,ibasf) = fetch_basin(igrif,ibasf) + basin_area(ib, ij)
+ ELSE
+ ! In the halo of the domain avoid double counting
+ fetch_basin(igrif,ibasf) = MAX(fetch_basin(igrif,ibasf), partial_sum(ib, ij))
ENDIF
- ENDIF
- ENDDO
+ !
+ it = outflow_grid(igrif, ibasf)
+ ibasf = outflow_basin(igrif, ibasf)
+ igrif = it
+ itt = itt + 1
+ IF ( itt .GT. maxiterations) THEN
+ WRITE(numout,&
+ "('Grid ',I5, ' and basin ',I5, ' did not converge after iteration ',I5)") ib, ij, itt
+ WRITE(numout, "('==> We flow through grid ',I5,' and basin ',I5)") igrif, ibasf
+ WRITE(numout,*) '==> Basin ID :', basin_id(igrif,ibasf), "Hierarchy :", basin_hierarchy(igrif,ibasf)
+ WRITE(numout,*) "==> Flow accumulation :", basin_fac(igrif,ibasf)
+ WRITE(numout,*) "==> Coordinates : ", lalo_g(igrif,2), lalo_g(igrif,1)
+ !
+ IF ( itt .GT. maxiterations+20) THEN
+ CALL ipslerr_p(3,'routing_reg_fetch','Problem...','','')
+ ENDIF
+ ENDIF
+ ENDDO
+ !
+ ENDIF
!
ENDDO
!
@@ -2763,40 +2813,40 @@ SUBROUTINE routing_reg_fetch(nbpt, gridarea, contfrac, nwbas, basin_count, basin
WRITE(numout,*) 'The smallest FETCH :', MINVAL(fetch_basin)
WRITE(numout,*) 'The largest FETCH :', MAXVAL(fetch_basin)
!
- ! Now we set for the 'num_largest' largest basins the outflow condition as stream flow
- ! (i.e. outflow_grid = -1) and all other outflow basins are set to coastal flow
- ! (i.e. outflow_grid = -2). The return flow is not touched.
+ ! We set all outflow points to coastal flow. This will be adjusted later once all procs have
+ ! exchanged their information. So we will have outflow_grid = -2.
!
nboutflow = 0
+ outflow_uparea(:) = zero
!
DO ib=1,nbpt
!
DO ij=1,basin_count(ib)
!
! We do not need any more the river flow flag as we are going to reset it.
+ ! We archive the upstream area of all outflow points so that we can sort them.
!
- IF ( outflow_grid(ib,ij) .EQ. -1) THEN
- outflow_grid(ib,ij) = -2
- ENDIF
- !
- IF ( outflow_grid(ib,ij) .EQ. -2) THEN
- !
+ IF ( outflow_grid(ib,ij) .LT. 0) THEN
+ IF ( outflow_grid(ib,ij) .EQ. -1) THEN
+ outflow_grid(ib,ij) = -2
+ ENDIF
nboutflow = nboutflow + 1
- tmp_area(nboutflow) = fetch_basin(ib,ij)
- tmpindex(nboutflow,1) = ib
- tmpindex(nboutflow,2) = ij
- !
+ IF ( nboutflow <= nbpt*nwbas) THEN
+ outflow_uparea(nboutflow) = fetch_basin(ib,ij)
+ ELSE
+ ! Not enought place so we replace a smaller one.
+ IF (fetch_basin(ib,ij) > MINVAL(outflow_uparea)) THEN
+ ff = MINLOC(outflow_uparea)
+ outflow_uparea(ff(1)) = fetch_basin(ib,ij)
+ ELSE
+ ! Ignore basin
+ ENDIF
+ ENDIF
ENDIF
!
ENDDO
ENDDO
!
- DO ib=1, num_largest
- ff = MAXLOC(tmp_area(1:nboutflow))
- outflow_grid(tmpindex(ff(1),1), tmpindex(ff(1),2)) = -1
- tmp_area(ff(1)) = zero
- ENDDO
- !
END SUBROUTINE routing_reg_fetch
!
!
@@ -2839,8 +2889,8 @@ SUBROUTINE routing_reg_fetch(nbpt, gridarea, contfrac, nwbas, basin_count, basin
!! \n
!_ ================================================================================================================================
-SUBROUTINE routing_reg_truncate(nbpt, nbasmax, gridarea, contfrac, nwbas, basin_count, basin_notrun, basin_area, basin_cg, &
- & basin_topoind, fetch_basin, basin_id, basin_coor, basin_type, basin_flowdir, outflow_grid, outflow_basin, &
+SUBROUTINE routing_reg_truncate(nbpt, nbasmax, gridarea, contfrac, nwbas, num_largest, basin_count, basin_notrun, basin_area, &
+ & basin_cg, basin_topoind, fetch_basin, basin_id, basin_coor, basin_type, basin_flowdir, outflow_grid, outflow_basin, &
& inflow_number, inflow_grid, inflow_basin)
!
IMPLICIT NONE
@@ -2857,6 +2907,7 @@ SUBROUTINE routing_reg_truncate(nbpt, nbasmax, gridarea, contfrac, nwbas, basin_
REAL(r_std), DIMENSION(nbpt), INTENT(in) :: contfrac !! Fraction of land in each grid box (unitless;0-1)
!
INTEGER(i_std), INTENT(in) :: nwbas !!
+ INTEGER(i_std), INTENT (in) :: num_largest
INTEGER(i_std), DIMENSION(nbpt), INTENT(inout) :: basin_count !!
INTEGER(i_std), DIMENSION(nbpt), INTENT(inout) :: basin_notrun !!
INTEGER(i_std), DIMENSION(nbpt,nwbas), INTENT(inout) :: basin_id !!
diff --git a/Interface.py b/Interface.py
index e264c61a4bfdea9ad29a16107580c3ecda17b0a2..6dffef3d5bb668ba1536561177ebc794977d1b58 100644
--- a/Interface.py
+++ b/Interface.py
@@ -160,40 +160,56 @@ class HydroSuper :
self.nbcoastal, self.coastal_basin, float(hydrodata.basinsmax))
return
#
- def fetch(self, part) :
+ def fetch(self, part, modelgrid) :
+ #
+ largest_pos = -50
#
fetch_basin = np.zeros(self.basin_area.shape, dtype=np.float32, order='F')
- self.fetch_basin = np.zeros(self.basin_area.shape, dtype=np.float32, order='F')
#
self.basin_area = routing_interface.areanorm(self.basin_count, self.basin_area)
- partial_sum = np.copy(self.basin_area)
+ partial_sum = np.zeros(self.basin_area.shape, dtype=np.float32, order='F')
#
for i in range(part.size) :
fetch_basin[:,:] = 0.0
- self.outflow_grid = routing_interface.fetch(self.basin_count, partial_sum, self.basin_id, self.basin_hierarchy, \
- self.basin_fac, self.outflow_grid, self.outflow_basin, fetch_basin)
- print("Rank : ", part.rank," iter =", i," Fetch_basin range :", np.min(fetch_basin), np.mean(fetch_basin), np.max(fetch_basin))
- print("Rank : ", part.rank," Partial sum range :", np.min(partial_sum), np.mean(partial_sum), np.max(partial_sum))
+ outflow_uparea = routing_interface.fetch(self.basin_count, self.basin_area, self.basin_id, self.basin_hierarchy, \
+ self.basin_fac, self.outflow_grid, self.outflow_basin, partial_sum, fetch_basin)
partial_sum = np.copy(fetch_basin)
part.landsendtohalo(partial_sum, order='F')
- partial_sum = part.copycore(partial_sum, self.basin_area, order='F')
- print("Rank :", part.rank," iter =", i," Copy core done")
+ partial_sum = part.zerocore(partial_sum, order='F')
+ #
+ yy=modelgrid.landscatter(np.sum(fetch_basin, axis=1)/np.sum(self.basin_area, axis=1))
- self.fetch_basin[:,:] = fetch_basin[:,:]
+ self.fetch_basin = np.copy(fetch_basin)
+ #
+ # Find area the largest basins need at least to have.
+ #
+ xtmp = np.hstack(part.comm.allgather(outflow_uparea[np.where(outflow_uparea > 0.0)]))
+ sorted_outareas = np.sort(np.array(xtmp))
+ largest_rivarea = sorted_outareas[largest_pos]
+ #
+ #
+ #
+ yy=modelgrid.landscatter(np.sum(self.fetch_basin, axis=1)/np.sum(self.basin_area, axis=1))
+ print("Rank :"+str(part.rank)+" OUT of fetch =+=+=+=+=+=+=+=+=+= \n"+str(yy)+"\n =+=+=+=+=+=+=+=+=+=")
+ self.num_largest = routing_interface.rivclassification(self.basin_count, self.outflow_grid, self.outflow_basin, \
+ self.fetch_basin, largest_rivarea)
return
#
#
#
class HydroGraph :
- def __init__(self, nbasmax, hydrosuper) :
+ def __init__(self, nbasmax, hydrosuper, part, modelgrid) :
self.nbasmax = nbasmax
self.routing_area, self.routing_cg, self.topo_resid, self.route_togrid, self.route_tobasin, self.route_nbintobas, self.global_basinid, \
self.route_outlet, self.route_type, self.origin_nbintobas, self.routing_fetch = \
- routing_interface.truncate(nbasmax, hydrosuper.basin_count, hydrosuper.basin_notrun, hydrosuper.basin_area, \
- hydrosuper.basin_cg, hydrosuper.basin_topoind, hydrosuper.fetch_basin, hydrosuper.basin_id, \
- hydrosuper.basin_outcoor, hydrosuper.basin_type, hydrosuper.basin_flowdir, \
- hydrosuper.outflow_grid, hydrosuper.outflow_basin, \
- hydrosuper.inflow_number,hydrosuper.inflow_grid,hydrosuper.inflow_basin)
+ routing_interface.truncate(nbasmax, hydrosuper.num_largest, hydrosuper.basin_count, hydrosuper.basin_notrun, \
+ hydrosuper.basin_area, hydrosuper.basin_cg, \
+ hydrosuper.basin_topoind, hydrosuper.fetch_basin, hydrosuper.basin_id, \
+ hydrosuper.basin_outcoor, hydrosuper.basin_type, hydrosuper.basin_flowdir, \
+ hydrosuper.outflow_grid, hydrosuper.outflow_basin, \
+ hydrosuper.inflow_number,hydrosuper.inflow_grid,hydrosuper.inflow_basin)
+ yy=modelgrid.landscatter(np.sum(self.routing_fetch, axis=1)/np.sum(self.routing_area, axis=1))
+ print ("Rank :"+str(part.rank)+" OUT of truncate =+=+=+=+=+=+=+=+=+= \n"+str(yy)+"\n =+=+=+=+=+=+=+=+=+=")
return
#
def dumpnetcdf(self, filename, globalgrid, procgrid, part) :
diff --git a/RoutingPreProc.py b/RoutingPreProc.py
index 0c6440a04232db012adbe88aabd26d0723bc69a3..0bc49e6a870f5606a4ad5a891013ce84222f5062 100644
--- a/RoutingPreProc.py
+++ b/RoutingPreProc.py
@@ -172,9 +172,9 @@ del hoverlap
gc.collect()
print("=================== Compute fetch ====================")
-hsuper.fetch(part)
+hsuper.fetch(part, modelgrid)
-hgraph = IF.HydroGraph(nbasmax, hsuper)
+hgraph = IF.HydroGraph(nbasmax, hsuper, part, modelgrid)
hgraph.dumpnetcdf(OutGraphFile, gg, modelgrid, part)