diff --git a/ChangeLog b/ChangeLog index 0e28ceac..70ee72ed 100644 --- a/ChangeLog +++ b/ChangeLog @@ -7,7 +7,7 @@ im_lintra(), im_lintra_vec(), im_black(), im_rot90, im_rot180(), im_rot270() im_sintra(), im_costra(), im_tantra(), im_asintra(), im_acostra(), im_atantra(), im_exptra(), im_exp10tra(), im_logtra(), im_log10tra(), - im_abs(), im_sign(), im_max(), im_maxpos() + im_abs(), im_sign(), im_max(), im_maxpos(), im_deviate() redone as classes - added argument priorites to help control arg ordering - generate has a 'stop' param to signal successful early termination diff --git a/libvips/arithmetic/Makefile.am b/libvips/arithmetic/Makefile.am index 08a1416e..c7401863 100644 --- a/libvips/arithmetic/Makefile.am +++ b/libvips/arithmetic/Makefile.am @@ -5,7 +5,7 @@ libarithmetic_la_SOURCES = \ abs.c \ im_bandmean.c \ im_cross_phase.c \ - im_deviate.c \ + deviate.c \ im_divide.c \ im_recomb.c \ im_linreg.c \ diff --git a/libvips/arithmetic/arithmetic.c b/libvips/arithmetic/arithmetic.c index 1419a877..c3bbda1b 100644 --- a/libvips/arithmetic/arithmetic.c +++ b/libvips/arithmetic/arithmetic.c @@ -507,9 +507,11 @@ vips_arithmetic_operation_init( void ) extern GType vips_avg_get_type( void ); extern GType vips_min_get_type( void ); extern GType vips_max_get_type( void ); + extern GType vips_deviate_get_type( void ); extern GType vips_linear_get_type( void ); extern GType vips_math_get_type( void ); extern GType vips_abs_get_type( void ); + extern GType vips_sign_get_type( void ); vips_add_get_type(); vips_invert_get_type(); @@ -517,8 +519,10 @@ vips_arithmetic_operation_init( void ) vips_avg_get_type(); vips_min_get_type(); vips_max_get_type(); + vips_deviate_get_type(); vips_linear_get_type(); vips_math_get_type(); vips_abs_get_type(); + vips_sign_get_type(); } diff --git a/libvips/arithmetic/avg.c b/libvips/arithmetic/avg.c index a75ac087..edab95dc 100644 --- a/libvips/arithmetic/avg.c +++ b/libvips/arithmetic/avg.c @@ -90,7 +90,7 @@ typedef struct _VipsAvg { VipsStatistic parent_instance; double sum; - double avg; + double out; } VipsAvg; typedef VipsStatisticClass VipsAvgClass; @@ -103,7 +103,7 @@ vips_avg_build( VipsObject *object ) VipsStatistic *statistic = VIPS_STATISTIC( object ); VipsAvg *avg = (VipsAvg *) object; - gint64 vals, pels; + gint64 vals; double average; if( VIPS_OBJECT_CLASS( vips_avg_parent_class )->build( object ) ) @@ -112,9 +112,10 @@ vips_avg_build( VipsObject *object ) /* Calculate average. For complex, we accumulate re*re + * im*im, so we need to sqrt. */ - pels = (gint64) vips_image_get_width( statistic->in ) * - vips_image_get_height( statistic->in ); - vals = pels * vips_image_get_bands( statistic->in ); + vals = (gint64) + vips_image_get_width( statistic->in ) * + vips_image_get_height( statistic->in ) * + vips_image_get_bands( statistic->in ); average = avg->sum / vals; if( vips_bandfmt_iscomplex( vips_image_get_format( statistic->in ) ) ) average = sqrt( average ); @@ -234,7 +235,7 @@ vips_avg_class_init( VipsAvgClass *class ) _( "Output" ), _( "Output value" ), VIPS_ARGUMENT_REQUIRED_OUTPUT, - G_STRUCT_OFFSET( VipsAvg, avg ), + G_STRUCT_OFFSET( VipsAvg, out ), -G_MAXDOUBLE, G_MAXDOUBLE, 0.0 ); } diff --git a/libvips/arithmetic/deviate.c b/libvips/arithmetic/deviate.c new file mode 100644 index 00000000..ed52c0c4 --- /dev/null +++ b/libvips/arithmetic/deviate.c @@ -0,0 +1,262 @@ +/* VipsDeviate + * + * Copyright: 1990, J. Cupitt + * + * Author: J. Cupitt + * Written on: 02/08/1990 + * Modified on: + * 5/5/93 JC + * - now does partial images + * - less likely to overflow + * - adapted from im_deviate + * 1/7/93 JC + * - adapted for partial v2 + * - ANSIfied + * 21/2/95 JC + * - modernised again + * 11/5/95 JC + * - oops! return( NULL ) in im_deviate(), instead of return( -1 ) + * 20/6/95 JC + * - now returns double, not float + * 13/1/05 + * - use 64 bit arithmetic + * 8/12/06 + * - add liboil support + * 2/9/09 + * - gtk-doc comment + * - minor reformatting + * 4/9/09 + * - use im__wrapscan() + * 31/7/10 + * - remove liboil + * 6/11/11 + * - rewrite as a class + */ + +/* + + This file is part of VIPS. + + VIPS is free software; you can redistribute it and/or modify + it under the terms of the GNU Lesser General Public License as published by + the Free Software Foundation; either version 2 of the License, or + (at your option) any later version. + + 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 Lesser General Public License for more details. + + You should have received a copy of the GNU Lesser General Public License + along with this program; if not, write to the Free Software + Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA + + */ + +/* + + These files are distributed with VIPS - http://www.vips.ecs.soton.ac.uk + + */ + +#ifdef HAVE_CONFIG_H +#include +#endif /*HAVE_CONFIG_H*/ +#include + +#include +#include +#include + +#include +#include + +#include "statistic.h" + +/** + * VipsDeviate: + * @in: input #IMAGE + * @out: output pixel standard deviation + * + * This operation finds the standard deviation of all pixels in @in. It + * operates on all bands of the input image: use im_stats() if you need + * to calculate an average for each band. + * + * Non-complex images only. + * + * See also: VipsDeviate. + */ + +typedef struct _VipsDeviate { + VipsStatistic parent_instance; + + double sum; + double sum2; + double out; +} VipsDeviate; + +typedef VipsStatisticClass VipsDeviateClass; + +G_DEFINE_TYPE( VipsDeviate, vips_deviate, VIPS_TYPE_STATISTIC ); + +static int +vips_deviate_build( VipsObject *object ) +{ + VipsStatistic *statistic = VIPS_STATISTIC( object ); + VipsDeviate *deviate = (VipsDeviate *) object; + + gint64 vals; + double s, s2; + + if( statistic->in && + vips_check_noncomplex( "VipsDeviate", statistic->in ) ) + return( -1 ); + + if( VIPS_OBJECT_CLASS( vips_deviate_parent_class )->build( object ) ) + return( -1 ); + + /* + + NOTE: NR suggests a two-pass algorithm to minimise roundoff. + But that's too expensive for us :-( so do it the old one-pass + way. + + */ + + /* Calculate and return deviation. Add a fabs to stop sqrt(<=0). + */ + vals = (gint64) + vips_image_get_width( statistic->in ) * + vips_image_get_height( statistic->in ) * + vips_image_get_bands( statistic->in ); + s = deviate->sum; + s2 = deviate->sum2; + + g_object_set( object, + "out", sqrt( fabs( s2 - (s * s / vals) ) / (vals - 1) ), + NULL ); + + return( 0 ); +} + +/* Start function: allocate space for an array in which we can accumulate the + * sum and sum of squares for this thread. + */ +static void * +vips_deviate_start( VipsStatistic *statistic ) +{ + double *ss2; + + ss2 = g_new( double, 2 ); + ss2[0] = 0.0; + ss2[1] = 0.0; + + return( (void *) ss2 ); +} + +/* Stop function. Add this little sum to the main sum. + */ +static int +vips_deviate_stop( VipsStatistic *statistic, void *seq ) +{ + VipsDeviate *deviate = (VipsDeviate *) statistic; + double *ss2 = (double *) seq; + + deviate->sum += ss2[0]; + deviate->sum2 += ss2[1]; + + g_free( ss2 ); + + return( 0 ); +} + +#define LOOP( TYPE ) { \ + TYPE *p = (TYPE *) in; \ + \ + for( x = 0; x < sz; x++ ) { \ + TYPE v = p[x]; \ + \ + sum += v; \ + sum2 += (double) v * (double) v; \ + } \ +} + +static int +vips_deviate_scan( VipsStatistic *statistic, void *seq, + int x, int y, void *in, int n ) +{ + const int sz = n * vips_image_get_bands( statistic->in ); + + double *ss2 = (double *) seq; + + double sum; + double sum2; + + sum = ss2[0]; + sum2 = ss2[1]; + + /* Now generate code for all types. + */ + switch( vips_image_get_format( statistic->in ) ) { + case VIPS_FORMAT_UCHAR: LOOP( unsigned char ); break; + case VIPS_FORMAT_CHAR: LOOP( signed char ); break; + case VIPS_FORMAT_USHORT: LOOP( unsigned short ); break; + case VIPS_FORMAT_SHORT: LOOP( signed short ); break; + case VIPS_FORMAT_UINT: LOOP( unsigned int ); break; + case VIPS_FORMAT_INT: LOOP( signed int ); break; + case VIPS_FORMAT_FLOAT: LOOP( float ); break; + case VIPS_FORMAT_DOUBLE: LOOP( double ); break; + + default: + g_assert( 0 ); + } + + ss2[0] = sum; + ss2[1] = sum2; + + return( 0 ); +} + +static void +vips_deviate_class_init( VipsDeviateClass *class ) +{ + GObjectClass *gobject_class = (GObjectClass *) class; + VipsObjectClass *object_class = (VipsObjectClass *) class; + VipsStatisticClass *sclass = VIPS_STATISTIC_CLASS( class ); + + gobject_class->set_property = vips_object_set_property; + gobject_class->get_property = vips_object_get_property; + + object_class->nickname = "deviate"; + object_class->description = _( "find image average" ); + object_class->build = vips_deviate_build; + + sclass->start = vips_deviate_start; + sclass->scan = vips_deviate_scan; + sclass->stop = vips_deviate_stop; + + VIPS_ARG_DOUBLE( class, "out", 2, + _( "Output" ), + _( "Output value" ), + VIPS_ARGUMENT_REQUIRED_OUTPUT, + G_STRUCT_OFFSET( VipsDeviate, out ), + -G_MAXDOUBLE, G_MAXDOUBLE, 0.0 ); +} + +static void +vips_deviate_init( VipsDeviate *deviate ) +{ +} + +int +vips_deviate( VipsImage *in, double *out, ... ) +{ + va_list ap; + int result; + + va_start( ap, out ); + result = vips_call_split( "deviate", ap, in, out ); + va_end( ap ); + + return( result ); +} diff --git a/libvips/arithmetic/im_deviate.c b/libvips/arithmetic/im_deviate.c deleted file mode 100644 index ae66e08a..00000000 --- a/libvips/arithmetic/im_deviate.c +++ /dev/null @@ -1,276 +0,0 @@ -/* im_deviate.c - * - * Copyright: 1990, J. Cupitt - * - * Author: J. Cupitt - * Written on: 02/08/1990 - * Modified on: - * 5/5/93 JC - * - now does partial images - * - less likely to overflow - * - adapted from im_avg - * 1/7/93 JC - * - adapted for partial v2 - * - ANSIfied - * 21/2/95 JC - * - modernised again - * 11/5/95 JC - * - oops! return( NULL ) in im_avg(), instead of return( -1 ) - * 20/6/95 JC - * - now returns double, not float - * 13/1/05 - * - use 64 bit arithmetic - * 8/12/06 - * - add liboil support - * 2/9/09 - * - gtk-doc comment - * - minor reformatting - * 4/9/09 - * - use im__wrapscan() - * 31/7/10 - * - remove liboil - */ - -/* - - This file is part of VIPS. - - VIPS is free software; you can redistribute it and/or modify - it under the terms of the GNU Lesser General Public License as published by - the Free Software Foundation; either version 2 of the License, or - (at your option) any later version. - - 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 Lesser General Public License for more details. - - You should have received a copy of the GNU Lesser General Public License - along with this program; if not, write to the Free Software - Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA - - */ - -/* - - These files are distributed with VIPS - http://www.vips.ecs.soton.ac.uk - - */ - -#ifdef HAVE_CONFIG_H -#include -#endif /*HAVE_CONFIG_H*/ -#include - -#include -#include -#include - -#include -#include - -typedef struct _Wrapscan { - IMAGE *in; - im_start_fn start; - im__wrapscan_fn scan; - im_stop_fn stop; - void *a; - void *b; -} Wrapscan; - -static void * -wrapscan_start( VipsImage *in, void *a, void *b ) -{ - Wrapscan *wrapscan = (Wrapscan *) a; - - return( wrapscan->start( in, wrapscan->a, wrapscan->b ) ); -} - -static int -wrapscan_stop( void *seq, void *a, void *b ) -{ - Wrapscan *wrapscan = (Wrapscan *) a; - - return( wrapscan->stop( seq, wrapscan->a, wrapscan->b ) ); -} - -static int -wrapscan_scan( REGION *reg, void *seq, void *a, void *b, gboolean *stop ) -{ - Wrapscan *wrapscan = (Wrapscan *) a; - Rect *r = ®->valid; - int lsk = IM_REGION_LSKIP( reg ); - - int y; - PEL *p; - - p = (PEL *) IM_REGION_ADDR( reg, r->left, r->top ); - for( y = 0; y < r->height; y++ ) { - if( wrapscan->scan( p, r->width, seq, - wrapscan->a, wrapscan->b ) ) - return( -1 ); - p += lsk; - } - - return( 0 ); -} - -/* Like vips_sink(), but the scan function works a line at a time, like - * im_wrap*(). Shared with im_min(), im_deviate() etc. - */ -int -im__wrapscan( VipsImage *in, - VipsStartFn start, im__wrapscan_fn scan, VipsStopFn stop, - void *a, void *b ) -{ - Wrapscan wrapscan; - - wrapscan.in = in; - wrapscan.start = start; - wrapscan.scan = scan; - wrapscan.stop = stop; - wrapscan.a = a; - wrapscan.b = b; - - return( vips_sink( in, - wrapscan_start, wrapscan_scan, wrapscan_stop, - &wrapscan, NULL ) ); -} - -/* Start function: allocate space for a pair of doubles in which we can - * accumulate the sum and the sum of squares. - */ -static void * -deviate_start( IMAGE *out, void *a, void *b ) -{ - double *ss2; - - if( !(ss2 = IM_ARRAY( NULL, 2, double )) ) - return( NULL ); - ss2[0] = 0.0; - ss2[1] = 0.0; - - return( (void *) ss2 ); -} - -/* Stop function. Add this little sum to the main sum. - */ -static int -deviate_stop( void *seq, void *a, void *b ) -{ - double *ss2 = (double *) seq; - double *global_ss2 = (double *) b; - - global_ss2[0] += ss2[0]; - global_ss2[1] += ss2[1]; - - im_free( seq ); - - return( 0 ); -} - -/* Sum pels in this section. - */ -#define LOOP( TYPE ) { \ - TYPE *p = (TYPE *) in; \ - \ - for( x = 0; x < sz; x++ ) { \ - TYPE v = p[x]; \ - \ - s += v; \ - s2 += (double) v * (double) v; \ - } \ -} - -/* Loop over region, accumulating a sum in *tmp. - */ -static int -deviate_scan( void *in, int n, void *seq, void *a, void *b ) -{ - const IMAGE *im = (IMAGE *) a; - const int sz = n * im->Bands; - - double *ss2 = (double *) seq; - - int x; - double s, s2; - - s = ss2[0]; - s2 = ss2[1]; - - /* Now generate code for all types. - */ - switch( im->BandFmt ) { - case IM_BANDFMT_UCHAR: LOOP( unsigned char ); break; - case IM_BANDFMT_CHAR: LOOP( signed char ); break; - case IM_BANDFMT_USHORT: LOOP( unsigned short ); break; - case IM_BANDFMT_SHORT: LOOP( signed short ); break; - case IM_BANDFMT_UINT: LOOP( unsigned int ); break; - case IM_BANDFMT_INT: LOOP( signed int ); break; - case IM_BANDFMT_FLOAT: LOOP( float ); break; - case IM_BANDFMT_DOUBLE: LOOP( double ); break; - - default: - g_assert( 0 ); - } - - ss2[0] = s; - ss2[1] = s2; - - return( 0 ); -} - -/** - * im_deviate: - * @in: input #IMAGE - * @out: output pixel standard deviation - * - * This operation finds the standard deviation of all pixels in @in. It - * operates on all bands of the input image: use im_stats() if you need - * to calculate an average for each band. - * - * Non-complex images only. - * - * See also: im_stats(), im_bandmean(), im_avg(), im_rank() - * - * Returns: 0 on success, -1 on error - */ -int -im_deviate( IMAGE *in, double *out ) -{ - double global_ss2[2]; - double s, s2; - gint64 vals; - - /* Check our args. - */ - if( im_pincheck( in ) || - im_check_uncoded( "im_deviate", in ) || - im_check_noncomplex( "im_deviate", in ) ) - return( -1 ); - - /* Loop over input, summing pixels. - */ - global_ss2[0] = 0.0; - global_ss2[1] = 0.0; - if( im__wrapscan( in, - deviate_start, deviate_scan, deviate_stop, in, global_ss2 ) ) - return( -1 ); - - /* - - NOTE: NR suggests a two-pass algorithm to minimise roundoff. - But that's too expensive for us :-( so do it the old one-pass - way. - - */ - - /* Calculate and return deviation. Add a fabs to stop sqrt(<=0). - */ - vals = (gint64) in->Xsize * (gint64) in->Ysize * (gint64) in->Bands; - s = global_ss2[0]; - s2 = global_ss2[1]; - *out = sqrt( fabs( s2 - (s * s / vals) ) / (vals - 1) ); - - return( 0 ); -} diff --git a/libvips/arithmetic/im_stats.c b/libvips/arithmetic/im_stats.c index cbf80150..26924ef4 100644 --- a/libvips/arithmetic/im_stats.c +++ b/libvips/arithmetic/im_stats.c @@ -68,6 +68,74 @@ #include #include +typedef struct _Wrapscan { + IMAGE *in; + im_start_fn start; + im__wrapscan_fn scan; + im_stop_fn stop; + void *a; + void *b; +} Wrapscan; + +static void * +wrapscan_start( VipsImage *in, void *a, void *b ) +{ + Wrapscan *wrapscan = (Wrapscan *) a; + + return( wrapscan->start( in, wrapscan->a, wrapscan->b ) ); +} + +static int +wrapscan_stop( void *seq, void *a, void *b ) +{ + Wrapscan *wrapscan = (Wrapscan *) a; + + return( wrapscan->stop( seq, wrapscan->a, wrapscan->b ) ); +} + +static int +wrapscan_scan( REGION *reg, void *seq, void *a, void *b, gboolean *stop ) +{ + Wrapscan *wrapscan = (Wrapscan *) a; + Rect *r = ®->valid; + int lsk = IM_REGION_LSKIP( reg ); + + int y; + PEL *p; + + p = (PEL *) IM_REGION_ADDR( reg, r->left, r->top ); + for( y = 0; y < r->height; y++ ) { + if( wrapscan->scan( p, r->width, seq, + wrapscan->a, wrapscan->b ) ) + return( -1 ); + p += lsk; + } + + return( 0 ); +} + +/* Like vips_sink(), but the scan function works a line at a time, like + * im_wrap*(). Shared with im_min(), im_deviate() etc. + */ +int +im__wrapscan( VipsImage *in, + VipsStartFn start, im__wrapscan_fn scan, VipsStopFn stop, + void *a, void *b ) +{ + Wrapscan wrapscan; + + wrapscan.in = in; + wrapscan.start = start; + wrapscan.scan = scan; + wrapscan.stop = stop; + wrapscan.a = a; + wrapscan.b = b; + + return( vips_sink( in, + wrapscan_start, wrapscan_scan, wrapscan_stop, + &wrapscan, NULL ) ); +} + /* Track min/max/sum/sum-of-squares for each thread during a scan. */ static void * diff --git a/libvips/deprecated/vips7compat.c b/libvips/deprecated/vips7compat.c index f1394f1c..f2a894f8 100644 --- a/libvips/deprecated/vips7compat.c +++ b/libvips/deprecated/vips7compat.c @@ -934,6 +934,12 @@ im_avg( IMAGE *in, double *out ) return( vips_avg( in, out, NULL ) ); } +int +im_deviate( IMAGE *in, double *out ) +{ + return( vips_deviate( in, out, NULL ) ); +} + int im_generate( VipsImage *im, im_start_fn start, im_generate_fn generate, im_stop_fn stop, void *a, void *b ) diff --git a/libvips/include/vips/arithmetic.h b/libvips/include/vips/arithmetic.h index 0698ddc0..3d91e16b 100644 --- a/libvips/include/vips/arithmetic.h +++ b/libvips/include/vips/arithmetic.h @@ -72,6 +72,8 @@ int vips_subtract( VipsImage *in1, VipsImage *in2, VipsImage **out, ... ) __attribute__((sentinel)); int vips_avg( VipsImage *in, double *out, ... ) __attribute__((sentinel)); +int vips_deviate( VipsImage *in, double *out, ... ) + __attribute__((sentinel)); int vips_min( VipsImage *in, double *out, ... ) __attribute__((sentinel)); int vips_max( VipsImage *in, double *out, ... ) @@ -99,7 +101,6 @@ DOUBLEMASK *im_measure_area( VipsImage *im, int h, int v, int *sel, int nsel, const char *name ); DOUBLEMASK *im_stats( VipsImage *in ); -int im_deviate( VipsImage *in, double *out ); int im_maxpos_avg( VipsImage *im, double *xpos, double *ypos, double *out ); int im_maxpos_vec( VipsImage *im, int *xpos, int *ypos, double *maxima, int n ); int im_minpos_vec( VipsImage *im, int *xpos, int *ypos, double *minima, int n ); diff --git a/libvips/include/vips/vips7compat.h b/libvips/include/vips/vips7compat.h index cf32198d..d68754cc 100644 --- a/libvips/include/vips/vips7compat.h +++ b/libvips/include/vips/vips7compat.h @@ -527,6 +527,7 @@ int im_minpos( VipsImage *in, int *xpos, int *ypos, double *out ); int im_max( VipsImage *in, double *out ); int im_maxpos( VipsImage *in, int *xpos, int *ypos, double *out ); int im_avg( VipsImage *in, double *out ); +int im_deviate( VipsImage *in, double *out ); int im_invert( VipsImage *in, VipsImage *out ); int im_lintra( double a, VipsImage *in, double b, VipsImage *out ); int im_lintra_vec( int n, double *a, VipsImage *in, double *b, VipsImage *out );