diff --git a/DocumentationInterface b/DocumentationInterface
index 0a431bb3357b8b4035f82d5e536d4a6d29583352..0f090224dfd7d7d9ef726aaab2438a37b1347f6a 100644
--- a/DocumentationInterface
+++ b/DocumentationInterface
@@ -194,7 +194,7 @@ Parameters
 ----------
 corepts : input rank-1 array('i') with bounds (nbcore)
 basin_count : input rank-1 array('i') with bounds (nbpt)
-basin_area : in/output rank-2 array('f') with bounds (nbpt,nwbas)
+basin_area : input rank-2 array('f') with bounds (nbpt,nwbas)
 basin_id : input rank-2 array('i') with bounds (nbpt,nwbas)
 basin_hierarchy : input rank-2 array('f') with bounds (nbpt,nwbas)
 basin_fac : input rank-2 array('f') with bounds (nbpt,nwbas)
diff --git a/F90subroutines/routing_interface.f90 b/F90subroutines/routing_interface.f90
index 59886d09cbac7e542237e3070123ed38b297c46d..4248e89cb9ccef1df591313eef9e45a9a3cddce7 100644
--- a/F90subroutines/routing_interface.f90
+++ b/F90subroutines/routing_interface.f90
@@ -344,13 +344,13 @@ SUBROUTINE fetch(nbpt, nwbas, nbcore, corepts, basin_count, basin_area, basin_id
   INTEGER(i_std), INTENT (in) :: nbcore
   INTEGER(i_std), DIMENSION(nbcore), INTENT(in)        :: corepts
   INTEGER(i_std), DIMENSION(nbpt), INTENT(in)          :: basin_count !!
-  REAL(r_std), DIMENSION(nbpt,nwbas), INTENT(inout)    :: basin_area !!
+  REAL(r_std), DIMENSION(nbpt,nwbas), INTENT(in)       :: basin_area !!
   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(in)    :: 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
+  REAL(r_std), DIMENSION(nbpt,nwbas), INTENT(in)       :: partial_sum
   !
   !! IN-OUTPUT VARIABLES
   REAL(r_std), DIMENSION(nbpt,nwbas), INTENT(inout) :: fetch_basin !!
diff --git a/Interface.py b/Interface.py
index d6b215d7f21ca2a0006925004cb6d2a3b3bb5535..12ee03da7f21eae3cbcb8fca87a7c22d0c876fc0 100644
--- a/Interface.py
+++ b/Interface.py
@@ -166,12 +166,12 @@ class HydroSuper :
         # Order of magnitude for the area precision in m^2.
         prec = 100.0
         #
-        fetch_basin = np.zeros(self.basin_area.shape, dtype=np.float64, order='F')
+        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, self.outflow_grid)
-        partial_sum = np.zeros(self.basin_area.shape, dtype=np.float64, order='F')
+        partial_sum = np.zeros(self.basin_area.shape, dtype=np.float32, order='F')
         #
-        old_sorted = np.zeros(largest_pos, dtype=np.float64, order='F')
+        old_sorted = np.zeros(largest_pos, dtype=np.float32, order='F')
         #
         maxdiff_sorted = prec*prec
         iter_count = 0
@@ -196,7 +196,7 @@ class HydroSuper :
                 maxdiff_sorted = np.max(np.abs(sorted_outareas[0:largest_pos]-old_sorted))
             old_sorted[:] = sorted_outareas[0:largest_pos]
             iter_count += 1
-            
+
         self.fetch_basin = np.copy(fetch_basin)
         #
         # Upstream area of the smalest river we call largest rivers.