subroutine gdt2gds(igds,igdstmpl,idefnum,ideflist,kgds, & igrid,iret) C$$$ SUBPROGRAM DOCUMENTATION BLOCK C . . . . C SUBPROGRAM: gdt2gds C PRGMMR: Gilbert ORG: W/NP11 DATE: 2003-06-17 C C ABSTRACT: This routine converts grid information from a GRIB2 C Grid Description Section as well as its C Grid Definition Template to GRIB1 GDS info. In addition, C a check is made to determine if the grid is an NCEP C predefined grid. C C PROGRAM HISTORY LOG: C 2003-06-17 Gilbert C 2004-04-27 Gilbert - Added support for gaussian grids. C C USAGE: CALL gdt2gds(igds,igdstmpl,idefnum,ideflist,kgds,igrid,iret) C INPUT ARGUMENT LIST: C igds() - Contains information read from the appropriate GRIB Grid C Definition Section 3 for the field being returned. C Must be dimensioned >= 5. C igds(1)=Source of grid definition (see Code Table 3.0) C igds(2)=Number of grid points in the defined grid. C igds(3)=Number of octets needed for each C additional grid points definition. C Used to define number of C points in each row ( or column ) for C non-regular grids. C = 0, if using regular grid. C igds(4)=Interpretation of list for optional points C definition. (Code Table 3.11) C igds(5)=Grid Definition Template Number (Code Table 3.1) C igdstmpl() - Grid Definition Template values for GDT 3.igds(5) C idefnum - The number of entries in array ideflist. C i.e. number of rows ( or columns ) C for which optional grid points are defined. C ideflist() - Optional integer array containing C the number of grid points contained in each row (or column). C C OUTPUT ARGUMENT LIST: (INCLUDING WORK ARRAYS) C kgds() - GRIB1 GDS as described in w3fi63 format. C igrid - NCEP predifined GRIB1 grid number C set to 255, if not NCEP grid C iret - Error return value: C 0 = Successful C 1 = Unrecognized GRIB2 GDT number 3.igds(5) C C REMARKS: LIST CAVEATS, OTHER HELPFUL HINTS OR INFORMATION C C ATTRIBUTES: C LANGUAGE: INDICATE EXTENSIONS, COMPILER OPTIONS C MACHINE: IBM SP C C$$$ ! integer,intent(in) :: idefnum integer,intent(in) :: igds(*),igdstmpl(*),ideflist(*) integer,intent(out) :: kgds(*),igrid,iret integer :: kgds72(200),kgds71(200),idum(200),jdum(200) iret=0 if (igds(5).eq.0) then ! Lat/Lon grid kgds(1)=0 kgds(2)=igdstmpl(8) ! Ni kgds(3)=igdstmpl(9) ! Nj kgds(4)=igdstmpl(12)/1000 ! Lat of 1st grid point kgds(5)=igdstmpl(13)/1000 ! Long of 1st grid point kgds(6)=0 ! resolution and component flags if (igdstmpl(1)==2 ) kgds(6)=64 if ( btest(igdstmpl(14),4).OR.btest(igdstmpl(14),5) ) & kgds(6)=kgds(6)+128 if ( btest(igdstmpl(14),3) ) kgds(6)=kgds(6)+8 kgds(7)=igdstmpl(15)/1000 ! Lat of last grid point kgds(8)=igdstmpl(16)/1000 ! Long of last grid point kgds(9)=igdstmpl(17)/1000 ! Di kgds(10)=igdstmpl(18)/1000 ! Dj kgds(11)=igdstmpl(19) ! Scanning mode kgds(12)=0 kgds(13)=0 kgds(14)=0 kgds(15)=0 kgds(16)=0 kgds(17)=0 kgds(18)=0 kgds(19)=0 kgds(20)=255 kgds(21)=0 kgds(22)=0 ! ! Process irreg grid stuff, if necessary ! if ( idefnum.ne.0 ) then if ( igdstmpl(8).eq.-1 ) then kgds(2)=65535 kgds(9)=65535 endif if ( igdstmpl(9).eq.-1 ) then kgds(3)=65535 kgds(10)=65535 endif kgds(19)=0 kgds(20)=33 if ( kgds(1).eq.1.OR.kgds(1).eq.3 ) kgds(20)=43 kgds(21)=igds(2) ! num of grid points do j=1,idefnum kgds(21+j)=ideflist(j) enddo endif elseif (igds(5).eq.10) then ! Mercator grid kgds(1)=1 ! Grid Definition Template number kgds(2)=igdstmpl(8) ! Ni kgds(3)=igdstmpl(9) ! Nj kgds(4)=igdstmpl(10)/1000 ! Lat of 1st grid point kgds(5)=igdstmpl(11)/1000 ! Long of 1st grid point kgds(6)=0 ! resolution and component flags if (igdstmpl(1)==2 ) kgds(6)=64 if ( btest(igdstmpl(12),4).OR.btest(igdstmpl(12),5) ) & kgds(6)=kgds(6)+128 if ( btest(igdstmpl(12),3) ) kgds(6)=kgds(6)+8 kgds(7)=igdstmpl(14)/1000 ! Lat of last grid point kgds(8)=igdstmpl(15)/1000 ! Long of last grid point kgds(9)=igdstmpl(13)/1000 ! Lat intersects earth kgds(10)=0 kgds(11)=igdstmpl(16) ! Scanning mode kgds(12)=igdstmpl(18)/1000 ! Di kgds(13)=igdstmpl(19)/1000 ! Dj kgds(14)=0 kgds(15)=0 kgds(16)=0 kgds(17)=0 kgds(18)=0 kgds(19)=0 kgds(20)=255 kgds(21)=0 kgds(22)=0 elseif (igds(5).eq.30) then ! Lambert Conformal Grid kgds(1)=3 kgds(2)=igdstmpl(8) ! Nx kgds(3)=igdstmpl(9) ! Ny kgds(4)=igdstmpl(10)/1000 ! Lat of 1st grid point kgds(5)=igdstmpl(11)/1000 ! Long of 1st grid point kgds(6)=0 ! resolution and component flags if (igdstmpl(1)==2 ) kgds(6)=64 if ( btest(igdstmpl(12),4).OR.btest(igdstmpl(12),5) ) & kgds(6)=kgds(6)+128 if ( btest(igdstmpl(12),3) ) kgds(6)=kgds(6)+8 kgds(7)=igdstmpl(14)/1000 ! Lon of orientation kgds(8)=igdstmpl(15)/1000 ! Dx kgds(9)=igdstmpl(16)/1000 ! Dy kgds(10)=igdstmpl(17) ! Projection Center Flag kgds(11)=igdstmpl(18) ! Scanning mode kgds(12)=igdstmpl(19)/1000 ! Lat in 1 kgds(13)=igdstmpl(20)/1000 ! Lat in 2 kgds(14)=igdstmpl(21)/1000 ! Lat of S. Pole of projection kgds(15)=igdstmpl(22)/1000 ! Lon of S. Pole of projection kgds(16)=0 kgds(17)=0 kgds(18)=0 kgds(19)=0 kgds(20)=255 kgds(21)=0 kgds(22)=0 elseif (igds(5).eq.40) then ! Gaussian Lat/Lon grid kgds(1)=4 kgds(2)=igdstmpl(8) ! Ni kgds(3)=igdstmpl(9) ! Nj kgds(4)=igdstmpl(12)/1000 ! Lat of 1st grid point kgds(5)=igdstmpl(13)/1000 ! Long of 1st grid point kgds(6)=0 ! resolution and component flags if (igdstmpl(1)==2 ) kgds(6)=64 if ( btest(igdstmpl(14),4).OR.btest(igdstmpl(14),5) ) & kgds(6)=kgds(6)+128 if ( btest(igdstmpl(14),3) ) kgds(6)=kgds(6)+8 kgds(7)=igdstmpl(15)/1000 ! Lat of last grid point kgds(8)=igdstmpl(16)/1000 ! Long of last grid point kgds(9)=igdstmpl(17)/1000 ! Di kgds(10)=igdstmpl(18) ! N - Number of parallels kgds(11)=igdstmpl(19) ! Scanning mode kgds(12)=0 kgds(13)=0 kgds(14)=0 kgds(15)=0 kgds(16)=0 kgds(17)=0 kgds(18)=0 kgds(19)=0 kgds(20)=255 kgds(21)=0 kgds(22)=0 elseif (igds(5).eq.20) then ! Polar Stereographic Grid kgds(1)=5 kgds(2)=igdstmpl(8) ! Nx kgds(3)=igdstmpl(9) ! Ny kgds(4)=igdstmpl(10)/1000 ! Lat of 1st grid point kgds(5)=igdstmpl(11)/1000 ! Long of 1st grid point kgds(6)=0 ! resolution and component flags if (igdstmpl(1)==2 ) kgds(6)=64 if ( btest(igdstmpl(12),4).OR.btest(igdstmpl(12),5) ) & kgds(6)=kgds(6)+128 if ( btest(igdstmpl(12),3) ) kgds(6)=kgds(6)+8 kgds(7)=igdstmpl(14)/1000 ! Lon of orientation kgds(8)=igdstmpl(15)/1000 ! Dx kgds(9)=igdstmpl(16)/1000 ! Dy kgds(10)=igdstmpl(17) ! Projection Center Flag kgds(11)=igdstmpl(18) ! Scanning mode kgds(12)=0 kgds(13)=0 kgds(14)=0 kgds(15)=0 kgds(16)=0 kgds(17)=0 kgds(18)=0 kgds(19)=0 kgds(20)=255 kgds(21)=0 kgds(22)=0 else Print *,'gdt2gds: Unrecognized GRIB2 GDT = 3.',igds(5) iret=1 kgds(1:22)=0 return endif ! ! Can we determine NCEP grid number ? ! igrid=255 do j=254,1,-1 !do j=225,225 kgds71=0 kgds72=0 call w3fi71(j,kgds71,ierr) if ( ierr.ne.0 ) cycle ! convert W to E for longitudes if ( kgds71(3).eq.0 ) then ! lat/lon if ( kgds71(7).lt.0 ) kgds71(7)=360000+kgds71(7) if ( kgds71(10).lt.0 ) kgds71(10)=360000+kgds71(10) elseif ( kgds71(3).eq.1 ) then ! mercator if ( kgds71(7).lt.0 ) kgds71(7)=360000+kgds71(7) if ( kgds71(10).lt.0 ) kgds71(10)=360000+kgds71(10) elseif ( kgds71(3).eq.3 ) then ! lambert conformal if ( kgds71(7).lt.0 ) kgds71(7)=360000+kgds71(7) if ( kgds71(9).lt.0 ) kgds71(9)=360000+kgds71(9) if ( kgds71(18).lt.0 ) kgds71(18)=360000+kgds71(18) elseif ( kgds71(3).eq.4 ) then ! Guassian lat/lon if ( kgds71(7).lt.0 ) kgds71(7)=360000+kgds71(7) if ( kgds71(10).lt.0 ) kgds71(10)=360000+kgds71(10) elseif ( kgds71(3).eq.5 ) then ! polar stereographic if ( kgds71(7).lt.0 ) kgds71(7)=360000+kgds71(7) if ( kgds71(9).lt.0 ) kgds71(9)=360000+kgds71(9) endif call r63w72(idum,kgds,jdum,kgds72) if ( kgds72(3).eq.3 ) kgds72(14)=0 ! lambert conformal fix if ( kgds72(3).eq.1 ) kgds72(15:18)=0 ! mercator fix if ( kgds72(3).eq.5 ) kgds72(14:18)=0 ! polar str fix !print *,'SAGT71:',(kgds71(k),k=1,30) !print *,'SAGT72:',(kgds72(k),k=1,30) if ( all(kgds71.eq.kgds72) ) then igrid=j return endif enddo return end