*IF DEF,A06_3A,OR,DEF,A06_3B ADR2F405.2
C ******************************COPYRIGHT****************************** GTS2F400.3673
C (c) CROWN COPYRIGHT 1995, METEOROLOGICAL OFFICE, All Rights Reserved. GTS2F400.3674
C GTS2F400.3675
C Use, duplication or disclosure of this code is subject to the GTS2F400.3676
C restrictions as set forth in the contract. GTS2F400.3677
C GTS2F400.3678
C Meteorological Office GTS2F400.3679
C London Road GTS2F400.3680
C BRACKNELL GTS2F400.3681
C Berkshire UK GTS2F400.3682
C RG12 2SZ GTS2F400.3683
C GTS2F400.3684
C If no contract has been raised with this copy of the code, the use, GTS2F400.3685
C duplication or disclosure of it is strictly prohibited. Permission GTS2F400.3686
C to do so must first be obtained in writing from the Head of Numerical GTS2F400.3687
C Modelling at the above address. GTS2F400.3688
C ******************************COPYRIGHT****************************** GTS2F400.3689
C GTS2F400.3690
! SUBROUTINE GW_SCOR TO DELINEATE VERT SCORER PROFILE FOR GW_LEE/GW_VERT GWSCOR3A.3
! GWSCOR3A.4
SUBROUTINE GW_SCOR 2GWSCOR3A.5
1 (PSTAR,PEXNER,THETA,U,V,LEVELS,START_L,H_JUMP,POINTS,AKH,BKH GWSCOR3A.6
2 ,UNIT_X,UNIT_Y,TRANS,K_LEE,H_LEE,LSQ_LEE,L_LEE) GWSCOR3A.7
GWSCOR3A.8
IMPLICIT NONE GWSCOR3A.9
! Description: TO CALCULATE THE PROFILE OF THE SCORER PARAMETER GWSCOR3A.10
! TO CALCULATE LEE HEIGHT AND TRANSMITTION COEFFICIENT GWSCOR3A.11
! GWSCOR3A.12
! Method: UNIFIED MODEL DOCUMENTATION PAPER NO. ? GWSCOR3A.13
! THE EQUATIONS USED ARE eq 11,30,49. GWSCOR3A.14
! GWSCOR3A.15
! Current code owner: S.Webster ASW1F403.49
ASW1F403.50
! GWSCOR3A.17
! History: GWSCOR3A.18
! Version Date Comment GWSCOR3A.19
! 3.4 18/10/94 Original Code. J.R.Mitchell GWSCOR3A.20
! 4.4 19/09/97 Remove *IF -DEF,CRAY compile options. S.Webster ASW1F404.9
!LL 4.5 14/05/98 Add Fujitsu directive to stop potential zero GRB0F405.144
!LL divide caused by optimizer. RBarnes@ecmwf.int GRB0F405.145
! GWSCOR3A.21
! Code Description: GWSCOR3A.22
! Language: Fortran 77 + common extensions GWSCOR3A.23
! This code is written to UMDP3 v6 programming standards. GWSCOR3A.24
! System component covered: ORIGINAL VERSION FOR CRAY Y-MP GWSCOR3A.25
! System task covered: PART OF P22 GWSCOR3A.26
! SUITABLE FOR SINGLE COLUMN USE,ROTATED GRIDS GWSCOR3A.27
! FURTHER ALTERATIONS MAY BE REQUIRED FOR AUTOTASKING EFFICIENCY GWSCOR3A.28
GWSCOR3A.29
! Global Variables GWSCOR3A.30
*CALL C_G 
GWSCOR3A.31
*CALL C_R_CP GWSCOR3A.32
! Local constants GWSCOR3A.33
*CALL C_GWAVE GWSCOR3A.34
! Subroutine arguements: GWSCOR3A.35
GWSCOR3A.36
INTEGER GWSCOR3A.37
* LEVELS !IN NUMBER OF MODEL LEVELS GWSCOR3A.38
*,START_L !IN START LEVEL FOR GWD CALCULATIONS GWSCOR3A.39
*,POINTS !IN NUMBER OF POINTS GWSCOR3A.40
*,K_LEE(POINTS,2) !OUT MODEL LEVEL OF LEE HEIGHT AND TOP GWSCOR3A.41
GWSCOR3A.42
REAL GWSCOR3A.43
* PSTAR(POINTS) !IN PSTAR FIELD GWSCOR3A.44
*,PEXNER(POINTS,LEVELS+1) !IN PEXNER GWSCOR3A.45
*,THETA(POINTS,LEVELS) !IN THETA FIELD GWSCOR3A.46
*,U(POINTS,LEVELS) !IN U FIELD GWSCOR3A.47
*,V(POINTS,LEVELS) !IN V FIELD GWSCOR3A.48
*,UNIT_X(POINTS) !IN X_COMPNT OF UNIT STRESS VECTOR GWSCOR3A.49
*,UNIT_Y(POINTS) !IN Y_COMPNT OF UNIT STRESS VECTOR GWSCOR3A.50
*,TRANS(POINTS) !OUT TRANSMITION COEFFICIENT GWSCOR3A.51
*,LSQ_LEE(POINTS,2) !OUT L1 SQUARED AND L2 SQUARED FOR GWSCOR3A.52
* ! LEE WAVE CALCULATIONS GWSCOR3A.53
*,H_LEE(POINTS) !OUT LEE HEIGHT TO TOP OF LEVEL GWSCOR3A.54
! AKH,BKH DEFINE HYBRID VERTICAL COORDINATES P=A+BP*-LAYER EDGES, GWSCOR3A.55
! DELTA_AK,DELTA_BK DEFINE PRESSURE DIFFERENCES ACROSS LAYERS GWSCOR3A.56
*,AKH(LEVELS+1) !IN VALUE AT LAYER BOUNDARY GWSCOR3A.57
*,BKH(LEVELS+1) !IN VALUE AT LAYER BOUMDARY GWSCOR3A.58
GWSCOR3A.59
LOGICAL GWSCOR3A.60
* H_JUMP(POINTS) !IN TRUE IF HYDROLIC JUMP REGIME GWSCOR3A.61
*,L_LEE(POINTS) !OUT TRUE IF LEE WAVE/TRANS POINT GWSCOR3A.62
GWSCOR3A.63
! Local parameters GWSCOR3A.64
REAL CPBYG GWSCOR3A.65
PARAMETER(CPBYG=CP/G) GWSCOR3A.66
! Local scalers GWSCOR3A.67
REAL GWSCOR3A.68
* DEL_EXNER,DEL_EXNER_KL,DEL_EXNER_KU, ! EXNER DIFFERENCE GWSCOR3A.69
* PU,PL GWSCOR3A.70
GWSCOR3A.71
LOGICAL FLAG GWSCOR3A.72
GWSCOR3A.73
INTEGER I,K,M ! LOOP COUNTER IN ROUTINE GWSCOR3A.74
INTEGER KK,KL,KU ! LEVEL COUNTERS IN ROUTINE GWSCOR3A.75
GWSCOR3A.76
! Local dynamic arrays GWSCOR3A.77
! LOCAL WORKSPACE ARRAYS: LEVELS+8 ARRAYS OF FULL FIELD LENGTH GWSCOR3A.78
! GWSCOR3A.79
INTEGER GWSCOR3A.85
* K_LEE_LIM(LEVELS-3) ! TOP MODEL LEVEL HEIGHT FOR CALCULATION GWSCOR3A.86
* ! OF L_SQ(I,2), FOR EACH POSSIBLE LEE LEVEL GWSCOR3A.87
*,K_LEE_MAX ! MAXIMUM LEVEL FOR FINDING LEE WAVE HEIGHT GWSCOR3A.88
*,K_TRANS ! LEVEL OF L_SQ SPLIT FOR CALCULATION OF GWSCOR3A.89
* ! TRANSMITION COEFFICIENT GWSCOR3A.90
GWSCOR3A.91
REAL GWSCOR3A.92
* ALPHA_L,PHASE ! LINKED LEE WAVE PARAMETER GWSCOR3A.93
*,DELTA_Z ! THICKNESS OF LAYER BETWEEN HALF LEVELS GWSCOR3A.94
*,P_EXNER_CENTRE(POINTS,2) ! EXNER PRESSURE AT LAYER CENTRES GWSCOR3A.95
*,N_SQ(POINTS,2) ! SQUARE OF BRUNT_VAISALA FREQUENCY AT GWSCOR3A.96
* ! LAYER BOUNDARIES GWSCOR3A.97
*,N_SQAV ! CENTRE LAYER AV. OF BOUYANCY FREQ SQ GWSCOR3A.98
*,ZH(POINTS) ! TOTAL HEIGHT OF JUMP CALCUALTION GWSCOR3A.99
*,P(LEVELS-3) ! PRESSURE AT LEVEL TOP BOUNDARIES GWSCOR3A.100
* ! FOR STANDARD POINT GWSCOR3A.101
*,P_LIM(LEVELS-3) ! PRESSURE TOP BOUNDARIES FOR GWSCOR3A.102
* ! K_LEE_LIM GWSCOR3A.103
*,UKP1 ! COMPONENT OF WIND AT LAYER K+1 IN GWSCOR3A.104
* ! DIRN. OF SURFACE STRESS GWSCOR3A.105
*,UK ! COMPONENT OF WIND AT LAYER K GWSCOR3A.106
*,UKM1 ! COMPONENT OF WIND AT LAYER K-1 GWSCOR3A.107
*,DU_DZ(POINTS,2) ! DELTAN OVER DELTAU - LAYER BOUNDARY GWSCOR3A.108
*,LSQ_SUM(POINTS,2) ! L1 SQUARED AND L2 SQUARED SUMMED GWSCOR3A.109
*,LSQ(POINTS,LEVELS-3) ! L SQUARED AT A LEVEL. GWSCOR3A.110
*,LSQ1, LSQ2, L1, L2 ! L CALCULATIONS FOR TRANSMITTION COEFF GWSCOR3A.111
*,CONST,CALC ! CALCULATIONS GWSCOR3A.112
GWSCOR3A.114
! Function and subroutine calls GWSCOR3A.115
*CALL P_EXNERC GWSCOR3A.116
GWSCOR3A.117
!------------------------------------------------------------------- GWSCOR3A.118
! 1. Parameter calculation. Tan(0.75*PI) = -1.0 GWSCOR3A.119
! Note ATAN produces answer in range -PI/2 to PI/2 GWSCOR3A.120
!-------------------------------------------------------------------- GWSCOR3A.121
PHASE=LEE_PHASE*3.14159 GWSCOR3A.122
ALPHA_L= TAN(PHASE) GWSCOR3A.123
CONST= ALPHA_L/( PHASE*SQRT( 1+ALPHA_L**2 ) ) GWSCOR3A.124
KL=1 GWSCOR3A.125
KU=2 GWSCOR3A.126
!--------------------------------------------------------------------- GWSCOR3A.127
! Find approximate level of 0.6 tropopause height - maximum lee limit GWSCOR3A.128
! Find approx half tropopause height - transmittion coefficient top GWSCOR3A.129
!--------------------------------------------------------------------- GWSCOR3A.130
FLAG = .TRUE. GWSCOR3A.131
K_LEE_MAX = LEVELS-4 GWSCOR3A.132
K_TRANS = K_LEE_MAX GWSCOR3A.133
DO K= 3,LEVELS-4 GWSCOR3A.134
IF (FLAG) THEN GWSCOR3A.135
PU=100000.*BKH(K+1) + AKH(K+1) GWSCOR3A.136
IF ( PU .GT. 50000. ) THEN GWSCOR3A.137
K_TRANS = K GWSCOR3A.138
ELSE IF ( PU .LT. 45000. ) THEN GWSCOR3A.139
K_LEE_MAX = K GWSCOR3A.140
FLAG = .FALSE. GWSCOR3A.141
ENDIF GWSCOR3A.142
ENDIF GWSCOR3A.143
ENDDO GWSCOR3A.144
!--------------------------------------------------------------------- GWSCOR3A.145
! Find array of top level scorer calculation limits K_LEE_LIM, for each GWSCOR3A.146
! possible lee level by using the same standard point as above. GWSCOR3A.147
!--------------------------------------------------------------------- GWSCOR3A.148
DO K= 2,LEVELS-3 GWSCOR3A.149
P(K)=100000.*BKH(K+1) + AKH(K+1) GWSCOR3A.150
IF ( K .LE. K_LEE_MAX ) THEN GWSCOR3A.151
P_LIM(K) = 2*P(K) - 100000. GWSCOR3A.152
K_LEE_LIM(K) = LEVELS-3 GWSCOR3A.153
ENDIF GWSCOR3A.154
ENDDO GWSCOR3A.155
GWSCOR3A.156
DO M= 2,K_LEE_MAX GWSCOR3A.157
DO K= M+1,LEVELS-3 GWSCOR3A.158
IF ( P(K) .LT. P_LIM(M) .AND. GWSCOR3A.159
* K_LEE_LIM(M).EQ.LEVELS-3 ) THEN GWSCOR3A.160
K_LEE_LIM(M) = K GWSCOR3A.161
ENDIF GWSCOR3A.162
ENDDO GWSCOR3A.163
ENDDO GWSCOR3A.164
K_LEE_LIM(1)=2 GWSCOR3A.165
!--------------------------------------------------------------------- GWSCOR3A.166
! Initialisation GWSCOR3A.167
!--------------------------------------------------------------------- GWSCOR3A.168
DO I=1,POINTS GWSCOR3A.169
K_LEE(I,1)=0 GWSCOR3A.170
K_LEE(I,2)=0 GWSCOR3A.171
IF( H_JUMP(I) ) THEN GWSCOR3A.172
L_LEE(I)=.FALSE. GWSCOR3A.173
ELSE GWSCOR3A.174
L_LEE(I)=.TRUE. GWSCOR3A.175
ENDIF GWSCOR3A.176
TRANS(I)=1.0 GWSCOR3A.177
ZH(I)=0.0 GWSCOR3A.178
ENDDO GWSCOR3A.179
GWSCOR3A.180
!--------------------------------------------------------------------- GWSCOR3A.181
! Bottom Level calculation - Calculates L_SQ at level 2 GWSCOR3A.182
!--------------------------------------------------------------------- GWSCOR3A.183
DO I=1,POINTS GWSCOR3A.184
IF( L_LEE(I) ) THEN GWSCOR3A.185
PU=PSTAR(I)*BKH(2) + AKH(2) GWSCOR3A.186
PL=PSTAR(I)*BKH(1) + AKH(1) GWSCOR3A.187
!-------------- KU / KL reversed ready for next stage------------ GWSCOR3A.188
P_EXNER_CENTRE(I,KU)= GWSCOR3A.189
& P_EXNER_C( PEXNER(I,2),PEXNER(I,1),PU,PL,KAPPA) GWSCOR3A.190
PL=PU GWSCOR3A.191
PU=PSTAR(I)*BKH(3) + AKH(3) GWSCOR3A.192
P_EXNER_CENTRE(I,KL)= GWSCOR3A.193
& P_EXNER_C( PEXNER(I,3),PEXNER(I,2),PU,PL,KAPPA) GWSCOR3A.194
!------------------------------------------------------------- GWSCOR3A.195
N_SQ(I,KL) = G*(THETA(I,2)-THETA(I,1))/(THETA(I,2)* GWSCOR3A.196
& THETA(I,1)*(P_EXNER_CENTRE(I,KU)-P_EXNER_CENTRE(I,KL)) GWSCOR3A.197
& *CPBYG) GWSCOR3A.198
PL=PU GWSCOR3A.199
PU=PSTAR(I)*BKH(4) + AKH(4) GWSCOR3A.200
P_EXNER_CENTRE(I,KU)= GWSCOR3A.201
& P_EXNER_C( PEXNER(I,4),PEXNER(I,3),PU,PL,KAPPA) GWSCOR3A.202
N_SQ(I,KU) = G*(THETA(I,3)-THETA(I,2))/(THETA(I,3)* GWSCOR3A.203
& THETA(I,2)*(P_EXNER_CENTRE(I,KL)-P_EXNER_CENTRE(I,KU)) GWSCOR3A.204
& *CPBYG) GWSCOR3A.205
!-------------------------------------------------------------------- GWSCOR3A.206
! N_SQAV at layer centre, interpolated from N_SQ(KL/KU) at boundaries GWSCOR3A.207
!-------------------------------------------------------------------- GWSCOR3A.208
N_SQAV = ( (PEXNER(I,2)-P_EXNER_CENTRE(I,KL))*N_SQ(I,KU) GWSCOR3A.209
& + (P_EXNER_CENTRE(I,KL) - PEXNER(I,3))*N_SQ(I,KL) ) GWSCOR3A.210
& / ( PEXNER(I,2) - PEXNER(I,3) ) GWSCOR3A.211
!-------------------------------------------------------------------- GWSCOR3A.212
! Note U is component parallel to stress vector GWSCOR3A.213
!-------------------------------------------------------------------- GWSCOR3A.214
UKP1 = U(I,3)*UNIT_X(I) + V(I,3)*UNIT_Y(I) GWSCOR3A.215
UK = U(I,2)*UNIT_X(I) + V(I,2)*UNIT_Y(I) GWSCOR3A.216
UKM1 = U(I,1)*UNIT_X(I) + V(I,1)*UNIT_Y(I) GWSCOR3A.217
IF( N_SQAV .LE. 0.0 .OR. UKP1.LE.0.3 GWSCOR3A.218
* .OR. UK.LE.0.2 .OR. UKM1.LE.0.1 ) THEN GWSCOR3A.219
L_LEE(I)=.FALSE. GWSCOR3A.220
ELSE GWSCOR3A.221
DEL_EXNER_KU = PEXNER(I,3) - P_EXNER_CENTRE(I,KU) GWSCOR3A.222
DEL_EXNER_KL = P_EXNER_CENTRE(I,KL) - PEXNER(I,3) GWSCOR3A.223
DELTA_Z= CPBYG*(THETA(I,2)*DEL_EXNER_KL + GWSCOR3A.224
* THETA(I,3)*DEL_EXNER_KU) GWSCOR3A.225
DU_DZ(I,KU)= ( UKP1-UK )/ DELTA_Z GWSCOR3A.226
PU=PSTAR(I)*BKH(2) + AKH(2) GWSCOR3A.227
PL=PSTAR(I)*BKH(1) + AKH(1) GWSCOR3A.228
CALC= GWSCOR3A.229
& P_EXNER_C( PEXNER(I,2),PEXNER(I,1),PU,PL,KAPPA) GWSCOR3A.230
DEL_EXNER_KU = PEXNER(I,2) - P_EXNER_CENTRE(I,KL) GWSCOR3A.231
DEL_EXNER_KL = CALC - PEXNER(I,2) GWSCOR3A.232
DELTA_Z= CPBYG*(THETA(I,1)*DEL_EXNER_KL + GWSCOR3A.233
* THETA(I,2)*DEL_EXNER_KU ) GWSCOR3A.234
GWSCOR3A.235
DU_DZ(I,KL)= ( UK-UKM1 )/ DELTA_Z GWSCOR3A.236
DEL_EXNER_KU = PEXNER(I,2) - PEXNER(I,3) GWSCOR3A.237
DEL_EXNER_KL = PEXNER(I,1) - PEXNER(I,2) GWSCOR3A.238
ZH(I) = CPBYG*(THETA(I,1)*DEL_EXNER_KL + GWSCOR3A.239
* THETA(I,2)*DEL_EXNER_KU ) GWSCOR3A.240
LSQ(I,2)=N_SQAV /UK**2 - ( DU_DZ(I,KU)-DU_DZ(I,KL) ) GWSCOR3A.241
* / (UK*CPBYG*THETA(I,2)*DEL_EXNER_KU) GWSCOR3A.242
LSQ_SUM(I,1)= 0.0 GWSCOR3A.243
LSQ_SUM(I,2)= LSQ(I,2) GWSCOR3A.244
ENDIF GWSCOR3A.245
ENDIF ! lee point GWSCOR3A.246
ENDDO ! Points GWSCOR3A.247
GWSCOR3A.248
!--------------------------------------------------------------------- GWSCOR3A.249
! Loop Levels GWSCOR3A.250
!--------------------------------------------------------------------- GWSCOR3A.251
DO M= 2,K_LEE_MAX GWSCOR3A.252
DO K= K_LEE_LIM(M-1), K_LEE_LIM(M) GWSCOR3A.253
GWSCOR3A.254
IF( K.NE.K_LEE_LIM(M-1) ) THEN GWSCOR3A.255
KK=KL GWSCOR3A.256
KL=KU GWSCOR3A.257
KU=KK GWSCOR3A.258
ENDIF GWSCOR3A.259
GWSCOR3A.260
! Optimizer precomputes values which can cause zero divide, so:- GRB0F405.146
!OCL NOPREEX,NOEVAL GRB0F405.147
DO I=1,POINTS GWSCOR3A.261
IF( L_LEE(I) .AND. ( K.NE.K_LEE_LIM(M-1) ) GWSCOR3A.262
& .AND. ( K_LEE(I,1).EQ.0 .OR. M.LE.K_TRANS ) ) THEN GWSCOR3A.263
! next level stage GWSCOR3A.264
PU=PSTAR(I)*BKH(K+2) + AKH(K+2) GWSCOR3A.265
PL=PSTAR(I)*BKH(K+1) + AKH(K+1) GWSCOR3A.266
P_EXNER_CENTRE(I,KU)= GWSCOR3A.267
& P_EXNER_C( PEXNER(I,K+2),PEXNER(I,K+1),PU,PL,KAPPA) GWSCOR3A.268
N_SQ(I,KU) = G*(THETA(I,K+1)-THETA(I,K))/(THETA(I,K+1)* GWSCOR3A.269
& THETA(I,K)*(P_EXNER_CENTRE(I,KL)-P_EXNER_CENTRE(I,KU))* GWSCOR3A.270
& CPBYG) GWSCOR3A.271
N_SQAV = ( (PEXNER(I,K)-P_EXNER_CENTRE(I,KL))*N_SQ(I,KU) GWSCOR3A.272
& + (P_EXNER_CENTRE(I,KL) - PEXNER(I,K+1))*N_SQ(I,KL) ) GWSCOR3A.273
& / ( PEXNER(I,K) - PEXNER(I,K+1) ) GWSCOR3A.274
!-------------------------------------------------------------------- GWSCOR3A.275
! Note U is component parallel to stress vector GWSCOR3A.276
!-------------------------------------------------------------------- GWSCOR3A.277
UKP1 = U(I,K+1)*UNIT_X(I) + V(I,K+1)*UNIT_Y(I) GWSCOR3A.278
UK = U(I,K) *UNIT_X(I) + V(I,K) *UNIT_Y(I) GWSCOR3A.279
UKM1 = U(I,K-1)*UNIT_X(I) + V(I,K-1)*UNIT_Y(I) GWSCOR3A.280
IF( N_SQAV .LE. 0.0 .OR. UKP1 .LE. 0.3) THEN GWSCOR3A.281
L_LEE(I)=.FALSE. GWSCOR3A.282
ELSE GWSCOR3A.283
DEL_EXNER_KU = PEXNER(I,K+1) - P_EXNER_CENTRE(I,KU) GWSCOR3A.284
DEL_EXNER_KL = P_EXNER_CENTRE(I,KL) - PEXNER(I,K+1) GWSCOR3A.285
DELTA_Z= CPBYG*(THETA(I,K)*DEL_EXNER_KL + GWSCOR3A.286
* THETA(I,K+1)*DEL_EXNER_KU) GWSCOR3A.287
DU_DZ(I,KU)= ( UKP1-UK )/ DELTA_Z GWSCOR3A.288
DEL_EXNER = PEXNER(I,K) - PEXNER(I,K+1) GWSCOR3A.289
LSQ(I,K)=N_SQAV/UK**2 - ( DU_DZ(I,KU)-DU_DZ(I,KL) ) GWSCOR3A.290
* / (UK*CPBYG*THETA(I,K)*DEL_EXNER) GWSCOR3A.291
LSQ_SUM(I,2)=LSQ_SUM(I,2)+LSQ(I,K) GWSCOR3A.292
ENDIF GWSCOR3A.293
ENDIF ! ...K.NE.K_Lee_Lim(M-1).. GWSCOR3A.294
GWSCOR3A.295
IF( L_LEE(I) .AND. GWSCOR3A.296
& ( K_LEE(I,1).EQ.0 .OR. M.LE.K_TRANS ) ) THEN GWSCOR3A.297
GWSCOR3A.298
IF( K .EQ. K_LEE_LIM(M) ) THEN GWSCOR3A.299
GWSCOR3A.300
DEL_EXNER = PEXNER(I,M) - PEXNER(I,M+1) GWSCOR3A.301
ZH(I) = ZH(I) + CPBYG*THETA(I,M)*DEL_EXNER GWSCOR3A.302
GWSCOR3A.303
LSQ_SUM(I,1)=LSQ_SUM(I,1)+LSQ(I,M) GWSCOR3A.304
LSQ_SUM(I,2)=LSQ_SUM(I,2)-LSQ(I,M) GWSCOR3A.305
IF( LSQ_SUM(I,2).LT.0.0 ) LSQ_SUM(I,2)=0.0 GWSCOR3A.306
CALC= LSQ_SUM(I,1)/(M-1) - LSQ_SUM(I,2)/(K-M) GWSCOR3A.307
IF( CALC .GT. 0.0 .AND. K_LEE(I,1).EQ.0 GWSCOR3A.308
& .AND. M.GE.START_L ) THEN GWSCOR3A.309
CALC=CONST*ZH(I)*SQRT(CALC) GWSCOR3A.310
IF( CALC.LT.-1.0 ) THEN GWSCOR3A.311
K_LEE(I,1)=M GWSCOR3A.312
K_LEE(I,2)=K GWSCOR3A.313
H_LEE(I)=ZH(I) GWSCOR3A.314
LSQ_LEE(I,1)=LSQ_SUM(I,1)/(M-1) GWSCOR3A.315
LSQ_LEE(I,2)=LSQ_SUM(I,2)/(K-M) GWSCOR3A.316
ENDIF GWSCOR3A.317
ENDIF GWSCOR3A.318
ENDIF GWSCOR3A.319
GWSCOR3A.320
IF( M .EQ. K_TRANS ) THEN GWSCOR3A.321
IF( LSQ_SUM(I,1).GT.0.0 .AND. LSQ_SUM(I,2).GT.0.0 )THEN GWSCOR3A.322
LSQ1=LSQ_SUM(I,1)/(M-1) GWSCOR3A.323
LSQ2=LSQ_SUM(I,2)/(K-M) GWSCOR3A.324
L1=SQRT(LSQ1) GWSCOR3A.325
L2=SQRT(LSQ2) GWSCOR3A.326
CALC= LSQ1+LSQ2+(LSQ1-LSQ2)*COS(2.*L1*ZH(I)) GWSCOR3A.327
IF( CALC.EQ.0.0) THEN GWSCOR3A.328
TRANS(I)=1.0 GWSCOR3A.329
ELSE GWSCOR3A.330
TRANS(I)=2*L1*L2/ CALC GWSCOR3A.331
ENDIF GWSCOR3A.332
IF( TRANS(I).GT.1.0 ) THEN GWSCOR3A.333
TRANS(I)=1.0 GWSCOR3A.334
ENDIF GWSCOR3A.335
ENDIF GWSCOR3A.336
ENDIF GWSCOR3A.337
GWSCOR3A.338
IF( K_LEE(I,1).EQ.0 .AND. M.EQ.K_LEE_MAX ) THEN GWSCOR3A.339
L_LEE(I)= .FALSE. GWSCOR3A.340
ENDIF GWSCOR3A.341
GWSCOR3A.342
ENDIF ! K_Lee or K_trans = 0 GWSCOR3A.343
ENDDO ! Points GWSCOR3A.344
ENDDO ! K-Loop GWSCOR3A.345
ENDDO ! M-Loop GWSCOR3A.346
GWSCOR3A.347
! FIND CRITICAL LEVELS? GWSCOR3A.348
RETURN GWSCOR3A.349
END GWSCOR3A.350
GWSCOR3A.351
*ENDIF GWSCOR3A.352