module fit_functions use constants, only: dp, pi implicit none private real(dp), parameter :: eps = 1e-15_dp public :: S_A, S_B contains elemental function S_A(angle) real(dp), intent(in) :: angle real(dp) :: S_A real(dp) :: angle_in_period real(dp), parameter :: fit_fac = 0.3_dp angle_in_period = modulo(angle + pi, 2.0_dp*pi) - pi S_A = fit_fac*(angle_in_period - sign(0.5_dp*pi, angle_in_period)) if (abs(mod(angle, pi)) < eps) S_A = 0.0_dp end function S_A elemental function S_B(angle) real(dp), intent(in) :: angle real(dp) :: S_B real(dp) :: angle_in_period real(dp), parameter :: fit_fac = 2.0_dp angle_in_period = modulo(angle + pi, 2.0_dp*pi) - pi S_B = sign(fit_fac, angle_in_period) if (abs(mod(angle, pi)) < eps) S_B = 0.0_dp end function S_B end module fit_functions