diff --git a/F90subroutines/routing_reg.f90 b/F90subroutines/routing_reg.f90
index 832fbfebbf86110e6cd8967d7b8bacc336bad0af..62aa6eef9addc9140041dd3585f618cdc9e623af 100644
--- a/F90subroutines/routing_reg.f90
+++ b/F90subroutines/routing_reg.f90
@@ -2073,10 +2073,6 @@ SUBROUTINE routing_reg_linkup(nbpt, neighbours, nwbas, ijdimmax, inflowmax, basi
! At this point any of the outflows is designated by a negative values in
! outflow_grid
!
- WRITE(numout,*) "*****"
- WRITE(numout,*) "Linkup 0 - sp, sb = ", sp, sb
- WRITE(numout,*) "Linkup 0 - outflow_grid = ", outflow_grid(sp,sb)
- !
found = 0
IF ( outflow_grid(sp,sb) == 0 ) THEN
found = 1
@@ -2130,14 +2126,14 @@ SUBROUTINE routing_reg_linkup(nbpt, neighbours, nwbas, ijdimmax, inflowmax, basi
& nbpt, nwbas, inp, basin_count, basin_id, basin_hierarchy, basin_fac, &
& basin_flowdir, nbcoastal, coastal_basin, bop, bopqual)
!
- IF ( bop .LT. undef_int ) THEN
+ IF ( bop .LT. undef_int .AND. bop .NE. sb) THEN
!
CALL routing_updateflow(sp, sb, inp, bop, nbpt, nwbas, inflowmax, outflow_grid, outflow_basin, &
& inflow_number, inflow_grid, inflow_basin)
IF ( outflow_basin(sp,sb) == bop ) THEN
solved(sp,1) = solved(sp,1) + 1
found = 1
- WRITE(numout,*) "Solution found in the original outflow_grid"
+ WRITE(numout,*) "Solution found in the original outflow_grid", sb, bop
ENDIF
!
ENDIF
@@ -2190,7 +2186,7 @@ SUBROUTINE routing_reg_linkup(nbpt, neighbours, nwbas, ijdimmax, inflowmax, basi
! we know it's wrong but it will serve to make decision later
angp1 = routing_anglediff_g(sp, dp1, basin_flowdir(sp,sb))
! Check if grid point
- IF (dp1 .GT. 0) THEN
+ IF (dp1 .GT. 0 .AND. dp1 .NE. sp) THEN
CALL routing_reg_bestsubbasin(sp, sb, basin_id(sp,sb), basin_hierarchy(sp,sb), basin_fac(sp,sb), &
& basin_flowdir(sp,sb), invented_basins, &
& nbpt, nwbas, dp1, basin_count, basin_id, basin_hierarchy, basin_fac, &
@@ -2205,7 +2201,7 @@ SUBROUTINE routing_reg_linkup(nbpt, neighbours, nwbas, ijdimmax, inflowmax, basi
dm1 = neighbours_g(sp,order_ref(nb*2+1))
! we know it's wrong but it will serve to make decision later
angm1 = routing_anglediff_g(sp, dm1, basin_flowdir(sp,sb))
- IF (dm1 .GT. 0) THEN
+ IF (dm1 .GT. 0 .AND. dm1 .NE. sp) THEN
CALL routing_reg_bestsubbasin(sp, sb, basin_id(sp,sb), basin_hierarchy(sp,sb), basin_fac(sp,sb), &
& basin_flowdir(sp,sb), invented_basins, &
& nbpt, nwbas, dm1, basin_count, basin_id, basin_hierarchy, basin_fac, &
@@ -2256,7 +2252,7 @@ SUBROUTINE routing_reg_linkup(nbpt, neighbours, nwbas, ijdimmax, inflowmax, basi
IF ( outflow_basin(sp,sb) == bop ) THEN
solved(sp,2) = solved(sp,2) + 1
found = 1
- WRITE (numout,*) "Neighbours, output found at:" , nb, "level"
+ WRITE (numout,*) "Neighbours, output found at:" , nb, "level. sp,sb=", sp, sb, " dop,bop=", dop, bop
ELSE
nb = nb+1
ENDIF
@@ -2340,7 +2336,8 @@ SUBROUTINE routing_reg_linkup(nbpt, neighbours, nwbas, ijdimmax, inflowmax, basi
IF (( neighbours_g(sp,order_ref(nb)) == minhloc(basin_id(sp,sb),1)) .AND. (found == 0)) THEN
dop = minhloc(basin_id(sp,sb),1)
bop = minhloc(basin_id(sp,sb),2)
- IF ((dop < undef_int ) .AND. (dop .GT. 0) .AND. (bop < undef_int) .AND. (bop .GT. 0)) THEN
+ IF ((dop < undef_int ) .AND. (dop .GT. 0) .AND. (bop < undef_int) .AND. (bop .GT. 0) &
+ & .AND. (dop .NE. sp) .AND. (bop .NE. sb) ) THEN
CALL routing_updateflow(sp, sb, dop, bop, nbpt,nwbas, inflowmax, outflow_grid, outflow_basin, &
& inflow_number, inflow_grid, inflow_basin)
! It is possible that the lowest hierarchy is in two grid cells
@@ -2352,6 +2349,7 @@ SUBROUTINE routing_reg_linkup(nbpt, neighbours, nwbas, ijdimmax, inflowmax, basi
ELSE
WRITE(numout,*) "Lowest hierarchy may be in two different grid points"
END IF
+ WRITE(numout,*) "sp,sb = ", sp,sb, " dop,bop = ", dop,bop
WRITE(numout,*) "hierarch of (sp,sb)", basin_hierarchy(sp,sb)
WRITE(numout,*) "Lowest basinid hierarch", basin_hierarchy(dop,bop)
END IF
diff --git a/ModelGrid.py b/ModelGrid.py
index 7661a43ff69ba8ea0ade4f465cb50999a3f04563..81b7a4c35bb01720219b9f2d0fac9a2af5184c1f 100644
--- a/ModelGrid.py
+++ b/ModelGrid.py
@@ -183,6 +183,8 @@ def getland (geo, ist, ni, jst, nj) :
vn=list(v.name for v in geo.variables.values())
if "LANDMASK" in vn :
land=geo.variables["LANDMASK"][0,jst:jst+nj,ist:ist+ni]
+ elif "Contfrac" in vn :
+ land=geo.variables["Contfrac"][jst:jst+nj,ist:ist+ni]
elif "elevation" in vn :
land=geo.variables["elevation"][jst:jst+nj,ist:ist+ni]
if "missing_value" in geo.variables["elevation"].ncattrs() :
diff --git a/Partition.py b/Partition.py
index 93165133451ee5dce7e52033e2a744c98f409f49..99e04a458a5722bf2dc58961fe40f13cfc641f4c 100644
--- a/Partition.py
+++ b/Partition.py
@@ -52,7 +52,7 @@ def fit_partition(partin, land):
for i, dom in enumerate(partin):
l = land[dom["jstart"]:dom["jstart"]+dom["nbj"],dom["istart"]:dom["istart"]+dom["nbi"]]
l1 = np.sum(l, axis = 0)
-
+
i0 = np.where(l1>0)[0][0]
i1 = np.where(l1>0)[0][-1]
diff --git a/tests/Mallorca/BuildHTUs_Mallorca.pbs b/tests/Mallorca/BuildHTUs_Mallorca.pbs
index 86964522b7c7b7912e9aeeedfe90dd77ba6271be..83428289aae04a0081ae2d0d637156847ee716b1 100644
--- a/tests/Mallorca/BuildHTUs_Mallorca.pbs
+++ b/tests/Mallorca/BuildHTUs_Mallorca.pbs
@@ -29,6 +29,7 @@ if [ $? -gt 0 ] ; then
echo "XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX"
echo "X Run MECORDEX on 2 Proc failed X"
echo "XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX"
+ exit
else
echo "======================================"
echo "= Run MEDCORDEX on 2 Proc successful ="
@@ -36,9 +37,27 @@ else
ls -lt
fi
#
+#
+# 4 Proc
+#
+/bin/rm -f run.def *.txt
+cp run_E2OFD.def run.def
+mpirun -n 2 python ../../RoutingPreProc.py
+if [ $? -gt 0 ] ; then
+ echo "XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX"
+ echo "X Run E2OFD on 4 Proc failed X"
+ echo "XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX"
+ exit
+else
+ echo "=================================="
+ echo "= Run E2OFD on 4 Proc successful ="
+ echo "=================================="
+ ls -lt
+fi
+#
# More points on 32 proc
#
-/bin/rm -f run.def
+/bin/rm -f run.def *.txt
cp run_EuroSW.def run.def
mpirun -n 32 python ../../RoutingPreProc.py
if [ $? -gt 0 ] ; then
diff --git a/tests/Mallorca/run_E2OFD.def b/tests/Mallorca/run_E2OFD.def
new file mode 100644
index 0000000000000000000000000000000000000000..544545cad9fcef02f7c037be7f28da6d86f37efa
--- /dev/null
+++ b/tests/Mallorca/run_E2OFD.def
@@ -0,0 +1,27 @@
+[OverAll]
+#
+#
+EarthRadius = 6370000.
+#
+ModelGridFile = /bdd/MEDI/workspaces/polcher/NewRouting/MED_E2OFD_Contfrac.nc
+HydroFile = /bdd/MEDI/workspaces/polcher/NewRouting/routing_MED.nc
+# Mallorca
+WEST_EAST = 2.3, 3.75
+SOUTH_NORTH = 39.00, 40.1
+#
+# FORTRAN interface parameters
+#
+Documentation = true
+#
+# Configuration for the graph to be generated
+#
+nbasmax = 35
+#
+# Number of operation of simplification performed together
+#
+numop = 200
+#
+# Output
+#
+GraphFile = E2OFD_test_graph.nc
+#
diff --git a/tests/Mallorca/run_EuroSW.def b/tests/Mallorca/run_EuroSW.def
index 2fba5f8f6b15f8a2e090529162caf636b4829277..d0186a0e682e47b296fe4ea7454eadb7debb49a3 100644
--- a/tests/Mallorca/run_EuroSW.def
+++ b/tests/Mallorca/run_EuroSW.def
@@ -10,7 +10,7 @@ HydroFile = /bdd/MEDI/workspaces/polcher/NewRouting/routing_MED.nc
WEST_EAST = 2.2, 3.6
SOUTH_NORTH = 39.20, 40.10
#
-##WeightFile = EuroSW_Weights.nc
+WeightFile = EuroSW_Weights.nc
#
# FORTRAN interface parameters
#