Removal of labelled statements (NOAA-EMC#1457)

Author: ukmo-juan-castillo

model/src/constants.F90

@@ -34,6 +34,7 @@ MODULE CONSTANTS
   !/    20-Jan-2017 : Add parameters for ESMF             ( version 6.02 )
   !/    01-Mar-2018 : Add UNDEF parameter                 ( version 6.02 )
   !/    05-Jun-2018 : Add PDLIB parameters                ( version 6.04 )
+  !/    04-Jul-2025 : Remove labelled statements          ( version X.XX )
   !/
   !/    Copyright 2009-2012 National Weather Service (NWS),
   !/       National Oceanic and Atmospheric Administration.  All rights
@@ -275,134 +276,131 @@ SUBROUTINE KZEONE(X, Y, RE0, IM0, RE1, IM1)
     ! ABSCISSAS AND THE WEIGHT FACTORS USED IN THE GAUSS-
     ! HERMITE QUADRATURE.
     R2 = X*X + Y*Y
-    IF (R2.GE.1.96D2) GO TO 50
-    IF (R2.GE.1.849D1) GO TO 30
-    ! THIS SECTION CALCULATES THE FUNCTIONS USING THE SERIES
-    ! EXPANSIONS
-    X2 = X/2.0D0
-    Y2 = Y/2.0D0
-    P1 = X2*X2
-    P2 = Y2*Y2
-    T1 = -(DLOG(P1+P2)/2.0D0+0.5772156649015329D0)
-    ! THE CONSTANT IN THE PRECEDING STATEMENT IS EULER*S
-    ! CONSTANT
-    T2 = -DATAN2(Y,X)
-    X2 = P1 - P2
-    Y2 = X*Y2
-    RTERM = 1.0D0
-    ITERM = 0.0D0
-    RE0 = T1
-    IM0 = T2
-    T1 = T1 + 0.5D0
-    RE1 = T1
-    IM1 = T2
-    P2 = DSQRT(R2)
-    L = 2.106D0*P2 + 4.4D0
-    IF (P2.LT.8.0D-1) L = 2.129D0*P2 + 4.0D0
-    DO N=1,INT(L)
-      P1 = N
-      P2 = N*N
-      R1 = RTERM
-      RTERM = (R1*X2-ITERM*Y2)/P2
-      ITERM = (R1*Y2+ITERM*X2)/P2
-      T1 = T1 + 0.5D0/P1
-      RE0 = RE0 + T1*RTERM - T2*ITERM
-      IM0 = IM0 + T1*ITERM + T2*RTERM
-      P1 = P1 + 1.0D0
-      T1 = T1 + 0.5D0/P1
-      RE1 = RE1 + (T1*RTERM-T2*ITERM)/P1
-      IM1 = IM1 + (T1*ITERM+T2*RTERM)/P1
-    END DO
-    R1 = X/R2 - 0.5D0*(X*RE1-Y*IM1)
-    R2 = -Y/R2 - 0.5D0*(X*IM1+Y*RE1)
-    P1 = DEXP(X)
-    RE0 = P1*RE0
-    IM0 = P1*IM0
-    RE1 = P1*R1
-    IM1 = P1*R2
-    RETURN
-    ! THIS SECTION CALCULATES THE FUNCTIONS USING THE INTEGRAL
-    ! REPRESENTATION, EQN 3, EVALUATED WITH 15 POINT GAUSS-
-    ! HERMITE QUADRATURE
-30  X2 = 2.0D0*X
-    Y2 = 2.0D0*Y
-    R1 = Y2*Y2
-    P1 = DSQRT(X2*X2+R1)
-    P2 = DSQRT(P1+X2)
-    T1 = EXSQ(1)/(2.0D0*P1)
-    RE0 = T1*P2
-    IM0 = T1/P2
-    RE1 = 0.0D0
-    IM1 = 0.0D0
-    DO N=2,8
-      T2 = X2 + TSQ(N)
-      P1 = DSQRT(T2*T2+R1)
-      P2 = DSQRT(P1+T2)
-      T1 = EXSQ(N)/P1
-      RE0 = RE0 + T1*P2
-      IM0 = IM0 + T1/P2
-      T1 = EXSQ(N)*TSQ(N)
-      RE1 = RE1 + T1*P2
-      IM1 = IM1 + T1/P2
-    END DO
-    T2 = -Y2*IM0
-    RE1 = RE1/R2
-    R2 = Y2*IM1/R2
-    RTERM = 1.41421356237309D0*DCOS(Y)
-    ITERM = -1.41421356237309D0*DSIN(Y)
-    ! THE CONSTANT IN THE PREVIOUS STATEMENTS IS,OF COURSE,
-    ! SQRT(2.0).
-    IM0 = RE0*ITERM + T2*RTERM
-    RE0 = RE0*RTERM - T2*ITERM
-    T1 = RE1*RTERM - R2*ITERM
-    T2 = RE1*ITERM + R2*RTERM
-    RE1 = T1*X + T2*Y
-    IM1 = -T1*Y + T2*X
-    RETURN
-    ! THIS SECTION CALCULATES THE FUNCTIONS USING THE
-    ! ASYMPTOTIC EXPANSIONS
-50  RTERM = 1.0D0
-    ITERM = 0.0D0
-    RE0 = 1.0D0
-    IM0 = 0.0D0
-    RE1 = 1.0D0
-    IM1 = 0.0D0
-    P1 = 8.0D0*R2
-    P2 = DSQRT(R2)
-    L = 3.91D0+8.12D1/P2
-    R1 = 1.0D0
-    R2 = 1.0D0
-    M = -8
-    K = 3
-    DO N=1,INT(L)
-      M = M + 8
-      K = K - M
-      R1 = FLOAT(K-4)*R1
-      R2 = FLOAT(K)*R2
-      T1 = FLOAT(N)*P1
-      T2 = RTERM
-      RTERM = (T2*X+ITERM*Y)/T1
-      ITERM = (-T2*Y+ITERM*X)/T1
-      RE0 = RE0 + R1*RTERM
-      IM0 = IM0 + R1*ITERM
-      RE1 = RE1 + R2*RTERM
-      IM1 = IM1 + R2*ITERM
-    END DO
-    T1 = DSQRT(P2+X)
-    T2 = -Y/T1
-    P1 = 8.86226925452758D-1/P2
-    ! THIS CONSTANT IS SQRT(PI)/2.0, WITH PI=3.14159...
-    RTERM = P1*DCOS(Y)
-    ITERM = -P1*DSIN(Y)
-    R1 = RE0*RTERM - IM0*ITERM
-    R2 = RE0*ITERM + IM0*RTERM
-    RE0 = T1*R1 - T2*R2
-    IM0 = T1*R2 + T2*R1
-    R1 = RE1*RTERM - IM1*ITERM
-    R2 = RE1*ITERM + IM1*RTERM
-    RE1 = T1*R1 - T2*R2
-    IM1 = T1*R2 + T2*R1
-    RETURN
+    IF (R2.GE.1.96D2) THEN
+      ! CALCULATE THE FUNCTIONS USING THE ASYMPTOTIC EXPANSIONS
+      RTERM = 1.0D0
+      ITERM = 0.0D0
+      RE0 = 1.0D0
+      IM0 = 0.0D0
+      RE1 = 1.0D0
+      IM1 = 0.0D0
+      P1 = 8.0D0*R2
+      P2 = DSQRT(R2)
+      L = 3.91D0+8.12D1/P2
+      R1 = 1.0D0
+      R2 = 1.0D0
+      M = -8
+      K = 3
+      DO N=1,INT(L)
+        M = M + 8
+        K = K - M
+        R1 = FLOAT(K-4)*R1
+        R2 = FLOAT(K)*R2
+        T1 = FLOAT(N)*P1
+        T2 = RTERM
+        RTERM = (T2*X+ITERM*Y)/T1
+        ITERM = (-T2*Y+ITERM*X)/T1
+        RE0 = RE0 + R1*RTERM
+        IM0 = IM0 + R1*ITERM
+        RE1 = RE1 + R2*RTERM
+        IM1 = IM1 + R2*ITERM
+      END DO
+      T1 = DSQRT(P2+X)
+      T2 = -Y/T1
+      P1 = 8.86226925452758D-1/P2
+      ! THIS CONSTANT IS SQRT(PI)/2.0, WITH PI=3.14159...
+      RTERM = P1*DCOS(Y)
+      ITERM = -P1*DSIN(Y)
+      R1 = RE0*RTERM - IM0*ITERM
+      R2 = RE0*ITERM + IM0*RTERM
+      RE0 = T1*R1 - T2*R2
+      IM0 = T1*R2 + T2*R1
+      R1 = RE1*RTERM - IM1*ITERM
+      R2 = RE1*ITERM + IM1*RTERM
+      RE1 = T1*R1 - T2*R2
+      IM1 = T1*R2 + T2*R1
+    ELSE IF (R2.GE.1.849D1) THEN
+      ! CALCULATE THE FUNCTIONS USING THE INTEGRAL
+      ! REPRESENTATION, EQN 3, EVALUATED WITH 15 POINT GAUSS-
+      ! HERMITE QUADRATURE
+      X2 = 2.0D0*X
+      Y2 = 2.0D0*Y
+      R1 = Y2*Y2
+      P1 = DSQRT(X2*X2+R1)
+      P2 = DSQRT(P1+X2)
+      T1 = EXSQ(1)/(2.0D0*P1)
+      RE0 = T1*P2
+      IM0 = T1/P2
+      RE1 = 0.0D0
+      IM1 = 0.0D0
+      DO N=2,8
+        T2 = X2 + TSQ(N)
+        P1 = DSQRT(T2*T2+R1)
+        P2 = DSQRT(P1+T2)
+        T1 = EXSQ(N)/P1
+        RE0 = RE0 + T1*P2
+        IM0 = IM0 + T1/P2
+        T1 = EXSQ(N)*TSQ(N)
+        RE1 = RE1 + T1*P2
+        IM1 = IM1 + T1/P2
+      END DO
+      T2 = -Y2*IM0
+      RE1 = RE1/R2
+      R2 = Y2*IM1/R2
+      RTERM = 1.41421356237309D0*DCOS(Y)
+      ITERM = -1.41421356237309D0*DSIN(Y)
+      ! THE CONSTANT IN THE PREVIOUS STATEMENTS IS,OF COURSE,
+      ! SQRT(2.0).
+      IM0 = RE0*ITERM + T2*RTERM
+      RE0 = RE0*RTERM - T2*ITERM
+      T1 = RE1*RTERM - R2*ITERM
+      T2 = RE1*ITERM + R2*RTERM
+      RE1 = T1*X + T2*Y
+      IM1 = -T1*Y + T2*X
+    ELSE
+      ! CALCULATE THE FUNCTIONS USING THE SERIES EXPANSIONS
+      X2 = X/2.0D0
+      Y2 = Y/2.0D0
+      P1 = X2*X2
+      P2 = Y2*Y2
+      T1 = -(DLOG(P1+P2)/2.0D0+0.5772156649015329D0)
+      ! THE CONSTANT IN THE PRECEDING STATEMENT IS EULER*S
+      ! CONSTANT
+      T2 = -DATAN2(Y,X)
+      X2 = P1 - P2
+      Y2 = X*Y2
+      RTERM = 1.0D0
+      ITERM = 0.0D0
+      RE0 = T1
+      IM0 = T2
+      T1 = T1 + 0.5D0
+      RE1 = T1
+      IM1 = T2
+      P2 = DSQRT(R2)
+      L = 2.106D0*P2 + 4.4D0
+      IF (P2.LT.8.0D-1) L = 2.129D0*P2 + 4.0D0
+      DO N=1,INT(L)
+        P1 = N
+        P2 = N*N
+        R1 = RTERM
+        RTERM = (R1*X2-ITERM*Y2)/P2
+        ITERM = (R1*Y2+ITERM*X2)/P2
+        T1 = T1 + 0.5D0/P1
+        RE0 = RE0 + T1*RTERM - T2*ITERM
+        IM0 = IM0 + T1*ITERM + T2*RTERM
+        P1 = P1 + 1.0D0
+        T1 = T1 + 0.5D0/P1
+        RE1 = RE1 + (T1*RTERM-T2*ITERM)/P1
+        IM1 = IM1 + (T1*ITERM+T2*RTERM)/P1
+      END DO
+      R1 = X/R2 - 0.5D0*(X*RE1-Y*IM1)
+      R2 = -Y/R2 - 0.5D0*(X*IM1+Y*RE1)
+      P1 = DEXP(X)
+      RE0 = P1*RE0
+      IM0 = P1*IM0
+      RE1 = P1*R1
+      IM1 = P1*R2
+    END IF
   END SUBROUTINE KZEONE
 
   ! ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++