updated top-level README and version_decl for V4.4.2 (#1795)
[WRF.git] / external / io_grib2 / g2lib / addlocal.F
blobfc1b31b01bf3e6e03bca52f325d8c17902230e28
1       subroutine addlocal(cgrib,lcgrib,csec2,lcsec2,ierr)
2 !$$$  SUBPROGRAM DOCUMENTATION BLOCK
3 !                .      .    .                                       .
4 ! SUBPROGRAM:    addlocal 
5 !   PRGMMR: Gilbert         ORG: W/NP11    DATE: 2000-05-01
7 ! ABSTRACT: This subroutine adds a Local Use Section (Section 2) to 
8 !   a GRIB2 message.
9 !   This routine is used with routines "gribcreate", "addgrid", "addfield",
10 !   and "gribend" to create a complete GRIB2 message.  Subroutine
11 !   gribcreate must be called first to initialize a new GRIB2 message.
13 ! PROGRAM HISTORY LOG:
14 ! 2000-05-01  Gilbert
16 ! USAGE:    CALL addlocal(cgrib,lcgrib,csec2,lcsec2,ierr)
17 !   INPUT ARGUMENT LIST:
18 !     cgrib    - Character array to contain the GRIB2 message
19 !     lcgrib   - Maximum length (bytes) of array cgrib.
20 !     csec2    - Character array containing information to be added to
21 !                Section 2.
22 !     lcsec2   - Number of bytes of character array csec2 to be added to
23 !                Section 2.
25 !   OUTPUT ARGUMENT LIST:      
26 !     cgrib    - Character array to contain the GRIB2 message
27 !     ierr     - Error return code.
28 !                0 = no error
29 !                1 = GRIB message was not initialized.  Need to call
30 !                    routine gribcreate first.
31 !                2 = GRIB message already complete.  Cannot add new section.
32 !                3 = Sum of Section byte counts doesn't add to total byte count.
33 !                4 = Previous Section was not 1 or 7.
35 ! REMARKS: Note that the Local Use Section ( Section 2 ) can only follow
36 !          Section 1 or Section 7 in a GRIB2 message.
38 ! ATTRIBUTES:
39 !   LANGUAGE: Fortran 90
40 !   MACHINE:  IBM SP
42 !$$$
44       character(len=1),intent(inout) :: cgrib(lcgrib)
45       character(len=1),intent(in) :: csec2(lcsec2)
46       integer,intent(in) :: lcgrib,lcsec2
47       integer,intent(out) :: ierr
48       
49       character(len=4),parameter :: grib='GRIB',c7777='7777'
50       character(len=4):: ctemp
51       integer,parameter :: two=2
52       integer lensec2,iofst,ibeg,lencurr,len
54       ierr=0
56 !  Check to see if beginning of GRIB message exists
58       ctemp=cgrib(1)//cgrib(2)//cgrib(3)//cgrib(4)
59       if ( ctemp.ne.grib ) then
60         print *,'addlocal: GRIB not found in given message.'
61         print *,'addlocal: Call to routine gribcreate required',
62      &          ' to initialize GRIB messge.'
63         ierr=1
64         return
65       endif
67 !  Get current length of GRIB message
68 !  
69       call g2lib_gbyte(cgrib,lencurr,96,32)
71 !  Check to see if GRIB message is already complete
72 !  
73       ctemp=cgrib(lencurr-3)//cgrib(lencurr-2)//cgrib(lencurr-1)
74      &      //cgrib(lencurr)
75       if ( ctemp.eq.c7777 ) then
76         print *,'addlocal: GRIB message already complete.  Cannot',
77      &          ' add new section.'
78         ierr=2
79         return
80       endif
82 !  Loop through all current sections of the GRIB message to
83 !  find the last section number.
85       len=16    ! length of Section 0
86       do 
87       !    Get section number and length of next section
88         iofst=len*8
89         call g2lib_gbyte(cgrib,ilen,iofst,32)
90         iofst=iofst+32
91         call g2lib_gbyte(cgrib,isecnum,iofst,8)
92         len=len+ilen
93       !    Exit loop if last section reached
94         if ( len.eq.lencurr ) exit
95       !    If byte count for each section doesn't match current
96       !    total length, then there is a problem.
97         if ( len.gt.lencurr ) then
98           print *,'addlocal: Section byte counts don''t add to total.'
99           print *,'addlocal: Sum of section byte counts = ',len
100           print *,'addlocal: Total byte count in Section 0 = ',lencurr
101           ierr=3
102           return
103         endif
104       enddo
106 !  Section 2 can only be added after sections 1 and 7.
108       if ( (isecnum.ne.1) .and. (isecnum.ne.7) ) then
109         print *,'addlocal: Section 2 can only be added after Section',
110      &          ' 1 or Section 7.'
111         print *,'addlocal: Section ',isecnum,' was the last found in',
112      &          ' given GRIB message.'
113         ierr=4
114         return
115       endif
117 !  Add Section 2  - Local Use Section
119       ibeg=lencurr*8        !   Calculate offset for beginning of section 2
120       iofst=ibeg+32         !   leave space for length of section
121       call g2lib_sbyte(cgrib,two,iofst,8)     ! Store section number ( 2 )
122       istart=lencurr+5
123       cgrib(istart+1:istart+lcsec2)=csec2(1:lcsec2)
124       !
125       !   Calculate length of section 2 and store it in octets
126       !   1-4 of section 2.
127       !
128       lensec2=lcsec2+5      ! bytes
129       call g2lib_sbyte(cgrib,lensec2,ibeg,32)
132 !  Update current byte total of message in Section 0
134       call g2lib_sbyte(cgrib,lencurr+lensec2,96,32)
136       return
137       end