!> @file data_output_netcdf4_module.f90
!--------------------------------------------------------------------------------------------------!
! This file is part of the PALM model system.
!
! PALM 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 3 of the License, or (at your option) any later
! version.
!
! PALM 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
! PALM. If not, see <http://www.gnu.org/licenses/>.
!
! Copyright 2019-2019 Leibniz Universitaet Hannover
!--------------------------------------------------------------------------------------------------!
!
! Current revisions:
! ------------------
!
!
! Former revisions:
! -----------------
! $Id: data_output_netcdf4_module.f90 4106 2019-07-19 08:54:42Z gronemeier $
! Initial revision
!
!
! Authors:
! --------
!> @author: Tobias Gronemeier
!
! Description:
! ------------
!> NetCDF output module to write data to NetCDF files.
!> This is either done in parallel mode via parallel NetCDF4 I/O or in serial mode only by PE0.
!--------------------------------------------------------------------------------------------------!
MODULE data_output_netcdf4_module

   USE kinds

#if defined( __parallel )
#if defined( __mpifh )
   INCLUDE "mpif.h"
#else
   USE MPI
#endif
#endif

#if defined( __netcdf4 )
   USE NETCDF
#endif

   IMPLICIT NONE

   CHARACTER(LEN=800) ::  internal_error_message = ''  !< string containing the last error message
   CHARACTER(LEN=800) ::  temp_string                  !< dummy string

   CHARACTER(LEN=*), PARAMETER ::  mode_parallel = 'parallel'  !< string selecting netcdf4 parallel mode
   CHARACTER(LEN=*), PARAMETER ::  mode_serial   = 'serial'    !< string selecting netcdf4 serial mode

   INTEGER(iwp) ::  debug_output_unit       !< Fortran Unit Number of the debug-output file
   INTEGER(iwp) ::  global_id_in_file = -1  !< value of global ID within a file

   LOGICAL ::  print_debug_output = .FALSE.  !< if true, debug output is printed

   SAVE

   PRIVATE

   INTERFACE netcdf4_init_module
      MODULE PROCEDURE netcdf4_init_module
   END INTERFACE netcdf4_init_module

   INTERFACE netcdf4_open_file
      MODULE PROCEDURE netcdf4_open_file
   END INTERFACE netcdf4_open_file

   INTERFACE netcdf4_init_dimension
      MODULE PROCEDURE netcdf4_init_dimension
   END INTERFACE netcdf4_init_dimension

   INTERFACE netcdf4_init_variable
      MODULE PROCEDURE netcdf4_init_variable
   END INTERFACE netcdf4_init_variable

   INTERFACE netcdf4_write_attribute
      MODULE PROCEDURE netcdf4_write_attribute
   END INTERFACE netcdf4_write_attribute

   INTERFACE netcdf4_init_end
      MODULE PROCEDURE netcdf4_init_end
   END INTERFACE netcdf4_init_end

   INTERFACE netcdf4_write_variable
      MODULE PROCEDURE netcdf4_write_variable
   END INTERFACE netcdf4_write_variable

   INTERFACE netcdf4_finalize
      MODULE PROCEDURE netcdf4_finalize
   END INTERFACE netcdf4_finalize

   INTERFACE netcdf4_get_error_message
      MODULE PROCEDURE netcdf4_get_error_message
   END INTERFACE netcdf4_get_error_message

   PUBLIC &
      netcdf4_finalize, &
      netcdf4_get_error_message, &
      netcdf4_init_dimension, &
      netcdf4_init_end, &
      netcdf4_init_module, &
      netcdf4_init_variable, &
      netcdf4_open_file, &
      netcdf4_write_attribute, &
      netcdf4_write_variable


CONTAINS


!--------------------------------------------------------------------------------------------------!
! Description:
! ------------
!> Initialize data-output module.
!--------------------------------------------------------------------------------------------------!
SUBROUTINE netcdf4_init_module( program_debug_output_unit, debug_output, dom_global_id )

   INTEGER(iwp), INTENT(IN) ::  dom_global_id              !< global id within a file defined by DOM
   INTEGER(iwp), INTENT(IN) ::  program_debug_output_unit  !< file unit number for debug output

   LOGICAL, INTENT(IN) ::  debug_output  !< if true, debug output is printed


   debug_output_unit = program_debug_output_unit

   print_debug_output = debug_output

   global_id_in_file = dom_global_id

END SUBROUTINE netcdf4_init_module

