Commit def325f4 authored by Uwe Schulzweida's avatar Uwe Schulzweida
Browse files

add support for NaN in DBL_IS_EQUAL

parent c31a5b36
2009-02-?? Uwe Schulzweida <Uwe.Schulzweida@zmaw.de>
2009-03-?? Uwe Schulzweida <Uwe.Schulzweida@zmaw.de>
* add support for NaN in DBL_IS_EQUAL
* add support for GRID type LCC2 (LCC PROJ.4 version)
* streamOpenAppen: set ncmode to 2 (bug fix)
* Version 1.3.0 released
......
......@@ -14,16 +14,26 @@
#ifndef DBL_IS_EQUAL
#if defined (HAVE_ISNAN)
# define DBL_IS_EQUAL(x,y) (isnan(x)||isnan(y)?(isnan(x)&&isnan(y)?1:0):!(fabs(x - y) > 0))
#if defined (XXX_HAVE_ISNAN)
# define DBL_IS_NAN(x) (isnan(x))
# define DBL_IS_EQUAL(x,y) (DBL_IS_NAN(x)||DBL_IS_NAN(y)?(DBL_IS_NAN(x)&&DBL_IS_NAN(y)?1:0):!(x < y || y < x))
#elif defined (XXX_FP_NAN)
# define DBL_IS_NAN(x) (fpclassify(x) == FP_NAN)
# define DBL_IS_EQUAL(x,y) (DBL_IS_NAN(x)||DBL_IS_NAN(y)?(DBL_IS_NAN(x)&&DBL_IS_NAN(y)?1:0):!(x < y || y < x))
#else
# define DBL_IS_NAN(x) ((x) != (x))
/*
# define DBL_IS_EQUAL(x,y) (fabs(x - y) <= 2.0*(y*DBL_EPSILON + DBL_MIN))
# define DBL_IS_EQUAL(x,y) (!(x < y || y < x))
*/
# define DBL_IS_EQUAL(x,y) (!(fabs(x - y) > 0))
# define DBL_IS_EQUAL(x,y) (DBL_IS_NAN(x)||DBL_IS_NAN(y)?(DBL_IS_NAN(x)&&DBL_IS_NAN(y)?1:0):!(x < y || y < x))
#endif
#endif
#ifndef IS_EQUAL
# define IS_NOT_EQUAL(x,y) (x < y || y < x)
# define IS_EQUAL(x,y) (!IS_NOT_EQUAL(x,y))
#endif
void decode_date(int date, int *year, int *month, int *day);
void decode_time(int time, int *hour, int *minute);
......@@ -135,7 +145,7 @@ static void printInfo(int gridtype, int date, int time, int code, double level,
arrmax = -1e50;
for ( i = 0; i < datasize; i++ )
{
if ( fabs(data[i]-missval) > 0 )
if ( !DBL_IS_EQUAL(data[i], missval) )
{
if ( data[i] < arrmin ) arrmin = data[i];
if ( data[i] > arrmax ) arrmax = data[i];
......@@ -294,7 +304,7 @@ static void printGridInfo(int vlistID)
fprintf(stdout, "size : dim = %d nlon = %d nlat = %d\n", gridsize, xsize, ysize);
fprintf(stdout, "%*s", nbyte0, "");
fprintf(stdout, "%-9s : first = %.9g last = %.9g", xname, xfirst, xlast);
if ( !DBL_IS_EQUAL(xinc, 0) )
if ( IS_NOT_EQUAL(xinc, 0) )
fprintf(stdout, " inc = %.9g", xinc);
fprintf(stdout, " %s", xunits);
if ( gridIsCircular(gridID) )
......@@ -303,7 +313,7 @@ static void printGridInfo(int vlistID)
}
fprintf(stdout, "%*s", nbyte0, "");
fprintf(stdout, "%-9s : first = %.9g last = %.9g", yname, yfirst, ylast);
if ( !DBL_IS_EQUAL(yinc, 0) &&
if ( IS_NOT_EQUAL(yinc, 0) &&
(gridtype == GRID_LONLAT || gridtype == GRID_SINUSOIDAL || gridtype == GRID_LAEA) )
fprintf(stdout, " inc = %.9g", yinc);
fprintf(stdout, " %s", yunits);
......
/* Generated automatically from m214003 on Thu Oct 23 13:39:51 CEST 2008 */
/* Generated automatically from m214003 on Tue Feb 17 11:26:33 CET 2009 */
/* GRIBLIB_VERSION="1.2.1" */
/* GRIBLIB_VERSION="1.2.2" */
#if defined (HAVE_CONFIG_H)
# include "config.h"
......@@ -20,6 +20,10 @@
#ifndef _GRIB_INT_H
#define _GRIB_INT_H
#if defined (HAVE_CONFIG_H)
# include "config.h"
#endif
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
......@@ -60,10 +64,24 @@
#ifndef DBL_IS_EQUAL
#if defined (XXX_HAVE_ISNAN)
# define DBL_IS_NAN(x) (isnan(x))
# define DBL_IS_EQUAL(x,y) (DBL_IS_NAN(x)||DBL_IS_NAN(y)?(DBL_IS_NAN(x)&&DBL_IS_NAN(y)?1:0):!(x < y || y < x))
#elif defined (XXX_FP_NAN)
# define DBL_IS_NAN(x) (fpclassify(x) == FP_NAN)
# define DBL_IS_EQUAL(x,y) (DBL_IS_NAN(x)||DBL_IS_NAN(y)?(DBL_IS_NAN(x)&&DBL_IS_NAN(y)?1:0):!(x < y || y < x))
#else
# define DBL_IS_NAN(x) ((x) != (x))
/*
#define DBL_IS_EQUAL(x,y) (fabs(x - y) <= 2.0*(y*DBL_EPSILON + DBL_MIN))
# define DBL_IS_EQUAL(x,y) (!(x < y || y < x))
*/
#define DBL_IS_EQUAL(x,y) (!(fabs(x - y) > 0))
# define DBL_IS_EQUAL(x,y) (DBL_IS_NAN(x)||DBL_IS_NAN(y)?(DBL_IS_NAN(x)&&DBL_IS_NAN(y)?1:0):!(x < y || y < x))
#endif
#endif
#ifndef IS_EQUAL
# define IS_NOT_EQUAL(x,y) (x < y || y < x)
# define IS_EQUAL(x,y) (!IS_NOT_EQUAL(x,y))
#endif
......@@ -3334,7 +3352,7 @@ int encodeBDS(GRIBPACK *lGrib, int *gribLen, int decscale, int *isec2, int *isec
{
binscale = 0;
}
else if ( (!DBL_IS_EQUAL(fmin, 0.0)) && (fabs(range/fmin) <= jpepsln) )
else if ( IS_NOT_EQUAL(fmin, 0.0) && (fabs(range/fmin) <= jpepsln) )
{
binscale = 0;
}
......@@ -8362,7 +8380,7 @@ int gribUnzip(unsigned char *dbuf, long dbufsize, unsigned char *sbuf, long sbu
return (gribLen);
}
static const char grb_libvers[] = "1.2.1" " of ""Oct 23 2008"" ""13:39:51";
static const char grb_libvers[] = "1.2.2" " of ""Feb 17 2009"" ""11:26:33";
const char *
gribLibraryVersion(void)
{
......
......@@ -387,7 +387,7 @@ void gridGenYvals(int gridtype, int ysize, double yfirst, double ylast, double y
{
calc_gaussaw(yvals, ysize, yfirst, ylast);
if ( ! (DBL_IS_EQUAL(yfirst, 0) && DBL_IS_EQUAL(ylast, 0)) )
if ( ! (IS_EQUAL(yfirst, 0) && IS_EQUAL(ylast, 0)) )
if ( fabs(yvals[0] - yfirst) > 0.001 || fabs(yvals[ysize-1] - ylast) > 0.001 )
{
double yinc = fabs(ylast-yfirst)/(ysize-1);
......@@ -432,7 +432,7 @@ void gridGenYvals(int gridtype, int ysize, double yfirst, double ylast, double y
{
if ( (! (fabs(yinc) > 0)) && ysize > 1 )
{
if ( DBL_IS_EQUAL(yfirst, ylast) && !DBL_IS_EQUAL(yfirst, 0) ) ylast *= -1;
if ( IS_EQUAL(yfirst, ylast) && IS_NOT_EQUAL(yfirst, 0) ) ylast *= -1;
if ( yfirst > ylast )
yinc = (yfirst-ylast)/(ysize-1);
......@@ -2142,9 +2142,9 @@ static void grid_check_cyclic(grid_t *gridptr)
{
xinc = xvals[1] - xvals[0];
x0 = 2*xvals[xsize-1]-xvals[xsize-2]-360;
if ( ! DBL_IS_EQUAL(xvals[0], xvals[xsize-1]) )
if ( IS_NOT_EQUAL(xvals[0], xvals[xsize-1]) )
if ( fabs(x0 - xvals[0]) < 1.e-4 ) gridptr->isCyclic = TRUE;
/* if ( DBL_IS_EQUAL(x0, xvals[0]) ) gridptr->isCyclic = TRUE; */
/* if ( IS_EQUAL(x0, xvals[0]) ) gridptr->isCyclic = TRUE; */
}
}
else if ( gridptr->type == GRID_CURVILINEAR )
......@@ -2246,12 +2246,12 @@ int gridCompare(int gridID, grid_t grid)
{
if ( grid.xdef == 2 && grid.ydef == 2 )
{
if ( ! (DBL_IS_EQUAL(grid.xfirst, 0) && DBL_IS_EQUAL(grid.xlast, 0) && DBL_IS_EQUAL(grid.xinc, 0)) &&
! (DBL_IS_EQUAL(grid.yfirst, 0) && DBL_IS_EQUAL(grid.ylast, 0) && DBL_IS_EQUAL(grid.yinc, 0)) &&
!DBL_IS_EQUAL(grid.xfirst, grid.xlast) && !DBL_IS_EQUAL(grid.yfirst, grid.ylast) )
if ( ! (IS_EQUAL(grid.xfirst, 0) && IS_EQUAL(grid.xlast, 0) && IS_EQUAL(grid.xinc, 0)) &&
! (IS_EQUAL(grid.yfirst, 0) && IS_EQUAL(grid.ylast, 0) && IS_EQUAL(grid.yinc, 0)) &&
IS_NOT_EQUAL(grid.xfirst, grid.xlast) && IS_NOT_EQUAL(grid.yfirst, grid.ylast) )
{
if ( !DBL_IS_EQUAL(grid.xfirst, gridInqXval(gridID, 0)) ||
!DBL_IS_EQUAL(grid.yfirst, gridInqYval(gridID, 0)))
if ( IS_NOT_EQUAL(grid.xfirst, gridInqXval(gridID, 0)) ||
IS_NOT_EQUAL(grid.yfirst, gridInqYval(gridID, 0)))
{
differ = 1;
}
......@@ -2305,8 +2305,8 @@ int gridCompare(int gridID, grid_t grid)
{
if ( grid.xdef == 2 && grid.ydef == 2 )
{
if ( ! (DBL_IS_EQUAL(grid.xfirst, 0) && DBL_IS_EQUAL(grid.xlast, 0) && DBL_IS_EQUAL(grid.xinc, 0)) &&
! (DBL_IS_EQUAL(grid.yfirst, 0) && DBL_IS_EQUAL(grid.ylast, 0)) )
if ( ! (IS_EQUAL(grid.xfirst, 0) && IS_EQUAL(grid.xlast, 0) && IS_EQUAL(grid.xinc, 0)) &&
! (IS_EQUAL(grid.yfirst, 0) && IS_EQUAL(grid.ylast, 0)) )
if ( fabs(grid.xfirst - gridInqXval(gridID, 0)) > 0.001 ||
fabs(grid.yfirst - gridInqYval(gridID, 0)) > 0.001 ||
fabs(fabs(grid.xinc) - fabs(gridInqXinc(gridID))) > fabs(grid.xinc/1000))
......@@ -2317,8 +2317,8 @@ int gridCompare(int gridID, grid_t grid)
else
{
if ( grid.xvals && grid.yvals )
if ( !DBL_IS_EQUAL(grid.xvals[0], gridInqXval(gridID, 0)) ||
!DBL_IS_EQUAL(grid.yvals[0], gridInqYval(gridID, 0)) ||
if ( IS_NOT_EQUAL(grid.xvals[0], gridInqXval(gridID, 0)) ||
IS_NOT_EQUAL(grid.yvals[0], gridInqYval(gridID, 0)) ||
fabs(fabs(grid.xvals[1]-grid.xvals[0]) - fabs(gridInqXinc(gridID)))
> fabs(grid.xinc/1000))
{
......@@ -3134,7 +3134,7 @@ void gridPrint(int gridID, int opt)
xinc = gridInqXinc(gridID);
}
if ( !DBL_IS_EQUAL(xinc, 0) && opt )
if ( IS_NOT_EQUAL(xinc, 0) && opt )
{
fprintf(fp, "xfirst = %g\n", xfirst);
fprintf(fp, "xinc = %g\n", xinc);
......@@ -3181,7 +3181,7 @@ void gridPrint(int gridID, int opt)
yinc = gridInqYinc(gridID);
}
if ( !DBL_IS_EQUAL(yinc, 0) && opt )
if ( IS_NOT_EQUAL(yinc, 0) && opt )
{
fprintf(fp, "yfirst = %g\n", yfirst);
fprintf(fp, "yinc = %g\n", yinc);
......
......@@ -1093,8 +1093,8 @@ void cdfDefXaxis(int streamID, int gridID)
{
dimlen0 = gridInqXsize(gridID0);
if ( dimlen == dimlen0 )
if ( DBL_IS_EQUAL(gridInqXval(gridID0, 0), gridInqXval(gridID, 0)) &&
DBL_IS_EQUAL(gridInqXval(gridID0, dimlen-1), gridInqXval(gridID, dimlen-1)) )
if ( IS_EQUAL(gridInqXval(gridID0, 0), gridInqXval(gridID, 0)) &&
IS_EQUAL(gridInqXval(gridID0, dimlen-1), gridInqXval(gridID, dimlen-1)) )
{
dimID = streamptr->xdimID[index];
break;
......@@ -1257,8 +1257,8 @@ void cdfDefYaxis(int streamID, int gridID)
{
dimlen0 = gridInqYsize(gridID0);
if ( dimlen == dimlen0 )
if ( DBL_IS_EQUAL(gridInqYval(gridID0, 0), gridInqYval(gridID, 0)) &&
DBL_IS_EQUAL(gridInqYval(gridID0, dimlen-1), gridInqYval(gridID, dimlen-1)) )
if ( IS_EQUAL(gridInqYval(gridID0, 0), gridInqYval(gridID, 0)) &&
IS_EQUAL(gridInqYval(gridID0, dimlen-1), gridInqYval(gridID, dimlen-1)) )
{
dimID = streamptr->ydimID[index];
break;
......@@ -1430,8 +1430,8 @@ void cdfDefLonLat2D(int streamID, int gridID)
{
dimlen0 = gridInqXsize(gridID0);
if ( xdimlen == dimlen0 )
if ( DBL_IS_EQUAL(gridInqXval(gridID0, 0), gridInqXval(gridID, 0)) &&
DBL_IS_EQUAL(gridInqXval(gridID0, xdimlen-1), gridInqXval(gridID, xdimlen-1)) )
if ( IS_EQUAL(gridInqXval(gridID0, 0), gridInqXval(gridID, 0)) &&
IS_EQUAL(gridInqXval(gridID0, xdimlen-1), gridInqXval(gridID, xdimlen-1)) )
{
xdimID = streamptr->xdimID[index];
break;
......@@ -2306,7 +2306,7 @@ void cdfDefMapping(int streamID, int gridID)
cdf_put_att_double(fileID, ncvarid, "earth_radius", NC_DOUBLE, 1L, &a);
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);
if ( DBL_IS_EQUAL(lat_1, lat_2) )
if ( IS_EQUAL(lat_1, lat_2) )
cdf_put_att_double(fileID, ncvarid, "standard_parallel", NC_DOUBLE, 1L, &lat_1);
else
{
......@@ -2662,13 +2662,13 @@ int cdfDefVar(int streamID, int varID)
addoffset = vlistInqVarAddoffset(vlistID, varID);
scalefactor = vlistInqVarScalefactor(vlistID, varID);
laddoffset = !DBL_IS_EQUAL(addoffset, 0);
lscalefactor = !DBL_IS_EQUAL(scalefactor, 1);
laddoffset = IS_NOT_EQUAL(addoffset, 0);
lscalefactor = IS_NOT_EQUAL(scalefactor, 1);
if ( laddoffset || lscalefactor )
{
if ( DBL_IS_EQUAL(addoffset, (double) ((float) addoffset)) &&
DBL_IS_EQUAL(scalefactor, (double) ((float) scalefactor)) )
if ( IS_EQUAL(addoffset, (double) ((float) addoffset)) &&
IS_EQUAL(scalefactor, (double) ((float) scalefactor)) )
{
astype = NC_FLOAT;
}
......@@ -2794,8 +2794,8 @@ void cdfReadVarDP(int streamID, int varID, double *data, int *nmiss)
addoffset = vlistInqVarAddoffset(vlistID, varID);
scalefactor = vlistInqVarScalefactor(vlistID, varID);
laddoffset = !DBL_IS_EQUAL(addoffset, 0);
lscalefactor = !DBL_IS_EQUAL(scalefactor, 1);
laddoffset = IS_NOT_EQUAL(addoffset, 0);
lscalefactor = IS_NOT_EQUAL(scalefactor, 1);
if ( laddoffset || lscalefactor )
{
......@@ -2933,8 +2933,8 @@ void cdfWriteVarDP(int streamID, int varID, double *data, int nmiss)
addoffset = vlistInqVarAddoffset(vlistID, varID);
scalefactor = vlistInqVarScalefactor(vlistID, varID);
laddoffset = !DBL_IS_EQUAL(addoffset, 0);
lscalefactor = !DBL_IS_EQUAL(scalefactor, 1);
laddoffset = IS_NOT_EQUAL(addoffset, 0);
lscalefactor = IS_NOT_EQUAL(scalefactor, 1);
nvals = gridInqSize(gridID)*zaxisInqSize(zaxisID);
......@@ -3177,8 +3177,8 @@ int cdfReadVarSliceDP(int streamID, int varID, int levelID, double *data, int *n
addoffset = vlistInqVarAddoffset(vlistID, varID);
scalefactor = vlistInqVarScalefactor(vlistID, varID);
laddoffset = !DBL_IS_EQUAL(addoffset, 0);
lscalefactor = !DBL_IS_EQUAL(scalefactor, 1);
laddoffset = IS_NOT_EQUAL(addoffset, 0);
lscalefactor = IS_NOT_EQUAL(scalefactor, 1);
if ( laddoffset || lscalefactor )
{
......@@ -3310,8 +3310,8 @@ int cdfWriteVarSliceDP(int streamID, int varID, int levelID, double *data, int n
addoffset = vlistInqVarAddoffset(vlistID, varID);
scalefactor = vlistInqVarScalefactor(vlistID, varID);
laddoffset = !DBL_IS_EQUAL(addoffset, 0);
lscalefactor = !DBL_IS_EQUAL(scalefactor, 1);
laddoffset = IS_NOT_EQUAL(addoffset, 0);
lscalefactor = IS_NOT_EQUAL(scalefactor, 1);
nvals = gridInqSize(gridID);
......@@ -4879,7 +4879,7 @@ int cdfInqContents(int streamID)
if ( islat && islon )
{
int lgauss = FALSE;
if ( DBL_IS_EQUAL(yinc, 0) && (int) ysize > 2 ) /* check if gaussian */
if ( IS_EQUAL(yinc, 0) && (int) ysize > 2 ) /* check if gaussian */
{
double *yvals, *yw;
yvals = (double *) malloc(ysize*sizeof(double));
......@@ -5298,9 +5298,9 @@ int cdfInqContents(int streamID)
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 ( !DBL_IS_EQUAL(ncvars[ncvarid].addoffset, 0) )
if ( IS_NOT_EQUAL(ncvars[ncvarid].addoffset, 0) )
vlistDefVarAddoffset(vlistID, varID, ncvars[ncvarid].addoffset);
if ( !DBL_IS_EQUAL(ncvars[ncvarid].scalefactor, 1) )
if ( IS_NOT_EQUAL(ncvars[ncvarid].scalefactor, 1) )
vlistDefVarScalefactor(vlistID, varID, ncvars[ncvarid].scalefactor);
vlistDefVarDatatype(vlistID, varID, cdfInqDatatype(ncvars[ncvarid].xtype));
......
......@@ -1463,8 +1463,8 @@ int grbScanTimestep2(int streamID)
gridID = vlistInqVarGrid(vlistID, varID);
if ( gridInqSize(gridID) == 1 && gridInqType(gridID) == GRID_LONLAT )
{
if ( !DBL_IS_EQUAL(gridInqXval(gridID, 0),ISEC2_FirstLon*0.001) ||
!DBL_IS_EQUAL(gridInqYval(gridID, 0),ISEC2_FirstLat*0.001) )
if ( IS_NOT_EQUAL(gridInqXval(gridID, 0),ISEC2_FirstLon*0.001) ||
IS_NOT_EQUAL(gridInqYval(gridID, 0),ISEC2_FirstLat*0.001) )
gridChangeType(gridID, GRID_TRAJECTORY);
}
}
......
......@@ -57,18 +57,29 @@ char *strdup(const char *s);
# include "ieg.h"
#endif
#ifndef DBL_IS_EQUAL
#if defined (XXX_HAVE_ISNAN)
# define DBL_IS_NAN(x) (isnan(x))
# define DBL_IS_NAN(x) (isnan(x))
# define DBL_IS_EQUAL(x,y) (DBL_IS_NAN(x)||DBL_IS_NAN(y)?(DBL_IS_NAN(x)&&DBL_IS_NAN(y)?1:0):!(x < y || y < x))
#elif defined (XXX_FP_NAN)
# define DBL_IS_NAN(x) (fpclassify(x) == FP_NAN)
# define DBL_IS_NAN(x) (fpclassify(x) == FP_NAN)
# define DBL_IS_EQUAL(x,y) (DBL_IS_NAN(x)||DBL_IS_NAN(y)?(DBL_IS_NAN(x)&&DBL_IS_NAN(y)?1:0):!(x < y || y < x))
#else
# define DBL_IS_NAN(x) ((x) != (x))
/*
# define DBL_IS_EQUAL(x,y) (!(x < y || y < x))
*/
# define DBL_IS_EQUAL(x,y) (DBL_IS_NAN(x)||DBL_IS_NAN(y)?(DBL_IS_NAN(x)&&DBL_IS_NAN(y)?1:0):!(x < y || y < x))
#endif
#endif
#ifndef IS_EQUAL
# define IS_NOT_EQUAL(x,y) (x < y || y < x)
# define IS_EQUAL(x,y) (!IS_NOT_EQUAL(x,y))
#endif
#ifndef INT
# define INT(x) ((int)(x))
#endif
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment