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

Last change on this file since 1032 was 668, checked in by suehring, 13 years ago

last commit documented

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