Home

Context Navigation

source: palm/trunk/UTIL/analyze_particle_netcdf_data.f90 @ 1026

Last change on this file since 1026 was 1, checked in by raasch, 18 years ago
Initial repository layout and content
File size: 8.0 KB

Line
1	PROGRAM analyze_particle_netcdf_data
2
3	!------------------------------------------------------------------------------!
4	! Actual revisions:
5	! -----------------
6	!
7	!
8	! Former revisions:
9	! -----------------
10	! Initial revision 12/07/05
11	!
12	!
13	! Description:
14	! ------------
15	! This is an EXAMPLE routine how to read NetCDF particle data output by PALM.
16	! As an example, the mean heigth and radius of all particles are calculated
17	! and output for each output time available on the file
18	!
19	!
20	! Any additonal analyzation requested has to be added by the user by following
21	! the steps given in this example!
22	!
23	!
24	!
25	! This routine must be compiled with:
26	! decalpha:
27	! f95 -fast -r8 -I/usr/local/netcdf-3.5.1/include
28	! -L/usr/local/netcdf-3.5.1/lib -lnetcdf
29	! IBM-Regatta:
30	! xlf95 -qsuffix=cpp=f90 -qrealsize=8 -q64 -qmaxmem=-1 -Q
31	! -I /aws/dataformats/netcdf-3.5.0/netcdf-64-32-3.5.0/include_F90_64
32	! -L/aws/dataformats/netcdf-3.5.0/netcdf-64-32-3.5.0/lib -lnetcdf -O3
33	! IBM-Regatta KISTI:
34	! xlf95 -qsuffix=cpp=f90 -qrealsize=8 -q64 -qmaxmem=-1 -Q
35	! -I /applic/netcdf64/src/f90
36	! -L/applic/lib/NETCDF64 -lnetcdf -O3
37	! IMUK:
38	! ifort analyze_particle_netcdf_data
39	! -I /muksoft/packages/netcdf/linux/include -axW -r8 -nbs -Vaxlib
40	! -L /muksoft/packages/netcdf/linux/lib -lnetcdf
41	! NEC-SX6:
42	! sxf90 analyze_particle_netcdf_data.f90
43	! -I /pool/SX-6/netcdf/netcdf-3.6.0-p1/include -C hopt -Wf '-A idbl4'
44	! -D__netcdf -L/pool/SX-6/netcdf/netcdf-3.6.0-p1/lib -lnetcdf
45	!------------------------------------------------------------------------------!
46
47	USE netcdf
48
49	IMPLICIT NONE
50
51	!
52	!-- Local variables
53	CHARACTER (LEN=7) :: id_char
54	CHARACTER (LEN=2000) :: string
55
56
57	INTEGER :: f, fn, i, id_dim_time, id_var_rnop, id_var_r, id_var_time, &
58	id_var_x, id_var_y, id_var_z, n, nc_stat, ntl, start
59
60	INTEGER, DIMENSION(1000) :: id_set
61	INTEGER, DIMENSION(1) :: id_dim_time_old
62
63	INTEGER, DIMENSION(:), ALLOCATABLE :: total_nop
64	INTEGER, DIMENSION(:,:), ALLOCATABLE :: nop
65
66	LOGICAL :: found
67
68	REAL :: mean_r, mean_z
69	REAL, DIMENSION(:), ALLOCATABLE :: prt_r, prt_x, prt_y, prt_z, tl
70
71	!
72	!-- Check, if file from PE0 exists. If it does not exist, PALM did not
73	!-- create any output for this cross-section.
74	fn = 0
75	WRITE (id_char,'(''_'',I4.4)') fn
76	INQUIRE ( FILE=id_char, EXIST=found )
77
78	!
79	!-- Find out the number of files (equal to the number of PEs which
80	!-- have been used in PALM) and open them
81	IF ( .NOT. found ) THEN
82	PRINT*, '+++ no file _0000 found in current working directory'
83	STOP
84	ENDIF
85
86	DO WHILE ( found )
87
88	!
89	!-- Open NetCDF dataset
90	nc_stat = NF90_OPEN( id_char, NF90_NOWRITE, id_set(fn) )
91	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( nc_stat, 1 )
92	fn = fn + 1
93	WRITE (id_char,'(''_'',I4.4)') fn
94	INQUIRE ( FILE=id_char, EXIST=found )
95
96	ENDDO
97	fn = fn - 1
98
99	PRINT, '** ', fn+1, ' files found'
100
101	!
102	!-- Get the run description header and print it out
103	string = ' '
104	nc_stat = NF90_GET_ATT( id_set(0), NF90_GLOBAL, 'title', string )
105	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( nc_stat, 2 )
106	PRINT, '** run: ', TRIM( string )
107
108	!
109	!-- Get the available time levels
110	nc_stat = NF90_INQ_VARID( id_set(0), 'time', id_var_time )
111	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( nc_stat, 3 )
112
113	nc_stat = NF90_INQUIRE_VARIABLE( id_set(0), id_var_time, &
114	dimids = id_dim_time_old )
115	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( nc_stat, 4 )
116	id_dim_time = id_dim_time_old(1)
117
118	nc_stat = NF90_INQUIRE_DIMENSION( id_set(0), id_dim_time, len = ntl )
119	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( nc_stat, 5 )
120	ALLOCATE( tl(1:ntl) )
121	print*, 'ntl=',ntl
122
123	nc_stat = NF90_GET_VAR( id_set(0), id_var_time, tl )
124	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( nc_stat, 6 )
125
126	DO n = 1, ntl
127	print, '** time_level(', n, ') =', tl(n)
128	ENDDO
129
130	!
131	!-- Get the number of particles used
132	nc_stat = NF90_INQ_VARID( id_set(0), 'real_num_of_prt', id_var_rnop )
133	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( nc_stat, 7 )
134
135	ALLOCATE( nop(1:ntl,0:fn), total_nop(1:ntl) )
136
137	DO f = 0, fn
138
139	nc_stat = NF90_GET_VAR( id_set(f), id_var_rnop, nop(1:ntl,f) )
140	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( nc_stat, 8 )
141
142	ENDDO
143
144	total_nop = 0
145	DO n = 1, ntl
146
147	DO f = 0, fn
148	total_nop(n) = total_nop(n) + nop(n,f)
149	ENDDO
150
151	PRINT, '** time = ', tl(n), ' total # of particles: ', total_nop(n)
152
153	ENDDO
154
155	!
156	!-- Get the particle x and y coordinates
157	nc_stat = NF90_INQ_VARID( id_set(0), 'pt_x', id_var_x )
158	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( nc_stat, 9 )
159
160	nc_stat = NF90_INQ_VARID( id_set(0), 'pt_y', id_var_y )
161	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( nc_stat, 10 )
162
163	nc_stat = NF90_INQ_VARID( id_set(0), 'pt_z', id_var_z )
164	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( nc_stat, 11 )
165
166	nc_stat = NF90_INQ_VARID( id_set(0), 'pt_radius', id_var_r )
167	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( nc_stat, 12 )
168
169	PRINT*, ' '
170	!
171	!-- Loop over all timelevels
172	DO n = 1, ntl
173
174	ALLOCATE( prt_x(total_nop(n)), prt_y(total_nop(n)), &
175	prt_z(total_nop(n)), prt_r(total_nop(n)) )
176
177	start = 1
178
179	!
180	!-- Read the data from the files (one file per processor)
181	DO f = 0, fn
182
183	nc_stat = NF90_GET_VAR( id_set(f), id_var_x, &
184	prt_x(start:start+nop(n,f)-1), &
185	start = (/ 1, n /), &
186	count = (/ nop(n,f), 1 /) )
187	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( nc_stat, 13 )
188
189	nc_stat = NF90_GET_VAR( id_set(f), id_var_y, &
190	prt_y(start:start+nop(n,f)-1), &
191	start = (/ 1, n /), &
192	count = (/ nop(n,f), 1 /) )
193	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( nc_stat, 14 )
194
195	nc_stat = NF90_GET_VAR( id_set(f), id_var_z, &
196	prt_z(start:start+nop(n,f)-1), &
197	start = (/ 1, n /), &
198	count = (/ nop(n,f), 1 /) )
199	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( nc_stat, 15 )
200
201	nc_stat = NF90_GET_VAR( id_set(f), id_var_r, &
202	prt_r(start:start+nop(n,f)-1), &
203	start = (/ 1, n /), &
204	count = (/ nop(n,f), 1 /) )
205	IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( nc_stat, 16 )
206
207	start = start + nop(n,f)
208
209	ENDDO
210
211	mean_z = 0.0
212	mean_r = 0.0
213	DO i = 1, total_nop(n)
214	mean_z = mean_z + prt_z(i)
215	mean_r = mean_r + prt_r(i)
216	ENDDO
217	mean_z = mean_z / total_nop(n)
218	mean_r = mean_r / total_nop(n)
219
220	PRINT, '** time = ', tl(n), ' mean height = ', mean_z, &
221	' mean radius = ', mean_r
222
223	!
224	!-- prt_x, prt_y, prt_z, and prt_r contain the particle coordinates and
225	!-- radii, respectively. Please output these data or carry out the
226	!-- requested analyzing in this program before you deallocate the arrays.
227
228	DEALLOCATE( prt_x, prt_y, prt_z, prt_r )
229
230	ENDDO
231
232
233	END PROGRAM analyze_particle_netcdf_data
234
235
236
237	SUBROUTINE handle_netcdf_error( nc_stat, position )
238	!
239	!-- Prints out a text message corresponding to the current NetCDF status
240
241	USE netcdf
242
243	IMPLICIT NONE
244
245	INTEGER, INTENT(IN) :: nc_stat, position
246
247	IF ( nc_stat /= NF90_NOERR ) THEN
248	PRINT*, '+++ analyze_particle_netcdf_data'
249	PRINT*, ' netcdf error: ', TRIM( nf90_strerror( nc_stat ) )
250	PRINT*, ' position = ', position
251	STOP
252	ENDIF
253
254	END SUBROUTINE handle_netcdf_error

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |