New operator: intlevel - Linear level interpolation

.gitattributes

@@ -236,6 +236,7 @@ src/Info.c -text
 src/Input.c -text
 src/Intgrid.c -text
 src/Intgridtraj.c -text
+src/Intlevel.c -text
 src/Intntime.c -text
 src/Inttime.c -text
 src/Intyear.c -text

ChangeLog

@@ -3,7 +3,8 @@
 	* using CDI library version 1.2.0
 	* add support for netCDF4 classic with deflate option (option -z zip)
 	* Arith: add filltype FILL_FILE [request: Rita Seiffert]
-	* New operator: invertlev
+	* New operator: intlevel - Linear level interpolation
+	* New operator: invertlev - Invert level [request: Marco Giorgetta]
 	* New module: Tests (normal, studentt, chisquare, beta, fisher)
 	* New module: Echam5ini (read_e5ml, write_e5ml)
 	* Version 1.2.0 released
@@ -324,7 +325,7 @@
 2006-02-01  Uwe Schulzweida  <schulzweida@dkrz.de>
 
 	* using CDI library version 0.9.5
-	* use marco DBL_IS_EQUAL to check floating-point equality
+	* use macro DBL_IS_EQUAL to check floating-point equality
 	* use taxisCopyTimestep if posible
 	* New operator: sp2gp2 (Spectral) for ERA40 data
 	* New operator: gp2sp2 (Spectral) for ERA40 data

doc/tex/Modules

@@ -79,6 +79,7 @@ Remap         Interpolation
 Intgrid       Interpolation
 Remapeta      Interpolation
 Vertint       Interpolation
+Intlevel      Interpolation
 Inttime       Interpolation
 Intyear       Interpolation
 Spectral      Transformation

doc/tex/mod/Vertint

@@ -12,7 +12,7 @@ Interpolate 3D variables on hybrid model level to pressure or height level.
 The input file must contain the log. surface pressure or the surface pressure.
 To interpolate the temperature, the orography is also needed.
 The pressure, temperature, and orography are identified by there code numbers.
-Supported code tables are: WMO standard table number 2 and ECMWF local table number 128.
+Supported parameter tables are: WMO standard table number 2 and ECMWF local table number 128.
 Use the alias  ml2plx/ml2hlx or the environment variable EXTRAPOLATE
 to extrapolate missing values.
 @EndDescription

src/Intlevel.c

