!####################################################################### ! #define TRACEBACK write (gol,'("in ",a," (",a,", line",i5,")")') rname, __FILE__, __LINE__; call goErr #define IF_NOTOK_RETURN(action) if (status/=0) then; TRACEBACK; action; return; end if #define IF_ERROR_RETURN(action) if (status >0) then; TRACEBACK; action; return; end if ! ! trap errors from GRIB_API: #define IF_GRIB_NOTOK_RETURN(action) if (status/=GRIB_SUCCESS) then; call GRIB_Get_Error_String(status,gol); call goErr; TRACEBACK; action; return; end if ! !####################################################################### MODULE FILE_GRIB_API USE GO, only : gol, goPr, goErr USE GRIB_API, only : HGRIB => kindOfInt ! handle kind USE GRIB_API, only : GRIB_SUCCESS, GRIB_Get_Error_String IMPLICIT NONE PRIVATE PUBLIC :: TGribFile PUBLIC :: Init, Done, Opened PUBLIC :: Get PUBLIC :: Check PUBLIC :: ReadRecord ! --- const ------------------------------ character(len=*), parameter :: mname = 'file_grib_api' ! --- types ------------------------------- TYPE TGribFile integer(kind=HGRIB) :: fu character(len=256) :: fname ! grib id and return value (end_of_file check) integer(kind=HGRIB) :: igrib integer :: iret logical :: opened END TYPE TGribFile ! --- interfaces --------------------------- interface Init module procedure gribfile_Init end interface interface Done module procedure gribfile_Done end interface interface Opened module procedure gribfile_Opened end interface interface ReadRecord module procedure gribfile_ReadRecord end interface interface Get module procedure gribfile_Get end interface interface Check module procedure gribfile_Check end interface CONTAINS ! ================================================================= ! Open gribfile. ! ! USAGE ! call Init( gribfile, 'input.gb', 'r'|'w' ) ! ! DESCRIPTION ! Interface around routine 'pbOpen'. ! In 'gribfile', space is allocated to store grib sections. ! subroutine gribfile_Init( F, file, mode, status ) use Grib_API, only : Grib_Open_File ! --- in/out ----------------------- type(TGribFile), intent(out) :: F character(len=*), intent(in) :: file character(len=*), intent(in) :: mode integer, intent(out) :: status ! --- const ------------------------ character(len=*), parameter :: rname = mname//'/gribfile_Init' ! --- begin ------------------------ F%fname = file ! filename select case ( mode ) case ( 'r' ) call grib_open_file ( F%fu, file, 'r' ) F%opened=.true. case ( 'w' ) write (gol,'("writing not yet supported for grib-api")') ; call goErr TRACEBACK; status=1; return case default write (gol,'("unknown mode `",a,"`")') mode; call goErr TRACEBACK; status=1; return end select ! no message loaded yet: F%igrib = -1 end subroutine gribfile_Init ! === ! Close gribfile, clear buffers ! ! USAGE ! call Done( gribfile ) ! ! DESCRIPTION ! Interface around routine 'pbClose'. ! subroutine gribfile_Done( F, status ) use Grib_API, only : GRIB_Release use Grib_API, only : Grib_Close_File ! --- in/out ----------------------- type(TGribFile), intent(inout) :: F integer, intent(inout) :: status ! --- const ------------------------ character(len=*), parameter :: rname = mname//'/gribfile_Done' ! --- begin ------------------------ ! old message loaded ? if ( F%igrib > 0 ) then ! remove from memory: call GRIB_Release( F%igrib, status=status ) IF_GRIB_NOTOK_RETURN(status=1) ! no message loaded yet: F%igrib = -1 end if ! write(*,'(1x,a,a)') 'Close grib file ',trim(F%fname) call grib_close_file ( F%fu, status ) if ( status /= 0 ) then write (gol,'("from closing grib file:")'); call goErr write (gol,'(" status : ",i6)') status; call goErr write (gol,'(" file : ",a)') trim(F%fname); call goErr TRACEBACK; status=1; return else F%opened=.false. end if end subroutine gribfile_Done ! === ! Grib file opened ? logical function gribfile_Opened( gribfile ) ! --- in/out ------------------------------ type(TGribFile), intent(in) :: gribfile ! --- begin ------------------------- !inquire( unit=gribfile%fu, opened=gribfile_Opened ) gribfile_Opened=gribfile%opened end function gribfile_Opened ! ================================================================== ! USAGE ! call ReadRecord( gribfile [,debug=0|1] [,status=status] ) ! ! DESCRIPTION ! Reads next grib record (a/k/a grib message) in buffers. ! ! RETURN STATUS ! 0 : no error ! 1 : eof ! 2 : some error ! Execution might stop if other errors are encountered. ! subroutine gribfile_ReadRecord( F, status ) use Grib_API, only : GRIB_END_OF_FILE use Grib_API, only : GRIB_SUCCESS use Grib_API, only : Grib_New_From_File use Grib_API, only : Grib_Get use GRIB_API, only : GRIB_Release ! --- in/out ------------------------- type(TGribFile), intent(inout) :: F integer, intent(out) :: status ! --- const -------------------------- character(len=*), parameter :: rname = mname//'/gribfile_ReadRecord' integer :: iret logical :: verbose ! --- begin -------------------------- ! old message loaded ? if ( F%igrib > 0 ) then ! remove from memory: call GRIB_Release( F%igrib, status=status ) IF_GRIB_NOTOK_RETURN(status=1) ! no message loaded yet: F%igrib = -1 end if ! a new grib message is loaded from file ! igrib is the grib id to be used in subsequent calls ! ret is the return value to check on end_of_file call grib_new_from_file ( F%fu, F%igrib, iret ) if ( iret == GRIB_END_OF_FILE ) then write (gol,'("end-of-file is hit before a GRIB product is read")') ; call goErr write (gol,'(" grib file : ",a)') F%fname ; call goErr TRACEBACK; status = 1; return end if if ( iret /= GRIB_SUCCESS ) then write (gol,'("problem reading GRIB data block")') ; call goErr write (gol,'(" grib file : ",a)') F%fname ; call goErr call grib_get_error_string(iret,gol) call goErr TRACEBACK; status = 2; return end if ! no error ... status = 0 end subroutine gribfile_ReadRecord ! ============================================================ ! Extract data from grib sections. ! Optional arguments are passed to read parameter id, date, etc. ! If the same option is preceded by 'check_', an error ! message is written if the grib file does not match ! ! pid : parameter identifier (integer number) ! ! level : index of hybride level ! ! reftime : 5 element integer array : yy, mm, dd, hh, min ! timerange : 4 element ingeger array : unit, val1, val2, indicator ! subroutine gribfile_Get( F, status, & model_id, pid, pidShortName, & levtype, level, nlev, hyb_a, hyb_b, & reftime, timerange, & gridtype, & ll, & lon_first, lon_last, lon_inc, lon_n, & lat_first, lat_last, lat_inc, lat_n, & T, sh, & N, gg ) use Grib_API, only : grib_get, grib_get_size use file_grib_base, only : gridtype_ll, gridtype_gg, gridtype_sh, gridtype_red_gg use file_grib_base, only : levtype_sfc, levtype_hyb, levtype_land, levtype_hyb2 ! --- in/out ------------------------- type(TGribFile), intent(in) :: F integer, intent(out) :: status integer, intent(out), optional :: model_id integer, intent(out), optional :: pid character(72), intent(out), optional :: pidShortName integer, intent(out), optional :: levtype, level, nlev real, intent(out), optional :: hyb_a(:), hyb_b(:) integer, intent(out), optional :: reftime(5) integer, intent(out), optional :: timerange(4) integer, intent(out), optional :: gridtype real, intent(out), optional :: lon_first, lon_last, lon_inc integer, intent(out), optional :: lon_n real, intent(out), optional :: lat_first, lat_last, lat_inc integer, intent(out), optional :: lat_n real, intent(out), optional :: ll(:,:) integer, intent(out), optional :: T complex, intent(out), optional :: sh(:) integer, intent(out), optional :: N real, intent(out), optional :: gg(:) ! --- const -------------------------------------- character(len=*), parameter :: rname = mname//'/gribfile_Get' ! --- local -------------------------- integer :: grib_reftime(5) integer :: grib_timerange(4) integer :: PVPresent, nb_pv, step, iret, ni, nj integer :: dataDate, dataTime integer :: numberOfValues real, allocatable :: pv(:) real, allocatable :: values(:) character(72) :: gridtypestring integer :: TypeOfLevel, edition ! --- begin -------------------------- ! this grib1/grib2 business is getting messy! call grib_get(F%igrib, 'editionNumber', edition) ! -- model id if ( present( model_id ) ) then call grib_get(F%igrib, 'generatingProcessIdentifier', model_id ) end if ! -- parameter if ( present( pid ) ) then call grib_get(F%igrib, 'paramId', pid ) end if ! -- parameter short name if ( present( pidShortName ) ) then call grib_get(F%igrib, 'shortName', pidShortName ) end if ! -- level if ( present(levtype) ) then call grib_get(F%igrib, 'levtype', levtype ) ! make TM5 code works with edition2: if ((edition==2).and.(levtype==levtype_hyb2)) then levtype=levtype_hyb endif endif if ( present(level ) ) & call grib_get(F%igrib, 'level', level ) ! -- Vertical Coordinate Parameters if ( present(hyb_a) .or. present(hyb_b) ) then if ( .not. (present(hyb_a) .and. present(hyb_b)) ) then write (gol,'("extract both hybrid params or none")'); call goErr TRACEBACK ; status = 1 ; return end if if (edition==1) then call grib_get(F%igrib,'indicatorOfTypeOfLevel', TypeOfLevel) select case (TypeOfLevel) case(levtype_hyb) ! get size, check, read, and split into 'a' and 'b' call grib_get_size(F%igrib,'pv',nb_pv) if ( nb_pv /= size(hyb_a) + size(hyb_b) ) then write (gol,'("numbers of hybride parameters do not match:")'); call goErr write (gol,'(" expected : ",i6)') size(hyb_a)+size(hyb_b); call goErr write (gol,'(" found : ",i6)') nb_pv; call goErr write (gol,'(" file : ",a)') trim(F%fname); call goErr TRACEBACK ; status = 1 ; return end if allocate(pv(nb_pv)) call grib_get(F%igrib,'pv',pv) hyb_a = pv(1:size(hyb_a)) hyb_b = pv(size(hyb_a)+1:nb_pv) deallocate(pv) ! case(levtype_sfc) ! !surface ! nlev=1 case default write (gol,'("there are no hybride parameters in the grib message")'); call goErr write (gol,'(" file : ",a)') trim(F%fname); call goErr TRACEBACK ; status = 1 ; return end select elseif (edition==2) then call grib_get(F%igrib,'PVPresent',PVPresent) if (PVPresent == 1) then ! get size, check, read, and split into 'a' and 'b' call grib_get_size(F%igrib,'pv',nb_pv) if ( nb_pv /= size(hyb_a) + size(hyb_b) ) then write (gol,'("numbers of hybride parameters do not match:")'); call goErr write (gol,'(" expected : ",i6)') size(hyb_a)+size(hyb_b); call goErr write (gol,'(" found : ",i6)') nb_pv; call goErr write (gol,'(" file : ",a)') trim(F%fname); call goErr TRACEBACK ; status = 1 ; return end if allocate(pv(nb_pv)) call grib_get(F%igrib,'pv',pv) hyb_a = pv(1:size(hyb_a)) hyb_b = pv(size(hyb_a)+1:nb_pv) deallocate(pv) else write (gol,'("there are no hybride parameters in the grib message")'); call goErr write (gol,'(" file : ",a)') trim(F%fname); call goErr TRACEBACK ; status = 1 ; return end if end if !PLS: PVPresent flag is not available anymore. (for ei in grib1 I think) !PLS !PLS ! set PVPresent as an integer !PLS call grib_get(F%igrib,'PVPresent',PVPresent) !PLS if (PVPresent == 1) then !PLS ! number of parameters: !PLS call grib_get_size(F%igrib,'pv',nb_pv) !PLS ! parameters are stored a(0),..,a(L),b(0),..,b(L), !PLS ! thus length should be combined size of 'a' and 'b' arrays: !PLS if ( nb_pv /= size(hyb_a) + size(hyb_b) ) then !PLS write (gol,'("numbers of hybride parameters do not match:")'); call goErr !PLS write (gol,'(" expected : ",i6)') size(hyb_a)+size(hyb_b); call goErr !PLS write (gol,'(" found : ",i6)') nb_pv; call goErr !PLS write (gol,'(" file : ",a)') trim(F%fname); call goErr !PLS TRACEBACK ; status = 1 ; return !PLS end if !PLS ! storage for all parameters: !PLS allocate(pv(nb_pv)) !PLS ! read: !PLS call grib_get(F%igrib,'pv',pv) !PLS ! split into 'a' and 'b' coeff: !PLS hyb_a = pv(1:size(hyb_a)) !PLS hyb_b = pv(size(hyb_a)+1:nb_pv) !PLS ! clear: !PLS deallocate(pv) !PLS else !PLS write (gol,'("there are no hybride parameters in the grib message")'); call goErr !PLS write (gol,'(" file : ",a)') trim(F%fname); call goErr !PLS TRACEBACK ; status = 1 ; return !PLS end if end if ! number of levels if ( present(nlev) ) then if (edition==1) then call grib_get(F%igrib,'indicatorOfTypeOfLevel', TypeOfLevel) select case (TypeOfLevel) case(levtype_hyb) call grib_get_size(F%igrib,'pv',nb_pv) ! half level coeffs are hyba(0:nlev) and hybb(0:nlev), compute number of (full) levels: nlev = nb_pv/2 - 1 case(levtype_sfc) !surface nlev=1 case(levtype_land) nlev=1 case default write (gol,'("do not know about that TypeOfLevel")'); call goErr TRACEBACK ; status = 1 ; return end select elseif (edition==2) then call grib_get(F%igrib,'PVPresent',PVPresent) if (PVPresent == 1) then ! get total number of half level hyb. coeffs: call grib_get_size(F%igrib,'pv',nb_pv) ! half level coeffs are hyba(0:nlev) and hybb(0:nlev), compute number of (full) levels: nlev = nb_pv/2 - 1 else write (gol,'("there are no hybrid parameters in the grib message")'); call goErr TRACEBACK ; status = 1 ; return end if end if endif ! -- reference time if ( present(reftime) ) then ! extract reference time: call grib_get(F%igrib,'dataDate',dataDate) call grib_get(F%igrib,'dataTime',dataTime) grib_reftime(1) = dataDate/10000 ! year grib_reftime(2) = Mod(dataDate/100,100) ! month grib_reftime(3) = Mod(datadate,100) ! day grib_reftime(4) = dataTime/100 ! hour grib_reftime(5) = Mod(dataTime,100) ! minutes reftime = grib_reftime end if ! -- time range if ( present(timerange) ) then ! write (*,'(" ",a," WARNING timerange argument not fully supported ")') trim(name) ! Available from Grib-api: stepType, stepUnits, startStep, endStep, stepRange, step call grib_get(F%igrib,'step',step) grib_timerange(1:4) = (/ 1, step, 0, 0 /) timerange = grib_timerange end if ! --- grid if ( present(gridtype) ) then call grib_get(F%igrib,'gridType',gridtypestring) select case ( trim(gridtypestring) ) case ( 'regular_ll' ) ; gridtype = gridtype_ll case ( 'regular_gg' ) ; gridtype = gridtype_gg ! gribex gives that for reduced. case ( 'regular_sh' ) ; gridtype = gridtype_sh case ( 'sh' ) ; gridtype = gridtype_sh ! used by grid_api case ( 'reduced_gg' ) ; gridtype = gridtype_red_gg ! used by grid_api case default write (gol,'("unsupported grid type string : ",a)') trim(gridtypestring); call goErr TRACEBACK ; status = 1 ; return end select end if if ( present(lon_n) ) call grib_get(F%igrib,'Ni',lon_n) if ( present(lat_n) ) call grib_get(F%igrib,'Nj',lat_n) if ( present(lat_first) ) call grib_get(F%igrib,'latitudeOfFirstGridPointInDegrees',lat_first) if ( present(lon_first) ) call grib_get(F%igrib,'longitudeOfFirstGridPointInDegrees',lon_first) if ( present(lat_last) ) call grib_get(F%igrib,'latitudeOfLastGridPointInDegrees',lat_last) if ( present(lon_last) ) call grib_get(F%igrib,'longitudeOfLastGridPointInDegrees',lon_last) if ( present(lon_inc) ) call grib_get(F%igrib,'iDirectionIncrementInDegrees',lon_inc) if ( present(lat_inc) ) call grib_get(F%igrib,'jDirectionIncrementInDegrees',lat_inc) ! -- lat/lon grid if ( present(ll) ) then !call Check( F, gridtype=gridtype_ll ) call grib_get(F%igrib,'gridType',gridtypestring) if ( trim(gridtypestring) /= 'regular_ll' ) then write (gol,'("grid types do not match:")'); call goErr write (gol,'(" expected : regular_ll")') ; call goErr write (gol,'(" found : ",i4)') gridtype; call goErr TRACEBACK ; status = 1 ; return end if call grib_get(F%igrib,'Ni',ni) call grib_get(F%igrib,'Nj',nj) if ( (size(ll,1)/=ni) .or. (size(ll,2)/=nj) ) then write (gol,'("grid sizes do not match:")'); call goErr write (gol,'(" expected : ",i4," x ",i4)') size(ll,1), size(ll,2); call goErr write (gol,'(" found : ",i4," x ",i4)') ni, nj; call goErr TRACEBACK ; status = 1 ; return end if ! get the size of the values array, and allocate storage call grib_get_size(F%igrib,'values',numberOfValues) allocate(values(numberOfValues), stat=iret) ! get the field call grib_get(F%igrib,'values',values) ll = reshape( values(1:size(ll)), shape(ll) ) ! release storage deallocate(values) end if ! -- gaussian lat/lon grid if ( present(N) ) then call grib_get(F%igrib,'N',N) end if if ( present(gg) ) then ! Get and check size of field call grib_get_size(F%igrib,'values',numberOfValues) if ( size(gg) /= numberOfValues ) then write (gol,'("gg grid sizes do not match:")') ; call goErr write (gol,'(" expected : ",i6)') size(gg) ; call goErr write (gol,'(" found : ",i6)') numberOfValues ; call goErr TRACEBACK; status=1; return end if ! Get the field call grib_get(F%igrib,'values',gg) end if ! -- spectral coef if ( present(T) ) then ! truncation nb call grib_get(F%igrib,'M',T) end if if ( present(sh) ) then ! Get and check size of field call grib_get_size(F%igrib,'values',numberOfValues) if ( size(sh) /= numberOfValues/2 ) then write (gol,'("numbers of spectral coefficients do not match:")') ; call goErr write (gol,'(" expected : ",i6)') size(sh) ; call goErr write (gol,'(" found : ",i6)') numberOfValues/2 ; call goErr TRACEBACK; status=1; return end if ! Get the field, convert from real to complex allocate(values(numberOfValues), stat=iret) call grib_get(F%igrib,'values',values) sh = transfer( values, (/(0.0,0.0)/) ) deallocate( values ) end if ! ok status = 0 end subroutine gribfile_Get ! === ! NOTE ! Do not check hybride coeff; different number of bits etc ! cause small differences ... subroutine gribfile_Check( F, status, debug, & model_id, pid, & levtype, level, & reftime, timerange, & gridtype, & lon_first, lon_last, lon_inc, lon_n, & lat_first, lat_last, lat_inc, lat_n, & T ) use Grib_API, only : grib_get use file_grib_base, only : gridtype_ll, gridtype_gg, gridtype_sh, gridtype_red_gg ! --- in/out ------------------------- type(TGribFile), intent(in) :: F integer, intent(out) :: status integer, intent(in), optional :: debug integer, intent(in), optional :: model_id integer, intent(in), optional :: pid integer, intent(in), optional :: levtype, level !real, intent(in), optional :: hyb_a(:), hyb_b(:) integer, intent(in), optional :: reftime(5) integer, intent(in), optional :: timerange(4) integer, intent(in), optional :: gridtype integer, intent(in), optional :: lon_first, lon_last, lon_inc, lon_n integer, intent(in), optional :: lat_first, lat_last, lat_inc, lat_n integer, intent(in), optional :: T ! --- const -------------------------------------- character(len=*), parameter :: rname = mname//'/gribfile_Check' ! --- local -------------------------- integer :: dataDate, dataTime integer :: grib_model_id integer :: grib_pid integer :: grib_levtype, grib_level integer :: grib_reftime(5) integer :: grib_timerange(4) integer :: grib_n ! integer :: grib_T ! ! !real, allocatable :: grib_hyb_a(:), grib_hyb_b(:) logical :: verbose integer :: step, recGridType character(255) :: gridtypestring real :: xdum ! --- begin -------------------------- ! write error messages ? verbose = .false. if ( present(debug) ) verbose = debug > 0 ! ok by default status = 0 ! -- check model id if ( present( model_id ) ) then call grib_get(F%igrib, 'generatingProcessIdentifier', grib_model_id ) if ( grib_model_id /= model_id ) then if ( verbose ) then write (gol,'("model id''s do not match:")') ; call goErr write (gol,'(" expected : ",i4)') model_id ; call goErr write (gol,'(" found : ",i4)') grib_model_id ; call goErr write (gol,'(" in ",a)') trim(F%fname) ; call goErr TRACEBACK end if status=status+1 end if end if ! -- check parameter if ( present( pid ) ) then call grib_get(F%igrib, 'paramId', grib_pid ) if ( grib_pid /= pid ) then if ( verbose ) then write (gol,'("parameter id''s do not match:")') ; call goErr write (gol,'(" expected : ",i6)') pid ; call goErr write (gol,'(" found : ",i6)') grib_pid ; call goErr write (gol,'(" in ",a)') trim(F%fname) ; call goErr TRACEBACK end if status = status+1 end if end if ! -- check level if ( present(levtype) ) then call grib_get(F%igrib, 'levtype', grib_levtype ) if ( grib_levtype /= levtype ) then if ( verbose ) then write (gol,'("level types do not match:")') ; call goErr write (gol,'(" expected : ",i6)') levtype ; call goErr write (gol,'(" found : ",i6)') grib_levtype ; call goErr write (gol,'(" in ",a)') trim(F%fname) ; call goErr TRACEBACK end if status=status+1 end if end if if ( present(level) ) then call grib_get(F%igrib, 'level', grib_level ) if ( grib_level /= level ) then if ( verbose ) then write (gol,'("levels do not match:")') ; call goErr write (gol,'(" expected : ",i6)') level ; call goErr write (gol,'(" found : ",i6)') grib_level ; call goErr write (gol,'(" in ",a)') trim(F%fname) ; call goErr TRACEBACK end if status=status+1 end if end if !if ( present(hyb_a) .or. present(hyb_b) ) then ! if ( .not. (present(hyb_a) .and. present(hyb_b)) ) then ! print *, 'gribsections_Extract : error : extract both hybride params or none' ! !stop ! end if ! allocate( grib_hyb_a(size(hyb_a)) ) ! allocate( grib_hyb_b(size(hyb_b)) ) ! call Get( F, hyb_a=grib_hyb_a, hyb_b=grib_hyb_b ) ! if ( any(grib_hyb_a/=hyb_a) .or. any(grib_hyb_b/=hyb_b) ) then ! if ( present(status) ) then ! status = 1 ! else ! print *, 'gribsections_Check : hybride coeff in grib file do not match' ! stop ! end if ! end if ! deallocate( grib_hyb_a ) ! deallocate( grib_hyb_b ) !end if ! -- check reference time ! extract reference time: call grib_get(F%igrib,'dataDate',dataDate) call grib_get(F%igrib,'dataTime',dataTime) grib_reftime(1) = dataDate/10000 ! year grib_reftime(2) = Mod(dataDate/100,100) ! month grib_reftime(3) = Mod(datadate,100) ! day grib_reftime(4) = dataTime/100 ! hour grib_reftime(5) = Mod(dataTime,100) ! minutes if ( present(reftime) ) then if ( any( grib_reftime /= reftime ) ) then if ( verbose ) then write (gol,'("reference times do not match:")') ; call goErr write (gol,'(" expected : ",i4,4i3)') reftime ; call goErr write (gol,'(" found : ",i4,4i3)') grib_reftime ; call goErr write (gol,'(" in ",a)') trim(F%fname) ; call goErr TRACEBACK end if status=status+1 end if end if ! -- check time range if ( present(timerange) ) then write (*,'(a," WARNING timerange argument not fully supported ")') rname write(*,*) ! Available from Grib-api: stepType, stepUnits, startStep, endStep, stepRange, step call grib_get(F%igrib,'step',step) grib_timerange(1:4) = (/ 1, step, 0, 0 /) if ( any( grib_timerange /= timerange ) ) then if ( verbose ) then write (gol,'("time ranges do not match:")') ; call goErr write (gol,'(" expected : ",4i3)') timerange ; call goErr write (gol,'(" found : ",4i3)') grib_timerange ; call goErr write (gol,'(" in ",a)') trim(F%fname) ; call goErr TRACEBACK end if status=status+1 end if end if ! --- check grid if ( present(gridtype) ) then call grib_get(F%igrib,'gridType',gridtypestring) ! See https://software.ecmwf.int/wiki/display/GRIB/Grib+API+keys for availability select case ( trim(gridtypestring) ) case ( 'regular_ll' ) ; recGridType = gridtype_ll case ( 'regular_gg' ) ; recGridType = gridtype_gg case ( 'regular_sh' ) ; recGridType = gridtype_sh case ( 'sh' ) ; recGridType = gridtype_sh ! used by grid_api case ( 'reduced_gg' ) ; recgridtype = gridtype_red_gg ! used by grid_api case default if ( verbose ) then write (gol,'("unsupported grid type string : ",a)') trim(gridtypestring); call goErr endif recGridType = -999 end select if ( recGridType /= gridtype ) then if ( verbose ) then write (gol,'("gridtypes do not match:")') ; call goErr write (gol,'(" expected : ",i6)') gridtype ; call goErr write (gol,'(" found : ",i6)') recGridType ; call goErr write (gol,'(" in ",a)') trim(F%fname) ; call goErr TRACEBACK end if status=status+1 end if end if if ( present(lon_n) ) then call grib_get(F%igrib,'Ni',grib_n) if ( grib_n /= lon_n ) then if ( verbose ) then write (gol,'("number of longitudes do not match:")') ; call goErr write (gol,'(" expected : ",i6)') lon_n ; call goErr write (gol,'(" found : ",i6)') grib_n ; call goErr write (gol,'(" in ",a)') trim(F%fname) ; call goErr TRACEBACK end if status=status+1 end if end if if ( present(lat_n) ) then call grib_get(F%igrib,'Nj',grib_n) if ( grib_n /= lat_n ) then if ( verbose ) then write (gol,'("number of latitudes do not match:")') ; call goErr write (gol,'(" expected : ",i6)') lat_n ; call goErr write (gol,'(" found : ",i6)') grib_n ; call goErr write (gol,'(" in ",a)') trim(F%fname) ; call goErr TRACEBACK end if status=status+1 end if end if if ( present(lat_first) ) then call grib_get(F%igrib,'latitudeOfFirstGridPointInDegrees',xdum) grib_n = nint(xdum*1000) if ( grib_n /= lat_first ) then if ( verbose ) then write (gol,'("first latitudes do not match:")') ; call goErr write (gol,'(" expected : ",i6)') lat_first ; call goErr write (gol,'(" found : ",i6)') grib_n ; call goErr write (gol,'(" in ",a)') trim(F%fname) ; call goErr TRACEBACK end if status=status+1 end if end if if ( present(lon_first) ) then call grib_get(F%igrib,'longitudeOfFirstGridPointInDegrees',xdum) grib_n = nint(xdum*1000) if ( grib_n /= lon_first ) then if ( verbose ) then write (gol,'("first longitudes do not match:")') ; call goErr write (gol,'(" expected : ",i6)') lon_first ; call goErr write (gol,'(" found : ",i6)') grib_n ; call goErr write (gol,'(" in ",a)') trim(F%fname) ; call goErr TRACEBACK end if status=status+1 end if end if if ( present(lat_last) ) then call grib_get(F%igrib,'latitudeOfLastGridPointInDegrees',xdum) grib_n = nint(xdum*1000) if ( grib_n /= lat_last ) then if ( verbose ) then write (gol,'("last latitudes do not match:")') ; call goErr write (gol,'(" expected : ",i6)') lat_last ; call goErr write (gol,'(" found : ",i6)') grib_n ; call goErr write (gol,'(" in ",a)') trim(F%fname) ; call goErr TRACEBACK end if status=status+1 end if end if if ( present(lon_last) ) then call grib_get(F%igrib,'longitudeOfLastGridPointInDegrees',xdum) grib_n = nint(xdum*1000) if ( grib_n /= lon_last ) then if ( verbose ) then write (gol,'("last longitudes do not match:")') ; call goErr write (gol,'(" expected : ",i6)') lon_last ; call goErr write (gol,'(" found : ",i6)') grib_n ; call goErr write (gol,'(" in ",a)') trim(F%fname) ; call goErr TRACEBACK end if status=status+1 end if end if if ( present(lon_inc) ) then call grib_get(F%igrib,'iDirectionIncrementInDegrees',xdum) grib_n = nint(xdum*1000) if ( grib_n /= lon_inc ) then if ( verbose ) then write (gol,'("longitude increments do not match:")') ; call goErr write (gol,'(" expected : ",i6)') lon_inc ; call goErr write (gol,'(" found : ",i6)') grib_n ; call goErr write (gol,'(" in ",a)') trim(F%fname) ; call goErr TRACEBACK end if status=status+1 end if end if if ( present(lat_inc) ) then call grib_get(F%igrib,'jDirectionIncrementInDegrees',xdum) grib_n = nint(xdum*1000) if ( grib_n /= lat_inc ) then if ( verbose ) then write (gol,'("latitude increments do not match:")') ; call goErr write (gol,'(" expected : ",i6)') lat_inc ; call goErr write (gol,'(" found : ",i6)') grib_n ; call goErr write (gol,'(" in ",a)') trim(F%fname) ; call goErr TRACEBACK end if status=status+1 end if end if ! --- check truncation if ( present(T) ) then write (*,'("ERROR, Spectral grid not yet supported for Grib-api")') stop end if end subroutine gribfile_Check end module file_grib_api