gridinfo_class.F90

! Copyright (C) 2010  ARPA-SIM <urpsim@smr.arpa.emr.it>
! authors:
! Davide Cesari <dcesari@arpa.emr.it>
! Paolo Patruno <ppatruno@arpa.emr.it>

! This program is free software; you can redistribute it and/or
! modify it under the terms of the GNU General Public License as
! published by the Free Software Foundation; either version 2 of 
! the License, or (at your option) any later version.

! This program is distributed in the hope that it will be useful,
! but WITHOUT ANY WARRANTY; without even the implied warranty of
! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
! GNU General Public License for more details.

! You should have received a copy of the GNU General Public License
! along with this program.  If not, see <http://www.gnu.org/licenses/>.
#include "config.h"

MODULE gridinfo_class

USE grid_class
USE datetime_class
USE vol7d_timerange_class
USE vol7d_level_class
USE volgrid6d_var_class
#ifdef HAVE_LIBGRIBAPI
USE grib_api
#endif
#ifdef HAVE_LIBGDAL
USE gdal
#endif
USE grid_id_class
use log4fortran
use optional_values


IMPLICIT NONE

character (len=255),parameter:: subcategory="gridinfo_class"


TYPE gridinfo_def
  TYPE(griddim_def) :: griddim
  TYPE(datetime) :: time
  TYPE(vol7d_timerange) :: timerange
  TYPE(vol7d_level) :: level
  TYPE(volgrid6d_var) :: var
  TYPE(grid_id) :: gaid
  TYPE(l4f_handle) :: category
END TYPE gridinfo_def

INTEGER, PARAMETER :: &
 cosmo_centre(3) = (/78,80,200/), & ! emission centres using COSMO coding
 ecmwf_centre(1) = (/98/) ! emission centres using ECMWF coding

INTERFACE init
  MODULE PROCEDURE gridinfo_init
END INTERFACE

INTERFACE delete
  MODULE PROCEDURE gridinfo_delete
END INTERFACE

INTERFACE clone
  MODULE PROCEDURE gridinfo_clone
END INTERFACE

INTERFACE import
  MODULE PROCEDURE gridinfo_import, gridinfo_import_from_file
!  MODULE PROCEDURE import_time,import_timerange,import_level,import_gridinfo, &
!   import_volgrid6d_var,gridinfo_import_from_file
END INTERFACE

INTERFACE export
  MODULE PROCEDURE gridinfo_export, gridinfo_export_to_file
!  MODULE PROCEDURE export_time,export_timerange,export_level,export_gridinfo, &
!   export_volgrid6d_var
END INTERFACE

INTERFACE display
  MODULE PROCEDURE gridinfo_display, gridinfov_display
END INTERFACE

INTERFACE decode_gridinfo
  MODULE PROCEDURE gridinfo_decode_data
END INTERFACE

INTERFACE encode_gridinfo
  MODULE PROCEDURE gridinfo_encode_data
END INTERFACE

#define ARRAYOF_ORIGTYPE TYPE(gridinfo_def)
#define ARRAYOF_TYPE arrayof_gridinfo
#define ARRAYOF_ORIGDESTRUCTOR(x) CALL delete(x)
#include "arrayof_pre.F90"
! from arrayof
PUBLIC insert, append, remove, packarray

PRIVATE
PUBLIC gridinfo_def,init,delete,import,export,clone
PUBLIC display,decode_gridinfo,encode_gridinfo

CONTAINS

#include "arrayof_post.F90"

SUBROUTINE gridinfo_init(this, gaid, griddim, time, timerange, level, var, &
 clone, categoryappend)
TYPE(gridinfo_def),intent(out) :: this
TYPE(grid_id),intent(in),optional :: gaid
type(griddim_def),intent(in),optional :: griddim
TYPE(datetime),intent(in),optional :: time
TYPE(vol7d_timerange),intent(in),optional :: timerange
TYPE(vol7d_level),intent(in),optional :: level
TYPE(volgrid6d_var),intent(in),optional :: var
logical , intent(in),optional :: clone
character(len=*),INTENT(in),OPTIONAL :: categoryappend

character(len=512) :: a_name

