| /* im_gradcor |
| * |
| * Copyright: 2007 Nottingham Trent University |
| * |
| * Author: Tom Vajzovic |
| * Written on: 2007-06-07 |
| * 3/2/10 |
| * - gtkdoc |
| * - cleanups |
| */ |
| |
| /* |
| |
| 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 |
| |
| */ |
| |
| |
| /** HEADERS **/ |
| |
| #ifdef HAVE_CONFIG_H |
| #include <config.h> |
| #endif /*HAVE_CONFIG_H*/ |
| #include <vips/intl.h> |
| |
| #include <stdlib.h> |
| |
| #include <vips/vips.h> |
| |
| #ifdef WITH_DMALLOC |
| #include <dmalloc.h> |
| #endif /*WITH_DMALLOC*/ |
| |
| |
| /** LOCAL TYPES **/ |
| |
| typedef struct { |
| REGION *reg; |
| int *region_xgrad; |
| int *region_ygrad; |
| size_t region_xgrad_area; |
| size_t region_ygrad_area; |
| } |
| gradcor_seq_t; |
| |
| |
| /** LOCAL FUNCTION DECLARATIONS **/ |
| |
| static void *gradcor_start( IMAGE *out, void *vptr_large, void *unrequired ); |
| static int gradcor_stop( void *vptr_seq, void *unrequired, void *unreq2 ); |
| static int gradcor_gen( REGION *to_make, void *vptr_seq, void *unrequired, void *vptr_grads ); |
| |
| #define XGRAD_GEN_DECLARATION( TYPE ) static int xgrad_gen_ ## TYPE( REGION *to_make, void *vptr_make_from, void *unrequired, void *unreq2 ) |
| #define YGRAD_GEN_DECLARATION( TYPE ) static int ygrad_gen_ ## TYPE( REGION *to_make, void *vptr_make_from, void *unrequired, void *unreq2 ) |
| |
| XGRAD_GEN_DECLARATION( guint8 ); |
| YGRAD_GEN_DECLARATION( guint8 ); |
| XGRAD_GEN_DECLARATION( gint8 ); |
| YGRAD_GEN_DECLARATION( gint8 ); |
| XGRAD_GEN_DECLARATION( guint16 ); |
| YGRAD_GEN_DECLARATION( guint16 ); |
| XGRAD_GEN_DECLARATION( gint16 ); |
| YGRAD_GEN_DECLARATION( gint16 ); |
| XGRAD_GEN_DECLARATION( guint32 ); |
| YGRAD_GEN_DECLARATION( guint32 ); |
| XGRAD_GEN_DECLARATION( gint32 ); |
| YGRAD_GEN_DECLARATION( gint32 ); |
| #if 0 |
| XGRAD_GEN_DECLARATION( float ); |
| YGRAD_GEN_DECLARATION( float ); |
| XGRAD_GEN_DECLARATION( double ); |
| YGRAD_GEN_DECLARATION( double ); |
| #endif |
| |
| |
| /** EXPORTED FUNCTION DEFINITIONS **/ |
| |
| int im_gradcor_raw( IMAGE *large, IMAGE *small, IMAGE *out ){ |
| #define FUNCTION_NAME "im_gradcor_raw" |
| |
| if( im_piocheck( large, out ) || im_pincheck( small ) ) |
| return -1; |
| |
| if( im_check_uncoded( "im_gradcor", large ) || |
| im_check_mono( "im_gradcor", large ) || |
| im_check_uncoded( "im_gradcor", small ) || |
| im_check_mono( "im_gradcor", small ) || |
| im_check_format_same( "im_gradcor", large, small ) || |
| im_check_int( "im_gradcor", large ) ) |
| return( -1 ); |
| |
| if( large-> Xsize < small-> Xsize || large-> Ysize < small-> Ysize ){ |
| im_error( FUNCTION_NAME, "second image must be smaller than first" ); |
| return -1; |
| } |
| if( im_cp_desc( out, large ) ) |
| return -1; |
| |
| out-> Xsize= 1 + large-> Xsize - small-> Xsize; |
| out-> Ysize= 1 + large-> Ysize - small-> Ysize; |
| out-> BandFmt= IM_BANDFMT_FLOAT; |
| |
| if( im_demand_hint( out, IM_FATSTRIP, large, NULL ) ) |
| return -1; |
| |
| { |
| IMAGE *xgrad= im_open_local( out, FUNCTION_NAME ": xgrad", "t" ); |
| IMAGE *ygrad= im_open_local( out, FUNCTION_NAME ": ygrad", "t" ); |
| IMAGE **grads= im_allocate_input_array( out, xgrad, ygrad, NULL ); |
| |
| return ! xgrad || ! ygrad || ! grads |
| || im_grad_x( small, xgrad ) |
| || im_grad_y( small, ygrad ) |
| || im_generate( out, gradcor_start, gradcor_gen, gradcor_stop, (void*) large, (void*) grads ); |
| } |
| #undef FUNCTION_NAME |
| } |
| |
| /** |
| * im_gradcor: |
| * @in: input image |
| * @ref: reference image |
| * @out: output image |
| * |
| * Calculate a correlation surface. |
| * |
| * @ref is placed at every position in @in and a correlation coefficient |
| * calculated. One-band, integer images only. @in and @ref must have the |
| * same #VipsBandFmt. The output |
| * image is always %IM_BANDFMT_FLOAT. @ref must be smaller than @in. The output |
| * image is the same size as the input. |
| * |
| * The method takes the gradient images of the two images then takes the |
| * dot-product correlation of the two vector images. |
| * The vector expression of this method is my (tcv) own creation. It is |
| * equivalent to the complex-number method of: |
| * |
| * ARGYRIOU, V. et al. 2003. Estimation of sub-pixel motion using |
| * gradient cross correlation. Electronics Letters, 39 (13). |
| * |
| * See also: im_spcor(). |
| * |
| * Returns: 0 on success, -1 on error |
| */ |
| int |
| im_gradcor( IMAGE *in, IMAGE *ref, IMAGE *out ) |
| { |
| #define FUNCTION_NAME "im_gradcor" |
| IMAGE *t1 = im_open_local( out, FUNCTION_NAME " intermediate", "p" ); |
| |
| if( !t1 || |
| im_embed( in, t1, 1, |
| ref->Xsize / 2, ref->Ysize / 2, |
| in->Xsize + ref->Xsize - 1, |
| in->Ysize + ref->Ysize - 1 ) || |
| im_gradcor_raw( t1, ref, out ) ) |
| return( -1 ); |
| |
| out->Xoffset = 0; |
| out->Yoffset = 0; |
| |
| return( 0 ); |
| #undef FUNCTION_NAME |
| } |
| |
| /** |
| * im_grad_x: |
| * @in: input image |
| * @out: output image |
| * |
| * Find horizontal differences between adjacent pixels. |
| * |
| * Generates an image where the value of each pixel is the difference between |
| * it and the pixel to its right. The output has the same height as the input |
| * and one pixel less width. One-band integer formats only. The result is |
| * always %IM_BANDFMT_INT. |
| * |
| * This operation is much faster than (though equivalent to) im_conv() with the |
| * mask [[-1, 1]]. |
| * |
| * See also: im_grad_y(), im_conv(). |
| * |
| * Returns: 0 on success, -1 on error |
| */ |
| int im_grad_x( IMAGE *in, IMAGE *out ){ |
| #define FUNCTION_NAME "im_grad_x" |
| |
| if( im_piocheck( in, out ) ) |
| return -1; |
| |
| if( im_check_uncoded( "im_grad_x", in ) || |
| im_check_mono( "im_grad_x", in ) || |
| im_check_int( "im_grad_x", in ) ) |
| return( -1 ); |
| if( im_cp_desc( out, in ) ) |
| return -1; |
| |
| -- out-> Xsize; |
| out-> BandFmt= IM_BANDFMT_INT; /* do not change without updating im_gradcor() */ |
| |
| if( im_demand_hint( out, IM_THINSTRIP, in, NULL ) ) |
| return -1; |
| |
| #define RETURN_GENERATE( TYPE ) return im_generate( out, im_start_one, xgrad_gen_ ## TYPE, im_stop_one, (void*) in, NULL ) |
| |
| switch( in-> BandFmt ){ |
| |
| case IM_BANDFMT_UCHAR: |
| RETURN_GENERATE( guint8 ); |
| |
| case IM_BANDFMT_CHAR: |
| RETURN_GENERATE( gint8 ); |
| |
| case IM_BANDFMT_USHORT: |
| RETURN_GENERATE( guint16 ); |
| |
| case IM_BANDFMT_SHORT: |
| RETURN_GENERATE( gint16 ); |
| |
| case IM_BANDFMT_UINT: |
| RETURN_GENERATE( guint32 ); |
| |
| case IM_BANDFMT_INT: |
| RETURN_GENERATE( gint32 ); |
| #if 0 |
| case IM_BANDFMT_FLOAT: |
| RETURN_GENERATE( float ); |
| case IM_BANDFMT_DOUBLE: |
| RETURN_GENERATE( double ); |
| #endif |
| #undef RETURN_GENERATE |
| default: |
| g_assert( 0 ); |
| } |
| |
| /* Keep gcc happy. |
| */ |
| return 0; |
| #undef FUNCTION_NAME |
| } |
| |
| |
| /** |
| * im_grad_y: |
| * @in: input image |
| * @out: output image |
| * |
| * Find vertical differences between adjacent pixels. |
| * |
| * Generates an image where the value of each pixel is the difference between |
| * it and the pixel below it. The output has the same width as the input |
| * and one pixel less height. One-band integer formats only. The result is |
| * always %IM_BANDFMT_INT. |
| * |
| * This operation is much faster than (though equivalent to) im_conv() with the |
| * mask [[-1], [1]]. |
| * |
| * See also: im_grad_x(), im_conv(). |
| * |
| * Returns: 0 on success, -1 on error |
| */ |
| int im_grad_y( IMAGE *in, IMAGE *out ){ |
| #define FUNCTION_NAME "im_grad_y" |
| |
| if( im_piocheck( in, out ) ) |
| return -1; |
| |
| if( im_check_uncoded( "im_grad_y", in ) || |
| im_check_mono( "im_grad_y", in ) || |
| im_check_int( "im_grad_y", in ) ) |
| return( -1 ); |
| |
| if( im_cp_desc( out, in ) ) |
| return -1; |
| |
| -- out-> Ysize; |
| out-> BandFmt= IM_BANDFMT_INT; /* do not change without updating im_gradcor() */ |
| |
| if( im_demand_hint( out, IM_FATSTRIP, in, NULL ) ) |
| return -1; |
| |
| #define RETURN_GENERATE( TYPE ) return im_generate( out, im_start_one, ygrad_gen_ ## TYPE, im_stop_one, (void*) in, NULL ) |
| |
| switch( in-> BandFmt ){ |
| |
| case IM_BANDFMT_UCHAR: |
| RETURN_GENERATE( guint8 ); |
| |
| case IM_BANDFMT_CHAR: |
| RETURN_GENERATE( gint8 ); |
| |
| case IM_BANDFMT_USHORT: |
| RETURN_GENERATE( guint16 ); |
| |
| case IM_BANDFMT_SHORT: |
| RETURN_GENERATE( gint16 ); |
| |
| case IM_BANDFMT_UINT: |
| RETURN_GENERATE( guint32 ); |
| |
| case IM_BANDFMT_INT: |
| RETURN_GENERATE( gint32 ); |
| #if 0 |
| case IM_BANDFMT_FLOAT: |
| RETURN_GENERATE( float ); |
| case IM_BANDFMT_DOUBLE: |
| RETURN_GENERATE( double ); |
| #endif |
| #undef RETURN_GENERATE |
| default: |
| g_assert( 0 ); |
| } |
| |
| /* Keep gcc happy. |
| */ |
| return 0; |
| #undef FUNCTION_NAME |
| } |
| |
| |
| /** LOCAL FUNCTION DEFINITIONS **/ |
| |
| static void *gradcor_start( IMAGE *out, void *vptr_large, void *unrequired ){ |
| |
| gradcor_seq_t *seq= IM_NEW( NULL, gradcor_seq_t ); |
| if( ! seq ) |
| return NULL; |
| |
| seq-> region_xgrad= (int*) NULL; |
| seq-> region_ygrad= (int*) NULL; |
| seq-> region_xgrad_area= 0; |
| seq-> region_ygrad_area= 0; |
| |
| seq-> reg= im_region_create( (IMAGE*) vptr_large ); |
| if( ! seq-> reg ){ |
| im_free( (void*) seq ); |
| return NULL; |
| } |
| return (void*) seq; |
| } |
| |
| static int gradcor_stop( void *vptr_seq, void *unrequired, void *unreq2 ){ |
| |
| gradcor_seq_t *seq= (gradcor_seq_t*) vptr_seq; |
| if( seq ){ |
| im_free( (void*) seq-> region_xgrad ); |
| im_free( (void*) seq-> region_ygrad ); |
| im_region_free( seq-> reg ); |
| seq-> region_xgrad= (int*) NULL; |
| seq-> region_ygrad= (int*) NULL; |
| seq-> reg= (REGION*) NULL; |
| im_free( (void*) seq ); |
| } |
| return 0; |
| } |
| |
| static int gradcor_gen( REGION *to_make, void *vptr_seq, void *unrequired, void *vptr_grads ){ |
| |
| gradcor_seq_t *seq= (gradcor_seq_t*) vptr_seq; |
| REGION *make_from= seq-> reg; |
| |
| IMAGE **grads= (IMAGE**) vptr_grads; |
| IMAGE *small_xgrad= grads[0]; |
| IMAGE *small_ygrad= grads[1]; |
| |
| Rect require= { |
| to_make-> valid. left, |
| to_make-> valid. top, |
| to_make-> valid. width + small_xgrad-> Xsize, |
| to_make-> valid. height + small_ygrad-> Ysize |
| }; |
| size_t region_xgrad_width= require. width - 1; |
| size_t region_ygrad_height= require. height - 1; |
| |
| if( im_prepare( make_from, &require ) ) |
| return -1; |
| |
| #define FILL_BUFFERS( TYPE ) /* fill region_xgrad */ \ |
| { \ |
| TYPE *reading= (TYPE*) IM_REGION_ADDR( make_from, require. left, require. top ); \ |
| size_t read_skip= ( IM_REGION_LSKIP( make_from ) / sizeof(TYPE) ) - region_xgrad_width; \ |
| size_t area_need= region_xgrad_width * require. height; \ |
| \ |
| if( seq-> region_xgrad_area < area_need ){ \ |
| free( seq-> region_xgrad ); \ |
| seq-> region_xgrad= malloc( area_need * sizeof(int) ); \ |
| if( ! seq-> region_xgrad ) \ |
| return -1; \ |
| seq-> region_xgrad_area= area_need; \ |
| } \ |
| { \ |
| int *writing= seq-> region_xgrad; \ |
| int *write_end= writing + area_need; \ |
| int *write_stop; \ |
| for( ; writing < write_end; reading+= read_skip ) \ |
| for( write_stop= writing + region_xgrad_width; writing < write_stop; ++reading, ++writing ) \ |
| *writing= reading[1] - reading[0]; \ |
| } \ |
| } \ |
| { /* fill region_ygrad */ \ |
| TYPE *reading= (TYPE*) IM_REGION_ADDR( make_from, require. left, require. top ); \ |
| size_t read_line= IM_REGION_LSKIP( make_from ) / sizeof(TYPE); \ |
| size_t read_skip= read_line - require. width; \ |
| size_t area_need= require. width * region_ygrad_height; \ |
| \ |
| if( seq-> region_ygrad_area < area_need ){ \ |
| free( seq-> region_ygrad ); \ |
| seq-> region_ygrad= malloc( area_need * sizeof(int) ); \ |
| if( ! seq-> region_ygrad ) \ |
| return -1; \ |
| seq-> region_ygrad_area= area_need; \ |
| } \ |
| { \ |
| int *writing= seq-> region_ygrad; \ |
| int *write_end= writing + area_need; \ |
| int *write_stop; \ |
| for( ; writing < write_end; reading+= read_skip ) \ |
| for( write_stop= writing + require. width; writing < write_stop; ++reading, ++writing ) \ |
| *writing= reading[ read_line ] - reading[0]; \ |
| } \ |
| } |
| switch( make_from-> im-> BandFmt ){ |
| case IM_BANDFMT_UCHAR: |
| FILL_BUFFERS( unsigned char ) |
| break; |
| case IM_BANDFMT_CHAR: |
| FILL_BUFFERS( signed char ) |
| break; |
| case IM_BANDFMT_USHORT: |
| FILL_BUFFERS( unsigned short int ) |
| break; |
| case IM_BANDFMT_SHORT: |
| FILL_BUFFERS( signed short int ) |
| break; |
| case IM_BANDFMT_UINT: |
| FILL_BUFFERS( unsigned int ) |
| break; |
| case IM_BANDFMT_INT: |
| FILL_BUFFERS( signed int ) |
| break; |
| default: |
| g_assert( 0 ); |
| } |
| { /* write to output */ |
| size_t write_skip= IM_REGION_LSKIP( to_make ) / sizeof( float ); |
| float *writing= (float*) IM_REGION_ADDR_TOPLEFT( to_make ); |
| float *write_end= writing + write_skip * to_make-> valid. height; |
| float *write_stop; |
| size_t write_width= to_make-> valid. width; |
| |
| size_t small_xgrad_width= small_xgrad-> Xsize; |
| size_t small_ygrad_width= small_ygrad-> Xsize; |
| int *small_xgrad_end= (int*) small_xgrad-> data + small_xgrad_width * small_xgrad-> Ysize; |
| int *small_ygrad_end= (int*) small_ygrad-> data + small_ygrad_width * small_ygrad-> Ysize; |
| |
| int *region_xgrad_start= seq-> region_xgrad; |
| int *region_ygrad_start= seq-> region_ygrad; |
| size_t region_xgrad_start_skip= region_xgrad_width - write_width; |
| size_t region_ygrad_start_skip= require. width - write_width; |
| |
| size_t region_xgrad_read_skip= region_xgrad_width - small_xgrad_width; |
| size_t region_ygrad_read_skip= require. width - small_ygrad_width; |
| |
| write_skip-= write_width; |
| |
| for( ; writing < write_end; writing+= write_skip, region_xgrad_start+= region_xgrad_start_skip, region_ygrad_start+= region_ygrad_start_skip ) |
| for( write_stop= writing + write_width; writing < write_stop; ++writing, ++region_xgrad_start, ++region_ygrad_start ){ |
| gint64 sum= 0; |
| { |
| int *small_xgrad_read= (int*) small_xgrad-> data; |
| int *small_xgrad_stop; |
| int *region_xgrad_read= region_xgrad_start; |
| |
| for( ; small_xgrad_read < small_xgrad_end; region_xgrad_read+= region_xgrad_read_skip ) |
| for( small_xgrad_stop= small_xgrad_read + small_xgrad_width; small_xgrad_read < small_xgrad_stop; ++small_xgrad_read, ++region_xgrad_read ) |
| sum+= *small_xgrad_read * *region_xgrad_read; |
| } |
| { |
| int *small_ygrad_read= (int*) small_ygrad-> data; |
| int *small_ygrad_stop; |
| int *region_ygrad_read= region_ygrad_start; |
| |
| for( ; small_ygrad_read < small_ygrad_end; region_ygrad_read+= region_ygrad_read_skip ) |
| for( small_ygrad_stop= small_ygrad_read + small_ygrad_width; small_ygrad_read < small_ygrad_stop; ++small_ygrad_read, ++region_ygrad_read ) |
| sum+= *small_ygrad_read * *region_ygrad_read; |
| } |
| *writing= sum; |
| } |
| } |
| return 0; |
| } |
| |
| #define XGRAD_GEN_DEFINITION( TYPE ) \ |
| static int xgrad_gen_ ## TYPE( REGION *to_make, void *vptr_make_from, void *unrequired, void *unreq2 ){ \ |
| \ |
| REGION *make_from= (REGION*) vptr_make_from; \ |
| Rect require= { \ |
| to_make-> valid. left, \ |
| to_make-> valid. top, \ |
| to_make-> valid. width + 1, \ |
| to_make-> valid. height \ |
| }; \ |
| if( im_prepare( make_from, &require ) ) \ |
| return -1; \ |
| \ |
| { \ |
| int *writing= (int*) IM_REGION_ADDR_TOPLEFT( to_make ); \ |
| size_t write_skip= IM_REGION_LSKIP( to_make ) / sizeof(int); \ |
| int *write_end= writing + write_skip * to_make-> valid. height; \ |
| size_t write_width= to_make-> valid. width; \ |
| int *write_stop; \ |
| \ |
| TYPE *reading= (TYPE*) IM_REGION_ADDR( make_from, require. left, require. top ); \ |
| size_t read_skip= ( IM_REGION_LSKIP( make_from ) / sizeof(TYPE) ) - write_width; \ |
| \ |
| write_skip-= write_width; \ |
| \ |
| for( ; writing < write_end; writing+= write_skip, reading+= read_skip ) \ |
| for( write_stop= writing + write_width; writing < write_stop; ++writing, ++reading ) \ |
| *writing= (int)( reading[1] - reading[0] ); \ |
| } \ |
| return 0; \ |
| } |
| |
| #define YGRAD_GEN_DEFINITION( TYPE ) \ |
| static int ygrad_gen_ ## TYPE( REGION *to_make, void *vptr_make_from, void *unrequired, void *unreq2 ){ \ |
| \ |
| REGION *make_from= (REGION*) vptr_make_from; \ |
| Rect require= { \ |
| to_make-> valid. left, \ |
| to_make-> valid. top, \ |
| to_make-> valid. width, \ |
| to_make-> valid. height + 1 \ |
| }; \ |
| if( im_prepare( make_from, &require ) ) \ |
| return -1; \ |
| \ |
| { \ |
| int *writing= (int*) IM_REGION_ADDR_TOPLEFT( to_make ); \ |
| size_t write_skip= IM_REGION_LSKIP( to_make ) / sizeof(int); \ |
| int *write_end= writing + write_skip * to_make-> valid. height; \ |
| size_t write_width= to_make-> valid. width; \ |
| int *write_stop; \ |
| \ |
| TYPE *reading= (TYPE*) IM_REGION_ADDR( make_from, require. left, require. top ); \ |
| size_t read_line= IM_REGION_LSKIP( make_from ) / sizeof(TYPE); \ |
| size_t read_skip= read_line - write_width; \ |
| \ |
| write_skip-= write_width; \ |
| \ |
| for( ; writing < write_end; writing+= write_skip, reading+= read_skip ) \ |
| for( write_stop= writing + write_width; writing < write_stop; ++writing, ++reading ) \ |
| *writing= (int)( reading[ read_line ] - reading[0] ); \ |
| } \ |
| return 0; \ |
| } |
| |
| XGRAD_GEN_DEFINITION( guint8 ) |
| YGRAD_GEN_DEFINITION( guint8 ) |
| XGRAD_GEN_DEFINITION( gint8 ) |
| YGRAD_GEN_DEFINITION( gint8 ) |
| XGRAD_GEN_DEFINITION( guint16 ) |
| YGRAD_GEN_DEFINITION( guint16 ) |
| XGRAD_GEN_DEFINITION( gint16 ) |
| YGRAD_GEN_DEFINITION( gint16 ) |
| XGRAD_GEN_DEFINITION( guint32 ) |
| YGRAD_GEN_DEFINITION( guint32 ) |
| XGRAD_GEN_DEFINITION( gint32 ) |
| YGRAD_GEN_DEFINITION( gint32 ) |
| #if 0 |
| XGRAD_GEN_DEFINITION( float ) |
| YGRAD_GEN_DEFINITION( float ) |
| XGRAD_GEN_DEFINITION( double ) |
| YGRAD_GEN_DEFINITION( double ) |
| #endif |