+/*
+  This file is part of CDO. CDO is a collection of Operators to
+  manipulate and analyse Climate model Data.
+
+  Copyright (C) 2003-2008 Uwe Schulzweida, Uwe.Schulzweida@zmaw.de
+  See COPYING file for copying and redistribution conditions.
+
+  This program is free software; you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation; version 2 of the License.
+
+  This program is distributed in the hope that it will be useful,
+  but WITHOUT ANY WARRANTY; without even the implied warranty of
+  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+  GNU General Public License for more details.
+*/
+
+/*
+   This module contains the following operators:
+
+      Intlevel   intlevel        Linear level interpolation
+*/
+
+
+#include <ctype.h>
+#include <string.h>
+#include <math.h>
+
+#include "cdi.h"
+#include "cdo.h"
+#include "cdo_int.h"
+#include "pstream.h"
+#include "list.h"
+
+
+static void interp_lev(int gridsize, double missval, double *vardata1, double *vardata2,
+		       int nlev2, int *lev_idx1, int *lev_idx2, double *lev_wgt1, double *lev_wgt2)
+{
+  int i, ilev;
+  int idx1, idx2;
+  double wgt1, wgt2;
+  double w1, w2;
+  double *var1L1, *var1L2, *var2;
+
+  for ( ilev = 0; ilev < nlev2; ++ilev )
+    {
+      idx1 = lev_idx1[ilev];
+      idx2 = lev_idx2[ilev];
+      wgt1 = lev_wgt1[ilev];
+      wgt2 = lev_wgt2[ilev];
+      var2 = vardata2+gridsize*ilev;
+
+      if ( cdoVerbose ) cdoPrint("level %d: idx1=%d idx2=%d wgt1=%g wgt2=%g", ilev, idx1, idx2, wgt1, wgt2);
+
+      var1L1 = vardata1+gridsize*idx1;
+      var1L2 = vardata1+gridsize*idx2;
+      for ( i = 0; i < gridsize; ++i )
+	{
+	  w1 = wgt1;
+	  w2 = wgt2;
+	  if ( DBL_IS_EQUAL(var1L1[i], missval) ) w1 = 0;
+	  if ( DBL_IS_EQUAL(var1L2[i], missval) ) w2 = 0;
+
+	  if ( DBL_IS_EQUAL(w1, 0) && DBL_IS_EQUAL(w2, 0) )
+	    {
+	      var2[i] = missval;
+	    }
+	  else if ( DBL_IS_EQUAL(w1, 0) )
+	    {
+	      if ( w2 >= 0.5 )
+		var2[i] = var1L2[i];
+	      else
+		var2[i] = missval;	      
+	    }
+	  else if ( DBL_IS_EQUAL(w2, 0) )
+	    {
+	      if ( w1 >= 0.5 )
+		var2[i] = var1L1[i];
+	      else
+		var2[i] = missval;	      
+	    }
+	  else
+	    {
+	      var2[i] = var1L1[i]*w1 + var1L2[i]*w2;
+	    }
+	}
+    }
+}
+
+
+static void gen_weights(int expol, int nlev1, double *lev1, int nlev2, double *lev2,
+			int *lev_idx1, int *lev_idx2, double *lev_wgt1, double *lev_wgt2)
+{
+  int i1, i2;
+  double val1, val2;
+  int    idx1 = 0, idx2 = 0;
+
+  for ( i2 = 0; i2 < nlev2; ++i2 )
+    {
+      for ( i1 = 1; i1 < nlev1; ++i1 )
+	{
+	  if ( lev1[i1-1] < lev1[i1] )
+	    {
+	      idx1 = i1 - 1;
+	      idx2 = i1;
+	    }
+	  else
+	    {
+	      idx1 = i1;
+	      idx2 = i1-1;
+	    }
+	  val1 = lev1[idx1];
+	  val2 = lev1[idx2];
+
+	  if ( lev2[i2] > val1 && lev2[i2] <= val2 ) break;
+	}
+
+      if ( i1 == nlev1 ) cdoAbort("Level %g not found!", lev2[i2]);
+      else
+	{
+	  if ( i1-1 == 0 )
+	    {
+	      lev_idx1[i2] = 1;
+	      lev_idx2[i2] = 1;
+	      lev_wgt1[i2] = 0;
+	      if ( expol || DBL_IS_EQUAL(lev2[i2], val2) )
+		lev_wgt2[i2] = 1;
+	      else
+		lev_wgt2[i2] = 0;
+	    }
+	  else if ( i1 == nlev1-1 )
+	    {
+	      lev_idx1[i2] = nlev1-2;
+	      lev_idx2[i2] = nlev1-2;
+	      if ( expol || DBL_IS_EQUAL(lev2[i2], val2) )
+		lev_wgt1[i2] = 1;
+	      else
+		lev_wgt1[i2] = 0;
+	      lev_wgt2[i2] = 0;
+	    }
+	  else
+	    {
+	      lev_idx1[i2] = idx1;
+	      lev_idx2[i2] = idx2;
+	      lev_wgt1[i2] = (lev1[lev_idx2[i2]] - lev2[i2]) / (lev1[lev_idx2[i2]] - lev1[lev_idx1[i2]]);
+	      lev_wgt2[i2] = (lev2[i2] - lev1[lev_idx1[i2]]) / (lev1[lev_idx2[i2]] - lev1[lev_idx1[i2]]);
+	    }
+	  lev_idx1[i2]--;
+	  lev_idx2[i2]--;
+	  /*
+	  printf("%d %g %d %d %g %g %d %d %g %g\n",
+	  i2, lev2[i2], idx1, idx2, lev1[idx1], lev1[idx2], 
+	  lev_idx1[i2], lev_idx2[i2], lev_wgt1[i2], lev_wgt2[i2]);
+	  */
+	}
+    }
+}
+
+
+void *Intlevel(void *argument)
+{
+  static char func[] = "Intlevel";
+  int INTLEVEL, INTLEVELX;
+  int operatorID;
+  int streamID1, streamID2;
+  int vlistID1, vlistID2;
+  int gridsize;
+  int recID, nrecs;
+  int i, offset;
+  int tsID, varID, levelID;
+  int nvars;
+  int nzaxis;
+  int zaxisID1 = -1, zaxisID2;
+  int gridID, zaxisID;
+  int nlev1, nlev2, nlevel = 0, maxlev;
+  int lup, ldown;
+  int **varnmiss = NULL;
+  int *varinterp = NULL;
+  int *vars = NULL;
+  int expol = FALSE;
+  double missval;
+  double *lev1 = NULL, *lev2 = NULL;
+  double *single1, *single2;
+  double **vardata1 = NULL, **vardata2 = NULL;
+  int taxisID1, taxisID2;
+  int *lev_idx1, *lev_idx2;
+  double *lev_wgt1, *lev_wgt2;
+  LIST *flist = listNew(FLT_LIST);
+
+  cdoInitialize(argument);
+
+  INTLEVEL  = cdoOperatorAdd("intlevel",  0, 0, NULL);
+  INTLEVELX = cdoOperatorAdd("intlevelx", 0, 0, NULL);
+
+  operatorID = cdoOperatorID();
+
+  if ( operatorID == INTLEVELX ) expol = TRUE;
+
+  operatorInputArg("levels");
+
+  nlev2 = args2fltlist(operatorArgc(), operatorArgv(), flist);
+  lev2  = (double *) listArrayPtr(flist);
+
+  if ( cdoVerbose ) for ( i = 0; i < nlev2; ++i ) printf("lev2 %d: %g\n", i, lev2[i]);
+
+  streamID1 = streamOpenRead(cdoStreamName(0));
+  if ( streamID1 < 0 ) cdiError(streamID1, "Open failed on %s", cdoStreamName(0));
+
+  vlistID1 = streamInqVlist(streamID1);
+  vlistID2 = vlistDuplicate(vlistID1);
+
+  taxisID1 = vlistInqTaxis(vlistID1);
+  taxisID2 = taxisDuplicate(taxisID1);
+  vlistDefTaxis(vlistID2, taxisID2);
+
+  nzaxis  = vlistNzaxis(vlistID1);
+  for ( i = 0; i < nzaxis; i++ )
+    {
+      zaxisID = vlistZaxis(vlistID1, i);
+      nlevel  = zaxisInqSize(zaxisID);
+      if ( zaxisInqType(zaxisID) != ZAXIS_HYBRID && zaxisInqType(zaxisID) != ZAXIS_HYBRID_HALF )
+	if ( nlevel > 1 )
+	  {
+	    zaxisID1 = zaxisID;
+	    break;
+	  }
+    }
+
+  if ( i == nzaxis ) cdoAbort("No processable variable found!");
+
+  nlev1 = nlevel;
+  lev1  = (double *) malloc((nlev1+2)*sizeof(double));
+  zaxisInqLevels(zaxisID1, lev1+1);
+
+  lup = FALSE;
+  ldown = FALSE;
+  for ( i = 1; i < nlev1; ++i )
+    {
+      if ( i == 1 )
+	{
+	  if ( lev1[i+1] > lev1[i] )
+	    lup = TRUE;
+	  else if ( lev1[i+1] < lev1[i] )
+	    ldown = TRUE;	
+	}
+      else
+	{
+	  if ( lup )
+	    {
+	      if ( !(lev1[i+1] > lev1[i]) ) lup = FALSE;
+	    }
+	  else if ( ldown )
+	    {
+	      if ( !(lev1[i+1] < lev1[i]) ) ldown = FALSE;
+	    }
+	}
+    }
+
+  if ( lup )
+    {
+      lev1[0]       = -1.e33; 
+      lev1[nlev1+1] =  1.e33;
+    }
+  else if ( ldown )
+    {
+      lev1[0]       =  1.e33; 
+      lev1[nlev1+1] = -1.e33;
+    }
+  else
+    cdoWarning("Non monotonic zaxis!");
+
+  if ( cdoVerbose ) for ( i = 0; i < nlev1+2; ++i ) printf("lev1 %d: %g\n", i, lev1[i]);
+
+  lev_idx1 = (int *) malloc(nlev2*sizeof(int));
+  lev_idx2 = (int *) malloc(nlev2*sizeof(int));
+  lev_wgt1 = (double *) malloc(nlev2*sizeof(double));
+  lev_wgt2 = (double *) malloc(nlev2*sizeof(double));
+
+  gen_weights(expol, nlev1+2, lev1, nlev2, lev2, lev_idx1, lev_idx2, lev_wgt1, lev_wgt2);
+
+  zaxisID2 = zaxisCreate(zaxisInqType(zaxisID1), nlev2);
+  zaxisDefLevels(zaxisID2, lev2);
+
+  for ( i = 0; i < nzaxis; i++ )
+    {
+      if ( zaxisID1 == vlistZaxis(vlistID1, i) )
+	{
+	  vlistChangeZaxisIndex(vlistID2, i, zaxisID2);
+	}
+    }
+
+  streamID2 = streamOpenWrite(cdoStreamName(1), cdoFiletype());
+  if ( streamID2 < 0 ) cdiError(streamID2, "Open failed on %s", cdoStreamName(1));
+
+  streamDefVlist(streamID2, vlistID2);
+
+  nvars = vlistNvars(vlistID1);
+
+  vars      = (int *) malloc(nvars*sizeof(int));
+  vardata1  = (double **) malloc(nvars*sizeof(double*));
+  vardata2  = (double **) malloc(nvars*sizeof(double*));
+  varnmiss  = (int **) malloc(nvars*sizeof(int*));
+  varinterp = (int *) malloc(nvars*sizeof(int));
+
+  maxlev   = nlev1 > nlev2 ? nlev1 : nlev2;
+
+  for ( varID = 0; varID < nvars; varID++ )
+    {
+      gridID   = vlistInqVarGrid(vlistID1, varID);
+      zaxisID  = vlistInqVarZaxis(vlistID1, varID);
+      gridsize = gridInqSize(gridID);
+      nlevel   = zaxisInqSize(zaxisID);
+
+      vardata1[varID] = (double *) malloc(gridsize*nlevel*sizeof(double));
+
+      if ( zaxisID == zaxisID1 )
+	{
+	  varinterp[varID] = TRUE;
+	  vardata2[varID]  = (double *) malloc(gridsize*nlev2*sizeof(double));
+	  varnmiss[varID]  = (int *) malloc(maxlev*sizeof(int));
+	  memset(varnmiss[varID], 0, maxlev*sizeof(int));
+	}
+      else
+	{
+	  varinterp[varID] = FALSE;
+	  vardata2[varID]  = vardata1[varID];
+	  varnmiss[varID]  = (int *) malloc(nlevel*sizeof(int));
+	}
+    }
+
+
+  tsID = 0;
+  while ( (nrecs = streamInqTimestep(streamID1, tsID)) )
+    {
+      for ( varID = 0; varID < nvars; ++varID ) vars[varID] = FALSE;
+
+      taxisCopyTimestep(taxisID2, taxisID1);
+
+      streamDefTimestep(streamID2, tsID);
+
+      for ( recID = 0; recID < nrecs; recID++ )
+	{
+	  streamInqRecord(streamID1, &varID, &levelID);
+	  gridsize = gridInqSize(vlistInqVarGrid(vlistID1, varID));
+	  nlevel   = zaxisInqSize(vlistInqVarZaxis(vlistID1, varID));
+	  offset   = gridsize*levelID;
+	  single1  = vardata1[varID] + offset;
+	  
+	  streamReadRecord(streamID1, single1, &varnmiss[varID][levelID]);
+	  vars[varID] = TRUE;
+	}
+
+      for ( varID = 0; varID < nvars; varID++ )
+	{
+	  if ( vars[varID] && varinterp[varID] )
+	    {
+	      gridID   = vlistInqVarGrid(vlistID1, varID);
+	      missval  = vlistInqVarMissval(vlistID1, varID);
+	      gridsize = gridInqSize(gridID);
+
+	      interp_lev(gridsize, missval, vardata1[varID], vardata2[varID],
+			 nlev2, lev_idx1, lev_idx2, lev_wgt1, lev_wgt2);
+	    }
+	}
+
+      for ( varID = 0; varID < nvars; varID++ )
+	{
+	  if ( vars[varID] )
+	    {
+	      nlevel = zaxisInqSize(vlistInqVarZaxis(vlistID2, varID));
+	      for ( levelID = 0; levelID < nlevel; levelID++ )
+		{
+		  gridsize = gridInqSize(vlistInqVarGrid(vlistID2, varID));
+		  offset   = gridsize*levelID;
+		  single2  = vardata2[varID] + offset;
+		  streamDefRecord(streamID2, varID, levelID);
+		  streamWriteRecord(streamID2, single2, varnmiss[varID][levelID]);
+		}
+	    }
+	}
+
+      tsID++;
+    }
+
+  streamClose(streamID2);
+  streamClose(streamID1);
+
+  for ( varID = 0; varID < nvars; varID++ )
+    {
+      free(varnmiss[varID]);
+      free(vardata1[varID]);
+      if ( varinterp[varID] ) free(vardata2[varID]);
+    }
+
+  free(varinterp);
+  free(varnmiss);
+  free(vardata2);
+  free(vardata1);
+  free(vars);
+
+  free(lev_idx1);
+  free(lev_idx2);
+  free(lev_wgt1);
+  free(lev_wgt2);
+
+  if (lev1) free(lev1);
+
+  listDelete(flist);
+
+  cdoFinish();
+
+  return (0);
+}

