NCEPLIBS-w3emc  2.11.0
w3fa12.f
Go to the documentation of this file.
1 C> @file
2 C> @brief Computes legendre polynomials at a given latitude.
3 C> @author Joe Sela @date 1980-10-28
4 C>
5 C> Subroutine computes legendre polynomials at a given latitude.
6 C>
7 C> Program history log:
8 C> - Joe Sela 1980-10-20
9 C> - Ralph Jones 1984-06-01 Change to ibm vs fortran.
10 C> - Ralph Jones 1993-04-12 Changes for cray, double precision to real.
11 C>
12 C> @param[out] PLN Real locations contain legendre
13 C> polynomials, size is (jcap+2)*(jcap+1)
14 C> @param[in] COLRAD Colatitude in radians of desired point.
15 C> @param[in] JCAP For rhomboiadal truncation of zonal wave
16 C> @param[in] EPS Coeff. used in recursion equation.
17 C> Dimension of eps is (jcap+2)*(jcap+1)
18 C>
19 C> @author Joe Sela @date 1980-10-28
20  SUBROUTINE w3fa12(PLN,COLRAD,JCAP,EPS)
21  REAL A
22  REAL B
23  REAL COLRAD
24  REAL COS2
25  REAL EPS(*)
26  REAL FL
27  REAL PROD
28  REAL P1
29  REAL P2
30  REAL P3
31  REAL SINLAT
32  REAL PLN(*)
33 C
34  SAVE
35 C
36  sinlat = cos(colrad)
37  cos2 = 1.0 - sinlat * sinlat
38  prod = 1.0
39  a = 1.0
40  b = 0.0
41  jcap1 = jcap+1
42  jcap2 = jcap+2
43 C
44  DO 300 ll = 1,jcap1
45  l = ll - 1
46  fl = l
47  jle = l * jcap2
48  IF (l.EQ.0) GO TO 100
49  a = a + 2.0
50  b = b + 2.0
51  prod = prod * cos2 * a / b
52  100 CONTINUE
53  p1 = sqrt(0.5 * prod)
54  pln(jle+1) = p1
55  p2 = sqrt(2.0 * fl + 3.0) * sinlat * p1
56  pln(jle+2) = p2
57 C
58  DO 200 n = 3,jcap2
59  lindex = jle + n
60  p3 = (sinlat*p2 - eps(lindex-1)*p1)/eps(lindex)
61  pln(lindex) = p3
62  p1 = p2
63  p2 = p3
64 200 CONTINUE
65 300 CONTINUE
66  RETURN
67 C
68  END