diff --git a/.gitattributes b/.gitattributes
index 0ffda24590f2f62020407b18bf4f396775ea1452..397e9f9562b546b54f460533b021161bf0cfaa43 100644
--- a/.gitattributes
+++ b/.gitattributes
@@ -298,3 +298,6 @@ tests/Makefile.am -text
 tests/Makefile.in -text
 tests/test_grib.c -text
 tests/test_grib.sh -text
+util/sunf95preproc-wrapper -text
+util/sxpreproc-wrapper -text
+util/xlfpreproc-wrapper -text
diff --git a/ChangeLog b/ChangeLog
index 2f5d6dcae5bc3332fa1495f95d0d12ba44fe9818..3272bd13ba896f92970bd7269122eb795c77e6fe 100644
--- a/ChangeLog
+++ b/ChangeLog
@@ -1,3 +1,58 @@
+2011-12-13 Uwe Schulzweida  <Uwe.Schulzweida@zmaw.de>
+
+	* verify_coordinate_vars: bug fix in check for units = "1" [report: Katharina Six]
+
+2011-11-11 Uwe Schulzweida  <Uwe.Schulzweida@zmaw.de>
+
+	* added support for netCDF attributes scale_factor and add_offset for lon/lat coordinates
+
+2011-11-04 Uwe Schulzweida  <Uwe.Schulzweida@zmaw.de>
+
+	* added support for GRIB1_LTYPE_SIGMA_LAYER
+
+2011-11-01 Uwe Schulzweida  <Uwe.Schulzweida@zmaw.de>
+
+	* added support for netcdf attribute valid_min/valid_max [request: Etienne Tourigny]
+
+2011-10-27 Uwe Schulzweida  <Uwe.Schulzweida@zmaw.de>
+
+	* cdilib.c::defineAttributes: bug fix atttxt [report: Florian Prill]
+
+2011-10-25 Uwe Schulzweida  <Uwe.Schulzweida@zmaw.de>
+
+	* added support for netcdf attribute valid_range [request: Etienne Tourigny]
+
+2011-10-17  Uwe Schulzweida  <Uwe.Schulzweida@zmaw.de>
+
+	* Version 1.5.3 released
+	* using CGRIBEX library version 1.5.1
+
+2011-10-11  Uwe Schulzweida  <Uwe.Schulzweida@zmaw.de>
+
+       * zaxisCompare: set epsilon from 0 to 1e-9 [request: Felicia Brisc]
+
+2011-10-06  Uwe Schulzweida  <Uwe.Schulzweida@zmaw.de>
+
+	* added level type ZAXIS_TOA, ZAXIS_SEA_BOTTOM, ZAXIS_ATMOSPHERE [request: Dörte Liermann]
+
+2011-10-05  Uwe Schulzweida  <Uwe.Schulzweida@zmaw.de>
+
+	* stream_cdf::cdfInqContents: check units of hybrid levels
+	* varAddRecord: used max number of bit_per_value for 3D GRIB data
+	* gribapiDefGrid: added parameter jScansPositively [report: Juan Jose Tasso]
+
+2011-10-02  Uwe Schulzweida  <Uwe.Schulzweida@zmaw.de>
+
+	* deflate compression with netCDF4 doesn't work (bug fix) [report: Geert Jan van Oldenborgh]
+
+2011-09-21  Uwe Schulzweida  <Uwe.Schulzweida@zmaw.de>
+
+	* correct netCDF dimension order of unstructured grids (bug fix) [report: Ralf Mueller]
+
+2011-08-22 Uwe Schulzweida  <Uwe.Schulzweida@zmaw.de>
+
+	* Version 1.5.2 released
+
 2011-08-15  Uwe Schulzweida  <Uwe.Schulzweida@zmaw.de>
 
 	* streamFilesuffix: added suffix for filetype NC4C (bug fix)
diff --git a/NEWS b/NEWS
index 6689445bdc5fd2ddc75dce086f38b584e62358e3..dd35b622b60ece1a31ad7e93fcfe3c04cb0258bd 100644
--- a/NEWS
+++ b/NEWS
@@ -1,6 +1,13 @@
 CDI NEWS
 --------
 
+Version 1.5.3 (17 October 2011):
+
+   New features:
+     * Added support for level type ZAXIS_TOA, ZAXIS_SEA_BOTTOM, ZAXIS_ATMOSPHERE
+   Fixed bugs:
+     * deflate compression with netCDF4 doesn't work
+
 Version 1.5.2 (22 August 2011):
 
    New features:
diff --git a/config/default b/config/default
index b098139f85056c7c12cdf5452915a41006e2ce59..af14564dc36ead7edee1bad7e9b030f020131b1a 100755
--- a/config/default
+++ b/config/default
@@ -87,6 +87,14 @@ case "${HOSTNAME}" in
 		    CC=sxc++ AR=sxar RANLIB=ls \
                     CFLAGS="-O -Onooverlap,restrict=all -pvctl,fullmsg,noassume,loopcnt=1000000"
 	;;
+    lxe0*)
+        echo 'Please choose compiler modules! Checkout with "module av"!'
+        ./configure --prefix=$(pwd)/build-SX  --host=sx9-nec-superux \
+              --with-netcdf=/usr/local/pkg-sx9 \
+              --with-grib_api=/usr/local/pkg-sx9/grib_api   CC=sxc++ FC=sxf90 \
+              LD=/SX/opt/crosskit/inst/bin/sxld AR=/SX/opt/crosskit/inst/bin/sxar \
+              RANLIB=echo
+	;;
 # mips-sgi-irix6.5
     ecgate1)
 	./configure --prefix=$HOME/local \
diff --git a/configure b/configure
index ce0990b041ce06d28dcf55512a1a65988fca5108..ee25571f6e3ac098106c2d45b64f3f673943f1cf 100755
--- a/configure
+++ b/configure
@@ -1,6 +1,6 @@
 #! /bin/sh
 # Guess values for system-dependent variables and create Makefiles.
-# Generated by GNU Autoconf 2.68 for cdi 1.5.2.
+# Generated by GNU Autoconf 2.68 for cdi 1.5.4.
 #
 # Report bugs to <http://code.zmaw.de/projects/cdi>.
 #
@@ -570,8 +570,8 @@ MAKEFLAGS=
 # Identity of this package.
 PACKAGE_NAME='cdi'
 PACKAGE_TARNAME='cdi'
-PACKAGE_VERSION='1.5.2'
-PACKAGE_STRING='cdi 1.5.2'
+PACKAGE_VERSION='1.5.4'
+PACKAGE_STRING='cdi 1.5.4'
 PACKAGE_BUGREPORT='http://code.zmaw.de/projects/cdi'
 PACKAGE_URL=''
 
@@ -1394,7 +1394,7 @@ if test "$ac_init_help" = "long"; then
   # Omit some internal or obsolete options to make the list less imposing.
   # This message is too long to be a string in the A/UX 3.1 sh.
   cat <<_ACEOF
