Versions Compared

Key

  • This line was added.
  • This line was removed.
  • Formatting was changed.

The following Fortran 77 program opens a GRIB file and calls INTF2 to interpolate fields to a 1.5/1.5 LatLon grid and extract a sub-area. The function INTIN is not used here because the input fields are in GRIB format and are self-defining.

Subroutines PBOPEN, PBCLOSE, PBGRIB and PBWRITE handle pure binary input and output files.

Depending on the parameter you wish to interpolate, there might be necessary to specify a path for the LSM files (using e.g. environment variable MARS_LSM_PATH.)

 

Code Block
titleFortran 77
linenumberstrue
C
C Copyright 2015 ECMWF.
C
C 

...

Code Block
      PROGRAM SAMPLE2
C
      IMPLICIT NONE
      INTEGER IPROD
      INTEGER INTV
      REAL REALV
      CHARACTER*20 CHARV
      DIMENSION INTV(4), REALV(4), CHARV(4)
C
      INTEGER JPGRIB, JPBYTES
C
      PARAMETER (JPGRIB = 2000000)
C
C     JPBYTES is the size in bytes on an 'INTEGER'
C     Set JPBYTES = 8 on a 64-bit machine.
CThis software is licensed under the terms of the Apache Licence
C Version 2.0 which can be obtained at http://www.apache.org/licenses/LICENSE-2.0.
C
C Unless required by applicable law or agreed to in writing, software
C distributed under the License is distributed on an "AS IS" BASIS,
C WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
C
C In applying this licence, ECMWF does not waive the privileges and immunities
C granted to it by virtue of its status as an intergovernmental organisation
C nor does it submit to any jurisdiction.
C

      PROGRAM EXAMPLE_INTF2
      IMPLICIT NONE


C     Parameters
      INTEGER JPGRIB   ! GRIB size (up to 1/16 deg)
      INTEGER JPBYTES  ! bytes/integer
      PARAMETER (JPGRIB = 33190420)
#ifdef INTEGER_8
      PARAMETER (JPBYTES = 8)
#else
      PARAMETER (JPBYTES = 4)
#endif


C
     Local INTEGERvariables
  INGRIB, NEWFLD
   INTEGER   DIMENSION  INGRIB(JPGRIB), NEWFLDINTV(JPGRIB4)
C
      REAL ZNFELDI, ZNFELDO
      DIMENSION ZNFELDIREALV(1), ZNFELDO(1)
C4)
      CHARACTER*20 CHARV

      INTEGER IUNIT1, IUNIT2, IREC, INLEN, NEWLENCHARACTER*128 INFILE, OUTFILE, ARG
      INTEGER INLEN, OUTLEN, U1, U2, IRET, NARGSN
      INTEGER*4 J
C INGRIB(JPGRIB), OUTGRIB(JPGRIB)


C     Externals
      INTEGER INTOUT, INTF2, IARGC


C      CHARACTER*128 INFILE, OUTFILE, CARG(4)
C
C **********************************************************************
C------------------------------------------------------------------


C     Pick up file names from command line.
C
      INFILE  = ' '
      NARGSOUTFILE = IARGC()' '
      IF( NARGSIARGC().LTEQ.4 ) THEN
        DO N = print*1,'Usage: interpolation_example -i inputfile -o outputfile'
 4, 2
           STOP
CALL GETARG(N,ARG)
        END  IF

   (ARG.EQ.'-i') THEN
    DO 101 J=1,NARGS
      CALL GETARG(J,CARG(J)N+1,INFILE)
 101  CONTINUE

      DO 102 J=1,NARGS,2
        IF(CARG(J) ELSEIF (ARG.EQ.'-io') THEN
            CALL INFILE=CARGGETARG(JN+1,OUTFILE)
          ENDIF
       ELSEIF(CARG(J).EQ.'-o') THEN ENDDO
      ENDIF
     OUTFILE=CARG(J+1)
        ELSE CALL CHECK(
     _  INDEX(INFILE,' ').EQ.1 .OR. INDEX(OUTFILE,' ').EQ.1,
     _   print*,'Usage: interpolationexample_exampleintf2 -i inputfileinfile -o outputfileoutfile' )

      INTV  = STOP0
      REALV = 0.
  END  IF
 102 CHARV CONTINUE= ''


C     Define the packing accuracy for the new field(s).
C
      INTV(1) = 24
      IRET = INTOUT('accuracy', INTV, REALV, CHARV)
      IFCALL CHECK( IRET.NE.0, ) THEN
        WRITE(*,*) ' First INTOUT failed.''INTOUT (accuracy) failed')
      INTV(1) = STOP
      ENDIF
