sergiotarxz/libvips
1
0
Fork 0
Code Issues Pull Requests Releases Wiki Activity

avoid NaN in max/min, better double ranges

Browse Source
This commit is contained in:
John Cupitt 2012-02-24 13:57:50 +00:00
parent 1fffc6a4de
commit 862cac1e4f
9 changed files with 248 additions and 126 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

8
TODO
View File

@ -7,8 +7,6 @@ blocking bugs
`gdouble' is invalid or out of range for property `yres' of
type `gdouble'
similarly for vips_max() of an image which contains NaN, -INFINITY etc.
- new binding is still missing constants
how do boxed types work? confusing
@ -141,6 +139,12 @@ arithmetic
- HAVE_HYPOT could define a hypot() macro?
- fix a better NaN policy
should we not generate images containing NaN (eg. divide tries to avoid /0),
or should vips_max() etc. try to avoid NaN in images (eg. vips_max() takes a
lot a care to skip NaN, though vips_stats() does not)?
iofuncs
=======

2
libvips/arithmetic/avg.c
View File

@ -219,7 +219,7 @@ vips_avg_class_init( VipsAvgClass *class )
_( "Output value" ),
VIPS_ARGUMENT_REQUIRED_OUTPUT,
G_STRUCT_OFFSET( VipsAvg, out ),
-G_MAXDOUBLE, G_MAXDOUBLE, 0.0 );
-INFINITY, INFINITY, 0.0 );
}
static void

2
libvips/arithmetic/deviate.c
View File

@ -220,7 +220,7 @@ vips_deviate_class_init( VipsDeviateClass *class )
_( "Output value" ),
VIPS_ARGUMENT_REQUIRED_OUTPUT,
G_STRUCT_OFFSET( VipsDeviate, out ),
-G_MAXDOUBLE, G_MAXDOUBLE, 0.0 );
-INFINITY, INFINITY, 0.0 );
}
static void

143
libvips/arithmetic/max.c
View File

@ -18,6 +18,9 @@
* 6/11/11
* - rewrite as a class
* - abandon scan if we find maximum possible value
* 24/2/12
* - avoid NaN in float/double/complex images
* - allow +/- INFINITY as a result
*/
/*
@ -90,25 +93,31 @@ vips_max_build( VipsObject *object )
VipsStatistic *statistic = VIPS_STATISTIC( object );
VipsMax *max = (VipsMax *) object;
double m;
if( VIPS_OBJECT_CLASS( vips_max_parent_class )->build( object ) )
return( -1 );
/* For speed we accumulate max^2 for complex.
/* Don't set if there's no value (eg. if every pixel is NaN). This
* will trigger an error later.
*/
m = max->max;
if( vips_bandfmt_iscomplex( vips_image_get_format( statistic->in ) ) )
m = sqrt( m );
if( max->set ) {
double m;
/* We have to set the props via g_object_set() to stop vips
* complaining they are unset.
*/
g_object_set( max,
"out", m,
"x", max->x,
"y", max->y,
NULL );
/* For speed we accumulate max^2 for complex.
*/
m = max->max;
if( vips_bandfmt_iscomplex(
vips_image_get_format( statistic->in ) ) )
m = sqrt( m );
/* We have to set the props via g_object_set() to stop vips
* complaining they are unset.
*/
g_object_set( max,
"out", m,
"x", max->x,
"y", max->y,
NULL );
}
return( 0 );
}
@ -147,32 +156,6 @@ vips_max_stop( VipsStatistic *statistic, void *seq )
return( 0 );
}
/* real max with no limits.
*/
#define LOOP( TYPE ) { \
TYPE *p = (TYPE *) in; \
TYPE m; \
\
if( max->set ) \
m = max->max; \
else { \
m = p[0]; \
max->x = x; \
max->y = y; \
} \
\
for( i = 0; i < sz; i++ ) { \
if( p[i] > m ) { \
m = p[i]; \
max->x = x + i / bands; \
max->y = y; \
} \
} \
\
max->max = m; \
max->set = TRUE; \
}
/* real max with an upper bound.
*/
#define LOOPU( TYPE, UPPER ) { \
@ -203,33 +186,77 @@ vips_max_stop( VipsStatistic *statistic, void *seq )
max->set = TRUE; \
}
#define CLOOP( TYPE ) { \
/* float/double max ... no limits, and we have to avoid NaN.
*
* NaN compares false to every float value, so if we were to take the first
* point in this buffer as our start max (as we do above) and it was NaN, we'd
* never replace it with a true value.
*/
#define LOOPF( TYPE ) { \
TYPE *p = (TYPE *) in; \
double m; \
TYPE m; \
gboolean set; \
\
if( max->set ) \
m = max->max; \
else { \
m = p[0] * p[0] + p[1] * p[1]; \
max->x = x; \
max->y = y; \
} \
set = max->set; \
m = max->max; \
\
for( i = 0; i < sz; i++ ) { \
double mod; \
if( set ) { \
if( p[i] > m ) { \
m = p[i]; \
max->x = x + i / bands; \
max->y = y; \
} \
} \
else if( !isnan( p[i] ) ) { \
m = p[i]; \
max->x = x + i / bands; \
max->y = y; \
set = TRUE; \
} \
} \
\
if( set ) { \
max->max = m; \
max->set = TRUE; \
} \
}
/* As LOOPF, but complex. Track max(mod) to avoid sqrt().
*/
#define LOOPC( TYPE ) { \
TYPE *p = (TYPE *) in; \
TYPE m; \
gboolean set; \
\
set = max->set; \
m = max->max; \
\
for( i = 0; i < sz; i++ ) { \
TYPE mod; \
\
mod = p[0] * p[0] + p[1] * p[1]; \
p += 2; \
\
if( mod > m ) { \
if( set ) { \
if( mod > m ) { \
m = mod; \
max->x = x + i / bands; \
max->y = y; \
} \
} \
else if( !isnan( mod ) ) { \
m = mod; \
max->x = x + i / bands; \
max->y = y; \
set = TRUE; \
} \
} \
\
max->max = m; \
max->set = TRUE; \
if( set ) { \
max->max = m; \
max->set = TRUE; \
} \
}
/* Loop over region, adding to seq.
@ -259,14 +286,14 @@ vips_max_scan( VipsStatistic *statistic, void *seq,
LOOPU( signed int, INT_MAX ); break;
case VIPS_FORMAT_FLOAT:
LOOP( float ); break;
LOOPF( float ); break;
case VIPS_FORMAT_DOUBLE:
LOOP( double ); break;
LOOPF( double ); break;
case VIPS_FORMAT_COMPLEX:
CLOOP( float ); break;
LOOPC( float ); break;
case VIPS_FORMAT_DPCOMPLEX:
CLOOP( double ); break;
LOOPC( double ); break;
default:
g_assert( 0 );
@ -298,7 +325,7 @@ vips_max_class_init( VipsMaxClass *class )
_( "Output value" ),
VIPS_ARGUMENT_REQUIRED_OUTPUT,
G_STRUCT_OFFSET( VipsMax, max ),
-G_MAXDOUBLE, G_MAXDOUBLE, 0.0 );
-INFINITY, INFINITY, 0.0 );
VIPS_ARG_INT( class, "x", 2,
_( "x" ),

143
libvips/arithmetic/min.c
View File

@ -19,6 +19,9 @@
* - rewrite as a class
* 5/9/11
* - abandon scan if we find minimum possible value
* 24/2/12
* - avoid NaN in float/double/complex images
* - allow +/- INFINITY as a result
*/
/*
@ -91,25 +94,31 @@ vips_min_build( VipsObject *object )
VipsStatistic *statistic = VIPS_STATISTIC( object );
VipsMin *min = (VipsMin *) object;
double m;
if( VIPS_OBJECT_CLASS( vips_min_parent_class )->build( object ) )
return( -1 );
/* For speed we accumulate min^2 for complex.
/* Don't set if there's no value (eg. if every pixel is NaN). This
* will trigger an error later.
*/
m = min->min;
if( vips_bandfmt_iscomplex( vips_image_get_format( statistic->in ) ) )
m = sqrt( m );
if( min->set ) {
double m;
/* We have to set the props via g_object_set() to stop vips
* complaining they are unset.
*/
g_object_set( min,
"out", m,
"x", min->x,
"y", min->y,
NULL );
/* For speed we accumulate min^2 for complex.
*/
m = min->min;
if( vips_bandfmt_iscomplex(
vips_image_get_format( statistic->in ) ) )
m = sqrt( m );
/* We have to set the props via g_object_set() to stop vips
* complaining they are unset.
*/
g_object_set( min,
"out", m,
"x", min->x,
"y", min->y,
NULL );
}
return( 0 );
}
@ -148,32 +157,6 @@ vips_min_stop( VipsStatistic *statistic, void *seq )
return( 0 );
}
/* real min with no limits.
*/
#define LOOP( TYPE ) { \
TYPE *p = (TYPE *) in; \
TYPE m; \
\
if( min->set ) \
m = min->min; \
else { \
m = p[0]; \
min->x = x; \
min->y = y; \
} \
\
for( i = 0; i < sz; i++ ) { \
if( p[i] < m ) { \
m = p[i]; \
min->x = x + i / bands; \
min->y = y; \
} \
} \
\
min->min = m; \
min->set = TRUE; \
}
/* real min with a lower bound.
*/
#define LOOPL( TYPE, LOWER ) { \
@ -204,33 +187,77 @@ vips_min_stop( VipsStatistic *statistic, void *seq )
min->set = TRUE; \
}
#define CLOOP( TYPE ) { \
/* float/double min ... no limits, and we have to avoid NaN.
*
* NaN compares false to every float value, so if we were to take the first
* point in this buffer as our start min (as we do above) and it was NaN, we'd
* never replace it with a true value.
*/
#define LOOPF( TYPE ) { \
TYPE *p = (TYPE *) in; \
double m; \
TYPE m; \
gboolean set; \
\
if( min->set ) \
m = min->min; \
else { \
m = p[0] * p[0] + p[1] * p[1]; \
min->x = x; \
min->y = y; \
} \
set = min->set; \
m = min->min; \
\
for( i = 0; i < sz; i++ ) { \
double mod; \
if( set ) { \
if( p[i] < m ) { \
m = p[i]; \
min->x = x + i / bands; \
min->y = y; \
} \
} \
else if( !isnan( p[i] ) ) { \
m = p[i]; \
min->x = x + i / bands; \
min->y = y; \
set = TRUE; \
} \
} \
\
if( set ) { \
min->min = m; \
min->set = TRUE; \
} \
}
/* As LOOPF, but complex. Track min(mod) to avoid sqrt().
*/
#define LOOPC( TYPE ) { \
TYPE *p = (TYPE *) in; \
TYPE m; \
gboolean set; \
\
set = min->set; \
m = min->min; \
\
for( i = 0; i < sz; i++ ) { \
TYPE mod; \
\
mod = p[0] * p[0] + p[1] * p[1]; \
p += 2; \
\
if( mod < m ) { \
if( set ) { \
if( mod > m ) { \
m = mod; \
min->x = x + i / bands; \
min->y = y; \
} \
} \
else if( !isnan( mod ) ) { \
m = mod; \
min->x = x + i / bands; \
min->y = y; \
set = TRUE; \
} \
} \
\
min->min = m; \
min->set = TRUE; \
if( set ) { \
min->min = m; \
min->set = TRUE; \
} \
}
/* Loop over region, adding to seq.
@ -260,14 +287,14 @@ vips_min_scan( VipsStatistic *statistic, void *seq,
LOOPL( signed int, INT_MIN ); break;
case VIPS_FORMAT_FLOAT:
LOOP( float ); break;
LOOPF( float ); break;
case VIPS_FORMAT_DOUBLE:
LOOP( double ); break;
LOOPF( double ); break;
case VIPS_FORMAT_COMPLEX:
CLOOP( float ); break;
LOOPC( float ); break;
case VIPS_FORMAT_DPCOMPLEX:
CLOOP( double ); break;
LOOPC( double ); break;
default:
g_assert( 0 );
@ -299,7 +326,7 @@ vips_min_class_init( VipsMinClass *class )
_( "Output value" ),
VIPS_ARGUMENT_REQUIRED_OUTPUT,
G_STRUCT_OFFSET( VipsMin, min ),
-G_MAXDOUBLE, G_MAXDOUBLE, 0.0 );
-INFINITY, INFINITY, 0.0 );
VIPS_ARG_INT( class, "x", 2,
_( "x" ),

64
libvips/arithmetic/stats.c
View File

@ -309,6 +309,70 @@ vips_stats_start( VipsStatistic *statistic )
local->set = TRUE; \
}
/* As above, but for float/double types where we have to avoid NaN.
*/
#define LOOPF( TYPE ) { \
for( b = 0; b < bands; b++ ) { \
TYPE *p = ((TYPE *) in) + b; \
double *q = ARY( local->out, 0, b + 1 ); \
TYPE small, big; \
double sum, sum2; \
int xmin, ymin; \
int xmax, ymax; \
\
if( local->set ) { \
small = q[COL_MIN]; \
big = q[COL_MAX]; \
sum = q[COL_SUM]; \
sum2 = q[COL_SUM2]; \
xmin = q[COL_XMIN]; \
ymin = q[COL_YMIN]; \
xmax = q[COL_XMAX]; \
ymax = q[COL_YMAX]; \
} \
else { \
small = p[0]; \
big = p[0]; \
sum = 0; \
sum2 = 0; \
xmin = x; \
ymin = y; \
xmax = x; \
ymax = y; \
} \
\
for( i = 0; i < n; i++ ) { \
TYPE value = *p; \
\
sum += value; \
sum2 += (double) value * (double) value; \
if( value > big ) { \
big = value; \
xmax = x + i; \
ymax = y; \
} \
else if( value < small ) { \
small = value; \
xmin = x + i; \
ymin = y; \
} \
\
p += bands; \
} \
\
q[COL_MIN] = small; \
q[COL_MAX] = big; \
q[COL_SUM] = sum; \
q[COL_SUM2] = sum2; \
q[COL_XMIN] = xmin; \
q[COL_YMIN] = ymin; \
q[COL_XMAX] = xmax; \
q[COL_YMAX] = ymax; \
} \
\
local->set = TRUE; \
}
/* Loop over region, accumulating a sum in *tmp.
*/
static int

4
libvips/conversion/copy.c
View File

@ -366,14 +366,14 @@ vips_copy_class_init( VipsCopyClass *class )
_( "Horizontal resolution in pixels/mm" ),
VIPS_ARGUMENT_OPTIONAL_INPUT,
G_STRUCT_OFFSET( VipsCopy, xres ),
0, 1000000, 0 );
-0.0, 1000000, 0 );
VIPS_ARG_DOUBLE( class, "yres", 10,
_( "Yres" ),
_( "Vertical resolution in pixels/mm" ),
VIPS_ARGUMENT_OPTIONAL_INPUT,
G_STRUCT_OFFSET( VipsCopy, yres ),
0, 1000000, 0 );
-0.0, 1000000, 0 );
VIPS_ARG_INT( class, "xoffset", 11,
_( "Xoffset" ),

4
libvips/foreign/tiffsave.c
View File

@ -250,14 +250,14 @@ vips_foreign_save_tiff_class_init( VipsForeignSaveTiffClass *class )
_( "Horizontal resolution in pixels/mm" ),
VIPS_ARGUMENT_NONE,
G_STRUCT_OFFSET( VipsImage, Xres ),
0, 1000000, 0 );
-0.0, 1000000, 0 );
VIPS_ARG_DOUBLE( class, "yres", 17,
_( "Yres" ),
_( "Vertical resolution in pixels/mm" ),
VIPS_ARGUMENT_NONE,
G_STRUCT_OFFSET( VipsImage, Yres ),
0, 1000000, 0 );
-0.0, 1000000, 0 );
VIPS_ARG_BOOL( class, "bigtiff", 18,
_( "Bigtiff" ),

4
libvips/iofuncs/image.c
View File

@ -888,14 +888,14 @@ vips_image_class_init( VipsImageClass *class )
_( "Horizontal resolution in pixels/mm" ),
VIPS_ARGUMENT_SET_ALWAYS,
G_STRUCT_OFFSET( VipsImage, Xres ),
0, 1000000, 0 );
-0.0, 1000000, 0 );
VIPS_ARG_DOUBLE( class, "yres", 9,
_( "Yres" ),
_( "Vertical resolution in pixels/mm" ),
VIPS_ARGUMENT_SET_ALWAYS,
G_STRUCT_OFFSET( VipsImage, Yres ),
0, 1000000, 0 );
-0.0, 1000000, 0 );
VIPS_ARG_INT( class, "xoffset", 10,
_( "Xoffset" ),