-\`configure' configures cdi 1.5.2 to adapt to many kinds of systems.
+\`configure' configures cdi 1.5.4 to adapt to many kinds of systems.
 
 Usage: $0 [OPTION]... [VAR=VALUE]...
 
@@ -1464,7 +1464,7 @@ fi
 
 if test -n "$ac_init_help"; then
   case $ac_init_help in
-     short | recursive ) echo "Configuration of cdi 1.5.2:";;
+     short | recursive ) echo "Configuration of cdi 1.5.4:";;
    esac
   cat <<\_ACEOF
 
@@ -1616,7 +1616,7 @@ fi
 test -n "$ac_init_help" && exit $ac_status
 if $ac_init_version; then
   cat <<\_ACEOF
-cdi configure 1.5.2
+cdi configure 1.5.4
 generated by GNU Autoconf 2.68
 
 Copyright (C) 2010 Free Software Foundation, Inc.
@@ -2331,7 +2331,7 @@ cat >config.log <<_ACEOF
 This file contains any messages produced by compilers while
 running configure, to aid debugging if configure makes a mistake.
 
-It was created by cdi $as_me 1.5.2, which was
+It was created by cdi $as_me 1.5.4, which was
 generated by GNU Autoconf 2.68.  Invocation command line was
 
   $ $0 $@
@@ -3263,7 +3263,7 @@ fi
 
 # Define the identity of the package.
  PACKAGE='cdi'
- VERSION='1.5.2'
+ VERSION='1.5.4'
 
 
 cat >>confdefs.h <<_ACEOF
@@ -24725,7 +24725,7 @@ fi
 
 done
 
-              { $as_echo "$as_me:${as_lineno-$LINENO}: checking for library containing deflate" >&5
+                      { $as_echo "$as_me:${as_lineno-$LINENO}: checking for library containing deflate" >&5
 $as_echo_n "checking for library containing deflate... " >&6; }
 if ${ac_cv_search_deflate+:} false; then :
   $as_echo_n "(cached) " >&6
@@ -25561,10 +25561,19 @@ $as_echo "no" >&6; }
 fi
 
 
-                                 $as_echo $NC_CONFIG
                                  if test "x$NC_CONFIG" != "x"; then :
-  { $as_echo "$as_me:${as_lineno-$LINENO}: result: found" >&5
-$as_echo "found" >&6; }
+  { $as_echo "$as_me:${as_lineno-$LINENO}: checking netcdf's OpenDAP support" >&5
+$as_echo_n "checking netcdf's OpenDAP support... " >&6; }
+                                    if test "x$($NC_CONFIG --has-dap)" = "xyes"; then :
+
+$as_echo "#define HAVE_LIBNC_DAP 1" >>confdefs.h
+
+                                           { $as_echo "$as_me:${as_lineno-$LINENO}: result: yes" >&5
+$as_echo "yes" >&6; }
+else
+  { $as_echo "$as_me:${as_lineno-$LINENO}: result: no" >&5
+$as_echo "no" >&6; }
+fi
                                    { $as_echo "$as_me:${as_lineno-$LINENO}: checking netcdf's nc2 support" >&5
 $as_echo_n "checking netcdf's nc2 support... " >&6; }
                                    if test "x$($NC_CONFIG --has-nc2)" = "xyes"; then :
@@ -26074,7 +26083,7 @@ ENABLE_IEG=$enable_ieg
 
 #  ----------------------------------------------------------------------
 # At the moment, there are two possible CDI bindings
-# (default for CDO) linking directly to CDI's object files in ./libcdi/src
+# (default for CDO) linking directly to CDI convenience library with libtool
 # (default for CDI) build and link to a shared CDI library
 if test "x$CDO_DISABLE_CDILIB" = "x1"; then :
   enable_cdi_lib=no
@@ -27180,7 +27189,7 @@ cat >>$CONFIG_STATUS <<\_ACEOF || ac_write_fail=1
 # report actual input values of CONFIG_FILES etc. instead of their
 # values after options handling.
 ac_log="
-This file was extended by cdi $as_me 1.5.2, which was
+This file was extended by cdi $as_me 1.5.4, which was
 generated by GNU Autoconf 2.68.  Invocation command line was
 
   CONFIG_FILES    = $CONFIG_FILES
@@ -27246,7 +27255,7 @@ _ACEOF
 cat >>$CONFIG_STATUS <<_ACEOF || ac_write_fail=1
 ac_cs_config="`$as_echo "$ac_configure_args" | sed 's/^ //; s/[\\""\`\$]/\\\\&/g'`"
 ac_cs_version="\\
-cdi config.status 1.5.2
+cdi config.status 1.5.4
 configured by $0, generated by GNU Autoconf 2.68,
   with options \\"\$ac_cs_config\\"
 
diff --git a/configure.ac b/configure.ac
index 70a09e23872f5fb5801db90a58227344d87f0023..b69e8efe1f3256af3cd82587e368bc53f6d9ed9c 100644
--- a/configure.ac
+++ b/configure.ac
@@ -1,6 +1,6 @@
 #  Process this file with autoconf to produce a configure script.
 
-AC_INIT([cdi], [1.5.2], [http://code.zmaw.de/projects/cdi])
+AC_INIT([cdi], [1.5.4], [http://code.zmaw.de/projects/cdi])
 
 echo "configuring ${PACKAGE_NAME} ${PACKAGE_VERSION}"
 
diff --git a/m4/acx_options.m4 b/m4/acx_options.m4
index 64c21b5dc98c4bfc59ca297ed61fd67787f02705..0305ea3edcf6072d641f9d6e8e5af1eecac85549 100644
--- a/m4/acx_options.m4
+++ b/m4/acx_options.m4
@@ -45,8 +45,8 @@ AC_ARG_WITH([zlib],
                           AC_SEARCH_LIBS([deflate],[z],[AC_DEFINE([HAVE_LIBZ],[1],[Define 1 for ZLIB support])])
                           ZLIB_INCLUDE=" -I$ZLIB_ROOT/include"
                           ZLIB_LIBS=" -L$ZLIB_ROOT/lib -lz"])],
-             [AC_CHECK_HEADERS(zlib.h)
-              AC_SEARCH_LIBS([deflate],[z],[AC_DEFINE([HAVE_LIBZ],[1],[Define 1 for ZLIB support])])
+                     [AC_CHECK_HEADERS(zlib.h)
+                      AC_SEARCH_LIBS([deflate],[z],[AC_DEFINE([HAVE_LIBZ],[1],[Define 1 for ZLIB support])])
               ZLIB_LIBS=" -lz"])
 AC_SUBST([ZLIB_INCLUDE])
 AC_SUBST([ZLIB_LIBS])
@@ -173,10 +173,12 @@ AC_ARG_WITH([netcdf],
                                  NETCDF_INCLUDE=" -I$NETCDF_ROOT/include"
                                  AC_MSG_CHECKING([nc-config script])
                                  AC_CHECK_PROG(NC_CONFIG,nc-config,[$NETCDF_ROOT/bin/nc-config],,["$NETCDF_ROOT/bin"])
-                                 AS_ECHO([$NC_CONFIG])
                                  AS_IF([test "x$NC_CONFIG" != "x"],
-                                   [AC_MSG_RESULT([found])
-                                   AC_MSG_CHECKING([netcdf's nc2 support])
+                                   [AC_MSG_CHECKING([netcdf's OpenDAP support])
+                                    AS_IF([test "x$($NC_CONFIG --has-dap)" = "xyes"],
+                                          [AC_DEFINE([HAVE_LIBNC_DAP],[1],[Define to 1 for NETCDF OpenDAP])
+                                           AC_MSG_RESULT([yes])],[AC_MSG_RESULT([no])])]
+                                   [AC_MSG_CHECKING([netcdf's nc2 support])
                                    AS_IF([test "x$($NC_CONFIG --has-nc2)" = "xyes"],
                                          [AC_DEFINE([HAVE_NETCDF2],[1],[Define to 1 for NETCDF2 support])
                                           AC_MSG_RESULT([yes])],[AC_MSG_RESULT([no])])
@@ -310,7 +312,7 @@ AC_MSG_RESULT([$enable_ieg])
 AC_SUBST([ENABLE_IEG],[$enable_ieg])
 #  ----------------------------------------------------------------------
 # At the moment, there are two possible CDI bindings
-# (default for CDO) linking directly to CDI's object files in ./libcdi/src
+# (default for CDO) linking directly to CDI convenience library with libtool
 # (default for CDI) build and link to a shared CDI library
 AS_IF([test "x$CDO_DISABLE_CDILIB" = "x1"],[enable_cdi_lib=no],[enable_cdi_lib=yes])
 # save CDI binding mode for later automake use
diff --git a/src/Makefile.am b/src/Makefile.am
index c25c38664e337f5680b534b6bbfb28d01ea2886e..f31b096491ceab05e9e29226535f0ab36ff98d9e 100644
--- a/src/Makefile.am
+++ b/src/Makefile.am
@@ -124,7 +124,11 @@ cdilib.c:
 cdilib.o: cdilib.c
 	$(COMPILE) -c $<
 
-LOCALTARGETS  = pkgconfig/cdi.pc
+LOCALTARGETS  = 
+if ENABLE_CDI_LIB
+LOCALTARGETS += pkgconfig/cdi.pc
+endif
+
 if CREATE_ISOC
 LOCALTARGETS += mo_cdi.o mo_cdi.$(FCMODEXT)
 endif
@@ -149,6 +153,9 @@ if CREATE_ISOC
 CLEANFILES += $(top_builddir)/src/mo_cdi.$(FCMODEXT) $(top_builddir)/src/mo_cdi.o
 endif
 
+if ENABLE_CDI_LIB
+CLEANFILES += pkgconfig/cdi.pc
+
 install-exec-local: pkgconfig/cdi.pc
 	$(mkinstalldirs) "$(DESTDIR)$(libdir)/pkgconfig"
 	$(install_sh_DATA) pkgconfig/cdi.pc "$(DESTDIR)$(libdir)/pkgconfig/cdi.pc"
@@ -157,3 +164,10 @@ uninstall-local:
 	rm -f "$(DESTDIR)$(libdir)/pkgconfig/cdi.pc"
 	rmdir "$(DESTDIR)$(libdir)/pkgconfig"
 
+endif
+
+install-exec-hook:
+	-@rmdir "$(DESTDIR)$(libdir)"
+install-data-hook:
+	-@rmdir "$(DESTDIR)$(includedir)"
+
diff --git a/src/Makefile.in b/src/Makefile.in
index b92009c9edfc0b165da93312918a44b6d43e6e2f..434d2a2528cf408a568133b761a894a018c235e9 100644
--- a/src/Makefile.in
+++ b/src/Makefile.in
@@ -35,8 +35,10 @@ PRE_UNINSTALL = :
 POST_UNINSTALL = :
 build_triplet = @build@
 host_triplet = @host@
-@CREATE_ISOC_TRUE@am__append_1 = mo_cdi.o mo_cdi.$(FCMODEXT)
-@CREATE_ISOC_TRUE@am__append_2 = $(top_builddir)/src/mo_cdi.$(FCMODEXT) $(top_builddir)/src/mo_cdi.o
+@ENABLE_CDI_LIB_TRUE@am__append_1 = pkgconfig/cdi.pc
+@CREATE_ISOC_TRUE@am__append_2 = mo_cdi.o mo_cdi.$(FCMODEXT)
+@CREATE_ISOC_TRUE@am__append_3 = $(top_builddir)/src/mo_cdi.$(FCMODEXT) $(top_builddir)/src/mo_cdi.o
+@ENABLE_CDI_LIB_TRUE@am__append_4 = pkgconfig/cdi.pc
 subdir = src
 DIST_COMMON = $(am__include_HEADERS_DIST) $(srcdir)/Makefile.am \
 	$(srcdir)/Makefile.in $(srcdir)/config.h.in
@@ -401,9 +403,9 @@ libcdi_la_SOURCES = \
         stream.c         \
         swap.c           
 
-LOCALTARGETS = pkgconfig/cdi.pc $(am__append_1)
+LOCALTARGETS = $(am__append_1) $(am__append_2)
 #
-CLEANFILES = `ls *~` cdilib.c $(am__append_2)
+CLEANFILES = `ls *~` cdilib.c $(am__append_3) $(am__append_4)
 all: config.h
 	$(MAKE) $(AM_MAKEFLAGS) all-am
 
@@ -726,6 +728,8 @@ distclean-generic:
 maintainer-clean-generic:
 	@echo "This command is intended for maintainers to use"
 	@echo "it deletes files that may require special tools to rebuild."
+@ENABLE_CDI_LIB_FALSE@uninstall-local:
+@ENABLE_CDI_LIB_FALSE@install-exec-local:
 clean: clean-am
 
 clean-am: clean-generic clean-libLTLIBRARIES clean-libtool \
@@ -750,13 +754,15 @@ info: info-am
 info-am:
 
 install-data-am: install-includeHEADERS
-
+	@$(NORMAL_INSTALL)
+	$(MAKE) $(AM_MAKEFLAGS) install-data-hook
 install-dvi: install-dvi-am
 
 install-dvi-am:
 
 install-exec-am: install-exec-local install-libLTLIBRARIES
-
+	@$(NORMAL_INSTALL)
+	$(MAKE) $(AM_MAKEFLAGS) install-exec-hook
 install-html: install-html-am
 
 install-html-am:
@@ -798,15 +804,16 @@ ps-am:
 uninstall-am: uninstall-includeHEADERS uninstall-libLTLIBRARIES \
 	uninstall-local
 
-.MAKE: all install-am install-strip
+.MAKE: all install-am install-data-am install-exec-am install-strip
 
 .PHONY: CTAGS GTAGS all all-am all-local check check-am clean \
 	clean-generic clean-libLTLIBRARIES clean-libtool \
 	clean-noinstLTLIBRARIES ctags distclean distclean-compile \
 	distclean-generic distclean-hdr distclean-libtool \
 	distclean-tags distdir dvi dvi-am html html-am info info-am \
-	install install-am install-data install-data-am install-dvi \
-	install-dvi-am install-exec install-exec-am install-exec-local \
+	install install-am install-data install-data-am \
+	install-data-hook install-dvi install-dvi-am install-exec \
+	install-exec-am install-exec-hook install-exec-local \
 	install-html install-html-am install-includeHEADERS \
 	install-info install-info-am install-libLTLIBRARIES \
 	install-man install-pdf install-pdf-am install-ps \
@@ -841,13 +848,18 @@ pkgconfig/cdi.pc: pkgconfig/cdi.pc.in ../config.status
 #
 all-local: $(LOCALTARGETS) 
 
-install-exec-local: pkgconfig/cdi.pc
-	$(mkinstalldirs) "$(DESTDIR)$(libdir)/pkgconfig"
-	$(install_sh_DATA) pkgconfig/cdi.pc "$(DESTDIR)$(libdir)/pkgconfig/cdi.pc"
+@ENABLE_CDI_LIB_TRUE@install-exec-local: pkgconfig/cdi.pc
+@ENABLE_CDI_LIB_TRUE@	$(mkinstalldirs) "$(DESTDIR)$(libdir)/pkgconfig"
+@ENABLE_CDI_LIB_TRUE@	$(install_sh_DATA) pkgconfig/cdi.pc "$(DESTDIR)$(libdir)/pkgconfig/cdi.pc"
+
+@ENABLE_CDI_LIB_TRUE@uninstall-local:
+@ENABLE_CDI_LIB_TRUE@	rm -f "$(DESTDIR)$(libdir)/pkgconfig/cdi.pc"
+@ENABLE_CDI_LIB_TRUE@	rmdir "$(DESTDIR)$(libdir)/pkgconfig"
 
-uninstall-local:
-	rm -f "$(DESTDIR)$(libdir)/pkgconfig/cdi.pc"
-	rmdir "$(DESTDIR)$(libdir)/pkgconfig"
+install-exec-hook:
+	-@rmdir "$(DESTDIR)$(libdir)"
+install-data-hook:
+	-@rmdir "$(DESTDIR)$(includedir)"
 
 # Tell versions [3.59,3.63) of GNU make to not export all variables.
 # Otherwise a system limit (for SysV at least) may be exceeded.
diff --git a/src/cdi.h b/src/cdi.h
index fcf80d283f036139b57036ac970361b255410d15..2efd2a4c3516f6f021efc23092474d31b8825fbf 100644
--- a/src/cdi.h
+++ b/src/cdi.h
@@ -147,6 +147,10 @@ extern "C" {
 #define  ZAXIS_ALTITUDE          10  /* Altitude above mean sea level in meters */
 #define  ZAXIS_SIGMA             11  /* Sigma level                             */
 #define  ZAXIS_MEANSEA           12  /* Mean sea level                          */
+#define  ZAXIS_TOA               13  /* Norminal top of atmosphere              */
+#define  ZAXIS_SEA_BOTTOM        14  /* Sea bottom                              */
+#define  ZAXIS_ATMOSPHERE        15  /* Entire atmosphere                       */
+#define  ZAXIS_REFERENCE         16  /* zaxis reference number                  */
 
 /* TAXIS types */
 
diff --git a/src/cdi.inc b/src/cdi.inc
index a75e20bb6c7aab06dfa037515983d17f93d8320c..45d968538b92c94596546e52e3d761a31fdb9497 100644
--- a/src/cdi.inc
+++ b/src/cdi.inc
@@ -3,7 +3,7 @@
 !
 ! Author:
 ! -------
-! Uwe Schulzweida, MPI-MET, Hamburg,   October 2011
+! Uwe Schulzweida, MPI-MET, Hamburg,   December 2011
 !
 
       INTEGER    CDI_MAX_NAME          
@@ -244,6 +244,14 @@
       PARAMETER (ZAXIS_SIGMA            = 11)
       INTEGER    ZAXIS_MEANSEA         
       PARAMETER (ZAXIS_MEANSEA          = 12)
+      INTEGER    ZAXIS_TOA             
+      PARAMETER (ZAXIS_TOA              = 13)
+      INTEGER    ZAXIS_SEA_BOTTOM      
+      PARAMETER (ZAXIS_SEA_BOTTOM       = 14)
+      INTEGER    ZAXIS_ATMOSPHERE      
+      PARAMETER (ZAXIS_ATMOSPHERE       = 15)
+      INTEGER    ZAXIS_REFERENCE       
+      PARAMETER (ZAXIS_REFERENCE        = 16)
 !
 !  TAXIS types
 !
@@ -327,28 +335,22 @@
       PARAMETER (PIO_ASYNCH             =  3)
       INTEGER    PIO_FPGUARD           
       PARAMETER (PIO_FPGUARD            =  4)
-      INTEGER         pioInit
-!                                    (INTEGER         ,
-!                                     INTEGER         ,
-!                                     INTEGER         ,
-!                                     INTEGER         ,
-!                                     INTEGER         )
-      EXTERNAL        pioInit
-
-!                     pioFinalize
-      EXTERNAL        pioFinalize
-
 !                     pioEndDef
       EXTERNAL        pioEndDef
 
 !                     pioEndTimestepping
       EXTERNAL        pioEndTimestepping
 
-!                     pioWriteTimestep
+!                     pioFinalize
+      EXTERNAL        pioFinalize
+
+      INTEGER         pioInit
 !                                    (INTEGER         ,
 !                                     INTEGER         ,
+!                                     INTEGER         ,
+!                                     INTEGER         ,
 !                                     INTEGER         )
-      EXTERNAL        pioWriteTimestep
+      EXTERNAL        pioInit
 
       INTEGER         pioInqVarDecoChunk
 !                                    (INTEGER         ,
@@ -364,6 +366,12 @@
 !                                    (INTEGER         )
       EXTERNAL        pioNamespaceSetActive
 
+!                     pioWriteTimestep
+!                                    (INTEGER         ,
+!                                     INTEGER         ,
+!                                     INTEGER         )
+      EXTERNAL        pioWriteTimestep
+
       CHARACTER*80    cdiStringError
 !                                    (INTEGER         cdiErrno)
       EXTERNAL        cdiStringError
diff --git a/src/cdiFortran.c b/src/cdiFortran.c
index 9c895e25b5d04285ae23824370c9021aab6f16e7..f08ea15036bb719db70d6578a77e6aef76e9506f 100644
--- a/src/cdiFortran.c
+++ b/src/cdiFortran.c
@@ -53,14 +53,14 @@
 
 /*  CALENDAR types  */
 
-FCALLSCFUN5 (INT, pioInit, PIOINIT, pioinit, INT, INT, INT, INT, PINT)
-FCALLSCSUB0 (pioFinalize, PIOFINALIZE, piofinalize)
 FCALLSCSUB0 (pioEndDef, PIOENDDEF, pioenddef)
 FCALLSCSUB0 (pioEndTimestepping, PIOENDTIMESTEPPING, pioendtimestepping)
-FCALLSCSUB3 (pioWriteTimestep, PIOWRITETIMESTEP, piowritetimestep, INT, INT, INT)
+FCALLSCSUB0 (pioFinalize, PIOFINALIZE, piofinalize)
+FCALLSCFUN5 (INT, pioInit, PIOINIT, pioinit, INT, INT, INT, INT, PINT)
 FCALLSCFUN2 (INT, pioInqVarDecoChunk, PIOINQVARDECOCHUNK, pioinqvardecochunk, INT, INT)
 FCALLSCFUN2 (INT, pioInqVarDecoOff, PIOINQVARDECOOFF, pioinqvardecooff, INT, INT)
 FCALLSCSUB1 (pioNamespaceSetActive, PIONAMESPACESETACTIVE, pionamespacesetactive, INT)
+FCALLSCSUB3 (pioWriteTimestep, PIOWRITETIMESTEP, piowritetimestep, INT, INT, INT)
 FCALLSCFUN1 (STRING, cdiStringError, CDISTRINGERROR, cdistringerror, INT)
 FCALLSCSUB1 (cdiDebug, CDIDEBUG, cdidebug, INT)
 FCALLSCFUN0 (STRING, cdiLibraryVersion, CDILIBRARYVERSION, cdilibraryversion)
diff --git a/src/cgribex.h b/src/cgribex.h
index 86d241ef495ba79603192814c3c5e0d231ab21d7..373cccb792f9b1ad5a7b10e4a37d987866d77218 100644
--- a/src/cgribex.h
+++ b/src/cgribex.h
@@ -8,12 +8,16 @@
 
 /* GRIB1 Level Types */
 #define  GRIB1_LTYPE_SURFACE               1
+#define  GRIB1_LTYPE_TOA           8
+#define  GRIB1_LTYPE_SEA_BOTTOM            9
+#define  GRIB1_LTYPE_ATMOSPHERE           10
 #define  GRIB1_LTYPE_99                   99
 #define  GRIB1_LTYPE_ISOBARIC            100
 #define  GRIB1_LTYPE_MEANSEA             102
 #define  GRIB1_LTYPE_ALTITUDE            103
 #define  GRIB1_LTYPE_HEIGHT              105
 #define  GRIB1_LTYPE_SIGMA               107
+#define  GRIB1_LTYPE_SIGMA_LAYER         108
 #define  GRIB1_LTYPE_HYBRID              109
 #define  GRIB1_LTYPE_HYBRID_LAYER        110
 #define  GRIB1_LTYPE_LANDDEPTH           111
@@ -240,12 +244,3 @@ double calculate_pfactor(const double* spectralField, long fieldTruncation, long
 
 #endif  /* _CGRIBEX_H */ 
 
-/*
- * Local Variables:
- * c-file-style: "Java"
- * c-basic-offset: 2
- * indent-tabs-mode: nil
- * show-trailing-whitespace: t
- * require-trailing-newline: t
- * End:
- */
diff --git a/src/cgribexlib.c b/src/cgribexlib.c
index baefaa7703569954a61194d473c3e7936615a4c9..9183c6240a0c9b3b98d1cf74940a9aedd98a51a7 100644
--- a/src/cgribexlib.c
+++ b/src/cgribexlib.c
@@ -1,7 +1,7 @@
 
-/* Automatically generated by m214003 at 2011-02-15, do not edit */
+/* Automatically generated by m214003 at 2011-08-29, do not edit */
 
-/* CGRIBEXLIB_VERSION="1.5.0" */
+/* CGRIBEXLIB_VERSION="1.5.1" */
 
 #if defined (HAVE_CONFIG_H)
 #  include "config.h"
@@ -3268,9 +3268,7 @@ void encodeBMS(GRIBPACK *lGrib, long *gribLen, double *fsec3, int *isec4, double
   long bmsLen, bmsUnusedBits;
   long fsec4size;
   long z = *gribLen;
-#ifdef VECTORCODE
   unsigned int *imask;
-#endif
   static int lmissvalinfo = 1;
   /*  unsigned int c, imask; */
 
@@ -3820,7 +3818,7 @@ void gribEncode(int *isec0, int *isec1, int *isec2, double *fsec2, int *isec3,
   bmsIncluded = ISEC1_Sec2Or3Flag & 64;
 
   /* set max header len */
-  len = 3000;
+  len = 16384;
 
   /* add data len */
   numBytes = (ISEC4_NumBits+7)>>3;
@@ -3830,7 +3828,6 @@ void gribEncode(int *isec0, int *isec1, int *isec2, double *fsec2, int *isec3,
   /* add bitmap len */
   if ( bmsIncluded ) len += (klenp+7)>>3;
 
-
 #if defined (VECTORCODE)
   lGrib = (GRIBPACK *) malloc(len*sizeof(GRIBPACK));
   if ( lGrib == NULL ) SysError("No Memory!");
@@ -3886,17 +3883,15 @@ void gribEncode(int *isec0, int *isec1, int *isec2, double *fsec2, int *isec3,
   bdsend = gribLen;
   encodeES(lGrib, &gribLen, bdsstart);
 
-  if ( (size_t) gribLen > len )
-    Error("lGrib buffer too small! len = %d  gribLen = %d", len, gribLen);
-
   if ( (size_t) gribLen > kleng*sizeof(int) )
     Error("kgrib buffer too small! kleng = %d  gribLen = %d", kleng, gribLen);
 
 #if defined (VECTORCODE)
+  if ( (size_t) gribLen > len )
+    Error("lGrib buffer too small! len = %d  gribLen = %d", len, gribLen);
+
   (void) PACK_GRIB(lGrib, (unsigned char *)CGrib, gribLen, -1L);
-#endif
 
-#if defined (VECTORCODE)
   free(lGrib);
 #endif
 
@@ -4421,13 +4416,20 @@ void decode_double_array_common(unsigned char *igrib, long jlend, int NumBits,
   unsigned char *bits = igrib;
   unsigned int jmask;
   long i;
-  int tbits = 0;
+  unsigned int tbits = 0;
   int n_bits = NumBits;
   int t_bits = 0;
       
   jmask = (1 << n_bits) - 1;
   for ( i = 0; i < jlend; i++ )
     {
+      if (n_bits - t_bits > 8)
+	{
+	  tbits = (tbits << 16) | (bits[0] << 8) | (bits[1]);
+	  bits += 2;
+	  t_bits += 16;
+	}
+
       while ( t_bits < n_bits )
 	{
 	  tbits = (tbits * 256) + *bits++;
@@ -4441,6 +4443,54 @@ void decode_double_array_common(unsigned char *igrib, long jlend, int NumBits,
     fpdata[i] = fmin + zscale*fpdata[i];
 }
 
+static unsigned int mask[] = {0,1,3,7,15,31,63,127,255};
+static double shift[9] = {1.0, 2.0, 4.0, 8.0, 16.0, 32.0, 64.0, 128.0, 256.0};
+
+static 
+void decode_double_array_common2(unsigned char *igrib, long jlend, int NumBits, 
+				 double fmin, double zscale, double *fpdata)
+{
+  /* code from wgrib routine BDS_unpack */
+  unsigned char *bits = igrib;
+  unsigned int map_mask;
+  long i;
+  int n_bits = NumBits;
+  int c_bits, j_bits;
+  double jj;
+    
+  /* older unoptimized code, not often used */
+  c_bits = 8;
+  map_mask = 128;
+  for ( i = 0; i < jlend; i++ )
+    {
+      jj = 0.0;
+      j_bits = n_bits;
+      while (c_bits <= j_bits)
+	{
+	  if (c_bits == 8)
+	    {
+	      jj = jj * 256.0  + (double) (*bits++);
+	      j_bits -= 8;
+	    }
+	  else
+	    {
+	      jj = (jj * shift[c_bits]) + (double) (*bits & mask[c_bits]);
+	      bits++;
+	      j_bits -= c_bits;
+	      c_bits = 8;
+	    }
+	}
+
+      if (j_bits)
+	{
+	  c_bits -= j_bits;
+	  jj = (jj * shift[j_bits]) + (double) ((*bits >> c_bits) & mask[j_bits]);
+	}
+      
+      fpdata[i] = fmin + zscale*jj;
+    }
+} 
+
 static 
 void decode_double_array(unsigned char *igrib, long jlend, int numBits, 
 			 double fmin, double zscale, double *fpdata)
@@ -4498,9 +4548,13 @@ void decode_double_array(unsigned char *igrib, long jlend, int numBits,
     {
       decode_double_array_common(igrib, jlend, numBits, fmin, zscale, fpdata);
     }
+  else if ( numBits > 25 && numBits < 32 )
+    {
+      decode_double_array_common2(igrib, jlend, numBits, fmin, zscale, fpdata);
+    }
   else
     {
-      fprintf(stderr," Unimplemented packing factor %d\n", numBits);
+      fprintf(stderr," Unimplemented packing factor %d!\n", numBits);
       exit(EXIT_FAILURE);
     }
 
@@ -4542,9 +4596,13 @@ void decode_double_array(unsigned char *igrib, long jlend, int numBits,
     {
       decode_double_array_common(igrib, jlend, numBits, fmin, zscale, fpdata);
     }
+  else if ( numBits > 25 && numBits < 32 )
+    {
+      decode_double_array_common2(igrib, jlend, numBits, fmin, zscale, fpdata);
+    }
   else
     {
-      fprintf(stderr," Unimplemented packing factor %d\n", numBits);
+      fprintf(stderr," Unimplemented packing factor %d!\n", numBits);
       exit(EXIT_FAILURE);
     }
 #endif
@@ -5472,57 +5530,85 @@ int gribReadSize(int fileID)
 
       pdssize = gribsize;
       fileSetPos(fileID, (off_t) 3, SEEK_CUR);
-      if ( CGRIBEX_Debug )
-	Message("pdssize     = %d", pdssize);
+      if ( CGRIBEX_Debug ) Message("pdssize     = %d", pdssize);
       flag = filePtrGetc(fileptr);
-      if ( CGRIBEX_Debug )
-	Message("flag        = %d", flag);
+      if ( CGRIBEX_Debug ) Message("flag        = %d", flag);
   
       fileSetPos(fileID, (off_t) pdssize-8, SEEK_CUR);
 
       if ( flag & 128 )
 	{
-	  b1 = filePtrGetc(fileptr);
-	  b2 = filePtrGetc(fileptr);
-	  b3 = filePtrGetc(fileptr);
+	  b1 = filePtrGetc(fileptr); b2 = filePtrGetc(fileptr); b3 = filePtrGetc(fileptr);
 	  gdssize = (b1 << 16) + (b2 << 8) + b3;
 	  fileSetPos(fileID, (off_t) gdssize-3, SEEK_CUR);
+	  if ( CGRIBEX_Debug ) Message("gdssize     = %d", gdssize);
 	}
-      if ( CGRIBEX_Debug )
-	Message("gdssize     = %d", gdssize);
 
       if ( flag & 64 )
 	{
-	  b1 = filePtrGetc(fileptr);
-	  b2 = filePtrGetc(fileptr);
-	  b3 = filePtrGetc(fileptr);
+	  b1 = filePtrGetc(fileptr); b2 = filePtrGetc(fileptr); b3 = filePtrGetc(fileptr);
 	  bmssize = (b1 << 16) + (b2 << 8) + b3;
 	  fileSetPos(fileID, (off_t) bmssize-3, SEEK_CUR);
+	  if ( CGRIBEX_Debug ) Message("bmssize     = %d", bmssize);
 	}
-      if ( CGRIBEX_Debug )
-	Message("bmssize     = %d", bmssize);
 
-      b1 = filePtrGetc(fileptr);
-      b2 = filePtrGetc(fileptr);
-      b3 = filePtrGetc(fileptr);
+      b1 = filePtrGetc(fileptr); b2 = filePtrGetc(fileptr); b3 = filePtrGetc(fileptr);
       bdssize = (b1 << 16) + (b2 << 8) + b3;
-
-      if ( CGRIBEX_Debug )
-	Message("bdssize     = %d", bdssize);
+      if ( CGRIBEX_Debug ) Message("bdssize     = %d", bdssize);
 
       gribsize = issize + pdssize + gdssize + bmssize + bdssize + essize;
     }
+  else if ( gribversion == 1 )
+    {
+      if ( gribsize > JP23SET ) /* Large GRIB record */
+	{
+	  int pdssize = 0, gdssize = 0, bmssize = 0, bdssize = 0;
+	  int issize = 4, essize = 4;
+	  int flag;
+
+	  b1 = filePtrGetc(fileptr); b2 = filePtrGetc(fileptr); b3 = filePtrGetc(fileptr);
+	  pdssize = (b1 << 16) + (b2 << 8) + b3;
+	  if ( CGRIBEX_Debug ) Message("pdssize     = %d", pdssize);
+
+	  for ( int i = 0; i < 5; ++i ) flag = filePtrGetc(fileptr);
+	  if ( CGRIBEX_Debug ) Message("flag        = %d", flag);
+  
+	  fileSetPos(fileID, (off_t) pdssize-8, SEEK_CUR);
+
+	  if ( flag & 128 )
+	    {
+	      b1 = filePtrGetc(fileptr); b2 = filePtrGetc(fileptr); b3 = filePtrGetc(fileptr);
+	      gdssize = (b1 << 16) + (b2 << 8) + b3;
+	      fileSetPos(fileID, (off_t) gdssize-3, SEEK_CUR);
+	      if ( CGRIBEX_Debug ) Message("gdssize     = %d", gdssize);
+	    }
+	  
+	  if ( flag & 64 )
+	    {
+	      b1 = filePtrGetc(fileptr); b2 = filePtrGetc(fileptr); b3 = filePtrGetc(fileptr);
+	      bmssize = (b1 << 16) + (b2 << 8) + b3;
+	      fileSetPos(fileID, (off_t) bmssize-3, SEEK_CUR);
+	      if ( CGRIBEX_Debug ) Message("bmssize     = %d", bmssize);
+	    }
+
+	  b1 = filePtrGetc(fileptr); b2 = filePtrGetc(fileptr); b3 = filePtrGetc(fileptr);
+	  bdssize = (b1 << 16) + (b2 << 8) + b3;
+	  bdssize = correct_bdslen(bdssize, gribsize, issize+pdssize+gdssize+bmssize);
+	  if ( CGRIBEX_Debug ) Message("bdssize     = %d", bdssize);
+
+	  gribsize = issize+pdssize+gdssize+bmssize+bdssize+essize;
+	}
+    }
   else if ( gribversion == 2 )
     {
       int i;
       /* we set gribsize the following way because it doesn't matter then
 	 whether int is 4 or 8 bytes long - we don't have to care if the size
 	 really fits: if it does not, the record can not be read at all */
-
       gribsize = 0;
       for ( i = 0; i < 8; i++ ) gribsize = (gribsize << 8) | filePtrGetc(fileptr);
     }
-  else if ( gribversion != 1 )
+  else
     {
       gribsize = 0;
       Warning("GRIB version %d unsupported!", gribversion);
@@ -9086,18 +9172,9 @@ int  gribUnzip(unsigned char *dbuf, long dbufsize, unsigned char *sbuf, long sbu
 
   return (gribLen);
 }
-static const char grb_libvers[] = "1.5.0" " of ""Feb 15 2011"" ""10:23:39";
+static const char grb_libvers[] = "1.5.1" " of ""Aug 29 2011"" ""20:30:27";
 const char *
 cgribexLibraryVersion(void)
 {
   return (grb_libvers);
 }
-/*
- * Local Variables:
- * c-file-style: "Java"
- * c-basic-offset: 2
- * indent-tabs-mode: nil
- * show-trailing-whitespace: t
- * require-trailing-newline: t
- * End:
- */
diff --git a/src/config.h.in b/src/config.h.in
index c1613ccc2026dd8503e4f147cd213d262557ea00..600b1ec8650f7ba0793a1d7110ee918b0ab13d7c 100644
--- a/src/config.h.in
+++ b/src/config.h.in
@@ -51,6 +51,9 @@
 /* Define to 1 if you have the `malloc' library (-lmalloc). */
 #undef HAVE_LIBMALLOC
 
