467 lines
11 KiB
C
467 lines
11 KiB
C
/* stats.c ... many image stats in a single pass
|
|
*
|
|
(C) Kirk Martinez 1993
|
|
23/4/93 J.Cupitt
|
|
- adapted to partial images
|
|
- special struct abandoned; now returns DOUBLEMASK.
|
|
1/7/93 JC
|
|
- adapted for partial v2
|
|
- ANSIfied
|
|
27/7/93 JC
|
|
- init of mask changed to use actual values, not IM_MAXDOUBLE and
|
|
(-IM_MAXDOUBLE). These caused problems when assigned to floats.
|
|
funny business with offset==42, yuk!
|
|
31/8/93 JC
|
|
- forgot to init global max/min properly! sorry.
|
|
21/6/95 JC
|
|
- still did not init max and min correctly --- now fixed for good
|
|
|
|
* 13/1/05
|
|
* - use 64 bit arithmetic
|
|
* 1/9/09
|
|
* - argh nope min/max was broken again for >1CPU in short pipelines on
|
|
* some architectures
|
|
* 7/9/09
|
|
* - rework based on new im__wrapscan() / im_max() ideas for a proper fix
|
|
* - gtkdoc comment
|
|
* 7/11/11
|
|
* - redone as a class
|
|
* - track maxpos / minpos too
|
|
*/
|
|
|
|
/*
|
|
|
|
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., 51 Franklin Street, Fifth Floor, Boston, MA
|
|
02110-1301 USA
|
|
|
|
*/
|
|
|
|
/*
|
|
|
|
These files are distributed with VIPS - http://www.vips.ecs.soton.ac.uk
|
|
|
|
*/
|
|
|
|
/*
|
|
#define VIPS_DEBUG
|
|
*/
|
|
|
|
#ifdef HAVE_CONFIG_H
|
|
#include <config.h>
|
|
#endif /*HAVE_CONFIG_H*/
|
|
#include <glib/gi18n-lib.h>
|
|
|
|
#include <stdio.h>
|
|
#include <stdlib.h>
|
|
#include <math.h>
|
|
|
|
#include <vips/vips.h>
|
|
#include <vips/internal.h>
|
|
|
|
#include "statistic.h"
|
|
|
|
typedef struct _VipsStats {
|
|
VipsStatistic parent_instance;
|
|
|
|
VipsImage *out;
|
|
|
|
gboolean set; /* FALSE means no value yet */
|
|
} VipsStats;
|
|
|
|
typedef VipsStatisticClass VipsStatsClass;
|
|
|
|
G_DEFINE_TYPE( VipsStats, vips_stats, VIPS_TYPE_STATISTIC );
|
|
|
|
/* Names for our columns.
|
|
*/
|
|
enum {
|
|
COL_MIN = 0,
|
|
COL_MAX = 1,
|
|
COL_SUM = 2,
|
|
COL_SUM2 = 3,
|
|
COL_AVG = 4,
|
|
COL_SD = 5,
|
|
COL_XMIN = 6,
|
|
COL_YMIN = 7,
|
|
COL_XMAX = 8,
|
|
COL_YMAX = 9,
|
|
COL_LAST = 10
|
|
};
|
|
|
|
static int
|
|
vips_stats_build( VipsObject *object )
|
|
{
|
|
VipsObjectClass *class = VIPS_OBJECT_GET_CLASS( object );
|
|
VipsStatistic *statistic = VIPS_STATISTIC( object );
|
|
VipsStats *stats = (VipsStats *) object;
|
|
|
|
gint64 vals, pels;
|
|
double *row0, *row;
|
|
int b, y, i;
|
|
|
|
if( vips_object_argument_isset( object, "in" ) ) {
|
|
int bands = vips_image_get_bands( statistic->in );
|
|
|
|
if( vips_check_noncomplex( class->nickname, statistic->in ) )
|
|
return( -1 );
|
|
|
|
g_object_set( object,
|
|
"out", vips_image_new_matrix( COL_LAST, bands + 1 ),
|
|
NULL );
|
|
}
|
|
|
|
if( VIPS_OBJECT_CLASS( vips_stats_parent_class )->build( object ) )
|
|
return( -1 );
|
|
|
|
pels = (gint64) vips_image_get_width( statistic->in ) *
|
|
vips_image_get_height( statistic->in );
|
|
vals = pels * vips_image_get_bands( statistic->in );
|
|
|
|
row0 = VIPS_MATRIX( stats->out, 0, 0 );
|
|
row = VIPS_MATRIX( stats->out, 0, 1 );
|
|
for( i = 0; i < COL_LAST; i++ )
|
|
row0[i] = row[i];
|
|
|
|
for( b = 1; b < vips_image_get_bands( statistic->in ); b++ ) {
|
|
row = VIPS_MATRIX( stats->out, 0, b + 1 );
|
|
|
|
if( row[COL_MIN] < row0[COL_MIN] ) {
|
|
row0[COL_MIN] = row[COL_MIN];
|
|
row0[COL_XMIN] = row[COL_XMIN];
|
|
row0[COL_YMIN] = row[COL_YMIN];
|
|
}
|
|
|
|
if( row[COL_MAX] > row0[COL_MAX] ) {
|
|
row0[COL_MAX] = row[COL_MAX];
|
|
row0[COL_XMAX] = row[COL_XMAX];
|
|
row0[COL_YMAX] = row[COL_YMAX];
|
|
}
|
|
|
|
row0[COL_SUM] += row[COL_SUM];
|
|
row0[COL_SUM2] += row[COL_SUM2];
|
|
}
|
|
|
|
for( y = 1; y < vips_image_get_height( stats->out ); y++ ) {
|
|
double *row = VIPS_MATRIX( stats->out, 0, y );
|
|
|
|
row[COL_AVG] = row[COL_SUM] / pels;
|
|
row[COL_SD] = sqrt( VIPS_FABS( row[COL_SUM2] -
|
|
(row[COL_SUM] * row[COL_SUM] / pels) ) / (pels - 1) );
|
|
}
|
|
|
|
row0[COL_AVG] = row0[COL_SUM] / vals;
|
|
row0[COL_SD] = sqrt( VIPS_FABS( row0[COL_SUM2] -
|
|
(row0[COL_SUM] * row0[COL_SUM] / vals) ) / (vals - 1) );
|
|
|
|
return( 0 );
|
|
}
|
|
|
|
/* Stop function. Add these little stats to the main set of stats.
|
|
*/
|
|
static int
|
|
vips_stats_stop( VipsStatistic *statistic, void *seq )
|
|
{
|
|
int bands = vips_image_get_bands( statistic->in );
|
|
VipsStats *global = (VipsStats *) statistic;
|
|
VipsStats *local = (VipsStats *) seq;
|
|
|
|
int b;
|
|
|
|
if( local->set && !global->set ) {
|
|
for( b = 0; b < bands; b++ ) {
|
|
double *p = VIPS_MATRIX( local->out, 0, b + 1 );
|
|
double *q = VIPS_MATRIX( global->out, 0, b + 1 );
|
|
|
|
int i;
|
|
|
|
for( i = 0; i < COL_LAST; i++ )
|
|
q[i] = p[i];
|
|
}
|
|
|
|
global->set = TRUE;
|
|
}
|
|
else if( local->set && global->set ) {
|
|
for( b = 0; b < bands; b++ ) {
|
|
double *p = VIPS_MATRIX( local->out, 0, b + 1 );
|
|
double *q = VIPS_MATRIX( global->out, 0, b + 1 );
|
|
|
|
if( p[COL_MIN] < q[COL_MIN] ) {
|
|
q[COL_MIN] = p[COL_MIN];
|
|
q[COL_XMIN] = p[COL_XMIN];
|
|
q[COL_YMIN] = p[COL_YMIN];
|
|
}
|
|
|
|
if( p[COL_MAX] > q[COL_MAX] ) {
|
|
q[COL_MAX] = p[COL_MAX];
|
|
q[COL_XMAX] = p[COL_XMAX];
|
|
q[COL_YMAX] = p[COL_YMAX];
|
|
}
|
|
|
|
q[COL_SUM] += p[COL_SUM];
|
|
q[COL_SUM2] += p[COL_SUM2];
|
|
}
|
|
}
|
|
|
|
VIPS_FREEF( g_object_unref, local->out );
|
|
VIPS_FREEF( g_free, seq );
|
|
|
|
return( 0 );
|
|
}
|
|
|
|
/* Start function: make a dummy local stats for the private use of this thread.
|
|
*/
|
|
static void *
|
|
vips_stats_start( VipsStatistic *statistic )
|
|
{
|
|
int bands = vips_image_get_bands( statistic->in );
|
|
|
|
VipsStats *stats;
|
|
|
|
stats = g_new( VipsStats, 1 );
|
|
if( !(stats->out = vips_image_new_matrix( COL_LAST, bands + 1 )) ) {
|
|
g_free( stats );
|
|
return( NULL );
|
|
}
|
|
stats->set = FALSE;
|
|
|
|
return( (void *) stats );
|
|
}
|
|
|
|
/* We scan lines bands times to avoid repeating band loops.
|
|
* Use temp variables of same type for min/max for faster comparisons.
|
|
*/
|
|
#define LOOP( TYPE ) { \
|
|
for( b = 0; b < bands; b++ ) { \
|
|
TYPE *p = ((TYPE *) in) + b; \
|
|
double *q = VIPS_MATRIX( 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; \
|
|
}
|
|
|
|
/* 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 = VIPS_MATRIX( 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
|
|
vips_stats_scan( VipsStatistic *statistic, void *seq,
|
|
int x, int y, void *in, int n )
|
|
{
|
|
const int bands = vips_image_get_bands( statistic->in );
|
|
VipsStats *local = (VipsStats *) seq;
|
|
|
|
int b, i;
|
|
|
|
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_not_reached();
|
|
}
|
|
|
|
return( 0 );
|
|
}
|
|
|
|
static void
|
|
vips_stats_class_init( VipsStatsClass *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 = "stats";
|
|
object_class->description = _( "find many image stats" );
|
|
object_class->build = vips_stats_build;
|
|
|
|
sclass->start = vips_stats_start;
|
|
sclass->scan = vips_stats_scan;
|
|
sclass->stop = vips_stats_stop;
|
|
|
|
VIPS_ARG_IMAGE( class, "out", 100,
|
|
_( "Output" ),
|
|
_( "Output array of statistics" ),
|
|
VIPS_ARGUMENT_REQUIRED_OUTPUT,
|
|
G_STRUCT_OFFSET( VipsStats, out ) );
|
|
}
|
|
|
|
static void
|
|
vips_stats_init( VipsStats *stats )
|
|
{
|
|
}
|
|
|
|
/**
|
|
* vips_stats: (method)
|
|
* @in: image to scan
|
|
* @out: (out): image of statistics
|
|
* @...: %NULL-terminated list of optional named arguments
|
|
*
|
|
* Find many image statistics in a single pass through the data. @out is a
|
|
* one-band #VIPS_FORMAT_DOUBLE image of at least 10 columns by n + 1
|
|
* (where n is number of bands in image @in)
|
|
* rows. Columns are statistics, and are, in order: minimum, maximum, sum,
|
|
* sum of squares, mean, standard deviation, x coordinate of minimum, y
|
|
* coordinate of minimum, x coordinate of maximum, y coordinate of maximum.
|
|
* Later versions of vips_stats() may add more columns.
|
|
*
|
|
* Row 0 has statistics for all
|
|
* bands together, row 1 has stats for band 1, and so on.
|
|
*
|
|
* If there is more than one maxima or minima, one of them will be chosen at
|
|
* random.
|
|
*
|
|
* See also: vips_avg(), vips_min().
|
|
*
|
|
* Returns: 0 on success, -1 on error
|
|
*/
|
|
int
|
|
vips_stats( VipsImage *in, VipsImage **out, ... )
|
|
{
|
|
va_list ap;
|
|
int result;
|
|
|
|
va_start( ap, out );
|
|
result = vips_call_split( "stats", ap, in, out );
|
|
va_end( ap );
|
|
|
|
return( result );
|
|
}
|