program getfield
implicit none

type selector
   character :: name*20
   integer   :: grib(3)
end type selector

type(selector), parameter :: selection(*) = [         &
   selector('orography'           , [  6, 105,   0]), &
   selector('z500'                , [  6, 100, 500]), &
   selector('pmsl'                , [  1, 103,   0]), &
   selector('landseamask'         , [ 81, 105,   0]), &
   selector('cloudcover'          , [ 71, 105,   0]), &
   selector('latentheatflux'      , [121, 105, 904]), &
   selector('sensibleheatflux'    , [122, 105, 904]), &
   selector('cape'                , [225, 105,  32]), &
   selector('sneeuwhoogte'        , [138, 105,   0]), &
   selector('zicht'               , [ 20, 105,   0]), &
   selector('temperatuur'         , [ 11, 105,   2]), &
   selector('windstoten'          , [228, 105,  10]), &
   selector('wolkentoptemperatuur', [131, 105,   0]), & ! cloud top temperature (infrared)
   selector('waterdamp'           , [132, 105,   0]), & ! water vapour brightness temperature
   selector('waterdamp+wolkencorr', [133, 105,   0]), & ! water vapour br. temp. + correction for clouds
   selector('wolkenwater'         , [134, 105,   0])  & ! cloud water reflectivity (visible)
]
integer, parameter :: jplb1 = 30, jplb2 = 400, jpsup = 40, jplf = 1600000

character :: filename*256, fieldname*20
real      :: zbs2(jplb2), zsup(jpsup), field(jplf)
integer   :: ibs1(jplb1), ibs2(jplb2)
integer   :: ibpw = 0, ibpwio = 0, imdi
integer   :: i, j, ilf, ifld, ierr, nlon, nlat

call get_command_argument(1, filename)
call get_command_argument(2, fieldname)

ifld = findloc(selection % name, fieldname, dim=1)
if (ifld <= 0 .or. ifld > size(selection)) then
   write(0,*) 'GETFIELD: UNKNOWN FIELD NAME: ',fieldname
   call abort
endif

call ASIMHM(ibpw, ibpwio, imdi)

call ASIMHO(1, filename, 0, ierr)
if (ierr /= 0) then
   write(0,*) 'GETFIELD: CANNOT OPEN ',trim(filename),' ASIMHO RETURNS:', ierr
   call abort
endif

ibs1 = imdi
ibs1(9:11) = selection(ifld) % grib

ilf = jplf
call GETFD(1, ibs1, jplb1, ibs2, jplb2, zbs2, jplb2, zsup, jpsup, field, ilf, 0, ierr)
if (ierr < 0) then
   write(0,*) 'GETFIELD: ERROR IN GETFD; IERR=',ierr
   write(0,'(30(1x,i3))') ibs1
   call abort
endif

call ASIMHC(1, ierr)

nlon = ibs2(7)
nlat = ibs2(9)

do i = 1, nlat
   do j = 1, nlon
      print*, j, i, field((i-1)*nlon+j)
   enddo
   print*
enddo

end program getfield