+/* Define to 1 for NETCDF OpenDAP */
+#undef HAVE_LIBNC_DAP
+
 /* Define to 1 for NETCDF support */
 #undef HAVE_LIBNETCDF
 
diff --git a/src/gribapi.h b/src/gribapi.h
index a662656087141cae16389462ab1a8cb011849a0c..ab1899f4a0dece01150d2c7c4563f82ee95f8748 100644
--- a/src/gribapi.h
+++ b/src/gribapi.h
@@ -5,6 +5,9 @@
 
 /* GRIB2 Level Types */
 #define  GRIB2_LTYPE_SURFACE               1
+#define  GRIB2_LTYPE_TOA           8
+#define  GRIB2_LTYPE_SEA_BOTTOM            9
+#define  GRIB2_LTYPE_ATMOSPHERE           10
 #define  GRIB2_LTYPE_ISOBARIC            100
 #define  GRIB2_LTYPE_MEANSEA             101
 #define  GRIB2_LTYPE_ALTITUDE            102
diff --git a/src/grid.c b/src/grid.c
index ef2085454c56246b5b862a3fcfe6e642065814c3..ad2b9197d870e417fb7c32a1e62b0aa97a4bd858 100644
--- a/src/grid.c
+++ b/src/grid.c
@@ -2222,41 +2222,42 @@ void grid_check_cyclic(grid_t *gridptr)
         }
 
       if ( xbounds && xsize > 1 )
-        {
-          double val1, val2;
+	{
+	  double val1, val2;
 
-          gridptr->isCyclic = TRUE;
-          for ( j = 0; j < ysize; ++j )
-            {
-              i1 = j*xsize*4;
-              i2 = j*xsize*4+(xsize-1)*4;
-              nc = 0;
-              for ( k1 = 0; k1 < 4; ++k1 )
-                {
-                  val1 = xbounds[i1+k1];
-                  for ( k2 = 0; k2 < 4; ++k2 )
-                    {
-                      val2 = xbounds[i2+k2];
-
-                      if ( val1 < 1 && val2 > 300 ) val1 += 360;
-                      if ( val2 < 1 && val1 > 300 ) val2 += 360;
-                      if ( val1 < -179 && val2 > 120 ) val1 += 360;
-                      if ( val2 < -179 && val1 > 120 ) val2 += 360;
-
-                      if ( fabs(val1-val2) < 0.001 )
-                        {
-                          nc++;
-                          break;
-                        }
-                    }
-                }
-              if ( nc < 1 )
-                {
-                  gridptr->isCyclic = FALSE;
-                  break;
-                }
-            }
-        }
+	  gridptr->isCyclic = TRUE;
+	  for ( j = 0; j < ysize; ++j )
+	    {
+	      i1 = j*xsize*4;
+	      i2 = j*xsize*4+(xsize-1)*4;
+	      nc = 0;
+	      for ( k1 = 0; k1 < 4; ++k1 )
+		{
+		  val1 = xbounds[i1+k1];
+		  for ( k2 = 0; k2 < 4; ++k2 )
+		    {
+		      val2 = xbounds[i2+k2];
+
+		      if ( val1 < 1 && val2 > 300 ) val1 += 360;
+		      if ( val2 < 1 && val1 > 300 ) val2 += 360;
+		      if ( val1 < -179 && val2 > 120 ) val1 += 360;
+		      if ( val2 < -179 && val1 > 120 ) val2 += 360;
+
+		      if ( fabs(val1-val2) < 0.001 )
+			{
+			  nc++;
+			  break;
+			}
+		    }
+		}
+
+	      if ( nc < 1 )
+		{
+		  gridptr->isCyclic = FALSE;
+		  break;
+		}
+	    }
+	}
     }
 }
 
diff --git a/src/stream.c b/src/stream.c
index f53dcb7e12e3055098bf796de2377a82c3780077..74e5a86804647175f1a993e1703ebcb63a277dec 100644
--- a/src/stream.c
+++ b/src/stream.c
@@ -2011,11 +2011,7 @@ void cdiPrintVersion(void)
   fprintf(stderr, "GRIB_API library version : %s\n", gribapiLibraryVersion());
 #endif
 #if  defined  (HAVE_LIBNETCDF)
-#  if  defined  (HAVE_LIBNC_DAP)
-  fprintf(stderr, "  nc_dap library version : %s\n", cdfLibraryVersion());
-#  else
   fprintf(stderr, "  netCDF library version : %s\n", cdfLibraryVersion());
-#  endif
 #endif
 #if  defined  (HAVE_LIBHDF5)
   fprintf(stderr, "    HDF5 library version : %s\n", hdfLibraryVersion());
