source: palm/trunk/UTIL/read_prt_data.f90 @ 2732

Last change on this file since 2732 was 2718, checked in by maronga, 7 years ago

deleting of deprecated files; headers updated where needed

File size: 7.0 KB
Line 
1 PROGRAM read_prt_data
2!------------------------------------------------------------------------------!
3! This file is part of the PALM model system.
4!
5! PALM is free software: you can redistribute it and/or modify it under the
6! terms of the GNU General Public License as published by the Free Software
7! Foundation, either version 3 of the License, or (at your option) any later
8! version.
9!
10! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
11! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
12! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
13!
14! You should have received a copy of the GNU General Public License along with
15! PALM. If not, see <http://www.gnu.org/licenses/>.
16!
17! Copyright 1997-2018  Leibniz Universitaet Hannover
18!------------------------------------------------------------------------------!
19!
20! Current revisions:
21! -----------------
22!
23!
24! Former revisions:
25! -----------------
26! $Id: read_prt_data.f90 2696 2017-12-14 17:12:51Z kanani $
27! Corrected "Former revisions" section
28!
29! 2696 2017-12-14 17:12:51Z kanani
30! Change in file header (GPL part)
31!
32! 2123 2017-01-18 12:49:59Z hoffmann
33! Particle variables renamed: tailpoints, tail_id, and dvrp_psize to id1, id2,
34! and user
35!
36! 1771 2016-02-29 10:57:56Z hoffmann
37! Initial revision.
38!
39! Description:
40! ------------
41! This routine reads the particle data generated by PALM, and enables user
42! analysis. Compile and execute this program in the folder where your particle
43! data (_000000, _000001, ...) is stored.
44!------------------------------------------------------------------------------!
45
46    IMPLICIT NONE
47
48    CHARACTER(LEN=7)    ::  i_proc_char
49    CHARACTER (LEN=80)  ::  rtext
50    CHARACTER (LEN=110) ::  run_description_header
51
52    REAL(KIND=8) ::  simulated_time
53    INTEGER      ::  ip, i_proc=0, i_proc_end, jp, kp,                         &
54                     max_number_of_particle_groups, number_of_particles,       &
55                     number_of_particle_groups, n, nxl,                        &
56                     nxr, nys, nyn, nzb, nzt, nbgp, status
57
58    INTEGER, DIMENSION(:,:,:), ALLOCATABLE ::  prt_count
59
60    TYPE particle_type
61       SEQUENCE
62       REAL(KIND=8) ::  radius, age, age_m, dt_sum, user, e_m,                 &
63                        origin_x, origin_y, origin_z, rvar1, rvar2, rvar3,     &
64                        speed_x, speed_y, speed_z, weight_factor, x, y, z
65      INTEGER       ::  class, group, id1, id2
66      LOGICAL       ::  particle_mask
67      INTEGER       ::  block_nr
68    END TYPE particle_type
69
70    TYPE(particle_type), DIMENSION(:), ALLOCATABLE ::  particles_temp
71    TYPE(particle_type), DIMENSION(:), POINTER     ::  particles
72
73    TYPE  grid_particle_def
74       TYPE(particle_type), DIMENSION(:), ALLOCATABLE ::  particles
75    END TYPE grid_particle_def
76
77    TYPE(grid_particle_def), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  grid_particles
78
79    TYPE particle_groups_type
80        SEQUENCE
81        REAL(KIND=8) ::  density_ratio, radius, exp_arg, exp_term
82    END TYPE particle_groups_type
83
84    TYPE(particle_groups_type), DIMENSION(:), ALLOCATABLE ::  particle_groups
85
86    LOGICAL ::  found
87!
88!-- Check if file from PE0 exists and terminate program if it doesn't.
89    WRITE (i_proc_char,'(''_'',I6.6)')  i_proc
90    INQUIRE ( FILE=i_proc_char, EXIST=found )
91!
92!-- Estimate the number of files (equal to the number of PEs which
93!-- have been used in PALM)
94    DO  WHILE ( found )
95       OPEN ( 85, FILE=i_proc_char, FORM='UNFORMATTED' )
96       i_proc = i_proc + 1
97       WRITE (i_proc_char,'(''_'',I6.6)')  i_proc
98       INQUIRE ( FILE=i_proc_char, EXIST=found )
99       CLOSE( 85 )
100    ENDDO
101!
102!-- Info-output
103    PRINT*, ''
104    PRINT*, '*** read_prt ***'
105    IF ( i_proc /= 0 )  THEN
106       PRINT*, '***', i_proc, ' file(s) found'
107    ELSE
108       PRINT*, '+++ file _000000 not found'
109       PRINT*, '    program terminated'
110       STOP
111    ENDIF
112!
113!-- Set number of files and opens them sequentially
114    i_proc_end = i_proc-1
115
116    DO  i_proc = 0, i_proc_end
117!
118!--    Open particle data file for each process
119       WRITE (i_proc_char,'(''_'',I6.6)')  i_proc
120       OPEN ( 85, FILE=i_proc_char, FORM='UNFORMATTED' ) 
121!
122!--    Read general description of file and inform user
123       READ ( 85 )  run_description_header
124       READ ( 85 )  rtext
125       IF ( i_proc == 0 )  THEN
126          PRINT*, ' '
127          PRINT*, '*** data description header:'
128          PRINT*, '*** ', run_description_header
129          PRINT*, '*** ', rtext
130          PRINT*, ' '
131       ENDIF
132       PRINT*, '*** data of file', i_proc, 'is analysed'
133       READ ( 85 )  number_of_particle_groups, max_number_of_particle_groups
134
135       IF ( .NOT. ALLOCATED(particle_groups) )  THEN
136          ALLOCATE( particle_groups(max_number_of_particle_groups) )
137       ENDIF
138   
139       READ ( 85 )  particle_groups
140       READ ( 85 )  nxl, nxr, nys, nyn, nzb, nzt, nbgp
141
142       IF ( ALLOCATED(prt_count)  .OR.  ALLOCATED(grid_particles) )  THEN
143          DEALLOCATE( prt_count, grid_particles )
144       ENDIF
145
146       ALLOCATE( prt_count(nzb:nzt+1,nys-nbgp:nyn+nbgp,nxl-nbgp:nxr+nbgp),     &
147                 grid_particles(nzb+1:nzt,nys:nyn,nxl:nxr) )
148
149!
150!--    Start individual time loop
151       DO
152          READ( 85, IOSTAT=status)  simulated_time
153          IF ( status < 0 )  THEN
154             PRINT*, '*** end of file reached'
155             EXIT
156          ENDIF
157          PRINT*, '*** time of analyzed data set:', simulated_time
158
159          READ ( 85 )  prt_count
160
161          DO  ip = nxl, nxr
162             DO  jp = nys, nyn
163                DO  kp = nzb+1, nzt
164                   number_of_particles = prt_count(kp,jp,ip)
165                   IF ( number_of_particles <= 0 )  CYCLE
166                   ALLOCATE( grid_particles(kp,jp,ip)%particles(1:number_of_particles) )
167                   ALLOCATE( particles_temp(1:number_of_particles) )
168                   READ ( 85 )  particles_temp
169                   grid_particles(kp,jp,ip)%particles = particles_temp
170                   DEALLOCATE( particles_temp )
171                ENDDO
172             ENDDO
173          ENDDO
174!
175!--       This part can be used for user analysis
176          DO  ip = nxl, nxr
177             DO  jp = nys, nyn
178                DO  kp = nzb+1, nzt
179                   number_of_particles = prt_count(kp,jp,ip)
180                   IF ( number_of_particles <= 0 )  CYCLE
181                   particles => grid_particles(kp,jp,ip)%particles
182!
183!--                Add your analysis here
184!                   DO  n = 1, number_of_particles
185!
186!                   ENDDO
187
188                ENDDO
189             ENDDO
190          ENDDO
191!
192!--       Deallocate grid_particles%particles in case of more than one timestep
193          DO  ip = nxl, nxr
194             DO  jp = nys, nyn
195                DO  kp = nzb+1, nzt
196                   number_of_particles = prt_count(kp,jp,ip)
197                   IF ( number_of_particles <= 0 )  CYCLE
198                   DEALLOCATE( grid_particles(kp,jp,ip)%particles )
199                ENDDO
200             ENDDO
201          ENDDO     
202
203       ENDDO ! end of time loop
204
205       CLOSE( 85 )
206
207    ENDDO ! end of file loop
208
209 END PROGRAM read_prt_data
Note: See TracBrowser for help on using the repository browser.