Fortran Calculation of Mask.

F90subroutines/diagst.f90

+MODULE diagst 
+
+CONTAINS
+
+SUBROUTINE upstmask(nbpt, nbasmax, routetogrid, routetobasin,&
+          & routenbbasin,basin_area, area,  g0,b0, mask)
+
+    IMPLICIT NONE
+    !! INPUT
+    INTEGER(4), INTENT(in) :: nbpt, nbasmax
+    !
+    INTEGER(4), DIMENSION(nbpt, nbasmax), INTENT(in)   :: routetogrid
+    INTEGER(4), DIMENSION(nbpt, nbasmax), INTENT(in)   :: routetobasin
+    INTEGER(4), DIMENSION(nbpt), INTENT(in)            :: routenbbasin
+    !
+    REAL(8), DIMENSION(nbpt, nbasmax), INTENT(in)      :: basin_area
+    REAL(8), DIMENSION(nbpt), INTENT(in)               :: area
+    !
+    INTEGER(4), INTENT(in)                             :: g0,b0
+    !
+    !! LOCAL
+    INTEGER(4)                                         :: ig,ib,jg,jb,jt,nbas
+    !
+    !! OUTPUT
+    REAL(8), DIMENSION(nbpt), INTENT(out)           ::  mask
+    !!
+    !!
+    mask(:) = 0
+
+    DO ig =1,nbpt
+      !WRITE(*,*) ig, "/",nbpt
+      nbas = routenbbasin(ig)
+      DO ib=1,nbas      
+         jg = ig
+         jb = ib
+         DO WHILE (jb .LE. nbasmax)
+            ! CHECK if mask > 0
+            ! valable si j'utilise mask ig,ib
+            IF (((jg .EQ. g0) .AND. (jb .EQ. b0))) THEN
+               ! add area as a variable
+               mask(ig) = mask(ig) + basin_area(ig,ib)/area(ig)
+               jb = nbasmax+1
+            ELSE
+               jt = routetogrid(jg,jb)
+               jb = routetobasin(jg,jb)
+               jg = jt
+            END IF
+         END DO
+      END DO
+    END DO
+    
+END SUBROUTINE upstmask
+
+END MODULE
+