1C Copyright 1981-2016 ECMWF. 2C 3C This software is licensed under the terms of the Apache Licence 4C Version 2.0 which can be obtained at http://www.apache.org/licenses/LICENSE-2.0. 5C 6C In applying this licence, ECMWF does not waive the privileges and immunities 7C granted to it by virtue of its status as an intergovernmental organisation 8C nor does it submit to any jurisdiction. 9C 10 11 SUBROUTINE JSPLEG1( PLEG, PLAT, KTRUNC) 12 IMPLICIT NONE 13C 14C----> 15C**** JSPLEG1 - Routine to calculate legendre functions 16C 17C Purpose 18C -------- 19C 20C This routine calculates the legendre functions for one latitude. 21C (but not their derivatives) 22C 23C 24C Interface 25C ---------- 26C 27C CALL JSPLEG1( PLEG, PLAT, KTRUNC) 28C 29C 30C Input parameters 31C ---------------- 32C 33C PLAT - Latitude in radians 34C KTRUNC - Spectral truncation 35C 36C 37C Output parameters 38C ----------------- 39C 40C PLEG - Array of legendre functions for one latitude. 41C The array must be at least (KTRUNC+1)*(KTRUNC+4)/2 42C words long. 43C 44C 45C Common block usage 46C ------------------ 47C 48C None 49C 50C 51C Method 52C ------ 53C 54C Recurrence relation with explicit relations for P(m,m) and 55C P(m,m+1) 56C 57C 58C Externals 59C --------- 60C 61C None 62C 63C 64C Reference 65C --------- 66C 67C None 68C 69C 70C Comments 71C -------- 72C 73C Rewritten from SPLEG1 to avoid using common blocks. 74C Work arrays ZHLPx currently dimensioned for T213 using 75C parameter JPTRP1. 76C 77C 78C AUTHOR 79C ------ 80C 81C J.D.Chambers ECMWF 9 November 1993 82C 83C 84C Modifications 85C ------------- 86C 87C None 88C 89C----< 90C _______________________________________________________ 91C 92C* Section 0. Definition of variables. 93C _______________________________________________________ 94C 95C* Prefix conventions for variable names 96C 97C Logical L (but not LP), global or common. 98C O, dummy argument 99C G, local variable 100C LP, parameter. 101C Character C, global or common. 102C H, dummy argument 103C Y (but not YP), local variable 104C YP, parameter. 105C Integer M and N, global or common. 106C K, dummy argument 107C I, local variable 108C J (but not JP), loop control 109C JP, parameter. 110C Real A to F and Q to X, global or common. 111C P (but not PP), dummy argument 112C Z, local variable 113C PP, parameter. 114#include "jparams.h" 115C 116C Subroutine arguments 117C 118 REAL PLEG, PLAT 119 DIMENSION PLEG(*) 120 INTEGER KTRUNC 121C 122C Local variables 123C 124 REAL ZHLP1, ZHLP2, ZHLP3 125 DIMENSION ZHLP1(JPTRP1), ZHLP2(JPTRP1), ZHLP3(JPTRP1) 126 INTEGER ITOUT1, I1M, ILM, JM, JCN, IM2 127 REAL ZSIN, ZCOS, ZF1M, ZRE1, ZF2M, ZM, ZN, ZE1, ZE2 128C 129C _______________________________________________________ 130C 131C 132C* Section 1. Initialization. 133C _______________________________________________________ 134C 135 100 CONTINUE 136C 137 ITOUT1 = KTRUNC+1 138 ZSIN = SIN(PLAT) 139 ZCOS = SQRT(1.-ZSIN*ZSIN) 140C 141C Step 1. M = 0, N = 0 and N = 1 142C 143 ILM = 2 144 PLEG(1) = 1.0 145 ZF1M = SQRT(3.0) 146 PLEG(2) = ZF1M*ZSIN 147C _______________________________________________________ 148C 149C 150C Step 2. Sum for M = 0 to T (T = truncation) 151C _______________________________________________________ 152C 153 200 CONTINUE 154C 155 DO 205 JM = 2,ITOUT1 156 ZM = JM-1 157 ZHLP1(JM) = SQRT(2.*ZM+3.) 158 ZHLP2(JM) = 1./SQRT(2.*ZM) 159 205 CONTINUE 160 ZHLP1(1) = SQRT(3.) 161C 162 DO 570 JM = 1,ITOUT1 163C 164 I1M = JM-1 165 ZM = I1M 166 ZRE1 = ZHLP1(JM) 167 ZE1 = 1./ZRE1 168C _______________________________________________________ 169C 170C 171C Step 3. M > 0 only 172C _______________________________________________________ 173C 174 300 CONTINUE 175C 176 IF (I1M.NE.0) THEN 177 ZF2M = ZF1M*ZCOS*ZHLP2(JM) 178 ZF1M = ZF2M*ZRE1 179C _______________________________________________________ 180C 181C 182C Step 4. N = M and N = M+1 183C _______________________________________________________ 184C 185 400 CONTINUE 186C 187 ILM = ILM+1 188 PLEG(ILM) = ZF2M 189 ILM = ILM+1 190 PLEG(ILM) = ZF1M*ZSIN 191C 192C When output truncation is reached, return to calling program 193 IF ( JM .EQ. ITOUT1 ) GOTO 570 194C 195 ENDIF 196C _______________________________________________________ 197C 198C 199C Step 5. Sum for N = M+2 to T+1 200C _______________________________________________________ 201C 202 500 CONTINUE 203C 204 IM2 = I1M+2 205C 206 DO 520 JCN = IM2,ITOUT1 207 ZN = JCN 208 ZHLP3(JCN) = SQRT((4.*ZN*ZN-1.)/(ZN*ZN-ZM*ZM)) 209 520 CONTINUE 210C 211 DO 560 JCN = IM2,ITOUT1 212C 213 ZE2 = ZHLP3(JCN) 214 ILM = ILM+1 215 PLEG(ILM) = ZE2*(ZSIN*PLEG(ILM-1)-ZE1*PLEG(ILM-2)) 216 ZE2 = 1./ZE2 217 ZE1 = ZE2 218C 219 560 CONTINUE 220C 221 570 CONTINUE 222C 223C _______________________________________________________ 224C 225C 226C* Section 9. Return to calling routine. 227C _______________________________________________________ 228C 229 900 CONTINUE 230C 231 RETURN 232 END 233