src/Makefile.am

@@ -36,6 +36,7 @@ cdo_SOURCES = Arith.c		\
               Input.c		\
               Intgrid.c 	\
               Intgridtraj.c 	\
+              Intlevel.c 	\
               Inttime.c 	\
               Intntime.c 	\
               Intyear.c 	\

src/Makefile.in

@@ -58,54 +58,54 @@ am_cdo_OBJECTS = Arith.$(OBJEXT) Arithc.$(OBJEXT) Arithdays.$(OBJEXT) \
 	Fldrms.$(OBJEXT) Fldstat.$(OBJEXT) Gradsdes.$(OBJEXT) \
 	Gridcell.$(OBJEXT) Harmonic.$(OBJEXT) Histogram.$(OBJEXT) \
 	Info.$(OBJEXT) Input.$(OBJEXT) Intgrid.$(OBJEXT) \
-	Intgridtraj.$(OBJEXT) Inttime.$(OBJEXT) Intntime.$(OBJEXT) \
-	Intyear.$(OBJEXT) Invert.$(OBJEXT) Invertlev.$(OBJEXT) \
-	Log.$(OBJEXT) Maskbox.$(OBJEXT) Mastrfu.$(OBJEXT) \
-	Math.$(OBJEXT) Merge.$(OBJEXT) Mergegrid.$(OBJEXT) \
-	Mergetime.$(OBJEXT) Merstat.$(OBJEXT) Monarith.$(OBJEXT) \
-	Mrotuv.$(OBJEXT) Mrotuvb.$(OBJEXT) Ninfo.$(OBJEXT) \
-	Nmltest.$(OBJEXT) Output.$(OBJEXT) Outputgmt.$(OBJEXT) \
-	Pinfo.$(OBJEXT) Regres.$(OBJEXT) Remap.$(OBJEXT) \
-	Remapeta.$(OBJEXT) Replace.$(OBJEXT) Rotuv.$(OBJEXT) \
-	Runpctl.$(OBJEXT) Runstat.$(OBJEXT) Seascount.$(OBJEXT) \
-	Seaspctl.$(OBJEXT) Seasstat.$(OBJEXT) Selbox.$(OBJEXT) \
-	Select.$(OBJEXT) Seloperator.$(OBJEXT) Selrec.$(OBJEXT) \
-	Seltime.$(OBJEXT) Set.$(OBJEXT) Setbox.$(OBJEXT) \
-	Setgatt.$(OBJEXT) Setgrid.$(OBJEXT) Sethalo.$(OBJEXT) \
-	Setmiss.$(OBJEXT) Setrange.$(OBJEXT) Setrcaname.$(OBJEXT) \
-	Settime.$(OBJEXT) Setzaxis.$(OBJEXT) Showinfo.$(OBJEXT) \
-	Sinfo.$(OBJEXT) Smooth9.$(OBJEXT) Sort.$(OBJEXT) \
-	Specinfo.$(OBJEXT) Spectral.$(OBJEXT) Spectrum.$(OBJEXT) \
-	Split.$(OBJEXT) Splitrec.$(OBJEXT) Splitsel.$(OBJEXT) \
-	Splittime.$(OBJEXT) Splityear.$(OBJEXT) Subtrend.$(OBJEXT) \
-	Templates.$(OBJEXT) Test.$(OBJEXT) Tests.$(OBJEXT) \
-	Timcount.$(OBJEXT) Timpctl.$(OBJEXT) Timsort.$(OBJEXT) \
-	Timselpctl.$(OBJEXT) Timselstat.$(OBJEXT) Timstat.$(OBJEXT) \
-	Timstat2.$(OBJEXT) Timstat3.$(OBJEXT) Tinfo.$(OBJEXT) \
-	Trend.$(OBJEXT) Trms.$(OBJEXT) Vardup.$(OBJEXT) \
-	Vargen.$(OBJEXT) Varrms.$(OBJEXT) Vertint.$(OBJEXT) \
-	Vertstat.$(OBJEXT) Wind.$(OBJEXT) Writegrid.$(OBJEXT) \
-	Writerandom.$(OBJEXT) Yhourstat.$(OBJEXT) Ydaypctl.$(OBJEXT) \
-	Ydaystat.$(OBJEXT) Ydrunpctl.$(OBJEXT) Ydrunstat.$(OBJEXT) \
-	Ymonarith.$(OBJEXT) Ymonpctl.$(OBJEXT) Ymonstat.$(OBJEXT) \
-	Yseaspctl.$(OBJEXT) Yseasstat.$(OBJEXT) Zonstat.$(OBJEXT) \
-	cdilib.$(OBJEXT) commandline.$(OBJEXT) exception.$(OBJEXT) \
-	expr_yacc.$(OBJEXT) expr_lex.$(OBJEXT) expr.$(OBJEXT) \
-	cdo.$(OBJEXT) cdo_pthread.$(OBJEXT) cdo_vlist.$(OBJEXT) \
-	field.$(OBJEXT) fieldc.$(OBJEXT) field2.$(OBJEXT) \
-	fieldmer.$(OBJEXT) fieldzon.$(OBJEXT) grid.$(OBJEXT) \
-	hetaeta.$(OBJEXT) history.$(OBJEXT) institution.$(OBJEXT) \
-	interpol.$(OBJEXT) job.$(OBJEXT) modules.$(OBJEXT) \
-	namelist.$(OBJEXT) normal.$(OBJEXT) pipe.$(OBJEXT) \
-	process.$(OBJEXT) remaplib.$(OBJEXT) timer.$(OBJEXT) \
-	realtime.$(OBJEXT) pstream.$(OBJEXT) table.$(OBJEXT) \
-	userlog.$(OBJEXT) util.$(OBJEXT) legendre.$(OBJEXT) \
-	fourier.$(OBJEXT) specspace.$(OBJEXT) readline.$(OBJEXT) \
-	julian.$(OBJEXT) vinterp.$(OBJEXT) zaxis.$(OBJEXT) \
-	pthread_debug.$(OBJEXT) color.$(OBJEXT) list.$(OBJEXT) \
-	percentiles.$(OBJEXT) nth_element.$(OBJEXT) ecacore.$(OBJEXT) \
-	ecautil.$(OBJEXT) EcaIndices.$(OBJEXT) Hi.$(OBJEXT) \
-	Wct.$(OBJEXT) statistic.$(OBJEXT)
+	Intgridtraj.$(OBJEXT) Intlevel.$(OBJEXT) Inttime.$(OBJEXT) \
+	Intntime.$(OBJEXT) Intyear.$(OBJEXT) Invert.$(OBJEXT) \
+	Invertlev.$(OBJEXT) Log.$(OBJEXT) Maskbox.$(OBJEXT) \
+	Mastrfu.$(OBJEXT) Math.$(OBJEXT) Merge.$(OBJEXT) \
+	Mergegrid.$(OBJEXT) Mergetime.$(OBJEXT) Merstat.$(OBJEXT) \
+	Monarith.$(OBJEXT) Mrotuv.$(OBJEXT) Mrotuvb.$(OBJEXT) \
+	Ninfo.$(OBJEXT) Nmltest.$(OBJEXT) Output.$(OBJEXT) \
+	Outputgmt.$(OBJEXT) Pinfo.$(OBJEXT) Regres.$(OBJEXT) \
+	Remap.$(OBJEXT) Remapeta.$(OBJEXT) Replace.$(OBJEXT) \
+	Rotuv.$(OBJEXT) Runpctl.$(OBJEXT) Runstat.$(OBJEXT) \
+	Seascount.$(OBJEXT) Seaspctl.$(OBJEXT) Seasstat.$(OBJEXT) \
+	Selbox.$(OBJEXT) Select.$(OBJEXT) Seloperator.$(OBJEXT) \
+	Selrec.$(OBJEXT) Seltime.$(OBJEXT) Set.$(OBJEXT) \
+	Setbox.$(OBJEXT) Setgatt.$(OBJEXT) Setgrid.$(OBJEXT) \
+	Sethalo.$(OBJEXT) Setmiss.$(OBJEXT) Setrange.$(OBJEXT) \
+	Setrcaname.$(OBJEXT) Settime.$(OBJEXT) Setzaxis.$(OBJEXT) \
+	Showinfo.$(OBJEXT) Sinfo.$(OBJEXT) Smooth9.$(OBJEXT) \
+	Sort.$(OBJEXT) Specinfo.$(OBJEXT) Spectral.$(OBJEXT) \
+	Spectrum.$(OBJEXT) Split.$(OBJEXT) Splitrec.$(OBJEXT) \
+	Splitsel.$(OBJEXT) Splittime.$(OBJEXT) Splityear.$(OBJEXT) \
+	Subtrend.$(OBJEXT) Templates.$(OBJEXT) Test.$(OBJEXT) \
+	Tests.$(OBJEXT) Timcount.$(OBJEXT) Timpctl.$(OBJEXT) \
+	Timsort.$(OBJEXT) Timselpctl.$(OBJEXT) Timselstat.$(OBJEXT) \
+	Timstat.$(OBJEXT) Timstat2.$(OBJEXT) Timstat3.$(OBJEXT) \
+	Tinfo.$(OBJEXT) Trend.$(OBJEXT) Trms.$(OBJEXT) \
+	Vardup.$(OBJEXT) Vargen.$(OBJEXT) Varrms.$(OBJEXT) \
+	Vertint.$(OBJEXT) Vertstat.$(OBJEXT) Wind.$(OBJEXT) \
+	Writegrid.$(OBJEXT) Writerandom.$(OBJEXT) Yhourstat.$(OBJEXT) \
+	Ydaypctl.$(OBJEXT) Ydaystat.$(OBJEXT) Ydrunpctl.$(OBJEXT) \
+	Ydrunstat.$(OBJEXT) Ymonarith.$(OBJEXT) Ymonpctl.$(OBJEXT) \
+	Ymonstat.$(OBJEXT) Yseaspctl.$(OBJEXT) Yseasstat.$(OBJEXT) \
+	Zonstat.$(OBJEXT) cdilib.$(OBJEXT) commandline.$(OBJEXT) \
+	exception.$(OBJEXT) expr_yacc.$(OBJEXT) expr_lex.$(OBJEXT) \
+	expr.$(OBJEXT) cdo.$(OBJEXT) cdo_pthread.$(OBJEXT) \
+	cdo_vlist.$(OBJEXT) field.$(OBJEXT) fieldc.$(OBJEXT) \
+	field2.$(OBJEXT) fieldmer.$(OBJEXT) fieldzon.$(OBJEXT) \
+	grid.$(OBJEXT) hetaeta.$(OBJEXT) history.$(OBJEXT) \
+	institution.$(OBJEXT) interpol.$(OBJEXT) job.$(OBJEXT) \
+	modules.$(OBJEXT) namelist.$(OBJEXT) normal.$(OBJEXT) \
+	pipe.$(OBJEXT) process.$(OBJEXT) remaplib.$(OBJEXT) \
+	timer.$(OBJEXT) realtime.$(OBJEXT) pstream.$(OBJEXT) \
+	table.$(OBJEXT) userlog.$(OBJEXT) util.$(OBJEXT) \
+	legendre.$(OBJEXT) fourier.$(OBJEXT) specspace.$(OBJEXT) \
+	readline.$(OBJEXT) julian.$(OBJEXT) vinterp.$(OBJEXT) \
+	zaxis.$(OBJEXT) pthread_debug.$(OBJEXT) color.$(OBJEXT) \
+	list.$(OBJEXT) percentiles.$(OBJEXT) nth_element.$(OBJEXT) \
+	ecacore.$(OBJEXT) ecautil.$(OBJEXT) EcaIndices.$(OBJEXT) \
+	Hi.$(OBJEXT) Wct.$(OBJEXT) statistic.$(OBJEXT)
 cdo_OBJECTS = $(am_cdo_OBJECTS)
 am__DEPENDENCIES_1 =
 cdo_DEPENDENCIES = $(am__DEPENDENCIES_1)