if (present(categoryappend))then
   call l4f_launcher(a_name,a_name_append=trim(subcategory)//"."//trim(categoryappend))
else
   call l4f_launcher(a_name,a_name_append=trim(subcategory))
end if
this%category=l4f_category_get_handle(a_name)

#ifdef DEBUG
call l4f_category_log(this%category,l4f_debug,"start init gridinfo")
#endif

if (present(gaid))then
  if (optio_log(clone))then
    CALL copy(gaid,this%gaid)
  else
    this%gaid = gaid
  endif
else
  this%gaid = grid_id_new()
end if

!#ifdef DEBUG
!call l4f_category_log(this%category,L4F_DEBUG,"gaid present: "&
! //to_char(c_e(this%gaid))//" value: "//to_char(this%gaid))
!#endif

if (present(griddim))then
  call copy(griddim,this%griddim)
else
  call init(this%griddim,categoryappend=categoryappend)
end if

if (present(time))then
  this%time=time
else
  call init(this%time)
end if

if (present(timerange))then
  this%timerange=timerange
else
  call init(this%timerange)
end if

if (present(level))then
  this%level=level
else
  call init(this%level)
end if

if (present(var))then
  this%var=var
else
  call init(this%var)
end if

END SUBROUTINE gridinfo_init


SUBROUTINE gridinfo_delete(this)
TYPE(gridinfo_def),intent(inout) :: this

#ifdef DEBUG
call l4f_category_log(this%category,l4f_debug,"start delete_gridinfo" )
#endif

call delete(this%griddim)
call delete(this%time)
call delete(this%timerange)
call delete(this%level)
call delete(this%var)

#ifdef DEBUG
call l4f_category_log(this%category,l4f_debug,"delete gaid" )
#endif
CALL delete(this%gaid)

!chiudo il logger
call l4f_category_delete(this%category)

END SUBROUTINE gridinfo_delete


SUBROUTINE gridinfo_display(this, namespace)
TYPE(gridinfo_def),intent(in) :: this
CHARACTER (len=*),OPTIONAL :: namespace

#ifdef DEBUG
call l4f_category_log(this%category,l4f_debug,"displaying gridinfo " )
#endif

print*,"----------------------- gridinfo display ---------------------"
call display(this%griddim)
call display(this%time)
call display(this%timerange)
call display(this%level)
call display(this%var)
call display(this%gaid, namespace=namespace)
print*,"--------------------------------------------------------------"

END SUBROUTINE gridinfo_display

SUBROUTINE gridinfov_display(this, namespace)
TYPE(arrayof_gridinfo),INTENT(in) :: this
CHARACTER (len=*),OPTIONAL :: namespace

INTEGER :: i

print*,"----------------------- gridinfo array -----------------------"

DO i = 1, this%arraysize

#ifdef DEBUG
  CALL l4f_category_log(this%array(i)%category,l4f_debug, &
   "displaying gridinfo array, element "//t2c(i))
#endif
  CALL display(this%array(i), namespace)

ENDDO
print*,"--------------------------------------------------------------"

END SUBROUTINE gridinfov_display


SUBROUTINE gridinfo_clone(this, that, categoryappend)
TYPE(gridinfo_def),INTENT(in) :: this
TYPE(gridinfo_def),INTENT(out) :: that
CHARACTER(len=*),INTENT(in),OPTIONAL :: categoryappend

CALL init(that, gaid=this%gaid, griddim=this%griddim, time=this%time, &
 timerange=this%timerange, level=this%level, var=this%var, clone=.true., &
 categoryappend=categoryappend)

END SUBROUTINE gridinfo_clone


SUBROUTINE gridinfo_import(this)
TYPE(gridinfo_def),INTENT(inout) :: this

#ifdef HAVE_LIBGRIBAPI
INTEGER :: gaid
#endif
#ifdef HAVE_LIBGDAL
TYPE(gdalrasterbandh) :: gdalid
#endif

#ifdef DEBUG
call l4f_category_log(this%category,l4f_debug,"import from gaid")
#endif

! griddim is imported separately in grid_class
CALL import(this%griddim, this%gaid)

#ifdef HAVE_LIBGRIBAPI
gaid = grid_id_get_gaid(this%gaid)
IF (c_e(gaid)) CALL gridinfo_import_gribapi(this, gaid)
#endif
#ifdef HAVE_LIBGDAL
gdalid = grid_id_get_gdalid(this%gaid)
IF (gdalassociated(gdalid)) CALL gridinfo_import_gdal(this, gdalid)
#endif

END SUBROUTINE gridinfo_import


SUBROUTINE gridinfo_import_from_file(this, filename, categoryappend)
TYPE(arrayof_gridinfo) :: this
CHARACTER(len=*),INTENT(in) :: filename
CHARACTER(len=*),INTENT(in),OPTIONAL :: categoryappend

type(gridinfo_def) :: gridinfol
INTEGER :: ngrid
TYPE(l4f_handle) :: category
CHARACTER(len=512) :: a_name
TYPE(grid_file_id) :: input_file
TYPE(grid_id) :: input_grid

IF (PRESENT(categoryappend)) THEN
  CALL l4f_launcher(a_name,a_name_append= &
   trim(subcategory)//"."//trim(categoryappend))
ELSE
  CALL l4f_launcher(a_name,a_name_append=trim(subcategory))
ENDIF
category=l4f_category_get_handle(a_name)

#ifdef DEBUG
CALL l4f_category_log(category,l4f_debug,"import from file")
#endif

input_file = grid_file_id_new(filename, 'r')

ngrid = 0
!$OMP PARALLEL DEFAULT(SHARED)
!$OMP MASTER
DO WHILE(.true.)
  input_grid = grid_id_new(input_file)
  IF (.NOT. c_e(input_grid)) EXIT

  CALL l4f_category_log(category,l4f_info,"import gridinfo")
  ngrid = ngrid + 1
!$OMP TASK FIRSTPRIVATE(input_grid, ngrid), PRIVATE(gridinfol)
  IF (PRESENT(categoryappend)) THEN
    CALL init(gridinfol, gaid=input_grid, &
     categoryappend=trim(categoryappend)//"-msg"//trim(to_char(ngrid)))
  ELSE
    CALL init(gridinfol, gaid=input_grid, &
     categoryappend="msg"//trim(to_char(ngrid)))
  ENDIF
  CALL import(gridinfol)
!$OMP CRITICAL
  CALL insert(this, gridinfol)
! gridinfol is intentionally not destroyed, since now it lives into this
!$OMP END CRITICAL
!$OMP END TASK
ENDDO
!$OMP END MASTER
!$OMP END PARALLEL

CALL packarray(this)

CALL l4f_category_log(category,l4f_info, &
 "gridinfo_import, "//t2c(ngrid)//" messages/bands imported from file "// &
 trim(filename))

! close file
CALL delete(input_file)
! close logger
CALL l4f_category_delete(category)

END SUBROUTINE gridinfo_import_from_file


SUBROUTINE gridinfo_export(this)
TYPE(gridinfo_def),INTENT(inout) :: this

#ifdef HAVE_LIBGRIBAPI
INTEGER :: gaid
#endif
#ifdef HAVE_LIBGDAL
!TYPE(gdalrasterbandh) :: gdalid
#endif

#ifdef DEBUG
call l4f_category_log(this%category,l4f_debug,"export to gaid" )
#endif

! griddim is exported separately in grid_class
CALL export(this%griddim, this%gaid)

#ifdef HAVE_LIBGRIBAPI
IF (grid_id_get_driver(this%gaid) == 'grib_api') THEN
  gaid = grid_id_get_gaid(this%gaid)
  IF (c_e(gaid)) CALL gridinfo_export_gribapi(this, gaid)
ENDIF
#endif
#ifdef HAVE_LIBGDAL
IF (grid_id_get_driver(this%gaid) == 'gdal') THEN
!gdalid = grid_id_get_gdalid(this%gaid)
  CALL l4f_category_log(this%category,l4f_warn,"export to gdal not implemented" )
ENDIF
#endif

END SUBROUTINE gridinfo_export


SUBROUTINE gridinfo_export_to_file(this, filename, categoryappend)
TYPE(arrayof_gridinfo) :: this
CHARACTER(len=*),INTENT(in) :: filename
CHARACTER(len=*),INTENT(in),OPTIONAL :: categoryappend

INTEGER :: i
TYPE(l4f_handle) :: category
CHARACTER(len=512) :: a_name
TYPE(grid_file_id) :: output_file
TYPE(grid_id) :: valid_grid_id

IF (PRESENT(categoryappend)) THEN
  CALL l4f_launcher(a_name,a_name_append= &
   trim(subcategory)//"."//trim(categoryappend))
ELSE
  CALL l4f_launcher(a_name,a_name_append=trim(subcategory))
ENDIF
category=l4f_category_get_handle(a_name)

#ifdef DEBUG
CALL l4f_category_log(category,l4f_debug, &
 "exporting to file "//trim(filename)//" "//t2c(this%arraysize)//" fields")
#endif

valid_grid_id = grid_id_new()
DO i = 1, this%arraysize ! find a valid grid_id in this
  IF (c_e(this%array(i)%gaid)) THEN
    valid_grid_id = this%array(i)%gaid
    EXIT
  ENDIF
ENDDO

IF (c_e(valid_grid_id)) THEN ! a valid grid_id has been found
! open file
  output_file = grid_file_id_new(filename, 'w', from_grid_id=valid_grid_id)
  IF (c_e(output_file)) THEN
!$OMP PARALLEL DEFAULT(SHARED)
!$OMP MASTER
    DO i = 1, this%arraysize
!$OMP TASK FIRSTPRIVATE(i)
      CALL export(this%array(i)) ! export information to gaid
!$OMP CRITICAL
      CALL export(this%array(i)%gaid, output_file) ! export gaid to file
!$OMP END CRITICAL
!$OMP END TASK
    ENDDO
!$OMP END MASTER
!$OMP END PARALLEL
! close file
    CALL delete(output_file)
  ELSE
    CALL l4f_category_log(category,l4f_error,"opening file "//trim(filename))
    CALL raise_error()
  ENDIF
ELSE ! no valid grid_id has been found
  CALL l4f_category_log(category,l4f_error, &
   "gridinfo object of size "//t2c(this%arraysize))
  CALL l4f_category_log(category,l4f_error, &
   "no valid grid id found when exporting to file "//trim(filename))
  CALL raise_error()
ENDIF

!chiudo il logger
CALL l4f_category_delete(category)

END SUBROUTINE gridinfo_export_to_file


FUNCTION gridinfo_decode_data(this) RESULT(field)
TYPE(gridinfo_def),INTENT(in) :: this
REAL :: field(this%griddim%dim%nx, this%griddim%dim%ny) ! array of decoded values

CALL grid_id_decode_data(this%gaid, field)

END FUNCTION gridinfo_decode_data


SUBROUTINE gridinfo_encode_data(this, field)
TYPE(gridinfo_def),INTENT(inout) :: this
REAL,intent(in) :: field(:,:)

IF (SIZE(field,1) /= this%griddim%dim%nx &
 .OR. SIZE(field,2) /= this%griddim%dim%ny) THEN
  CALL l4f_category_log(this%category,l4f_error, &
   'gridinfo_encode: field and gridinfo object non conformal, field: ' &
   //trim(to_char(SIZE(field,1)))//'X'//trim(to_char(SIZE(field,2)))//', nx,ny:' &
   //trim(to_char(this%griddim%dim%nx))//'X'//trim(to_char(this%griddim%dim%ny)))
  CALL raise_error()
  RETURN
ENDIF

CALL grid_id_encode_data(this%gaid, field)

END SUBROUTINE gridinfo_encode_data


! =========================================
! grib_api driver specific code
! could this be moved to a separate module?
! =========================================
#ifdef HAVE_LIBGRIBAPI
SUBROUTINE gridinfo_import_gribapi(this, gaid)
TYPE(gridinfo_def),INTENT(inout) :: this ! gridinfo object
INTEGER, INTENT(in) :: gaid ! grib_api id of the grib loaded in memory to import

call time_import_gribapi(this%time, gaid)
call timerange_import_gribapi(this%timerange,gaid)
call level_import_gribapi(this%level, gaid)
call var_import_gribapi(this%var, gaid)

call normalize_gridinfo(this)

END SUBROUTINE gridinfo_import_gribapi


! grib_api driver
SUBROUTINE gridinfo_export_gribapi(this, gaid)
TYPE(gridinfo_def),INTENT(inout) :: this ! gridinfo object
INTEGER, INTENT(in) :: gaid ! grib_api id of the grib loaded in memory to export

TYPE(conv_func) :: c_func
REAL,ALLOCATABLE :: tmparr(:,:)

! convert variable and values to the correct edition if required
CALL volgrid6d_var_normalize(this%var, c_func, grid_id_new(grib_api_id=gaid))
IF (this%var == volgrid6d_var_miss) THEN
  CALL l4f_log(l4f_error, &
   'A suitable variable has not been found in table when converting template')
  CALL raise_error()
ENDIF
IF (c_func /= conv_func_miss) THEN ! convert values as well
  tmparr = decode_gridinfo(this) ! f2003 implicit allocation
  CALL compute(c_func, tmparr)
  CALL encode_gridinfo(this, tmparr)
ENDIF

CALL unnormalize_gridinfo(this)

CALL time_export_gribapi(this%time, gaid, this%timerange)
CALL timerange_export_gribapi(this%timerange, gaid, this%time)
CALL level_export_gribapi(this%level, gaid)
CALL var_export_gribapi(this%var, gaid)

END SUBROUTINE gridinfo_export_gribapi


SUBROUTINE time_import_gribapi(this,gaid)
TYPE(datetime),INTENT(out) :: this ! datetime object
INTEGER,INTENT(in) :: gaid ! grib_api id of the grib loaded in memory to import 

INTEGER :: editionnumber, ttimeincr, tprocdata, centre, p2g, p2, unit, status
CHARACTER(len=9) :: date
CHARACTER(len=10) :: time

CALL grib_get(gaid,'GRIBEditionNumber',editionnumber)

IF (editionnumber == 1 .OR. editionnumber == 2) THEN

  CALL grib_get(gaid,'dataDate',date )
  CALL grib_get(gaid,'dataTime',time(:5) )

  CALL init(this,simpledate=date(:8)//time(:4))

  IF (editionnumber == 2) THEN

    CALL grib_get(gaid,'typeOfProcessedData',tprocdata,status)
    CALL grib_get(gaid,'typeOfTimeIncrement',ttimeincr,status)
    IF (ttimeincr == 255) ttimeincr = 2 ! fix some MeteosWiss data
! if analysis-like statistically processed data is encountered, the
! reference time must be shifted to the end of the processing period
    IF (status == grib_success .AND. ttimeincr == 1) THEN
! old libsim convention, to be removed sometime in the future
      CALL grib_get(gaid,'lengthOfTimeRange',p2g)
      CALL grib_get(gaid,'indicatorOfUnitForTimeRange',unit)
      CALL g2_interval_to_second(unit, p2g, p2)
      this = this + timedelta_new(sec=p2)
    ELSE IF (status == grib_success .AND. ttimeincr == 2 .AND. tprocdata == 0) THEN
! generally accepted grib2 convention, DWD exception for cosmo
! "accumulated" analysis is such that reftime points to the end of the
! interval, so no time shift in that case
      CALL grib_get(gaid,'lengthOfTimeRange',p2g)
      CALL grib_get(gaid,'indicatorOfUnitForTimeRange',unit)
      CALL g2_interval_to_second(unit, p2g, p2)
      CALL grib_get(gaid,'centre',centre)
      IF (centre /= 78) THEN
        this = this + timedelta_new(sec=p2)
      ENDIF
    ELSE IF ((status == grib_success .AND. ttimeincr == 2) .OR. &
     status /= grib_success) THEN ! usual case
! do nothing
    ELSE ! valid but unsupported typeOfTimeIncrement
      CALL l4f_log(l4f_error,'typeOfTimeIncrement '//t2c(ttimeincr)// &
       ' not supported')
      CALL raise_error()
    ENDIF
  ENDIF

ELSE
  CALL l4f_log(l4f_error,'GribEditionNumber '//t2c(editionnumber)//' not supported')
  CALL raise_error()

ENDIF

END SUBROUTINE time_import_gribapi


SUBROUTINE time_export_gribapi(this, gaid, timerange)
TYPE(datetime),INTENT(in) :: this ! datetime object
INTEGER,INTENT(in) :: gaid ! grib_api id of the grib loaded in memory to export
TYPE(vol7d_timerange) :: timerange ! timerange, used for grib2 coding of statistically processed analysed data

INTEGER :: editionnumber, centre
CHARACTER(len=8) :: env_var
LOGICAL :: g2cosmo_behavior

CALL grib_get(gaid,'GRIBEditionNumber',editionnumber)

IF (editionnumber == 1) THEN

  CALL code_referencetime(this)

ELSE IF (editionnumber == 2 )THEN

  IF (timerange%p1 >= timerange%p2) THEN ! forecast-like
    CALL code_referencetime(this)
  ELSE IF (timerange%p1 == 0) THEN ! analysis-like
! ready for coding with general convention
    CALL getenv('LIBSIM_G2COSMO_BEHAVIOR', env_var)
    g2cosmo_behavior = len_trim(env_var) > 0
    CALL grib_get(gaid,'centre',centre)
    IF (g2cosmo_behavior .AND. centre == 78) THEN ! DWD analysis exception 
      CALL code_referencetime(this)
    ELSE ! cosmo or old simc convention
      CALL code_referencetime(this-timedelta_new(sec=timerange%p2))
    ENDIF
  ELSE ! bad timerange
    CALL l4f_log( l4f_error, 'Timerange with 0>p1>p2 cannot be exported in grib2')
    CALL raise_error()
  ENDIF

ELSE

  CALL l4f_log(l4f_error,'GribEditionNumber '//t2c(editionnumber)//' not supported')
  CALL raise_error()

ENDIF

CONTAINS

SUBROUTINE code_referencetime(reftime)
TYPE(datetime),INTENT(in) :: reftime

INTEGER :: date,time
CHARACTER(len=17) :: date_time

! datetime is AAAAMMGGhhmmssmsc
CALL getval(reftime, simpledate=date_time)
READ(date_time(:8),'(I8)')date
READ(date_time(9:12),'(I4)')time
CALL grib_set(gaid,'dataDate',date)
CALL grib_set(gaid,'dataTime',time)

END SUBROUTINE code_referencetime

END SUBROUTINE time_export_gribapi


SUBROUTINE level_import_gribapi(this, gaid)
TYPE(vol7d_level),INTENT(out) :: this ! vol7d_level object
INTEGER,INTENT(in) :: gaid ! grib_api id of the grib loaded in memory to import

INTEGER :: EditionNumber,level1,l1,level2,l2
INTEGER :: ltype,ltype1,scalef1,scalev1,ltype2,scalef2,scalev2

call grib_get(gaid,'GRIBEditionNumber',editionnumber)

if (editionnumber == 1)then
  
  call grib_get(gaid,'indicatorOfTypeOfLevel',ltype)
  call grib_get(gaid,'topLevel',l1)
  call grib_get(gaid,'bottomLevel',l2)

  call level_g1_to_g2(ltype,l1,l2,ltype1,scalef1,scalev1,ltype2,scalef2,scalev2)

else if (editionnumber == 2)then

  call grib_get(gaid,'typeOfFirstFixedSurface',ltype1)
  call grib_get(gaid,'scaleFactorOfFirstFixedSurface',scalef1)
  call grib_get(gaid,'scaledValueOfFirstFixedSurface',scalev1)
  IF (scalef1 == -1 .OR. scalev1 == -1) THEN
    scalef1 = imiss; scalev1 = imiss
  ENDIF

  call grib_get(gaid,'typeOfSecondFixedSurface',ltype2)
  call grib_get(gaid,'scaleFactorOfSecondFixedSurface',scalef2)
  call grib_get(gaid,'scaledValueOfSecondFixedSurface',scalev2)
  IF (scalef2 == -1 .OR. scalev2 == -1) THEN
    scalef2 = imiss; scalev2 = imiss
  ENDIF

else

  CALL l4f_log(l4f_error,'GribEditionNumber '//t2c(editionnumber)//' not supported')
  CALL raise_error()

end if

! Convert missing levels and units m -> mm
call level_g2_to_dballe(ltype1,scalef1,scalev1,ltype2,scalef2,scalev2, &
 level1,l1,level2,l2)

call init (this,level1,l1,level2,l2)

END SUBROUTINE level_import_gribapi


SUBROUTINE level_export_gribapi(this, gaid)
TYPE(vol7d_level),INTENT(in) :: this ! vol7d_level object
INTEGER,INTENT(in) :: gaid ! grib_api id of the grib loaded in memory to export

INTEGER :: EditionNumber, ltype1, scalef1, scalev1, ltype2, scalef2, scalev2, &
 ltype, l1, l2

CALL level_dballe_to_g2(this%level1, this%l1, this%level2, this%l2, &
 ltype1, scalef1, scalev1, ltype2, scalef2, scalev2)

call grib_get(gaid,'GRIBEditionNumber',editionnumber)

if (editionnumber == 1)then

  CALL level_g2_to_g1(ltype1,scalef1,scalev1,ltype2,scalef2,scalev2,ltype,l1,l2)

  call grib_set(gaid,'indicatorOfTypeOfLevel',ltype)
! it is important to set topLevel after, otherwise, in case of single levels
! bottomLevel=0 overwrites topLevel (aliases in grib_api)
  call grib_set(gaid,'bottomLevel',l2)
  call grib_set(gaid,'topLevel',l1)

else if (editionnumber == 2)then

  CALL grib_set(gaid,'typeOfFirstFixedSurface',ltype1)
  IF (.NOT.c_e(scalef1) .OR. .NOT.c_e(scalev1)) THEN ! code missing values correctly
    CALL grib_set_missing(gaid,'scaleFactorOfFirstFixedSurface')
    CALL grib_set_missing(gaid,'scaledValueOfFirstFixedSurface')
  ELSE
    CALL grib_set(gaid,'scaleFactorOfFirstFixedSurface',scalef1)
    CALL grib_set(gaid,'scaledValueOfFirstFixedSurface',scalev1)
  ENDIF

  CALL grib_set(gaid,'typeOfSecondFixedSurface',ltype2)
  IF (ltype2 == 255 .OR. .NOT.c_e(scalef2) .OR. .NOT.c_e(scalev2)) THEN ! code missing values correctly
    CALL grib_set_missing(gaid,'scaleFactorOfSecondFixedSurface')
    CALL grib_set_missing(gaid,'scaledValueOfSecondFixedSurface')
  ELSE
    CALL grib_set(gaid,'scaleFactorOfSecondFixedSurface',scalef2)
    CALL grib_set(gaid,'scaledValueOfSecondFixedSurface',scalev2)
  ENDIF

else

  CALL l4f_log(l4f_error,'GribEditionNumber '//t2c(editionnumber)//' not supported')
  CALL raise_error()

end if

END SUBROUTINE level_export_gribapi


SUBROUTINE timerange_import_gribapi(this, gaid)
TYPE(vol7d_timerange),INTENT(out) :: this ! vol7d_timerange object
INTEGER,INTENT(in) :: gaid ! grib_api id of the grib loaded in memory to import

INTEGER :: EditionNumber, tri, unit, p1g, p2g, p1, p2, statproc, &
 ttimeincr, tprocdata, status

call grib_get(gaid,'GRIBEditionNumber',editionnumber)

IF (editionnumber == 1) THEN

  CALL grib_get(gaid,'timeRangeIndicator',tri)
  CALL grib_get(gaid,'P1',p1g)
  CALL grib_get(gaid,'P2',p2g)
  CALL grib_get(gaid,'indicatorOfUnitOfTimeRange',unit)
  CALL timerange_g1_to_v7d(tri, p1g, p2g, unit, statproc, p1, p2)

ELSE IF (editionnumber == 2) THEN
  
  CALL grib_get(gaid,'forecastTime',p1g)
  CALL grib_get(gaid,'indicatorOfUnitOfTimeRange',unit)
  CALL g2_interval_to_second(unit, p1g, p1)
  call grib_get(gaid,'typeOfStatisticalProcessing',statproc,status)

  IF (status == grib_success .AND. statproc >= 0 .AND. statproc < 254) THEN ! statistically processed
    CALL grib_get(gaid,'lengthOfTimeRange',p2g)
    CALL grib_get(gaid,'indicatorOfUnitForTimeRange',unit)
    CALL g2_interval_to_second(unit, p2g, p2)

! for forecast-like timeranges p1 has to be shifted to the end of interval
    CALL grib_get(gaid,'typeOfProcessedData',tprocdata,status)
    CALL grib_get(gaid,'typeOfTimeIncrement',ttimeincr)
    IF (ttimeincr == 2 .AND. tprocdata /= 0) THEN
      p1 = p1 + p2
    ELSE
      IF (p1 > 0) THEN
        CALL l4f_log(l4f_warn,'Found p1>0 in grib2 analysis data, strange things may happen')
      ENDIF
    ENDIF
  ELSE ! point in time
    statproc = 254
    p2 = 0
  
  ENDIF

ELSE

  CALL l4f_log(l4f_error,'GribEditionNumber '//t2c(editionnumber)//' not supported')
  CALL raise_error()

ENDIF

CALL init(this, statproc, p1, p2)

END SUBROUTINE timerange_import_gribapi


SUBROUTINE timerange_export_gribapi(this, gaid, reftime)
TYPE(vol7d_timerange),INTENT(in) :: this ! vol7d_timerange object
INTEGER,INTENT(in) :: gaid ! grib_api id of the grib loaded in memory to export
TYPE(datetime) :: reftime ! reference time of data, used for coding correct end of statistical processing period in grib2

INTEGER :: EditionNumber, centre, tri, currentunit, unit, p1_g1, p2_g1, p1, p2, pdtn
CHARACTER(len=8) :: env_var
LOGICAL :: g2cosmo_behavior

CALL grib_get(gaid,'GRIBEditionNumber',editionnumber)

IF (editionnumber == 1 ) THEN
! Convert vol7d_timerange members to grib1 with reasonable time unit
  CALL grib_get(gaid,'indicatorOfUnitOfTimeRange',currentunit)
  CALL timerange_v7d_to_g1(this%timerange, this%p1, this%p2, &
   tri, p1_g1, p2_g1, unit)
! Set the native keys
  CALL grib_set(gaid,'timeRangeIndicator',tri)
  CALL grib_set(gaid,'P1',p1_g1)
  CALL grib_set(gaid,'P2',p2_g1)
  CALL grib_set(gaid,'indicatorOfUnitOfTimeRange',unit)

ELSE IF (editionnumber == 2) THEN
  CALL grib_get(gaid,'productDefinitionTemplateNumber', pdtn)

  IF (this%timerange == 254) THEN ! point in time -> template 4.0
    IF (pdtn < 0 .OR. pdtn > 7) &
     CALL grib_set(gaid,'productDefinitionTemplateNumber', 0)
! Set reasonable time unit
    CALL timerange_v7d_to_g2(this%p1,p1,unit)
! Set the native keys
    CALL grib_set(gaid,'indicatorOfUnitOfTimeRange',unit)
    CALL grib_set(gaid,'forecastTime',p1)

  ELSE IF (this%timerange >= 0 .AND. this%timerange < 254) THEN
! statistically processed -> template 4.8
    IF (pdtn < 8 .OR. pdtn > 14) &
     CALL grib_set(gaid,'productDefinitionTemplateNumber', 8)

    IF (this%p1 >= this%p2) THEN ! forecast-like
! Set reasonable time unit
      CALL timerange_v7d_to_g2(this%p1-this%p2,p1,unit)
      CALL grib_set(gaid,'indicatorOfUnitOfTimeRange',unit)
      CALL grib_set(gaid,'forecastTime',p1)
      CALL code_endoftimeinterval(reftime+timedelta_new(sec=this%p1))
! Successive times processed have same start time of forecast,
! forecast time is incremented
      CALL grib_set(gaid,'typeOfStatisticalProcessing',this%timerange)
! typeOfTimeIncrement to be replaced with a check that typeOfProcessedData /= 0
      CALL grib_set(gaid,'typeOfTimeIncrement',2)
      CALL timerange_v7d_to_g2(this%p2,p2,unit)
      CALL grib_set(gaid,'indicatorOfUnitForTimeRange',unit)
      CALL grib_set(gaid,'lengthOfTimeRange',p2)

    ELSE IF (this%p1 == 0) THEN ! analysis-like
! Set reasonable time unit
      CALL timerange_v7d_to_g2(this%p2,p2,unit)
      CALL grib_set(gaid,'indicatorOfUnitOfTimeRange',unit)
      CALL grib_set(gaid,'forecastTime',0)
      CALL code_endoftimeinterval(reftime)
! Successive times processed have same forecast time, start time of
! forecast is incremented
      CALL grib_set(gaid,'typeOfStatisticalProcessing',this%timerange)
! typeOfTimeIncrement to be replaced with typeOfProcessedData
      CALL getenv('LIBSIM_G2COSMO_BEHAVIOR', env_var)
      g2cosmo_behavior = len_trim(env_var) > 0
      IF (g2cosmo_behavior) THEN
        CALL grib_set(gaid,'typeOfProcessedData',0)
      ELSE
        CALL grib_set(gaid,'typeOfTimeIncrement',1)
      ENDIF
      CALL grib_set(gaid,'indicatorOfUnitForTimeRange',unit)
      CALL grib_set(gaid,'lengthOfTimeRange',p2)

! warn about local use
      IF (this%timerange >= 192) THEN
        CALL l4f_log(l4f_warn, &
         'coding in grib2 a nonstandard typeOfStatisticalProcessing '// &
         t2c(this%timerange))
      ENDIF
    ELSE ! bad timerange
      CALL l4f_log(l4f_error, &
       'Timerange with 0>p1>p2 cannot be exported in grib2')
      CALL raise_fatal_error()
    ENDIF
  ELSE
    CALL l4f_log(l4f_error, &
     'typeOfStatisticalProcessing not supported: '//trim(to_char(this%timerange)))
    CALL raise_fatal_error()
  ENDIF

ELSE
  CALL l4f_log(l4f_error,'GribEditionNumber '//t2c(editionnumber)//' not supported')
  CALL raise_fatal_error()
ENDIF

CONTAINS

! Explicitely compute and code in grib2 template 4.8 the end of
! overalltimeinterval which is not done automatically by grib_api
SUBROUTINE code_endoftimeinterval(endtime)
TYPE(datetime),INTENT(in) :: endtime

INTEGER :: year, month, day, hour, minute, msec

CALL getval(endtime, year=year, month=month, day=day, &
 hour=hour, minute=minute, msec=msec)
  CALL grib_set(gaid,'yearOfEndOfOverallTimeInterval',year)
  CALL grib_set(gaid,'monthOfEndOfOverallTimeInterval',month)
  CALL grib_set(gaid,'dayOfEndOfOverallTimeInterval',day)
  CALL grib_set(gaid,'hourOfEndOfOverallTimeInterval',hour)
  CALL grib_set(gaid,'minuteOfEndOfOverallTimeInterval',minute)
  CALL grib_set(gaid,'secondOfEndOfOverallTimeInterval',msec/1000)

END SUBROUTINE code_endoftimeinterval

END SUBROUTINE timerange_export_gribapi


SUBROUTINE var_import_gribapi(this, gaid)
TYPE(volgrid6d_var),INTENT(out) :: this ! volgrid6d_var object
INTEGER,INTENT(in) :: gaid ! grib_api id of the grib loaded in memory to import

INTEGER :: EditionNumber, centre, discipline, category, number

call grib_get(gaid,'GRIBEditionNumber',editionnumber)

if (editionnumber == 1) then

  call grib_get(gaid,'centre',centre)
  call grib_get(gaid,'gribTablesVersionNo',category)
  call grib_get(gaid,'indicatorOfParameter',number)

  call init(this, centre, category, number)

else if (editionnumber == 2) then

  call grib_get(gaid,'centre',centre)
  call grib_get(gaid,'discipline',discipline)
  call grib_get(gaid,'parameterCategory',category)
  call grib_get(gaid,'parameterNumber',number)

  call init(this, centre, category, number, discipline)
  
else

  CALL l4f_log(l4f_error,'GribEditionNumber '//t2c(editionnumber)//' not supported')
  CALL raise_error()

endif

END SUBROUTINE var_import_gribapi


SUBROUTINE var_export_gribapi(this, gaid)
TYPE(volgrid6d_var),INTENT(in) :: this ! volgrid6d_var object
INTEGER,INTENT(in) :: gaid ! grib_api id of the grib loaded in memory to export

INTEGER ::EditionNumber

call grib_get(gaid,'GRIBEditionNumber',editionnumber)

if (editionnumber == 1) then

  IF (this%centre /= 255) & ! if centre missing (coming from bufr), keep template
   CALL grib_set(gaid,'centre',this%centre)
  CALL grib_set(gaid,'gribTablesVersionNo',this%category)
  CALL grib_set(gaid,'indicatorOfParameter',this%number)

else if (editionnumber == 2) then

! this must be changed to 65535!!!!
  IF (this%centre /= 255) & ! if centre missing (coming from bufr), keep template
   CALL grib_set(gaid,'centre',this%centre)
  CALL grib_set(gaid,'discipline',this%discipline)
  CALL grib_set(gaid,'parameterCategory',this%category)
  CALL grib_set(gaid,'parameterNumber',this%number)

else

  CALL l4f_log(l4f_error,'GribEditionNumber '//t2c(editionnumber)//' not supported')
  CALL raise_error()

end if

END SUBROUTINE var_export_gribapi


SUBROUTINE level_g2_to_dballe(ltype1,scalef1,scalev1,ltype2,scalef2,scalev2, lt1,l1,lt2,l2)
integer,intent(in) :: ltype1,scalef1,scalev1,ltype2,scalef2,scalev2
integer,intent(out) ::lt1,l1,lt2,l2


CALL g2_to_dballe(ltype1, scalef1, scalev1, lt1, l1)
CALL g2_to_dballe(ltype2, scalef2, scalev2, lt2, l2)

CONTAINS

SUBROUTINE g2_to_dballe(ltype, scalef, scalev, lt, l)
integer,intent(in) :: ltype,scalef,scalev
integer,intent(out) :: lt,l

doubleprecision :: sl


IF (ltype == 255 .OR. ltype == -1) THEN
  lt = imiss
  l = imiss
ELSE IF (ltype <= 10 .OR. ltype == 101 .OR. (ltype >= 162 .AND. ltype <= 184)) THEN
  lt = ltype
  l = imiss
ELSE
  lt = ltype
  IF (c_e(scalef) .AND. c_e(scalev)) THEN
    sl = scalev*(10.d0**(-scalef))
! apply unit conversions
    IF (any(ltype == height_level)) THEN
      l = nint(sl*1000.d0)
    ELSE IF (any(ltype == thermo_level)) THEN
      l = nint(sl*10.d0)
    ELSE IF (any(ltype == sigma_level)) THEN
      l = nint(sl*10000.d0)
    ELSE
      l = nint(sl)
    ENDIF
  ELSE
    l = imiss
  ENDIF
ENDIF

END SUBROUTINE g2_to_dballe

END SUBROUTINE level_g2_to_dballe


SUBROUTINE level_dballe_to_g2(lt1,l1,lt2,l2, ltype1,scalef1,scalev1,ltype2,scalef2,scalev2)
integer,intent(in) :: lt1,l1,lt2,l2
integer,intent(out) :: ltype1,scalef1,scalev1,ltype2,scalef2,scalev2

CALL dballe_to_g2(lt1, l1, ltype1, scalef1, scalev1)
CALL dballe_to_g2(lt2, l2, ltype2, scalef2, scalev2)

CONTAINS

SUBROUTINE dballe_to_g2(lt, l, ltype, scalef, scalev)
INTEGER,INTENT(in) :: lt,l
INTEGER,INTENT(out) :: ltype,scalef,scalev


IF (lt == imiss) THEN
  ltype = 255
  scalev = 0
  scalef = 0
ELSE IF (lt <= 10 .OR. lt == 101 .OR. (lt >= 162 .AND. lt <= 184)) THEN
  ltype = lt
  scalev = 0
  scalef = 0
ELSE IF (lt == 256 .AND. l == imiss) THEN ! special case for cloud level -> surface
  ltype = 1
  scalev = 0
  scalef = 0
ELSE
  ltype = lt
  scalev = l
  IF (any(ltype == height_level)) THEN
    scalef = 3
  ELSE IF (any(ltype == thermo_level)) THEN
    scalef = 1
  ELSE IF (any(ltype == sigma_level)) THEN
    scalef = 4
  ELSE
    scalef = 0
  ENDIF
ENDIF

!Caso generale reale
!IF (ANY(ltype == height_level)) THEN
!  sl=l/1000.D0
!ELSE
!  sl=l
!END IF
!IF (ABS(sl) < PRECISION) THEN
!  scalef = 0
!  scalev = 0
!ELSE
!  magn = LOG10(ABS(sl))
!  DO scalef = magn, MAX(magn, 20)
!    IF (ABS((sl*10.D0**(scalef) - ANINT(sl*10.D0**(scalef))) / &
!     sl*10.D0**(scalef)) < PRECISION) EXIT
!  ENDDO
!  sl = scalev*(10.D0**(-scalef))
!ENDIF

END SUBROUTINE dballe_to_g2

END SUBROUTINE level_dballe_to_g2


SUBROUTINE level_g1_to_g2(ltype,l1,l2,ltype1,scalef1,scalev1,ltype2,scalef2,scalev2)
integer,intent(in) :: ltype,l1,l2
integer,intent(out) :: ltype1,scalef1,scalev1,ltype2,scalef2,scalev2

ltype1=255
scalef1=0
scalev1=0
ltype2=255
scalef2=0
scalev2=0

if (ltype > 0 .and. ltype <= 9)then 
   ltype1=ltype
else if (ltype == 20) then
  ltype1=20
  scalev1=l1
  scalef1=2
else if (ltype == 100) then
  ltype1=100
  scalev1=l1*100
else if (ltype == 101) then
  ltype1=100
  scalev1=l1*1000
  ltype2=100
  scalev2=l2*1000
else if (ltype == 102) then
  ltype1=101
else if (ltype == 103) then
  ltype1=102
  scalev1=l1
else if (ltype == 104) then
  ltype1=102
  scalev1=l1*100
  ltype2=102
  scalev2=l2*100
else if (ltype == 105) then
  ltype1=103
  scalev1=l1
else if (ltype == 106) then
  ltype1=103
  scalev1=l1*100
  ltype2=103
  scalev2=l2*100
else if (ltype == 107) then
  ltype1=104
  scalef1=4
  scalev1=l1
else if (ltype == 108) then
  ltype1=104
  scalef1=2
  scalev1=l1
  ltype2=104
  scalef2=2
  scalev2=l2
else if (ltype == 109) then
  ltype1=105
  scalev1=l1
else if (ltype == 110) then
  ltype1=105
  scalev1=l1
  ltype2=105
  scalev2=l2
else if (ltype == 111) then
  ltype1=106
  scalef1=2
  scalev1=l1
else if (ltype == 112) then
  ltype1=106
  scalef1=2
  scalev1=l1
  ltype2=106
  scalef2=2
  scalev2=l2
else if (ltype == 113) then
  ltype1=107
  scalev1=l1
else if (ltype == 114) then
  ltype1=107
  scalev1=475+l1
  ltype2=107
  scalev2=475+l2
else if (ltype == 115) then
  ltype1=108
  scalev1=l1*100
else if (ltype == 116) then
  ltype1=108
  scalev1=l1*100
  ltype2=108
  scalev2=l2*100
else if (ltype == 117) then
  ltype1=109
  scalef1=9
  scalev1=l1
  if ( btest(l1,15) ) then
    scalev1=-1*mod(l1,32768)
  endif
else if (ltype == 119) then
  ltype1=111
  scalef1=4
  scalev1=l1
else if (ltype == 120) then
  ltype1=111
  scalef1=2
  scalev1=l1
  ltype2=111
  scalef2=2
  scalev2=l2
else if (ltype == 121) then
  ltype1=100
  scalev1=(1100+l1)*100
  ltype2=100
  scalev2=(1100+l2)*100
else if (ltype == 125) then
  ltype1=103
  scalef1=2
  scalev1=l1
else if (ltype == 128) then
  ltype1=104
  scalef1=3
  scalev1=1100+l1
  ltype2=104
  scalef2=3
  scalev2=1100+l2
else if (ltype == 141) then
  ltype1=100
  scalev1=l1*100
  ltype2=100
  scalev2=(1100+l2)*100
else if (ltype == 160) then
  ltype1=160
  scalev1=l1
else if (ltype == 200) then ! entire atmosphere
  ltype1=1
  ltype2=8
else if (ltype == 210) then ! entire ocean
  ltype1=1
  ltype2=9
else

  call l4f_log(l4f_error,'level_g1_to_g2: GRIB1 level '//trim(to_char(ltype)) &
   //' cannot be converted to GRIB2.')
  call raise_error()
endif

END SUBROUTINE level_g1_to_g2


SUBROUTINE level_g2_to_g1(ltype1,scalef1,scalev1,ltype2,scalef2,scalev2,ltype,l1,l2)
integer,intent(in) :: ltype1,scalef1,scalev1,ltype2,scalef2,scalev2
integer,intent(out) :: ltype,l1,l2

if (ltype1 > 0 .and. ltype1 <= 9 .and. ltype2 == 255) then ! simple
  ltype = ltype1
  l1 = 0
  l2 = 0
else if (ltype1 == 20 .and. ltype2 == 255) then ! isothermal
  ltype = 2
  l1 = rescale2(scalef1-2,scalev1)!*100
  l2 = 0
else if (ltype1 == 100 .and. ltype2 == 255) then ! isobaric
  ltype = 100
  l1 = rescale2(scalef1+2,scalev1)!/100
  l2 = 0
else if (ltype1 == 100 .and. ltype2 == 100) then
  ltype = 101
  l1 = rescale1(scalef1+3,scalev1)!/1000
  l2 = rescale1(scalef2+3,scalev2)!/1000
else if (ltype1 == 101 .and. ltype2 == 255) then
  ltype = 102
  l1 = 0
  l2 = 0
else if (ltype1 == 102 .and. ltype2 == 255) then ! altitude over sea level
  ltype = 103
  l1 = rescale2(scalef1,scalev1)
  l2 = 0
else if (ltype1 == 102 .and. ltype2 == 102) then
  ltype = 104
  l1 = rescale1(scalef1+2,scalev1)!/100
  l2 = rescale1(scalef2+2,scalev2)!/100
else if (ltype1 == 103 .and. ltype2 == 255) then ! height over ground
  ltype = 105
  l1 = rescale2(scalef1,scalev1)
  l2 = 0
else if (ltype1 == 103 .and. ltype2 == 103) then
  ltype = 106
  l1 = rescale1(scalef1+2,scalev1)!/100
  l2 = rescale1(scalef2+2,scalev2)!/100
else if (ltype1 == 104 .and. ltype2 == 255) then ! sigma
  ltype = 107
  l1 = rescale2(scalef1,scalev1-4)!*10000
  l2 = 0
else if (ltype1 == 104 .and. ltype2 == 104) then
  ltype = 108
  l1 = rescale1(scalef1-2,scalev1)!*100
  l2 = rescale1(scalef2-2,scalev2)!*100
else if (ltype1 == 105 .and. ltype2 == 255) then ! hybrid
  ltype = 109
  l1 = rescale2(scalef1,scalev1)
  l2 = 0
else if (ltype1 == 105 .and. ltype2 == 105) then
  ltype = 110
  l1 = rescale1(scalef1,scalev1)
  l2 = rescale1(scalef2,scalev2)
else if (ltype1 == 106 .and. ltype2 == 255) then ! depth
  ltype = 111
  l1 = rescale2(scalef1-2,scalev1)!*100
  l2 = 0
else if (ltype1 == 106 .and. ltype2 == 106) then
  ltype = 112
  l1 = rescale1(scalef1-2,scalev1)!*100
  l2 = rescale1(scalef2-2,scalev2)!*100
else if (ltype1 == 107 .and. ltype2 == 255) then ! isentropic
  ltype = 113
  l1 = rescale2(scalef1,scalev1)
  l2 = 0
else if (ltype1 == 107 .and. ltype2 == 107) then
  ltype = 114
  l1 = rescale1(scalef1,scalev1)
  l2 = rescale1(scalef2,scalev2)
else if (ltype1 == 108 .and. ltype2 == 255) then ! pressure diff to ground
  ltype = 115
  l1 = rescale2(scalef1+2,scalev1)!/100
  l2 = 0
else if (ltype1 == 108 .and. ltype2 == 108) then
  ltype = 116
  l1 = rescale1(scalef1+2,scalev1)!/100
  l2 = rescale1(scalef2+2,scalev2)!/100
else if (ltype1 == 109 .and. ltype2 == 255) then ! potential vorticity surf
  ltype = 117
  l1 = rescale2(scalef1+9,scalev1)!/10**9
  l2 = 0
else if (ltype1 == 111 .and. ltype2 == 255) then ! eta level
  ltype = 119
  l1 = rescale2(scalef1-2,scalev1)!*100
  l2 = 0
else if (ltype1 == 111 .and. ltype2 == 111) then
  ltype = 120
  l1 = rescale1(scalef1-4,scalev1)!*10000
  l2 = rescale1(scalef2-4,scalev2)!*10000
else if (ltype1 == 160 .and. ltype2 == 255) then ! depth below sea
  ltype = 160
  l1 = rescale2(scalef1,scalev1)
  l2 = 0
else if (ltype1 == 1 .and. ltype2 == 8) then ! entire atmosphere
  ltype = 200
else if (ltype1 == 1 .and. ltype2 == 9) then ! entire ocean
  ltype = 201
else ! mi sono rotto per ora

  ltype = 255
  l1 = 0
  l2 = 0
  call l4f_log(l4f_error,'level_g2_to_g1: GRIB2 levels '//trim(to_char(ltype1))//' ' &
   //trim(to_char(ltype2))//' cannot be converted to GRIB1.')
  call raise_error()
endif

CONTAINS

FUNCTION rescale1(scalef, scalev) RESULT(rescale)
INTEGER,INTENT(in) :: scalef, scalev
INTEGER :: rescale

rescale = min(255, nint(scalev*10.0d0**(-scalef)))

END FUNCTION rescale1

FUNCTION rescale2(scalef, scalev) RESULT(rescale)
INTEGER,INTENT(in) :: scalef, scalev
INTEGER :: rescale

rescale = min(65535, nint(scalev*10.0d0**(-scalef)))

END FUNCTION rescale2

END SUBROUTINE level_g2_to_g1

! Convert timerange from grib1 to grib2-like way:
!
! Tri 2 (point in time) gives (hopefully temporarily) statproc 205.
!
! Tri 13 (COSMO-nudging) gives p1 (forecast time) 0 and a temporary
! 257 statproc.
!
! Further processing and correction of the values returned is
! performed in normalize_gridinfo.
SUBROUTINE timerange_g1_to_v7d(tri, p1_g1, p2_g1, unit, statproc, p1, p2)
INTEGER,INTENT(in) :: tri, p1_g1, p2_g1, unit
INTEGER,INTENT(out) :: statproc, p1, p2

IF (tri == 0 .OR. tri == 1) THEN ! point in time
  statproc = 254
  CALL g1_interval_to_second(unit, p1_g1, p1)
  p2 = 0
ELSE IF (tri == 10) THEN ! point in time
  statproc = 254
  CALL g1_interval_to_second(unit, p1_g1*256+p2_g1, p1)
  p2 = 0
ELSE IF (tri == 2) THEN ! somewhere between p1 and p2, not good for grib2 standard
  statproc = 205
  CALL g1_interval_to_second(unit, p2_g1, p1)
  CALL g1_interval_to_second(unit, p2_g1-p1_g1, p2)
ELSE IF (tri == 3) THEN ! average
  statproc = 0
  CALL g1_interval_to_second(unit, p2_g1, p1)
  CALL g1_interval_to_second(unit, p2_g1-p1_g1, p2)
ELSE IF (tri == 4) THEN ! accumulation
  statproc = 1
  CALL g1_interval_to_second(unit, p2_g1, p1)
  CALL g1_interval_to_second(unit, p2_g1-p1_g1, p2)
ELSE IF (tri == 5) THEN ! difference
  statproc = 4
  CALL g1_interval_to_second(unit, p2_g1, p1)
  CALL g1_interval_to_second(unit, p2_g1-p1_g1, p2)
ELSE IF (tri == 13) THEN ! COSMO-nudging, use a temporary value then normalize
  statproc = 257 ! check if 257 is legal!
  p1 = 0 ! analysis regardless of p2_g1
  CALL g1_interval_to_second(unit, p2_g1-p1_g1, p2)
ELSE
  call l4f_log(l4f_error,'timerange_g1_to_g2: GRIB1 timerange '//trim(to_char(tri)) &
   //' cannot be converted to GRIB2.')
  CALL raise_error()
ENDIF

if (statproc == 254 .and. p2 /= 0 ) then
  call l4f_log(l4f_warn,"inconsistence in timerange:254,"//trim(to_char(p1))//","//trim(to_char(p2)))
end if

END SUBROUTINE timerange_g1_to_v7d


!0       Minute
!1       Hour
!2       Day
!3       Month
!4       Year
!5       Decade (10 years)
!6       Normal (30 years)
!7       Century(100 years)
!8-9     Reserved
!10      3 hours
!11      6 hours
!12      12 hours
! in COSMO, grib1:
!13 = 15 minuti
!14 = 30 minuti
! in grib2:
!13      Second

SUBROUTINE g1_interval_to_second(unit, valuein, valueout)
INTEGER,INTENT(in) :: unit, valuein
INTEGER,INTENT(out) :: valueout

INTEGER,PARAMETER :: unitlist(0:14)=(/ 60,3600,86400,2592000, &
 31536000,315360000,946080000,imiss,imiss,imiss,10800,21600,43200,900,1800/)

valueout = imiss
IF (unit >= lbound(unitlist,1) .AND. unit <= ubound(unitlist,1)) THEN
  IF (c_e(unitlist(unit))) THEN
    valueout = valuein*unitlist(unit)
  ENDIF
ENDIF

END SUBROUTINE g1_interval_to_second


SUBROUTINE g2_interval_to_second(unit, valuein, valueout)
INTEGER,INTENT(in) :: unit, valuein
INTEGER,INTENT(out) :: valueout

INTEGER,PARAMETER :: unitlist(0:13)=(/ 60,3600,86400,2592000, &
 31536000,315360000,946080000,imiss,imiss,imiss,10800,21600,43200,1/)

valueout = imiss
IF (unit >= lbound(unitlist,1) .AND. unit <= ubound(unitlist,1)) THEN
  IF (c_e(unitlist(unit))) THEN
    valueout = valuein*unitlist(unit)
  ENDIF
ENDIF

END SUBROUTINE g2_interval_to_second


! Convert timerange from grib2-like to grib1 way:
!
! Statproc 205 (point in time, non standard, not good in grib2) is
! correctly converted to tri 2.
!
! Statproc 257 (COSMO nudging-like, non standard, not good in grib2)
! should not appear here, but is anyway converted to tri 13 (real
! COSMO-nudging).
!
! In case p1_g1<0 (i.e. statistically processed analysis data, with
! p1=0 and p2>0), COSMO-nudging tri 13 is (mis-)used.
SUBROUTINE timerange_v7d_to_g1(statproc, p1, p2, tri, p1_g1, p2_g1, unit)
INTEGER,INTENT(in) :: statproc, p1, p2
INTEGER,INTENT(out) :: tri, p1_g1, p2_g1, unit

INTEGER :: ptmp, pdl

pdl = p1 - p2
IF (statproc == 254) pdl = p1 ! avoid unexpected situations (necessary?)

CALL timerange_choose_unit_g1(p1, pdl, p2_g1, p1_g1, unit)
IF (statproc == 0) THEN ! average
  tri = 3
ELSE IF (statproc == 1) THEN ! accumulation
  tri = 4
ELSE IF (statproc == 4) THEN ! difference
  tri = 5
ELSE IF (statproc == 205) THEN ! point in time interval, not good for grib2 standard
  tri = 2
ELSE IF (statproc == 257) THEN ! COSMO-nudging (statistical processing in the past)
! this should never happen (at least from COSMO grib1), since 257 is
! converted to something else in normalize_gridinfo; now a negative
! p1_g1 is set, it will be corrected in the next section
  tri = 13
!  CALL second_to_gribtr(p1, p2_g1, unit, 1)
!  CALL second_to_gribtr(p1-p2, p1_g1, unit, 1)
ELSE IF (statproc == 254) THEN ! point in time
  tri = 0
  p2_g1 = 0
ELSE
  CALL l4f_log(l4f_error,'timerange_v7d_to_g1: GRIB2 statisticalprocessing ' &
   //trim(to_char(statproc))//' cannot be converted to GRIB1.')
  CALL raise_error()
ENDIF

IF (p1_g1 > 255 .OR. p2_g1 > 255) THEN
  ptmp = max(p1_g1,p2_g1)
  p2_g1 = mod(ptmp,256)
  p1_g1 = ptmp/256
  IF (tri /= 0) THEN ! if not instantaneous give warning
    CALL l4f_log(l4f_warn,'timerange_v7d_to_g1: timerange too long for grib1 ' &
     //trim(to_char(ptmp))//', forcing time range indicator to 10.')
  ENDIF
  tri = 10
ENDIF
  

! p1 < 0 is not allowed, use COSMO trick
IF (p1_g1 < 0) THEN
  ptmp = p1_g1
  p1_g1 = 0
  p2_g1 = p2_g1 - ptmp
  tri = 13
ENDIF

END SUBROUTINE timerange_v7d_to_g1


SUBROUTINE timerange_v7d_to_g2(valuein, valueout, unit)
INTEGER,INTENT(in) :: valuein
INTEGER,INTENT(out) :: valueout, unit

IF (valuein == imiss) THEN
  valueout = imiss
  unit = 1
ELSE IF (mod(valuein,3600) == 0) THEN ! prefer hours
  valueout = valuein/3600
  unit = 1
ELSE IF (mod(valuein,60) == 0) THEN ! then minutes
  valueout = valuein/60
  unit = 0
ELSE ! seconds
  valueout = valuein
  unit = 13
ENDIF

END SUBROUTINE timerange_v7d_to_g2


! These units are tested for applicability:
! 0  Minute
! 1  Hour
! 2  Day
! 10 3 hours
! 11 6 hours
! 12 12 hours
SUBROUTINE timerange_choose_unit_g1(valuein1, valuein2, valueout1, valueout2, unit)
INTEGER,INTENT(in) :: valuein1, valuein2
INTEGER,INTENT(out) :: valueout1, valueout2, unit

INTEGER :: i
TYPE unitchecker
  INTEGER :: unit
  INTEGER :: sectounit
END TYPE unitchecker

TYPE(unitchecker),PARAMETER :: hunit(5) = (/ &
 unitchecker(1, 3600), unitchecker(10, 10800), unitchecker(11, 21600), &
 unitchecker(12, 43200), unitchecker(2, 86400) /)
TYPE(unitchecker),PARAMETER :: munit(3) = (/ & ! 13 14 COSMO extensions
 unitchecker(0, 60), unitchecker(13, 900), unitchecker(14, 1800) /)

unit = imiss
IF (.NOT.c_e(valuein1) .OR. .NOT.c_e(valuein2)) THEN
  valueout1 = imiss
  valueout2 = imiss
  unit = 1
ELSE IF (mod(valuein1,3600) == 0 .AND. mod(valuein2,3600) == 0) THEN ! prefer hours
  DO i = 1, SIZE(hunit)
    IF (mod(valuein1, hunit(i)%sectounit) == 0 &
     .AND. mod(valuein2, hunit(i)%sectounit) == 0 &
     .AND. valuein1/hunit(i)%sectounit < 255 &
     .AND. valuein2/hunit(i)%sectounit < 255) THEN
      valueout1 = valuein1/hunit(i)%sectounit
      valueout2 = valuein2/hunit(i)%sectounit
      unit = hunit(i)%unit
      EXIT
    ENDIF
  ENDDO
  IF (.NOT.c_e(unit)) THEN
! last chance, disable overflow check and start from longest unit,
    DO i = SIZE(hunit), 1, -1
      IF (mod(valuein1, hunit(i)%sectounit) == 0 &
       .AND. mod(valuein2, hunit(i)%sectounit) == 0) THEN
        valueout1 = valuein1/hunit(i)%sectounit
        valueout2 = valuein2/hunit(i)%sectounit
        unit = hunit(i)%unit
        EXIT
      ENDIF
    ENDDO
  ENDIF
ELSE IF (mod(valuein1,60) == 0. .AND. mod(valuein2,60) == 0) THEN ! then minutes
  DO i = 1, SIZE(munit)
    IF (mod(valuein1, munit(i)%sectounit) == 0 &
     .AND. mod(valuein2, munit(i)%sectounit) == 0 &
     .AND. valuein1/munit(i)%sectounit < 255 &
     .AND. valuein2/munit(i)%sectounit < 255) THEN
      valueout1 = valuein1/munit(i)%sectounit
      valueout2 = valuein2/munit(i)%sectounit
      unit = munit(i)%unit
      EXIT
    ENDIF
  ENDDO
  IF (.NOT.c_e(unit)) THEN
! last chance, disable overflow check and start from longest unit,
    DO i = SIZE(munit), 1, -1
      IF (mod(valuein1, munit(i)%sectounit) == 0 &
       .AND. mod(valuein2, munit(i)%sectounit) == 0) THEN
        valueout1 = valuein1/munit(i)%sectounit
        valueout2 = valuein2/munit(i)%sectounit
        unit = munit(i)%unit
        EXIT
      ENDIF
    ENDDO
  ENDIF
ENDIF

IF (.NOT.c_e(unit)) THEN
  CALL l4f_log(l4f_error,'timerange_second_to_g1: cannot find a grib1 timerange unit for coding ' &
   //t2c(valuein1)//','//t2c(valuein2)//'s intervals' )
    CALL raise_error()
ENDIF

END SUBROUTINE timerange_choose_unit_g1


! Standardize variables and timerange in DB-all.e thinking:
!
! Timerange 205 (point in time interval) is converted to maximum or
! minimum if parameter is recognized as such and parameter is
! corrected as well, otherwise 205 is kept (with possible error
! conditions later).
!
! Timerange 257 (COSMO nudging) is converted to point in time if
! interval length is 0, or to a proper timerange if parameter is
! recognized as a COSMO statistically processed parameter (and in case
! of maximum or minimum the parameter is corrected as well); if
! parameter is not recognized, it is converted to instantaneous with a
! warning.
SUBROUTINE normalize_gridinfo(this)
TYPE(gridinfo_def),intent(inout) :: this

IF (this%timerange%timerange == 254) THEN ! instantaneous

!tmin
  IF (this%var == volgrid6d_var_new(255,2,16,255)) THEN
    this%var%number=11
    RETURN
  ENDIF

!tmax
  IF (this%var == volgrid6d_var_new(255,2,15,255)) THEN
    this%var%number=11
    RETURN
  ENDIF

ELSE IF (this%timerange%timerange == 205) THEN ! point in time interval

!tmin
  IF (this%var == volgrid6d_var_new(255,2,16,255)) THEN
    this%var%number=11
    this%timerange%timerange=3
    RETURN
  ENDIF

!tmax
  IF (this%var == volgrid6d_var_new(255,2,15,255)) THEN
    this%var%number=11
    this%timerange%timerange=2
    RETURN
  ENDIF

! it is accepted to keep 187 since it is wind gust, not max wind
  IF (this%var%discipline == 255 .AND. &
   any(this%var%centre == cosmo_centre)) THEN ! grib1 & COSMO

    IF (this%var%category == 201) THEN ! table 201

      IF (this%var%number == 187) THEN ! wind gust
!        this%var%category=2
!        this%var%number=32
        this%timerange%timerange=2
      ENDIF
    ENDIF
  ENDIF

ELSE IF (this%timerange%timerange == 257) THEN ! COSMO-nudging

  IF (this%timerange%p2 == 0) THEN ! point in time

    this%timerange%timerange=254

  ELSE ! guess timerange according to parameter

    IF (this%var%discipline == 255 .AND. &
     any(this%var%centre == cosmo_centre)) THEN ! grib1 & COSMO

      IF (this%var%category >= 1 .AND. this%var%category <= 3) THEN ! WMO table 2

        if (this%var%number == 11) then ! T
          this%timerange%timerange=0 ! average

        else if (this%var%number == 15) then ! T
          this%timerange%timerange=2 ! maximum
          this%var%number=11 ! reset also parameter

        else if (this%var%number == 16) then ! T
          this%timerange%timerange=3 ! minimum
          this%var%number=11 ! reset also parameter

        else if (this%var%number == 17) then ! TD
          this%timerange%timerange=0 ! average

        else if (this%var%number == 33) then ! U
          this%timerange%timerange=0 ! average

        else if (this%var%number == 34) then ! V
          this%timerange%timerange=0 ! average

        else if (this%var%number == 57) then ! evaporation
          this%timerange%timerange=1 ! accumulated

        else if (this%var%number == 61) then ! TOT_PREC
          this%timerange%timerange=1 ! accumulated

        else if (this%var%number == 78) then ! SNOW_CON
          this%timerange%timerange=1 ! accumulated

        else if (this%var%number == 79) then ! SNOW_GSP
          this%timerange%timerange=1 ! accumulated

        else if (this%var%number == 90) then ! RUNOFF
          this%timerange%timerange=1 ! accumulated

        else if (this%var%number == 111) then ! fluxes
          this%timerange%timerange=0 ! average
        else if (this%var%number == 112) then
          this%timerange%timerange=0 ! average
        else if (this%var%number == 113) then
          this%timerange%timerange=0 ! average
        else if (this%var%number == 114) then
          this%timerange%timerange=0 ! average
        else if (this%var%number == 121) then
          this%timerange%timerange=0 ! average
        else if (this%var%number == 122) then
          this%timerange%timerange=0 ! average
        else if (this%var%number == 124) then
          this%timerange%timerange=0 ! average
        else if (this%var%number == 125) then
          this%timerange%timerange=0 ! average
        else if (this%var%number == 126) then
          this%timerange%timerange=0 ! average
        else if (this%var%number == 127) then
          this%timerange%timerange=0 ! average

        endif

      ELSE IF (this%var%category == 201) THEN ! table 201

        if (this%var%number == 5) then ! photosynthetic flux
          this%timerange%timerange=0 ! average

        else if (this%var%number == 20) then ! SUN_DUR
          this%timerange%timerange=1 ! accumulated

        else if (this%var%number == 22) then ! radiation fluxes
          this%timerange%timerange=0 ! average
        else if (this%var%number == 23) then
          this%timerange%timerange=0 ! average
        else if (this%var%number == 24) then
          this%timerange%timerange=0 ! average
        else if (this%var%number == 25) then
          this%timerange%timerange=0 ! average
        else if (this%var%number == 26) then
          this%timerange%timerange=0 ! average
        else if (this%var%number == 27) then
          this%timerange%timerange=0 ! average

        else if (this%var%number == 42) then ! water divergence
          this%timerange%timerange=1 ! accumulated

        else if (this%var%number == 102) then ! RAIN_GSP
          this%timerange%timerange=1 ! accumulated

        else if (this%var%number == 113) then ! RAIN_CON
          this%timerange%timerange=1 ! accumulated

        else if (this%var%number == 132) then ! GRAU_GSP
          this%timerange%timerange=1 ! accumulated

        else if (this%var%number == 135) then ! HAIL_GSP
          this%timerange%timerange=1 ! accumulated

        else if (this%var%number == 187) then ! UVMAX
!          this%var%category=2 ! reset also parameter
!          this%var%number=32
          this%timerange%timerange=2 ! maximum

        else if (this%var%number == 218) then ! maximum 10m dynamical gust
          this%timerange%timerange=2 ! maximum

        else if (this%var%number == 219) then ! maximum 10m convective gust
          this%timerange%timerange=2 ! maximum

        endif

      ELSE IF (this%var%category == 202) THEN ! table 202

        if (this%var%number == 231) then ! sso parameters
          this%timerange%timerange=0 ! average
        else if (this%var%number == 232) then
          this%timerange%timerange=0 ! average
        else if (this%var%number == 233) then
          this%timerange%timerange=0 ! average
        endif

      ELSE ! parameter not recognized, set instantaneous?

        call l4f_category_log(this%category,l4f_warn, &
         'normalize_gridinfo: found COSMO non instantaneous analysis 13,0,'//&
         trim(to_char(this%timerange%p2)))
        call l4f_category_log(this%category,l4f_warn, &
         'associated to an apparently instantaneous parameter '//&
         trim(to_char(this%var%centre))//','//trim(to_char(this%var%category))//','//&
         trim(to_char(this%var%number))//','//trim(to_char(this%var%discipline)))
        call l4f_category_log(this%category,l4f_warn, 'forcing to instantaneous')

        this%timerange%p2 = 0
        this%timerange%timerange = 254

      ENDIF ! category
    ENDIF ! grib1 & COSMO
  ENDIF ! p2
ENDIF ! timerange 

IF (this%var%discipline == 255 .AND. &
 any(this%var%centre == ecmwf_centre)) THEN ! grib1 & ECMWF
! useful references:
! http://www.ecmwf.int/publications/manuals/libraries/tables/tables_index.html
! http://www.ecmwf.int/products/data/technical/soil/discret_soil_lay.html

  IF (this%var%category == 128) THEN ! table 128

    IF ((this%var%number == 142 .OR. & ! large scale precipitation
     this%var%number == 143 .OR. & ! convective precipitation
     this%var%number == 144 .OR. & ! total snow
     this%var%number == 228 .OR. & ! total precipitation
     this%var%number == 145 .OR. & ! boundary layer dissipation
     this%var%number == 146 .OR. & ! surface sensible heat flux
     this%var%number == 147 .OR. & ! surface latent heat flux
     this%var%number == 169) .AND. & ! surface solar radiation downwards
     this%timerange%timerange == 254) THEN
      this%timerange%timerange = 1 ! accumulated
      this%timerange%p2 = this%timerange%p1 ! length of period = forecast time

    ELSE IF ((this%var%number == 165 .OR. & ! 10m U
     this%var%number == 166) .AND. & ! 10m V
     this%level%level1 == 1) THEN
      this%level%level1 = 103
      this%level%l1 = 10000 ! 10m

    ELSE IF ((this%var%number == 167 .OR. & ! 2m T
     this%var%number == 168) .AND. & ! 2m Td
     this%level%level1 == 1) THEN
      this%level%level1 = 103
      this%level%l1 = 2000 ! 2m

    ELSE IF (this%var%number == 39 .OR. this%var%number == 139 .OR. this%var%number == 140) THEN ! SWVL1,STL1,SWL1
      this%level%level1 = 106 ! below surface
      this%level%l1 = 0
      this%level%l2 = 70 ! 7cm

    ELSE IF (this%var%number == 40 .OR. this%var%number == 170) THEN ! SWVL2,STL2 (STL2 wrong before 2000)
      this%level%level1 = 106 ! below surface
      this%level%l1 = 70
      this%level%l2 = 280

    ELSE IF (this%var%number == 171) THEN ! SWL2
      this%level%level1 = 106 ! below surface
      this%level%l1 = 70
      this%level%l2 = 210

    ELSE IF (this%var%number == 41 .OR. this%var%number == 183) THEN ! SWVL3,STL3 (STL3 wrong before 2000)
      this%level%level1 = 106 ! below surface
      this%level%l1 = 280
      this%level%l2 = 1000

    ELSE IF (this%var%number == 184) THEN ! SWL3
      this%level%level1 = 106 ! below surface
      this%level%l1 = 210
      this%level%l2 = 1000

    ELSE IF (this%var%number == 42 .OR. this%var%number == 236 .OR. this%var%number == 237) THEN ! SWVL4,STL4,SWL4
      this%level%level1 = 106 ! below surface
      this%level%l1 = 1000
      this%level%l2 = 2890

    ELSE IF (this%var%number == 121 .AND. &
     (this%timerange%timerange == 254 .OR. this%timerange%timerange == 205)) THEN ! MX2T6
      this%timerange%timerange = 2 ! max
      this%timerange%p2 = 21600 ! length of period = 6 hours
      this%var%number=167 ! set to T2m, it could be 130 T as well
      this%level%level1 = 103
      this%level%l1 = 2000 ! 2m

    ELSE IF (this%var%number == 122 .AND. &
     (this%timerange%timerange == 254 .OR. this%timerange%timerange == 205)) THEN ! MN2T6
      this%timerange%timerange = 3 ! min
      this%timerange%p2 = 21600 ! length of period = 6 hours
      this%var%number=1
      this%var%number=167 ! set to T2m, it could be 130 T as well
      this%level%level1 = 103
      this%level%l1 = 2000 ! 2m

    ELSE IF (this%var%number == 123 .AND. &
     (this%timerange%timerange == 254 .OR. this%timerange%timerange == 205)) THEN ! 10FG6
      this%timerange%timerange = 2 ! max
      this%timerange%p2 = 21600 ! length of period = 6 hours
      this%level%level1 = 103
      this%level%l1 = 10000 ! 10m

! set cloud cover to bufr style
    ELSE IF (this%var%number == 186) THEN ! low cloud cover
      this%var%number = 248
      this%level = vol7d_level_new(level1=256, level2=258, l2=1)
    ELSE IF (this%var%number == 187) THEN ! medium cloud cover
      this%var%number = 248
      this%level = vol7d_level_new(level1=256, level2=258, l2=2)
    ELSE IF (this%var%number == 188) THEN ! high cloud cover
      this%var%number = 248
      this%level = vol7d_level_new(level1=256, level2=258, l2=3)

    ENDIF
  ELSE IF (this%var%category == 228) THEN ! table 228

    IF (this%var%number == 24) THEN
      this%level%level1 = 4 ! Level of 0C Isotherm
      this%level%l1 = 0
      this%level%level2 = 255
      this%level%l2 = 0

    ELSE IF (this%var%number == 26 .AND. &
     (this%timerange%timerange == 254 .OR. this%timerange%timerange == 205)) THEN ! MX2T3
      this%timerange%timerange = 2 ! max
      this%timerange%p2 = 10800 ! length of period = 3 hours
      this%var%category = 128 ! reset to table version 128
      this%var%number=167 ! set to T2m, it could be 130 T as well
      this%level%level1 = 103
      this%level%l1 = 2000 ! 2m

    ELSE IF (this%var%number == 27 .AND. &
     (this%timerange%timerange == 254 .OR. this%timerange%timerange == 205)) THEN ! MN2T3
      this%timerange%timerange = 3 ! min
      this%timerange%p2 = 10800 ! length of period = 3 hours
      this%var%category = 128 ! reset to table version 128
      this%var%number=167 ! set to T2m, it could be 130 T as well
      this%level%level1 = 103
      this%level%l1 = 2000 ! 2m

    ELSE IF (this%var%number == 28 .AND. &
     (this%timerange%timerange == 254 .OR. this%timerange%timerange == 205)) THEN ! 10FG3
      this%timerange%timerange = 2 ! max
      this%timerange%p2 = 10800 ! length of period = 3 hours
      this%level%level1 = 103
      this%level%l1 = 10000 ! 10m

    ENDIF

  ENDIF ! table 128
ENDIF ! grib1 & ECMWF

IF (this%var%discipline == 255 .AND. &
 this%var%category >= 1 .AND. this%var%category <= 3) THEN ! grib1 WMO table 2

! set cloud cover to bufr style
  IF (this%var%number == 73) THEN ! low cloud cover
    this%var%number = 71
    this%level = vol7d_level_new(level1=256, level2=258, l2=1)
  ELSE IF (this%var%number == 74) THEN ! medium cloud cover
    this%var%number = 71
    this%level = vol7d_level_new(level1=256, level2=258, l2=2)
  ELSE IF (this%var%number == 75) THEN ! high cloud cover
    this%var%number = 71
    this%level = vol7d_level_new(level1=256, level2=258, l2=3)

  ENDIF

ENDIF


END SUBROUTINE normalize_gridinfo


! Destandardize variables and timerange from DB-all.e thinking:
!
! Parameters having maximum or minimum statistical processing and
! having a corresponding extreme parameter in grib table, are
! temporarily converted to timerange 205 and to the corresponding
! extreme parameter; if parameter is not recognized, the max or min
! statistical processing is kept (with possible error conditions
! later).
SUBROUTINE unnormalize_gridinfo(this)
TYPE(gridinfo_def),intent(inout) :: this

IF (this%timerange%timerange == 3) THEN ! min

  IF (this%var == volgrid6d_var_new(255,2,11,255)) THEN ! tmin
    this%var%number=16
    this%timerange%timerange=205

  ELSE IF (any(this%var%centre == ecmwf_centre)) THEN ! ECMWF
    IF (this%var == volgrid6d_var_new(this%var%centre,128,167,255)) THEN ! tmin
      this%var%number=122
      this%timerange%timerange=205

    ENDIF
  ENDIF
ELSE IF (this%timerange%timerange == 2) THEN ! max

  IF (this%var == volgrid6d_var_new(255,2,11,255)) THEN ! tmax
    this%var%number=15
    this%timerange%timerange=205

  ELSE IF (any(this%var%centre == ecmwf_centre)) THEN ! ECMWF
    IF (this%var == volgrid6d_var_new(this%var%centre,128,167,255)) THEN ! tmax
      this%var%number=121
      this%timerange%timerange=205

    ELSE IF(this%var == volgrid6d_var_new(this%var%centre,128,123,255)) THEN ! uvmax
      this%timerange%timerange=205

    ELSE IF(this%var == volgrid6d_var_new(this%var%centre,228,28,255)) THEN ! uvmax
      this%timerange%timerange=205

    ENDIF
  ELSE IF (any(this%var%centre == cosmo_centre)) THEN ! grib1 & COSMO

! wind
! it is accepted to keep 187 since it is wind gust, not max wind
!    IF (this%var == volgrid6d_var_new(255,2,32,255)) THEN
!      this%var%category=201
!      this%var%number=187
!      this%timerange%timerange=205
!    ENDIF
    IF (this%var == volgrid6d_var_new(this%var%centre,201,187,255)) THEN
      this%timerange%timerange=205
    ENDIF
  ENDIF
ENDIF

! reset cloud cover to grib1 style
IF (this%var%discipline == 255 .AND. this%var%category == 2) THEN ! grib1 table 2
  IF (this%var%number == 71 .AND. &
   this%level%level1 == 256 .AND. this%level%level2 == 258) THEN ! l/m/h cloud cover
    IF (this%level%l2 == 1) THEN ! l
      this%var%number = 73
    ELSE IF (this%level%l2 == 2) THEN ! m
      this%var%number = 74
    ELSE IF (this%level%l2 == 3) THEN ! h
      this%var%number = 75
    ENDIF
    this%level = vol7d_level_new(level1=1) ! reset to surface
  ENDIF
ENDIF

IF (any(this%var%centre == ecmwf_centre)) THEN ! ECMWF
! reset cloud cover to grib1 style
  IF (this%var%discipline == 255 .AND. this%var%category == 128) THEN ! grib1 table 128
    IF ((this%var%number == 248 .OR. this%var%number == 164) .AND. &
     this%level%level1 == 256 .AND. this%level%level2 == 258) THEN ! l/m/h cloud cover
      IF (this%level%l2 == 1) THEN ! l
        this%var%number = 186
      ELSE IF (this%level%l2 == 2) THEN ! m
        this%var%number = 187
      ELSE IF (this%level%l2 == 3) THEN ! h
        this%var%number = 188
      ENDIF
      this%level = vol7d_level_new(level1=1) ! reset to surface
    ENDIF
  ENDIF
ENDIF

END SUBROUTINE unnormalize_gridinfo
#endif


! =========================================
! gdal driver specific code
! could this be moved to a separate module?
! =========================================
#ifdef HAVE_LIBGDAL
SUBROUTINE gridinfo_import_gdal(this, gdalid)
TYPE(gridinfo_def),INTENT(inout) :: this ! gridinfo object
TYPE(gdalrasterbandh),INTENT(in) :: gdalid ! gdal rasterband pointer

TYPE(gdaldataseth) :: hds


!call time_import_gdal(this%time, gaid)
this%time = datetime_new(year=2010, month=1, day=1) ! conventional year

!call timerange_import_gdal(this%timerange,gaid)
this%timerange = vol7d_timerange_new(254, 0, 0) ! instantaneous data

!call level_import_gdal(this%level, gaid)
hds = gdalgetbanddataset(gdalid) ! go back to dataset
IF (gdalgetrastercount(hds) == 1) THEN ! single level dataset
  this%level = vol7d_level_new(1, 0) ! surface
ELSE
  this%level = vol7d_level_new(105, gdalgetbandnumber(gdalid)) ! hybrid level
ENDIF

!call var_import_gdal(this%var, gaid)
this%var = volgrid6d_var_new(centre=255, category=2, number=8) ! topography height

END SUBROUTINE gridinfo_import_gdal
#endif


END MODULE gridinfo_class