Index: /trunk/tools/module_scientific.f90
===================================================================
--- /trunk/tools/module_scientific.f90	(revision 2281)
+++ /trunk/tools/module_scientific.f90	(revision 2282)
@@ -28,4 +28,7 @@
 ! look_clockwise_borders: Subroutine to look clock-wise for a next point within a collection of borders 
 !   (limits of a region)
+! multi_index_mat2DRK: Subroutine to provide the indices of the different locations of a value inside a 2D RK matrix
+! multi_index_mat3DRK: Subroutine to provide the indices of the different locations of a value inside a 3D RK matrix
+! multi_index_mat4DRK: Subroutine to provide the indices of the different locations of a value inside a 4D RK matrix
 ! NcountR: Subroutine to count real values
 ! paths_border: Subroutine to search the paths of a border field.
@@ -5208,3 +5211,178 @@
   END SUBROUTINE spaceweightstats
 
+
+  SUBROUTINE multi_index_mat2DRK(d1, d2, d12, mat, value, Nindices, indices)
+  ! Subroutine to provide the indices of the different locations of a value inside a 2D RK matrix
+
+    IMPLICIT NONE
+
+    INTEGER, INTENT(in)                                  :: d1, d2, d12
+    REAL(r_k), DIMENSION(d1,d2), INTENT(in)              :: mat
+    REAL(r_k),INTENT(in)                                 :: value
+    INTEGER, INTENT(out)                                 :: Nindices
+    INTEGER, DIMENSION(2,d12), INTENT(out)               :: indices
+
+! Local
+    INTEGER                                              :: i,j
+    INTEGER                                              :: Ncounts1D, icount1D
+    REAL(r_k), DIMENSION(d2)                             :: diffmat1D
+
+    !!!!!!! Variables
+    ! d1, d2: shape of the 2D matrix
+    ! mat: 2D matrix
+    ! value: value to be looking for
+    ! Nindices: number of times value found within matrix
+    ! indices: indices of the found values
+
+    fname = 'multi_index_mat2DRK'
+
+    Nindices = 0
+    indices = 0
+    DO i=1, d1
+      diffmat1D = mat(i,:) - value
+      IF (ANY(diffmat1D == zeroRK)) THEN
+        Ncounts1D = COUNT(diffmat1D == zeroRK)
+        icount1D = 0
+        DO j=1, d2
+          IF (diffmat1D(j) == zeroRK) THEN
+            Nindices = Nindices + 1
+            indices(1,Nindices) = i
+            indices(2,Nindices) = j
+            icount1D = icount1D + 1
+            IF (icount1D == Ncounts1D) EXIT
+          END IF
+        END DO
+      END IF
+    END DO
+
+  END SUBROUTINE multi_index_mat2DRK
+
+  SUBROUTINE multi_index_mat3DRK(d1, d2, d3, d123, mat, value, Nindices, indices)
+  ! Subroutine to provide the indices of the different locations of a value inside a 3D RK matrix
+
+    IMPLICIT NONE
+
+    INTEGER, INTENT(in)                                  :: d1, d2, d3, d123
+    REAL(r_k), DIMENSION(d1,d2,d3), INTENT(in)           :: mat
+    REAL(r_k),INTENT(in)                                 :: value
+    INTEGER, INTENT(out)                                 :: Nindices
+    INTEGER, DIMENSION(3,d123), INTENT(out)              :: indices
+
+! Local
+    INTEGER                                              :: i,j,k
+    INTEGER                                              :: Ncounts1D, icount1D
+    INTEGER                                              :: Ncounts2D, icount2D
+    REAL(r_k), DIMENSION(d2,d3)                          :: diffmat2D
+    REAL(r_k), DIMENSION(d3)                             :: diffmat1D
+
+    !!!!!!! Variables
+    ! d1, d2, d3: shape of the 3D matrix
+    ! mat: 3D matrix
+    ! value: value to be looking for
+    ! Nindices: number of times value found within matrix
+    ! indices: indices of the found values
+
+    fname = 'multi_index_mat3DRK'
+
+    Nindices = 0
+    indices = 0
+    DO i=1, d1
+      diffmat2D = mat(i,:,:) - value
+      IF (ANY(diffmat2D == zeroRK)) THEN
+        Ncounts2D = COUNT(diffmat2D == zeroRK)
+        icount2D = 0
+        DO j=1, d2
+          diffmat1D = mat(i,j,:) - value
+          IF (ANY(diffmat1D == zeroRK)) THEN
+            Ncounts1D = COUNT(diffmat1D == zeroRK)
+            icount1D = 0
+            DO k=1, d3
+              IF (diffmat1D(k) == zeroRK) THEN
+                Nindices = Nindices + 1
+                indices(1,Nindices) = i
+                indices(2,Nindices) = j
+                indices(3,Nindices) = k
+                icount1D = icount1D + 1
+                IF (icount1D == Ncounts1D) EXIT
+              END IF
+            END DO
+            icount2D = icount2D + icount1D
+            IF (icount2D == Ncounts2D) EXIT
+          END IF
+        END DO
+      END IF
+    END DO
+
+  END SUBROUTINE multi_index_mat3DRK
+
+  SUBROUTINE multi_index_mat4DRK(d1, d2, d3, d4, d1234, mat, value, Nindices, indices)
+  ! Subroutine to provide the indices of the different locations of a value inside a 4D RK matrix
+
+    IMPLICIT NONE
+
+    INTEGER, INTENT(in)                                  :: d1, d2, d3, d4, d1234
+    REAL(r_k), DIMENSION(d1,d2,d3,d4), INTENT(in)        :: mat
+    REAL(r_k),INTENT(in)                                 :: value
+    INTEGER, INTENT(out)                                 :: Nindices
+    INTEGER, DIMENSION(4,d1234), INTENT(out)             :: indices
+
+! Local
+    INTEGER                                              :: i,j,k,l
+    INTEGER                                              :: Ncounts1D, icount1D
+    INTEGER                                              :: Ncounts2D, icount2D
+    INTEGER                                              :: Ncounts3D, icount3D
+    REAL(r_k), DIMENSION(d2,d3,d4)                       :: diffmat3D
+    REAL(r_k), DIMENSION(d3,d4)                          :: diffmat2D
+    REAL(r_k), DIMENSION(d4)                             :: diffmat1D
+
+    !!!!!!! Variables
+    ! d1, d2, d3, d4: shape of the 4D matrix
+    ! mat: 4D matrix
+    ! value: value to be looking for
+    ! Nindices: number of times value found within matrix
+    ! indices: indices of the found values
+
+    fname = 'multi_index_mat4DRK'
+
+    Nindices = 0
+    indices = 0
+    DO i=1, d1
+      diffmat3D = mat(i,:,:,:) - value
+      IF (ANY(diffmat3D == zeroRK)) THEN
+        Ncounts3D = COUNT(diffmat3D == zeroRK)
+        icount3D = 0
+        DO j=1, d2
+          diffmat2D = mat(i,j,:,:) - value
+          IF (ANY(diffmat2D == zeroRK)) THEN
+            Ncounts2D = COUNT(diffmat2D == zeroRK)
+            icount2D = 0
+            DO k=1, d3
+              diffmat1D = mat(i,j,k,:) - value
+              IF (ANY(diffmat1D == zeroRK)) THEN
+                Ncounts1D = COUNT(diffmat1D == zeroRK)
+                icount1D = 0
+                DO l=1, d4
+                  IF (diffmat1D(l) == zeroRK) THEN
+                    Nindices = Nindices + 1
+                    indices(1,Nindices) = i
+                    indices(2,Nindices) = j
+                    indices(3,Nindices) = k
+                    indices(4,Nindices) = l
+                    icount1D = icount1D + 1
+                    IF (icount1D == Ncounts1D) EXIT
+                  END IF
+                END DO
+              icount2D = icount2D + icount1D
+              IF (icount2D == Ncounts2D) EXIT
+              END IF
+            END DO
+            icount3D = icount3D + icount1D
+            IF (icount3D == Ncounts3D) EXIT
+          END IF
+        END DO
+      END IF
+    END DO
+
+  END SUBROUTINE multi_index_mat4DRK
+
 END MODULE module_scientific