C0


C     Define the geographical area for the new field(s).
C
      REALV(1) =  60.0
      REALV(2) = -10.0
      REALV(3) =  40.0
      REALV(4) =  15.0
      IRET = INTOUT('area', INTV, REALV, CHARV)
      IFCALL CHECK( IRET.NE.0 ) THEN
        WRITE(*,*) ' Second INTOUT failed.', 'INTOUT (area) failed')
      REALV = STOP
      ENDIF
C0.


C     Define the grid interval for the new field(s).
C
      REALV(1) = 1.5
      REALV(2) = 1.5
      IRET = INTOUT('grid', INTV, REALV, CHARV)
      IFCALL CHECK( IRET.NE.0 ) THEN
        WRITE(*,*) ' Third INTOUT failed.', 'INTOUT (grid) failed')
      REALV = STOP
      ENDIF
C0.


C     Open input and output files.
C
      CALL PBOPEN(IUNIT1U1, INFILE, 'r', IRET)
      IFCALL CHECK( IRET.NE.0, ) STOP ' PBOPEN failed (r)')
      CALL PBOPEN(IUNIT2U2, OUTFILE, 'w', IRET)
      IFCALL CHECK( IRET.NE.0 ) STOP, ' PBOPEN failed (w)')


C     Start of IPRODloop on =input 0
C
Cfields
      PRINT *, 'Start of loop through input GRIB-coded fields
C
 200interpolation...'
      N = 0
 220  CONTINUE
        IPRODN = IPRODN + 1
C
C       Read next product.
Cfield
        CALL PBGRIB(IUNIT1U1, INGRIB, JPGRIB*JPBYTES, IRECINLEN, IRET)
        IF ( IRET.EQ.-1 )) THEN
        IRET = 0
        GOTO 900290
      ENDIF
      IFCALL CHECK( IRET.NE.0 ), STOP ' PBGRIB failed')
C
C       Interpolate.
C
       PRINT WRITE(*,*) ' Interpolate productfield number #', IPRODN
       OUTLEN NEWLEN = JPGRIB
        INLEN  = IREC
        IRET = INTF2(INGRIB,INLEN,NEWFLDOUTGRIB,NEWLENOUTLEN)
      CALL  IF CHECK( IRET.NE.0, 'INTF failed') THEN

C     Write the new field to  WRITE(*,*) ' INTF failed.'file
      IF (OUTLEN.GT.0) THEN
        CALL  STOP
  PBWRITE(U2, OUTGRIB, OUTLEN, IRET)
      ENDIFELSE
C
C        WritePRINT the*, new'Output productsame toas file.
C
input'
        CALL PBWRITE( IUNIT2U2, NEWFLDINGRIB, NEWLEN*JPBYTESINLEN, IRET)
      ENDIF
      IFCALL CHECK( IRET.LT.(NEWLEN*JPBYTES) ) STOP ' 0, 'PBWRITE failed')

C
C       Loop back for next product.
Cfield
      GOTO 200220

C
C     Closedown.Close
C
 900290  CONTINUE
C
      IPROD = IPROD - 1
CALL PBCLOSE(U1, IRET)
      CALL WRITEPBCLOSE(*U2,*) ' All done after IRET)

      PRINT *, 'Interpolated ', IPROD(N-1), ' productsfield(s).'
      END


C
C     ------------------------------------------------------------------


      SUBROUTINE CHECK(OOPS,MSG)
      IMPLICIT CloseNONE
    input and outputLOGICAL files.
C OOPS
      CHARACTER MSG*(*)
      CALLIF PBCLOSE(IUNIT1, IRET)
(OOPS) THEN
        PRINT *, MSG
        CALL PBCLOSE(IUNIT2, IRET)
CEXIT(3)
      STOPENDIF
      END