source: flexpart.git/preproc/src/grphreal.f90 @ 5763793

10.4.1_peseiFPv9.3.1FPv9.3.1b_testingFPv9.3.2GFS_025bugfixes+enhancementsdevfp9.3.1-20161214-nc4grib2nc4_repairrelease-10release-10.4.1scaling-bugunivie
Last change on this file since 5763793 was 5763793, checked in by Anne Fouilloux <annefou@…>, 9 years ago

Add pre-processing programs and examples. This is still under development!

  • Property mode set to 100644
File size: 5.1 KB
Line 
1MODULE GRTOPH
2
3  USE PHTOGR
4
5CONTAINS
6  !
7  SUBROUTINE GRPH213(CXMN,FELD,WSAVE,IFAX,Z,W,MLAT, &
8       MNAUF,MAXL,MAXB,MLEVEL)
9
10    !     DIE ROUTINE F]HRT EINE TRANSFORMATION EINER
11    !     FELDVARIABLEN VOM PHASENRAUM IN  DEN PHYSIKALISCHEN
12    !     RAUM AUF KUGELKOORDINATEN DURCH
13    !
14    !     CXMN   = SPEKTRALKOEFFIZIENTEN IN DER REIHENFOLGE
15    !              CX00,CX01,CX11,CX02,....CXMNAUFMNAUF
16    !                   CXM              = FOURIERKOEFFIZIENTEN - nur ein Hilfsfeld
17    !     FELD   = FELD DER METEOROLOGISCHEN VARIABLEN
18    !                   WSAVE  = Working Array fuer Fouriertransformation
19    !     Z      = LEGENDREFUNKTIONSWERTE
20    !
21    !     MNAUF    ANZAHL DER FOURIERKOEFFIZIENTEN
22    !     MAXL     ANZAHL DER FUER DAS GITTER BENUTZTEN LAENGEN
23    !     MAXB     ANZAHL DER FUER DAS GITTER BENOETIGTEN BREITEN
24    !     MLEVEL   ANZAHL DER LEVELS, DIE TRANSFORMIERT WERDEN
25    !
26    IMPLICIT REAL (A-H,O-Z)
27
28
29    !                   Anzahl der Gitterpunkte pro Breitenkreis des reduzierten
30    !                   Gauss'schen Gitters
31    INTEGER MLAT(MAXB),ISIZE,IFAX(10,MAXB)
32
33    !    FELD DER LEGENDREPOLYNOME FUER EINE BREITE
34    REAL*8 Z(MAXB/2,0:((MNAUF+3)*(MNAUF+4))/2)
35
36    !      LOGICAL*1 USED(((216*217)/2+1)*160)
37
38    DIMENSION CXMN(0:(MNAUF+1)*(MNAUF+2)-1,MLEVEL)
39    REAL FELD(MAXL,MLEVEL)
40    DIMENSION WSAVE(8*MAXB+15,MAXB/2)
41    REAL*8 W(MAXB)
42    DIMENSION IND(MAXB)
43
44
45    IND(1)=0
46    DO J=2,MAXB/2
47       IND(j)=IND(J-1)+MLAT(J-1)
48    ENDDO
49    !$OMP PARALLEL DO SCHEDULE(DYNAMIC)
50    DO  L=1,MLEVEL
51       CALL GRPHSUB(L,IND,CXMN,FELD,WSAVE,IFAX,Z,W,MLAT, &
52            &MNAUF,MAXL,MAXB,MLEVEL)
53    ENDDO
54    !$omp end parallel do
55
56
57    RETURN
58  END SUBROUTINE GRPH213
59  !
60  SUBROUTINE GRPHSUB(L,IND,CXMN,FELD,WSAVE,IFAX,Z,W,MLAT, &
61       &MNAUF,MAXL,MAXB,MLEVEL)
62
63    !     DIE ROUTINE F]HRT EINE TRANSFORMATION EINER
64    !     FELDVARIABLEN VOM PHASENRAUM IN  DEN PHYSIKALISCHEN
65    !     RAUM AUF KUGELKOORDINATEN DURCH
66    !
67    !     CXMN   = SPEKTRALKOEFFIZIENTEN IN DER REIHENFOLGE
68    !              CX00,CX01,CX11,CX02,....CXMNAUFMNAUF
69    !                   CXM              = FOURIERKOEFFIZIENTEN - nur ein Hilfsfeld
70    !     FELD   = FELD DER METEOROLOGISCHEN VARIABLEN
71    !                   WSAVE  = Working Array fuer Fouriertransformation
72    !     Z      = LEGENDREFUNKTIONSWERTE
73    !
74    !     MNAUF    ANZAHL DER FOURIERKOEFFIZIENTEN
75    !     MAXL     ANZAHL DER FUER DAS GITTER BENUTZTEN LAENGEN
76    !     MAXB     ANZAHL DER FUER DAS GITTER BENOETIGTEN BREITEN
77    !     MLEVEL   ANZAHL DER LEVELS, DIE TRANSFORMIERT WERDEN
78    !
79    IMPLICIT REAL (A-H,O-Z)
80
81    !    FELD DER FOURIERKOEFFIZIENTEN
82    REAL CXMS(4*(MNAUF+1))
83    REAL CXMA(4*(MNAUF+1))
84    REAL,ALLOCATABLE :: CXM(:,:)
85
86    !                   Anzahl der Gitterpunkte pro Breitenkreis des reduzierten
87    !                   Gauss'schen Gitters
88    INTEGER MLAT(MAXB),ISIZE
89
90    !    FELD DER LEGENDREPOLYNOME FUER EINE BREITE
91    REAL Z(MAXB/2,0:((MNAUF+3)*(MNAUF+4))/2)
92
93    !      LOGICAL*1 USED(((216*217)/2+1)*160)
94
95    REAL CXMN(0:(MNAUF+1)*(MNAUF+2)-1,MLEVEL)
96    REAL FELD(MAXL,MLEVEL)
97    REAL WSAVE(8*MAXB+15,MAXB/2)
98    INTEGER IFAX(10,MAXB)
99    REAL W(MAXB)
100    INTEGER IND(MAXB)
101
102    ALLOCATE(CXM( 4*MAXB,MAXB))
103    DO  J=1,MAXB/2
104       CXMS(1:MLAT(J))=FELD(IND(J)+1:IND(J)+MLAT(J),L)
105       CALL RFOUFTR(CXMS,WSAVE(1,J),IFAX(:,J),MNAUF,MLAT(J),1)
106       CXMA(1:MLAT(J))=FELD(MAXL-IND(J)-MLAT(J)+1:MAXL-IND(J),L)
107       CALL RFOUFTR(CXMA, &
108            &WSAVE(1,J),IFAX(:,J),MNAUF,MLAT(J),1)
109       DO  I=1,2*(MNAUF+1)
110          CXM(I,J)=CXMS(I)+CXMA(I)
111          CXM(I,MAXB+1-J)=CXMS(I)-CXMA(I)
112       ENDDO
113    ENDDO
114    CALL LGTR213(CXMN(0,L),CXM,Z,W,MLAT,MNAUF,MAXB)
115
116    DEALLOCATE(CXM)
117
118    RETURN
119  END SUBROUTINE GRPHSUB
120  !
121  SUBROUTINE LGTR213(CXMN,CXM,Z,W,MLAT,MNAUF,MAXB)
122    IMPLICIT REAL (A-H,O-Z)
123    INTEGER MLAT(MAXB)
124    DIMENSION CXM(0:4*MAXB-1,MAXB)
125    DIMENSION CXMN(0:2*(((MNAUF+1)*MNAUF)/2+MNAUF)+1)
126    REAL*8 Z(MAXB/2,0:((MNAUF+3)*(MNAUF+4))/2)
127    REAL*8 W(MAXB),CR,CI,HILF
128    LOGICAL EVEN
129    !
130    !     DIESE ROUTINE BERECHNET DIE KFFKs CXMN
131    !
132    LL=0
133    LLP=0
134    DO  I=0,MNAUF
135       KM=0
136       DO
137          KM=KM+1
138          IF(MLAT(KM).GT.2*I) THEN
139             EXIT
140          ENDIF
141       ENDDO
142       DO  J=I,MNAUF
143          CR=0
144          CI=0
145          EVEN=MOD(I+J,2).EQ.0
146          IF(EVEN) THEN
147             DO  K=KM,MAXB/2
148                HILF=W(K)*Z(K,LLP)
149                CR=CR+CXM(2*I,K)*HILF
150                CI=CI+CXM(2*I+1,K)*HILF
151             ENDDO
152          ELSE
153             DO  K=KM,MAXB/2
154                HILF=W(K)*Z(K,LLP)
155                CR=CR+CXM(2*I,MAXB+1-K)*HILF
156                CI=CI+CXM(2*I+1,MAXB+1-K)*HILF
157             ENDDO
158          ENDIF
159          CXMN(2*LL)=CR
160          CXMN(2*LL+1)=CI
161          LL=LL+1
162          LLP=LLP+1
163       ENDDO
164       LLP=LLP+2
165    ENDDO
166    RETURN
167  END SUBROUTINE LGTR213
168  !
169
170  !
171  SUBROUTINE RFOUFTR(CXM,TRIGS,IFAX,MNAUF,MAXL,ISIGN)
172    !     BERECHNET DIE FOURIERSUMME MIT EINEM FFT-ALGORITHMUS
173    IMPLICIT REAL (A-H,O-Z)
174    DIMENSION CXM(0:2*MAXL-1)
175    DIMENSION FELD(MAXL),TRIGS(2*MAXL)
176    DIMENSION WSAVE(MAXAUF)
177    INTEGER IFAX(10)
178
179
180    ! NORMIERUNG...
181    WSAVE(1)=CXM(MAXL-1)
182
183    CXM(1:MAXL)=CXM(0:MAXL-1)/2
184    CXM(0)=WSAVE(1)/2
185    CALL FFT99(CXM,WSAVE,TRIGS,IFAX,1,1,MAXL,1,-1)
186    RETURN
187  END SUBROUTINE RFOUFTR
188
189END MODULE GRTOPH
Note: See TracBrowser for help on using the repository browser.
hosted by ZAMG