source: palm/trunk/SOURCE/write_compressed.f90 @ 1074

Last change on this file since 1074 was 1037, checked in by raasch, 12 years ago

last commit documented

  • Property svn:keywords set to Id
File size: 6.6 KB
Line 
1 SUBROUTINE write_compressed( field, fid_avs, fid_fld, my_id, nxl, nxr, nyn, &
2                              nys, nzb, nz_do3d, prec, nbgp )
3
4!--------------------------------------------------------------------------------!
5! This file is part of PALM.
6!
7! PALM is free software: you can redistribute it and/or modify it under the terms
8! of the GNU General Public License as published by the Free Software Foundation,
9! either version 3 of the License, or (at your option) any later version.
10!
11! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
12! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
13! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
14!
15! You should have received a copy of the GNU General Public License along with
16! PALM. If not, see <http://www.gnu.org/licenses/>.
17!
18! Copyright 1997-2012  Leibniz University Hannover
19!--------------------------------------------------------------------------------!
20!
21! Current revisions:
22! -----------------
23!
24! Former revisions:
25! ---------------------
26! $Id: write_compressed.f90 1037 2012-10-22 14:10:22Z hoffmann $
27!
28! 1036 2012-10-22 13:43:42Z raasch
29! code put under GPL (PALM 3.9)
30!
31! 667 2010-12-23 12:06:00Z suehring/gryschka
32! Array bounds and nx, ny adapted with nbgp
33!
34! 622 2010-12-10 08:08:13Z raasch
35! optional barriers included in order to speed up collective operations
36!
37! Feb. 2007
38! RCS Log replace by Id keyword, revision history cleaned up
39!
40! Revision 1.4  2006/02/23 13:15:09  raasch
41! nz_plot3d renamed nz_do3d
42!
43! Revision 1.1  1999/03/02 09:25:21  raasch
44! Initial revision
45!
46!
47! Description:
48! ------------
49! In this routine, 3D-data (to be plotted) are scaled and compressed by
50! the method of bit shifting. It is designed for the use outside of PALM
51! also, which is the reason why most of the data is passed by subroutine
52! arguments. Nevertheless, the module pegrid is needed by MPI calls.
53!
54! Arguments:
55! field         = data array to be compressed
56! fid_avs       = file-ID of AVS-data-file
57! fid_fld       = file-ID of AVS-header-file
58! my_id         = ID of the calling PE
59! nxl, nxr      = index bounds of the subdomain along x
60! nyn, nys      = index bounds of the subdomain along y
61! nzb,nz_do3d   = index bounds of the domain along z (can be smaller than
62!                 the total domain)
63! prec          = precision of packed data (number of digits after decimal
64!                 point)
65!------------------------------------------------------------------------------!
66
67    USE pegrid         !  needed for MPI_ALLREDUCE
68
69    IMPLICIT NONE
70
71    INTEGER, PARAMETER   :: ip4 = SELECTED_INT_KIND ( 9 )
72    INTEGER, PARAMETER   :: spk = SELECTED_REAL_KIND( 6 )
73
74    INTEGER ::  ampl, dummy1, dummy2, factor, i, ifieldmax, ifieldmax_l, &
75                ifieldmin, ifieldmin_l, ii, j, k, length, nfree, npack, nsize, &
76                nx, ny, nz, pos, startpos
77    INTEGER(ip4) ::  imask (32)
78    INTEGER(ip4), DIMENSION(:), ALLOCATABLE ::  ifield, packed_ifield
79
80    INTEGER, INTENT(IN) ::  fid_avs, fid_fld, my_id, nxl, nxr, nyn, nys, nzb, &
81                            nz_do3d, prec, nbgp
82
83    REAL(spk), INTENT(IN) :: field(1:((nxr-nxl+1+2*nbgp)*(nyn-nys+1+2*nbgp)*(nz_do3d-nzb+1)))
84
85!
86!-- Initialise local variables
87    ampl      = 0
88    ifieldmax = 0
89    ifieldmin = 0
90    npack = 0
91    nsize = 0
92    DO  i = 1,32
93       imask(i) = (2**i) - 1
94    ENDDO
95
96    nx     = nxr - nxl + 2*nbgp
97    ny     = nyn - nys + 2*nbgp
98    nz     = nz_do3d - nzb
99    length = (nx+1) * (ny+1) * (nz+1)
100
101!
102!-- Allocate memory for integer array
103    ALLOCATE ( ifield(1:length) )
104
105!
106!-- Store data on integer (in desired precision)
107    factor = 10**prec
108    DO  i = 1, length
109       ifield(i) = NINT( field(i) * factor )
110    ENDDO
111
112!
113!-- Find minimum and maximum
114    ifieldmax_l = MAXVAL( ifield )
115    ifieldmin_l = MINVAL( ifield )
116
117#if defined( __parallel )
118    IF ( collective_wait )  CALL MPI_BARRIER( comm2d, ierr )
119    CALL MPI_ALLREDUCE( ifieldmax_l, ifieldmax, 1, MPI_INTEGER, MPI_MAX, &
120                        comm2d, ierr )
121    IF ( collective_wait )  CALL MPI_BARRIER( comm2d, ierr )
122    CALL MPI_ALLREDUCE( ifieldmin_l, ifieldmin, 1, MPI_INTEGER, MPI_MIN, &
123                        comm2d, ierr )
124#else
125    ifieldmax = ifieldmax_l
126    ifieldmin = ifieldmin_l
127#endif
128
129!
130!-- Minimum scaling
131    ifield = ifield - ifieldmin
132
133!
134!-- Calculate the number of bits needed for each value
135    ampl  = ifieldmax - ifieldmin
136    nsize = 1
137
138    DO WHILE ( imask(nsize) < ampl )
139       nsize = nsize + 1
140    ENDDO
141
142!
143!-- Calculate size of the packed array
144    npack = length * nsize
145    IF ( MOD( npack, 32 ) /= 0 )  npack = npack + 32
146    npack = npack / 32
147
148!
149!-- Start packing the data
150    ALLOCATE ( packed_ifield(1:npack) )
151    packed_ifield = 0
152
153!
154!-- Starting position of a word
155    startpos = 0
156
157!
158!-- Starting position of the word to which data are actually written
159    ii = 1
160
161!
162!-- Compress all data
163    DO  i = 1, length
164
165!
166!--    Cut the significant bits from the actual grid point value (GPV)
167       dummy1 = IAND( ifield(i), imask(nsize) )
168
169!
170!--    Calculate number of free bits of the actual word after packing the GPV
171       nfree = 32 - startpos - nsize
172
173       IF ( nfree > 0 )  THEN 
174!
175!--       GPV fits to the actual word (ii), additional bits are still free.
176!--       Shift GPV to the new position
177          dummy2 = ISHFT( dummy1 ,nfree )
178
179!
180!--       Append bits to the already packed data
181          packed_ifield(ii) = packed_ifield(ii) + dummy2
182
183!
184!--       Calculate new starting position
185          startpos = startpos + nsize
186
187       ELSEIF ( nfree .EQ. 0 )  THEN 
188!
189!--       GPV fills the actual word (ii) exactly
190          packed_ifield(ii) = packed_ifield(ii) + dummy1
191
192!
193!--       Activate next (new) word
194          ii = ii + 1
195
196!
197!--       Reset starting position of the new word
198          startpos = 0
199
200       ELSE 
201!
202!--       GPV must be split up to the actual (ii) and the following (ii+1)
203!--       word. Shift first part of GPV to its position.
204          dummy2 = ISHFT( dummy1, nfree )
205
206!
207!--       Append bits
208          packed_ifield(ii) = packed_ifield(ii) + dummy2 
209
210!
211!--       Store rest of GPV on the next word
212          ii = ii + 1
213          packed_ifield(ii) = ISHFT( dummy1, 32+nfree )
214!
215!--       Calculate starting position of the next GPV
216          startpos = -nfree
217
218       ENDIF
219
220    ENDDO
221
222!
223!-- Write the compressed 3D array
224    WRITE ( fid_avs )  packed_ifield
225
226!
227!-- Write additional informations on  FLD-file
228    IF ( my_id == 0 )  WRITE ( fid_fld, 100 )  prec, ifieldmin, nsize, length
229
230    DEALLOCATE( ifield, packed_ifield )
231
232!
233!-- Formats
234100 FORMAT ('# precision = ',I4/ &
235            '# feldmin   = ',I8/ &
236            '# nbits     = ',I2/ &
237            '# nskip     = ',I8)
238
239END SUBROUTINE write_compressed
Note: See TracBrowser for help on using the repository browser.