diff --git a/src/stream_cdf.c b/src/stream_cdf.c
index 19794b46dfc8d98e4d5ce22e19eedd1a092fcea9..13e3a41371663df85b293408d6ff944f0c4a02d8 100644
--- a/src/stream_cdf.c
+++ b/src/stream_cdf.c
@@ -88,11 +88,14 @@ typedef struct {
   int      natts;
   int     *atts;
   int      deflate;
+  int      lunsigned;
+  int      lvalidrange;
   size_t   vlen;
   double  *vdata;
   double   missval;
   double   addoffset;
   double   scalefactor;
+  double   validrange[2];
   char     name[MAXNAMELEN];
   char     longname[MAXNAMELEN];
   char     stdname[MAXNAMELEN];
@@ -160,7 +163,7 @@ int isTimeUnits(const char *timeunits)
       if ( *ptu )
         {
           while ( isspace(*ptu) ) ptu++;
-          
+
           if ( memcmp(ptu, "as", 2) == 0 )
             timetype = TAXIS_ABSOLUTE;
           else if ( memcmp(ptu, "since", 5) == 0 )
@@ -275,7 +278,7 @@ int splitBasetime(const char *timeunits, taxis_t *taxis)
               (*taxis).rtime = rtime;
 
               if ( CDI_Debug )
-                Message("rdate = %d  rtime = %d", rdate, rtime);        
+                Message("rdate = %d  rtime = %d", rdate, rtime);
             }
         }
     }
@@ -366,11 +369,14 @@ void cdfGetAttText(int fileID, int ncvarid, char *attname, int attlen, char *att
 
 #if  defined  (HAVE_LIBNETCDF)
 static
-int cdfInqDatatype(int xtype)
+int cdfInqDatatype(int xtype, int lunsigned)
 {
   int datatype = -1;
 
+  if ( xtype == NC_BYTE && lunsigned ) xtype = NC_UBYTE;
+
   if      ( xtype == NC_BYTE   )  datatype = DATATYPE_INT8;
+  else if ( xtype == NC_UBYTE  )  datatype = DATATYPE_UINT8;
   /* else if ( xtype == NC_CHAR   )  datatype = DATATYPE_UINT8; */
   else if ( xtype == NC_SHORT  )  datatype = DATATYPE_INT16;
   else if ( xtype == NC_INT    )  datatype = DATATYPE_INT32;
@@ -378,7 +384,6 @@ int cdfInqDatatype(int xtype)
   else if ( xtype == NC_DOUBLE )  datatype = DATATYPE_FLT64;
 #if  defined  (HAVE_NETCDF4)
   else if ( xtype == NC_LONG   )  datatype = DATATYPE_INT32;
-  else if ( xtype == NC_UBYTE  )  datatype = DATATYPE_UINT8;
   else if ( xtype == NC_USHORT )  datatype = DATATYPE_UINT16;
   else if ( xtype == NC_UINT   )  datatype = DATATYPE_UINT32;
   else if ( xtype == NC_INT64  )  datatype = DATATYPE_FLT64;
@@ -408,7 +413,7 @@ int cdfDefDatatype(int datatype, int filetype)
       else if ( datatype == DATATYPE_UINT16 ) xtype = NC_USHORT;
       else if ( datatype == DATATYPE_UINT32 ) xtype = NC_UINT;
 #else
-      else if ( datatype == DATATYPE_UINT8  ) xtype = NC_SHORT;
+      else if ( datatype == DATATYPE_UINT8  ) xtype = NC_BYTE;
       else if ( datatype == DATATYPE_UINT16 ) xtype = NC_INT;
       else if ( datatype == DATATYPE_UINT32 ) xtype = NC_INT;
 #endif
@@ -420,7 +425,7 @@ int cdfDefDatatype(int datatype, int filetype)
       if      ( datatype == DATATYPE_INT8   ) xtype = NC_BYTE;
       else if ( datatype == DATATYPE_INT16  ) xtype = NC_SHORT;
       else if ( datatype == DATATYPE_INT32  ) xtype = NC_INT;
-      else if ( datatype == DATATYPE_UINT8  ) xtype = NC_SHORT;
+      else if ( datatype == DATATYPE_UINT8  ) xtype = NC_BYTE;
       else if ( datatype == DATATYPE_UINT16 ) xtype = NC_INT;
       else if ( datatype == DATATYPE_UINT32 ) xtype = NC_INT;
       else if ( datatype == DATATYPE_FLT64  ) xtype = NC_DOUBLE;
@@ -441,7 +446,6 @@ void defineAttributes(int vlistID, int varID, int fileID, int ncvarID)
   int atttype, attlen;
   size_t len;
   char attname[1024];
-  char atttxt[8192];
 
   vlistInqNatts(vlistID, varID, &natts);
 
@@ -451,9 +455,12 @@ void defineAttributes(int vlistID, int varID, int fileID, int ncvarID)
 
       if ( atttype == DATATYPE_TXT )
         {
-          vlistInqAttTxt(vlistID, varID, attname, sizeof(atttxt), atttxt);
+          char *atttxt;
+          atttxt = (char *) malloc(attlen*sizeof(char));
+          vlistInqAttTxt(vlistID, varID, attname, attlen, atttxt);
           len = attlen;
           cdf_put_att_text(fileID, ncvarID, attname, len, atttxt);
+          free(atttxt);
         }
       else if ( atttype == DATATYPE_INT16 || atttype == DATATYPE_INT32 )
         {
@@ -515,7 +522,7 @@ int cdfCopyRecord(int streamID2, int streamID1)
 
   gridID = vlistInqVarGrid(vlistID1, ivarID);
 
-  datasize = gridInqSize(gridID);      
+  datasize = gridInqSize(gridID);
   /* bug fix for constant netCDF fields */
   if ( datasize < 1048576 ) datasize = 1048576;
 
@@ -684,7 +691,7 @@ void cdfDefVarSzip(int ncid, int ncvarid)
             {
               lwarn = FALSE;
               Warning("netCDF4/Szip compression not compiled in!");
-            }     
+            }
         }
       else
         Error("nc_def_var_szip failed, status = %d", retval);
@@ -716,10 +723,10 @@ void cdfDefVarMissval(int streamID, int varID, int dtype, int lcheck)
 
       xtype = cdfDefDatatype(dtype, streamptr->filetype);
 
-      cdf_put_att_double(fileID, ncvarid, "_FillValue", (nc_type) xtype, 1L, &missval);
+      cdf_put_att_double(fileID, ncvarid, "_FillValue", (nc_type) xtype, 1, &missval);
 
       if ( cdiNcMissingValue == 1 )
-        cdf_put_att_double(fileID, ncvarid, "missing_value", (nc_type) xtype, 1L, &missval);
+        cdf_put_att_double(fileID, ncvarid, "missing_value", (nc_type) xtype, 1, &missval);
 
       if ( lcheck && streamptr->ncmode == 2 ) cdf_enddef(fileID);
 
@@ -1309,10 +1316,10 @@ int checkGridName(int type, char *axisname, int fileID, int vlistID, int gridID,
         }
 
       if ( checkname ) iz++;
-      
+
       if ( iz > 99 ) break;
     }
-  
+
   if ( iz ) sprintf(&axisname[strlen(axisname)], "_%d", iz+1);
 
   return (iz);
@@ -1432,7 +1439,7 @@ void cdfDefXaxis(int streamID, int gridID)
           if ( gridIsRotated(gridID) )
             {
               double north_pole = gridInqXpole(gridID);
-              cdf_put_att_double(fileID, ncvarid, "north_pole", NC_DOUBLE, 1L, &north_pole);
+              cdf_put_att_double(fileID, ncvarid, "north_pole", NC_DOUBLE, 1, &north_pole);
             }
           */
         }
@@ -1561,7 +1568,7 @@ void cdfDefYaxis(int streamID, int gridID)
           if ( gridIsRotated(gridID) )
             {
               double north_pole = gridInqYpole(gridID);
-              cdf_put_att_double(fileID, ncvarid, "north_pole", NC_DOUBLE, 1L, &north_pole);
+              cdf_put_att_double(fileID, ncvarid, "north_pole", NC_DOUBLE, 1, &north_pole);
             }
           */
         }
@@ -1996,7 +2003,7 @@ void cdfDefUnstructured(int streamID, int gridID)
 
       nvertex = gridInqNvertex(gridID);
       if ( nvertex > 0 ) cdf_def_dim(fileID, vertname, nvertex, &nvdimID);
-        
+
       if ( gridInqXvalsPtr(gridID) )
         {
           cdf_def_var(fileID, xaxisname, (nc_type) xtype, 1, &dimID, &ncxvarid);
@@ -2092,7 +2099,7 @@ void cdfDefVCT(int streamID, int zaxisID)
       int ilev = zaxisInqVctSize(zaxisID)/2;
       int mlev = ilev - 1;
       size_t start;
-      size_t count = 1;      
+      size_t count = 1;
       int ncdimid, ncdimid2;
       int hyaiid, hybiid, hyamid, hybmid;
       double mval;
@@ -2207,8 +2214,11 @@ void cdfDefZaxis(int streamID, int zaxisID)
   dimlen = zaxisInqSize(zaxisID);
   type   = zaxisInqType(zaxisID);
 
-  if ( dimlen == 1 && type == ZAXIS_SURFACE ) return;
-  if ( dimlen == 1 && type == ZAXIS_MEANSEA ) return;
+  if ( dimlen == 1 && type == ZAXIS_SURFACE     ) return;
+  if ( dimlen == 1 && type == ZAXIS_TOA         ) return;
+  if ( dimlen == 1 && type == ZAXIS_SEA_BOTTOM  ) return;
+  if ( dimlen == 1 && type == ZAXIS_ATMOSPHERE  ) return;
+  if ( dimlen == 1 && type == ZAXIS_MEANSEA     ) return;
 
   zaxisInqName(zaxisID, axisname);
   /*
@@ -2318,7 +2328,7 @@ void cdfDefZaxis(int streamID, int zaxisID)
 	      cdf_put_att_text(fileID, ncvarid, "formula", strlen(tmpname), tmpname);
 	      strcpy(tmpname, "ap: hyai b: hybi ps: aps");
 	      cdf_put_att_text(fileID, ncvarid, "formula_terms", strlen(tmpname), tmpname);
-              
+
 	      cdf_enddef(fileID);
 	      streamptr->ncmode = 2;
 
@@ -2367,8 +2377,8 @@ void cdfDefZaxis(int streamID, int zaxisID)
 	    }
 
           cdf_put_att_text(fileID, ncvarid, "axis", 1, "Z");
-              
-	  if ( zaxisInqLbounds(zaxisID, NULL) && zaxisInqUbounds(zaxisID, NULL) ) 
+
+	  if ( zaxisInqLbounds(zaxisID, NULL) && zaxisInqUbounds(zaxisID, NULL) )
             {
 	      if ( nc_inq_dimid(fileID, "nb2", &nvdimID) != NC_NOERR )
 		cdf_def_dim(fileID, "nb2", nvertex, &nvdimID);
@@ -2387,7 +2397,7 @@ void cdfDefZaxis(int streamID, int zaxisID)
           streamptr->ncmode = 2;
 
           cdf_put_var_double(fileID, ncvarid, zaxisInqLevelsPtr(zaxisID));
-	  
+
           if ( ncbvarid != UNDEFID )
 	    {
 	      int i;
@@ -2443,10 +2453,10 @@ void cdfDefPole(int streamID, int gridID)
   if ( ncerr == NC_NOERR )
     {
       cdf_put_att_text(fileID, ncvarid, "grid_mapping_name", strlen(mapname), mapname);
-      cdf_put_att_double(fileID, ncvarid, "grid_north_pole_latitude", NC_DOUBLE, 1L, &ypole);
-      cdf_put_att_double(fileID, ncvarid, "grid_north_pole_longitude", NC_DOUBLE, 1L, &xpole);
+      cdf_put_att_double(fileID, ncvarid, "grid_north_pole_latitude", NC_DOUBLE, 1, &ypole);
+      cdf_put_att_double(fileID, ncvarid, "grid_north_pole_longitude", NC_DOUBLE, 1, &xpole);
       if ( angle > 0 )
-        cdf_put_att_double(fileID, ncvarid, "north_pole_grid_longitude", NC_DOUBLE, 1L, &angle);
+        cdf_put_att_double(fileID, ncvarid, "north_pole_grid_longitude", NC_DOUBLE, 1, &angle);
     }
 
   cdf_enddef(fileID);
@@ -2475,8 +2485,8 @@ void cdfDefMapping(int streamID, int gridID)
         {
           cdf_put_att_text(fileID, ncvarid, "grid_mapping_name", strlen(mapname), mapname);
           /*
-          cdf_put_att_double(fileID, ncvarid, "grid_north_pole_latitude", NC_DOUBLE, 1L, &ypole);
-          cdf_put_att_double(fileID, ncvarid, "grid_north_pole_longitude", NC_DOUBLE, 1L, &xpole);
+          cdf_put_att_double(fileID, ncvarid, "grid_north_pole_latitude", NC_DOUBLE, 1, &ypole);
+          cdf_put_att_double(fileID, ncvarid, "grid_north_pole_longitude", NC_DOUBLE, 1, &xpole);
           */
         }
 
@@ -2499,9 +2509,9 @@ void cdfDefMapping(int streamID, int gridID)
           gridInqLaea(gridID, &a, &lon_0, &lat_0);
 
           cdf_put_att_text(fileID, ncvarid, "grid_mapping_name", strlen(mapname), mapname);
-          cdf_put_att_double(fileID, ncvarid, "earth_radius", NC_DOUBLE, 1L, &a);
-          cdf_put_att_double(fileID, ncvarid, "longitude_of_projection_origin", NC_DOUBLE, 1L, &lon_0);
-          cdf_put_att_double(fileID, ncvarid, "latitude_of_projection_origin", NC_DOUBLE, 1L, &lat_0);
+          cdf_put_att_double(fileID, ncvarid, "earth_radius", NC_DOUBLE, 1, &a);
+          cdf_put_att_double(fileID, ncvarid, "longitude_of_projection_origin", NC_DOUBLE, 1, &lon_0);
+          cdf_put_att_double(fileID, ncvarid, "latitude_of_projection_origin", NC_DOUBLE, 1, &lat_0);
         }
 
       cdf_enddef(fileID);
@@ -2524,17 +2534,17 @@ void cdfDefMapping(int streamID, int gridID)
 
           cdf_put_att_text(fileID, ncvarid, "grid_mapping_name", strlen(mapname), mapname);
           if ( radius > 0 )
-            cdf_put_att_double(fileID, ncvarid, "earth_radius", NC_DOUBLE, 1L, &radius);
-          cdf_put_att_double(fileID, ncvarid, "longitude_of_central_meridian", NC_DOUBLE, 1L, &lon_0);
-          cdf_put_att_double(fileID, ncvarid, "latitude_of_projection_origin", NC_DOUBLE, 1L, &lat_0);
+            cdf_put_att_double(fileID, ncvarid, "earth_radius", NC_DOUBLE, 1, &radius);
+          cdf_put_att_double(fileID, ncvarid, "longitude_of_central_meridian", NC_DOUBLE, 1, &lon_0);
+          cdf_put_att_double(fileID, ncvarid, "latitude_of_projection_origin", NC_DOUBLE, 1, &lat_0);
           if ( IS_EQUAL(lat_1, lat_2) )
-            cdf_put_att_double(fileID, ncvarid, "standard_parallel", NC_DOUBLE, 1L, &lat_1);
+            cdf_put_att_double(fileID, ncvarid, "standard_parallel", NC_DOUBLE, 1, &lat_1);
           else
             {
               double lat_1_2[2];
               lat_1_2[0] = lat_1;
               lat_1_2[1] = lat_2;
-              cdf_put_att_double(fileID, ncvarid, "standard_parallel", NC_DOUBLE, 2L, lat_1_2);
+              cdf_put_att_double(fileID, ncvarid, "standard_parallel", NC_DOUBLE, 2, lat_1_2);
             }
         }
 
@@ -2654,7 +2664,7 @@ int cdfDefVar(int streamID, int varID)
   int ndims = 0;
   int len;
   int timeID;
-  int xtype;
+  int xtype, dtype;
   int gridtype, gridsize;
   int gridindex, zaxisindex;
   int tablenum;
@@ -2705,8 +2715,7 @@ int cdfDefVar(int streamID, int varID)
   dimorder[0] = ixyz/100;
   dimorder[1] = (ixyz-dimorder[0]*100)/10;
   dimorder[2] = (ixyz-dimorder[0]*100-dimorder[1]*10);
