SUBROUTINE ReadRayElevationAngles( freq, Depth, TopOpt, RunType )
REAL (KIND=8), INTENT( IN ) :: freq, Depth
CHARACTER (LEN= 6), INTENT( IN ) :: TopOpt, RunType
REAL (KIND=8) :: d_theta_recommended
IF ( TopOpt( 6 : 6 ) == 'I' ) THEN
READ( ENVFile, * ) Angles%Nalpha, Angles%iSingle_alpha ! option to trace a single beam
ELSE
READ( ENVFile, * ) Angles%Nalpha
END IF
IF ( Angles%Nalpha == 0 ) THEN ! automatically estimate Nalpha to use
IF ( RunType( 1 : 1 ) == 'R' ) THEN
Angles%Nalpha = 50 ! For a ray trace plot, we don't want too many rays ...
ELSE
! you're letting ME choose? OK: ideas based on an isospeed ocean
! limit based on phase of adjacent beams at maximum range
Angles%Nalpha = MAX( INT( ( ( 0.3 * Pos%Rr( Pos%NRr ) ) * freq ) / c0 ), 300 )
! limit based on having beams that are thin with respect to the water depth
! assumes also a full 360 degree angular spread of rays
! Should check which Depth is used here, in case where there is a variable bathymetry
d_theta_recommended = ATAN( Depth / ( 10.0 * Pos%Rr( Pos%NRr ) ) )
Angles%Nalpha = MAX( INT( pi / d_theta_recommended ), Angles%Nalpha )
END IF
END IF
ALLOCATE( Angles%alpha( MAX( 3, Angles%Nalpha ) ), STAT = AllocateStatus )
IF ( AllocateStatus /= 0 ) CALL ERROUT( 'ReadRayElevationAngles', 'Insufficient memory to store beam angles' )
IF ( Angles%Nalpha > 2 ) Angles%alpha( 3 ) = -999.9
READ( ENVFile, * ) Angles%alpha
CALL SubTab( Angles%alpha, Angles%Nalpha )
CALL Sort( Angles%alpha, Angles%Nalpha )
! full 360-degree sweep? remove duplicate beam
! LP: Changed from TINY( ), see README.md.
IF ( Angles%Nalpha > 1 .AND. ABS( MOD( Angles%alpha( Angles%Nalpha ) - Angles%alpha( 1 ), &
360.0D0 ) ) < 10.0 * SPACING( 360.0D0 ) ) &
Angles%Nalpha = Angles%Nalpha - 1
WRITE( PRTFile, * ) '__________________________________________________________________________'
WRITE( PRTFile, * )
WRITE( PRTFile, * ) ' Number of beams in elevation = ', Angles%Nalpha
IF ( Angles%iSingle_alpha > 0 ) WRITE( PRTFile, * ) 'Trace only beam number ', Angles%iSingle_alpha
WRITE( PRTFile, * ) ' Beam take-off angles (degrees)'
IF ( Angles%Nalpha >= 1 ) WRITE( PRTFile, "( 5G14.6 )" ) Angles%alpha( 1 : MIN( Angles%Nalpha, Number_to_Echo ) )
IF ( Angles%Nalpha > Number_to_Echo ) WRITE( PRTFile, "( G14.6 )" ) ' ... ', Angles%alpha( Angles%Nalpha )
IF ( Angles%Nalpha > 1 .AND. Angles%alpha( Angles%Nalpha ) == Angles%alpha( 1 ) ) &
CALL ERROUT( 'ReadRayElevationAngles', 'First and last beam take-off angle are identical' )
IF ( TopOpt( 6 : 6 ) == 'I' ) THEN
IF ( Angles%iSingle_alpha < 1 .OR. Angles%iSingle_alpha > Angles%Nalpha ) &
CALL ERROUT( 'ReadRayElevationAngles', 'Selected beam, iSingle_alpha not in [ 1, Angles%Nalpha ]' )
END IF
END SUBROUTINE ReadRayElevationAngles