| 1 |
C $Header: /u/gcmpack/MITgcm/pkg/fizhi/fizhi_utils.F,v 1.14 2005/05/20 23:30:11 molod Exp $ |
| 2 |
C $Name: $ |
| 3 |
|
| 4 |
#include "FIZHI_OPTIONS.h" |
| 5 |
function minval (q,im) |
| 6 |
implicit none |
| 7 |
integer im, i |
| 8 |
_RL q(im), minval |
| 9 |
minval = 1.e15 |
| 10 |
do i=1,im |
| 11 |
if( q(i).lt.minval ) minval = q(i) |
| 12 |
enddo |
| 13 |
return |
| 14 |
end |
| 15 |
FUNCTION ERRF (ARG) |
| 16 |
C*********************************************************************** |
| 17 |
C FUNCTION ERRF |
| 18 |
C PURPOSE |
| 19 |
C COMPUTES ERROR FUNCTION OF ARGUMENT |
| 20 |
C USAGE |
| 21 |
C CALLED BY TRBFLX |
| 22 |
C DESCRIPTION OF PARAMETERS |
| 23 |
C ARG - INPUTED ARGUMENT |
| 24 |
C REMARKS: |
| 25 |
C USED TO COMPUTE FRACTIONAL CLOUD COVER AND LIQUID WATER CONTENT |
| 26 |
C FROM TURBULENCE STATISTICS |
| 27 |
C ********************************************************************** |
| 28 |
implicit none |
| 29 |
_RL arg,errf |
| 30 |
|
| 31 |
_RL aa1,aa2,aa3,aa4,aa5,pp,x2,x3,x4,x5 |
| 32 |
PARAMETER ( AA1 = 0.254829592 ) |
| 33 |
PARAMETER ( AA2 = -0.284496736 ) |
| 34 |
PARAMETER ( AA3 = 1.421413741 ) |
| 35 |
PARAMETER ( AA4 = -1.453152027 ) |
| 36 |
PARAMETER ( AA5 = 1.061405429 ) |
| 37 |
PARAMETER ( PP = 0.3275911 ) |
| 38 |
PARAMETER ( X2 = AA2 / AA1 ) |
| 39 |
PARAMETER ( X3 = AA3 / AA2 ) |
| 40 |
PARAMETER ( X4 = AA5 / AA3 ) |
| 41 |
PARAMETER ( X5 = AA5 / AA4 ) |
| 42 |
|
| 43 |
_RL aarg,tt |
| 44 |
|
| 45 |
ERRF = 1. |
| 46 |
AARG=ABS(ARG) |
| 47 |
|
| 48 |
IF ( AARG .LT. 4.0 ) THEN |
| 49 |
TT = 1./(1.+PP*AARG) |
| 50 |
ERRF = 1. - |
| 51 |
1 (AA1*TT*(1.+X2*TT*(1.+X3*TT*(1.+X4*TT*(1.+X5*TT))))) |
| 52 |
2 * EXP(-AARG*AARG) |
| 53 |
ENDIF |
| 54 |
|
| 55 |
IF ( ARG .LT. 0.0 ) ERRF = -ERRF |
| 56 |
|
| 57 |
RETURN |
| 58 |
END |
| 59 |
|
| 60 |
SUBROUTINE STRIP(A,B,IA,IB,L,K) |
| 61 |
implicit none |
| 62 |
integer ia,ib,L,K |
| 63 |
_RL A(IA,L), B(IB,L) |
| 64 |
|
| 65 |
INTEGER OFFSET,Lena,i,j |
| 66 |
|
| 67 |
OFFSET = IB*(K-1) |
| 68 |
Lena = MIN(IB,IA-OFFSET) |
| 69 |
OFFSET = OFFSET+1 |
| 70 |
|
| 71 |
IF(Lena.EQ.IB) THEN |
| 72 |
DO 100 J=1,L |
| 73 |
DO 100 I=1,Lena |
| 74 |
B(I,J) = A(I+OFFSET-1,J) |
| 75 |
100 CONTINUE |
| 76 |
ELSE |
| 77 |
DO 200 J=1,L |
| 78 |
DO 300 I=1,Lena |
| 79 |
B(I,J) = A(I+OFFSET-1,J) |
| 80 |
300 CONTINUE |
| 81 |
DO 400 I=1,IB-Lena |
| 82 |
B(Lena+I,J) = A(Lena+OFFSET-1,J) |
| 83 |
400 CONTINUE |
| 84 |
200 CONTINUE |
| 85 |
ENDIF |
| 86 |
|
| 87 |
RETURN |
| 88 |
END |
| 89 |
SUBROUTINE STRIPINT(A,B,IA,IB,L,K) |
| 90 |
implicit none |
| 91 |
integer ia,ib,L,K |
| 92 |
INTEGER A(IA,L), B(IB,L) |
| 93 |
|
| 94 |
INTEGER OFFSET,Lena,i,j |
| 95 |
|
| 96 |
OFFSET = IB*(K-1) |
| 97 |
Lena = MIN(IB,IA-OFFSET) |
| 98 |
OFFSET = OFFSET+1 |
| 99 |
|
| 100 |
IF(Lena.EQ.IB) THEN |
| 101 |
DO 100 J=1,L |
| 102 |
DO 100 I=1,Lena |
| 103 |
B(I,J) = A(I+OFFSET-1,J) |
| 104 |
100 CONTINUE |
| 105 |
ELSE |
| 106 |
DO 200 J=1,L |
| 107 |
DO 300 I=1,Lena |
| 108 |
B(I,J) = A(I+OFFSET-1,J) |
| 109 |
300 CONTINUE |
| 110 |
DO 400 I=1,IB-Lena |
| 111 |
B(Lena+I,J) = A(Lena+OFFSET-1,J) |
| 112 |
400 CONTINUE |
| 113 |
200 CONTINUE |
| 114 |
ENDIF |
| 115 |
|
| 116 |
RETURN |
| 117 |
END |
| 118 |
SUBROUTINE PASTE(B,A,IB,IA,L,K) |
| 119 |
implicit none |
| 120 |
integer ia,ib,L,K |
| 121 |
_RL A(IA,L), B(IB,L) |
| 122 |
|
| 123 |
INTEGER OFFSET,Lena,i,j |
| 124 |
|
| 125 |
OFFSET = IB*(K-1) |
| 126 |
Lena = MIN(IB,IA-OFFSET) |
| 127 |
OFFSET = OFFSET+1 |
| 128 |
|
| 129 |
DO 100 J=1,L |
| 130 |
DO 100 I=1,Lena |
| 131 |
A(I+OFFSET-1,J) = B(I,J) |
| 132 |
100 CONTINUE |
| 133 |
|
| 134 |
RETURN |
| 135 |
END |
| 136 |
SUBROUTINE PSTBMP(B,A,IB,IA,L,K) |
| 137 |
implicit none |
| 138 |
integer ia,ib,L,K |
| 139 |
_RL A(IA,L), B(IB,L) |
| 140 |
|
| 141 |
INTEGER OFFSET,Lena,i,j |
| 142 |
|
| 143 |
OFFSET = IB*(K-1) |
| 144 |
Lena = MIN(IB,IA-OFFSET) |
| 145 |
OFFSET = OFFSET+1 |
| 146 |
|
| 147 |
DO 100 J=1,L |
| 148 |
DO 100 I=1,Lena |
| 149 |
A(I+OFFSET-1,J) = A(I+OFFSET-1,J) + B(I,J) |
| 150 |
100 CONTINUE |
| 151 |
C |
| 152 |
RETURN |
| 153 |
END |
| 154 |
SUBROUTINE STRINT(A,B,IA,IB,L,K) |
| 155 |
implicit none |
| 156 |
integer ia,ib,L,K |
| 157 |
INTEGER A(IA,L), B(IB,L) |
| 158 |
|
| 159 |
INTEGER OFFSET,Lena,i,j |
| 160 |
|
| 161 |
OFFSET = IB*(K-1) |
| 162 |
Lena = MIN(IB,IA-OFFSET) |
| 163 |
OFFSET = OFFSET+1 |
| 164 |
|
| 165 |
IF(Lena.EQ.IB) THEN |
| 166 |
DO 100 J=1,L |
| 167 |
DO 100 I=1,Lena |
| 168 |
B(I,J) = A(I+OFFSET-1,J) |
| 169 |
100 CONTINUE |
| 170 |
ELSE |
| 171 |
DO 200 J=1,L |
| 172 |
DO 300 I=1,Lena |
| 173 |
B(I,J) = A(I+OFFSET-1,J) |
| 174 |
300 CONTINUE |
| 175 |
DO 400 I=1,IB-Lena |
| 176 |
B(Lena+I,J) = A(Lena+OFFSET-1,J) |
| 177 |
400 CONTINUE |
| 178 |
200 CONTINUE |
| 179 |
ENDIF |
| 180 |
|
| 181 |
RETURN |
| 182 |
END |
| 183 |
SUBROUTINE QSAT (TT,P,Q,DQDT,LDQDT) |
| 184 |
C*********************************************************************** |
| 185 |
C |
| 186 |
C PURPOSE: |
| 187 |
C ======== |
| 188 |
C Compute Saturation Specific Humidity |
| 189 |
C |
| 190 |
C INPUT: |
| 191 |
C ====== |
| 192 |
C TT ......... Temperature (Kelvin) |
| 193 |
C P .......... Pressure (mb) |
| 194 |
C LDQDT ...... Logical Flag to compute QSAT Derivative |
| 195 |
C |
| 196 |
C OUTPUT: |
| 197 |
C ======= |
| 198 |
C Q .......... Saturation Specific Humidity |
| 199 |
C DQDT ....... Saturation Specific Humidity Derivative wrt Temperature |
| 200 |
C |
| 201 |
C |
| 202 |
C*********************************************************************** |
| 203 |
|
| 204 |
IMPLICIT NONE |
| 205 |
_RL TT, P, Q, DQDT |
| 206 |
LOGICAL LDQDT |
| 207 |
|
| 208 |
_RL AIRMW, H2OMW |
| 209 |
|
| 210 |
PARAMETER ( AIRMW = 28.97 ) |
| 211 |
PARAMETER ( H2OMW = 18.01 ) |
| 212 |
|
| 213 |
_RL ESFAC, ERFAC |
| 214 |
PARAMETER ( ESFAC = H2OMW/AIRMW ) |
| 215 |
PARAMETER ( ERFAC = (1.0-ESFAC)/ESFAC ) |
| 216 |
|
| 217 |
_RL aw0, aw1, aw2, aw3, aw4, aw5, aw6 |
| 218 |
_RL bw0, bw1, bw2, bw3, bw4, bw5, bw6 |
| 219 |
_RL ai0, ai1, ai2, ai3, ai4, ai5, ai6 |
| 220 |
_RL bi0, bi1, bi2, bi3, bi4, bi5, bi6 |
| 221 |
|
| 222 |
_RL d0, d1, d2, d3, d4, d5, d6 |
| 223 |
_RL e0, e1, e2, e3, e4, e5, e6 |
| 224 |
_RL f0, f1, f2, f3, f4, f5, f6 |
| 225 |
_RL g0, g1, g2, g3, g4, g5, g6 |
| 226 |
|
| 227 |
c ******************************************************** |
| 228 |
c *** Polynomial Coefficients WRT Water (Lowe, 1977) **** |
| 229 |
c *** (Valid +50 C to -50 C) **** |
| 230 |
c ******************************************************** |
| 231 |
|
| 232 |
parameter ( aw0 = 6.107799961e+00 * esfac ) |
| 233 |
parameter ( aw1 = 4.436518521e-01 * esfac ) |
| 234 |
parameter ( aw2 = 1.428945805e-02 * esfac ) |
| 235 |
parameter ( aw3 = 2.650648471e-04 * esfac ) |
| 236 |
parameter ( aw4 = 3.031240396e-06 * esfac ) |
| 237 |
parameter ( aw5 = 2.034080948e-08 * esfac ) |
| 238 |
parameter ( aw6 = 6.136820929e-11 * esfac ) |
| 239 |
|
| 240 |
parameter ( bw0 = +4.438099984e-01 * esfac ) |
| 241 |
parameter ( bw1 = +2.857002636e-02 * esfac ) |
| 242 |
parameter ( bw2 = +7.938054040e-04 * esfac ) |
| 243 |
parameter ( bw3 = +1.215215065e-05 * esfac ) |
| 244 |
parameter ( bw4 = +1.036561403e-07 * esfac ) |
| 245 |
parameter ( bw5 = +3.532421810e-10 * esfac ) |
| 246 |
parameter ( bw6 = -7.090244804e-13 * esfac ) |
| 247 |
|
| 248 |
|
| 249 |
c ******************************************************** |
| 250 |
c *** Polynomial Coefficients WRT Ice (Lowe, 1977) **** |
| 251 |
c *** (Valid +0 C to -50 C) **** |
| 252 |
c ******************************************************** |
| 253 |
|
| 254 |
parameter ( ai0 = +6.109177956e+00 * esfac ) |
| 255 |
parameter ( ai1 = +5.034698970e-01 * esfac ) |
| 256 |
parameter ( ai2 = +1.886013408e-02 * esfac ) |
| 257 |
parameter ( ai3 = +4.176223716e-04 * esfac ) |
| 258 |
parameter ( ai4 = +5.824720280e-06 * esfac ) |
| 259 |
parameter ( ai5 = +4.838803174e-08 * esfac ) |
| 260 |
parameter ( ai6 = +1.838826904e-10 * esfac ) |
| 261 |
|
| 262 |
parameter ( bi0 = +5.030305237e-01 * esfac ) |
| 263 |
parameter ( bi1 = +3.773255020e-02 * esfac ) |
| 264 |
parameter ( bi2 = +1.267995369e-03 * esfac ) |
| 265 |
parameter ( bi3 = +2.477563108e-05 * esfac ) |
| 266 |
parameter ( bi4 = +3.005693132e-07 * esfac ) |
| 267 |
parameter ( bi5 = +2.158542548e-09 * esfac ) |
| 268 |
parameter ( bi6 = +7.131097725e-12 * esfac ) |
| 269 |
|
| 270 |
|
| 271 |
c ******************************************************** |
| 272 |
c *** Polynomial Coefficients WRT Ice **** |
| 273 |
c *** Starr and Cox (1985) (Valid -40 C to -70 C) **** |
| 274 |
c ******************************************************** |
| 275 |
|
| 276 |
|
| 277 |
parameter ( d0 = 0.535098336e+01 * esfac ) |
| 278 |
parameter ( d1 = 0.401390832e+00 * esfac ) |
| 279 |
parameter ( d2 = 0.129690326e-01 * esfac ) |
| 280 |
parameter ( d3 = 0.230325039e-03 * esfac ) |
| 281 |
parameter ( d4 = 0.236279781e-05 * esfac ) |
| 282 |
parameter ( d5 = 0.132243858e-07 * esfac ) |
| 283 |
parameter ( d6 = 0.314296723e-10 * esfac ) |
| 284 |
|
| 285 |
parameter ( e0 = 0.469290530e+00 * esfac ) |
| 286 |
parameter ( e1 = 0.333092511e-01 * esfac ) |
| 287 |
parameter ( e2 = 0.102164528e-02 * esfac ) |
| 288 |
parameter ( e3 = 0.172979242e-04 * esfac ) |
| 289 |
parameter ( e4 = 0.170017544e-06 * esfac ) |
| 290 |
parameter ( e5 = 0.916466531e-09 * esfac ) |
| 291 |
parameter ( e6 = 0.210844486e-11 * esfac ) |
| 292 |
|
| 293 |
|
| 294 |
c ******************************************************** |
| 295 |
c *** Polynomial Coefficients WRT Ice **** |
| 296 |
c *** Starr and Cox (1985) (Valid -65 C to -95 C) **** |
| 297 |
c ******************************************************** |
| 298 |
|
| 299 |
parameter ( f0 = 0.298152339e+01 * esfac ) |
| 300 |
parameter ( f1 = 0.191372282e+00 * esfac ) |
| 301 |
parameter ( f2 = 0.517609116e-02 * esfac ) |
| 302 |
parameter ( f3 = 0.754129933e-04 * esfac ) |
| 303 |
parameter ( f4 = 0.623439266e-06 * esfac ) |
| 304 |
parameter ( f5 = 0.276961083e-08 * esfac ) |
| 305 |
parameter ( f6 = 0.516000335e-11 * esfac ) |
| 306 |
|
| 307 |
parameter ( g0 = 0.312654072e+00 * esfac ) |
| 308 |
parameter ( g1 = 0.195789002e-01 * esfac ) |
| 309 |
parameter ( g2 = 0.517837908e-03 * esfac ) |
| 310 |
parameter ( g3 = 0.739410547e-05 * esfac ) |
| 311 |
parameter ( g4 = 0.600331350e-07 * esfac ) |
| 312 |
parameter ( g5 = 0.262430726e-09 * esfac ) |
| 313 |
parameter ( g6 = 0.481960676e-12 * esfac ) |
| 314 |
|
| 315 |
_RL TMAX, TICE |
| 316 |
PARAMETER ( TMAX=323.15, TICE=273.16) |
| 317 |
|
| 318 |
_RL T, D, W, QX, DQX |
| 319 |
T = MIN(TT,TMAX) - TICE |
| 320 |
DQX = 0. |
| 321 |
QX = 0. |
| 322 |
|
| 323 |
c Fitting for temperatures above 0 degrees centigrade |
| 324 |
c --------------------------------------------------- |
| 325 |
if(t.gt.0.) then |
| 326 |
qx = aw0+T*(aw1+T*(aw2+T*(aw3+T*(aw4+T*(aw5+T*aw6))))) |
| 327 |
if (ldqdt) then |
| 328 |
dqx = bw0+T*(bw1+T*(bw2+T*(bw3+T*(bw4+T*(bw5+T*bw6))))) |
| 329 |
endif |
| 330 |
endif |
| 331 |
|
| 332 |
c Fitting for temperatures between 0 and -40 |
| 333 |
c ------------------------------------------ |
| 334 |
if( t.le.0. .and. t.gt.-40.0 ) then |
| 335 |
w = (40.0 + t)/40.0 |
| 336 |
qx = w *(aw0+T*(aw1+T*(aw2+T*(aw3+T*(aw4+T*(aw5+T*aw6)))))) |
| 337 |
. + (1.-w)*(ai0+T*(ai1+T*(ai2+T*(ai3+T*(ai4+T*(ai5+T*ai6)))))) |
| 338 |
if (ldqdt) then |
| 339 |
dqx = w *(bw0+T*(bw1+T*(bw2+T*(bw3+T*(bw4+T*(bw5+T*bw6)))))) |
| 340 |
. + (1.-w)*(bi0+T*(bi1+T*(bi2+T*(bi3+T*(bi4+T*(bi5+T*bi6)))))) |
| 341 |
endif |
| 342 |
endif |
| 343 |
|
| 344 |
c Fitting for temperatures between -40 and -70 |
| 345 |
c -------------------------------------------- |
| 346 |
if( t.le.-40.0 .and. t.ge.-70.0 ) then |
| 347 |
qx = d0+T*(d1+T*(d2+T*(d3+T*(d4+T*(d5+T*d6))))) |
| 348 |
if (ldqdt) then |
| 349 |
dqx = e0+T*(e1+T*(e2+T*(e3+T*(e4+T*(e5+T*e6))))) |
| 350 |
endif |
| 351 |
endif |
| 352 |
|
| 353 |
c Fitting for temperatures less than -70 |
| 354 |
c -------------------------------------- |
| 355 |
if(t.lt.-70.0) then |
| 356 |
qx = f0+t*(f1+t*(f2+t*(f3+t*(f4+t*(f5+t*f6))))) |
| 357 |
if (ldqdt) then |
| 358 |
dqx = g0+t*(g1+t*(g2+t*(g3+t*(g4+t*(g5+t*g6))))) |
| 359 |
endif |
| 360 |
endif |
| 361 |
|
| 362 |
c Compute Saturation Specific Humidity |
| 363 |
c ------------------------------------ |
| 364 |
D = (P-ERFAC*QX) |
| 365 |
IF(D.LT.0.) THEN |
| 366 |
Q = 1.0 |
| 367 |
IF (LDQDT) DQDT = 0. |
| 368 |
ELSE |
| 369 |
D = 1.0 / D |
| 370 |
Q = MIN(QX * D,1.0 _d 0) |
| 371 |
IF (LDQDT) DQDT = (1.0 + ERFAC*Q) * D * DQX |
| 372 |
ENDIF |
| 373 |
RETURN |
| 374 |
END |
| 375 |
subroutine vqsat (tt,p,q,dqdt,ldqdt,n) |
| 376 |
implicit none |
| 377 |
integer i,n |
| 378 |
logical ldqdt |
| 379 |
_RL tt(n), p(n), q(n), dqdt(n) |
| 380 |
#ifdef CRAY |
| 381 |
#ifdef f77 |
| 382 |
cfpp$ expand (qsat) |
| 383 |
#endif |
| 384 |
#endif |
| 385 |
do i=1,n |
| 386 |
call qsat ( tt(i),p(i),q(i),dqdt(i),ldqdt ) |
| 387 |
enddo |
| 388 |
return |
| 389 |
end |
| 390 |
|
| 391 |
subroutine stripit(a,b,irun,ia,ib,l,k) |
| 392 |
implicit none |
| 393 |
integer ia,ib,irun,l,k |
| 394 |
_RL a(ia,l), b(ib,l) |
| 395 |
integer i,j,Lena,offset |
| 396 |
|
| 397 |
offset = ib*(k-1) |
| 398 |
Lena = min(ib,irun-offset) |
| 399 |
offset = offset+1 |
| 400 |
|
| 401 |
if(Lena.eq.ib) then |
| 402 |
do 100 j=1,l |
| 403 |
do 100 i=1,Lena |
| 404 |
b(i,j) = a(i+offset-1,j) |
| 405 |
100 continue |
| 406 |
else |
| 407 |
do 200 j=1,l |
| 408 |
do 300 i=1,Lena |
| 409 |
b(i,j) = a(i+offset-1,j) |
| 410 |
300 continue |
| 411 |
do 400 i=1,ib-Lena |
| 412 |
b(Lena+i,j) = a(Lena+offset-1,j) |
| 413 |
400 continue |
| 414 |
200 continue |
| 415 |
endif |
| 416 |
return |
| 417 |
end |
| 418 |
|
| 419 |
subroutine stripitint(a,b,irun,ia,ib,l,k) |
| 420 |
implicit none |
| 421 |
integer ia,ib,irun,l,k,a(ia,l),b(ib,l) |
| 422 |
integer i,j,Lena,offset |
| 423 |
|
| 424 |
offset = ib*(k-1) |
| 425 |
Lena = min(ib,irun-offset) |
| 426 |
offset = offset+1 |
| 427 |
|
| 428 |
if(Lena.eq.ib) then |
| 429 |
do 100 j=1,l |
| 430 |
do 100 i=1,Lena |
| 431 |
b(i,j) = a(i+offset-1,j) |
| 432 |
100 continue |
| 433 |
else |
| 434 |
do 200 j=1,l |
| 435 |
do 300 i=1,Lena |
| 436 |
b(i,j) = a(i+offset-1,j) |
| 437 |
300 continue |
| 438 |
do 400 i=1,ib-Lena |
| 439 |
b(Lena+i,j) = a(Lena+offset-1,j) |
| 440 |
400 continue |
| 441 |
200 continue |
| 442 |
endif |
| 443 |
return |
| 444 |
end |
| 445 |
|
| 446 |
subroutine pastit(b,a,ib,ia,irun,L,k) |
| 447 |
implicit none |
| 448 |
integer ib,ia,L,k,irun,Lena,offset |
| 449 |
integer i,j |
| 450 |
_RL a(ia,l), b(ib,l) |
| 451 |
|
| 452 |
offset = ib*(k-1) |
| 453 |
Lena = min(ib,irun-offset) |
| 454 |
offset = offset+1 |
| 455 |
|
| 456 |
do 100 j=1,L |
| 457 |
do 100 i=1,Lena |
| 458 |
a(i+offset-1,j) = b(i,j) |
| 459 |
100 continue |
| 460 |
return |
| 461 |
end |
| 462 |
|
| 463 |
subroutine pstbitint(b,a,ib,ia,irun,l,k) |
| 464 |
implicit none |
| 465 |
integer ib,ia,L,k,irun,Lena,offset |
| 466 |
_RL a(ia,l) |
| 467 |
integer b(ib,l) |
| 468 |
integer i,j |
| 469 |
|
| 470 |
offset = ib*(k-1) |
| 471 |
Lena = min(ib,irun-offset) |
| 472 |
offset = offset+1 |
| 473 |
|
| 474 |
do 100 j=1,L |
| 475 |
do 100 i=1,Lena |
| 476 |
a(i+offset-1,j) = a(i+offset-1,j) + float(b(i,j)) |
| 477 |
100 continue |
| 478 |
return |
| 479 |
end |
| 480 |
|
| 481 |
|
| 482 |
subroutine pstbmpit(b,a,ib,ia,irun,l,k) |
| 483 |
implicit none |
| 484 |
integer ib,ia,L,k,irun,Lena,offset |
| 485 |
_RL a(ia,l), b(ib,l) |
| 486 |
integer i,j |
| 487 |
|
| 488 |
offset = ib*(k-1) |
| 489 |
Lena = min(ib,irun-offset) |
| 490 |
offset = offset+1 |
| 491 |
|
| 492 |
do 100 j=1,L |
| 493 |
do 100 i=1,Lena |
| 494 |
a(i+offset-1,j) = a(i+offset-1,j) + b(i,j) |
| 495 |
100 continue |
| 496 |
return |
| 497 |
end |
| 498 |
|
| 499 |
subroutine strip2tile(a,indx,b,irun,ia,ib,levs,npeice) |
| 500 |
c----------------------------------------------------------------------- |
| 501 |
c subroutine strip2tile - extract one processors worth of grid points |
| 502 |
c from a grid space array to a stripped tile |
| 503 |
c space array |
| 504 |
c |
| 505 |
c input: |
| 506 |
c a - array to be stripped FROM [ia,levs] |
| 507 |
c indx - array of horizontal indeces of grid points to convert to |
| 508 |
c tile space |
| 509 |
c irun - number of points in array a that need to be stripped |
| 510 |
c ia - inner of dimension of source array |
| 511 |
c ib - inner dimension of target array AND the number of points |
| 512 |
c in the target array to be filled |
| 513 |
c levs - number of vertical levels AND outer dimension of arrays |
| 514 |
c npeice - the current strip number to be filled |
| 515 |
c output: |
| 516 |
c b - array to be filled, ie, one processors field [ib,levs] |
| 517 |
c----------------------------------------------------------------------- |
| 518 |
implicit none |
| 519 |
integer ia,ib,irun,levs,npeice |
| 520 |
_RL a(ia,levs), b(ib,levs) |
| 521 |
integer indx(irun) |
| 522 |
integer i,k,Lena,offset |
| 523 |
|
| 524 |
offset = ib*(npeice-1) |
| 525 |
Lena = min(ib,irun-offset) |
| 526 |
offset = offset+1 |
| 527 |
|
| 528 |
if(Lena.eq.ib) then |
| 529 |
do 100 k=1,levs |
| 530 |
do 100 i=1,Lena |
| 531 |
b(i,k) = a(indx(i+offset-1),k) |
| 532 |
100 continue |
| 533 |
else |
| 534 |
do 200 k=1,levs |
| 535 |
do 300 i=1,Lena |
| 536 |
b(i,k) = a(indx(i+offset-1),k) |
| 537 |
300 continue |
| 538 |
do 400 i=1,ib-Lena |
| 539 |
b(Lena+i,k) = a(indx(Lena+offset-1),k) |
| 540 |
400 continue |
| 541 |
200 continue |
| 542 |
endif |
| 543 |
return |
| 544 |
end |
| 545 |
|
| 546 |
subroutine paste2grd_old(b,indx,chfr,ib,numpts,a,ia,levs,npeice) |
| 547 |
c----------------------------------------------------------------------- |
| 548 |
c subroutine paste2grd - paste one processors worth of grid points |
| 549 |
c from a stripped tile array to a grid |
| 550 |
c space array |
| 551 |
c |
| 552 |
c input: |
| 553 |
c b - array to be pasted back into grid space [ib,levs] |
| 554 |
c indx - array of horizontal indeces of grid points to convert to |
| 555 |
c tile space[numpts] |
| 556 |
c chfr - fractional area covered by the tile [ib] |
| 557 |
c ib - inner dimension of source array AND number of points in |
| 558 |
c array a that need to be pasted |
| 559 |
c numpts - total number of points which were stripped |
| 560 |
c ia - inner of dimension of target array |
| 561 |
c levs - number of vertical levels AND outer dimension of arrays |
| 562 |
c npeice - the current strip number to be filled |
| 563 |
c output: |
| 564 |
c a - grid space array to be filled [ia,levs] |
| 565 |
c |
| 566 |
c IMPORTANT NOTE: |
| 567 |
c |
| 568 |
c This routine will result in roundoff differences if called from |
| 569 |
c within a parallel region. |
| 570 |
c----------------------------------------------------------------------- |
| 571 |
|
| 572 |
implicit none |
| 573 |
integer ia,ib,levs,numpts,npeice |
| 574 |
integer indx(numpts) |
| 575 |
_RL a(ia,levs), b(ib,levs), chfr(ib) |
| 576 |
|
| 577 |
integer i,L,offset,Lena |
| 578 |
|
| 579 |
offset = ib*(npeice-1) |
| 580 |
Lena = min(ib,numpts-offset) |
| 581 |
offset = offset+1 |
| 582 |
|
| 583 |
do L = 1,levs |
| 584 |
do i=1,Lena |
| 585 |
a(indx(i+offset-1),L) = a(indx(i+offset-1),L) + b(i,L)*chfr(i) |
| 586 |
enddo |
| 587 |
enddo |
| 588 |
return |
| 589 |
end |
| 590 |
subroutine paste2grd (b,indx,chfr,ib,numpts,a,ia,levs,npeice, |
| 591 |
. check) |
| 592 |
c----------------------------------------------------------------------- |
| 593 |
c subroutine paste2grd - paste one processors worth of grid points |
| 594 |
c from a stripped tile array to a grid |
| 595 |
c space array |
| 596 |
c |
| 597 |
c input: |
| 598 |
c b - array to be pasted back into grid space [ib,levs] |
| 599 |
c indx - array of horizontal indeces of grid points to convert to |
| 600 |
c tile space[numpts] |
| 601 |
c chfr - fractional area covered by the tile [ib] |
| 602 |
c ib - inner dimension of source array AND number of points in |
| 603 |
c array a that need to be pasted |
| 604 |
c numpts - total number of points which were stripped |
| 605 |
c ia - inner of dimension of target array |
| 606 |
c levs - number of vertical levels AND outer dimension of arrays |
| 607 |
c npeice - the current strip number to be filled |
| 608 |
c check - logical to check for undefined values |
| 609 |
c output: |
| 610 |
c a - grid space array to be filled [ia,levs] |
| 611 |
c |
| 612 |
c IMPORTANT NOTE: |
| 613 |
c |
| 614 |
c This routine will result in roundoff differences if called from |
| 615 |
c within a parallel region. |
| 616 |
c----------------------------------------------------------------------- |
| 617 |
|
| 618 |
implicit none |
| 619 |
integer ia,ib,levs,numpts,npeice |
| 620 |
integer indx(numpts) |
| 621 |
_RL a(ia,levs), b(ib,levs), chfr(ib) |
| 622 |
logical check |
| 623 |
|
| 624 |
integer i,L,offset,Lena |
| 625 |
_RL undef,getcon |
| 626 |
|
| 627 |
offset = ib*(npeice-1) |
| 628 |
Lena = min(ib,numpts-offset) |
| 629 |
offset = offset+1 |
| 630 |
|
| 631 |
if( check ) then |
| 632 |
undef = getcon('UNDEF') |
| 633 |
do L= 1,levs |
| 634 |
do i= 1,Lena |
| 635 |
if( a(indx(i+offset-1),L).eq.undef .or. b(i,L).eq.undef ) then |
| 636 |
a(indx(i+offset-1),L) = undef |
| 637 |
else |
| 638 |
a(indx(i+offset-1),L)=a(indx(i+offset-1),L) + b(i,L)*chfr(i) |
| 639 |
endif |
| 640 |
enddo |
| 641 |
enddo |
| 642 |
else |
| 643 |
do L= 1,levs |
| 644 |
do i= 1,Lena |
| 645 |
a(indx(i+offset-1),L)=a(indx(i+offset-1),L) + b(i,L)*chfr(i) |
| 646 |
enddo |
| 647 |
enddo |
| 648 |
endif |
| 649 |
|
| 650 |
return |
| 651 |
end |
| 652 |
SUBROUTINE GRD2MSC(A,IM,JM,IGRD,B,MXCHPS,NCHP) |
| 653 |
|
| 654 |
implicit none |
| 655 |
integer im,jm,mxchps,nchp |
| 656 |
integer igrd(mxchps) |
| 657 |
_RL A(IM,JM), B(MXCHPS) |
| 658 |
|
| 659 |
integer i |
| 660 |
|
| 661 |
IF(NCHP.GE.0) THEN |
| 662 |
DO I=1,NCHP |
| 663 |
B(I) = A(IGRD(I),1) |
| 664 |
ENDDO |
| 665 |
ELSE |
| 666 |
PRINT *, 'ERROR IN GRD2MSC' |
| 667 |
ENDIF |
| 668 |
|
| 669 |
RETURN |
| 670 |
END |
| 671 |
|
| 672 |
SUBROUTINE MSC2GRD(IGRD,CHFR,B,MXCHPS,NCHP,FRACG,A,IM,JM) |
| 673 |
|
| 674 |
implicit none |
| 675 |
_RL zero,one |
| 676 |
parameter ( one = 1.) |
| 677 |
parameter (zero = 0.) |
| 678 |
integer im,jm,mxchps,nchp |
| 679 |
integer igrd(mxchps) |
| 680 |
_RL A(IM,JM), B(MXCHPS), CHFR(MXCHPS), FRACG(IM,JM) |
| 681 |
|
| 682 |
_RL VT1(IM,JM) |
| 683 |
integer i |
| 684 |
|
| 685 |
IF(NCHP.GE.0) THEN |
| 686 |
DO I=1,IM*JM |
| 687 |
VT1(I,1) = ZERO |
| 688 |
ENDDO |
| 689 |
|
| 690 |
DO I=1,NCHP |
| 691 |
VT1(IGRD(I),1) = VT1(IGRD(I),1) + B(I)*CHFR(I) |
| 692 |
ENDDO |
| 693 |
|
| 694 |
DO I=1,IM*JM |
| 695 |
A(I,1) = A(I,1)*(ONE-FRACG(I,1)) + VT1(I,1) |
| 696 |
ENDDO |
| 697 |
ELSE |
| 698 |
PRINT *, 'ERROR IN MSC2GRD' |
| 699 |
ENDIF |
| 700 |
|
| 701 |
RETURN |
| 702 |
END |
| 703 |
|
| 704 |
subroutine chpprm(nymd,nhms,mxchps,nchp,chlt,ityp,alai, |
| 705 |
1 agrn,zoch,z2ch,cdrc,cdsc,sqsc,ufac,rsl1,rsl2,rdcs) |
| 706 |
|
| 707 |
implicit none |
| 708 |
|
| 709 |
integer nymd,nhms,nchp,mxchps,ityp(mxchps) |
| 710 |
_RL chlt(mxchps) |
| 711 |
_RL alai(mxchps),agrn(mxchps) |
| 712 |
_RL zoch(mxchps), z2ch(mxchps), cdrc(mxchps), cdsc(mxchps) |
| 713 |
_RL sqsc(mxchps), ufac(mxchps), rsl1(mxchps), rsl2(mxchps) |
| 714 |
_RL rdcs(mxchps) |
| 715 |
|
| 716 |
C********************************************************************* |
| 717 |
C********************* SUBROUTINE CHPPRM **************************** |
| 718 |
C********************** 14 JUNE 1991 ****************************** |
| 719 |
C********************************************************************* |
| 720 |
|
| 721 |
integer ntyps |
| 722 |
parameter (ntyps=10) |
| 723 |
|
| 724 |
_RL pblzet |
| 725 |
parameter (pblzet = 50.) |
| 726 |
integer k1,k2,nymd1,nhms1,nymd2,nhms2,i |
| 727 |
_RL getcon,vkrm,rootl,vroot,dum1,dum2,alphaf |
| 728 |
_RL facm,facp |
| 729 |
_RL scat,d |
| 730 |
|
| 731 |
_RL |
| 732 |
& vgdd(12, ntyps), vgz0(12, ntyps), |
| 733 |
& vgrd(12, ntyps), vgrt(12, ntyps), |
| 734 |
|
| 735 |
& vgrf11(ntyps), vgrf12(ntyps), |
| 736 |
& vgtr11(ntyps), vgtr12(ntyps), |
| 737 |
& vgroca(ntyps), vgrotd(ntyps), |
| 738 |
& vgrdrs(ntyps), vgz2 (ntyps) |
| 739 |
|
| 740 |
|
| 741 |
data vgz0 / |
| 742 |
1 2.6530, 2.6530, 2.6530, 2.6530, 2.6530, 2.6530, 2.6530, |
| 743 |
1 2.6530, 2.6530, 2.6530, 2.6530, 2.6530, |
| 744 |
2 0.5200, 0.5200, 0.6660, 0.9100, 1.0310, 1.0440, 1.0420, |
| 745 |
2 1.0370, 1.0360, 0.9170, 0.6660, 0.5200, |
| 746 |
3 1.1120, 1.1030, 1.0880, 1.0820, 1.0760, 1.0680, 1.0730, |
| 747 |
3 1.0790, 1.0820, 1.0880, 1.1030, 1.1120, |
| 748 |
4 0.0777, 0.0778, 0.0778, 0.0779, 0.0778, 0.0771, 0.0759, |
| 749 |
4 0.0766, 0.0778, 0.0779, 0.0778, 0.0778, |
| 750 |
5 0.2450, 0.2450, 0.2270, 0.2000, 0.2000, 0.2000, 0.2000, |
| 751 |
5 0.267, 0.292, 0.280, 0.258, 0.2450, |
| 752 |
6 0.0752, 0.0752, 0.0752, 0.0752, 0.0752, 0.0757, 0.0777, |
| 753 |
6 0.0778, 0.0774, 0.0752, 0.0752, 0.0752, |
| 754 |
7 0.0112, 0.0112, 0.0112, 0.0112, 0.0112, 0.0112, 0.0112, |
| 755 |
7 0.0112, 0.0112, 0.0112, 0.0112, 0.0112, |
| 756 |
8 0.0112, 0.0112, 0.0112, 0.0112, 0.0112, 0.0112, 0.0112, |
| 757 |
8 0.0112, 0.0112, 0.0112, 0.0112, 0.0112, |
| 758 |
9 0.0112, 0.0112, 0.0112, 0.0112, 0.0112, 0.0112, 0.0112, |
| 759 |
9 0.0112, 0.0112, 0.0112, 0.0112, 0.0112, |
| 760 |
1 0.0112, 0.0112, 0.0112, 0.0112, 0.0112, 0.0112, 0.0112, |
| 761 |
1 0.0112, 0.0112, 0.0112, 0.0112, 0.0112 |
| 762 |
& / |
| 763 |
|
| 764 |
|
| 765 |
data vgrd / |
| 766 |
1 285.87, 285.87, 285.87, 285.87, 285.87, 285.87, 285.87, |
| 767 |
1 285.87, 285.87, 285.87, 285.87, 285.87, |
| 768 |
2 211.32, 211.32, 218.78, 243.40, 294.87, 345.90, 355.18, |
| 769 |
2 341.84, 307.22, 244.84, 218.78, 211.32, |
| 770 |
3 565.41, 587.05, 623.46, 638.13, 652.86, 675.04, 660.24, |
| 771 |
3 645.49, 638.13, 623.46, 587.05, 565.41, |
| 772 |
4 24.43, 24.63, 24.80, 24.96, 25.72, 27.74, 30.06, |
| 773 |
4 28.86, 25.90, 25.11, 24.80, 24.63, |
| 774 |
5 103.60, 103.60, 102.35, 100.72, 100.72, 100.72, 100.72, |
| 775 |
5 105.30, 107.94, 106.59, 104.49, 103.60, |
| 776 |
6 22.86, 22.86, 22.86, 22.86, 22.86, 23.01, 24.36, |
| 777 |
6 24.69, 24.04, 22.86, 22.86, 22.86, |
| 778 |
7 23.76, 23.76, 23.76, 23.76, 23.76, 23.76, 23.76, |
| 779 |
7 23.76, 23.76, 23.76, 23.76, 23.76, |
| 780 |
8 23.76, 23.76, 23.76, 23.76, 23.76, 23.76, 23.76, |
| 781 |
8 23.76, 23.76, 23.76, 23.76, 23.76, |
| 782 |
9 23.76, 23.76, 23.76, 23.76, 23.76, 23.76, 23.76, |
| 783 |
9 23.76, 23.76, 23.76, 23.76, 23.76, |
| 784 |
1 23.76, 23.76, 23.76, 23.76, 23.76, 23.76, 23.76, |
| 785 |
1 23.76, 23.76, 23.76, 23.76, 23.76 |
| 786 |
& / |
| 787 |
|
| 788 |
data vgrt / |
| 789 |
1 19737.8, 19737.8, 19737.8, 19737.8, 19737.8, 19737.8, 19737.8, |
| 790 |
1 19737.8, 19737.8, 19737.8, 19737.8, 19737.8, |
| 791 |
2 5010.0, 5010.0, 5270.0, 6200.0, 8000.0, 9700.0, 9500.0, |
| 792 |
2 8400.0, 6250.0, 5270.0, 5010.0, 5010.0, |
| 793 |
3 9000.0, 9200.0, 9533.3, 9666.7, 9800.0, 9866.7, 9733.3, |
| 794 |
3 9666.7, 9533.3, 9200.0, 9000.0, 9000.0, |
| 795 |
4 5500.0, 5625.0, 5750.0, 5875.0, 6625.0, 8750.0, 9375.0, |
| 796 |
4 6875.0, 6000.0, 5750.0, 5625.0, 5500.0, |
| 797 |
5 6500.0, 6000.0, 5500.0, 5500.0, 5500.0, 5500.0, 5500.0, |
| 798 |
5 7500.0, 8500.0, 7000.0, 6500.0, 6500.0, |
| 799 |
6 10625.0, 10625.0, 10625.0, 10625.0, 10625.0, 11250.0, 18750.0, |
| 800 |
6 17500.0, 10625.0, 10625.0, 10625.0, 10625.0, |
| 801 |
7 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, |
| 802 |
7 1.0, 1.0, 1.0, 1.0, 1.0, |
| 803 |
8 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, |
| 804 |
8 1.0, 1.0, 1.0, 1.0, 1.0, |
| 805 |
9 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, |
| 806 |
9 1.0, 1.0, 1.0, 1.0, 1.0, |
| 807 |
1 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, |
| 808 |
1 1.0, 1.0, 1.0, 1.0, 1.0 |
| 809 |
& / |
| 810 |
|
| 811 |
|
| 812 |
data vgdd / |
| 813 |
1 27.37, 27.37, 27.37, 27.37, 27.37, 27.37, 27.37, |
| 814 |
1 27.37, 27.37, 27.37, 27.37, 27.37, |
| 815 |
2 13.66, 13.66, 14.62, 15.70, 16.33, 16.62, 16.66, |
| 816 |
2 16.60, 16.41, 15.73, 14.62, 13.66, |
| 817 |
3 13.76, 13.80, 13.86, 13.88, 13.90, 13.93, 13.91, |
| 818 |
3 13.89, 13.88, 13.86, 13.80, 13.76, |
| 819 |
4 0.218, 0.227, 0.233, 0.239, 0.260, 0.299, 0.325, |
| 820 |
4 0.313, 0.265, 0.244, 0.233, 0.227, |
| 821 |
5 2.813, 2.813, 2.662, 2.391, 2.391, 2.391, 2.391, |
| 822 |
5 2.975, 3.138, 3.062, 2.907, 2.813, |
| 823 |
6 0.10629, 0.10629, 0.10629, 0.10629, 0.10629, 0.12299, 0.21521, |
| 824 |
6 0.22897, 0.19961, 0.10629, 0.10629, 0.10629, |
| 825 |
7 0.0001, 0.0001, 0.0001, 0.0001, 0.0001, 0.0001, 0.0001, |
| 826 |
7 0.0001, 0.0001, 0.0001, 0.0001, 0.0001, |
| 827 |
8 0.0001, 0.0001, 0.0001, 0.0001, 0.0001, 0.0001, 0.0001, |
| 828 |
8 0.0001, 0.0001, 0.0001, 0.0001, 0.0001, |
| 829 |
9 0.0001, 0.0001, 0.0001, 0.0001, 0.0001, 0.0001, 0.0001, |
| 830 |
9 0.0001, 0.0001, 0.0001, 0.0001, 0.0001, |
| 831 |
1 0.0001, 0.0001, 0.0001, 0.0001, 0.0001, 0.0001, 0.0001, |
| 832 |
1 0.0001, 0.0001, 0.0001, 0.0001, 0.0001 |
| 833 |
& / |
| 834 |
|
| 835 |
|
| 836 |
data vgrf11 /0.10,0.10,0.07,0.105,0.10,0.10,.001,.001,.001,.001/ |
| 837 |
|
| 838 |
data vgrf12 /0.16,0.16,0.16,0.360,0.16,0.16,.001,.001,.001,.001/ |
| 839 |
|
| 840 |
data vgtr11 /0.05,0.05,0.05,0.070,0.05,0.05,.001,.001,.001,.001/ |
| 841 |
|
| 842 |
data vgtr12 /.001,.001,.001, .220,.001,.001,.001,.001,.001,.001/ |
| 843 |
|
| 844 |
data vgroca / |
| 845 |
& 0.384E-6, 0.384E-6, 0.384E-6, 0.384E-6, 0.384E-6, 0.384E-6, |
| 846 |
& .1E-6, .1E-6, .1E-6, .1E-6 / |
| 847 |
|
| 848 |
data vgrotd /1.00,1.00,0.50,0.50,0.50,0.20,0.10,0.10,0.10,0.10/ |
| 849 |
|
| 850 |
data vgrdrs / |
| 851 |
& 0.75E13, 0.75E13, 0.75E13, 0.40E13, 0.75E13, 0.75E13, |
| 852 |
& 0.10E13, 0.10E13, 0.10E13, 0.10E13 / |
| 853 |
|
| 854 |
data vgz2 /35.0, 20.0, 17.0, 0.6, 5.0, 0.6, 0.1, 0.1, 0.1, 0.1/ |
| 855 |
|
| 856 |
vkrm = GETCON('VON KARMAN') |
| 857 |
|
| 858 |
call time_bound ( nymd,nhms, nymd1,nhms1, nymd2,nhms2, k1,k2 ) |
| 859 |
call interp_time ( nymd,nhms, nymd1,nhms1, nymd2,nhms2, facm,facp) |
| 860 |
|
| 861 |
do i=1,nchp |
| 862 |
|
| 863 |
zoch(i) = vgz0(k2,ityp(i))*facp + vgz0(k1,ityp(i))*facm |
| 864 |
rdcs(i) = vgrd(k2,ityp(i))*facp + vgrd(k1,ityp(i))*facm |
| 865 |
|
| 866 |
rootl = vgrt(k2,ityp(i))*facp + vgrt(k1,ityp(i))*facm |
| 867 |
|
| 868 |
vroot = rootl * vgroca(ityp (i)) |
| 869 |
dum1 = log (vroot / (1. - vroot)) |
| 870 |
dum2 = 1. / (8. * 3.14159 * rootl) |
| 871 |
alphaf = dum2 * (vroot - 3. -2. * dum1) |
| 872 |
|
| 873 |
rsl1(i) = vgrdrs (ityp (i)) / (rootl * vgrotd (ityp (i))) |
| 874 |
rsl2(i) = alphaf / vgrotd (ityp (i)) |
| 875 |
|
| 876 |
scat = agrn(i) *(vgtr11(ityp(i)) + vgrf11(ityp(i))) |
| 877 |
& + (1. - agrn(i))*(vgtr12(ityp(i)) + vgrf12(ityp(i))) |
| 878 |
sqsc(i) = sqrt (1. - scat) |
| 879 |
|
| 880 |
d = vgdd(k2,ityp(i))*facp + vgdd(k1,ityp(i))*facm |
| 881 |
ufac(i) = log( (vgz2(ityp(i)) - d) / zoch(i) ) |
| 882 |
* / log( pblzet / zoch(i) ) |
| 883 |
|
| 884 |
z2ch(i) = vgz2(ityp (i)) |
| 885 |
|
| 886 |
cdsc(i) = pblzet/zoch(i)+1. |
| 887 |
cdrc(i) = vkrm/log(cdsc(i)) |
| 888 |
cdrc(i) = cdrc(i)*cdrc(i) |
| 889 |
cdsc(i) = sqrt(cdsc(i)) |
| 890 |
cdsc(i) = cdrc(i)*cdsc(i) |
| 891 |
|
| 892 |
enddo |
| 893 |
|
| 894 |
return |
| 895 |
end |
| 896 |
|
| 897 |
subroutine pkappa (im,jm,lm,ple,pkle,pkz) |
| 898 |
C*********************************************************************** |
| 899 |
C Purpose |
| 900 |
C Calculate Phillips P**Kappa |
| 901 |
C |
| 902 |
C Arguments |
| 903 |
C PLE .... edge-level pressure |
| 904 |
C PKLE ... edge-level pressure**kappa |
| 905 |
C IM ..... longitude dimension |
| 906 |
C JM ..... latitude dimension |
| 907 |
C LM ..... vertical dimension |
| 908 |
C PKZ .... mid-level pressure**kappa |
| 909 |
C*********************************************************************** |
| 910 |
implicit none |
| 911 |
|
| 912 |
integer im,jm,lm |
| 913 |
_RL ple(im,jm,lm+1) |
| 914 |
_RL pkle(im,jm,lm+1) |
| 915 |
_RL pkz(im,jm,lm) |
| 916 |
|
| 917 |
_RL akap1,getcon |
| 918 |
integer i,j,L |
| 919 |
|
| 920 |
akap1 = 1.0 + getcon('KAPPA') |
| 921 |
|
| 922 |
do L = 1,lm |
| 923 |
do j = 1,jm |
| 924 |
do i = 1,im |
| 925 |
pkz(i,j,L) = ( ple (i,j,l+1)*pkle(i,j,l+1) |
| 926 |
. - ple (i,j,l)*pkle(i,j,l) ) |
| 927 |
. / ( akap1* (ple (i,j,l+1)-ple (i,j,l)) ) |
| 928 |
enddo |
| 929 |
enddo |
| 930 |
enddo |
| 931 |
|
| 932 |
return |
| 933 |
end |
Parent Directory
| 
Revision Log
| 
Revision Graph