NCEPLIBS-bufr  12.1.0
cidecode.F90
Go to the documentation of this file.
1 
5 
25 subroutine upc(chr,nchr,ibay,ibit,cnvnull)
26 
27  use modv_vars, only: nbytw, iordle, iordbe
28 
29  implicit none
30 
31  character*(*), intent(out) :: chr
32  character*1 cval(8)
33 
34  integer, intent(in) :: nchr, ibay(*)
35  integer, intent(inout) :: ibit
36  integer ival(2), lb, i, numchr
37 
38  logical, intent(in) :: cnvnull
39 
40  equivalence(cval,ival)
41 
42  ! Set lb to point to the "low-order" (i.e. least significant) byte within a machine word.
43 
44 #ifdef BIG_ENDIAN
45  lb = iordbe(nbytw)
46 #else
47  lb = iordle(nbytw)
48 #endif
49 
50  cval = ' '
51 
52  numchr = min(nchr,len(chr))
53  do i=1,numchr
54  call upb(ival(1),8,ibay,ibit)
55  if((ival(1)==0).and.(cnvnull)) then
56  chr(i:i) = ' '
57  else
58  chr(i:i) = cval(lb)
59  endif
60  enddo
61 
62  return
63 end subroutine upc
64 
79 subroutine upb8(nval,nbits,ibit,ibay)
80 
81  use modv_vars, only: nbitw
82 
83  implicit none
84 
85  integer, intent(in) :: nbits,ibit,ibay(*)
86  integer*8, intent(out) :: nval
87 
88  integer :: nvals(2), jbit, ival
89  integer*8 :: nval8
90 
91  equivalence(nval8,nvals)
92 
93  if(nbits<0) then
94  call bort('BUFRLIB: UPB8 - nbits < zero !!!!!')
95  elseif(nbits<=32) then
96  jbit=ibit; ival=0
97  call upb(ival,nbits,ibay,jbit)
98  nval=ival
99  elseif(nbits<=64) then
100  jbit=ibit; nvals=0
101  call upb(nvals(2),max(nbits-nbitw,0),ibay,jbit)
102  call upb(nvals(1),min(nbitw,nbits ),ibay,jbit)
103  nval=nval8
104  else
105  nval=0
106  endif
107 
108  return
109 end subroutine upb8
110 
127 subroutine up8(nval,nbits,ibay,ibit)
128 
129  implicit none
130 
131  integer, intent(in) :: nbits, ibay(*)
132  integer, intent(inout) :: ibit
133  integer*8, intent(out) :: nval
134 
135  call upb8(nval,nbits,ibit,ibay)
136  ibit = ibit+nbits
137 
138  return
139 end subroutine up8
140 
153 subroutine upbb(nval,nbits,ibit,ibay)
154 
155  use modv_vars, only: nbitw
156 
157  implicit none
158 
159  integer, intent(in) :: ibay(*), ibit, nbits
160  integer, intent(out) :: nval
161  integer nwd, nbt, int, jnt, irev, lbt
162 
163  ! If nbits=0, then just set nval=0 and return
164 
165  if(nbits==0) then
166  nval=0
167  return
168  endif
169 
170  nwd = ibit/nbitw + 1
171  nbt = mod(ibit,nbitw)
172  int = ishft(irev(ibay(nwd)),nbt)
173  int = ishft(int,nbits-nbitw)
174  lbt = nbt+nbits
175  if(lbt>nbitw) then
176  jnt = irev(ibay(nwd+1))
177  int = ior(int,ishft(jnt,lbt-2*nbitw))
178  endif
179  nval = int
180 
181  return
182 end subroutine upbb
183 
201 subroutine upb(nval,nbits,ibay,ibit)
202 
203  implicit none
204 
205  integer, intent(in) :: ibay(*), nbits
206  integer, intent(out) :: nval
207  integer, intent(inout) :: ibit
208 
209  call upbb(nval,nbits,ibit,ibay)
210  ibit = ibit+nbits
211 
212  return
213 end subroutine upb
214 
225 recursive integer function iupb(mbay,nbyt,nbit) result(iret)
226 
227  use modv_vars, only: im8b
228 
229  implicit none
230 
231  integer, intent(in) :: mbay(*), nbit, nbyt
232  integer my_nbit, my_nbyt, mbit
233 
234  ! Check for I8 integers.
235 
236  if(im8b) then
237  im8b=.false.
238 
239  call x84(nbyt,my_nbyt,1)
240  call x84(nbit,my_nbit,1)
241  iret = iupb(mbay,my_nbyt,my_nbit)
242 
243  im8b=.true.
244  return
245  endif
246 
247  mbit = (nbyt-1)*8
248  call upb(iret,nbit,mbay,mbit)
249 
250  return
251 end function iupb
252 
264 recursive integer function iupm(cbay,nbits) result(iret)
265 
266  use modv_vars, only: im8b, nbitw
267 
268  implicit none
269 
270  character*4, intent(in) :: cbay
271  character*4 cint
272  character*128 bort_str
273 
274  integer, intent(in) :: nbits
275  integer my_nbits, int, irev
276 
277  equivalence(cint,int)
278 
279  ! Check for I8 integers.
280 
281  if(im8b) then
282  im8b=.false.
283 
284  call x84(nbits,my_nbits,1)
285  iret = iupm(cbay,my_nbits)
286 
287  im8b=.true.
288  return
289  endif
290 
291  iret = 0
292  if(nbits>nbitw) then
293  write(bort_str,'("BUFRLIB: IUPM - NUMBER OF BITS BEING UNPACKED'// &
294  ', NBITS (",I4,"), IS > THE INTEGER WORD LENGTH ON THIS MACHINE, NBITW (",I3,")")') nbits,nbitw
295  call bort(bort_str)
296  endif
297  cint = cbay
298  int = irev(int)
299  iret = ishft(int,nbits-nbitw)
300 
301  return
302 end function iupm
303 
318 real*8 function ups(ival,node) result(r8ret)
319 
320  use moda_tables
321  use moda_nrv203
322 
323  implicit none
324 
325  integer*8, intent(in) :: ival
326  integer*8 imask
327  integer, intent(in) :: node
328  integer jj
329 
330  real*8, parameter :: ten = 10.
331 
332  r8ret = ( ival + irf(node) ) * ten**(-isc(node))
333 
334  if ( nnrv > 0 ) then
335  ! There are redefined reference values in the jump/link table, so we need to check if this node is affected by any of them.
336  do jj = 1, nnrv
337  if ( node == inodnrv(jj) ) then
338  ! This node contains a redefined reference value. Per the rules of BUFR, negative values may be encoded as positive
339  ! integers with the left-most bit set to 1.
340  imask = 2_8**(ibt(node)-1)
341  if ( iand(ival,imask) > 0 ) then
342  nrv(jj) = (-1) * ( ival - imask )
343  else
344  nrv(jj) = ival
345  end if
346  r8ret = nrv(jj)
347  return
348  else if ( ( tag(node)(1:8) == tagnrv(jj) ) .and. ( node >= isnrv(jj) ) .and. ( node <= ienrv(jj) ) ) then
349  ! The corresponding redefinded reference value needs to be used when decoding this value.
350  r8ret = ( ival + nrv(jj) ) * ten**(-isc(node))
351  return
352  end if
353  end do
354  end if
355 
356  return
357 end function ups
subroutine bort(str)
Log an error message, then abort the application program.
Definition: borts.F90:15
subroutine upb(nval, nbits, ibay, ibit)
Decode an integer value from within a specified number of bits of an integer array,...
Definition: cidecode.F90:202
subroutine upbb(nval, nbits, ibit, ibay)
Decode an integer value from within a specified number of bits of an integer array,...
Definition: cidecode.F90:154
subroutine upb8(nval, nbits, ibit, ibay)
Decode an 8-byte integer value from within a specified number of bits of an integer array,...
Definition: cidecode.F90:80
real *8 function ups(ival, node)
Unpack a real*8 value from an integer by applying the proper scale and reference values.
Definition: cidecode.F90:319
recursive integer function iupb(mbay, nbyt, nbit)
Decode an integer value from within a specified number of bits of an integer array,...
Definition: cidecode.F90:226
recursive integer function iupm(cbay, nbits)
Decode an integer value from within a specified number of bits of a character string,...
Definition: cidecode.F90:265
subroutine up8(nval, nbits, ibay, ibit)
Decode an 8-byte integer value from within a specified number of bits of an integer array,...
Definition: cidecode.F90:128
subroutine upc(chr, nchr, ibay, ibit, cnvnull)
Decode a character string from within a specified number of bytes of an integer array,...
Definition: cidecode.F90:26
integer function irev(n)
Return a copy of an integer value with the bytes possibly reversed.
Definition: misc.F90:257
Declare arrays and variables for use with any 2-03-YYY (change reference value) operators present wit...
integer, dimension(:), allocatable ienrv
End of entry range in jump/link table, within which the corresponding new reference value in nrv will...
character *8, dimension(:), allocatable tagnrv
Table B mnemonic to which the corresponding new reference value in nrv applies.
integer, dimension(:), allocatable isnrv
Start of entry range in jump/link table, within which the corresponding new reference value in nrv wi...
integer nnrv
Number of entries in the jump/link table which contain new reference values (up to a maximum of mxnrv...
integer *8, dimension(:), allocatable nrv
New reference values corresponding to inodnrv.
integer, dimension(:), allocatable inodnrv
Entries within jump/link table which contain new reference values.
Declare arrays and variables used to store the internal jump/link table.
integer, dimension(:), allocatable irf
Reference values corresponding to tag and typ:
integer, dimension(:), allocatable isc
Scale factors corresponding to tag and typ:
integer, dimension(:), allocatable ibt
Bit widths corresponding to tag and typ:
character *10, dimension(:), allocatable tag
Mnemonics in the jump/link table.
subroutine x84(iin8, iout4, nval)
Encode one or more 8-byte integer values as 4-byte integer values.
Definition: x4884.F90:65