Note
Function to calculate the gradient of a 2D array
It implements the details given in ISO 25178.
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| real(kind=R8), | intent(in), | dimension(1:nx, 1:ny) | :: | tab |
Input 2D array |
|
| integer(kind=I4), | intent(in) | :: | nx |
number of pixels along x |
||
| integer(kind=I4), | intent(in) | :: | ny |
number of pixels along x |
||
| real(kind=R8), | intent(in) | :: | dx |
x lag |
||
| real(kind=R8), | intent(in) | :: | dy |
y lag |
||
| real(kind=R8), | intent(out), | dimension(1:nx, 1:ny) | :: | gradx |
derivative along x 2D array |
|
| real(kind=R8), | intent(out), | dimension(1:nx, 1:ny) | :: | grady |
derivative along y 2D array |
subroutine gradient(tab, nx, ny, dx, dy, gradx, grady) !================================================================================================ !< @note Function to calculate the gradient of a 2D array !< !< It implements the details given in ISO 25178. !< !< @endnote !------------------------------------------------------------------------------------------------ implicit none integer(kind=I4), intent(in ) :: nx !! *number of pixels along x* integer(kind=I4), intent(in ) :: ny !! *number of pixels along x* real (kind=R8), intent(in ) :: dx !! *x lag* real (kind=R8), intent(in ) :: dy !! *y lag* real (kind=R8), intent(in ), dimension(1:nx, 1:ny) :: tab !! *Input 2D array* real (kind=R8), intent(out), dimension(1:nx, 1:ny) :: gradx !! *derivative along x 2D array* real (kind=R8), intent(out), dimension(1:nx, 1:ny) :: grady !! *derivative along y 2D array* integer(kind=I4) :: i, j !------------------------------------------------------------------ GRADX i = 1 gradx(1, 1:ny) = (UN / (60 * dx) ) * ( -147 * tab(i+0, 1:ny) & ! +360 * tab(i+1, 1:ny) & ! -450 * tab(i+2, 1:ny) & ! +400 * tab(i+3, 1:ny) & ! -225 * tab(i+4, 1:ny) & ! +072 * tab(i+5, 1:ny) & ! -010 * tab(i+6, 1:ny) ) ! i = 1 gradx(2, 1:ny) = (UN / (60 * dx) ) * ( -010 * tab(i+0, 1:ny) & ! -077 * tab(i+1, 1:ny) & ! +150 * tab(i+2, 1:ny) & ! -100 * tab(i+3, 1:ny) & ! +050 * tab(i+4, 1:ny) & ! -015 * tab(i+5, 1:ny) & ! +002 * tab(i+6, 1:ny) ) ! i = 1 gradx(3, 1:ny) = (UN / (60 * dx) ) * ( +002 * tab(i+0, 1:ny) & ! -024 * tab(i+1, 1:ny) & ! -035 * tab(i+2, 1:ny) & ! +080 * tab(i+3, 1:ny) & ! -030 * tab(i+4, 1:ny) & ! +008 * tab(i+5, 1:ny) & ! -001 * tab(i+6, 1:ny) ) ! do i = 4, nx - 3 gradx(i, 1:ny) = (UN / (60 * dx) ) * ( +01 * tab(i+3, 1:ny) & ! -09 * tab(i+2, 1:ny) & ! +45 * tab(i+1, 1:ny) & ! -45 * tab(i-1, 1:ny) & ! +09 * tab(i-2, 1:ny) & ! -01 * tab(i-3, 1:ny) ) ! enddo i = nx gradx(nx , 1:ny) = -(UN / (60 * dx) ) * ( -147 * tab(i-0, 1:ny) & ! +360 * tab(i-1, 1:ny) & ! -450 * tab(i-2, 1:ny) & ! +400 * tab(i-3, 1:ny) & ! -225 * tab(i-4, 1:ny) & ! +072 * tab(i-5, 1:ny) & ! -010 * tab(i-6, 1:ny) ) ! i = nx gradx(nx-1, 1:ny) = -(UN / (60 * dx) ) * ( -010 * tab(i-0, 1:ny) & ! -077 * tab(i-1, 1:ny) & ! +150 * tab(i-2, 1:ny) & ! -100 * tab(i-3, 1:ny) & ! +050 * tab(i-4, 1:ny) & ! -015 * tab(i-5, 1:ny) & ! +002 * tab(i-6, 1:ny) ) ! i = nx gradx(nx-2, 1:ny) = -(UN / (60 * dx) ) * ( +002 * tab(i-0, 1:ny) & ! -024 * tab(i-1, 1:ny) & ! -035 * tab(i-2, 1:ny) & ! +080 * tab(i-3, 1:ny) & ! -030 * tab(i-4, 1:ny) & ! +008 * tab(i-5, 1:ny) & ! -001 * tab(i-6, 1:ny) ) ! !------------------------------------------------------------------ GRADY j = 1 grady(1:nx, 1) = (UN / (60 * dy) ) * ( -147 * tab(1:nx, j+0) & ! +360 * tab(1:nx, j+1) & ! -450 * tab(1:nx, j+2) & ! +400 * tab(1:nx, j+3) & ! -225 * tab(1:nx, j+4) & ! +072 * tab(1:nx, j+5) & ! -010 * tab(1:nx, j+6) ) ! j = 1 grady(1:nx, 2) = (UN / (60 * dy) ) * ( -010 * tab(1:nx, j+0) & ! -077 * tab(1:nx, j+1) & ! +150 * tab(1:nx, j+2) & ! -100 * tab(1:nx, j+3) & ! +050 * tab(1:nx, j+4) & ! -015 * tab(1:nx, j+5) & ! +002 * tab(1:nx, j+6) ) ! j = 1 grady(1:nx, 3) = (UN / (60 * dy) ) * ( +002 * tab(1:nx, j+0) & ! -024 * tab(1:nx, j+1) & ! -035 * tab(1:nx, j+2) & ! +080 * tab(1:nx, j+3) & ! -030 * tab(1:nx, j+4) & ! +008 * tab(1:nx, j+5) & ! -001 * tab(1:nx, j+6) ) ! do j = 4, ny - 3 grady(1:nx, j) = (UN / (60 * dy) ) * ( +01 * tab(1:nx, j+3) & ! -09 * tab(1:nx, j+2) & ! +45 * tab(1:nx, j+1) & ! -45 * tab(1:nx, j-1) & ! +09 * tab(1:nx, j-2) & ! -01 * tab(1:nx, j-3) ) ! enddo j = ny grady(1:nx, ny ) = -(UN / (60 * dy) ) * ( -147 * tab(1:nx, j-0) & ! +360 * tab(1:nx, j-1) & ! -450 * tab(1:nx, j-2) & ! +400 * tab(1:nx, j-3) & ! -225 * tab(1:nx, j-4) & ! +072 * tab(1:nx, j-5) & ! -010 * tab(1:nx, j-6) ) ! j = ny grady(1:nx, ny-1) = -(UN / (60 * dy) ) * ( -010 * tab(1:nx, j-0) & ! -077 * tab(1:nx, j-1) & ! +150 * tab(1:nx, j-2) & ! -100 * tab(1:nx, j-3) & ! +050 * tab(1:nx, j-4) & ! -015 * tab(1:nx, j-5) & ! +002 * tab(1:nx, j-6) ) ! j = ny grady(1:nx, ny-2) = -(UN / (60 * dy) ) * ( +002 * tab(1:nx, j-0) & ! -024 * tab(1:nx, j-1) & ! -035 * tab(1:nx, j-2) & ! +080 * tab(1:nx, j-3) & ! -030 * tab(1:nx, j-4) & ! +008 * tab(1:nx, j-5) & ! -001 * tab(1:nx, j-6) ) ! return endsubroutine gradient