-
-  if ( dimorder[0] == 3 ) lchunk = FALSE; /* ZYX and ZXY */
+  if ( dimorder[0] != 3 ) lchunk = FALSE; /* ZYX and ZXY */
 
   if ( ((dimorder[0]>0)+(dimorder[1]>0)+(dimorder[2]>0)) < ((xid!=UNDEFID)+(yid!=UNDEFID)+(zid!=UNDEFID)) )
     {
@@ -2847,7 +2856,8 @@ int cdfDefVar(int streamID, int varID)
 
   /* if ( streamptr->ncmode == 2 ) cdf_redef(fileID); */
 
-  xtype = cdfDefDatatype(vlistInqVarDatatype(vlistID, varID), streamptr->filetype);
+  dtype = vlistInqVarDatatype(vlistID, varID);
+  xtype = cdfDefDatatype(dtype, streamptr->filetype);
 
   cdf_def_var(fileID, name, (nc_type) xtype, ndims, dims, &ncvarid);
 
@@ -2867,12 +2877,15 @@ int cdfDefVar(int streamID, int varID)
         }
       else
         {
-          static int lwarn = TRUE;
-
-          if ( lwarn )
+          if ( lchunk )
             {
-              lwarn = FALSE;
-              Warning("Deflate compression is only available for netCDF4!");
+              static int lwarn = TRUE;
+
+              if ( lwarn )
+                {
+                  lwarn = FALSE;
+                  Warning("Deflate compression is only available for netCDF4!");
+                }
             }
         }
     }
@@ -2915,7 +2928,7 @@ int cdfDefVar(int streamID, int varID)
     cdf_put_att_text(fileID, ncvarid, "units", strlen(units), units);
 
   if ( code > 0 && pdis == 255 )
-    cdf_put_att_int(fileID, ncvarid, "code", NC_INT, 1L, &code);
+    cdf_put_att_int(fileID, ncvarid, "code", NC_INT, 1, &code);
 
   if ( pdis != 255 )
     {
@@ -2928,7 +2941,7 @@ int cdfDefVar(int streamID, int varID)
     {
       tablenum = tableInqNum(tableID);
       if ( tablenum > 0 )
-        cdf_put_att_int(fileID, ncvarid, "table", NC_INT, 1L, &tablenum);
+        cdf_put_att_int(fileID, ncvarid, "table", NC_INT, 1, &tablenum);
     }
 
   if ( gridtype != GRID_GENERIC && gridtype != GRID_LONLAT  &&
@@ -3003,7 +3016,7 @@ int cdfDefVar(int streamID, int varID)
       axis[iax++] = '-';
       axis[iax++] = '-';
       cdf_put_att_text(fileID, ncvarid, "axis", iax, axis);
-      cdf_put_att_int(fileID, ncvarid, "truncation", NC_INT, 1L, &gridTruncation);
+      cdf_put_att_int(fileID, ncvarid, "truncation", NC_INT, 1, &gridTruncation);
     }
 
   /*  if ( xtype == NC_BYTE || xtype == NC_SHORT || xtype == NC_INT ) */
@@ -3027,16 +3040,34 @@ int cdfDefVar(int streamID, int varID)
 
           if ( xtype == (int) NC_FLOAT ) astype = NC_FLOAT;
 
-          cdf_put_att_double(fileID, ncvarid, "add_offset",   (nc_type) astype, 1L, &addoffset);
-          cdf_put_att_double(fileID, ncvarid, "scale_factor", (nc_type) astype, 1L, &scalefactor);
+          cdf_put_att_double(fileID, ncvarid, "add_offset",   (nc_type) astype, 1, &addoffset);
+          cdf_put_att_double(fileID, ncvarid, "scale_factor", (nc_type) astype, 1, &scalefactor);
         }
     }
 
+  if ( dtype == DATATYPE_UINT8 && xtype == NC_BYTE )
+    {
+      int validrange[2] = {0, 255};
+      cdf_put_att_int(fileID, ncvarid, "valid_range", NC_SHORT, 2, validrange);
+      cdf_put_att_text(fileID, ncvarid, "_Unsigned", 4, "true");
+    }
+
   streamptr->vars[varID].ncvarid = ncvarid;
 
   if ( vlistInqVarMissvalUsed(vlistID, varID) )
     cdfDefVarMissval(streamID, varID, vlistInqVarDatatype(vlistID, varID), 0);
 
+  if ( zid == -1 )
+    {
+      if ( zaxisInqType(zaxisID) == ZAXIS_TOA         || 
+           zaxisInqType(zaxisID) == ZAXIS_SEA_BOTTOM  ||
+           zaxisInqType(zaxisID) == ZAXIS_ATMOSPHERE )
+        {
+          zaxisInqName(zaxisID, varname);
+          cdf_put_att_text(fileID, ncvarid, "level_type", strlen(varname), varname);
+        }
+    }
+
   /* Attributes */
   defineAttributes(vlistID, varID, fileID, ncvarid);
 
@@ -3046,6 +3077,27 @@ int cdfDefVar(int streamID, int varID)
 }
 #endif
 