@@ -257,6 +257,7 @@ cdo_SOURCES = Arith.c		\
               Input.c		\
               Intgrid.c 	\
               Intgridtraj.c 	\
+              Intlevel.c 	\
               Inttime.c 	\
               Intntime.c 	\
               Intyear.c 	\
@@ -568,6 +569,7 @@ distclean-compile:
 @AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/Input.Po@am__quote@
 @AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/Intgrid.Po@am__quote@
 @AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/Intgridtraj.Po@am__quote@
+@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/Intlevel.Po@am__quote@
 @AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/Intntime.Po@am__quote@
 @AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/Inttime.Po@am__quote@
 @AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/Intyear.Po@am__quote@

src/Vertint.c

@@ -31,7 +31,6 @@
 #include "cdo.h"
 #include "cdo_int.h"
 #include "pstream.h"
-#include "functs.h"
 #include "vinterp.h"
 #include "list.h"
 
@@ -453,6 +452,7 @@ void *Vertint(void *argument)
   while ( (nrecs = streamInqTimestep(streamID1, tsID)) )
     {
       for ( varID = 0; varID < nvars; ++varID ) vars[varID] = FALSE;
+
       taxisCopyTimestep(taxisID2, taxisID1);
 
       streamDefTimestep(streamID2, tsID);

src/modules.c

@@ -68,6 +68,7 @@ void *Info(void *argument);
 void *Input(void *argument);
 void *Intgrid(void *argument);
 void *Intgridtraj(void *argument);
+void *Intlevel(void *argument);
 void *Inttime(void *argument);
 void *Intntime(void *argument);
 void *Intyear(void *argument);
@@ -259,6 +260,7 @@ void *Wct(void *argument);
 #define  InputOperators         {"input", "inputsrv", "inputext"}
 #define  IntgridOperators       {"intgridbil", "intpoint", "interpolate", "intarea", "boxavg", "thinout"}
 #define  IntgridtrajOperators   {"intgridtraj"}
+#define  IntlevelOperators      {"intlevel", "intlevelx"}
 #define  InttimeOperators       {"inttime"}
 #define  IntntimeOperators      {"intntime"}