source: palm/trunk/SOURCE/sum_up_3d_data.f90 @ 72

Last change on this file since 72 was 72, checked in by raasch, 17 years ago

preliminary changes for precipitation output

  • Property svn:keywords set to Id
File size: 12.2 KB
Line 
1 SUBROUTINE sum_up_3d_data
2
3!------------------------------------------------------------------------------!
4! Actual revisions:
5! -----------------
6! +sum-up of precipitation rate and roughness length (prr*, z0*)
7!
8! Former revisions:
9! -----------------
10! $Id: sum_up_3d_data.f90 72 2007-03-19 08:20:46Z raasch $
11! RCS Log replace by Id keyword, revision history cleaned up
12!
13! Revision 1.1  2006/02/23 12:55:23  raasch
14! Initial revision
15!
16!
17! Description:
18! ------------
19! Sum-up the values of 3d-arrays. The real averaging is later done in routine
20! average_3d_data.
21!------------------------------------------------------------------------------!
22
23    USE arrays_3d
24    USE averaging
25    USE cloud_parameters
26    USE control_parameters
27    USE cpulog
28    USE indices
29    USE interfaces
30    USE particle_attributes
31
32    IMPLICIT NONE
33
34    INTEGER ::  i, ii, j, k, n, psi
35
36    REAL    ::  mean_r, s_r3, s_r4
37
38
39    CALL cpu_log (log_point(34),'sum_up_3d_data','start')
40
41!
42!-- Allocate and initialize the summation arrays if called for the very first
43!-- time or the first time after average_3d_data has been called
44!-- (some or all of the arrays may have been already allocated
45!-- in read_3d_binary)
46    IF ( average_count_3d == 0 )  THEN
47
48       DO  ii = 1, doav_n
49
50          SELECT CASE ( TRIM( doav(ii) ) )
51
52             CASE ( 'e' )
53                IF ( .NOT. ALLOCATED( e_av ) )  THEN
54                   ALLOCATE( e_av(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
55                ENDIF
56                e_av = 0.0
57
58             CASE ( 'lwp*' )
59                IF ( .NOT. ALLOCATED( lwp_av ) )  THEN
60                   ALLOCATE( lwp_av(nys-1:nyn+1,nxl-1:nxr+1) )
61                ENDIF
62                lwp_av = 0.0
63
64             CASE ( 'p' )
65                IF ( .NOT. ALLOCATED( p_av ) )  THEN
66                   ALLOCATE( p_av(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
67                ENDIF
68                p_av = 0.0
69
70             CASE ( 'pc' )
71                IF ( .NOT. ALLOCATED( pc_av ) )  THEN
72                   ALLOCATE( pc_av(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
73                ENDIF
74                pc_av = 0.0
75
76             CASE ( 'pr' )
77                IF ( .NOT. ALLOCATED( pr_av ) )  THEN
78                   ALLOCATE( pr_av(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
79                ENDIF
80                pr_av = 0.0
81
82             CASE ( 'prr*' )
83                IF ( .NOT. ALLOCATED( precipitation_rate_av ) )  THEN
84                   ALLOCATE( precipitation_rate_av(nys-1:nyn+1,nxl-1:nxr+1) )
85                ENDIF
86                precipitation_rate_av = 0.0
87
88             CASE ( 'pt' )
89                IF ( .NOT. ALLOCATED( pt_av ) )  THEN
90                   ALLOCATE( pt_av(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
91                ENDIF
92                pt_av = 0.0
93
94             CASE ( 'q' )
95                IF ( .NOT. ALLOCATED( q_av ) )  THEN
96                   ALLOCATE( q_av(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
97                ENDIF
98                q_av = 0.0
99
100             CASE ( 'ql' )
101                IF ( .NOT. ALLOCATED( ql_av ) )  THEN
102                   ALLOCATE( ql_av(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
103                ENDIF
104                ql_av = 0.0
105
106             CASE ( 'ql_c' )
107                IF ( .NOT. ALLOCATED( ql_c_av ) )  THEN
108                   ALLOCATE( ql_c_av(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
109                ENDIF
110                ql_c_av = 0.0
111
112             CASE ( 'ql_v' )
113                IF ( .NOT. ALLOCATED( ql_v_av ) )  THEN
114                   ALLOCATE( ql_v_av(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
115                ENDIF
116                ql_v_av = 0.0
117
118             CASE ( 'ql_vp' )
119                IF ( .NOT. ALLOCATED( ql_vp_av ) )  THEN
120                   ALLOCATE( ql_vp_av(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
121                ENDIF
122                ql_vp_av = 0.0
123
124             CASE ( 'qv' )
125                IF ( .NOT. ALLOCATED( qv_av ) )  THEN
126                   ALLOCATE( qv_av(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
127                ENDIF
128                qv_av = 0.0
129
130             CASE ( 's' )
131                IF ( .NOT. ALLOCATED( s_av ) )  THEN
132                   ALLOCATE( s_av(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
133                ENDIF
134                s_av = 0.0
135
136             CASE ( 't*' )
137                IF ( .NOT. ALLOCATED( ts_av ) )  THEN
138                   ALLOCATE( ts_av(nys-1:nyn+1,nxl-1:nxr+1) )
139                ENDIF
140                ts_av = 0.0
141
142             CASE ( 'u' )
143                IF ( .NOT. ALLOCATED( u_av ) )  THEN
144                   ALLOCATE( u_av(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
145                ENDIF
146                u_av = 0.0
147
148             CASE ( 'u*' )
149                IF ( .NOT. ALLOCATED( us_av ) )  THEN
150                   ALLOCATE( us_av(nys-1:nyn+1,nxl-1:nxr+1) )
151                ENDIF
152                us_av = 0.0
153
154             CASE ( 'v' )
155                IF ( .NOT. ALLOCATED( v_av ) )  THEN
156                   ALLOCATE( v_av(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
157                ENDIF
158                v_av = 0.0
159
160             CASE ( 'vpt' )
161                IF ( .NOT. ALLOCATED( vpt_av ) )  THEN
162                   ALLOCATE( vpt_av(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
163                ENDIF
164                vpt_av = 0.0
165
166             CASE ( 'w' )
167                IF ( .NOT. ALLOCATED( w_av ) )  THEN
168                   ALLOCATE( w_av(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
169                ENDIF
170                w_av = 0.0
171
172             CASE ( 'z0*' )
173                IF ( .NOT. ALLOCATED( z0_av ) )  THEN
174                   ALLOCATE( z0_av(nys-1:nyn+1,nxl-1:nxr+1) )
175                ENDIF
176                z0_av = 0.0
177
178             CASE DEFAULT
179!
180!--             User-defined quantity
181                CALL user_3d_data_averaging( 'allocate', doav(ii) )
182
183          END SELECT
184
185       ENDDO
186
187    ENDIF
188
189!
190!-- Loop of all variables to be averaged.
191    DO  ii = 1, doav_n
192
193!
194!--    Store the array chosen on the temporary array.
195       SELECT CASE ( TRIM( doav(ii) ) )
196
197          CASE ( 'e' )
198             DO  i = nxl-1, nxr+1
199                DO  j = nys-1, nyn+1
200                   DO  k = nzb, nzt+1
201                      e_av(k,j,i) = e_av(k,j,i) + e(k,j,i)
202                   ENDDO
203                ENDDO
204             ENDDO
205
206          CASE ( 'lwp*' )
207             DO  i = nxl-1, nxr+1
208                DO  j = nys-1, nyn+1
209                   lwp_av(j,i) = lwp_av(j,i) + SUM( ql(nzb:nzt,j,i) * &
210                                                    dzw(1:nzt+1) )
211                ENDDO
212             ENDDO
213
214          CASE ( 'p' )
215             DO  i = nxl-1, nxr+1
216                DO  j = nys-1, nyn+1
217                   DO  k = nzb, nzt+1
218                      p_av(k,j,i) = p_av(k,j,i) + p(k,j,i)
219                   ENDDO
220                ENDDO
221             ENDDO
222
223          CASE ( 'pc' )
224             DO  i = nxl, nxr
225                DO  j = nys, nyn
226                   DO  k = nzb, nzt+1
227                      pc_av(k,j,i) = pc_av(k,j,i) + prt_count(k,j,i)
228                   ENDDO
229                ENDDO
230             ENDDO
231
232          CASE ( 'pr' )
233             DO  i = nxl, nxr
234                DO  j = nys, nyn
235                   DO  k = nzb, nzt+1
236                      psi = prt_start_index(k,j,i)
237                      s_r3 = 0.0
238                      s_r4 = 0.0
239                      DO  n = psi, psi+prt_count(k,j,i)-1
240                         s_r3 = s_r3 + particles(n)%radius**3
241                         s_r4 = s_r4 + particles(n)%radius**4
242                      ENDDO
243                      IF ( s_r3 /= 0.0 )  THEN
244                         mean_r = s_r4 / s_r3
245                      ELSE
246                         mean_r = 0.0
247                      ENDIF
248                      pr_av(k,j,i) = pr_av(k,j,i) + mean_r
249                   ENDDO
250                ENDDO
251             ENDDO
252
253          CASE ( 'pr*' )
254             DO  i = nxl-1, nxr+1
255                DO  j = nys-1, nyn+1
256                   precipitation_rate_av(j,i) = precipitation_rate_av(j,i) + &
257                                                precipitation_rate(j,i)
258                ENDDO
259             ENDDO
260
261          CASE ( 'pt' )
262             IF ( .NOT. cloud_physics ) THEN
263                DO  i = nxl-1, nxr+1
264                   DO  j = nys-1, nyn+1
265                      DO  k = nzb, nzt+1
266                         pt_av(k,j,i) = pt_av(k,j,i) + pt(k,j,i)
267                      ENDDO
268                   ENDDO
269                ENDDO
270             ELSE
271                DO  i = nxl-1, nxr+1
272                   DO  j = nys-1, nyn+1
273                      DO  k = nzb, nzt+1
274                         pt_av(k,j,i) = pt_av(k,j,i) + pt(k,j,i) + l_d_cp * &
275                                                       pt_d_t(k) * ql(k,j,i)
276                      ENDDO
277                   ENDDO
278                ENDDO
279             ENDIF
280
281          CASE ( 'q' )
282             DO  i = nxl-1, nxr+1
283                DO  j = nys-1, nyn+1
284                   DO  k = nzb, nzt+1
285                      q_av(k,j,i) = q_av(k,j,i) + q(k,j,i)
286                   ENDDO
287                ENDDO
288             ENDDO
289             
290          CASE ( 'ql' )
291             DO  i = nxl-1, nxr+1
292                DO  j = nys-1, nyn+1
293                   DO  k = nzb, nzt+1
294                      ql_av(k,j,i) = ql_av(k,j,i) + ql(k,j,i)
295                   ENDDO
296                ENDDO
297             ENDDO
298
299          CASE ( 'ql_c' )
300             DO  i = nxl-1, nxr+1
301                DO  j = nys-1, nyn+1
302                   DO  k = nzb, nzt+1
303                      ql_c_av(k,j,i) = ql_c_av(k,j,i) + ql_c(k,j,i)
304                   ENDDO
305                ENDDO
306             ENDDO
307
308          CASE ( 'ql_v' )
309             DO  i = nxl-1, nxr+1
310                DO  j = nys-1, nyn+1
311                   DO  k = nzb, nzt+1
312                      ql_v_av(k,j,i) = ql_v_av(k,j,i) + ql_v(k,j,i)
313                   ENDDO
314                ENDDO
315             ENDDO
316
317          CASE ( 'ql_vp' )
318             DO  i = nxl-1, nxr+1
319                DO  j = nys-1, nyn+1
320                   DO  k = nzb, nzt+1
321                      ql_vp_av(k,j,i) = ql_vp_av(k,j,i) + ql_vp(k,j,i)
322                   ENDDO
323                ENDDO
324             ENDDO
325
326          CASE ( 'qv' )
327             DO  i = nxl-1, nxr+1
328                DO  j = nys-1, nyn+1
329                   DO  k = nzb, nzt+1
330                      qv_av(k,j,i) = qv_av(k,j,i) + q(k,j,i) - ql(k,j,i)
331                   ENDDO
332                ENDDO
333             ENDDO
334
335          CASE ( 's' )
336             DO  i = nxl-1, nxr+1
337                DO  j = nys-1, nyn+1
338                   DO  k = nzb, nzt+1
339                      s_av(k,j,i) = s_av(k,j,i) + q(k,j,i)
340                   ENDDO
341                ENDDO
342             ENDDO
343             
344          CASE ( 't*' )
345             DO  i = nxl-1, nxr+1
346                DO  j = nys-1, nyn+1
347                   ts_av(j,i) = ts_av(j,i) + ts(j,i)
348                ENDDO
349             ENDDO
350
351          CASE ( 'u' )
352             DO  i = nxl-1, nxr+1
353                DO  j = nys-1, nyn+1
354                   DO  k = nzb, nzt+1
355                      u_av(k,j,i) = u_av(k,j,i) + u(k,j,i)
356                   ENDDO
357                ENDDO
358             ENDDO
359
360          CASE ( 'u*' )
361             DO  i = nxl-1, nxr+1
362                DO  j = nys-1, nyn+1
363                   us_av(j,i) = us_av(j,i) + us(j,i)
364                ENDDO
365             ENDDO
366
367          CASE ( 'v' )
368             DO  i = nxl-1, nxr+1
369                DO  j = nys-1, nyn+1
370                   DO  k = nzb, nzt+1
371                      v_av(k,j,i) = v_av(k,j,i) + v(k,j,i)
372                   ENDDO
373                ENDDO
374             ENDDO
375
376          CASE ( 'vpt' )
377             DO  i = nxl-1, nxr+1
378                DO  j = nys-1, nyn+1
379                   DO  k = nzb, nzt+1
380                      vpt_av(k,j,i) = vpt_av(k,j,i) + vpt(k,j,i)
381                   ENDDO
382                ENDDO
383             ENDDO
384
385          CASE ( 'w' )
386             DO  i = nxl-1, nxr+1
387                DO  j = nys-1, nyn+1
388                   DO  k = nzb, nzt+1
389                      w_av(k,j,i) = w_av(k,j,i) + w(k,j,i)
390                   ENDDO
391                ENDDO
392             ENDDO
393
394          CASE ( 'z0*' )
395             DO  i = nxl-1, nxr+1
396                DO  j = nys-1, nyn+1
397                   z0_av(j,i) = z0_av(j,i) + z0(j,i)
398                ENDDO
399             ENDDO
400
401          CASE DEFAULT
402!
403!--          User-defined quantity
404             CALL user_3d_data_averaging( 'sum', doav(ii) )
405
406       END SELECT
407
408    ENDDO
409
410    CALL cpu_log (log_point(34),'sum_up_3d_data','stop','nobarrier')
411
412
413 END SUBROUTINE sum_up_3d_data
Note: See TracBrowser for help on using the repository browser.