+static
+void scale_add(long size, double *data, double addoffset, double scalefactor)
+{
+  long i;
+  int laddoffset;
+  int lscalefactor;
+
+  laddoffset   = IS_NOT_EQUAL(addoffset, 0);
+  lscalefactor = IS_NOT_EQUAL(scalefactor, 1);
+
+  if ( laddoffset || lscalefactor )
+    {
+      for ( i = 0; i < size; ++i )
+        {
+          if ( lscalefactor ) data[i] *= scalefactor;
+          if ( laddoffset )   data[i] += addoffset;
+        }
+    }
+}
+
+
 void cdfReadVarDP(int streamID, int varID, double *data, int *nmiss)
 {
 #if  defined  (HAVE_LIBNETCDF)
@@ -3365,6 +3417,8 @@ int cdfReadVarSliceDP(int streamID, int varID, int levelID, double *data, int *n
   int dimorder[3];
   int ixyz;
   int swapxy = FALSE;
+  int lvalidrange;
+  double validrange[2];
   double missval;
   int laddoffset, lscalefactor;
   double addoffset, scalefactor;
@@ -3475,10 +3529,31 @@ int cdfReadVarSliceDP(int streamID, int varID, int levelID, double *data, int *n
       free(tdata);
     }
 
+  if ( vlistInqVarDatatype(vlistID, varID) == DATATYPE_UINT8 )
+    {
+      nc_type xtype;
+      cdf_inq_vartype(fileID, ncvarid, &xtype);
+      if ( xtype == NC_BYTE )
+        {
+          for ( i = 0; i < gridsize; i++ )
+            if ( data[i] < 0 ) data[i] += 256;
+        }
+    }
+
   *nmiss = 0;
-  if ( vlistInqVarMissvalUsed(vlistID, varID) == TRUE  )
+  if ( vlistInqVarMissvalUsed(vlistID, varID) == TRUE )
     {
       missval = vlistInqVarMissval(vlistID, varID);
+
+      lvalidrange = vlistInqVarValidrange(vlistID, varID, validrange);
+      // printf("readvarslice: validrange %d %g %g\n", lvalidrange, validrange[0], validrange[1]);
+      if ( lvalidrange )
+        for ( i = 0; i < gridsize; i++ )
+          {
+            if ( IS_NOT_EQUAL(validrange[0], VALIDMISS) && data[i] < validrange[0] ) data[i] = missval;
+            if ( IS_NOT_EQUAL(validrange[1], VALIDMISS) && data[i] > validrange[1] ) data[i] = missval;
+          }
+
       // printf("XXX %31.0f %31.0f %31.0f %31.0f\n", missval, (float)data[0]);
       for ( i = 0; i < gridsize; i++ )
         if ( DBL_IS_EQUAL(data[i], missval) ) *nmiss += 1;
@@ -3531,6 +3606,7 @@ int cdfWriteVarSliceDP(int streamID, int varID, int levelID, const double *data,
   size_t size, xsize, ysize;
   size_t start[4];
   size_t count[4];
+  int nvals;
   int ndims = 0;
   int idim;
   int timeID;
@@ -3632,9 +3708,10 @@ int cdfWriteVarSliceDP(int streamID, int varID, int levelID, const double *data,
 
   if ( nmiss > 0 ) cdfDefVarMissval(streamID, varID, dtype, 1);
 
+  nvals = gridInqSize(gridID);
+
   /*  if ( dtype == DATATYPE_INT8 || dtype == DATATYPE_INT16 || dtype == DATATYPE_INT32 ) */
     {
-      int nvals;
       int laddoffset, lscalefactor;
       double addoffset, scalefactor;
       double missval;
@@ -3646,8 +3723,6 @@ int cdfWriteVarSliceDP(int streamID, int varID, int levelID, const double *data,
 
       missval      = vlistInqVarMissval(vlistID, varID);
 
-      nvals = gridInqSize(gridID);
-
       if ( laddoffset || lscalefactor )
         {
           mdata = (double *) malloc(nvals*sizeof(double));
@@ -3686,6 +3761,17 @@ int cdfWriteVarSliceDP(int streamID, int varID, int levelID, const double *data,
               mdata = (double *) malloc(nvals*sizeof(double));
               for ( i = 0; i < nvals; i++ ) mdata[i] = NINT(data[i]);
             }
+
+          if ( dtype == DATATYPE_UINT8 )
+            {
+              nc_type xtype;
+              cdf_inq_vartype(fileID, ncvarid, &xtype);
+              if ( xtype == NC_BYTE )
+                {
+                  for ( i = 0; i < nvals; ++i )
+                    if ( mdata[i] > 127 ) mdata[i] -= 256;
+                }
+            }
         }
 
       if ( CDF_Debug )
@@ -3710,8 +3796,6 @@ int cdfWriteVarSliceDP(int streamID, int varID, int levelID, const double *data,
 
   if ( swapxy )
     {
-      int nvals;
-      nvals = gridInqSize(gridID);
       sdata = (double *) malloc(nvals*sizeof(double));
       for ( int j = 0; j < (int)ysize; ++j )
         for ( int i = 0; i < (int)xsize; ++i )
@@ -3927,7 +4011,7 @@ void init_ncvars(long nvars, ncvar_t *ncvars)
   long ncvarid;
 
   for ( ncvarid = 0; ncvarid < nvars; ncvarid++ )
-    {      
+    {
       ncvars[ncvarid].ignore          = FALSE;
       ncvars[ncvarid].isvar           = UNDEFID;
       ncvars[ncvarid].islon           = FALSE;
@@ -3976,6 +4060,10 @@ void init_ncvars(long nvars, ncvar_t *ncvars)
       ncvars[ncvarid].natts           = 0;
       ncvars[ncvarid].atts            = NULL;
       ncvars[ncvarid].deflate         = 0;
+      ncvars[ncvarid].lunsigned       = 0;
+      ncvars[ncvarid].lvalidrange     = 0;
+      ncvars[ncvarid].validrange[0]   = VALIDMISS;
+      ncvars[ncvarid].validrange[1]   = VALIDMISS;
     }
 }
 
@@ -4124,11 +4212,11 @@ int isGaussGrid(long ysize, double yinc, double *yvals)
 
       /* check S->N */
       if ( lgauss == FALSE )
-        {                 
+        {
           for ( i = 0; i < ysize; i++ )
             if ( fabs(yv[i] - yvals[ysize-i-1]) >
                  ((yv[0] - yv[1])/500) ) break;
-      
+
           if ( i == ysize ) lgauss = TRUE;
         }
 
@@ -4152,7 +4240,7 @@ void cdfSetVar(ncvar_t *ncvars, int ncvarid, int isvar)
     {
       if ( ! ncvars[ncvarid].ignore )
         Warning("Inconsistent variable definition for %s!", ncvars[ncvarid].name);
-        
+
       ncvars[ncvarid].warn = TRUE;
       isvar = FALSE;
     }
@@ -4373,7 +4461,30 @@ void scanVarAttributes(int fileID, int nvars, ncvar_t *ncvars, ncdim_t *ncdims,
                   if ( warn )
                     {
                       warn = FALSE;
-                      Warning("Gridtype %s unsupported!", attstring);
+                      Warning("Grid type %s unsupported!", attstring);
+                    }
+                }
+
+              cdfSetVar(ncvars, ncvarid, TRUE);
+            }
+          else if ( strcmp(attname, "level_type") == 0 && atttype == NC_CHAR )
+            {
+              cdfGetAttText(fileID, ncvarid, attname, attstringlen-1, attstring);
+              strtolower(attstring);
+
+              if      ( strcmp(attstring, "toa") == 0 )
+                ncvars[ncvarid].zaxistype = ZAXIS_TOA;
+              else if ( strcmp(attstring, "seabottom") == 0 )
+                ncvars[ncvarid].zaxistype = ZAXIS_SEA_BOTTOM;
+              else if ( strcmp(attstring, "atmosphere") == 0 )
+                ncvars[ncvarid].zaxistype = ZAXIS_ATMOSPHERE;
+              else
+                { 
+                  static int warn = TRUE;
+                  if ( warn )
+                    {
+                      warn = FALSE;
+                      Warning("Zaxis type %s unsupported!", attstring);
                     }
                 }
 
@@ -4395,7 +4506,7 @@ void scanVarAttributes(int fileID, int nvars, ncvar_t *ncvars, ncdim_t *ncdims,
 		if ( ncvars[ncvarid].addoffset != 0 )
 		Warning("attribute add_offset not supported for atttype %d", atttype);
 	      */
-	      cdfSetVar(ncvars, ncvarid, TRUE);
+	      /* (also used for lon/lat) cdfSetVar(ncvars, ncvarid, TRUE); */
             }
           else if ( strcmp(attname, "scale_factor") == 0 && atttype != NC_CHAR )
             {
@@ -4405,7 +4516,7 @@ void scanVarAttributes(int fileID, int nvars, ncvar_t *ncvars, ncdim_t *ncdims,
 		if ( ncvars[ncvarid].scalefactor != 1 )
 		Warning("attribute scale_factor not supported for atttype %d", atttype);
 	      */
-	      cdfSetVar(ncvars, ncvarid, TRUE);
+	      /* (also used for lon/lat) cdfSetVar(ncvars, ncvarid, TRUE); */
             }
           else if ( strcmp(attname, "bounds") == 0 && atttype == NC_CHAR )
             {
@@ -4568,6 +4679,43 @@ void scanVarAttributes(int fileID, int nvars, ncvar_t *ncvars, ncdim_t *ncdims,
 	      ncvars[ncvarid].defmiss = TRUE;
 	      /* cdfSetVar(ncvars, ncvarid, TRUE); */
             }
+          else if ( strcmp(attname, "valid_range") == 0 && attlen == 2 )
+            {
+              if ( ncvars[ncvarid].lvalidrange == FALSE )
+                {
+                  cdfGetAttDouble(fileID, ncvarid, attname, 2, ncvars[ncvarid].validrange);
+                  ncvars[ncvarid].lvalidrange = TRUE;
+                  if ( ((int)ncvars[ncvarid].validrange[0]) == 0 && ((int)ncvars[ncvarid].validrange[1]) == 255 )
+                    ncvars[ncvarid].lunsigned = TRUE;
+                  /* cdfSetVar(ncvars, ncvarid, TRUE); */
+                }
+            }
+          else if ( strcmp(attname, "valid_min") == 0 && attlen == 1 )
+            {
+              cdfGetAttDouble(fileID, ncvarid, attname, 1, &(ncvars[ncvarid].validrange)[0]);
+              ncvars[ncvarid].lvalidrange = TRUE;
+            }
+          else if ( strcmp(attname, "valid_max") == 0 && attlen == 1 )
+            {
+              cdfGetAttDouble(fileID, ncvarid, attname, 1, &(ncvars[ncvarid].validrange)[1]);
+              ncvars[ncvarid].lvalidrange = TRUE;
+            }
+          else if ( strcmp(attname, "_Unsigned") == 0 && atttype == NC_CHAR )
+            {
+              cdfGetAttText(fileID, ncvarid, attname, attstringlen-1, attstring);
+              strtolower(attstring);
+
+              if ( memcmp(attstring, "true", 4) == 0 )
+                {
+                  ncvars[ncvarid].lunsigned = TRUE;
+                  /*
+                  ncvars[ncvarid].lvalidrange = TRUE;
+                  ncvars[ncvarid].validrange[0] = 0;
+                  ncvars[ncvarid].validrange[1] = 255;
+                  */
+                }
+	      /* cdfSetVar(ncvars, ncvarid, TRUE); */
+            }
           else if ( strcmp(attname, "cdi") == 0 && atttype == NC_CHAR )
             {
 	      cdfGetAttText(fileID, ncvarid, attname, attstringlen-1, attstring);
@@ -4797,7 +4945,7 @@ void verify_coordinate_vars_1(int ndims, ncdim_t *ncdims, ncvar_t *ncvars, int t
 		{
 		  ncvars[ncvarid].zaxistype = ZAXIS_PRESSURE;
 		}
-	      else if ( strcmp(ncvars[ncvarid].units, "level") == 0 )
+	      else if ( strcmp(ncvars[ncvarid].units, "level") == 0 || strcmp(ncvars[ncvarid].units, "1") == 0 )
 		{
 		  if      ( strcmp(ncvars[ncvarid].longname, "hybrid level at layer midpoints") == 0 )
 		    ncvars[ncvarid].zaxistype = ZAXIS_HYBRID;
@@ -4807,7 +4955,7 @@ void verify_coordinate_vars_1(int ndims, ncdim_t *ncdims, ncvar_t *ncvars, int t
 		    ncvars[ncvarid].zaxistype = ZAXIS_HYBRID_HALF;
 		  else if ( memcmp(ncvars[ncvarid].longname, "hybrid level at interfaces", 26) == 0 )
 		    ncvars[ncvarid].zaxistype = ZAXIS_HYBRID_HALF;
-		  else
+		  else if ( strcmp(ncvars[ncvarid].units, "level") == 0 )
 		    ncvars[ncvarid].zaxistype = ZAXIS_GENERIC;
 		}
 	      else if ( strcmp(ncvars[ncvarid].units, "cm")  == 0 )
@@ -4878,7 +5026,7 @@ void verify_coordinate_vars_2(int nvars, ncvar_t *ncvars)
 		  ncvars[ncvarid].zaxistype = ZAXIS_PRESSURE;
 		  continue;
 		}
-	      else if ( strcmp(ncvars[ncvarid].units, "level") == 0 )
+	      else if ( strcmp(ncvars[ncvarid].units, "level") == 0 || strcmp(ncvars[ncvarid].units, "1") == 0 )
 		{
 		  if      ( strcmp(ncvars[ncvarid].longname, "hybrid level at layer midpoints") == 0 )
 		    ncvars[ncvarid].zaxistype = ZAXIS_HYBRID;
@@ -4888,7 +5036,7 @@ void verify_coordinate_vars_2(int nvars, ncvar_t *ncvars)
 		    ncvars[ncvarid].zaxistype = ZAXIS_HYBRID_HALF;
 		  else if ( memcmp(ncvars[ncvarid].longname, "hybrid level at interfaces", 26) == 0 )
 		    ncvars[ncvarid].zaxistype = ZAXIS_HYBRID_HALF;
-		  else
+		  else if ( strcmp(ncvars[ncvarid].units, "level") == 0 )
 		    ncvars[ncvarid].zaxistype = ZAXIS_GENERIC;
 		  continue;
 		}
@@ -5084,8 +5232,7 @@ void define_all_grids(stream_t *streamptr, int fileID, int vlistID, ncdim_t *ncd
 			}
 		      else
 			{
-			  Warning("Unsupported grid structure for variable %s (grid dims > 2)!",
-				  ncvars[ncvarid].name);
+			  Warning("Unsupported grid structure for variable %s (grid dims > 2)!", ncvars[ncvarid].name);
 			  ncvars[ncvarid].xvarid = UNDEFID;
 			  ncvars[ncvarid].yvarid = UNDEFID;
 			  xvarid = UNDEFID;
@@ -5112,8 +5259,7 @@ void define_all_grids(stream_t *streamptr, int fileID, int vlistID, ncdim_t *ncd
 			dimsize2 = ncdims[dimid].len;
 			if ( dimsize1*dimsize2 != size )
 			  {
-			    Warning("Unsupported array structure, skipped variable %s!",
-				    ncvars[ncvarid].name);
+			    Warning("Unsupported array structure, skipped variable %s!", ncvars[ncvarid].name);
 			    ncvars[ncvarid].isvar = -1;
 			    continue;
 			  }
@@ -5130,8 +5276,7 @@ void define_all_grids(stream_t *streamptr, int fileID, int vlistID, ncdim_t *ncd
 			dimsize = ncdims[dimid].len;
 			if ( dimsize != size )
 			  {
-			    Warning("Unsupported array structure, skipped variable %s!",
-				    ncvars[ncvarid].name);
+			    Warning("Unsupported array structure, skipped variable %s!", ncvars[ncvarid].name);
 			    ncvars[ncvarid].isvar = -1;
 			    continue;
 			  }
@@ -5146,6 +5291,8 @@ void define_all_grids(stream_t *streamptr, int fileID, int vlistID, ncdim_t *ncd
 		  else
 		    cdf_get_var_double(fileID, xvarid, grid.xvals);
 
+                  scale_add(size, grid.xvals, ncvars[xvarid].addoffset, ncvars[xvarid].scalefactor);
+
 		  strcpy(grid.xname, ncvars[xvarid].name);
 		  strcpy(grid.xlongname, ncvars[xvarid].longname);
 		  strcpy(grid.xunits, ncvars[xvarid].units);
@@ -5183,8 +5330,7 @@ void define_all_grids(stream_t *streamptr, int fileID, int vlistID, ncdim_t *ncd
 			dimsize2 = ncdims[dimid].len;
 			if ( dimsize1*dimsize2 != size )
 			  {
-			    Warning("Unsupported array structure, skipped variable %s!",
-				    ncvars[ncvarid].name);
+			    Warning("Unsupported array structure, skipped variable %s!", ncvars[ncvarid].name);
 			    ncvars[ncvarid].isvar = -1;
 			    continue;
 			  }
@@ -5203,8 +5349,7 @@ void define_all_grids(stream_t *streamptr, int fileID, int vlistID, ncdim_t *ncd
 			dimsize = ncdims[dimid].len;
 			if ( dimsize != size )
 			  {
-			    Warning("Unsupported array structure, skipped variable %s!",
-				    ncvars[ncvarid].name);
+			    Warning("Unsupported array structure, skipped variable %s!", ncvars[ncvarid].name);
 			    ncvars[ncvarid].isvar = -1;
 			    continue;
 			  }
@@ -5219,6 +5364,8 @@ void define_all_grids(stream_t *streamptr, int fileID, int vlistID, ncdim_t *ncd
 		  else
 		    cdf_get_var_double(fileID, yvarid, grid.yvals);
 
+                  scale_add(size, grid.yvals, ncvars[xvarid].addoffset, ncvars[xvarid].scalefactor);
+
 		  strcpy(grid.yname, ncvars[yvarid].name);
 		  strcpy(grid.ylongname, ncvars[yvarid].longname);
 		  strcpy(grid.yunits, ncvars[yvarid].units);
@@ -5624,7 +5771,11 @@ void define_all_zaxes(stream_t *streamptr, int fileID, int vlistID, ncdim_t *ncd
 
 	      if ( zsize == 1 )
 		{
-		  zaxisType = ZAXIS_SURFACE;
+                  if ( ncvars[ncvarid].zaxistype != UNDEFID )
+                    zaxisType = ncvars[ncvarid].zaxistype;
+                  else
+                    zaxisType = ZAXIS_SURFACE;
+
 		  zvar[0] = 0;
 		  /*
 		  if ( zdimid == UNDEFID )
@@ -5746,12 +5897,16 @@ void define_all_vars(int fileID, int streamID, int vlistID, int instID, int mode
       if ( ncvars[ncvarid].longname[0] )      vlistDefVarLongname(vlistID, varID, ncvars[ncvarid].longname);
       if ( ncvars[ncvarid].stdname[0] )       vlistDefVarStdname(vlistID, varID, ncvars[ncvarid].stdname);
       if ( ncvars[ncvarid].units[0] )         vlistDefVarUnits(vlistID, varID, ncvars[ncvarid].units);
+
+      if ( ncvars[ncvarid].lvalidrange )
+        vlistDefVarValidrange(vlistID, varID, ncvars[ncvarid].validrange);
+
       if ( IS_NOT_EQUAL(ncvars[ncvarid].addoffset, 0) )
 	vlistDefVarAddoffset(vlistID, varID, ncvars[ncvarid].addoffset);
       if ( IS_NOT_EQUAL(ncvars[ncvarid].scalefactor, 1) )
 	vlistDefVarScalefactor(vlistID, varID, ncvars[ncvarid].scalefactor);
 
-      vlistDefVarDatatype(vlistID, varID, cdfInqDatatype(ncvars[ncvarid].xtype));
+      vlistDefVarDatatype(vlistID, varID, cdfInqDatatype(ncvars[ncvarid].xtype, ncvars[ncvarid].lunsigned));
 
       vlistDefVarInstitut(vlistID, varID, instID);
       vlistDefVarModel(vlistID, varID, modelID);
@@ -5778,14 +5933,28 @@ void define_all_vars(int fileID, int streamID, int vlistID, int instID, int mode
 
       if ( ncvars[ncvarid].timeID == TIME_VARIABLE ) iodim++;
 
-      for ( int idim = iodim; idim < ndims; idim++ )
+      if ( gridInqType(gridID) == GRID_UNSTRUCTURED && ndims-iodim <= 2 && ydimid == xdimid )
+        {
+          if ( xdimid == ncvars[ncvarid].dimids[ndims-1] )
+            {
+              ixyz = 321;
+            }
+          else
+            {
+              ixyz = 213;
+            }
+        }
+      else
         {
-          if      ( xdimid == ncvars[ncvarid].dimids[idim] )
-            ixyz += 1*ipow10[(ndims-idim-1)];
-          else if ( ydimid == ncvars[ncvarid].dimids[idim] )
-            ixyz += 2*ipow10[(ndims-idim-1)];
-          else if ( zdimid == ncvars[ncvarid].dimids[idim] )
-            ixyz += 3*ipow10[(ndims-idim-1)];
+          for ( int idim = iodim; idim < ndims; idim++ )
+            {
+              if      ( xdimid == ncvars[ncvarid].dimids[idim] )
+                ixyz += 1*ipow10[ndims-idim-1];
+              else if ( ydimid == ncvars[ncvarid].dimids[idim] )
+                ixyz += 2*ipow10[ndims-idim-1];
+              else if ( zdimid == ncvars[ncvarid].dimids[idim] )
+                ixyz += 3*ipow10[ndims-idim-1];
+            }
         }
 
       vlistDefVarXYZ(vlistID, varID, ixyz);
@@ -6258,7 +6427,7 @@ int cdfInqContents(int streamID)
 	  continue;
 	}
 
-      if ( cdfInqDatatype(ncvars[ncvarid].xtype) == -1 )
+      if ( cdfInqDatatype(ncvars[ncvarid].xtype, ncvars[ncvarid].lunsigned) == -1 )
 	{
 	  ncvars[ncvarid].isvar = 0;
 	  Warning("Variable %s has an unsupported data type, skipped!", ncvars[ncvarid].name);
diff --git a/src/stream_cgribex.c b/src/stream_cgribex.c
index efb1910f950f1d0c5fe2fce2100ca4755874fe5c..d36657b37af2f53755f6d258ebbf8727b8293b1f 100644
--- a/src/stream_cgribex.c
+++ b/src/stream_cgribex.c
@@ -77,6 +77,7 @@ int cgribexGetZaxisHasBounds(int grb_ltype)
 
   switch (grb_ltype)
     {
+    case GRIB1_LTYPE_SIGMA_LAYER:
     case GRIB1_LTYPE_HYBRID_LAYER:
     case GRIB1_LTYPE_LANDDEPTH_LAYER:
       {
@@ -1589,13 +1590,13 @@ void cgribexDefTime(int *isec1, int vdate, int vtime, int tsteptype, int numavg,
       int rdate, rtime;
       int ip1 = 0, ip2 = 0;
       int calendar;
-      
+
       calendar = taxisInqCalendar(taxisID);
       rdate    = taxisInqRdate(taxisID);
       rtime    = taxisInqRtime(taxisID);
 
       factor = cgribexDefDateTime(isec1, timeunit, rdate, rtime);
-      
+
       timerange = cgribexDefTimerange(tsteptype, factor, calendar,
 				      rdate, rtime, vdate, vtime, &ip1, &ip2);
 
@@ -1615,7 +1616,7 @@ void cgribexDefTime(int *isec1, int vdate, int vtime, int tsteptype, int numavg,
 	  if ( ip1 < 0 || ip1 > 0xFF   ) timetype = TAXIS_ABSOLUTE;
 	  if ( ip2 < 0 || ip2 > 0xFF   ) timetype = TAXIS_ABSOLUTE;
 	}
-      
+
       ISEC1_TimeRange   = timerange;
       ISEC1_TimePeriod1 = ip1;
       ISEC1_TimePeriod2 = ip2;
@@ -1921,6 +1922,27 @@ void cgribexDefLevel(int *isec1, int *isec2, double *fsec2, int zaxisID, int lev
 	ISEC1_Level2    = 0;
 	break;
       }
+    case ZAXIS_TOA:
+      {
+	ISEC1_LevelType = GRIB1_LTYPE_TOA;
+	ISEC1_Level1    = 0;
+	ISEC1_Level2    = 0;
+	break;
+      }
+    case ZAXIS_SEA_BOTTOM:
+      {
+	ISEC1_LevelType = GRIB1_LTYPE_SEA_BOTTOM;
+	ISEC1_Level1    = 0;
+	ISEC1_Level2    = 0;
+	break;
+      }
+    case ZAXIS_ATMOSPHERE:
+      {
+	ISEC1_LevelType = GRIB1_LTYPE_ATMOSPHERE;
+	ISEC1_Level1    = 0;
+	ISEC1_Level2    = 0;
+	break;
+      }
     case ZAXIS_MEANSEA:
       {
 	ISEC1_LevelType = GRIB1_LTYPE_MEANSEA;
@@ -2020,12 +2042,21 @@ void cgribexDefLevel(int *isec1, int *isec2, double *fsec2, int zaxisID, int lev
       }
     case ZAXIS_SIGMA:
       {
-	level = zaxisInqLevel(zaxisID, levelID);
+	if ( zaxisInqLbounds(zaxisID, NULL) && zaxisInqUbounds(zaxisID, NULL) )
+	  {
+	    ISEC1_LevelType = GRIB1_LTYPE_SIGMA_LAYER;
+	    ISEC1_Level1    = (int) zaxisInqLbound(zaxisID, levelID);
+	    ISEC1_Level2    = (int) zaxisInqUbound(zaxisID, levelID);
+	  }
+	else
+	  {
+            level = zaxisInqLevel(zaxisID, levelID);
 
-	ilevel = (int) level;
-	ISEC1_LevelType = GRIB1_LTYPE_SIGMA;
-	ISEC1_Level1    = ilevel;
-	ISEC1_Level2    = 0;
+            ilevel = (int) level;
+            ISEC1_LevelType = GRIB1_LTYPE_SIGMA;
+            ISEC1_Level1    = ilevel;
+            ISEC1_Level2    = 0;
+          }
 
 	break;
       }
diff --git a/src/stream_grb.c b/src/stream_grb.c
index cb55a7abd3d341d7ed468e7b2fc79b97a264c386..12a8fcd5bef8329e95ce0dc47a483f75c8e7ca34 100644
--- a/src/stream_grb.c
+++ b/src/stream_grb.c
@@ -33,12 +33,16 @@ int grib1ltypeToZaxisType(int grib_ltype)
   switch ( grib_ltype )
     {
     case GRIB1_LTYPE_SURFACE:         { zaxistype = ZAXIS_SURFACE;           break; }
+    case GRIB1_LTYPE_TOA:             { zaxistype = ZAXIS_TOA;               break; }
+    case GRIB1_LTYPE_SEA_BOTTOM:      { zaxistype = ZAXIS_SEA_BOTTOM;        break; }
+    case GRIB1_LTYPE_ATMOSPHERE:      { zaxistype = ZAXIS_ATMOSPHERE;        break; }
     case GRIB1_LTYPE_MEANSEA:         { zaxistype = ZAXIS_MEANSEA;           break; }
     case GRIB1_LTYPE_99:
     case GRIB1_LTYPE_ISOBARIC:        { zaxistype = ZAXIS_PRESSURE;          break; }
     case GRIB1_LTYPE_HEIGHT:          { zaxistype = ZAXIS_HEIGHT;            break; }
     case GRIB1_LTYPE_ALTITUDE:        { zaxistype = ZAXIS_ALTITUDE;	     break; }
-    case GRIB1_LTYPE_SIGMA:           { zaxistype = ZAXIS_SIGMA;	     break; }
+    case GRIB1_LTYPE_SIGMA:
+    case GRIB1_LTYPE_SIGMA_LAYER:     { zaxistype = ZAXIS_SIGMA;	     break; }
     case GRIB1_LTYPE_HYBRID:
     case GRIB1_LTYPE_HYBRID_LAYER:    { zaxistype = ZAXIS_HYBRID;	     break; }
     case GRIB1_LTYPE_LANDDEPTH:
@@ -58,6 +62,9 @@ int grib2ltypeToZaxisType(int grib_ltype)
   switch ( grib_ltype )
     {
     case GRIB2_LTYPE_SURFACE:            { zaxistype = ZAXIS_SURFACE;           break; }
+    case GRIB2_LTYPE_TOA:                { zaxistype = ZAXIS_TOA;               break; }
+    case GRIB2_LTYPE_SEA_BOTTOM:         { zaxistype = ZAXIS_SEA_BOTTOM;        break; }
+    case GRIB2_LTYPE_ATMOSPHERE:         { zaxistype = ZAXIS_ATMOSPHERE;        break; }
     case GRIB2_LTYPE_MEANSEA:            { zaxistype = ZAXIS_MEANSEA;           break; }
     case GRIB2_LTYPE_ISOBARIC:           { zaxistype = ZAXIS_PRESSURE;          break; }
     case GRIB2_LTYPE_HEIGHT:             { zaxistype = ZAXIS_HEIGHT;            break; }
diff --git a/src/stream_gribapi.c b/src/stream_gribapi.c
index 423b29056f58fa150a04ee92f6d718d5950cb47e..cca697bfc17e4a5f06386994173171686b721511 100644
--- a/src/stream_gribapi.c
+++ b/src/stream_gribapi.c
@@ -293,7 +293,7 @@ void gribapiGetGrid(grib_handle *gh, grid_t *grid)
 	  {
 	    if ( grid->xsize > 1 )
 	      {
-		if ( (grid->xfirst > grid->xlast) && (grid->xfirst >= 180) ) grid->xfirst -= 360;
+		if ( (grid->xfirst >= grid->xlast) && (grid->xfirst >= 180) ) grid->xfirst -= 360;
 
 		if ( editionNumber <= 1 )
 		  {
@@ -301,7 +301,7 @@ void gribapiGetGrid(grib_handle *gh, grid_t *grid)
 		    if ( IS_EQUAL(grid->xfirst, 0) && grid->xlast > 354 )
 		      {
 			double xinc = 360. / grid->xsize;
-		      
+
 			if ( fabs(grid->xinc-xinc) > 0.0 )
 			  {
 			    grid->xinc = xinc;
@@ -310,9 +310,9 @@ void gribapiGetGrid(grib_handle *gh, grid_t *grid)
 		      }
 		  }
 	      }
-	    grid->xdef   = 2;	    
+	    grid->xdef = 2;
 	  }
-	grid->ydef  = 0;
+	grid->ydef = 0;
         /* if ( IS_NOT_EQUAL(grid->yfirst, 0) || IS_NOT_EQUAL(grid->ylast, 0) ) */
 	  {
 	    if ( grid->ysize > 1 )
@@ -321,7 +321,7 @@ void gribapiGetGrid(grib_handle *gh, grid_t *grid)
 		  {
 		  }
 	      }
-	    grid->ydef   = 2;	    
+	    grid->ydef = 2;
 	  }
 	break;
       }
@@ -1973,7 +1973,7 @@ void gribapiDefTime(grib_handle *gh , int vdate, int vtime, int tsteptype, int n
       rtime    = taxisInqRtime(taxisID);
 
       factor = gribapiDefDateTime(gh, timeunit, rdate, rtime);
-      
+
       timerange = gribapiDefTimerange(tsteptype, factor, calendar,
 				      rdate, rtime, vdate, vtime, &ip);
       // printf("timerange: %d %d\n", timerange, ip);
@@ -2005,7 +2005,7 @@ void gribapiDefGrid(grib_handle *gh, int gridID, int ljpeg)
   static short lwarn = TRUE;
   size_t len;
   char *mesg;
-  
+
   gridtype = gridInqType(gridID);
 
   GRIB_CHECK(grib_get_long(gh, "editionNumber", &editionNumber), 0);
@@ -2024,7 +2024,7 @@ void gribapiDefGrid(grib_handle *gh, int gridID, int ljpeg)
       if ( (ysize ==  32 || ysize ==  48 || ysize ==  64 ||
 	    ysize ==  96 || ysize == 160 || ysize == 192 ||
 	    ysize == 240 || ysize == 320 || ysize == 384 ||
-	    ysize == 480 || ysize == 768 ) && 
+	    ysize == 480 || ysize == 768 ) &&
 	   (xsize == 2*ysize || xsize == 1) )
 	{
 	  gridtype = GRID_GAUSSIAN;
@@ -2136,6 +2136,12 @@ void gribapiDefGrid(grib_handle *gh, int gridID, int ljpeg)
 	GRIB_CHECK(grib_set_double(gh, "latitudeOfFirstGridPointInDegrees",  yfirst), 0);
 	GRIB_CHECK(grib_set_double(gh, "latitudeOfLastGridPointInDegrees",   ylast), 0);
 	GRIB_CHECK(grib_set_double(gh, "iDirectionIncrementInDegrees", xinc), 0);
+
+        {
+          long jscan = 0;
+          if ( yfirst < ylast ) jscan = 1;
+          GRIB_CHECK(grib_set_long(gh, "jScansPositively", jscan), 0);
+        }
 	/*
 	if ( fabs(xinc*1000 - ISEC2_LonIncr) > FLT_EPSILON )
 	  ISEC2_LonIncr = 0;
@@ -2338,6 +2344,30 @@ void gribapiDefLevel(grib_handle *gh, int param, int zaxisID, int levelID)
 	GRIB_CHECK(grib_set_long(gh, "level", (long) zaxisInqLevel(zaxisID, levelID)), 0);
 	break;
       }
+    case ZAXIS_TOA:
+      {
+	if ( editionNumber <= 1 )
+	  GRIB_CHECK(grib_set_long(gh, "indicatorOfTypeOfLevel", GRIB1_LTYPE_TOA), 0);
+	else
+	  GRIB_CHECK(grib_set_long(gh, "typeOfFirstFixedSurface", GRIB2_LTYPE_TOA), 0);
+	break;
+      }
+    case ZAXIS_SEA_BOTTOM:
+      {
+	if ( editionNumber <= 1 )
+	  GRIB_CHECK(grib_set_long(gh, "indicatorOfTypeOfLevel", GRIB1_LTYPE_SEA_BOTTOM), 0);
+	else
+	  GRIB_CHECK(grib_set_long(gh, "typeOfFirstFixedSurface", GRIB2_LTYPE_SEA_BOTTOM), 0);
+	break;
+      }
+    case ZAXIS_ATMOSPHERE:
+      {
+	if ( editionNumber <= 1 )
+	  GRIB_CHECK(grib_set_long(gh, "indicatorOfTypeOfLevel", GRIB1_LTYPE_ATMOSPHERE), 0);
+	else
+	  GRIB_CHECK(grib_set_long(gh, "typeOfFirstFixedSurface", GRIB2_LTYPE_ATMOSPHERE), 0);
+	break;
+      }
     case ZAXIS_MEANSEA:
       {
 	level = zaxisInqLevel(zaxisID, levelID);
diff --git a/src/varscan.c b/src/varscan.c
index e95d61b5a57bb72cb034968657dcc60a89330864..a3b610cdc3a70e4f3bba435243959ff21b2b7718 100644
--- a/src/varscan.c
+++ b/src/varscan.c
@@ -279,7 +279,6 @@ void varAddRecord(int recID, int param, int gridID, int zaxistype, int lbounds,
     {
       nvars++;
       varID = paramNewEntry(param);
-      if ( prec > vartable[varID].prec ) vartable[varID].prec = prec;
       vartable[varID].gridID    = gridID;
       vartable[varID].zaxistype = zaxistype;
       vartable[varID].ltype     = ltype;
@@ -296,19 +295,21 @@ void varAddRecord(int recID, int param, int gridID, int zaxistype, int lbounds,
       if ( vartable[varID].gridID != gridID )
 	{
 	  char paramstr[32];
-	  cdiParamToString(param, paramstr, sizeof(paramstr));	  
+	  cdiParamToString(param, paramstr, sizeof(paramstr));
 	  Message("param = %s gridID = %d", paramstr, gridID);
 	  Error("horizontal grid must not change for same param!");
 	}
       if ( vartable[varID].zaxistype != zaxistype )
 	{
 	  char paramstr[32];
-	  cdiParamToString(param, paramstr, sizeof(paramstr));	  
+	  cdiParamToString(param, paramstr, sizeof(paramstr));
 	  Message("param = %s zaxistype = %d", paramstr, zaxistype);
 	  Error("zaxistype must not change for same param!");
 	}
     }
 
+  if ( prec > vartable[varID].prec ) vartable[varID].prec = prec;
+
   levelID = levelNewEntry(varID, level1, level2);
   vartable[varID].levelTable[levelID].recID = recID;
 
@@ -713,11 +714,11 @@ int zaxisCompare(int zaxisID, int zaxistype, int nlevels, int lbounds, double *l
 	  const double *dlevels;
 	  char zlongname[256];
 	  char zunits[256];
-	  
+
 	  dlevels = zaxisInqLevelsPtr(zaxisID);
 	  for ( levelID = 0; levelID < nlevels; levelID++ )
 	    {
-	      if ( fabs(dlevels[levelID] - levels[levelID]) > 0 )
+	      if ( fabs(dlevels[levelID] - levels[levelID]) > 1.e-9 )
 		break;
 	    }
 
@@ -743,8 +744,8 @@ int zaxisCompare(int zaxisID, int zaxistype, int nlevels, int lbounds, double *l
 }
 
 
-int varDefZaxis(int vlistID, int zaxistype, int nlevels, double *levels, int lbounds, 
-		double *levels1, double *levels2, int vctsize, double *vct, char *name, 
+int varDefZaxis(int vlistID, int zaxistype, int nlevels, double *levels, int lbounds,
+		double *levels1, double *levels2, int vctsize, double *vct, char *name,
 		char *longname, char *units, int prec, int mode, int ltype)
 {
   /*
diff --git a/src/vlist.h b/src/vlist.h
index dfb4ae6b64fd630a36ff406ffc68480b05e21431..1607ddfbd858888a1bd0cf36940b31ae06c8c097 100644
--- a/src/vlist.h
+++ b/src/vlist.h
@@ -15,6 +15,7 @@
 #  include "cdi_limits.h"
 #endif
 
+#define VALIDMISS 1.e+303
 
 /*
  * CDI attribute
@@ -80,6 +81,7 @@ typedef struct
   int         timaccu;
   int         xyz;
   int         missvalused; /* TRUE if missval is defined */
+  int         lvalidrange;
   char       *name;
   char       *longname;
   char       *stdname;
@@ -87,6 +89,7 @@ typedef struct
   double      missval;
   double      scalefactor;
   double      addoffset;
+  double      validrange[2];
   levinfo_t  *levinfo;
   int         comptype;     // compression type
   int         complevel;    // compression level
@@ -140,6 +143,12 @@ void     vlistDefVarDeco ( int vlistID, int varID, int decoSize,
                               deco_t * deco );
 #endif
 
+/*      vlistDefVarValidrange: Define the valid range of a Variable */
+void    vlistDefVarValidrange(int vlistID, int varID, const double *validrange);
+
+/*      vlistInqVarValidrange: Get the valid range of a Variable */
+int     vlistInqVarValidrange(int vlistID, int varID, double *validrange);
+
 #endif  /* _VLIST_H */
 /*
  * Local Variables:
diff --git a/src/vlist_var.c b/src/vlist_var.c
index cfb558d3a45c1ddbecef04a2add4bf5b2f6c18c5..a011e9a1db2dea8de1a639c0df4b7cbf5c5562c1 100644
--- a/src/vlist_var.c
+++ b/src/vlist_var.c
@@ -27,39 +27,41 @@ void vlistvarInitEntry(int vlistID, int varID)
 
   vlistptr = vlist_to_pointer(vlistID);
 
-  vlistptr->vars[varID].fvarID       = varID;
-  vlistptr->vars[varID].mvarID       = varID;
-  vlistptr->vars[varID].flag         = 0;
-  vlistptr->vars[varID].param        = 0;
-  vlistptr->vars[varID].timeID       = CDI_UNDEFID;
-  vlistptr->vars[varID].datatype     = CDI_UNDEFID;
-  vlistptr->vars[varID].tsteptype    = TSTEP_INSTANT;
-  vlistptr->vars[varID].timave       = 0;
-  vlistptr->vars[varID].timaccu      = 0;
-  vlistptr->vars[varID].xyz          = 0;
-  vlistptr->vars[varID].gridID       = CDI_UNDEFID;
-  vlistptr->vars[varID].zaxisID      = CDI_UNDEFID;
-  vlistptr->vars[varID].instID       = CDI_UNDEFID;
-  vlistptr->vars[varID].modelID      = CDI_UNDEFID;
-  vlistptr->vars[varID].tableID      = CDI_UNDEFID;
-  vlistptr->vars[varID].missvalused  = FALSE;
-  vlistptr->vars[varID].missval      = cdiDefaultMissval;
-  vlistptr->vars[varID].addoffset    = 0.0;
-  vlistptr->vars[varID].scalefactor  = 1.0;
-  vlistptr->vars[varID].name         = NULL;
-  vlistptr->vars[varID].longname     = NULL;
-  vlistptr->vars[varID].stdname      = NULL;
-  vlistptr->vars[varID].units        = NULL;
-  vlistptr->vars[varID].nlevs        = 0;
-  vlistptr->vars[varID].levinfo      = NULL;
-  vlistptr->vars[varID].comptype     = COMPRESS_NONE;
-  vlistptr->vars[varID].complevel    = 1;
-  vlistptr->vars[varID].atts.nalloc  = MAX_ATTRIBUTES;
-  vlistptr->vars[varID].atts.nelems  = 0;
+  vlistptr->vars[varID].fvarID        = varID;
+  vlistptr->vars[varID].mvarID        = varID;
+  vlistptr->vars[varID].flag          = 0;
+  vlistptr->vars[varID].param         = 0;
+  vlistptr->vars[varID].timeID        = CDI_UNDEFID;
+  vlistptr->vars[varID].datatype      = CDI_UNDEFID;
+  vlistptr->vars[varID].tsteptype     = TSTEP_INSTANT;
+  vlistptr->vars[varID].timave        = 0;
+  vlistptr->vars[varID].timaccu       = 0;
+  vlistptr->vars[varID].xyz           = 0;
+  vlistptr->vars[varID].gridID        = CDI_UNDEFID;
+  vlistptr->vars[varID].zaxisID       = CDI_UNDEFID;
+  vlistptr->vars[varID].instID        = CDI_UNDEFID;
+  vlistptr->vars[varID].modelID       = CDI_UNDEFID;
+  vlistptr->vars[varID].tableID       = CDI_UNDEFID;
+  vlistptr->vars[varID].missvalused   = FALSE;
+  vlistptr->vars[varID].missval       = cdiDefaultMissval;
+  vlistptr->vars[varID].addoffset     = 0.0;
+  vlistptr->vars[varID].scalefactor   = 1.0;
+  vlistptr->vars[varID].name          = NULL;
+  vlistptr->vars[varID].longname      = NULL;
+  vlistptr->vars[varID].stdname       = NULL;
+  vlistptr->vars[varID].units         = NULL;
+  vlistptr->vars[varID].nlevs         = 0;
+  vlistptr->vars[varID].levinfo       = NULL;
+  vlistptr->vars[varID].comptype      = COMPRESS_NONE;
+  vlistptr->vars[varID].complevel     = 1;
+  vlistptr->vars[varID].atts.nalloc   = MAX_ATTRIBUTES;
+  vlistptr->vars[varID].atts.nelems   = 0;
+  vlistptr->vars[varID].lvalidrange   = 0;
+  vlistptr->vars[varID].validrange[0] = VALIDMISS;
+  vlistptr->vars[varID].validrange[1] = VALIDMISS;
   vlistptr->vars[varID].iorank       = CDI_UNDEFID;
   vlistptr->vars[varID].decoSize     = 0;
   vlistptr->vars[varID].deco         = NULL;
-
 }
 
 static
@@ -1152,8 +1154,22 @@ double vlistInqVarMissval(int vlistID, int varID)
   return (vlistptr->vars[varID].missval);
 }
 
+/*
+@Function  vlistDefVarMissval
+@Title     Define the missing value of a Variable
 
-double vlistInqVarScalefactor(int vlistID, int varID)
+@Prototype void vlistDefVarMissval(int vlistID, int varID, double missval)
+@Parameter
+    @Item  vlistID  Variable list ID, from a previous call to @fref{vlistCreate}.
+    @Item  varID    Variable identifier.
+    @Item  missval  Missing value.
+
+@Description
+The function @func{vlistDefVarMissval} defines the missing value of a variable.
+
+@EndFunction
+*/
+void vlistDefVarMissval(int vlistID, int varID, double missval)
 {
   vlist_t *vlistptr;
 
@@ -1161,11 +1177,12 @@ double vlistInqVarScalefactor(int vlistID, int varID)
 
   vlistCheckVarID(__func__, vlistID, varID);
 
-  return (vlistptr->vars[varID].scalefactor);
+  vlistptr->vars[varID].missval = missval;
+  vlistptr->vars[varID].missvalused = TRUE;
 }
 
 
-double vlistInqVarAddoffset(int vlistID, int varID)
+int vlistInqVarValidrange(int vlistID, int varID, double *validrange)
 {
   vlist_t *vlistptr;
 
@@ -1173,25 +1190,31 @@ double vlistInqVarAddoffset(int vlistID, int varID)
 
   vlistCheckVarID(__func__, vlistID, varID);
 
-  return (vlistptr->vars[varID].addoffset);
+  if ( validrange != NULL && vlistptr->vars[varID].lvalidrange )
+    {
+      validrange[0] = vlistptr->vars[varID].validrange[0];
+      validrange[1] = vlistptr->vars[varID].validrange[1];
+    }
+
+  return (vlistptr->vars[varID].lvalidrange);
 }
 
-/*
-@Function  vlistDefVarMissval
-@Title     Define the missing value of a Variable
 
-@Prototype void vlistDefVarMissval(int vlistID, int varID, double missval)
-@Parameter
-    @Item  vlistID  Variable list ID, from a previous call to @fref{vlistCreate}.
-    @Item  varID    Variable identifier.
-    @Item  missval  Missing value.
+void vlistDefVarValidrange(int vlistID, int varID, const double *validrange)
+{
+  vlist_t *vlistptr;
 
-@Description
-The function @func{vlistDefVarMissval} defines the missing value of a variable.
+  vlistptr = vlist_to_pointer(vlistID);
 
-@EndFunction
-*/
-void vlistDefVarMissval(int vlistID, int varID, double missval)
+  vlistCheckVarID(__func__, vlistID, varID);
+
+  vlistptr->vars[varID].validrange[0] = validrange[0];
+  vlistptr->vars[varID].validrange[1] = validrange[1];
+  vlistptr->vars[varID].lvalidrange = TRUE;
+}
+
+
+double vlistInqVarScalefactor(int vlistID, int varID)
 {
   vlist_t *vlistptr;
 
@@ -1205,11 +1228,21 @@ void vlistDefVarMissval(int vlistID, int varID, double missval)
 
   vlistCheckVarID(__func__, vlistID, varID);
 
-  vlistptr->vars[varID].missval = missval;
-  vlistptr->vars[varID].missvalused = TRUE;
+  return (vlistptr->vars[varID].scalefactor);
 }
 
 
+double vlistInqVarAddoffset(int vlistID, int varID)
+{
+  vlist_t *vlistptr;
+
+  vlistptr = vlist_to_pointer(vlistID);
+
+  vlistCheckVarID(__func__, vlistID, varID);
+
+  return (vlistptr->vars[varID].addoffset);
+}
+
 void vlistDefVarScalefactor(int vlistID, int varID, double scalefactor)
 {
   vlist_t *vlistptr;
diff --git a/src/zaxis.c b/src/zaxis.c
index 70e2458e1226264749d592bca9bff280f372db8a..c358abdd55a6ee97515862e8b3deab0ea8340e0b 100644
--- a/src/zaxis.c
+++ b/src/zaxis.c
@@ -28,19 +28,22 @@ static struct {
   char *units;    // 1: up;  2: down
 }
 ZaxistypeEntry[] = {
-  {0, "sfc",     "surface",           "",               ""},
-  {0, "lev",     "generic",           "",               "level"},
-  {2, "lev",     "hybrid",            "",               "level"},
-  {2, "lev",     "hybrid_half",       "",               "level"},
-  {2, "lev",     "pressure",          "air_pressure",   "Pa"},
-  {1, "height",  "height",            "height",         "m"},
-  {2, "depth",   "depth_below_sea",   "depth",          "m"},
-  {2, "depth",   "depth_below_land",  "",               "cm"},
-  {0, "lev",     "isentropic",        "",               "K"},
-  {0, "lev",     "trajectory",        "",               ""},
-  {1, "alt",     "altitude",          "",               "m"},
-  {0, "lev",     "sigma",             "",               "level"},
-  {0, "lev",     "meansea",           "",               "level"},
+  { /*  0 */ 0, "sfc",        "surface",           "",               ""},
+  { /*  1 */ 0, "lev",        "generic",           "",               "level"},
+  { /*  2 */ 2, "lev",        "hybrid",            "",               "level"},
+  { /*  3 */ 2, "lev",        "hybrid_half",       "",               "level"},
+  { /*  4 */ 2, "lev",        "pressure",          "air_pressure",   "Pa"},
+  { /*  5 */ 1, "height",     "height",            "height",         "m"},
+  { /*  6 */ 2, "depth",      "depth_below_sea",   "depth",          "m"},
+  { /*  7 */ 2, "depth",      "depth_below_land",  "",               "cm"},
+  { /*  8 */ 0, "lev",        "isentropic",        "",               "K"},
+  { /*  9 */ 0, "lev",        "trajectory",        "",               ""},
+  { /* 10 */ 1, "alt",        "altitude",          "",               "m"},
+  { /* 11 */ 0, "lev",        "sigma",             "",               "level"},
+  { /* 12 */ 0, "lev",        "meansea",           "",               "level"},
+  { /* 13 */ 0, "toa",        "top_of_atmosphere", "",               ""},
+  { /* 14 */ 0, "seabottom",  "sea_bottom",        "",               ""},
+  { /* 15 */ 0, "atmosphere", "atmosphere",        "",               ""},
 };
 
 
@@ -85,7 +88,6 @@ resOps zaxisOps = { zaxisCompareP, zaxisDestroyP, zaxisPrintP
 
 static int  ZAXIS_Debug = 0;   /* If set to 1, debugging */
 
-
 static
 void zaxisDefaultValue ( zaxis_t *zaxisptr )
 {
@@ -408,7 +410,7 @@ void zaxisDefUnits(int zaxisID, const char *units)
 @Prototype void zaxisInqName(int zaxisID, char *name)
 @Parameter
     @Item  zaxisID  Z-axis ID, from a previous call to @fref{zaxisCreate}.
-    @Item  name     Name of the Z-axis. The caller must allocate space for the 
+    @Item  name     Name of the Z-axis. The caller must allocate space for the
                     returned string. The maximum possible length, in characters, of
                     the string is given by the predefined constant CDI_MAX_NAME.
 
@@ -438,7 +440,7 @@ void zaxisInqName(int zaxisID, char *name)
 @Prototype void zaxisInqLongname(int zaxisID, char *longname)
 @Parameter
     @Item  zaxisID  Z-axis ID, from a previous call to @fref{zaxisCreate}.
-    @Item  longname Longname of the Z-axis. The caller must allocate space for the 
+    @Item  longname Longname of the Z-axis. The caller must allocate space for the
                     returned string. The maximum possible length, in characters, of
                     the string is given by the predefined constant CDI_MAX_NAME.
 
@@ -468,7 +470,7 @@ void zaxisInqLongname(int zaxisID, char *longname)
 @Prototype void zaxisInqUnits(int zaxisID, char *units)
 @Parameter
     @Item  zaxisID  Z-axis ID, from a previous call to @fref{zaxisCreate}
-    @Item  units    Units of the Z-axis. The caller must allocate space for the 
+    @Item  units    Units of the Z-axis. The caller must allocate space for the
                     returned string. The maximum possible length, in characters, of
                     the string is given by the predefined constant CDI_MAX_NAME.
 
@@ -733,7 +735,7 @@ const double *zaxisInqLevelsPtr(int zaxisID)
     @Item  zaxisID  Z-axis ID, from a previous call to @fref{zaxisCreate}.
     @Item  levels   Pointer to the location into which the levels are read.
                     The caller must allocate space for the returned values.
-    
+
 @Description
 The function @func{zaxisInqLevels} returns all levels of a Z-axis.
 
diff --git a/util/sunf95preproc-wrapper b/util/sunf95preproc-wrapper
new file mode 100755
index 0000000000000000000000000000000000000000..49a1d649f1f208cc531483251ed7601f296de1ab
--- /dev/null
+++ b/util/sunf95preproc-wrapper
@@ -0,0 +1,76 @@
+#! /bin/sh
+if [ "${DEBUG+set}" = set ]; then
+  set -x
+  outputRedir=">&2"
+else
+  outputRedir=">/dev/null 2>&1"
+fi
+while echo "$1" | grep '^-' >/dev/null 2>&1; do
+  FPPFLAGS="${FPPFLAGS+${FPPFLAGS} }\"$1\""
+  shift
+done
+TRAPCMD=':'
+trap 'eval $TRAPCMD' 0
+set -e
+test "${DEBUG+set}" = set && echo "$FPPFLAGS" >&2
+FCFLAGS=${FCFLAGS--EP}
+if [ "${FC+set}" != set ]; then
+  for F90C in sunf95 '' ; do
+    test -n "$F90C" || exit 1
+    set +e
+    F90BIN=`which $F90C 2>/dev/null`
+    set -e
+    test ! -x "$F90BIN" || break
+  done
+fi
+FC=${FC-$F90C}
+# append -fpp if necessary
+IFStr=`echo "$IFS" | sed -n '$!s/$/\\\\n/
+H
+$x
+$s/\n//g
+$P'`
+if echo "$FCFLAGS" \
+  | grep -v '\('"[$IFStr]\\|^\\)-[cf]pp\\([$IFStr]\\|\$\\)" >/dev/null
+then
+  FCFLAGS="${FCFLAGS+$FCFLAGS }-fpp"
+fi
+TMPDIR="${TMPDIR-/tmp}"
+{
+  tmp=`
+  (umask 077 && mktemp -d "$TMPDIR/fooXXXXXX") 2>/dev/null
+  ` &&
+  test -n "$tmp" && test -d "$tmp"
+} || {
+  tmp="$TMPDIR/foo$$-$RANDOM"
+  (umask 077 && mkdir "$tmp")
+} || exit $?
+TRAPCMD="$TRAPCMD ; rm -rf \"$tmp\""
+#echo \"$FCFLAGS\"
+for FortranFile in "$@" ; do
+  fppInput=`echo "$FortranFile" | sed 's:.*/::
+s:\.[^./]*:.F90:'`
+  cp "$FortranFile" "$tmp/$fppInput"
+  fppOutput=`echo "$fppInput" | sed 's:.*/::
+s:\.F90:.f90:'`
+  for i in fppOutput ; do
+    if eval test -e \$$i ; then
+      eval backup$i=true
+      eval mv \"\$$i\" \"\$$i.bak\"
+      eval TRAPCMD${i}Save=\"$TRAPCMD\"
+      TRAPCMD="$TRAPCMD ; "`eval echo mv \"\\$$i.bak\" \"\\$$i\"`
+    fi
+  done
+  set +e
+  eval \$FC \$FCFLAGS -F $FPPFLAGS \"\$tmp/\$fppInput\" $outputRedir
+  set -e
+  grep -v '^#' "$fppOutput"
+  test "${DEBUG+set}" = set && cat "$fppOutput" >&2
+  rm "$fppOutput" "$tmp/$fppInput"
+  for i in fppOutput ; do
+    if eval test \"\$backup$i\" = true ; then
+      eval mv \"\$$i.bak\" \"\$$i\"
+      TRAPCMD=`eval echo \\$TRAPCMD\${i}Save`
+    fi
+  done
+done
diff --git a/util/sxpreproc-wrapper b/util/sxpreproc-wrapper
new file mode 100755
index 0000000000000000000000000000000000000000..ceae5c9a43e6b181bf6a08c2bd3c1053f897b788
--- /dev/null
+++ b/util/sxpreproc-wrapper
@@ -0,0 +1,59 @@
+#! /bin/sh
+if [ "${DEBUG+set}" = set ]; then
+  set -x
+  outputRedir=">&2"
+else
+  outputRedir=">/dev/null 2>&1"
+fi
+while echo "$1" | grep '^-' >/dev/null 2>&1; do
+  FPPFLAGS="${FPPFLAGS+${FPPFLAGS} }$1"
+  shift
+done
+set -e
+test "${DEBUG+set}" = set && echo "$FPPFLAGS" >&2
+FCFLAGS=${FCFLAGS--EP}
+if [ "${FC+set}" != set ]; then
+  for F90C in sxf90 f90 '' ; do
+    test -n "$F90C" || exit 1
+    set +e
+    F90BIN=`which $F90C 2>/dev/null`
+    set -e
+    test ! -x "$F90BIN" || break
+  done
+fi
+FC=${FC-$F90C}
+IFStr=`echo "$IFS" | sed -n '$!s/$/\\\\n/
+H
+$x
+$s/\n//g
+$P'`
+#translate -Ep to -EP in FCFLAGS
+FCFLAGS=`echo "$FCFLAGS" | sed -e 's/\('"[$IFStr]\\|^\\)-Ep\\([$IFStr]\\|\$\\)"'/\1-EP\2/'`
+# append -EP if necessary
+if echo "$FCFLAGS" \
+  | grep -v '\('"[$IFStr]\\|^\\)-EP\\([$IFStr]\\|\$\\)" >/dev/null
+then
+  FCFLAGS="${FCFLAGS+$FCFLAGS }-EP"
+fi
+TMPDIR="${TMPDIR-/tmp}"
+{
+  tmp=`
+  (umask 077 && mktemp -d "$TMPDIR/fooXXXXXX") 2>/dev/null
+  ` &&
+  test -n "$tmp" && test -d "$tmp"
+} || {
+  tmp=$TMPDIR/foo$$-$RANDOM
+  (umask 077 && mkdir "$tmp")
+} || exit $?
+#echo \"$FCFLAGS\"
+for FortranFile in "$@" ; do
+  fppOutput=`echo "$FortranFile" | sed 's:.*/::
+s:^:'"$tmp/"'i.:'`
+  set +e
+  eval \$FC \$FCFLAGS -ts \"'$tmp'\" \$FPPFLAGS \"\$FortranFile\" $outputRedir
+  set -e
+  cat "$fppOutput" 2>/dev/null
+  test "${DEBUG+set}" = set && cat "$fppOutput" >&2
+  rm "$fppOutput"
+done
+rm -rf "$tmp"
\ No newline at end of file
diff --git a/util/xlfpreproc-wrapper b/util/xlfpreproc-wrapper
new file mode 100755
index 0000000000000000000000000000000000000000..647c2341bf60480a4a9a1548b6aed46380c29cf6
--- /dev/null
+++ b/util/xlfpreproc-wrapper
@@ -0,0 +1,22 @@
+#! /bin/sh
+#
+if [ "${DEBUG+set}" = set ]; then
+  set -x
+  outputRedir=">&2"
+else
+  outputRedir=">/dev/null 2>&1"
+fi
+while echo "$1" | grep '^-' >/dev/null 2>&1; do
+  FPPFLAGS="${FPPFLAGS+${FPPFLAGS} }$1"
+  shift
+done
+set -e
+for srcfile in "$@" ; do
+  set +e
+  eval \$FC \$FCFLAGS \$FPPFLAGS -d -qnoobject \"\$srcfile\" $outputRedir
+  set -e
+  FPPOUTNAME=`echo $srcfile | sed -e 's:\(.*/\)*\([^/]*\)\.[^./]*$:F\2.f:'`
+  cat "$FPPOUTNAME" 2>/dev/null
+  test "${DEBUG+set}" = set && cat "$FPPOUTNAME" >&2
+  rm "$FPPOUTNAME"
+done