!--------------------------------------------------------------------------------------------------!
! Description:
! ------------
!> Open netcdf file.
!--------------------------------------------------------------------------------------------------!
SUBROUTINE netcdf4_open_file( mode, filename, file_id, return_value )

   CHARACTER(LEN=*), INTENT(IN) ::  filename  !< name of file
   CHARACTER(LEN=*), INTENT(IN) ::  mode      !< operation mode (either parallel or serial)

   CHARACTER(LEN=*), PARAMETER ::  routine_name = 'netcdf4_open_file'  !< name of this routine

   INTEGER(iwp), INTENT(OUT) ::  file_id       !< file ID
   INTEGER(iwp)              ::  nc_stat       !< netcdf return value
   INTEGER(iwp), INTENT(OUT) ::  return_value  !< return value


   return_value = 0

   !-- Open new file
   CALL internal_message( 'debug', routine_name // ': create file "' // TRIM( filename ) // '"' )

   IF ( TRIM( mode ) == mode_serial )  THEN

#if defined( __netcdf4 )
      nc_stat = NF90_CREATE( TRIM( filename ),                    &
                             IOR( NF90_NOCLOBBER, NF90_NETCDF4 ), &
                             file_id )
#else
      file_id = -1
      nc_stat = 0
      return_value = 1
      CALL internal_message( 'error', routine_name //                               &
                             ': pre-processor directive "__netcdf4" not given. ' // &
                             'Using NetCDF4 output not possible' )
#endif

   ELSEIF ( TRIM( mode ) == mode_parallel )  THEN

#if defined( __parallel ) && defined( __netcdf4 ) && defined( __netcdf4_parallel )
      nc_stat = NF90_CREATE( TRIM( filename ),                                       &
                             IOR( NF90_NOCLOBBER, IOR( NF90_NETCDF4, NF90_MPIIO ) ), &
                             file_id, COMM = MPI_COMM_WORLD, INFO = MPI_INFO_NULL )
#else
      file_id = -1
      nc_stat = 0
      return_value = 1
      CALL internal_message( 'error', routine_name //                                 &
                             ': pre-processor directives "__parallel" and/or ' //     &
                             '"__netcdf4" and/or "__netcdf4_parallel" not given. ' // &
                             'Using parallel NetCDF4 output not possible' )
#endif

   ELSE
      file_id = -1
      nc_stat = 0
      return_value = 1
      CALL internal_message( 'error', routine_name // ': selected mode "' // &
                                      TRIM( mode ) // '" must be either "' // &
                                      mode_serial // '" or "' // mode_parallel // '"' )
   ENDIF

#if defined( __netcdf4 )
   IF ( nc_stat /= NF90_NOERR )  THEN
      return_value = 1
      CALL internal_message( 'error', routine_name // ': NetCDF error while opening file "' // &
                                      TRIM( filename ) // '": ' // NF90_STRERROR( nc_stat ) )
   ENDIF
#endif

END SUBROUTINE netcdf4_open_file

!--------------------------------------------------------------------------------------------------!
! Description:
! ------------
!> Write attribute to netcdf file.
!--------------------------------------------------------------------------------------------------!
SUBROUTINE netcdf4_write_attribute( file_id, var_id, att_name, att_value_char, &
                 att_value_int8, att_value_int16, att_value_int32,             &
                 att_value_real32, att_value_real64, return_value )

   CHARACTER(LEN=*), INTENT(IN)           ::  att_name        !< name of attribute
   CHARACTER(LEN=*), INTENT(IN), OPTIONAL ::  att_value_char  !< value of attribute

   CHARACTER(LEN=*), PARAMETER ::  routine_name = 'netcdf4_write_attribute'  !< name of this routine

   INTEGER(iwp) ::  nc_stat    !< netcdf return value
   INTEGER(iwp) ::  target_id  !< ID of target which gets attribute (either global or var_id)

   INTEGER(iwp), INTENT(IN)  ::  file_id       !< file ID
   INTEGER(iwp), INTENT(OUT) ::  return_value  !< return value
   INTEGER(iwp), INTENT(IN)  ::  var_id        !< variable ID

   INTEGER(KIND=1), INTENT(IN), OPTIONAL ::  att_value_int8   !< value of attribute
   INTEGER(KIND=2), INTENT(IN), OPTIONAL ::  att_value_int16  !< value of attribute
   INTEGER(KIND=4), INTENT(IN), OPTIONAL ::  att_value_int32  !< value of attribute

   REAL(KIND=4), INTENT(IN), OPTIONAL ::  att_value_real32  !< value of attribute
   REAL(KIND=8), INTENT(IN), OPTIONAL ::  att_value_real64  !< value of attribute


#if defined( __netcdf4 )
   return_value = 0

   IF ( var_id == global_id_in_file )  THEN
      target_id = NF90_GLOBAL
   ELSE
      target_id = var_id
   ENDIF

   CALL internal_message( 'debug', &
                          routine_name // ': write attribute "' // TRIM( att_name ) // '"' )

   IF ( PRESENT( att_value_char ) )  THEN
      nc_stat = NF90_PUT_ATT( file_id, target_id, TRIM( att_name ), TRIM( att_value_char ) )
   ELSEIF ( PRESENT( att_value_int8 ) )  THEN
      nc_stat = NF90_PUT_ATT( file_id, target_id, TRIM( att_name ), att_value_int8 )
   ELSEIF ( PRESENT( att_value_int16 ) )  THEN
      nc_stat = NF90_PUT_ATT( file_id, target_id, TRIM( att_name ), att_value_int16 )
   ELSEIF ( PRESENT( att_value_int32 ) )  THEN
      nc_stat = NF90_PUT_ATT( file_id, target_id, TRIM( att_name ), att_value_int32 )
   ELSEIF ( PRESENT( att_value_real32 ) )  THEN
      nc_stat = NF90_PUT_ATT( file_id, target_id, TRIM( att_name ), att_value_real32 )
   ELSEIF ( PRESENT( att_value_real64 ) )  THEN
      nc_stat = NF90_PUT_ATT( file_id, target_id, TRIM( att_name ), att_value_real64 )
   ELSE
      return_value = 1
      CALL internal_message( 'error', TRIM( routine_name ) // &
                                      ': attribute "' // TRIM( att_name ) // '": no value given' )
   ENDIF

   IF ( return_value == 0 )  THEN
      IF ( nc_stat /= NF90_NOERR )  THEN
         return_value = 1
         CALL internal_message( 'error',                                       &
                 routine_name // ': NetCDF error while writing attribute "' // &
                 TRIM( att_name ) // '": ' // NF90_STRERROR( nc_stat ) )
      ENDIF
   ENDIF
#else
   return_value = 1
#endif

END SUBROUTINE netcdf4_write_attribute

!--------------------------------------------------------------------------------------------------!
! Description:
! ------------
!> Initialize dimension.
!--------------------------------------------------------------------------------------------------!
SUBROUTINE netcdf4_init_dimension( mode, file_id, dim_id, var_id, &
              dim_name, dim_type, dim_length, return_value )

   CHARACTER(LEN=*), INTENT(IN) ::  dim_name  !< name of dimension
   CHARACTER(LEN=*), INTENT(IN) ::  dim_type  !< data type of dimension
   CHARACTER(LEN=*), INTENT(IN) ::  mode      !< operation mode (either parallel or serial)

   CHARACTER(LEN=*), PARAMETER ::  routine_name = 'netcdf4_init_dimension'  !< name of this routine

   INTEGER(iwp), INTENT(OUT) ::  dim_id         !< dimension ID
   INTEGER(iwp), INTENT(IN)  ::  dim_length     !< length of dimension
   INTEGER(iwp), INTENT(IN)  ::  file_id        !< file ID
   INTEGER(iwp)              ::  nc_dim_length  !< length of dimension
   INTEGER(iwp)              ::  nc_stat        !< netcdf return value
   INTEGER(iwp), INTENT(OUT) ::  return_value   !< return value
   INTEGER(iwp), INTENT(OUT) ::  var_id         !< variable ID


#if defined( __netcdf4 )
   return_value = 0
   var_id = -1

   CALL internal_message( 'debug', &
                          routine_name // ': init dimension "' // TRIM( dim_name ) // '"' )

   !-- Check if dimension is unlimited
   IF ( dim_length < 0 )  THEN
      nc_dim_length = NF90_UNLIMITED
   ELSE
      nc_dim_length = dim_length
   ENDIF

   !-- Define dimension in file
   nc_stat = NF90_DEF_DIM( file_id, dim_name, nc_dim_length, dim_id )

   IF ( nc_stat == NF90_NOERR )  THEN

      !-- Define variable holding dimension values in file
      CALL netcdf4_init_variable( mode, file_id, var_id, dim_name, dim_type, (/dim_id/), &
                                           is_global=.TRUE., return_value=return_value )

   ELSE
      return_value = 1
      CALL internal_message( 'error', routine_name //                                    &
                                      ': NetCDF error while initializing dimension "' // &
                                      TRIM( dim_name ) // '": ' // NF90_STRERROR( nc_stat ) )
   ENDIF
#else
   return_value = 1
   var_id = -1
   dim_id = -1
#endif

END SUBROUTINE netcdf4_init_dimension

!--------------------------------------------------------------------------------------------------!
! Description:
! ------------
!> Initialize variable.
!--------------------------------------------------------------------------------------------------!
SUBROUTINE netcdf4_init_variable( mode, file_id, var_id, var_name, var_type, var_dim_ids, &
                                  is_global, return_value )

   CHARACTER(LEN=*), INTENT(IN) ::  mode      !< operation mode (either parallel or serial)
   CHARACTER(LEN=*), INTENT(IN) ::  var_name  !< name of variable
   CHARACTER(LEN=*), INTENT(IN) ::  var_type  !< data type of variable

   CHARACTER(LEN=*), PARAMETER ::  routine_name = 'netcdf4_init_variable'  !< name of this routine

   INTEGER(iwp), INTENT(IN)  ::  file_id       !< file ID
   INTEGER(iwp)              ::  nc_stat       !< netcdf return value
   INTEGER(iwp)              ::  nc_var_type   !< netcdf data type
   INTEGER(iwp), INTENT(OUT) ::  return_value  !< return value
   INTEGER(iwp), INTENT(OUT) ::  var_id        !< variable ID

   INTEGER(iwp), DIMENSION(:), INTENT(IN) ::  var_dim_ids  !< list of dimension IDs used by variable

   LOGICAL, INTENT(IN)  ::  is_global  !< true if variable is global (same on all PE)


#if defined( __netcdf4 )
   return_value = 0

   WRITE( temp_string, * ) is_global
   CALL internal_message( 'debug', routine_name // ': init variable "' // TRIM( var_name ) // &
                                   '" ( is_global = ' // TRIM( temp_string ) // ')' )

   nc_var_type = get_netcdf_data_type( var_type )

   IF ( nc_var_type /= -1_iwp )  THEN

      !-- Define variable in file
      nc_stat = NF90_DEF_VAR( file_id, var_name, nc_var_type, var_dim_ids, var_id )

#if defined( __netcdf4_parallel )
      !-- Define how variable can be accessed by PEs in parallel netcdf file
      IF ( nc_stat == NF90_NOERR  .AND.  TRIM( mode ) == mode_parallel )  THEN
         IF ( is_global )  THEN
            nc_stat = NF90_VAR_PAR_ACCESS( file_id, var_id, NF90_INDEPENDENT )
         ELSE
            nc_stat = NF90_VAR_PAR_ACCESS( file_id, var_id, NF90_COLLECTIVE )
         ENDIF
      ENDIF
#endif

      IF ( nc_stat /= NF90_NOERR)  THEN
         return_value = 1
         CALL internal_message( 'error', routine_name //                                   &
                                         ': NetCDF error while initializing variable "' // &
                                         TRIM( var_name ) // '": ' // NF90_STRERROR( nc_stat ) )
      ENDIF

   ELSE
      return_value = 1
   ENDIF

#else
   return_value = 1
   var_id = -1
#endif

END SUBROUTINE netcdf4_init_variable

!--------------------------------------------------------------------------------------------------!
! Description:
! ------------
!> Leave file definition state.
!--------------------------------------------------------------------------------------------------!
SUBROUTINE netcdf4_init_end( file_id, return_value )

   CHARACTER(LEN=*), PARAMETER ::  routine_name = 'netcdf4_init_end'  !< name of this routine

   INTEGER(iwp), INTENT(IN)  ::  file_id       !< file ID
   INTEGER(iwp)              ::  nc_stat       !< netcdf return value
   INTEGER(iwp)              ::  old_mode      !< previous netcdf fill mode
   INTEGER(iwp), INTENT(OUT) ::  return_value  !< return value


#if defined( __netcdf4 )
   return_value = 0

   WRITE( temp_string, * ) file_id
   CALL internal_message( 'debug',        &
                          routine_name // &
                          ': finalize file definition (file_id=' // TRIM( temp_string ) // ')' )

   !-- Set general no fill, otherwise the performance drops significantly
   nc_stat = NF90_SET_FILL( file_id, NF90_NOFILL, old_mode )

   IF ( nc_stat == NF90_NOERR )  THEN
      nc_stat = NF90_ENDDEF( file_id )
   ENDIF

   IF ( nc_stat /= NF90_NOERR )  THEN
      return_value = 1
      CALL internal_message( 'error', routine_name // ': NetCDF error: ' // &
                                      NF90_STRERROR( nc_stat ) )
   ENDIF
#else
   return_value = 1
#endif

END SUBROUTINE netcdf4_init_end

!--------------------------------------------------------------------------------------------------!
! Description:
! ------------
!> Write variable of different kind into netcdf file.
!--------------------------------------------------------------------------------------------------!
SUBROUTINE netcdf4_write_variable(                                        &
              file_id, var_id, bounds_start, bounds_end, bounds_origin,   &
              do_output, is_global,                                       &
              var_int8_0d,   var_int8_1d,   var_int8_2d,   var_int8_3d,   &
              var_int16_0d,  var_int16_1d,  var_int16_2d,  var_int16_3d,  &
              var_int32_0d,  var_int32_1d,  var_int32_2d,  var_int32_3d,  &
              var_intwp_0d,  var_intwp_1d,  var_intwp_2d,  var_intwp_3d,  &
              var_real32_0d, var_real32_1d, var_real32_2d, var_real32_3d, &
              var_real64_0d, var_real64_1d, var_real64_2d, var_real64_3d, &
              var_realwp_0d, var_realwp_1d, var_realwp_2d, var_realwp_3d, &
              return_value )

   CHARACTER(LEN=*), PARAMETER ::  routine_name = 'netcdf4_write_variable'  !< name of this routine

   INTEGER(iwp)              ::  d             !< loop index
   INTEGER(iwp), INTENT(IN)  ::  file_id       !< file ID
   INTEGER(iwp)              ::  myid          !< id number of processor element
   INTEGER(iwp)              ::  nc_stat       !< netcdf return value
   INTEGER(iwp)              ::  ndim          !< number of dimensions of variable in file
   INTEGER(iwp), INTENT(OUT) ::  return_value  !< return value
   INTEGER(iwp), INTENT(IN)  ::  var_id        !< variable ID

   INTEGER(iwp), DIMENSION(:), INTENT(IN)  ::  bounds_origin  !< starting index of each dimension
   INTEGER(iwp), DIMENSION(:), INTENT(IN)  ::  bounds_end     !< ending index of variable
   INTEGER(iwp), DIMENSION(:), INTENT(IN)  ::  bounds_start   !< starting index of variable
   INTEGER(iwp), DIMENSION(:), ALLOCATABLE ::  dim_ids        !< IDs of dimensions of variable in file
   INTEGER(iwp), DIMENSION(:), ALLOCATABLE ::  dim_lengths    !< length of dimensions of variable in file
   INTEGER(iwp), DIMENSION(:), ALLOCATABLE ::  value_count    !< count of values along each dimension to be written

   INTEGER(KIND=1), POINTER, INTENT(IN), OPTIONAL                               ::  var_int8_0d  !< output variable
   INTEGER(KIND=1), POINTER, CONTIGUOUS, INTENT(IN), OPTIONAL, DIMENSION(:)     ::  var_int8_1d  !< output variable
   INTEGER(KIND=1), POINTER, CONTIGUOUS, INTENT(IN), OPTIONAL, DIMENSION(:,:)   ::  var_int8_2d  !< output variable
   INTEGER(KIND=1), POINTER, CONTIGUOUS, INTENT(IN), OPTIONAL, DIMENSION(:,:,:) ::  var_int8_3d  !< output variable

   INTEGER(KIND=2), POINTER, INTENT(IN), OPTIONAL                               ::  var_int16_0d  !< output variable
   INTEGER(KIND=2), POINTER, CONTIGUOUS, INTENT(IN), OPTIONAL, DIMENSION(:)     ::  var_int16_1d  !< output variable
   INTEGER(KIND=2), POINTER, CONTIGUOUS, INTENT(IN), OPTIONAL, DIMENSION(:,:)   ::  var_int16_2d  !< output variable
   INTEGER(KIND=2), POINTER, CONTIGUOUS, INTENT(IN), OPTIONAL, DIMENSION(:,:,:) ::  var_int16_3d  !< output variable

   INTEGER(KIND=4), POINTER, INTENT(IN), OPTIONAL                               ::  var_int32_0d  !< output variable
   INTEGER(KIND=4), POINTER, CONTIGUOUS, INTENT(IN), OPTIONAL, DIMENSION(:)     ::  var_int32_1d  !< output variable
   INTEGER(KIND=4), POINTER, CONTIGUOUS, INTENT(IN), OPTIONAL, DIMENSION(:,:)   ::  var_int32_2d  !< output variable
   INTEGER(KIND=4), POINTER, CONTIGUOUS, INTENT(IN), OPTIONAL, DIMENSION(:,:,:) ::  var_int32_3d  !< output variable

   INTEGER(iwp), POINTER, INTENT(IN), OPTIONAL                               ::  var_intwp_0d  !< output variable
   INTEGER(iwp), POINTER, CONTIGUOUS, INTENT(IN), OPTIONAL, DIMENSION(:)     ::  var_intwp_1d  !< output variable
   INTEGER(iwp), POINTER, CONTIGUOUS, INTENT(IN), OPTIONAL, DIMENSION(:,:)   ::  var_intwp_2d  !< output variable
   INTEGER(iwp), POINTER, CONTIGUOUS, INTENT(IN), OPTIONAL, DIMENSION(:,:,:) ::  var_intwp_3d  !< output variable

   LOGICAL, INTENT(IN) ::  do_output  !< if false, set count to 0 and do no output
   LOGICAL, INTENT(IN) ::  is_global  !< true if variable is global (same on all PE)

   REAL(KIND=4), POINTER, INTENT(IN), OPTIONAL                               ::  var_real32_0d  !< output variable
   REAL(KIND=4), POINTER, CONTIGUOUS, INTENT(IN), OPTIONAL, DIMENSION(:)     ::  var_real32_1d  !< output variable
   REAL(KIND=4), POINTER, CONTIGUOUS, INTENT(IN), OPTIONAL, DIMENSION(:,:)   ::  var_real32_2d  !< output variable
   REAL(KIND=4), POINTER, CONTIGUOUS, INTENT(IN), OPTIONAL, DIMENSION(:,:,:) ::  var_real32_3d  !< output variable

   REAL(KIND=8), POINTER, INTENT(IN), OPTIONAL                               ::  var_real64_0d  !< output variable
   REAL(KIND=8), POINTER, CONTIGUOUS, INTENT(IN), OPTIONAL, DIMENSION(:)     ::  var_real64_1d  !< output variable
   REAL(KIND=8), POINTER, CONTIGUOUS, INTENT(IN), OPTIONAL, DIMENSION(:,:)   ::  var_real64_2d  !< output variable
   REAL(KIND=8), POINTER, CONTIGUOUS, INTENT(IN), OPTIONAL, DIMENSION(:,:,:) ::  var_real64_3d  !< output variable

   REAL(wp), POINTER, INTENT(IN), OPTIONAL                               ::  var_realwp_0d  !< output variable
   REAL(wp), POINTER, CONTIGUOUS, INTENT(IN), OPTIONAL, DIMENSION(:)     ::  var_realwp_1d  !< output variable
   REAL(wp), POINTER, CONTIGUOUS, INTENT(IN), OPTIONAL, DIMENSION(:,:)   ::  var_realwp_2d  !< output variable
   REAL(wp), POINTER, CONTIGUOUS, INTENT(IN), OPTIONAL, DIMENSION(:,:,:) ::  var_realwp_3d  !< output variable


#if defined( __netcdf4 )

#if defined( __parallel )
   CALL MPI_COMM_RANK( MPI_COMM_WORLD, myid, return_value )
   IF ( return_value /= 0 )  THEN
      CALL internal_message( 'error', routine_name // ': MPI_COMM_RANK error' )
   ENDIF
#else
   myid = 0
   return_value = 0
#endif

   IF ( return_value == 0  .AND.  ( .NOT. is_global  .OR.  myid == 0 ) )  THEN

      WRITE( temp_string, * ) var_id
      CALL internal_message( 'debug', routine_name // ': write variable ' // TRIM( temp_string ) )

      ndim = SIZE( bounds_start )

      ALLOCATE( value_count(ndim) )

      IF ( do_output ) THEN
         value_count = bounds_end - bounds_start + 1
      ELSE
         value_count = 0
      END IF

      !-- 8bit integer output
      IF ( PRESENT( var_int8_0d ) )  THEN
         nc_stat = NF90_PUT_VAR( file_id, var_id, (/ var_int8_0d /),       &
                                 start = bounds_start - bounds_origin + 1, &
                                 count = value_count )
      ELSEIF ( PRESENT( var_int8_1d ) )  THEN
         nc_stat = NF90_PUT_VAR( file_id, var_id, var_int8_1d,             &
                                 start = bounds_start - bounds_origin + 1, &
                                 count = value_count )
      ELSEIF ( PRESENT( var_int8_2d ) )  THEN
         nc_stat = NF90_PUT_VAR( file_id, var_id, var_int8_2d,             &
                                 start = bounds_start - bounds_origin + 1, &
                                 count = value_count )
      ELSEIF ( PRESENT( var_int8_3d ) )  THEN
         nc_stat = NF90_PUT_VAR( file_id, var_id, var_int8_3d,             &
                                 start = bounds_start - bounds_origin + 1, &
                                 count = value_count )
      !-- 16bit integer output
      ELSEIF ( PRESENT( var_int16_0d ) )  THEN
         nc_stat = NF90_PUT_VAR( file_id, var_id, (/ var_int16_0d /),      &
                                 start = bounds_start - bounds_origin + 1, &
                                 count = value_count )
      ELSEIF ( PRESENT( var_int16_1d ) )  THEN
         nc_stat = NF90_PUT_VAR( file_id, var_id, var_int16_1d,            &
                                 start = bounds_start - bounds_origin + 1, &
                                 count = value_count )
      ELSEIF ( PRESENT( var_int16_2d ) )  THEN
         nc_stat = NF90_PUT_VAR( file_id, var_id, var_int16_2d,            &
                                 start = bounds_start - bounds_origin + 1, &
                                 count = value_count )
      ELSEIF ( PRESENT( var_int16_3d ) )  THEN
         nc_stat = NF90_PUT_VAR( file_id, var_id, var_int16_3d,            &
                                 start = bounds_start - bounds_origin + 1, &
                                 count = value_count )
      !-- 32bit integer output
      ELSEIF ( PRESENT( var_int32_0d ) )  THEN
         nc_stat = NF90_PUT_VAR( file_id, var_id, (/ var_int32_0d /),      &
                                 start = bounds_start - bounds_origin + 1, &
                                 count = value_count )
      ELSEIF ( PRESENT( var_int32_1d ) )  THEN
         nc_stat = NF90_PUT_VAR( file_id, var_id, var_int32_1d,            &
                                 start = bounds_start - bounds_origin + 1, &
                                 count = value_count )
      ELSEIF ( PRESENT( var_int32_2d ) )  THEN
         nc_stat = NF90_PUT_VAR( file_id, var_id, var_int32_2d,            &
                                 start = bounds_start - bounds_origin + 1, &
                                 count = value_count )
      ELSEIF ( PRESENT( var_int32_3d ) )  THEN
         nc_stat = NF90_PUT_VAR( file_id, var_id, var_int32_3d,            &
                                 start = bounds_start - bounds_origin + 1, &
                                 count = value_count )
      !-- working-precision integer output
      ELSEIF ( PRESENT( var_intwp_0d ) )  THEN
         nc_stat = NF90_PUT_VAR( file_id, var_id, (/ var_intwp_0d /),      &
                                 start = bounds_start - bounds_origin + 1, &
                                 count = value_count )
      ELSEIF ( PRESENT( var_intwp_1d ) )  THEN
         nc_stat = NF90_PUT_VAR( file_id, var_id, var_intwp_1d,            &
                                 start = bounds_start - bounds_origin + 1, &
                                 count = value_count )
      ELSEIF ( PRESENT( var_intwp_2d ) )  THEN
         nc_stat = NF90_PUT_VAR( file_id, var_id, var_intwp_2d,            &
                                 start = bounds_start - bounds_origin + 1, &
                                 count = value_count )
      ELSEIF ( PRESENT( var_intwp_3d ) )  THEN
         nc_stat = NF90_PUT_VAR( file_id, var_id, var_intwp_3d,            &
                                 start = bounds_start - bounds_origin + 1, &
                                 count = value_count )
      !-- 32bit real output
      ELSEIF ( PRESENT( var_real32_0d ) )  THEN
         nc_stat = NF90_PUT_VAR( file_id, var_id, (/ var_real32_0d /),     &
                                 start = bounds_start - bounds_origin + 1, &
                                 count = value_count )
      ELSEIF ( PRESENT( var_real32_1d ) )  THEN
         nc_stat = NF90_PUT_VAR( file_id, var_id, var_real32_1d,           &
                                 start = bounds_start - bounds_origin + 1, &
                                 count = value_count )
      ELSEIF ( PRESENT( var_real32_2d ) )  THEN
         nc_stat = NF90_PUT_VAR( file_id, var_id, var_real32_2d,           &
                                 start = bounds_start - bounds_origin + 1, &
                                 count = value_count )
      ELSEIF ( PRESENT( var_real32_3d ) )  THEN
         nc_stat = NF90_PUT_VAR( file_id, var_id, var_real32_3d,           &
                                 start = bounds_start - bounds_origin + 1, &
                                 count = value_count )
      !-- 64bit real output
      ELSEIF ( PRESENT( var_real64_0d ) )  THEN
         nc_stat = NF90_PUT_VAR( file_id, var_id, (/ var_real64_0d /),     &
                                 start = bounds_start - bounds_origin + 1, &
                                 count = value_count )
      ELSEIF ( PRESENT( var_real64_1d ) )  THEN
         nc_stat = NF90_PUT_VAR( file_id, var_id, var_real64_1d,           &
                                 start = bounds_start - bounds_origin + 1, &
                                 count = value_count )
      ELSEIF ( PRESENT( var_real64_2d ) )  THEN
         nc_stat = NF90_PUT_VAR( file_id, var_id, var_real64_2d,           &
                                 start = bounds_start - bounds_origin + 1, &
                                 count = value_count )
      ELSEIF ( PRESENT( var_real64_3d ) )  THEN
         nc_stat = NF90_PUT_VAR( file_id, var_id, var_real64_3d,           &
                                 start = bounds_start - bounds_origin + 1, &
                                 count = value_count )
      !-- working-precision real output
      ELSEIF ( PRESENT( var_realwp_0d ) )  THEN
         nc_stat = NF90_PUT_VAR( file_id, var_id, (/ var_realwp_0d /),     &
                                 start = bounds_start - bounds_origin + 1, &
                                 count = value_count )
      ELSEIF ( PRESENT( var_realwp_1d ) )  THEN
         nc_stat = NF90_PUT_VAR( file_id, var_id, var_realwp_1d,           &
                                 start = bounds_start - bounds_origin + 1, &
                                 count = value_count )
      ELSEIF ( PRESENT( var_realwp_2d ) )  THEN
         nc_stat = NF90_PUT_VAR( file_id, var_id, var_realwp_2d,           &
                                 start = bounds_start - bounds_origin + 1, &
                                 count = value_count )
      ELSEIF ( PRESENT( var_realwp_3d ) )  THEN
         nc_stat = NF90_PUT_VAR( file_id, var_id, var_realwp_3d,           &
                                 start = bounds_start - bounds_origin + 1, &
                                 count = value_count )
      ELSE
         return_value = 1
         nc_stat = NF90_NOERR
         WRITE( temp_string, '(": variable_id=",I6,", file_id=",I6)' ) var_id, file_id
         CALL internal_message( 'error', routine_name // &
                                TRIM( temp_string ) //   &
                                ': no output values given' )
      ENDIF

      !-- Check for errors
      IF ( nc_stat /= NF90_NOERR )  THEN
         return_value = 1

         IF ( nc_stat == NF90_EEDGE .OR. nc_stat == NF90_EINVALCOORDS )  THEN

            !-- If given bounds exceed dimension bounds, get information of bounds in file
            WRITE( temp_string, * )  NF90_STRERROR( nc_stat )

            ALLOCATE( dim_ids(ndim) )
            ALLOCATE( dim_lengths(ndim) )

            nc_stat = NF90_INQUIRE_VARIABLE( file_id, var_id, dimids=dim_ids )

            d = 1
            DO WHILE ( d <= ndim .AND. nc_stat == NF90_NOERR )
               nc_stat = NF90_INQUIRE_DIMENSION( file_id, dim_ids(d), len=dim_lengths(d) )
               d = d + 1
            ENDDO

            IF ( nc_stat == NF90_NOERR )  THEN
               WRITE( temp_string, * )  TRIM( temp_string ) // '; given variable bounds: ' //  &
                  'start=', bounds_start, ', end=', bounds_end, '; file dimension bounds: ' // &
                  'start=', bounds_origin, ', end=', bounds_origin + dim_lengths - 1
               CALL internal_message( 'error', routine_name //     &
                                      ': error while writing: ' // &
                                      TRIM( temp_string ) )
            ELSE
               !-- Error occured during NF90_INQUIRE_VARIABLE or NF90_INQUIRE_DIMENSION
               CALL internal_message( 'error', routine_name //            &
                                      ': error while accessing file: ' // &
                                      NF90_STRERROR( nc_stat ) )
            ENDIF

         ELSE
            !-- Other NetCDF error
            CALL internal_message( 'error', routine_name //     &
                                   ': error while writing: ' // &
                                   NF90_STRERROR( nc_stat ) )
         ENDIF
      ENDIF

   ENDIF
#else
   return_value = 1
#endif

END SUBROUTINE netcdf4_write_variable

!--------------------------------------------------------------------------------------------------!
! Description:
! ------------
!> Close netcdf file.
!--------------------------------------------------------------------------------------------------!
SUBROUTINE netcdf4_finalize( file_id, return_value )

   CHARACTER(LEN=*), PARAMETER ::  routine_name = 'netcdf4_finalize'  !< name of routine

   INTEGER(iwp), INTENT(IN)  ::  file_id       !< file ID
   INTEGER(iwp)              ::  nc_stat       !< netcdf return value
   INTEGER(iwp), INTENT(OUT) ::  return_value  !< return value


#if defined( __netcdf4 )
   WRITE( temp_string, * ) file_id
   CALL internal_message( 'debug', routine_name // &
                                   ': close file (file_id=' // TRIM( temp_string ) // ')' )

   nc_stat = NF90_CLOSE( file_id )
   IF ( nc_stat == NF90_NOERR )  THEN
      return_value = 0
   ELSE
      return_value = 1
      CALL internal_message( 'error', routine_name // &
                                      ': NetCDF error: ' // NF90_STRERROR( nc_stat ) )
   ENDIF
#else
   return_value = 1
#endif

END SUBROUTINE netcdf4_finalize

!--------------------------------------------------------------------------------------------------!
! Description:
! ------------
!> Convert data_type string into netcdf data type value.
!--------------------------------------------------------------------------------------------------!
FUNCTION get_netcdf_data_type( data_type ) RESULT( return_value )

   CHARACTER(LEN=*), INTENT(IN) ::  data_type  !< requested data type

   CHARACTER(LEN=*), PARAMETER ::  routine_name = 'get_netcdf_data_type'  !< name of this routine

   INTEGER(iwp) ::  return_value  !< netcdf data type


   SELECT CASE ( TRIM( data_type ) )

#if defined( __netcdf4 )
      CASE ( 'char' )
         return_value = NF90_CHAR

      CASE ( 'int8' )
         return_value = NF90_BYTE

      CASE ( 'int16' )
         return_value = NF90_SHORT

      CASE ( 'int32' )
         return_value = NF90_INT

      CASE ( 'real32' )
         return_value = NF90_FLOAT

      CASE ( 'real64' )
         return_value = NF90_DOUBLE
#endif

      CASE DEFAULT
         CALL internal_message( 'error', routine_name // &
                                         ': data type unknown (' // TRIM( data_type ) // ')' )
         return_value = -1_iwp

   END SELECT

END FUNCTION get_netcdf_data_type

!--------------------------------------------------------------------------------------------------!
! Description:
! ------------
!> Message routine writing debug information into the debug file
!> or creating the error message string.
!--------------------------------------------------------------------------------------------------!
SUBROUTINE internal_message( level, string )

   CHARACTER(LEN=*), INTENT(IN) ::  level   !< message importance level
   CHARACTER(LEN=*), INTENT(IN) ::  string  !< message string


   IF ( TRIM( level ) == 'error' )  THEN

      WRITE( internal_error_message, '(A,A)' ) ': ', string

   ELSEIF ( TRIM( level ) == 'debug'  .AND.  print_debug_output )  THEN

      WRITE( debug_output_unit, '(A,A)' ) 'DOM DEBUG: ', string
      FLUSH( debug_output_unit )

   ENDIF

END SUBROUTINE internal_message

!--------------------------------------------------------------------------------------------------!
! Description:
! ------------
!> Return the last created error message.
!--------------------------------------------------------------------------------------------------!
SUBROUTINE netcdf4_get_error_message( error_message )

   CHARACTER(LEN=800), INTENT(OUT) ::  error_message  !< return error message to main program


   error_message = internal_error_message

END SUBROUTINE netcdf4_get_error_message


END MODULE data_output_netcdf4_module