InterpolateReflectionCoefficient Subroutine

public subroutine InterpolateReflectionCoefficient(RInt, R, NPts, PRTFile)

Arguments

Type IntentOptional Attributes Name
type(ReflectionCoef), intent(inout) :: RInt
type(ReflectionCoef), intent(in) :: R(NPts)
integer, intent(in) :: NPts
integer, intent(in) :: PRTFile

Called by

proc~~interpolatereflectioncoefficient~~CalledByGraph proc~interpolatereflectioncoefficient InterpolateReflectionCoefficient proc~reflect2d Reflect2D proc~reflect2d->proc~interpolatereflectioncoefficient proc~reflect2d~2 Reflect2D proc~reflect2d~2->proc~interpolatereflectioncoefficient proc~reflect3d Reflect3D proc~reflect3d->proc~interpolatereflectioncoefficient proc~traceray2d TraceRay2D proc~traceray2d->proc~reflect2d proc~traceray2d~2 TraceRay2D proc~traceray2d~2->proc~reflect2d~2 proc~traceray3d TraceRay3D proc~traceray3d->proc~reflect3d proc~bellhopcore BellhopCore proc~bellhopcore->proc~traceray2d proc~bellhopcore->proc~traceray3d proc~bellhopcore~2 BellhopCore proc~bellhopcore~2->proc~traceray2d~2 program~bellhop BELLHOP program~bellhop->proc~bellhopcore~2 program~bellhop3d BELLHOP3D program~bellhop3d->proc~bellhopcore

Source Code

  SUBROUTINE InterpolateReflectionCoefficient( RInt, R, NPts, PRTFile )

    ! Given an angle RInt%ThetaInt, returns the magnitude and
    ! phase of the reflection coefficient (RInt%R, RInt%phi).

    ! Uses linear interpolation using the two nearest abscissas
    ! Assumes phi has been unwrapped so that it varies smoothly.
    ! I tried modifying it to allow a complex angle of incidence but
    ! stopped when I realized I needed to fuss with a complex atan2 routine

    IMPLICIT NONE
    INTEGER,              INTENT( IN    ) :: PRTFile, NPts   ! unit number of print file, # pts in refl. coef.
    TYPE(ReflectionCoef), INTENT( IN    ) :: R( NPts )       ! Reflection coefficient table
    TYPE(ReflectionCoef), INTENT( INOUT ) :: RInt            ! interpolated value of refl. coef.
    INTEGER         :: iLeft, iRight, iMid
    REAL ( KIND=8 ) :: alpha, Thetaintr

    iLeft  = 1
    iRight = NPts

    thetaIntr = REAL( RInt%Theta )   ! This should be unnecessary? probably used when I was doing complex angles

    ! Three cases: ThetaInt left, in, or right of tabulated interval

    IF     ( thetaIntr < R( iLeft  )%theta ) THEN
       !iRight = 2
       RInt%R   = 0.0     ! R( iLeft  )%R
       RInt%phi = 0.0     ! R( iLeft  )%phi
       WRITE( PRTFile, * ) 'Warning in InterpolateReflectionCoefficient : Refl. Coef. being set to 0 outside tabulated domain'
       WRITE( PRTFile, * ) 'angle = ', thetaintr, 'lower limit = ', R( iLeft)%theta

    ELSE IF( thetaIntr > R( iRight )%theta ) THEN
       !iLeft = NPts - 1
       RInt%R   = 0.0     ! R( iRight )%R
       RInt%phi = 0.0     ! R( iRight )%phi
       ! WRITE( PRTFile, * ) 'Warning in InterpolateReflectionCoefficient : Refl. Coef. being set to 0 outside tabulated domain'
       ! WRITE( PRTFile, * ) 'angle = ', ThetaIntr, 'upper limit = ', R( iRight )%theta

    ELSE
       ! Search for bracketing abscissas: Log2( NPts ) stabs required for a bracket

       DO WHILE ( iLeft /= iRight - 1 )
          iMid = ( iLeft + iRight ) / 2
          IF ( R( iMid )%theta > thetaIntr ) THEN
             iRight = iMid
          ELSE
             iLeft  = iMid
          ENDIF
       END DO

       ! Linear interpolation for reflection coef

       alpha    = ( RInt%theta - R( iLeft )%theta ) / ( R( iRight )%theta - R( iLeft )%theta )
       RInt%R   = ( 1 - alpha ) * R( iLeft )%R   + alpha * R( iRight )%R
       RInt%phi = ( 1 - alpha ) * R( iLeft )%phi + alpha * R( iRight )%phi

    ENDIF

  END SUBROUTINE InterpolateReflectionCoefficient