NCEPLIBS-bufr  12.1.0
readwritesb.F90
Go to the documentation of this file.
1 
5 
31 recursive subroutine readsb(lunit,iret)
32 
33  use modv_vars, only: im8b
34 
35  use moda_msgcwd
36  use moda_unptyp
37  use moda_bitbuf
38  use moda_bitmaps
39  use moda_stcode
40 
41  implicit none
42 
43  integer, intent(in) :: lunit
44  integer, intent(out) :: iret
45  integer my_lunit, lun, il, im, ier, nbyt
46 
47  ! Check for I8 integers
48 
49  if(im8b) then
50  im8b=.false.
51 
52  call x84(lunit,my_lunit,1)
53  call readsb(my_lunit,iret)
54  call x48(iret,iret,1)
55 
56  im8b=.true.
57  return
58  endif
59 
60  iret = -1
61 
62  ! Check the file status
63 
64  call status(lunit,lun,il,im)
65  if(il==0) call bort('BUFRLIB: READSB - INPUT BUFR FILE IS CLOSED, IT MUST BE OPEN FOR INPUT')
66  if(il>0) call bort('BUFRLIB: READSB - INPUT BUFR FILE IS OPEN FOR OUTPUT, IT MUST BE OPEN FOR INPUT')
67  if(im==0) return
68 
69  ! See if there is another subset in the message
70 
71  if(nsub(lun)==msub(lun)) return
72  nsub(lun) = nsub(lun) + 1
73 
74  ! Read the next subset and reset the pointers
75 
76  nbtm = 0
77  lstnod = 0
78  lstnodct = 0
79  iscodes(lun) = 0
80  linbtm = .false.
81 
82  if(msgunp(lun)==0) then
83  ibit = mbyt(lun)*8
84  call upb(nbyt,16,mbay(1,lun),ibit)
85  call rdtree(lun,ier)
86  if(ier/=0) return
87  mbyt(lun) = mbyt(lun) + nbyt
88  elseif(msgunp(lun)==1) then
89  ! message with "standard" Section 3
90  ibit = mbyt(lun)
91  call rdtree(lun,ier)
92  if(ier/=0) return
93  mbyt(lun) = ibit
94  else
95  ! compressed message
96  call rdcmps(lun)
97  if (iscodes(lun) /= 0) return
98  endif
99 
100  iret = 0
101 
102  return
103 end subroutine readsb
104 
116 recursive integer function ireadsb(lunit) result(iret)
117 
118  use modv_vars, only: im8b
119 
120  implicit none
121 
122  integer, intent(in) :: lunit
123  integer my_lunit
124 
125  ! Check for I8 integers.
126 
127  if(im8b) then
128  im8b=.false.
129 
130  call x84(lunit,my_lunit,1)
131  iret=ireadsb(my_lunit)
132 
133  im8b=.true.
134  return
135  endif
136 
137  call readsb(lunit,iret)
138 
139  return
140 end function ireadsb
141 
167 recursive subroutine readns(lunit,subset,jdate,iret)
168 
169  use modv_vars, only: im8b
170 
171  use moda_msgcwd
172  use moda_tables
173 
174  implicit none
175 
176  integer, intent(in) :: lunit
177  integer, intent(out) :: jdate, iret
178  integer my_lunit, lun, il, im
179 
180  character*8, intent(out) :: subset
181 
182  ! Check for I8 integers
183 
184  if(im8b) then
185  im8b=.false.
186 
187  call x84(lunit,my_lunit,1)
188  call readns(my_lunit,subset,jdate,iret)
189  call x48(jdate,jdate,1)
190  call x48(iret,iret,1)
191 
192  im8b=.true.
193  return
194  endif
195 
196  ! Refresh the subset and jdate parameters
197 
198  call status(lunit,lun,il,im)
199  if(il==0) call bort('BUFRLIB: READNS - INPUT BUFR FILE IS CLOSED, IT MUST BE OPEN FOR INPUT')
200  if(il>0) call bort('BUFRLIB: READNS - INPUT BUFR FILE IS OPEN FOR OUTPUT, IT MUST BE OPEN FOR INPUT')
201  if(inode(lun)==0) then
202  subset = ' '
203  else
204  subset = tag(inode(lun))(1:8)
205  endif
206  jdate = idate(lun)
207 
208  ! Read the next subset in the BUFR file
209 
210  do while (.true.)
211  call readsb(lunit,iret)
212  if (iret==0) exit
213  call readmg(lunit,subset,jdate,iret)
214  if (iret/=0) exit
215  enddo
216 
217  return
218 end subroutine readns
219 
238 recursive integer function ireadns(lunit,subset,idate) result(iret)
239 
240  use modv_vars, only: im8b
241 
242  implicit none
243 
244  integer, intent(in) :: lunit
245  integer, intent(out) :: idate
246  integer my_lunit
247 
248  character*8, intent(out) :: subset
249 
250  ! Check for I8 integers.
251 
252  if(im8b) then
253  im8b=.false.
254 
255  call x84(lunit,my_lunit,1)
256  iret=ireadns(my_lunit,subset,idate)
257  call x48(idate,idate,1)
258 
259  im8b=.true.
260  return
261  endif
262 
263  call readns(lunit,subset,idate,iret)
264 
265  return
266 end function ireadns
267 
302 recursive subroutine writsb(lunit)
303 
304  use modv_vars, only: im8b
305 
306  use moda_msgcmp
307 
308  implicit none
309 
310  integer, intent(in) :: lunit
311  integer my_lunit, lun, il, im
312 
313  ! Check for I8 integers
314 
315  if(im8b) then
316  im8b=.false.
317 
318  call x84 ( lunit, my_lunit, 1 )
319  call writsb ( my_lunit )
320 
321  im8b=.true.
322  return
323  endif
324 
325  ! Check the file status
326 
327  call status(lunit,lun,il,im)
328  if(il==0) call bort('BUFRLIB: WRITSB - OUTPUT BUFR FILE IS CLOSED, IT MUST BE OPEN FOR OUTPUT')
329  if(il<0) call bort('BUFRLIB: WRITSB - OUTPUT BUFR FILE IS OPEN FOR INPUT, IT MUST BE OPEN FOR OUTPUT')
330  if(im==0) call bort('BUFRLIB: WRITSB - A MESSAGE MUST BE OPEN IN OUTPUT BUFR FILE, NONE ARE')
331 
332  ! Pack up the subset and put it into the message
333 
334  call wrtree(lun)
335  if( ccmf=='Y' ) then
336  call wrcmps(lunit)
337  else
338  call msgupd(lunit,lun)
339  endif
340 
341  return
342 end subroutine writsb
343 
422 recursive subroutine writsa(lunxx,lmsgt,msgt,msgl)
423 
424  use modv_vars, only: im8b
425 
426  use moda_bufrmg
427  use moda_msgcmp
428 
429  implicit none
430 
431  integer, intent(in) :: lunxx, lmsgt
432  integer, intent(out) :: msgt(*), msgl
433  integer my_lunxx, my_lmsgt, lunit, lun, il, im, n
434 
435  ! Check for I8 integers
436 
437  if(im8b) then
438  im8b=.false.
439 
440  call x84 ( lunxx, my_lunxx, 1 )
441  call x84 ( lmsgt, my_lmsgt, 1 )
442  call writsa ( my_lunxx, my_lmsgt*2, msgt, msgl )
443  msgl = msgl/2
444  call x48 ( msgl, msgl, 1 )
445 
446  im8b=.true.
447  return
448  endif
449 
450  lunit = abs(lunxx)
451 
452  ! Check the file status
453 
454  call status(lunit,lun,il,im)
455  if(il==0) call bort('BUFRLIB: WRITSA - OUTPUT BUFR FILE IS CLOSED, IT MUST BE OPEN FOR OUTPUT')
456  if(il<0) call bort('BUFRLIB: WRITSA - OUTPUT BUFR FILE IS OPEN FOR INPUT, IT MUST BE OPEN FOR OUTPUT')
457  if(im==0) call bort('BUFRLIB: WRITSA - A MESSAGE MUST BE OPEN IN OUTPUT BUFR FILE, NONE ARE')
458 
459  ! If lunxx < 0, force memory msg to be written (w/o any current subset)
460 
461  if(lunxx<0) call closmg(lunit)
462 
463  ! Is there a completed BUFR message to be returned?
464 
465  if(msglen(lun)>0) then
466  if(msglen(lun)>lmsgt) call bort('BUFRLIB: WRITSA - OVERFLOW OF OUTPUT BUFR MESSAGE ARRAY; TRY A LARGER '// &
467  'DIMENSION FOR THIS ARRAY')
468  msgl = msglen(lun)
469  do n=1,msgl
470  msgt(n) = msgtxt(n,lun)
471  enddo
472  msglen(lun) = 0
473  else
474  msgl = 0
475  endif
476 
477  if(lunxx<0) return
478 
479  ! Pack up the subset and put it into the message
480 
481  call wrtree(lun)
482  if( ccmf=='Y' ) then
483  call wrcmps(lunit)
484  else
485  call msgupd(lunit,lun)
486  endif
487 
488  ! If the just-completed call to wrcmps() or msgupd() for this subset caused a message to be flushed to abs(lunxx), then
489  ! attempt to retrieve and return that message now. Otherwise, we run the risk that the next call to openmb() or openmg()
490  ! might cause another message to be flushed, and thus overwrite the current message within array msgtxt before we
491  ! had the chance to retrieve it during the next call to writsa().
492 
493  ! Also note that, in rare instances (e.g. if the byte count of the most recent subset is > 65530), we could end up with
494  ! two BUFR messages available to be returned from this one call to writsa(). If sufficient space is available in the
495  ! msgt array, then go ahead and return both messages now.
496 
497  if( (msglen(lun)>0) .and. (msgl+msglen(lun)<=lmsgt) ) then
498  do n = 1,msglen(lun)
499  msgt(msgl+n) = msgtxt(n,lun)
500  enddo
501  msgl = msgl+msglen(lun)
502  msglen(lun) = 0
503  endif
504 
505  return
506 end subroutine writsa
507 
535 recursive subroutine rdmgsb(lunit,imsg,isub)
536 
537  use modv_vars, only: im8b
538 
539  use moda_msgcwd
540  use moda_bitbuf
541 
542  implicit none
543 
544  integer, intent(in) :: lunit, imsg, isub
545  integer my_lunit, my_imsg, my_isub, lun, il, im, i, jdate, iret
546 
547  character*128 bort_str
548  character*8 subset
549 
550  ! Check for I8 integers
551 
552  if(im8b) then
553  im8b=.false.
554 
555  call x84(lunit,my_lunit,1)
556  call x84(imsg,my_imsg,1)
557  call x84(isub,my_isub,1)
558  call rdmgsb(my_lunit,my_imsg,my_isub)
559 
560  im8b=.true.
561  return
562  endif
563 
564  ! Open the file and skip to message #imsg
565 
566  call openbf(lunit,'IN',lunit)
567  call status(lunit,lun,il,im)
568 
569  ! Note that we need to use subroutine readmg() to actually read in all of the messages (including the
570  ! first (imsg-1) messages!), just in case there are any embedded dictionary messages in the file.
571 
572  do i=1,imsg
573  call readmg(lunit,subset,jdate,iret)
574  if(iret<0) then
575  write(bort_str,'("BUFRLIB: RDMGSB - HIT END OF FILE BEFORE READING REQUESTED MESSAGE NO.",I5," IN '//&
576  'BUFR FILE CONNECTED TO UNIT",I4)') imsg,lunit
577  call bort(bort_str)
578  endif
579  enddo
580 
581  ! Position at subset #isub
582 
583  do i=1,isub
584  call readsb(lunit,iret)
585  if(iret<0) then
586  write(bort_str,'("BUFRLIB: RDMGSB - ALL SUBSETS READ BEFORE READING REQ. SUBSET NO.",I3," IN '// &
587  'REQ. MSG NO.",I5," IN BUFR FILE CONNECTED TO UNIT",I4)') isub,imsg,lunit
588  call bort(bort_str)
589  endif
590  enddo
591 
592  return
593 end subroutine rdmgsb
594 
612 subroutine msgupd(lunit,lun)
613 
614  use moda_msgcwd
615  use moda_bitbuf
616  use moda_h4wlc
617 
618  implicit none
619 
620  integer, intent(in) :: lunit, lun
621  integer nby0, nby1, nby2, nby3, nby4, nby5, iprt, ibyt, lbyt, lbit, nbyt, ii, iupb
622 
623  logical msgfull
624 
625  character*128 errstr
626 
627  common /msgptr/ nby0, nby1, nby2, nby3, nby4, nby5
628  common /quiet/ iprt
629 
630  ! Pad the subset buffer
631 
632  call pad(ibay,ibit,ibyt,8)
633 
634  ! Check whether the new subset should be written into the currently open message
635 
636  if(msgfull(mbyt(lun),ibyt,maxbyt) .or. ((ibyt>65530).and.(nsub(lun)>0))) then
637  ! No it should not, either because it doesn't fit
638  ! OR
639  ! It has byte count > 65530 (sufficiently close to the upper limit for the 16 bit byte counter placed at the beginning
640  ! of each subset), and the current message has at least one subset in it
641  !
642  ! In either of these cases, we need to write out the current message and then create a new one to hold the current subset
643  call msgwrt(lunit,mbay(1,lun),mbyt(lun))
644  call msgini(lun)
645  endif
646 
647  if(msgfull(mbyt(lun),ibyt,maxbyt)) then
648  ! This is an overlarge subset that won't fit in any message given the current value of maxbyt, so discard the subset
649  ! and exit gracefully.
650  if(iprt>=0) then
651  call errwrt('+++++++++++++++++++++WARNING+++++++++++++++++++++++')
652  write ( unit=errstr, fmt='(A,A,I7,A)') 'BUFRLIB: MSGUPD - SUBSET LONGER THAN ANY POSSIBLE MESSAGE ', &
653  '{MAXIMUM MESSAGE LENGTH = ', maxbyt, '}'
654  call errwrt(errstr)
655  call errwrt('>>>>>>>OVERLARGE SUBSET DISCARDED FROM FILE<<<<<<<<')
656  call errwrt('+++++++++++++++++++++WARNING+++++++++++++++++++++++')
657  call errwrt(' ')
658  endif
659  call usrtpl(lun,1,1)
660  return
661  endif
662 
663  ! Set a byte count and transfer the subset buffer into the message
664 
665  lbit = 0
666  call pkb(ibyt,16,ibay,lbit)
667 
668  ! Note that we want to append the data for this subset to the end of Section 4, but the value in mbyt(lun) already includes
669  ! the length of Section 5 (i.e. 4 bytes). Therefore, we need to begin writing at the point 3 bytes prior to the byte
670  ! currently pointed to by mbyt(lun).
671 
672  call mvb(ibay,1,mbay(1,lun),mbyt(lun)-3,ibyt)
673 
674  ! Update the subset and byte counters
675 
676  mbyt(lun) = mbyt(lun) + ibyt
677  nsub(lun) = nsub(lun) + 1
678 
679  lbit = (nby0+nby1+nby2+4)*8
680  call pkb(nsub(lun),16,mbay(1,lun),lbit)
681 
682  lbyt = nby0+nby1+nby2+nby3
683  nbyt = iupb(mbay(1,lun),lbyt+1,24)
684  lbit = lbyt*8
685  call pkb(nbyt+ibyt,24,mbay(1,lun),lbit)
686 
687  ! If any long character strings are being held internally for storage into this subset, store them now
688 
689  if(nh4wlc>0) then
690  do ii = 1, nh4wlc
691  call writlc(luh4wlc(ii),chh4wlc(ii),sth4wlc(ii))
692  enddo
693  nh4wlc = 0
694  endif
695 
696  ! If the subset byte count is > 65530, then give it its own one-subset message (cannot have any other subsets in this
697  ! message because their beginning would be beyond the upper limit of 65535 in the 16-bit byte counter, meaning they
698  ! could not be located!)
699 
700  if(ibyt>65530) then
701  if(iprt>=1) then
702  call errwrt('+++++++++++++++++++++WARNING+++++++++++++++++++++++')
703  write ( unit=errstr, fmt='(A,I7,A,A)') 'BUFRLIB: MSGUPD - SUBSET HAS BYTE COUNT = ',ibyt,' > UPPER LIMIT OF 65535'
704  call errwrt(errstr)
705  call errwrt('>>>>>>>WILL BE WRITTEN INTO ITS OWN MESSAGE<<<<<<<<')
706  call errwrt('+++++++++++++++++++++WARNING+++++++++++++++++++++++')
707  call errwrt(' ')
708  endif
709  call msgwrt(lunit,mbay(1,lun),mbyt(lun))
710  call msgini(lun)
711  endif
712 
713  ! Reset the user arrays
714 
715  call usrtpl(lun,1,1)
716 
717  return
718 end subroutine msgupd
719 
749 subroutine pad(ibay,ibit,ibyt,ipadb)
750 
751  implicit none
752 
753  integer, intent(inout) :: ibay(*), ibit
754  integer, intent(in) :: ipadb
755  integer, intent(out) :: ibyt
756  integer ipad
757 
758  character*128 bort_str
759 
760  ! Pad the subset to an ipadb bit boundary
761 
762  ipad = ipadb - mod(ibit+8,ipadb)
763  ! First pack the # of bits being padded (this is a delayed replication factor)
764  call pkb(ipad,8,ibay,ibit)
765  ! Now pad with zeroes to the byte boundary
766  call pkb(0,ipad,ibay,ibit)
767  ibyt = ibit/8
768 
769  if(mod(ibit,8)/=0) then
770  write(bort_str,'("BUFRLIB: PAD - THE NUMBER OF BITS IN A PACKED'// &
771  ' SUBSET AFTER PADDING (",I8,") IS NOT A MULTIPLE OF 8")') ibit
772  call bort(bort_str)
773  endif
774 
775  return
776 end subroutine pad
777 
803 recursive integer function lcmgdf(lunit,subset) result(iret)
804 
805  use modv_vars, only: im8b
806 
807  use moda_tables
808 
809  implicit none
810 
811  integer, intent(in) :: lunit
812  integer my_lunit, lun, il, im, mtyp, msbt, inod, nte, i
813 
814  character*8, intent(in) :: subset
815 
816  ! Check for I8 integers.
817 
818  if(im8b) then
819  im8b=.false.
820 
821  call x84(lunit,my_lunit,1)
822  iret=lcmgdf(my_lunit,subset)
823 
824  im8b=.true.
825  return
826  endif
827 
828  iret = 0
829 
830  ! Get lun from lunit.
831 
832  call status(lunit,lun,il,im)
833  if (il==0) call bort('BUFRLIB: LCMGDF - INPUT BUFR FILE IS CLOSED, IT MUST BE OPEN')
834 
835  ! Confirm that subset is defined for this logical unit.
836 
837  call nemtba(lun,subset,mtyp,msbt,inod)
838 
839  ! Check if there's a long character string in the definition.
840 
841  nte = isc(inod)-inod
842 
843  do i = 1, nte
844  if ( (typ(inod+i)=='CHR') .and. (ibt(inod+i)>64) ) then
845  iret = 1
846  return
847  endif
848  enddo
849 
850  iret = 0
851 
852  return
853 end function lcmgdf
854 
879 recursive subroutine ufbpos(lunit,irec,isub,subset,jdate)
880 
881  use bufrlib
882 
883  use modv_vars, only: im8b
884 
885  use moda_msgcwd
886  use moda_bitbuf
887 
888  implicit none
889 
890  integer, intent(in) :: lunit, irec, isub
891  integer, intent(out) :: jdate
892  integer my_lunit, my_irec, my_isub, lun, il, im, jrec, jsub, iret
893 
894  character*128 bort_str
895  character*8, intent(out) :: subset
896 
897  ! Check for I8 integers
898 
899  if(im8b) then
900  im8b=.false.
901  call x84(lunit,my_lunit,1)
902  call x84(irec,my_irec,1)
903  call x84(isub,my_isub,1)
904  call ufbpos(my_lunit,my_irec,my_isub,subset,jdate)
905  call x48(jdate,jdate,1)
906  im8b=.true.
907  return
908  endif
909 
910  ! Make sure a file is open for input
911 
912  call status(lunit,lun,il,im)
913  if(il==0) call bort('BUFRLIB: UFBPOS - INPUT BUFR FILE IS CLOSED, IT MUST BE OPEN FOR INPUT')
914  if(il>0) call bort('BUFRLIB: UFBPOS - INPUT BUFR FILE IS OPEN FOR OUTPUT, IT MUST BE OPEN FOR INPUT')
915 
916  if(irec<=0) then
917  write(bort_str,'("BUFRLIB: UFBPOS - REQUESTED MESSAGE NUMBER TO READ IN (",I5,") IS NOT VALID")') irec
918  call bort(bort_str)
919  endif
920  if(isub<=0) then
921  write(bort_str,'("BUFRLIB: UFBPOS - REQUESTED SUBSET NUMBER TO READ IN (",I5,") IS NOT VALID")') isub
922  call bort(bort_str)
923  endif
924 
925  ! See where pointers are currently located
926 
927  call ufbcnt(lunit,jrec,jsub)
928 
929  ! Rewind file if requested pointers are behind current pointers
930 
931  if(irec<jrec .or. (irec==jrec.and.isub<jsub)) then
932  call cewind_c(lun)
933  nmsg(lun) = 0
934  nsub(lun) = 0
935  call ufbcnt(lunit,jrec,jsub)
936  endif
937 
938  ! Read subset #isub from message #irec from file
939 
940  do while (irec>jrec)
941  call readmg(lunit,subset,jdate,iret)
942  if(iret<0) then
943  write(bort_str,'("BUFRLIB: UFBPOS - REQUESTED MESSAGE NUMBER '// &
944  'TO READ IN (",I5,") EXCEEDS THE NUMBER OF MESSAGES IN THE FILE (",I5,")")') irec, jrec
945  call bort(bort_str)
946  endif
947  call ufbcnt(lunit,jrec,jsub)
948  enddo
949 
950  do while (isub>jsub)
951  call readsb(lunit,iret)
952  if(iret/=0) then
953  write(bort_str,'("BUFRLIB: UFBPOS - REQ. SUBSET NUMBER TO READ'// &
954  ' IN (",I5,") EXCEEDS THE NUMBER OF SUBSETS (",I5,") IN THE REQ. MESSAGE (",I5,")")') isub, jsub, irec
955  call bort(bort_str)
956  endif
957  call ufbcnt(lunit,jrec,jsub)
958  enddo
959 
960  return
961 end subroutine ufbpos
962 
974 subroutine rdtree(lun,iret)
975 
976  use modv_vars, only: bmiss
977 
978  use moda_usrint
979  use moda_usrbit
980  use moda_ival
981  use moda_bitbuf
982  use moda_tables
983 
984  implicit none
985 
986  integer, intent(in) :: lun
987  integer, intent(out) :: iret
988  integer ier, n, node, kbit, nbt, icbfms
989 
990  character*8 cval
991 
992  real*8 rval, ups
993 
994  equivalence(cval,rval)
995 
996  iret = 0
997 
998  ! Cycle through a subset setting up the template
999 
1000  mbit(1) = ibit
1001  nbit(1) = 0
1002  call rcstpl(lun,ier)
1003  if(ier/=0) then
1004  iret = -1
1005  return
1006  endif
1007 
1008  ! Unpack a subset into the user array ival
1009 
1010  do n=1,nval(lun)
1011  call upb8(ival(n),nbit(n),mbit(n),mbay(1,lun))
1012  enddo
1013 
1014  ! Loop through each element of the subset, converting the unpacked values to the proper types
1015 
1016  do n=1,nval(lun)
1017  node = inv(n,lun)
1018  if(itp(node)==1) then
1019 
1020  ! The unpacked value is a delayed descriptor replication factor.
1021 
1022  val(n,lun) = ival(n)
1023  elseif(itp(node)==2) then
1024 
1025  ! The unpacked value is a real.
1026 
1027  if (ival(n)<2_8**ibt(node)-1) then
1028  val(n,lun) = ups(ival(n),node)
1029  else
1030  val(n,lun) = bmiss
1031  endif
1032  elseif(itp(node)==3) then
1033 
1034  ! The value is a character string, so unpack it using an equivalenced real*8 value. Note that a maximum of 8 characters
1035  ! will be unpacked here, so a separate subsequent call to subroutine readlc() will be needed to fully unpack any string
1036  ! longer than 8 characters.
1037 
1038  cval = ' '
1039  kbit = mbit(n)
1040  nbt = min(8,nbit(n)/8)
1041  call upc(cval,nbt,mbay(1,lun),kbit,.true.)
1042  if (nbit(n)<=64 .and. icbfms(cval,nbt)/=0) then
1043  val(n,lun) = bmiss
1044  else
1045  val(n,lun) = rval
1046  endif
1047  endif
1048  enddo
1049 
1050  ibit = nbit(nval(lun))+mbit(nval(lun))
1051 
1052  return
1053 end subroutine rdtree
1054 
1063 subroutine wrtree(lun)
1064 
1065  use moda_usrint
1066  use moda_ival
1067  use moda_ufbcpl
1068  use moda_bitbuf
1069  use moda_tables
1070 
1071  implicit none
1072 
1073  integer, intent(in) :: lun
1074  integer*8 ipks
1075  integer n, node, ncr, numchr, jj, ibfms
1076 
1077  character*120 lstr
1078  character*8 cval
1079 
1080  real*8 rval
1081 
1082  equivalence(cval,rval)
1083 
1084  ! Convert user numbers into scaled integers
1085 
1086  do n=1,nval(lun)
1087  node = inv(n,lun)
1088  if(itp(node)==1) then
1089  ival(n) = nint(val(n,lun))
1090  elseif(typ(node)=='NUM') then
1091  if( (ibfms(val(n,lun))==1) .or. (val(n,lun)/=val(n,lun)) ) then
1092  ! The user number is either "missing" or NaN.
1093  ival(n) = -1
1094  else
1095  ival(n) = ipks(val(n,lun),node)
1096  endif
1097  endif
1098  enddo
1099 
1100  ! Pack the user array into the subset buffer
1101 
1102  ibit = 16
1103 
1104  do n=1,nval(lun)
1105  node = inv(n,lun)
1106  if(itp(node)<3) then
1107  ! The value to be packed is numeric.
1108  call pkb8(ival(n),ibt(node),ibay,ibit)
1109  else
1110  ! The value to be packed is a character string.
1111  ncr=ibt(node)/8
1112  if ( ncr>8 .and. luncpy(lun)/=0 ) then
1113  ! The string is longer than 8 characters and there was a preceeding call to ufbcpy() involving this output unit,
1114  ! so read the long string with readlc() and then write it into the output buffer using pkc().
1115  call readlc(luncpy(lun),lstr,tag(node))
1116  call pkc(lstr,ncr,ibay,ibit)
1117  else
1118  rval = val(n,lun)
1119  if(ibfms(rval)/=0) then
1120  ! The value is "missing", so set all bits to 1 before packing the field as a character string.
1121  numchr = min(ncr,len(lstr))
1122  do jj = 1, numchr
1123  call ipkm(lstr(jj:jj),1,255)
1124  enddo
1125  call pkc(lstr,numchr,ibay,ibit)
1126  else
1127  ! The value is not "missing", so pack the equivalenced character string. Note that a maximum of 8 characters
1128  ! will be packed here, so a separate subsequent call to subroutine writlc() will be needed to fully encode any
1129  ! string longer than 8 characters.
1130  call pkc(cval,ncr,ibay,ibit)
1131  endif
1132  endif
1133  endif
1134  enddo
1135 
1136  ! Reset ufbcpy() file pointer
1137 
1138  luncpy(lun)=0
1139 
1140  return
1141 end subroutine wrtree
1142 
1155 subroutine rcstpl(lun,iret)
1156 
1157  use modv_vars, only: bmiss, maxjl, maxss, maxrcr
1158 
1159  use moda_usrint
1160  use moda_usrbit
1161  use moda_msgcwd
1162  use moda_bitbuf
1163  use moda_tables
1164  use moda_usrtmp
1165 
1166  implicit none
1167 
1168  character*128 bort_str
1169 
1170  integer, intent(in) :: lun
1171  integer, intent(out) :: iret
1172  integer nbmp(2,maxrcr), newn(2,maxrcr), knx(maxrcr), iprt, nodi, node, mbmp, knvn, nr, i, j, n, nn, n1, n2, new, &
1173  idpri, igetrfel
1174 
1175  common /quiet/ iprt
1176 
1177  iret = 0
1178 
1179  ! Set the initial values for the template
1180 
1181  inv(1,lun) = inode(lun)
1182  val(1,lun) = 0
1183  nbmp(1,1) = 1
1184  nbmp(2,1) = 1
1185  nodi = inode(lun)
1186  node = inode(lun)
1187  mbmp = 1
1188  knvn = 1
1189  nr = 0
1190 
1191  do i=1,maxrcr
1192  knx(i) = 0
1193  enddo
1194 
1195  outer: do while (.true.)
1196 
1197  ! Set up the parameters for a level of recursion
1198 
1199  nr = nr+1
1200  if(nr>maxrcr) then
1201  write(bort_str,'("BUFRLIB: RCSTPL - THE NUMBER OF RECURSION LEVELS EXCEEDS THE LIMIT (",I3,")")') maxrcr
1202  call bort(bort_str)
1203  endif
1204  nbmp(1,nr) = 1
1205  nbmp(2,nr) = mbmp
1206 
1207  n1 = iseq(node,1)
1208  n2 = iseq(node,2)
1209  if(n1==0) then
1210  write(bort_str,'("BUFRLIB: RCSTPL - UNSET EXPANSION SEGMENT ",A)') tag(nodi)
1211  call bort(bort_str)
1212  endif
1213  if(n2-n1+1>maxjl) then
1214  if(iprt>=0) then
1215  call errwrt('++++++++++++++BUFR ARCHIVE LIBRARY+++++++++++++++++')
1216  call errwrt('BUFRLIB: RCSTPL - MAXJL OVERFLOW; SUBSET SKIPPED')
1217  call errwrt('++++++++++++++BUFR ARCHIVE LIBRARY+++++++++++++++++')
1218  endif
1219  iret = -1
1220  return
1221  endif
1222  newn(1,nr) = 1
1223  newn(2,nr) = n2-n1+1
1224 
1225  do n=1,newn(2,nr)
1226  nn = jseq(n+n1-1)
1227  iutmp(n,nr) = nn
1228  vutmp(n,nr) = vali(nn)
1229  enddo
1230 
1231  do while (.true.)
1232 
1233  ! Store nodes at some recursion level
1234 
1235  do i=nbmp(1,nr),nbmp(2,nr)
1236  if(knx(nr)==0) knx(nr) = knvn
1237  if(i>nbmp(1,nr)) newn(1,nr) = 1
1238  do j=newn(1,nr),newn(2,nr)
1239  if(knvn+1>maxss) then
1240  if(iprt>=0) then
1241  call errwrt('++++++++++++++BUFR ARCHIVE LIBRARY+++++++++++++++++')
1242  call errwrt('BUFRLIB: RCSTPL - MAXSS OVERFLOW; SUBSET SKIPPED')
1243  call errwrt('++++++++++++++BUFR ARCHIVE LIBRARY+++++++++++++++++')
1244  endif
1245  iret = -1
1246  return
1247  endif
1248  knvn = knvn+1
1249  node = iutmp(j,nr)
1250  ! inv is positional index in internal jump/link table for packed subset element knvn in mbay
1251  inv(knvn,lun) = node
1252  ! mbit is the bit in mbay pointing to where the packed subset element knvn begins
1253  mbit(knvn) = mbit(knvn-1)+nbit(knvn-1)
1254  ! nbit is the number of bits in mbay occupied by packed subset element knvn
1255  nrfelm(knvn,lun) = igetrfel(knvn,lun)
1256  nbit(knvn) = ibt(node)
1257  if(tag(node)(1:5)=='DPRI ') then
1258  ! This is a bitmap entry, so get and store the corresponding value
1259  call upbb(idpri,nbit(knvn),mbit(knvn),mbay(1,lun))
1260  if(idpri==0) then
1261  val(knvn,lun) = 0.0
1262  else
1263  val(knvn,lun) = bmiss
1264  endif
1265  call strbtm(knvn,lun)
1266  endif
1267  ! Actual unpacked subset values are initialized here
1268  val(knvn,lun) = vutmp(j,nr)
1269  if(itp(node)==1) then
1270  call upbb(mbmp,nbit(knvn),mbit(knvn),mbay(1,lun))
1271  newn(1,nr) = j+1
1272  nbmp(1,nr) = i
1273  cycle outer
1274  endif
1275  enddo
1276  new = knvn-knx(nr)
1277  val(knx(nr)+1,lun) = val(knx(nr)+1,lun) + new
1278  knx(nr) = 0
1279  enddo
1280 
1281  ! Check if we need to continue one recursion level back
1282 
1283  if(nr-1 == 0) exit outer
1284  nr = nr-1
1285  enddo
1286 
1287  enddo outer
1288 
1289  ! Finally store the length of (i.e. the number of elements in) the subset template
1290 
1291  nval(lun) = knvn
1292 
1293  return
1294 end subroutine rcstpl
1295 
1306 subroutine usrtpl(lun,invn,nbmp)
1307 
1308  use modv_vars, only: maxjl, maxss
1309 
1310  use moda_usrint
1311  use moda_msgcwd
1312  use moda_tables
1313  use moda_ivttmp
1314  use moda_stcode
1315 
1316  implicit none
1317 
1318  integer, intent(in) :: lun, invn, nbmp
1319  integer iprt, i, j, ival, jval, n, n1, n2, nodi, node, newn, invr, knvn
1320 
1321  character*128 bort_str, errstr
1322 
1323  logical drp, drs, drb, drx
1324 
1325  common /quiet/ iprt
1326 
1327  if(iprt>=2) then
1328  call errwrt('++++++++++++++BUFR ARCHIVE LIBRARY+++++++++++++++++')
1329  write ( unit=errstr, fmt='(A,I3,A,I7,A,I5,A,A10)' ) &
1330  'BUFRLIB: USRTPL - LUN:INVN:NBMP:TAG(INODE(LUN)) = ', lun, ':', invn, ':', nbmp, ':', tag(inode(lun))
1331  call errwrt(errstr)
1332  call errwrt('++++++++++++++BUFR ARCHIVE LIBRARY+++++++++++++++++')
1333  call errwrt(' ')
1334  endif
1335 
1336  if(nbmp<=0) then
1337  if(iprt>=1) then
1338  call errwrt('+++++++++++++++++++++WARNING+++++++++++++++++++++++')
1339  call errwrt(.LE.'BUFRLIB: USRTPL - NBMP 0 - IMMEDIATE RETURN')
1340  call errwrt('+++++++++++++++++++++WARNING+++++++++++++++++++++++')
1341  call errwrt(' ')
1342  endif
1343  return
1344  endif
1345 
1346  drp = .false.
1347  drs = .false.
1348  drx = .false.
1349 
1350  ! Set up a node expansion
1351 
1352  if(invn==1) then
1353  ! The node is a Table A mnemonic
1354  nodi = inode(lun)
1355  inv(1,lun) = nodi
1356  nval(lun) = 1
1357  if(nbmp/=1) then
1358  write(bort_str,'("BUFRLIB: USRTPL - THIRD ARGUMENT (INPUT) = ",'// &
1359  'I4,", MUST BE 1 WHEN SECOND ARGUMENT (INPUT) IS 1 (SUBSET NODE) (",A,")")') nbmp, tag(nodi)
1360  call bort(bort_str)
1361  endif
1362  elseif(invn>0 .and. invn<=nval(lun)) then
1363  ! The node is (hopefully) a delayed replication factor
1364  nodi = inv(invn,lun)
1365  drp = typ(nodi) == 'DRP'
1366  drs = typ(nodi) == 'DRS'
1367  drb = typ(nodi) == 'DRB'
1368  drx = drp .or. drs .or. drb
1369  ival = nint(val(invn,lun))
1370  jval = 2**ibt(nodi)-1
1371  val(invn,lun) = ival+nbmp
1372  if(drb.and.nbmp/=1) then
1373  write(bort_str,'("BUFRLIB: USRTPL - THIRD ARGUMENT (INPUT) = ",'// &
1374  'I4,", MUST BE 1 WHEN NODE IS DRB (1-BIT DELAYED REPL. FACTOR) (",A,")")') nbmp, tag(nodi)
1375  call bort(bort_str)
1376  endif
1377  if(.not.drx) then
1378  write(bort_str,'("BUFRLIB: USRTPL - NODE IS OF TYPE ",A," - IT '// &
1379  'MUST BE EITHER A SUBSET OR DELAYED REPL. FACTOR (",A,")")') typ(nodi), tag(nodi)
1380  call bort(bort_str)
1381  endif
1382  if(ival<0) then
1383  write(bort_str,'("BUFRLIB: USRTPL - REPLICATION FACTOR IS NEGATIVE (=",I5,") (",A,")")') ival, tag(nodi)
1384  call bort(bort_str)
1385  endif
1386  if(ival+nbmp>jval) then
1387  write(bort_str,'("BUFRLIB: USRTPL - REPLICATION FACTOR OVERFLOW (EXCEEDS MAXIMUM OF",I6," (",A,")")') jval, tag(nodi)
1388  call errwrt(bort_str)
1389  iscodes(lun) = 1
1390  return
1391  endif
1392  else
1393  write(bort_str,'("BUFRLIB: USRTPL - INVENTORY INDEX {FIRST '// &
1394  'ARGUMENT (INPUT)} OUT OF BOUNDS (=",I5,", RANGE IS 1 TO",I6,") ")') invn, nval(lun)
1395  call bort(bort_str)
1396  endif
1397 
1398  ! Recall a pre-fab node expansion segment
1399 
1400  newn = 0
1401  n1 = iseq(nodi,1)
1402  n2 = iseq(nodi,2)
1403 
1404  if(n1==0) then
1405  write(bort_str,'("BUFRLIB: USRTPL - UNSET EXPANSION SEGMENT (",A,")")') tag(nodi)
1406  call bort(bort_str)
1407  endif
1408  if(n2-n1+1>maxjl) then
1409  write(bort_str,'("BUFRLIB: USRTPL - TEMPLATE ARRAY OVERFLOW, EXCEEDS THE LIMIT (",I6,") (",A,")")') maxjl, tag(nodi)
1410  call bort(bort_str)
1411  endif
1412 
1413  do n=n1,n2
1414  newn = newn+1
1415  itmp(newn) = jseq(n)
1416  vtmp(newn) = vali(jseq(n))
1417  enddo
1418 
1419  ! Move old nodes and store new ones
1420 
1421  if(nval(lun)+newn*nbmp>maxss) then
1422  write(bort_str,'("BUFRLIB: USRTPL - INVENTORY OVERFLOW (",I6,"), EXCEEDS THE LIMIT (",I6,") (",A,")")') &
1423  nval(lun)+newn*nbmp, maxss, tag(nodi)
1424  call bort(bort_str)
1425  endif
1426 
1427  do j=nval(lun),invn+1,-1
1428  inv(j+newn*nbmp,lun) = inv(j,lun)
1429  val(j+newn*nbmp,lun) = val(j,lun)
1430  enddo
1431 
1432  if(drp.or.drs) vtmp(1) = newn
1433  knvn = invn
1434 
1435  do i=1,nbmp
1436  do j=1,newn
1437  knvn = knvn+1
1438  inv(knvn,lun) = itmp(j)
1439  val(knvn,lun) = vtmp(j)
1440  enddo
1441  enddo
1442 
1443  ! Reset pointers and counters
1444 
1445  nval(lun) = nval(lun) + newn*nbmp
1446 
1447  if(iprt>=2) then
1448  call errwrt('++++++++++++++BUFR ARCHIVE LIBRARY+++++++++++++++++')
1449  write ( unit=errstr, fmt='(A,A,A10,2(A,I5),A,I7)' ) 'BUFRLIB: USRTPL - TAG(INV(INVN,LUN)):NEWN:NBMP:', &
1450  'NVAL(LUN) = ', tag(inv(invn,lun)), ':', newn, ':', nbmp, ':', nval(lun)
1451  call errwrt(errstr)
1452  do i=1,newn
1453  write ( unit=errstr, fmt='(2(A,I5),A,A10)' ) 'For I = ', i, ', ITMP(I) = ', itmp(i), ', TAG(ITMP(I)) = ', tag(itmp(i))
1454  call errwrt(errstr)
1455  enddo
1456  call errwrt('++++++++++++++BUFR ARCHIVE LIBRARY+++++++++++++++++')
1457  call errwrt(' ')
1458  endif
1459 
1460  if(drx) then
1461  node = nodi
1462  invr = invn
1463  outer: do while (.true.)
1464  node = jmpb(node)
1465  if(node<=0) exit
1466  if(itp(node)==0) then
1467  do invr=invr-1,1,-1
1468  if(inv(invr,lun)==node) then
1469  val(invr,lun) = val(invr,lun)+newn*nbmp
1470  cycle outer
1471  endif
1472  enddo
1473  write(bort_str,'("BUFRLIB: USRTPL - BAD BACKUP STRATEGY (",A,")")') tag(nodi)
1474  call bort(bort_str)
1475  else
1476  cycle
1477  endif
1478  enddo outer
1479  endif
1480 
1481  return
1482 end subroutine usrtpl
1483 
1497 recursive subroutine invmrg(lubfi,lubfj)
1498 
1499  use modv_vars, only: im8b
1500 
1501  use moda_usrint
1502  use moda_tables
1503 
1504  implicit none
1505 
1506  integer, intent(in) :: lubfi, lubfj
1507  integer nrpl, nmrg, namb, ntot, my_lubfi, my_lubfj, luni, il, im, lunj, jl, jm, is, js, node, nodj, ityp, iwrds, jwrds, &
1508  n, ioff, nwords, ibfms
1509 
1510  character*128 bort_str
1511 
1512  logical herei, herej, missi, missj, samei
1513 
1514  common /mrgcom/ nrpl, nmrg, namb, ntot
1515 
1516  ! Check for I8 integers
1517 
1518  if(im8b) then
1519  im8b=.false.
1520  call x84(lubfi,my_lubfi,1)
1521  call x84(lubfj,my_lubfj,1)
1522  call invmrg(my_lubfi,my_lubfj)
1523  im8b=.true.
1524  return
1525  endif
1526 
1527  is = 1
1528  js = 1
1529 
1530  ! Get the unit pointers
1531 
1532  call status(lubfi,luni,il,im)
1533  call status(lubfj,lunj,jl,jm)
1534 
1535  ! Step through the buffers comparing the inventory and merging data
1536 
1537  do while(is<=nval(luni))
1538  ! Confirm we're at the same node in each buffer
1539  node = inv(is,luni)
1540  nodj = inv(js,lunj)
1541  if(node/=nodj) then
1542  write(bort_str,'("BUFRLIB: INVMRG - NODE FROM INPUT BUFR FILE '// &
1543  '(",I7,") DOES NOT EQUAL NODE FROM OUTPUT BUFR FILE (",I7,"), TABULAR MISMATCH")') node, nodj
1544  call bort(bort_str)
1545  endif
1546 
1547  ityp = itp(node)
1548  if(ityp==1) then
1549  ! Do an entire sequence replacement
1550  if(typ(node)=='DRB') then
1551  ioff = 0
1552  else
1553  ioff = 1
1554  endif
1555  iwrds = nwords(is,luni)+ioff
1556  jwrds = nwords(js,lunj)+ioff
1557  if(iwrds>ioff .and. jwrds==ioff) then
1558  do n=nval(lunj),js+1,-1
1559  inv(n+iwrds-jwrds,lunj) = inv(n,lunj)
1560  val(n+iwrds-jwrds,lunj) = val(n,lunj)
1561  enddo
1562  do n=0,iwrds
1563  inv(js+n,lunj) = inv(is+n,luni)
1564  val(js+n,lunj) = val(is+n,luni)
1565  enddo
1566  nval(lunj) = nval(lunj)+iwrds-jwrds
1567  jwrds = iwrds
1568  nrpl = nrpl+1
1569  endif
1570  is = is+iwrds
1571  js = js+jwrds
1572  elseif((ityp==2).or.(ityp==3)) then
1573  ! Fill missing values
1574  herei = ibfms(val(is,luni))==0
1575  herej = ibfms(val(js,lunj))==0
1576  missi = .not.(herei)
1577  missj = .not.(herej)
1578  samei = val(is,luni)==val(js,lunj)
1579  if(herei.and.missj) then
1580  val(js,lunj) = val(is,luni)
1581  nmrg = nmrg+1
1582  elseif(herei.and.herej.and..not.samei) then
1583  namb = namb+1
1584  endif
1585  endif
1586 
1587  ! Bump the counters and go check the next pair
1588  is = is + 1
1589  js = js + 1
1590  enddo
1591 
1592  ntot = ntot+1
1593 
1594  return
1595 end subroutine invmrg
1596 
1605 integer function nwords(n,lun) result(iret)
1606 
1607  use moda_usrint
1608 
1609  implicit none
1610 
1611  integer, intent(in) :: n, lun
1612  integer k
1613 
1614  iret = 0
1615 
1616  do k=1,nint(val(n,lun))
1617  iret = iret + nint(val(iret+n+1,lun))
1618  enddo
1619 
1620  return
1621 end function nwords
subroutine strbtm(n, lun)
Store internal information in module moda_bitmaps if the input element is part of a bitmap.
Definition: bitmaps.F90:13
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
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
subroutine pkc(chr, nchr, ibay, ibit)
Encode a character string within a specified number of bytes of an integer array, starting at the bit...
Definition: ciencode.F90:25
recursive subroutine ipkm(cbay, nbyt, n)
Encode an integer value within a specified number of bytes of a character string, up to a maximum of ...
Definition: ciencode.F90:194
subroutine pkb(nval, nbits, ibay, ibit)
Encode an integer value within a specified number of bits of an integer array, starting at the bit im...
Definition: ciencode.F90:140
subroutine pkb8(nval, nbits, ibay, ibit)
Encode an 8-byte integer value within a specified number of bits of an integer array,...
Definition: ciencode.F90:97
subroutine rdcmps(lun)
Read the next compressed BUFR data subset into internal arrays.
Definition: compress.F90:110
subroutine wrcmps(lunix)
Write a compressed BUFR data subset.
Definition: compress.F90:387
subroutine mvb(ib1, nb1, ib2, nb2, nbm)
Copy a specified number of bytes from one packed binary array to another.
Definition: copydata.F90:731
subroutine nemtba(lun, nemo, mtyp, msbt, inod)
Get information about a Table A descriptor from the internal DX BUFR tables.
Definition: dxtable.F90:1247
subroutine errwrt(str)
Specify a custom location for the logging of error and diagnostic messages generated by the NCEPLIBS-...
Definition: errwrt.F90:32
integer function ibfms(r8val)
Check whether a real*8 data value returned from a previous call to any of the NCEPLIBS-bufr values-re...
Definition: missing.F90:25
Wrap C NCEPLIBS-bufr functions so they can be called from within the Fortran part of the library.
Definition: bufrlib.F90:11
Declare arrays and variables used to store BUFR messages internally for multiple file IDs.
integer, dimension(:), allocatable ibay
Current data subset.
integer ibit
Bit pointer within ibay.
integer, dimension(:,:), allocatable mbay
Current BUFR message for each file ID.
integer, dimension(:), allocatable mbyt
Length (in bytes) of current BUFR message for each file ID.
integer maxbyt
Maximum length of an output BUFR message.
Declare arrays and variables used to store bitmaps internally within a data subset definition.
integer lstnod
Most recent jump/link table entry that was processed by function igetrfel() and whose corresponding v...
integer nbtm
Number of stored bitmaps for the current data subset (up to a maximum of mxbtm).
integer lstnodct
Current count of consecutive occurrences of lstnod.
logical linbtm
true if a bitmap is in the process of being read for the current data subset; false otherwise.
Declare arrays used to store, for each output file ID, a copy of the BUFR message that was most recen...
integer, dimension(:), allocatable msglen
Length (in integers) of BUFR message most recently written to each output file ID.
integer, dimension(:,:), allocatable msgtxt
BUFR message most recently written to each output file ID.
Declare arrays and variables needed to store long character strings (greater than 8 bytes) via subrou...
integer nh4wlc
Number of long character strings being stored.
character *14, dimension(:), allocatable sth4wlc
Table B mnemonics associated with long character strings.
integer, dimension(:), allocatable luh4wlc
File ID for associated output file.
character *120, dimension(:), allocatable chh4wlc
Long character strings.
Declare an array used to pack or unpack all of the values of a BUFR data subset.
integer *8, dimension(:), allocatable ival
BUFR data subset values.
Declare arrays which provide working space in several subprograms (usrtpl() and ufbcup()) which manip...
real *8, dimension(:), allocatable vtmp
val array elements for new sections of a growing subset buffer.
integer, dimension(:), allocatable itmp
inv array elements for new sections of a growing subset buffer.
Declare a variable used to indicate whether output BUFR messages should be compressed.
character ccmf
Flag indicating whether BUFR output messages are to be compressed; this variable is initialized to a ...
Declare arrays used to store information about the current BUFR message that is in the process of bei...
integer, dimension(:), allocatable inode
Table A mnemonic for type of BUFR message.
integer, dimension(:), allocatable idate
Section 1 date-time of message.
integer, dimension(:), allocatable nmsg
Current message pointer within logical unit.
integer, dimension(:), allocatable msub
Total number of data subsets in message.
integer, dimension(:), allocatable nsub
Current subset pointer within message.
Declare an array used to store a status code for each file ID if an error or other abnormal result oc...
integer, dimension(:), allocatable iscodes
Abnormal status codes.
Declare arrays and variables used to store the internal jump/link table.
integer, dimension(:), allocatable jseq
Temporary storage used in expanding sequences.
integer, dimension(:,:), allocatable iseq
Temporary storage used in expanding sequences.
integer, dimension(:), allocatable isc
Scale factors corresponding to tag and typ:
integer, dimension(:), allocatable ibt
Bit widths corresponding to tag and typ:
real *8, dimension(:), allocatable vali
Initialized data values corresponding to typ:
character *3, dimension(:), allocatable typ
Type indicators corresponding to tag:
integer, dimension(:), allocatable jmpb
Jump backward indices corresponding to tag and typ:
character *10, dimension(:), allocatable tag
Mnemonics in the jump/link table.
integer, dimension(:), allocatable itp
Integer type values corresponding to typ:
Declare an array used to store, for each file ID, the logical unit number corresponding to a separate...
integer, dimension(:), allocatable luncpy
Logical unit numbers used to copy long character strings between BUFR data subsets.
Declare an array used to store, for each file ID from which a BUFR message is currently being read as...
integer, dimension(:), allocatable msgunp
Flag indicating how to unpack data subsets from BUFR message:
Declare arrays for internal storage of pointers to BUFR data subset values.
integer, dimension(:), allocatable nbit
Length (in bits) of each packed data value in data subset.
integer, dimension(:), allocatable mbit
Pointer in data subset to first bit of each packed data value.
Declare arrays used to store data values and associated metadata for the current BUFR data subset in ...
integer, dimension(:), allocatable nval
Number of data values in BUFR data subset.
real *8, dimension(:,:), allocatable, target val
Data values.
integer, dimension(:,:), allocatable, target inv
Inventory pointer which links each data value to its corresponding node in the internal jump/link tab...
integer, dimension(:,:), allocatable nrfelm
Referenced data value, for data values which refer to a previous data value in the BUFR data subset v...
Declare arrays used in subroutine rcstpl() to store subset segments that are being copied from a subs...
integer, dimension(:,:), allocatable iutmp
inv array elements for new sections of a growing subset buffer.
real *8, dimension(:,:), allocatable vutmp
val array elements for new sections of a growing subset buffer.
recursive subroutine openbf(lunit, io, lundx)
Connect a new file to the NCEPLIBS-bufr software for input or output operations, or initialize the li...
recursive subroutine ufbcnt(lunit, kmsg, ksub)
Get the current location of the file pointer within a BUFR file, in terms of a message number countin...
recursive subroutine status(lunit, lun, il, im)
Check whether a specified Fortran logical unit number is currently connected to the NCEPLIBS-bufr sof...
recursive subroutine closmg(lunin)
Close the BUFR message that is currently open for writing within internal arrays associated with logi...
recursive subroutine readmg(lunxx, subset, jdate, iret)
Read the next BUFR message from logical unit abs(lunxx) into internal arrays.
Definition: readwritemg.F90:44
subroutine msgini(lun)
Initialize, within the internal arrays, a new uncompressed BUFR message for output.
subroutine msgwrt(lunit, mesg, mgbyt)
Perform final checks and updates on a BUFR message before writing it to a specified Fortran logical u...
recursive subroutine ufbpos(lunit, irec, isub, subset, jdate)
Jump forwards or backwards to a specified data subset within a BUFR file.
subroutine pad(ibay, ibit, ibyt, ipadb)
Pad a BUFR data subset with zeroed-out bits up to the next byte boundary.
subroutine rdtree(lun, iret)
Read the next uncompressed BUFR data subset into internal arrays.
subroutine wrtree(lun)
Pack a BUFR data subset.
subroutine msgupd(lunit, lun)
Write an uncompressed BUFR data subset.
subroutine usrtpl(lun, invn, nbmp)
Expand a subset template within internal arrays.
recursive subroutine writsa(lunxx, lmsgt, msgt, msgl)
Write a complete data subset into a BUFR message, and return each completed message within a memory a...
recursive integer function ireadns(lunit, subset, idate)
Call subroutine readns() and pass back its return code as the function value.
integer function nwords(n, lun)
Compute the length of a specified delayed replication sequence within a data subset.
recursive subroutine writsb(lunit)
Write a complete data subset into a BUFR message, for eventual output to logical unit lunit.
recursive subroutine invmrg(lubfi, lubfj)
Merge parts of data subsets which have duplicate space and time coordinates but different or unique o...
recursive integer function ireadsb(lunit)
Call subroutine readsb() and pass back its return code as the function value.
recursive subroutine readsb(lunit, iret)
Read the next data subset from a BUFR message.
Definition: readwritesb.F90:32
subroutine rcstpl(lun, iret)
Initialize a subset template within internal arrays.
recursive integer function lcmgdf(lunit, subset)
Check whether the subset definition for a given message type contains any long character strings (gre...
recursive subroutine rdmgsb(lunit, imsg, isub)
Read a specified data subset from a BUFR file.
recursive subroutine readns(lunit, subset, jdate, iret)
Read the next data subset from a BUFR file.
recursive subroutine readlc(lunit, chr, str)
Read a long character string (greater than 8 bytes) from a data subset.
recursive subroutine writlc(lunit, chr, str)
Write a long character string (greater than 8 bytes) to a data subset.
subroutine x48(iin4, iout8, nval)
Encode one or more 4-byte integer values as 8-byte integer values.
Definition: x4884.F90:18
subroutine x84(iin8, iout4, nval)
Encode one or more 8-byte integer values as 4-byte integer values.
Definition: x4884.F90:65