[2fb99de] | 1 | MODULE FTRAFO |
---|
| 2 | |
---|
| 3 | CONTAINS |
---|
| 4 | |
---|
| 5 | |
---|
| 6 | |
---|
| 7 | C |
---|
| 8 | C Implementierung der spektralen Transformationsmethode unter Verwendung |
---|
| 9 | C des reduzierten Gauss'schen Gitters |
---|
| 10 | C |
---|
| 11 | C Berechnung der scale winds aus Vorticity und Divergenz |
---|
| 12 | C uebergibt man in XMN die Divergenz, so wird der divergente Anteil des |
---|
| 13 | C Windes (XPHI=Ud,XPHI=Vd) zurueckgegeben, uebergibt man die Vorticity, so |
---|
| 14 | C erhaelt man den rotationellen Wind (XLAM=Vrot,XPHI=-Urot). |
---|
| 15 | C Summiert man beide, erhaelt man den gesamten Scale wind |
---|
| 16 | C GWSAVE ist ein Hilfsfeld fuer die FFT |
---|
| 17 | C P enthaelt die assoziierten Legendrepolynome, H deren Ableitung |
---|
| 18 | C MLAT enthaelt die Anzahl der Gitterpunkte pro Breitenkreis |
---|
| 19 | C MNAUF gibt die spektrale Aufloesung an, |
---|
| 20 | C NI = Anzahl der Gauss'schen Gitterpunkte pro Flaeche |
---|
| 21 | C NJ = Anzahl der Gauss'schen Breiten, |
---|
| 22 | C NK = Anzahl der Niveaus |
---|
| 23 | |
---|
| 24 | SUBROUTINE VDTOUV(XMN,XLAM,XPHI,GWSAVE,IFAX,P,MLAT,MNAUF,NI,NJ,NK) |
---|
| 25 | |
---|
| 26 | |
---|
| 27 | USE PHTOGR |
---|
| 28 | |
---|
| 29 | IMPLICIT NONE |
---|
| 30 | INTEGER J,N,NI,NJ,NK,MNAUF,GGIND(NJ/2) |
---|
| 31 | INTEGER MLAT(NJ),IFAX(10,NJ) |
---|
| 32 | REAL XMN(0:(MNAUF+1)*(MNAUF+2)-1,NK) |
---|
| 33 | REAL P(0:(MNAUF+3)*(MNAUF+4)/2,NJ/2) |
---|
| 34 | REAL H(0:(MNAUF+2)*(MNAUF+3)/2) |
---|
| 35 | REAL XLAM(NI,NK),XPHI(NI,NK) |
---|
| 36 | REAL GWSAVE(8*NJ+15,NJ/2) |
---|
| 37 | REAL SCR,SCI,ACR,ACI,MUSCR,MUSCI,MUACR,MUACI |
---|
| 38 | REAL RT,IT |
---|
| 39 | |
---|
| 40 | GGIND(1)=0 |
---|
| 41 | DO 4 J = 2,NJ/2 |
---|
| 42 | GGIND(J)=GGIND(J-1)+MLAT(J-1) |
---|
| 43 | 4 CONTINUE |
---|
| 44 | !$OMP PARALLEL DO SCHEDULE(DYNAMIC) |
---|
| 45 | DO 5 J = 1,NJ/2 |
---|
| 46 | CALL VDUVSUB(J,XMN,XLAM,XPHI,GWSAVE,IFAX,P,GGIND(J), |
---|
| 47 | *MLAT,MNAUF,NI,NJ,NK) |
---|
| 48 | 5 CONTINUE |
---|
| 49 | !$OMP END PARALLEL DO |
---|
| 50 | RETURN |
---|
| 51 | END SUBROUTINE VDTOUV |
---|
| 52 | |
---|
| 53 | SUBROUTINE VDUVSUB(J,XMN,XLAM,XPHI,GWSAVE,IFAX,P, |
---|
| 54 | *GGIND,MLAT,MNAUF,NI,NJ,NK) |
---|
| 55 | |
---|
| 56 | USE PHTOGR |
---|
| 57 | |
---|
| 58 | IMPLICIT NONE |
---|
| 59 | INTEGER J,K,M,N,NI,NJ,NK,MNAUF,GGIND,LL,LLP,LLH,LLS,LLPS,LLHS |
---|
| 60 | INTEGER MLAT(NJ),IFAX(10,NJ) |
---|
| 61 | REAL UFOUC(0:MAXAUF),MUFOUC(0:MAXAUF) |
---|
| 62 | REAL VFOUC(0:MAXAUF),MVFOUC(0:MAXAUF) |
---|
| 63 | REAL XMN(0:(MNAUF+1)*(MNAUF+2)-1,NK) |
---|
| 64 | REAL P(0:(MNAUF+3)*(MNAUF+4)/2,NJ/2) |
---|
| 65 | REAL H(0:(MNAUF+2)*(MNAUF+3)/2) |
---|
| 66 | REAL XLAM(NI,NK),XPHI(NI,NK) |
---|
| 67 | REAL GWSAVE(8*NJ+15,NJ/2) |
---|
| 68 | REAL ERAD,SCR,SCI,ACR,ACI,MUSCR,MUSCI,MUACR,MUACI |
---|
| 69 | REAL FAC(0:MNAUF),RT,IT |
---|
| 70 | |
---|
| 71 | |
---|
| 72 | ERAD = 6367470.D0 |
---|
| 73 | |
---|
| 74 | FAC(0)=0.D0 |
---|
| 75 | DO 12 N=1,MNAUF |
---|
| 76 | FAC(N)=-ERAD/DBLE(N)/DBLE(N+1) |
---|
| 77 | 12 CONTINUE |
---|
| 78 | |
---|
| 79 | CALL DPLGND(MNAUF,P(0,J),H) |
---|
| 80 | DO 3 K = 1,NK |
---|
| 81 | LL=0 |
---|
| 82 | LLP=0 |
---|
| 83 | LLH=0 |
---|
| 84 | DO 2 M = 0,MNAUF |
---|
| 85 | SCR=0.D0 |
---|
| 86 | SCI=0.D0 |
---|
| 87 | ACR=0.D0 |
---|
| 88 | ACI=0.D0 |
---|
| 89 | MUSCR=0.D0 |
---|
| 90 | MUSCI=0.D0 |
---|
| 91 | MUACR=0.D0 |
---|
| 92 | MUACI=0.D0 |
---|
| 93 | LLS=LL |
---|
| 94 | LLPS=LLP |
---|
| 95 | LLHS=LLH |
---|
| 96 | IF(2*M+1.LT.MLAT(J)) THEN |
---|
| 97 | DO 1 N = M,MNAUF,2 |
---|
| 98 | RT=XMN(2*LL,K)*FAC(N) |
---|
| 99 | IT=XMN(2*LL+1,K)*FAC(N) |
---|
| 100 | SCR =SCR+ RT*P(LLP,J) |
---|
| 101 | SCI =SCI+ IT*P(LLP,J) |
---|
| 102 | MUACR =MUACR+ RT*H(LLH) |
---|
| 103 | MUACI =MUACI+ IT*H(LLH) |
---|
| 104 | LL=LL+2 |
---|
| 105 | LLP=LLP+2 |
---|
| 106 | LLH=LLH+2 |
---|
| 107 | 1 CONTINUE |
---|
| 108 | LL=LLS+1 |
---|
| 109 | LLP=LLPS+1 |
---|
| 110 | LLH=LLHS+1 |
---|
| 111 | DO 11 N = M+1,MNAUF,2 |
---|
| 112 | RT=XMN(2*LL,K)*FAC(N) |
---|
| 113 | IT=XMN(2*LL+1,K)*FAC(N) |
---|
| 114 | ACR =ACR+ RT*P(LLP,J) |
---|
| 115 | ACI =ACI+ IT*P(LLP,J) |
---|
| 116 | MUSCR =MUSCR+ RT*H(LLH) |
---|
| 117 | MUSCI =MUSCI+ IT*H(LLH) |
---|
| 118 | LL=LL+2 |
---|
| 119 | LLP=LLP+2 |
---|
| 120 | LLH=LLH+2 |
---|
| 121 | 11 CONTINUE |
---|
| 122 | ENDIF |
---|
| 123 | LL=LLS+(MNAUF-M+1) |
---|
| 124 | LLP=LLPS+(MNAUF-M+3) |
---|
| 125 | LLH=LLHS+(MNAUF-M+2) |
---|
| 126 | |
---|
| 127 | UFOUC(2*M)=-M*(SCI-ACI) |
---|
| 128 | UFOUC(2*M+1)=M*(SCR-ACR) |
---|
| 129 | VFOUC(2*M)=-M*(SCI+ACI) |
---|
| 130 | VFOUC(2*M+1)=M*(SCR+ACR) |
---|
| 131 | |
---|
| 132 | MUFOUC(2*M)=-(MUSCR-MUACR) |
---|
| 133 | MUFOUC(2*M+1)=-(MUSCI-MUACI) |
---|
| 134 | MVFOUC(2*M)=-(MUSCR+MUACR) |
---|
| 135 | MVFOUC(2*M+1)=-(MUSCI+MUACI) |
---|
| 136 | 2 CONTINUE |
---|
| 137 | |
---|
| 138 | CALL RFOURTR(VFOUC, |
---|
| 139 | *GWSAVE(:,J),IFAX(:,J),MNAUF,MLAT(J),1) |
---|
| 140 | XLAM(GGIND+1:GGIND+MLAT(J),K)=VFOUC(0:MLAT(J)-1) |
---|
| 141 | CALL RFOURTR(UFOUC, |
---|
| 142 | *GWSAVE(:,J),IFAX(:,J),MNAUF,MLAT(J),1) |
---|
| 143 | XLAM(NI-GGIND-MLAT(J)+1:NI-GGIND,K)=UFOUC(0:MLAT(J)-1) |
---|
| 144 | |
---|
| 145 | CALL RFOURTR(MVFOUC, |
---|
| 146 | *GWSAVE(:,J),IFAX(:,J),MNAUF,MLAT(J),1) |
---|
| 147 | XPHI(GGIND+1:GGIND+MLAT(J),K)=MVFOUC(0:MLAT(J)-1) |
---|
| 148 | CALL RFOURTR(MUFOUC, |
---|
| 149 | *GWSAVE(:,J),IFAX(:,J),MNAUF,MLAT(J),1) |
---|
| 150 | XPHI(NI-GGIND-MLAT(J)+1:NI-GGIND,K)=MUFOUC(0:MLAT(J)-1) |
---|
| 151 | |
---|
| 152 | 3 CONTINUE |
---|
| 153 | |
---|
| 154 | RETURN |
---|
| 155 | END SUBROUTINE VDUVSUB |
---|
| 156 | |
---|
| 157 | C Berechnung des Gradienten eines Skalars aus dem Feld des |
---|
| 158 | C Skalars XMN im Phasenraum. Zurueckgegeben werden die Felder der |
---|
| 159 | C Komponenten des horizontalen Gradienten XLAM,XPHI auf dem Gauss'schen Gitter. |
---|
| 160 | C GWSAVE ist ein Hilfsfeld fuer die FFT |
---|
| 161 | C P enthaelt die assoziierten Legendrepolynome, H deren Ableitung |
---|
| 162 | C MLAT enthaelt die Anzahl der Gitterpunkte pro Breitenkreis |
---|
| 163 | C MNAUF gibt die spektrale Aufloesung an, |
---|
| 164 | C NI = Anzahl der Gauss'schen Gitterpunkte, |
---|
| 165 | C NJ = Anzahl der Gauss'schen Breiten, |
---|
| 166 | C NK = Anzahl der Niveaus |
---|
| 167 | |
---|
| 168 | SUBROUTINE PHGRAD(XMN,XLAM,XPHI,GWSAVE,IFAX,P,H,MLAT, |
---|
| 169 | *MNAUF,NI,NJ,NK) |
---|
| 170 | |
---|
| 171 | USE PHTOGR |
---|
| 172 | IMPLICIT NONE |
---|
| 173 | INTEGER J,K,M,N,NI,NJ,NK,MNAUF,GGIND,LL,LLP,LLH,LLS,LLPS,LLHS |
---|
| 174 | INTEGER MLAT(NJ),IFAX(10,NJ) |
---|
| 175 | REAL UFOUC(0:MAXAUF),MUFOUC(0:MAXAUF) |
---|
| 176 | REAL VFOUC(0:MAXAUF),MVFOUC(0:MAXAUF) |
---|
| 177 | REAL XMN(0:(MNAUF+1)*(MNAUF+2)-1,NK) |
---|
| 178 | REAL P(0:(MNAUF+3)*(MNAUF+4)/2,NJ/2) |
---|
| 179 | REAL H(0:(MNAUF+2)*(MNAUF+3)/2) |
---|
| 180 | REAL XLAM(NI,NK),XPHI(NI,NK) |
---|
| 181 | REAL GWSAVE(8*NJ+15,NJ/2) |
---|
| 182 | REAL ERAD |
---|
| 183 | REAL SCR,SCI,ACR,ACI,MUSCR,MUSCI,MUACR,MUACI,RT,IT |
---|
| 184 | |
---|
| 185 | ERAD = 6367470.0 |
---|
| 186 | |
---|
| 187 | GGIND=0 |
---|
| 188 | DO 4 J = 1,NJ/2 |
---|
| 189 | CALL DPLGND(MNAUF,P(0,J),H) |
---|
| 190 | DO 3 K = 1,NK |
---|
| 191 | LL=0 |
---|
| 192 | LLP=0 |
---|
| 193 | LLH=0 |
---|
| 194 | DO 2 M = 0,MNAUF |
---|
| 195 | SCR=0.D0 |
---|
| 196 | SCI=0.D0 |
---|
| 197 | ACR=0.D0 |
---|
| 198 | ACI=0.D0 |
---|
| 199 | MUSCR=0.D0 |
---|
| 200 | MUSCI=0.D0 |
---|
| 201 | MUACR=0.D0 |
---|
| 202 | MUACI=0.D0 |
---|
| 203 | LLS=LL |
---|
| 204 | LLPS=LLP |
---|
| 205 | LLHS=LLH |
---|
| 206 | IF(2*M+1.LT.MLAT(J)) THEN |
---|
| 207 | DO 1 N = M,MNAUF,2 |
---|
| 208 | RT=XMN(2*LL,K) |
---|
| 209 | IT=XMN(2*LL+1,K) |
---|
| 210 | SCR =SCR+ RT*P(LLP,J) |
---|
| 211 | SCI =SCI+ IT*P(LLP,J) |
---|
| 212 | MUACR =MUACR+RT*H(LLH) |
---|
| 213 | MUACI =MUACI+ IT*H(LLH) |
---|
| 214 | LL=LL+2 |
---|
| 215 | LLP=LLP+2 |
---|
| 216 | LLH=LLH+2 |
---|
| 217 | 1 CONTINUE |
---|
| 218 | LL=LLS+1 |
---|
| 219 | LLP=LLPS+1 |
---|
| 220 | LLH=LLHS+1 |
---|
| 221 | DO 11 N = M+1,MNAUF,2 |
---|
| 222 | RT=XMN(2*LL,K) |
---|
| 223 | IT=XMN(2*LL+1,K) |
---|
| 224 | ACR =ACR+ RT*P(LLP,J) |
---|
| 225 | ACI =ACI+ IT*P(LLP,J) |
---|
| 226 | MUSCR =MUSCR+ RT*H(LLH) |
---|
| 227 | MUSCI =MUSCI+ IT*H(LLH) |
---|
| 228 | LL=LL+2 |
---|
| 229 | LLP=LLP+2 |
---|
| 230 | LLH=LLH+2 |
---|
| 231 | 11 CONTINUE |
---|
| 232 | ENDIF |
---|
| 233 | LL=LLS+(MNAUF-M+1) |
---|
| 234 | LLP=LLPS+(MNAUF-M+3) |
---|
| 235 | LLH=LLHS+(MNAUF-M+2) |
---|
| 236 | |
---|
| 237 | UFOUC(2*M)=-M*(SCI-ACI)/ERAD |
---|
| 238 | UFOUC(2*M+1)=M*(SCR-ACR)/ERAD |
---|
| 239 | VFOUC(2*M)=-M*(SCI+ACI)/ERAD |
---|
| 240 | VFOUC(2*M+1)=M*(SCR+ACR)/ERAD |
---|
| 241 | |
---|
| 242 | MUFOUC(2*M)=-(MUSCR-MUACR)/ERAD |
---|
| 243 | MUFOUC(2*M+1)=-(MUSCI-MUACI)/ERAD |
---|
| 244 | MVFOUC(2*M)=-(MUSCR+MUACR)/ERAD |
---|
| 245 | MVFOUC(2*M+1)=-(MUSCI+MUACI)/ERAD |
---|
| 246 | 2 CONTINUE |
---|
| 247 | |
---|
| 248 | CALL RFOURTR(VFOUC, |
---|
| 249 | *GWSAVE(:,J),IFAX(:,J),MNAUF,MLAT(J),1) |
---|
| 250 | XLAM(GGIND+1:GGIND+MLAT(J),K)=VFOUC(0:MLAT(J)-1) |
---|
| 251 | CALL RFOURTR(UFOUC, |
---|
| 252 | *GWSAVE(:,J),IFAX(:,J),MNAUF,MLAT(J),1) |
---|
| 253 | XLAM(NI-GGIND-MLAT(J)+1:NI-GGIND,K)=UFOUC(0:MLAT(J)-1) |
---|
| 254 | |
---|
| 255 | CALL RFOURTR(MVFOUC, |
---|
| 256 | *GWSAVE(:,J),IFAX(:,J),MNAUF,MLAT(J),1) |
---|
| 257 | XPHI(GGIND+1:GGIND+MLAT(J),K)=MVFOUC(0:MLAT(J)-1) |
---|
| 258 | CALL RFOURTR(MUFOUC, |
---|
| 259 | *GWSAVE(:,J),IFAX(:,J),MNAUF,MLAT(J),1) |
---|
| 260 | XPHI(NI-GGIND-MLAT(J)+1:NI-GGIND,K)=MUFOUC(0:MLAT(J)-1) |
---|
| 261 | |
---|
| 262 | 3 CONTINUE |
---|
| 263 | GGIND=GGIND+MLAT(J) |
---|
| 264 | 4 CONTINUE |
---|
| 265 | |
---|
| 266 | |
---|
| 267 | RETURN |
---|
| 268 | END SUBROUTINE PHGRAD |
---|
| 269 | |
---|
| 270 | C Berechnung des Gradienten eines Skalars aus dem Feld des |
---|
| 271 | C Skalars XMN im Phasenraum. Zurueckgegeben werden die Felder der |
---|
| 272 | C Komponenten des horizontalen Gradienten XLAM,XPHI auf dem Gauss'schen Gitter. |
---|
| 273 | C GWSAVE ist ein Hilfsfeld fuer die FFT |
---|
| 274 | C P enthaelt die assoziierten Legendrepolynome, H deren Ableitung |
---|
| 275 | C MLAT enthaelt die Anzahl der Gitterpunkte pro Breitenkreis |
---|
| 276 | C MNAUF gibt die spektrale Aufloesung an, |
---|
| 277 | C NI = Anzahl der Gauss'schen Gitterpunkte, |
---|
| 278 | C NJ = Anzahl der Gauss'schen Breiten, |
---|
| 279 | C NK = Anzahl der Niveaus |
---|
| 280 | |
---|
| 281 | SUBROUTINE PHGRACUT(XMN,XLAM,XPHI,GWSAVE,IFAX,P,H,MAUF, |
---|
| 282 | *MNAUF,NI,NJ,MANF,NK) |
---|
| 283 | |
---|
| 284 | USE PHTOGR |
---|
| 285 | IMPLICIT NONE |
---|
| 286 | INTEGER J,K,M,N,NI,NJ,NK,MNAUF,GGIND,LL,LLP,LLH,LLS,LLPS,LLHS |
---|
| 287 | INTEGER MAUF,MANF,I,IFAX(10) |
---|
| 288 | REAL UFOUC(0:MAXAUF),MUFOUC(0:MAXAUF) |
---|
| 289 | REAL VFOUC(0:MAXAUF),MVFOUC(0:MAXAUF) |
---|
| 290 | REAL XMN(0:(MNAUF+1)*(MNAUF+2)-1,NK) |
---|
| 291 | REAL P(0:(MNAUF+3)*(MNAUF+4)/2,NJ) |
---|
| 292 | REAL H(0:(MNAUF+2)*(MNAUF+3)/2) |
---|
| 293 | REAL XLAM(NI,NJ,NK),XPHI(NI,NJ,NK) |
---|
| 294 | REAL HLAM(MAXAUF,2),HPHI(MAXAUF,2) |
---|
| 295 | REAL GWSAVE(4*MAUF+15) |
---|
| 296 | REAL ERAD |
---|
| 297 | REAL SCR,SCI,ACR,ACI,MUSCR,MUSCI,MUACR,MUACI,RT,IT |
---|
| 298 | |
---|
| 299 | ERAD = 6367470.0 |
---|
| 300 | |
---|
| 301 | GGIND=0 |
---|
| 302 | DO 4 J = 1,NJ |
---|
| 303 | CALL DPLGND(MNAUF,P(0,J),H) |
---|
| 304 | DO 3 K = 1,NK |
---|
| 305 | LL=0 |
---|
| 306 | LLP=0 |
---|
| 307 | LLH=0 |
---|
| 308 | DO 2 M = 0,MNAUF |
---|
| 309 | SCR=0.D0 |
---|
| 310 | SCI=0.D0 |
---|
| 311 | ACR=0.D0 |
---|
| 312 | ACI=0.D0 |
---|
| 313 | MUSCR=0.D0 |
---|
| 314 | MUSCI=0.D0 |
---|
| 315 | MUACR=0.D0 |
---|
| 316 | MUACI=0.D0 |
---|
| 317 | LLS=LL |
---|
| 318 | LLPS=LLP |
---|
| 319 | LLHS=LLH |
---|
| 320 | IF(2*M+1.LT.MAUF) THEN |
---|
| 321 | DO 1 N = M,MNAUF,2 |
---|
| 322 | RT=XMN(2*LL,K) |
---|
| 323 | IT=XMN(2*LL+1,K) |
---|
| 324 | SCR =SCR+ RT*P(LLP,J) |
---|
| 325 | SCI =SCI+ IT*P(LLP,J) |
---|
| 326 | MUACR =MUACR+RT*H(LLH) |
---|
| 327 | MUACI =MUACI+ IT*H(LLH) |
---|
| 328 | LL=LL+2 |
---|
| 329 | LLP=LLP+2 |
---|
| 330 | LLH=LLH+2 |
---|
| 331 | 1 CONTINUE |
---|
| 332 | LL=LLS+1 |
---|
| 333 | LLP=LLPS+1 |
---|
| 334 | LLH=LLHS+1 |
---|
| 335 | DO 11 N = M+1,MNAUF,2 |
---|
| 336 | RT=XMN(2*LL,K) |
---|
| 337 | IT=XMN(2*LL+1,K) |
---|
| 338 | ACR =ACR+ RT*P(LLP,J) |
---|
| 339 | ACI =ACI+ IT*P(LLP,J) |
---|
| 340 | MUSCR =MUSCR+ RT*H(LLH) |
---|
| 341 | MUSCI =MUSCI+ IT*H(LLH) |
---|
| 342 | LL=LL+2 |
---|
| 343 | LLP=LLP+2 |
---|
| 344 | LLH=LLH+2 |
---|
| 345 | 11 CONTINUE |
---|
| 346 | ENDIF |
---|
| 347 | LL=LLS+(MNAUF-M+1) |
---|
| 348 | LLP=LLPS+(MNAUF-M+3) |
---|
| 349 | LLH=LLHS+(MNAUF-M+2) |
---|
| 350 | |
---|
| 351 | UFOUC(2*M)=-M*(SCI-ACI)/ERAD |
---|
| 352 | UFOUC(2*M+1)=M*(SCR-ACR)/ERAD |
---|
| 353 | VFOUC(2*M)=-M*(SCI+ACI)/ERAD |
---|
| 354 | VFOUC(2*M+1)=M*(SCR+ACR)/ERAD |
---|
| 355 | |
---|
| 356 | MUFOUC(2*M)=-(MUSCR-MUACR)/ERAD |
---|
| 357 | MUFOUC(2*M+1)=-(MUSCI-MUACI)/ERAD |
---|
| 358 | MVFOUC(2*M)=-(MUSCR+MUACR)/ERAD |
---|
| 359 | MVFOUC(2*M+1)=-(MUSCI+MUACI)/ERAD |
---|
| 360 | 2 CONTINUE |
---|
| 361 | |
---|
| 362 | CALL RFOURTR(VFOUC, |
---|
| 363 | *GWSAVE,IFAX,MNAUF,MAUF,1) |
---|
| 364 | |
---|
| 365 | CALL RFOURTR(MVFOUC, |
---|
| 366 | *GWSAVE,IFAX,MNAUF,MAUF,1) |
---|
| 367 | |
---|
| 368 | DO 6 I=0,NI-1 |
---|
| 369 | IF(MANF+I.LE. MAUF) THEN |
---|
| 370 | XLAM(I+1,J,K)=VFOUC(MANF+I-1) |
---|
| 371 | XPHI(I+1,J,K)=MVFOUC(MANF+I-1) |
---|
| 372 | ELSE |
---|
| 373 | XLAM(I+1,J,K)=VFOUC(MANF-MAUF+I-1) |
---|
| 374 | XPHI(I+1,J,K)=MVFOUC(MANF-MAUF+I-1) |
---|
| 375 | ENDIF |
---|
| 376 | 6 CONTINUE |
---|
| 377 | 3 CONTINUE |
---|
| 378 | GGIND=GGIND+MAUF |
---|
| 379 | 4 CONTINUE |
---|
| 380 | |
---|
| 381 | RETURN |
---|
| 382 | END SUBROUTINE PHGRACUT |
---|
| 383 | |
---|
| 384 | C Berechnung der Divergenz aus dem Windfeld (U,V) |
---|
| 385 | C im Phasenraum. Zurueckgegeben werden die Felder der |
---|
| 386 | C Komponenten des horizontalen Gradienten XLAM,XPHI auf dem Gauss'schen Gitter. |
---|
| 387 | C GWSAVE ist ein Hilfsfeld fuer die FFT |
---|
| 388 | C P enthaelt die assoziierten Legendrepolynome, H deren Ableitung |
---|
| 389 | C MLAT enthaelt die Anzahl der Gitterpunkte pro Breitenkreis |
---|
| 390 | C MNAUF gibt die spektrale Aufloesung an, |
---|
| 391 | C NI = Anzahl der Gauss'schen Gitterpunkte, |
---|
| 392 | C NJ = Anzahl der Gauss'schen Breiten, |
---|
| 393 | C NK = Anzahl der Niveaus |
---|
| 394 | C Beachte, dass das Windfeld eine um 1 erhoehte Aufloesung in mu-Richtung hat. |
---|
| 395 | |
---|
| 396 | SUBROUTINE CONTGL(PS,DPSDL,DPSDM,DIV,U,V,BREITE,ETA, |
---|
| 397 | *MLAT,A,B,NI,NJ,NK) |
---|
| 398 | |
---|
| 399 | IMPLICIT NONE |
---|
| 400 | |
---|
| 401 | INTEGER NI,NJ,NK,I,J,K,MLAT(NJ),L |
---|
| 402 | |
---|
| 403 | REAL A(NK+1),B(NK+1) |
---|
| 404 | REAL PS(NI),DPSDL(NI),DPSDM(NI) |
---|
| 405 | REAL DIV(NI,NK),U(NI,NK),V(NI,NK),ETA(NI,NK) |
---|
| 406 | REAL BREITE(NJ) |
---|
| 407 | |
---|
| 408 | REAL DIVT1,DIVT2,POB,PUN,DPSDT,COSB |
---|
| 409 | |
---|
| 410 | L=0 |
---|
| 411 | DO 4 J=1,NJ |
---|
| 412 | COSB=(1.0-BREITE(J)*BREITE(J)) |
---|
| 413 | DO 3 I=1,MLAT(J) |
---|
| 414 | L=L+1 |
---|
| 415 | DIVT1=0.0 |
---|
| 416 | DIVT2=0.0 |
---|
| 417 | DO 1 K=1,NK |
---|
| 418 | POB=A(K)+B(K)*PS(L) |
---|
| 419 | PUN=A(K+1)+B(K+1)*PS(L) |
---|
| 420 | |
---|
| 421 | DIVT1=DIVT1+DIV(L,K)*(PUN-POB) |
---|
| 422 | if(cosb .gt. 0.) then |
---|
| 423 | DIVT2=DIVT2+(B(K+1)-B(K))*PS(L)* |
---|
| 424 | *(U(L,K)*DPSDL(L)+V(L,K)*DPSDM(L))/COSB |
---|
| 425 | endif |
---|
| 426 | |
---|
| 427 | ETA(L,K)=-DIVT1-DIVT2 |
---|
| 428 | 1 CONTINUE |
---|
| 429 | |
---|
| 430 | DPSDT=(-DIVT1-DIVT2)/PS(L) |
---|
| 431 | |
---|
| 432 | DO 2 K=1,NK |
---|
| 433 | ETA(L,K)=ETA(L,K)-DPSDT*B(K+1)*PS(L) |
---|
| 434 | 2 CONTINUE |
---|
| 435 | PS(L)=DPSDT*PS(L) |
---|
| 436 | 3 CONTINUE |
---|
| 437 | 4 CONTINUE |
---|
| 438 | RETURN |
---|
| 439 | END SUBROUTINE CONTGL |
---|
| 440 | |
---|
| 441 | C OMEGA berechnet omega im Hybridkoordinatensystem |
---|
| 442 | C PS ist der Bodendruck, |
---|
| 443 | C DPSDL,DPSDM sind die Komponenten des Gradienten des Logarithmus des |
---|
| 444 | C Bodendrucks |
---|
| 445 | C DIV,U,V sind die horizontale Divergenz und das horizontale Windfeld |
---|
| 446 | C BREITE ist das Feld der Gauss'schen Breiten |
---|
| 447 | C E ist omega, |
---|
| 448 | |
---|
| 449 | SUBROUTINE OMEGA(PS,DPSDL,DPSDM,DIV,U,V,BREITE,E,MLAT,A,B,NGI |
---|
| 450 | * ,NGJ,MKK) |
---|
| 451 | |
---|
| 452 | IMPLICIT NONE |
---|
| 453 | |
---|
| 454 | INTEGER I,J,K,L,NGI,NGJ,MKK,MLAT(NGJ) |
---|
| 455 | |
---|
| 456 | REAL PS(NGI),DPSDL(NGI),DPSDM(NGI),A(MKK+1),B(MKK+1) |
---|
| 457 | REAL DIV(NGI,MKK),U(NGI,MKK),V(NGI,MKK),E(NGI,MKK) |
---|
| 458 | REAL BREITE(NGJ) |
---|
| 459 | |
---|
| 460 | REAL DIVT1,DIVT2,POB,PUN,DP,X,Y,COSB |
---|
| 461 | REAL DIVT3(MKK+2) |
---|
| 462 | |
---|
| 463 | L=0 |
---|
| 464 | DO 4 J=1,NGJ |
---|
| 465 | COSB=(1.0-BREITE(J)*BREITE(J)) |
---|
| 466 | DO 3 I=1,MLAT(J) |
---|
| 467 | L=L+1 |
---|
| 468 | DIVT1=0.0 |
---|
| 469 | DIVT2=0.0 |
---|
| 470 | DIVT3(1)=0.0 |
---|
| 471 | DO 1 K=1,MKK |
---|
| 472 | POB=A(K)+B(K)*PS(L) |
---|
| 473 | PUN=A(K+1)+B(K+1)*PS(L) |
---|
| 474 | DP=PUN-POB |
---|
| 475 | |
---|
| 476 | Y=PS(L)*(U(L,K)*DPSDL(L)+V(L,K)*DPSDM(L))/COSB |
---|
| 477 | IF(K.LT.3) THEN |
---|
| 478 | X=0.0 |
---|
| 479 | ELSE |
---|
| 480 | X=(B(K+1)-B(K))*Y |
---|
| 481 | ENDIF |
---|
| 482 | |
---|
| 483 | DIVT1=DIVT1+DIV(L,K)*DP |
---|
| 484 | DIVT2=DIVT2+X |
---|
| 485 | |
---|
| 486 | DIVT3(K+1)=-DIVT1-DIVT2 |
---|
| 487 | |
---|
| 488 | IF(K.GT.1) THEN |
---|
| 489 | E(L,K) = 0.5*(POB+PUN)/DP*Y* |
---|
| 490 | *((B(K+1)-B(K))+(A(K+1)*B(K)-A(K)*B(K+1))/ |
---|
| 491 | *DP*LOG(PUN/POB)) |
---|
| 492 | ELSE |
---|
| 493 | E(L,K) = 0.0 |
---|
| 494 | ENDIF |
---|
| 495 | |
---|
| 496 | E(L,K) = E(L,K)+0.5*(DIVT3(K)+DIVT3(K+1)) |
---|
| 497 | |
---|
| 498 | 1 CONTINUE |
---|
| 499 | 3 CONTINUE |
---|
| 500 | 4 CONTINUE |
---|
| 501 | RETURN |
---|
| 502 | END SUBROUTINE OMEGA |
---|
| 503 | |
---|
| 504 | END MODULE FTRAFO |
---|