Note
Function to calculate the double derivatives 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) | :: | gradxx |
double derivative along x, x 2D array |
|
| real(kind=R8), | intent(out), | dimension(1:nx, 1:ny) | :: | gradyy |
double derivative along y, y 2D array |
subroutine curv2(tab, nx, ny, dx, dy, gradxx, gradyy) !================================================================================================ !< @note Function to calculate the double derivatives 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) :: gradxx !! *double derivative along x, x 2D array* real (kind=R8), intent(out), dimension(1:nx, 1:ny) :: gradyy !! *double derivative along y, y 2D array* integer(kind=I4) :: i, j !------------------------------------------------------------------ GRADXX i = 1 gradxx(1, 1:ny) = (UN / (180 * dx**2) ) * ( +0812 * tab(i+0, 1:ny) & ! -3132 * tab(i+1, 1:ny) & ! +5265 * tab(i+2, 1:ny) & ! -5080 * tab(i+3, 1:ny) & ! +2970 * tab(i+4, 1:ny) & ! -0972 * tab(i+5, 1:ny) & ! +0137 * tab(i+6, 1:ny) ) ! i = 1 gradxx(2, 1:ny) = (UN / (180 * dx**2) ) * ( +0137 * tab(i+0, 1:ny) & ! -0147 * tab(i+1, 1:ny) & ! -0255 * tab(i+2, 1:ny) & ! +0470 * tab(i+3, 1:ny) & ! -0285 * tab(i+4, 1:ny) & ! +0093 * tab(i+5, 1:ny) & ! -0013 * tab(i+6, 1:ny) ) ! i = 1 gradxx(3, 1:ny) = (UN / (180 * dx**2) ) * ( -0013 * tab(i+0, 1:ny) & ! +0228 * tab(i+1, 1:ny) & ! -0420 * tab(i+2, 1:ny) & ! +0200 * tab(i+3, 1:ny) & ! +0015 * tab(i+4, 1:ny) & ! -0012 * tab(i+5, 1:ny) & ! +0002 * tab(i+6, 1:ny) ) ! do i = 4, nx - 3 gradxx(i, 1:ny) = (UN / (180 * dx**2) ) * ( +002 * tab(i+3, 1:ny) & ! -027 * tab(i+2, 1:ny) & ! +270 * tab(i+1, 1:ny) & ! -490 * tab(i , 1:ny) & ! +270 * tab(i-1, 1:ny) & ! -027 * tab(i-2, 1:ny) & ! +002 * tab(i-3, 1:ny) ) ! enddo i = nx gradxx(nx , 1:ny) = (UN / (180 * dx**2) ) * ( +0812 * tab(i-0, 1:ny) & ! -3132 * tab(i-1, 1:ny) & ! +5265 * tab(i-2, 1:ny) & ! -5080 * tab(i-3, 1:ny) & ! +2970 * tab(i-4, 1:ny) & ! -0972 * tab(i-5, 1:ny) & ! +0137 * tab(i-6, 1:ny) ) ! i = nx gradxx(nx-1, 1:ny) = (UN / (180 * dx**2) ) * ( +0137 * tab(i-0, 1:ny) & ! -0147 * tab(i-1, 1:ny) & ! -0255 * tab(i-2, 1:ny) & ! +0470 * tab(i-3, 1:ny) & ! -0285 * tab(i-4, 1:ny) & ! +0093 * tab(i-5, 1:ny) & ! -0013 * tab(i-6, 1:ny) ) ! i = nx gradxx(nx-2, 1:ny) = (UN / (180 * dx**2) ) * ( -0013 * tab(i-0, 1:ny) & ! +0228 * tab(i-1, 1:ny) & ! -0420 * tab(i-2, 1:ny) & ! +0200 * tab(i-3, 1:ny) & ! +0015 * tab(i-4, 1:ny) & ! -0012 * tab(i-5, 1:ny) & ! +0002 * tab(i-6, 1:ny) ) ! !------------------------------------------------------------------ GRADYY j = 1 gradyy(1:nx, 1) = (UN / (180 * dy**2) ) * ( +0812 * tab(1:nx, j+0) & ! -3132 * tab(1:nx, j+1) & ! +5265 * tab(1:nx, j+2) & ! -5080 * tab(1:nx, j+3) & ! +2970 * tab(1:nx, j+4) & ! -0972 * tab(1:nx, j+5) & ! +0137 * tab(1:nx, j+6) ) ! j = 1 gradyy(1:nx, 2) = (UN / (180 * dy**2) ) * ( +0137 * tab(1:nx, j+0) & ! -0147 * tab(1:nx, j+1) & ! -0255 * tab(1:nx, j+2) & ! +0470 * tab(1:nx, j+3) & ! -0285 * tab(1:nx, j+4) & ! +0093 * tab(1:nx, j+5) & ! -0013 * tab(1:nx, j+6) ) ! j = 1 gradyy(1:nx, 3) = (UN / (180 * dy**2) ) * ( -0013 * tab(1:nx, j+0) & ! +0228 * tab(1:nx, j+1) & ! -0420 * tab(1:nx, j+2) & ! +0200 * tab(1:nx, j+3) & ! +0015 * tab(1:nx, j+4) & ! -0012 * tab(1:nx, j+5) & ! +0002 * tab(1:nx, j+6) ) ! do j = 4, ny - 3 gradyy(1:nx, j) = (UN / (180 * dy**2) ) * ( +002 * tab(1:nx, j+3) & ! -027 * tab(1:nx, j+2) & ! +270 * tab(1:nx, j+1) & ! -490 * tab(1:nx, j ) & ! +270 * tab(1:nx, j-1) & ! -027 * tab(1:nx, j-2) & ! +002 * tab(1:nx, j-3) ) ! enddo j = ny gradyy(1:nx, ny ) = (UN / (180 * dy**2) ) * ( +0812 * tab(1:nx, j-0) & ! -3132 * tab(1:nx, j-1) & ! +5265 * tab(1:nx, j-2) & ! -5080 * tab(1:nx, j-3) & ! +2970 * tab(1:nx, j-4) & ! -0972 * tab(1:nx, j-5) & ! +0137 * tab(1:nx, j-6) ) ! j = ny gradyy(1:nx, ny-1) = (UN / (180 * dy**2) ) * ( +0137 * tab(1:nx, j-0) & ! -0147 * tab(1:nx, j-1) & ! -0255 * tab(1:nx, j-2) & ! +0470 * tab(1:nx, j-3) & ! -0285 * tab(1:nx, j-4) & ! +0093 * tab(1:nx, j-5) & ! -0013 * tab(1:nx, j-6) ) ! j = ny gradyy(1:nx, ny-2) = (UN / (180 * dy**2) ) * ( -0013 * tab(1:nx, j-0) & ! +0228 * tab(1:nx, j-1) & ! -0420 * tab(1:nx, j-2) & ! +0200 * tab(1:nx, j-3) & ! +0015 * tab(1:nx, j-4) & ! -0012 * tab(1:nx, j-5) & ! +0002 * tab(1:nx, j-6) ) ! return endsubroutine curv2