FLEXPART v9.3.0
 All Classes Files Functions Variables
processmetfields.F90
Go to the documentation of this file.
1 subroutine processmetfields(ind,metdata_format)
2  ! i o
3  !*****************************************************************************
4  ! *
5 
6  ! This subrotuine produces all the meteorological data in binary (fp) *
7  ! format. It includes all the subroutines that the traditional getfields *
8  ! would use to generate the meteorological data *
9  ! *
10  ! This routine is essential part of the grib2flexpart utility *
11  !
12  ! D. Morton (Boreal Scienctific Computing) *
13  ! D. Arnold (ZAMG) *
14  ! 15 October 2015 *
15  ! *
16  !*****************************************************************************
17  !*****************************************************************************
18 
19  ! Changes Arnold, D. and Morton, D. (2015): *
20  ! -- description of local and common variables *
21  !*****************************************************************************
22 
23 
24  use par_mod
25  use com_mod
26 
27  implicit none
28 
29  !****************************************************************************
30  ! Subroutine Parameters: *
31  ! input *
32  ! nstop > 0, if trajectory has to be terminated *
33  ! output *
34  ! metdata_format 0 = unknown, 1 = ecmwf, 2 = gfs meteorological data *
35  !
36  integer :: ind, metdata_format, i, lastslash
37  integer :: dumpdata
38  character(len=512):: fpfname, dumppath, filename
39 
40 
41  !****************************************************************************
42 
43  !****************************************************************************
44  ! Local variables *
45  !
46  ! indj indicates the number of the wind field to be read in *
47  ! memaux auxiliary variable to shuffle winds in memory *
48  ! indmin remembers the number of wind fields already treated *
49  ! uu(0:nxmax,0:nymax,nuvzmax,2) wind components in x-direction [m/s] *
50  ! vv(0:nxmax,0:nymax,nuvzmax,2) wind components in y-direction [m/s] *
51  ! ww(0:nxmax,0:nymax,nwzmax,2) wind components in z-direction [deltaeta/s]*
52  ! tt(0:nxmax,0:nymax,nuvzmax,2) temperature [K] *
53  ! ps(0:nxmax,0:nymax,2) surface pressure [Pa] *
54  !
55  integer :: memaux
56  real :: uuh(0:nxmax-1,0:nymax-1,nuvzmax)
57  real :: vvh(0:nxmax-1,0:nymax-1,nuvzmax)
58  real :: pvh(0:nxmax-1,0:nymax-1,nuvzmax)
59  real :: wwh(0:nxmax-1,0:nymax-1,nwzmax)
60  real :: uuhn(0:nxmaxn-1,0:nymaxn-1,nuvzmax,maxnests)
61  real :: vvhn(0:nxmaxn-1,0:nymaxn-1,nuvzmax,maxnests)
62  real :: pvhn(0:nxmaxn-1,0:nymaxn-1,nuvzmax,maxnests)
63  real :: wwhn(0:nxmaxn-1,0:nymaxn-1,nwzmax,maxnests)
64  !****************************************************************************
65 
66 
67  !****************************************************************************
68  ! Global variables *
69  ! from par_mod and com_mod *
70  ! nx,ny,nuvz,nwz field dimensions in x,y and z direction *
71  ! indmin remembers the number of wind fields already treated *
72  ! memind(2) pointer, on which place the wind fields are stored *
73  ! memtime(2) [s] times of the wind fields, which are kept in memory *
74  ! ndinterval [s] time difference between the two wind fields read in *
75  ! ldirect 1 for forward, -1 for backward simulation *
76  ! numbwf actual number of wind fields *
77  ! wftime(maxwf) [s] times relative to beginning time of wind fields *
78  ! wfname(maxwf) name of met file (used for performance timing out)*
79  !****************************************************************************
80 
81 
82 #ifdef PERFTIMER
83  INTEGER millisecs_start, millisecs_stop, count_rate, count_max
84 #endif
85 
86 !-----------------------------------------------------------------------------
87 
88 
89  ! Current time is after 2nd wind field
90  ! -> Resort wind field pointers, so that current time is between 1st and 2nd
91  !***************************************************************************
92 
93  !memaux=memind(1)
94  !memind(1)=memind(2)
95  !memind(2)=memaux
96  !memtime(1)=memtime(2)
97 
98 
99  ! Read a new wind field and store it on place memind(2)
100  !******************************************************
101 #ifdef PERFTIMER
102  CALL system_clock(millisecs_start, count_rate, count_max)
103 #endif
104  if (metdata_format.eq.ecmwf_metdata) then
105  call readwind_ecmwf(ind,memind(1),uuh,vvh,wwh)
106  call readwind_nests(ind,memind(1),uuhn,vvhn,wwhn)
107  call calcpar_ecmwf(memind(1),uuh,vvh,pvh)
108  call calcpar_nests(memind(1),uuhn,vvhn,pvhn,metdata_format)
109  call verttransform_ecmwf(memind(1),uuh,vvh,wwh,pvh)
110  call verttransform_nests(memind(1),uuhn,vvhn,wwhn,pvhn)
111  memtime(1)=wftime(ind)
112  endif
113  if (metdata_format.eq.gfs_metdata) then
114  call readwind_gfs(ind,memind(1),uuh,vvh,wwh)
115  call readwind_nests(ind,memind(1),uuhn,vvhn,wwhn)
116  call calcpar_gfs(memind(1),uuh,vvh,pvh)
117  call calcpar_nests(memind(1),uuhn,vvhn,pvhn,metdata_format)
118  call verttransform_gfs(memind(1),uuh,vvh,wwh,pvh)
119  call verttransform_nests(memind(1),uuhn,vvhn,wwhn,pvhn)
120  memtime(1)=wftime(ind)
121  endif
122 
123 #ifdef PERFTIMER
124  CALL system_clock(millisecs_stop, count_rate, count_max)
125  print *, 'Wall time to process: ', trim(wfname(ind)), &
126  ': ', (millisecs_stop-millisecs_start)/1000.0, ' seconds'
127 #endif
128 
129 
130  !lastSlash = 0
131  !do i=1,len(wfname(ind))
132  ! if (wfname(ind)(i:i).eq.'/') then
133  ! lastSlash = i
134  ! end if
135  !end do
136  filename = wfname(ind)(lastslash+1:len(wfname(ind)))
137 
138  !if ( ldirect.eq.1 ) then
139  ! fpfname = TRIM(filename) // '_fwd.fp'
140  !else
141  ! fpfname = TRIM(filename) // '_bwd.fp'
142  !endif
143  !print *, 'writing ', TRIM(dumpPath) // '/' // TRIM(fpfname)
144 
145 #ifdef PERFTIMER
146  CALL system_clock(millisecs_start, count_rate, count_max)
147 #endif
148 
149  !CALL fpmetbinary_dump( TRIM(dumpPath) // '/' // TRIM(fpfname), memind(1))
150 
151 #ifdef PERFTIMER
152  !CALL SYSTEM_CLOCK(millisecs_stop, count_rate, count_max)
153  !PRINT *, 'Wall time to process: ', &
154  ! ': ', (millisecs_stop-millisecs_start)/1000.0, ' seconds'
155 #endif
156 
157 end subroutine processmetfields
158 
159 
subroutine readwind_nests(indj, n, uuhn, vvhn, wwhn)
subroutine verttransform_ecmwf(n, uuh, vvh, wwh, pvh)
subroutine verttransform_nests(n, uuhn, vvhn, wwhn, pvhn)
subroutine calcpar_nests(n, uuhn, vvhn, pvhn, metdata_format)
subroutine readwind_gfs(indj, n, uuh, vvh, wwh)
subroutine readwind_ecmwf(indj, n, uuh, vvh, wwh)
Definition: readwind.F90:22
subroutine processmetfields(ind, metdata_format)
subroutine calcpar_gfs(n, uuh, vvh, pvh)
Definition: calcpar_gfs.f90:22
subroutine calcpar_ecmwf(n, uuh, vvh, pvh)
Definition: calcpar.f90:22
subroutine verttransform_gfs(n, uuh, vvh, wwh, pvh)