From 616bb9ab7f45558f2a2ca6dece9b84fd5b939035 Mon Sep 17 00:00:00 2001 From: John Cupitt Date: Mon, 7 Sep 2009 08:06:53 +0000 Subject: [PATCH] stuff --- ChangeLog | 2 + TODO | 38 +++---- libvips/arithmetic/im_avg.c | 26 ++--- libvips/arithmetic/im_linreg.c | 22 ++-- libvips/arithmetic/im_max.c | 188 +++++++++++++++----------------- libvips/arithmetic/im_maxpos.c | 36 +++--- libvips/arithmetic/im_min.c | 30 ++--- libvips/arithmetic/im_minpos.c | 36 +++--- libvips/include/vips/internal.h | 1 + 9 files changed, 186 insertions(+), 193 deletions(-) diff --git a/ChangeLog b/ChangeLog index 6470cb51..f6c3e50e 100644 --- a/ChangeLog +++ b/ChangeLog @@ -26,6 +26,8 @@ - im_remainder() has a float result for float types - im_measure() allows sel == NULL, meaning all patches - added "deprecated" package +- faster, simpler, better im_max() +- im_max() returns true modulus, not square of modulus, for complex images 25/3/09 started 7.18.0 - revised version numbers diff --git a/TODO b/TODO index 1c6e0aa2..f5a27593 100644 --- a/TODO +++ b/TODO @@ -1,4 +1,17 @@ -- im_stats() needs a rewrite ... scrap the stupid 42 thing +- new im_max() needs testing + + im_avg() needs to work for complex: return modulus + + use im_max() model to rewrite im_min(), im_stats() + +- im_minpos(), im_maxpos() need partialing + +- im_maxpos_vec needs gtkdoc + + + + + - 1-bit PNG read is broken? @@ -10,29 +23,6 @@ -- im__cast_and_call() no longer does - - out->Bbits = im_bits_of_fmt( out->BandFmt ); - - make sure we do this for each use of im__cast_and_call() - -- im__cast_and_call() no longer casts input to output type, it casts inputs - up to the smallest common type - - so we need to change all the buffer processors from - - T op T -> T - - to - - T1 op T1 -> T2 - - done add/sub/mul, doing div - - others? - - - - rename vipsCC in SWIG as pyvips? diff --git a/libvips/arithmetic/im_avg.c b/libvips/arithmetic/im_avg.c index 04629ab6..1b735645 100644 --- a/libvips/arithmetic/im_avg.c +++ b/libvips/arithmetic/im_avg.c @@ -100,6 +100,19 @@ stop_fn( void *seq, void *a, void *b ) return( 0 ); } +/* Sum pels in this section. + */ +#define LOOP( TYPE ) { \ + TYPE *p; \ + \ + for( y = to; y < bo; y++ ) { \ + p = (TYPE *) IM_REGION_ADDR( reg, le, y ); \ + \ + for( x = 0; x < sz; x++ ) \ + sum += *p++; \ + } \ +} + /* Loop over region, accumulating a sum in *tmp. */ static int @@ -115,19 +128,6 @@ scan_fn( REGION *reg, void *seq, void *a, void *b ) double sum = 0.0; int x, y; -/* Sum pels in this section. - */ -#define LOOP( TYPE ) { \ - TYPE *p; \ - \ - for( y = to; y < bo; y++ ) { \ - p = (TYPE *) IM_REGION_ADDR( reg, le, y ); \ - \ - for( x = 0; x < sz; x++ ) \ - sum += *p++; \ - } \ -} - /* Now generate code for all types. */ switch( im->BandFmt ) { diff --git a/libvips/arithmetic/im_linreg.c b/libvips/arithmetic/im_linreg.c index bf6aba26..95bcc5f8 100644 --- a/libvips/arithmetic/im_linreg.c +++ b/libvips/arithmetic/im_linreg.c @@ -1,12 +1,4 @@ -/* @(#) Function to find perform pixelwise linear regression on an array of - * @(#) single band images. - * @(#) - * @(#) int im_linreg( - * @(#) IMAGE **ins, - * @(#) IMAGE *out, - * @(#) double *xs - * @(#) ); - * @(#) +/* im_linreg.c * * Copyright: 2006, The Nottingham Trent University * @@ -16,6 +8,8 @@ * * 12/5/09 * - make x_anal() static, fix some signed/unsigned warnings + * 3/9/09 + * - gtkdoc comment */ /* @@ -148,6 +142,16 @@ SKIP_ALL_DECL( double ); /** EXPORTED FUNCTION DEFINITION **/ +/** im_linreg: + * @ins: NULL-terminated array of input images + * @out: results of analysis + * @xs: X position of each image (pixel value is Y) + * + * Function to find perform pixelwise linear regression on an array of + * single band images. The output is a seven-band douuble image + * + * TODO: figure out how this works and fix up these docs! + */ int im_linreg( IMAGE **ins, IMAGE *out, double *xs ){ #define FUNCTION_NAME "im_linreg" int n; diff --git a/libvips/arithmetic/im_max.c b/libvips/arithmetic/im_max.c index 2713432d..81e8b8bd 100644 --- a/libvips/arithmetic/im_max.c +++ b/libvips/arithmetic/im_max.c @@ -18,6 +18,9 @@ * - random wrong result for >1 thread :-( (thanks Joe) * 15/10/07 * - oh, heh, seq->inf was not being set correctly, not that it mattered + * 4/9/09 + * - rewrite with im__value(), much simpler and fixes a race condition + * - gtkdoc comment */ /* @@ -58,7 +61,6 @@ #include #include #include -#include #include @@ -66,39 +68,15 @@ #include #endif /*WITH_DMALLOC*/ -/* Per-call state. - */ -typedef struct _MaxInfo { - /* Parameters. - */ - IMAGE *in; - double *out; - - /* Global max so far. - */ - double value; - int valid; /* zero means value is unset */ -} MaxInfo; - -/* Per thread state. - */ -typedef struct _Seq { - MaxInfo *inf; - - double value; - int valid; /* zero means value is unset */ -} Seq; - /* New sequence value. */ static void * max_start( IMAGE *in, void *a, void *b ) { - MaxInfo *inf = (MaxInfo *) a; - Seq *seq = IM_NEW( NULL, Seq ); + double *value = (double *) a; + double *seq = IM_NEW( NULL, double ); - seq->inf = inf; - seq->valid = 0; + *seq = *value; return( (void *) seq ); } @@ -108,47 +86,19 @@ max_start( IMAGE *in, void *a, void *b ) static int max_stop( void *vseq, void *a, void *b ) { - Seq *seq = (Seq *) vseq; - MaxInfo *inf = (MaxInfo *) a; + double *seq = (double *) vseq; + double *value = (double *) a; - if( seq->valid ) { - if( !inf->valid ) - /* Just copy. - */ - inf->value = seq->value; - else - /* Merge. - */ - inf->value = IM_MAX( inf->value, seq->value ); - - inf->valid = 1; - } + /* Merge. + */ + *value = IM_MAX( *value, *seq ); im_free( seq ); return( 0 ); } -/* Loop over region, adding to seq. - */ -static int -max_scan( REGION *reg, void *vseq, void *a, void *b ) -{ - Seq *seq = (Seq *) vseq; - Rect *r = ®->valid; - IMAGE *im = reg->im; - int le = r->left; - int to = r->top; - int bo = IM_RECT_BOTTOM(r); - int nel = IM_REGION_N_ELEMENTS( reg ); - - int x, y; - - double m; - -#define loop(TYPE) { \ - m = *((TYPE *) IM_REGION_ADDR( reg, le, to )); \ - \ +#define LOOP( TYPE ) { \ for( y = to; y < bo; y++ ) { \ TYPE *p = (TYPE *) IM_REGION_ADDR( reg, le, y ); \ \ @@ -161,22 +111,17 @@ max_scan( REGION *reg, void *vseq, void *a, void *b ) } \ } -#define complex_loop(TYPE) { \ - TYPE *p = (TYPE *) IM_REGION_ADDR( reg, le, to ); \ - double real = p[0]; \ - double imag = p[1]; \ - \ - m = real * real + imag * imag; \ - \ +#define CLOOP( TYPE ) { \ for( y = to; y < bo; y++ ) { \ TYPE *p = (TYPE *) IM_REGION_ADDR( reg, le, y ); \ \ - for( x = 0; x < nel * 2; x += 2 ) { \ - double mod; \ + for( x = 0; x < nel; x++ ) { \ + double mod, re, im; \ \ - real = p[x]; \ - imag = p[x + 1]; \ - mod = real * real + imag * imag; \ + re = p[0]; \ + im = p[1]; \ + p += 2; \ + mod = re * re + im * im; \ \ if( mod > m ) \ m = mod; \ @@ -184,47 +129,76 @@ max_scan( REGION *reg, void *vseq, void *a, void *b ) } \ } +/* Loop over region, adding to seq. + */ +static int +max_scan( REGION *reg, void *vseq, void *a, void *b ) +{ + double *seq = (double *) vseq; + Rect *r = ®->valid; + IMAGE *im = reg->im; + int le = r->left; + int to = r->top; + int bo = IM_RECT_BOTTOM(r); + int nel = IM_REGION_N_ELEMENTS( reg ); + + int x, y; + double m; + + m = *seq; + 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; - case IM_BANDFMT_COMPLEX: complex_loop( float ); break; - case IM_BANDFMT_DPCOMPLEX: complex_loop( double ); break; + 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; + case IM_BANDFMT_COMPLEX: CLOOP( float ); break; + case IM_BANDFMT_DPCOMPLEX: CLOOP( double ); break; default: - assert( 0 ); + g_assert( 0 ); } - if( seq->valid ) { - seq->value = IM_MAX( seq->value, m ); - } - else { - seq->value = m; - seq->valid = 1; - } + *seq = m; #ifdef DEBUG printf( "im_max: left = %d, top = %d, width = %d, height = %d\n", r->left, r->top, r->width, r->height ); - printf( " (max = %g)\n", seq->value ); + printf( " (max = %g)\n", *seq ); #endif /*DEBUG*/ return( 0 ); } +/* Get the value of pixel (0, 0). Use this to init the min/max value for + * threads. Shared with im_min(), im_stats() etc. This will return mod for + * complex. + */ +int +im__value( IMAGE *im, double *value ) +{ + IMAGE *t; + + if( !(t = im_open_local( im, "im__value", "p" )) || + im_extract_area( im, t, 0, 0, 1, 1 ) || + im_avg( t, value ) ) + return( -1 ); + + return( 0 ); +} + /** * im_max: * @in: input #IMAGE * @out: output double * * Finds the the maximum value of image #in and returns it at the - * location pointed by out. If input is complex, the max square amplitude - * (re*re+im*im) is returned. im_max() finds the maximum of all bands: if you + * location pointed by out. If input is complex, the max modulus + * is returned. im_max() finds the maximum of all bands: if you * want to find the maximum of each band separately, use im_stats(). * * See also: im_maxpos(), im_min(), im_stats(). @@ -233,21 +207,29 @@ max_scan( REGION *reg, void *vseq, void *a, void *b ) */ int im_max( IMAGE *in, double *out ) -{ - MaxInfo inf; - - inf.in = in; - inf.out = out; - inf.valid = 0; +{ + double value; if( im_pincheck( in ) || im_check_uncoded( "im_max", in ) ) return( -1 ); - if( im_iterate( in, max_start, max_scan, max_stop, &inf, NULL ) ) + if( im__value( in, &value ) ) + return( -1 ); + /* We use square mod for scanning, for speed. + */ + if( im_iscomplex( in ) ) + value *= value; + + if( im_iterate( in, max_start, max_scan, max_stop, &value, NULL ) ) return( -1 ); - *out = inf.value; + /* Back to modulus. + */ + if( im_iscomplex( in ) ) + value = sqrt( value ); + + *out = value; return( 0 ); } diff --git a/libvips/arithmetic/im_maxpos.c b/libvips/arithmetic/im_maxpos.c index 3a25214f..dec838e8 100644 --- a/libvips/arithmetic/im_maxpos.c +++ b/libvips/arithmetic/im_maxpos.c @@ -1,15 +1,4 @@ -/* @(#) Function to find the maximum of an image. Works for any - * @(#) image type. Returns a double and the location of max. - * @(#) - * @(#) Function im_maxpos() assumes that input - * @(#) is either memory mapped or in a buffer. - * @(#) - * @(#) int im_maxpos(in, xpos, ypos, max) - * @(#) IMAGE *in; - * @(#) int *xpos, *ypos; - * @(#) double *max; - * @(#) - * @(#) Returns 0 on success and -1 on error +/* im_maxpos.c * * Copyright: 1990, J. Cupitt * @@ -21,6 +10,8 @@ * - im_incheck() call added * 20/6/95 JC * - now returns double for value, like im_max() + * 4/9/09 + * - gtkdoc comment */ /* @@ -55,7 +46,6 @@ #include #include -#include #include @@ -77,8 +67,22 @@ case IM_BANDFMT_COMPLEX: loopcmplx(float); break; \ case IM_BANDFMT_DPCOMPLEX: loopcmplx(double); break; -/* Find the position of the maximum of an image. Take any format, returns a - * float and its position. +/** + * im_maxpos: + * @in: image to search + * @xpos: returned x position of maximum + * @ypos: returned y position of maximum + * @out: returned pixel value at that position + * + * Function to find the maximum of an image. Works for any + * image type. Returns a double and the location of max. + * + * This is not a PIO operation! It may use a lot of memory and take a while. + * Needs a rewrite. + * + * See also: im_minpos(), im_min(), im_stats(), im_maxpos_avg(). + * + * Returns: 0 on success, -1 on error */ int im_maxpos( IMAGE *in, int *xpos, int *ypos, double *out ) @@ -135,7 +139,7 @@ im_maxpos( IMAGE *in, int *xpos, int *ypos, double *out ) switch( in->BandFmt ) { im_for_all_types(); default: { - assert( 0 ); + g_assert( 0 ); return( -1 ); } } diff --git a/libvips/arithmetic/im_min.c b/libvips/arithmetic/im_min.c index 55a76c70..75cb095f 100644 --- a/libvips/arithmetic/im_min.c +++ b/libvips/arithmetic/im_min.c @@ -1,11 +1,4 @@ -/* @(#) Function to find the minimum of an image. Works for any - * @(#) image type. Returns a double. - * @(#) - * @(#) int im_min(in, min) - * @(#) IMAGE *in; - * @(#) double *min; - * @(#) - * @(#) Returns 0 on success and -1 on error +/* im_min.c * * Copyright: 1990, J. Cupitt * @@ -221,6 +214,20 @@ scan_fn( REGION *reg, void *vseq, void *a, void *b ) return( 0 ); } +/** + * im_min: + * @in: input #IMAGE + * @out: output double + * + * Finds the the minimum value of image @in and returns it at the + * location pointed by @out. If input is complex, the min square amplitude + * (re*re+im*im) is returned. im_min() finds the maximum of all bands: if you + * want to find the maximum of each band separately, use im_stats(). + * + * See also: im_minpos(), im_max(), im_stats(). + * + * Returns: 0 on success, -1 on error + */ int im_min( IMAGE *in, double *out ) { @@ -230,12 +237,9 @@ im_min( IMAGE *in, double *out ) inf.out = out; inf.valid = 0; - if( im_pincheck( in ) ) + if( im_pincheck( in ) || + im_check_uncoded( "im_min", in ) ) return( -1 ); - if( in->Coding != IM_CODING_NONE ) { - im_error( "im_min", "%s", _( "not uncoded" ) ); - return( -1 ); - } if( im_iterate( in, start_fn, scan_fn, stop_fn, &inf, NULL ) ) return( -1 ); diff --git a/libvips/arithmetic/im_minpos.c b/libvips/arithmetic/im_minpos.c index 6861d69a..72ccc3e5 100644 --- a/libvips/arithmetic/im_minpos.c +++ b/libvips/arithmetic/im_minpos.c @@ -1,15 +1,4 @@ -/* @(#) Function to find the minimim of an image. Works for any - * @(#) image type. Returns a double and the location of min - * @(#) - * @(#) Function im_minpos() assumes that input - * @(#) is either memory mapped or in a buffer. - * @(#) - * @(#) int im_minpos(in, xpos, ypos, out) - * @(#) IMAGE *in; - * @(#) int *xpos, *ypos; - * @(#) double *out; - * @(#) - * @(#) Returns 0 on success and -1 on error +/* im_minpos.c * * Copyright: 1990, J. Cupitt * @@ -22,6 +11,8 @@ * - im_incheck() added * 20/6/95 JC * - now returns double for value, like im_max() + * 4/9/09 + * - gtkdoc comment */ /* @@ -56,7 +47,6 @@ #include #include -#include #include @@ -64,7 +54,23 @@ #include #endif /*WITH_DMALLOC*/ -/* Find the minimum of an image. Take any format, returns a double. */ +/** + * im_minpos: + * @in: image to search + * @xpos: returned x position of minimum + * @ypos: returned y position of minimum + * @out: returned pixel value at that position + * + * Function to find the minimum of an image. Works for any + * image type. Returns a double and the location of min. + * + * This is not a PIO operation! It may use a lot of memory and take a while. + * Needs a rewrite. + * + * See also: im_maxpos(), im_max(), im_stats(), im_maxpos_avg(). + * + * Returns: 0 on success, -1 on error + */ int im_minpos( IMAGE *in, int *xpos, int *ypos, double *out ) { @@ -131,7 +137,7 @@ im_minpos( IMAGE *in, int *xpos, int *ypos, double *out ) case IM_BANDFMT_DPCOMPLEX: loopcmplx(double); break; default: - assert( 0 ); + g_assert( 0 ); } /* Take out bands on x. diff --git a/libvips/include/vips/internal.h b/libvips/include/vips/internal.h index 1321de1d..44041e34 100644 --- a/libvips/include/vips/internal.h +++ b/libvips/include/vips/internal.h @@ -135,6 +135,7 @@ int im__cast_and_call( IMAGE *in1, IMAGE *in2, IMAGE *out, im_wrapmany_fn fn, void *a ); VipsBandFmt im__format_common( IMAGE *in1, IMAGE *in2 ); int im__math( const char *name, IMAGE *in, IMAGE *out, im_wrapone_fn gen ); +int im__value( IMAGE *im, double *value ); int im__test_kill( IMAGE *im ); void *im__mmap( int fd, int writeable, size_t length, gint64 offset );