| 1 | subroutine addfield(cgrib,lcgrib,ipdsnum,ipdstmpl,ipdstmplen, |
|---|
| 2 | & coordlist,numcoord,idrsnum,idrstmpl, |
|---|
| 3 | & idrstmplen,fld,ngrdpts,ibmap,bmap,ierr) |
|---|
| 4 | !$$$ SUBPROGRAM DOCUMENTATION BLOCK |
|---|
| 5 | ! . . . . |
|---|
| 6 | ! SUBPROGRAM: addfield |
|---|
| 7 | ! PRGMMR: Gilbert ORG: W/NP11 DATE: 2000-05-02 |
|---|
| 8 | ! |
|---|
| 9 | ! ABSTRACT: This subroutine packs up Sections 4 through 7 for a given field |
|---|
| 10 | ! and adds them to a GRIB2 message. They are Product Definition Section, |
|---|
| 11 | ! Data Representation Section, Bit-Map Section and Data Section, |
|---|
| 12 | ! respectively. |
|---|
| 13 | ! This routine is used with routines "gribcreate", "addlocal", "addgrid", |
|---|
| 14 | ! and "gribend" to create a complete GRIB2 message. Subroutine |
|---|
| 15 | ! gribcreate must be called first to initialize a new GRIB2 message. |
|---|
| 16 | ! Also, subroutine addgrid must be called after gribcreate and |
|---|
| 17 | ! before this routine to add the appropriate grid description to |
|---|
| 18 | ! the GRIB2 message. Also, a call to gribend is required to complete |
|---|
| 19 | ! GRIB2 message after all fields have been added. |
|---|
| 20 | ! |
|---|
| 21 | ! PROGRAM HISTORY LOG: |
|---|
| 22 | ! 2000-05-02 Gilbert |
|---|
| 23 | ! 2002-12-17 Gilbert - Added support for new templates using |
|---|
| 24 | ! PNG and JPEG2000 algorithms/templates. |
|---|
| 25 | ! 2004-06-22 Gilbert - Added check to determine if packing algorithm failed. |
|---|
| 26 | ! |
|---|
| 27 | ! USAGE: CALL addfield(cgrib,lcgrib,ipdsnum,ipdstmpl,ipdstmplen, |
|---|
| 28 | ! coordlist,numcoord,idrsnum,idrstmpl, |
|---|
| 29 | ! idrstmplen,fld,ngrdpts,ibmap,bmap,ierr) |
|---|
| 30 | ! INPUT ARGUMENT LIST: |
|---|
| 31 | ! cgrib - Character array to contain the GRIB2 message |
|---|
| 32 | ! lcgrib - Maximum length (bytes) of array cgrib. |
|---|
| 33 | ! ipdsnum - Product Definition Template Number ( see Code Table 4.0) |
|---|
| 34 | ! ipdstmpl - Contains the data values for the specified Product Definition |
|---|
| 35 | ! Template ( N=ipdsnum ). Each element of this integer |
|---|
| 36 | ! array contains an entry (in the order specified) of Product |
|---|
| 37 | ! Defintion Template 4.N |
|---|
| 38 | ! ipdstmplen - Max dimension of ipdstmpl() |
|---|
| 39 | ! coordlist- Array containg floating point values intended to document |
|---|
| 40 | ! the vertical discretisation associated to model data |
|---|
| 41 | ! on hybrid coordinate vertical levels. |
|---|
| 42 | ! numcoord - number of values in array coordlist. |
|---|
| 43 | ! idrsnum - Data Representation Template Number ( see Code Table 5.0 ) |
|---|
| 44 | ! idrstmpl - Contains the data values for the specified Data Representation |
|---|
| 45 | ! Template ( N=idrsnum ). Each element of this integer |
|---|
| 46 | ! array contains an entry (in the order specified) of Data |
|---|
| 47 | ! Representation Template 5.N |
|---|
| 48 | ! Note that some values in this template (eg. reference |
|---|
| 49 | ! values, number of bits, etc...) may be changed by the |
|---|
| 50 | ! data packing algorithms. |
|---|
| 51 | ! Use this to specify scaling factors and order of |
|---|
| 52 | ! spatial differencing, if desired. |
|---|
| 53 | ! idrstmplen - Max dimension of idrstmpl() |
|---|
| 54 | ! fld() - Array of data points to pack. |
|---|
| 55 | ! ngrdpts - Number of data points in grid. |
|---|
| 56 | ! i.e. size of fld and bmap. |
|---|
| 57 | ! ibmap - Bitmap indicator ( see Code Table 6.0 ) |
|---|
| 58 | ! 0 = bitmap applies and is included in Section 6. |
|---|
| 59 | ! 1-253 = Predefined bitmap applies |
|---|
| 60 | ! 254 = Previously defined bitmap applies to this field |
|---|
| 61 | ! 255 = Bit map does not apply to this product. |
|---|
| 62 | ! bmap() - Logical*1 array containing bitmap to be added. |
|---|
| 63 | ! ( if ibmap=0 or ibmap=254) |
|---|
| 64 | ! |
|---|
| 65 | ! OUTPUT ARGUMENT LIST: |
|---|
| 66 | ! cgrib - Character array to contain the GRIB2 message |
|---|
| 67 | ! ierr - Error return code. |
|---|
| 68 | ! 0 = no error |
|---|
| 69 | ! 1 = GRIB message was not initialized. Need to call |
|---|
| 70 | ! routine gribcreate first. |
|---|
| 71 | ! 2 = GRIB message already complete. Cannot add new section. |
|---|
| 72 | ! 3 = Sum of Section byte counts does not add to total |
|---|
| 73 | ! byte count. |
|---|
| 74 | ! 4 = Previous Section was not 3 or 7. |
|---|
| 75 | ! 5 = Could not find requested Product Definition Template. |
|---|
| 76 | ! 6 = Section 3 (GDS) not previously defined in message |
|---|
| 77 | ! 7 = Tried to use unsupported Data Representationi Template |
|---|
| 78 | ! 8 = Specified use of a previously defined bitmap, but one |
|---|
| 79 | ! does not exist in the GRIB message. |
|---|
| 80 | ! 9 = GDT of one of 5.50 through 5.53 required to pack |
|---|
| 81 | ! using DRT 5.51 |
|---|
| 82 | ! 10 = Error packing data field. |
|---|
| 83 | ! |
|---|
| 84 | ! REMARKS: Note that the Local Use Section ( Section 2 ) can only follow |
|---|
| 85 | ! Section 1 or Section 7 in a GRIB2 message. |
|---|
| 86 | ! |
|---|
| 87 | ! ATTRIBUTES: |
|---|
| 88 | ! LANGUAGE: Fortran 90 |
|---|
| 89 | ! MACHINE: IBM SP |
|---|
| 90 | ! |
|---|
| 91 | !$$$ |
|---|
| 92 | |
|---|
| 93 | use pdstemplates |
|---|
| 94 | use drstemplates |
|---|
| 95 | |
|---|
| 96 | character(len=1),intent(inout) :: cgrib(lcgrib) |
|---|
| 97 | integer,intent(in) :: ipdsnum,ipdstmpl(*) |
|---|
| 98 | integer,intent(in) :: idrsnum,numcoord,ipdstmplen,idrstmplen |
|---|
| 99 | integer,intent(in) :: lcgrib,ngrdpts,ibmap |
|---|
| 100 | real,intent(in) :: coordlist(numcoord) |
|---|
| 101 | real,target,intent(in) :: fld(ngrdpts) |
|---|
| 102 | integer,intent(out) :: ierr |
|---|
| 103 | integer,intent(inout) :: idrstmpl(*) |
|---|
| 104 | logical*1,intent(in) :: bmap(ngrdpts) |
|---|
| 105 | |
|---|
| 106 | character(len=4),parameter :: grib='GRIB',c7777='7777' |
|---|
| 107 | character(len=4):: ctemp |
|---|
| 108 | character(len=1),allocatable :: cpack(:) |
|---|
| 109 | real,pointer,dimension(:) :: pfld |
|---|
| 110 | real(4) :: coordieee(numcoord),re00 |
|---|
| 111 | integer(4) :: ire00,allones |
|---|
| 112 | integer :: mappds(ipdstmplen),intbmap(ngrdpts),mapdrs(idrstmplen) |
|---|
| 113 | integer,parameter :: zero=0,one=1,four=4,five=5,six=6,seven=7 |
|---|
| 114 | integer,parameter :: minsize=50000 |
|---|
| 115 | integer iofst,ibeg,lencurr,len,mappdslen,mapdrslen,lpos3 |
|---|
| 116 | integer width,height,ndpts |
|---|
| 117 | integer lensec3,lensec4,lensec5,lensec6,lensec7 |
|---|
| 118 | logical issec3,needext,isprevbmap |
|---|
| 119 | |
|---|
| 120 | ierr=0 |
|---|
| 121 | do jj=0,31 |
|---|
| 122 | allones=ibset(allones,jj) |
|---|
| 123 | enddo |
|---|
| 124 | ! |
|---|
| 125 | ! Check to see if beginning of GRIB message exists |
|---|
| 126 | ! |
|---|
| 127 | ctemp=cgrib(1)//cgrib(2)//cgrib(3)//cgrib(4) |
|---|
| 128 | if ( ctemp.ne.grib ) then |
|---|
| 129 | print *,'addfield: GRIB not found in given message.' |
|---|
| 130 | print *,'addfield: Call to routine gribcreate required', |
|---|
| 131 | & ' to initialize GRIB messge.' |
|---|
| 132 | ierr=1 |
|---|
| 133 | return |
|---|
| 134 | endif |
|---|
| 135 | ! |
|---|
| 136 | ! Get current length of GRIB message |
|---|
| 137 | ! |
|---|
| 138 | call g2lib_gbyte(cgrib,lencurr,96,32) |
|---|
| 139 | ! |
|---|
| 140 | ! Check to see if GRIB message is already complete |
|---|
| 141 | ! |
|---|
| 142 | ctemp=cgrib(lencurr-3)//cgrib(lencurr-2)//cgrib(lencurr-1) |
|---|
| 143 | & //cgrib(lencurr) |
|---|
| 144 | if ( ctemp.eq.c7777 ) then |
|---|
| 145 | print *,'addfield: GRIB message already complete. Cannot', |
|---|
| 146 | & ' add new section.' |
|---|
| 147 | ierr=2 |
|---|
| 148 | return |
|---|
| 149 | endif |
|---|
| 150 | ! |
|---|
| 151 | ! Loop through all current sections of the GRIB message to |
|---|
| 152 | ! find the last section number. |
|---|
| 153 | ! |
|---|
| 154 | issec3=.false. |
|---|
| 155 | isprevbmap=.false. |
|---|
| 156 | len=16 ! length of Section 0 |
|---|
| 157 | do |
|---|
| 158 | ! Get number and length of next section |
|---|
| 159 | iofst=len*8 |
|---|
| 160 | call g2lib_gbyte(cgrib,ilen,iofst,32) |
|---|
| 161 | iofst=iofst+32 |
|---|
| 162 | call g2lib_gbyte(cgrib,isecnum,iofst,8) |
|---|
| 163 | iofst=iofst+8 |
|---|
| 164 | ! Check if previous Section 3 exists and save location of |
|---|
| 165 | ! the section 3 in case needed later. |
|---|
| 166 | if (isecnum.eq.3) then |
|---|
| 167 | issec3=.true. |
|---|
| 168 | lpos3=len+1 |
|---|
| 169 | lensec3=ilen |
|---|
| 170 | endif |
|---|
| 171 | ! Check if a previous defined bitmap exists |
|---|
| 172 | if (isecnum.eq.6) then |
|---|
| 173 | call g2lib_gbyte(cgrib,ibmprev,iofst,8) |
|---|
| 174 | iofst=iofst+8 |
|---|
| 175 | if ((ibmprev.ge.0).and.(ibmprev.le.253)) isprevbmap=.true. |
|---|
| 176 | endif |
|---|
| 177 | len=len+ilen |
|---|
| 178 | ! Exit loop if last section reached |
|---|
| 179 | if ( len.eq.lencurr ) exit |
|---|
| 180 | ! If byte count for each section does not match current |
|---|
| 181 | ! total length, then there is a problem. |
|---|
| 182 | if ( len.gt.lencurr ) then |
|---|
| 183 | print *,'addfield: Section byte counts don''t add to total.' |
|---|
| 184 | print *,'addfield: Sum of section byte counts = ',len |
|---|
| 185 | print *,'addfield: Total byte count in Section 0 = ',lencurr |
|---|
| 186 | ierr=3 |
|---|
| 187 | return |
|---|
| 188 | endif |
|---|
| 189 | enddo |
|---|
| 190 | ! |
|---|
| 191 | ! Sections 4 through 7 can only be added after section 3 or 7. |
|---|
| 192 | ! |
|---|
| 193 | if ( (isecnum.ne.3) .and. (isecnum.ne.7) ) then |
|---|
| 194 | print *,'addfield: Sections 4-7 can only be added after', |
|---|
| 195 | & ' Section 3 or 7.' |
|---|
| 196 | print *,'addfield: Section ',isecnum,' was the last found in', |
|---|
| 197 | & ' given GRIB message.' |
|---|
| 198 | ierr=4 |
|---|
| 199 | return |
|---|
| 200 | ! |
|---|
| 201 | ! Sections 4 through 7 can only be added if section 3 was previously defined. |
|---|
| 202 | ! |
|---|
| 203 | elseif (.not.issec3) then |
|---|
| 204 | print *,'addfield: Sections 4-7 can only be added if Section', |
|---|
| 205 | & ' 3 was previously included.' |
|---|
| 206 | print *,'addfield: Section 3 was not found in', |
|---|
| 207 | & ' given GRIB message.' |
|---|
| 208 | print *,'addfield: Call to routine addgrid required', |
|---|
| 209 | & ' to specify Grid definition.' |
|---|
| 210 | ierr=6 |
|---|
| 211 | return |
|---|
| 212 | endif |
|---|
| 213 | ! |
|---|
| 214 | ! Add Section 4 - Product Definition Section |
|---|
| 215 | ! |
|---|
| 216 | ibeg=lencurr*8 ! Calculate offset for beginning of section 4 |
|---|
| 217 | iofst=ibeg+32 ! leave space for length of section |
|---|
| 218 | call g2lib_sbyte(cgrib,four,iofst,8) ! Store section number ( 4 ) |
|---|
| 219 | iofst=iofst+8 |
|---|
| 220 | call g2lib_sbyte(cgrib,numcoord,iofst,16) ! Store num of coordinate values |
|---|
| 221 | iofst=iofst+16 |
|---|
| 222 | call g2lib_sbyte(cgrib,ipdsnum,iofst,16) ! Store Prod Def Template num. |
|---|
| 223 | iofst=iofst+16 |
|---|
| 224 | ! |
|---|
| 225 | ! Get Product Definition Template |
|---|
| 226 | ! |
|---|
| 227 | call getpdstemplate(ipdsnum,mappdslen,mappds,needext,iret) |
|---|
| 228 | if (iret.ne.0) then |
|---|
| 229 | ierr=5 |
|---|
| 230 | return |
|---|
| 231 | endif |
|---|
| 232 | ! |
|---|
| 233 | ! Extend the Product Definition Template, if necessary. |
|---|
| 234 | ! The number of values in a specific template may vary |
|---|
| 235 | ! depending on data specified in the "static" part of the |
|---|
| 236 | ! template. |
|---|
| 237 | ! |
|---|
| 238 | if ( needext ) then |
|---|
| 239 | call extpdstemplate(ipdsnum,ipdstmpl,mappdslen,mappds) |
|---|
| 240 | endif |
|---|
| 241 | ! |
|---|
| 242 | ! Pack up each input value in array ipdstmpl into the |
|---|
| 243 | ! the appropriate number of octets, which are specified in |
|---|
| 244 | ! corresponding entries in array mappds. |
|---|
| 245 | ! |
|---|
| 246 | do i=1,mappdslen |
|---|
| 247 | nbits=iabs(mappds(i))*8 |
|---|
| 248 | if ( (mappds(i).ge.0).or.(ipdstmpl(i).ge.0) ) then |
|---|
| 249 | call g2lib_sbyte(cgrib,ipdstmpl(i),iofst,nbits) |
|---|
| 250 | else |
|---|
| 251 | call g2lib_sbyte(cgrib,one,iofst,1) |
|---|
| 252 | call g2lib_sbyte(cgrib,iabs(ipdstmpl(i)),iofst+1,nbits-1) |
|---|
| 253 | endif |
|---|
| 254 | iofst=iofst+nbits |
|---|
| 255 | enddo |
|---|
| 256 | ! |
|---|
| 257 | ! Add Optional list of vertical coordinate values |
|---|
| 258 | ! after the Product Definition Template, if necessary. |
|---|
| 259 | ! |
|---|
| 260 | if ( numcoord .ne. 0 ) then |
|---|
| 261 | call mkieee(coordlist,coordieee,numcoord) |
|---|
| 262 | call g2lib_sbytes(cgrib,coordieee,iofst,32,0,numcoord) |
|---|
| 263 | iofst=iofst+(32*numcoord) |
|---|
| 264 | endif |
|---|
| 265 | ! |
|---|
| 266 | ! Calculate length of section 4 and store it in octets |
|---|
| 267 | ! 1-4 of section 4. |
|---|
| 268 | ! |
|---|
| 269 | lensec4=(iofst-ibeg)/8 |
|---|
| 270 | call g2lib_sbyte(cgrib,lensec4,ibeg,32) |
|---|
| 271 | ! |
|---|
| 272 | ! Pack Data using appropriate algorithm |
|---|
| 273 | ! |
|---|
| 274 | ! |
|---|
| 275 | ! Get Data Representation Template |
|---|
| 276 | ! |
|---|
| 277 | call getdrstemplate(idrsnum,mapdrslen,mapdrs,needext,iret) |
|---|
| 278 | if (iret.ne.0) then |
|---|
| 279 | ierr=5 |
|---|
| 280 | return |
|---|
| 281 | endif |
|---|
| 282 | ! |
|---|
| 283 | ! contract data field, removing data at invalid grid points, |
|---|
| 284 | ! if bit-map is provided with field. |
|---|
| 285 | ! |
|---|
| 286 | if ( ibmap.eq.0 .OR. ibmap.eq.254 ) then |
|---|
| 287 | allocate(pfld(ngrdpts)) |
|---|
| 288 | ndpts=0; |
|---|
| 289 | do jj=1,ngrdpts |
|---|
| 290 | intbmap(jj)=0 |
|---|
| 291 | if ( bmap(jj) ) then |
|---|
| 292 | intbmap(jj)=1 |
|---|
| 293 | ndpts=ndpts+1 |
|---|
| 294 | pfld(ndpts)=fld(jj); |
|---|
| 295 | endif |
|---|
| 296 | enddo |
|---|
| 297 | else |
|---|
| 298 | ndpts=ngrdpts; |
|---|
| 299 | pfld=>fld; |
|---|
| 300 | endif |
|---|
| 301 | lcpack=0 |
|---|
| 302 | nsize=ndpts*4 |
|---|
| 303 | if (nsize .lt. minsize) nsize=minsize |
|---|
| 304 | allocate(cpack(nsize),stat=istat) |
|---|
| 305 | if (idrsnum.eq.0) then ! Simple Packing |
|---|
| 306 | call simpack(pfld,ndpts,idrstmpl,cpack,lcpack) |
|---|
| 307 | elseif (idrsnum.eq.2.or.idrsnum.eq.3) then ! Complex Packing |
|---|
| 308 | call cmplxpack(pfld,ndpts,idrsnum,idrstmpl,cpack,lcpack) |
|---|
| 309 | elseif (idrsnum.eq.50) then ! Sperical Harmonic Simple Packing |
|---|
| 310 | call simpack(pfld(2),ndpts-1,idrstmpl,cpack,lcpack) |
|---|
| 311 | call mkieee(real(pfld(1)),re00,1) ! ensure RE(0,0) value is IEEE format |
|---|
| 312 | !call g2lib_gbyte(re00,idrstmpl(5),0,32) |
|---|
| 313 | ire00=transfer(re00,ire00) |
|---|
| 314 | idrstmpl(5)=ire00 |
|---|
| 315 | elseif (idrsnum.eq.51) then ! Sperical Harmonic Complex Packing |
|---|
| 316 | call getpoly(cgrib(lpos3),lensec3,jj,kk,mm) |
|---|
| 317 | if (jj.ne.0 .AND. kk.ne.0 .AND. mm.ne.0) then |
|---|
| 318 | call specpack(pfld,ndpts,jj,kk,mm,idrstmpl,cpack,lcpack) |
|---|
| 319 | else |
|---|
| 320 | print *,'addfield: Cannot pack DRT 5.51.' |
|---|
| 321 | ierr=9 |
|---|
| 322 | return |
|---|
| 323 | endif |
|---|
| 324 | #ifdef USE_JPEG2000 |
|---|
| 325 | elseif (idrsnum.eq.40 .OR. idrsnum.eq.40000) then ! JPEG2000 encoding |
|---|
| 326 | if (ibmap.eq.255) then |
|---|
| 327 | call getdim(cgrib(lpos3),lensec3,width,height,iscan) |
|---|
| 328 | if ( width.eq.0 .OR. height.eq.0 ) then |
|---|
| 329 | width=ndpts |
|---|
| 330 | height=1 |
|---|
| 331 | elseif ( width.eq.allones .OR. height.eq.allones ) then |
|---|
| 332 | width=ndpts |
|---|
| 333 | height=1 |
|---|
| 334 | elseif ( ibits(iscan,5,1) .eq. 1) then ! Scanning mode: bit 3 |
|---|
| 335 | itemp=width |
|---|
| 336 | width=height |
|---|
| 337 | height=itemp |
|---|
| 338 | endif |
|---|
| 339 | else |
|---|
| 340 | width=ndpts |
|---|
| 341 | height=1 |
|---|
| 342 | endif |
|---|
| 343 | lcpack=nsize |
|---|
| 344 | call jpcpack(pfld,width,height,idrstmpl,cpack,lcpack,ierr) |
|---|
| 345 | #endif /* USE_JPEG2000 */ |
|---|
| 346 | #ifdef USE_PNG |
|---|
| 347 | elseif (idrsnum.eq.41 .OR. idrsnum.eq.40010) then ! PNG encoding |
|---|
| 348 | if (ibmap.eq.255) then |
|---|
| 349 | call getdim(cgrib(lpos3),lensec3,width,height,iscan) |
|---|
| 350 | if ( width.eq.0 .OR. height.eq.0 ) then |
|---|
| 351 | width=ndpts |
|---|
| 352 | height=1 |
|---|
| 353 | elseif ( width.eq.allones .OR. height.eq.allones ) then |
|---|
| 354 | width=ndpts |
|---|
| 355 | height=1 |
|---|
| 356 | elseif ( ibits(iscan,5,1) .eq. 1) then ! Scanning mode: bit 3 |
|---|
| 357 | itemp=width |
|---|
| 358 | width=height |
|---|
| 359 | height=itemp |
|---|
| 360 | endif |
|---|
| 361 | else |
|---|
| 362 | width=ndpts |
|---|
| 363 | height=1 |
|---|
| 364 | endif |
|---|
| 365 | call pngpack(pfld,width,height,idrstmpl,cpack,lcpack) |
|---|
| 366 | #endif /* USE_PNG */ |
|---|
| 367 | else |
|---|
| 368 | print *,'addfield: Data Representation Template 5.',idrsnum, |
|---|
| 369 | * ' not yet implemented.' |
|---|
| 370 | ierr=7 |
|---|
| 371 | return |
|---|
| 372 | endif |
|---|
| 373 | if ( ibmap.eq.0 .OR. ibmap.eq.254 ) then |
|---|
| 374 | deallocate(pfld) |
|---|
| 375 | endif |
|---|
| 376 | if ( lcpack .lt. 0 ) then |
|---|
| 377 | if( allocated(cpack) )deallocate(cpack) |
|---|
| 378 | ierr=10 |
|---|
| 379 | return |
|---|
| 380 | endif |
|---|
| 381 | |
|---|
| 382 | ! |
|---|
| 383 | ! Add Section 5 - Data Representation Section |
|---|
| 384 | ! |
|---|
| 385 | ibeg=iofst ! Calculate offset for beginning of section 5 |
|---|
| 386 | iofst=ibeg+32 ! leave space for length of section |
|---|
| 387 | call g2lib_sbyte(cgrib,five,iofst,8) ! Store section number ( 5 ) |
|---|
| 388 | iofst=iofst+8 |
|---|
| 389 | call g2lib_sbyte(cgrib,ndpts,iofst,32) ! Store num of actual data points |
|---|
| 390 | iofst=iofst+32 |
|---|
| 391 | call g2lib_sbyte(cgrib,idrsnum,iofst,16) ! Store Data Repr. Template num. |
|---|
| 392 | iofst=iofst+16 |
|---|
| 393 | ! |
|---|
| 394 | ! Pack up each input value in array idrstmpl into the |
|---|
| 395 | ! the appropriate number of octets, which are specified in |
|---|
| 396 | ! corresponding entries in array mapdrs. |
|---|
| 397 | ! |
|---|
| 398 | do i=1,mapdrslen |
|---|
| 399 | nbits=iabs(mapdrs(i))*8 |
|---|
| 400 | if ( (mapdrs(i).ge.0).or.(idrstmpl(i).ge.0) ) then |
|---|
| 401 | call g2lib_sbyte(cgrib,idrstmpl(i),iofst,nbits) |
|---|
| 402 | else |
|---|
| 403 | call g2lib_sbyte(cgrib,one,iofst,1) |
|---|
| 404 | call g2lib_sbyte(cgrib,iabs(idrstmpl(i)),iofst+1,nbits-1) |
|---|
| 405 | endif |
|---|
| 406 | iofst=iofst+nbits |
|---|
| 407 | enddo |
|---|
| 408 | ! |
|---|
| 409 | ! Calculate length of section 5 and store it in octets |
|---|
| 410 | ! 1-4 of section 5. |
|---|
| 411 | ! |
|---|
| 412 | lensec5=(iofst-ibeg)/8 |
|---|
| 413 | call g2lib_sbyte(cgrib,lensec5,ibeg,32) |
|---|
| 414 | |
|---|
| 415 | ! |
|---|
| 416 | ! Add Section 6 - Bit-Map Section |
|---|
| 417 | ! |
|---|
| 418 | ibeg=iofst ! Calculate offset for beginning of section 6 |
|---|
| 419 | iofst=ibeg+32 ! leave space for length of section |
|---|
| 420 | call g2lib_sbyte(cgrib,six,iofst,8) ! Store section number ( 6 ) |
|---|
| 421 | iofst=iofst+8 |
|---|
| 422 | call g2lib_sbyte(cgrib,ibmap,iofst,8) ! Store Bit Map indicator |
|---|
| 423 | iofst=iofst+8 |
|---|
| 424 | ! |
|---|
| 425 | ! Store bitmap, if supplied |
|---|
| 426 | ! |
|---|
| 427 | if (ibmap.eq.0) then |
|---|
| 428 | call g2lib_sbytes(cgrib,intbmap,iofst,1,0,ngrdpts) ! Store BitMap |
|---|
| 429 | iofst=iofst+ngrdpts |
|---|
| 430 | endif |
|---|
| 431 | ! |
|---|
| 432 | ! If specifying a previously defined bit-map, make sure |
|---|
| 433 | ! one already exists in the current GRIB message. |
|---|
| 434 | ! |
|---|
| 435 | if ((ibmap.eq.254).and.(.not.isprevbmap)) then |
|---|
| 436 | print *,'addfield: Requested previously defined bitmap, ', |
|---|
| 437 | & ' but one does not exist in the current GRIB message.' |
|---|
| 438 | ierr=8 |
|---|
| 439 | return |
|---|
| 440 | endif |
|---|
| 441 | ! |
|---|
| 442 | ! Calculate length of section 6 and store it in octets |
|---|
| 443 | ! 1-4 of section 6. Pad to end of octect, if necessary. |
|---|
| 444 | ! |
|---|
| 445 | left=8-mod(iofst,8) |
|---|
| 446 | if (left.ne.8) then |
|---|
| 447 | call g2lib_sbyte(cgrib,zero,iofst,left) ! Pad with zeros to fill Octet |
|---|
| 448 | iofst=iofst+left |
|---|
| 449 | endif |
|---|
| 450 | lensec6=(iofst-ibeg)/8 |
|---|
| 451 | call g2lib_sbyte(cgrib,lensec6,ibeg,32) |
|---|
| 452 | |
|---|
| 453 | ! |
|---|
| 454 | ! Add Section 7 - Data Section |
|---|
| 455 | ! |
|---|
| 456 | ibeg=iofst ! Calculate offset for beginning of section 7 |
|---|
| 457 | iofst=ibeg+32 ! leave space for length of section |
|---|
| 458 | call g2lib_sbyte(cgrib,seven,iofst,8) ! Store section number ( 7 ) |
|---|
| 459 | iofst=iofst+8 |
|---|
| 460 | ! Store Packed Binary Data values, if non-constant field |
|---|
| 461 | if (lcpack.ne.0) then |
|---|
| 462 | ioctet=iofst/8 |
|---|
| 463 | cgrib(ioctet+1:ioctet+lcpack)=cpack(1:lcpack) |
|---|
| 464 | iofst=iofst+(8*lcpack) |
|---|
| 465 | endif |
|---|
| 466 | ! |
|---|
| 467 | ! Calculate length of section 7 and store it in octets |
|---|
| 468 | ! 1-4 of section 7. |
|---|
| 469 | ! |
|---|
| 470 | lensec7=(iofst-ibeg)/8 |
|---|
| 471 | call g2lib_sbyte(cgrib,lensec7,ibeg,32) |
|---|
| 472 | |
|---|
| 473 | if( allocated(cpack) )deallocate(cpack) |
|---|
| 474 | ! |
|---|
| 475 | ! Update current byte total of message in Section 0 |
|---|
| 476 | ! |
|---|
| 477 | newlen=lencurr+lensec4+lensec5+lensec6+lensec7 |
|---|
| 478 | call g2lib_sbyte(cgrib,newlen,96,32) |
|---|
| 479 | |
|---|
| 480 | return |
|---|
| 481 | end |
|---|
| 482 | |
|---|
