ReorderKernelForGFT.f90

Go to the documentation of this file.

!==============================================================================
!==============================================================================
Subroutine shiftkerneltoeven(ispec, y_dim_kernel_orig,x_dim_kernel_orig, y_dim_kernel,x_dim_kernel,kernel_orig,kernel_shifted) 
    USE gft, only: gft_prec
    IMPLICIT NONE
    
    INTEGER, INTENT(IN):: y_dim_kernel_orig,x_dim_kernel_orig, y_dim_kernel,x_dim_kernel, ispec
    REAL(kind=gft_prec) , DIMENSION(1: y_dim_kernel_orig ,1: x_dim_kernel_orig), INTENT(IN):: kernel_orig
    
    REAL(kind=gft_prec) , DIMENSION(0: y_dim_kernel-1 ,x_dim_kernel), INTENT(OUT):: kernel_shifted
    
    
    REAL(kind=gft_prec),  DIMENSION(1: y_dim_kernel ,x_dim_kernel)   :: kernel_help   
    INTEGER:: y_dim_kernel2, x_dim_kernel2 
    
    
    kernel_shifted=0
    
    if (mod(y_dim_kernel_orig, 2)==0 .and. mod(x_dim_kernel_orig,2) == 0)  THEN ! if equal then exact the same kernel is stored, just shifted by -1 in y dir.
        kernel_shifted(0: y_dim_kernel-1,1:x_dim_kernel)=kernel_orig
    else  ! the maximum value in the center is distributed to the four center cells.
        y_dim_kernel2=int(y_dim_kernel_orig/2) ! half the original dimension
        x_dim_kernel2=int(x_dim_kernel_orig/2) ! half the original dimension
        kernel_help=0
        kernel_help( 1:y_dim_kernel2     , 1:x_dim_kernel_orig)=kernel_orig(1:(y_dim_kernel2 ),:) 
        
        kernel_help((y_dim_kernel2+1)    , 1:x_dim_kernel_orig)=kernel_orig((y_dim_kernel2  +1),:)!/2
        kernel_help((y_dim_kernel2+2)    , 1:x_dim_kernel_orig)=kernel_orig((y_dim_kernel2  +1),:)!/2
        kernel_help((y_dim_kernel2+3):y_dim_kernel, 1:x_dim_kernel_orig)=kernel_orig((y_dim_kernel2+2):y_dim_kernel_orig,:)
        kernel_shifted(0:x_dim_kernel-1, 1:x_dim_kernel2)                = kernel_help(:, 1      :x_dim_kernel2 )
        kernel_shifted(0:x_dim_kernel-1,(x_dim_kernel2+1))               = kernel_help(:,(x_dim_kernel2+1) ) ! /2
        kernel_shifted(0:x_dim_kernel-1,(x_dim_kernel2+2))               = kernel_help(:,(x_dim_kernel2+1) ) !/2
        kernel_shifted(0:x_dim_kernel-1,(x_dim_kernel2+3):x_dim_kernel ) = kernel_help(:,(x_dim_kernel2+2): x_dim_kernel_orig )
    END IF
    return 
end subroutine shiftkerneltoeven
 
!==============================================================================
!==============================================================================
SUBROUTINE reorderkernelforgft(dim_y,dim_x, shifted_kernel, reorderedkernel)
    USE gft
    IMPLICIT NONE
    INTEGER, INTENT(IN):: dim_y,dim_x 
    REAL(kind=gft_prec), DIMENSION(0:(dim_y-1),dim_x), INTENT(IN)::  shifted_kernel
    
    REAL(kind=gft_prec), DIMENSION(0:(dim_y-1),dim_x), INTENT(OUT):: reorderedkernel
    INTEGER:: center_dim_y, center_dim_x,i1,j1,i,j
    REAL:: dim_y2,dim_x2,sumReOrderedKAr
    
    
    !   call PrintMatrixASXY2File(shifted_kernel,    0,dim_y-1,1,dim_x,"kernel_orig_shifted")
    
    ! Order the kernel to the four corners
    dim_y2=dim_y/2.
    dim_x2=dim_x/2. 
    
    center_dim_y=anint( dim_y2)
    center_dim_x=anint( dim_x2)
 
    reorderedkernel(:,:)=0.
    do i= 0,dim_y-1
        i1= mod((i+ 1 +center_dim_y),dim_y)         ! this is the symmetrical version   
        do j= 1,dim_x
            j1 = mod((j-1+center_dim_x),dim_x)  + 1   ! the -1 makes the reordered kernel symmetrical   
            reorderedkernel(i ,j ) = shifted_kernel(i1,j1)
        end do ! j
    end do ! i
    sumreorderedkar = sum(reorderedkernel)
    reorderedkernel=reorderedkernel/sumreorderedkar ! to make sure that the kernel has a sum of 1
  
    !   call PrintMatrixASXY2File (reorderedkernel,0,dim_y-1,1,dim_x,"kernel_reordered" )
    return
END SUBROUTINE reorderkernelforgft
 
!===================================================================== 
!===============================================================
REAL FUNCTION distance (x,x0,y,y0 )
IMPLICIT NONE
REAL, INTENT(in) :: x,x0,y,y0 
distance=((  x- x0)**2+( y-y0)**2)**0.5     
END FUNCTION distance