1 /* find image overlaps
2  *
3  * Copyright: 1990, N. Dessipris.
4  *
5  * Author: Nicos Dessipris
6  * Written on: 02/05/1990
7  * Modified on : 18/04/1991
8  * 8/7/93 JC
9  *	- allows IM_CODING_LABQ coding
10  *	- now calls im_incheck()
11  * 13/7/95 JC
12  *	- rewritten
13  *	- now uses im_spcor()
14  * 13/8/96 JC
15  *	- order of args changed to help C++ API
16  * 24/1/11
17  * 	- gtk-doc
18  * 18/6/20 kleisauke
19  * 	- convert to vips8
20  */
21 
22 /*
23 
24     This file is part of VIPS.
25 
26     VIPS is free software; you can redistribute it and/or modify
27     it under the terms of the GNU Lesser General Public License as published by
28     the Free Software Foundation; either version 2 of the License, or
29     (at your option) any later version.
30 
31     This program is distributed in the hope that it will be useful,
32     but WITHOUT ANY WARRANTY; without even the implied warranty of
33     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
34     GNU Lesser General Public License for more details.
35 
36     You should have received a copy of the GNU Lesser General Public License
37     along with this program; if not, write to the Free Software
38     Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
39     02110-1301  USA
40 
41  */
42 
43 /*
44 
45     These files are distributed with VIPS - http://www.vips.ecs.soton.ac.uk
46 
47  */
48 
49 #ifdef HAVE_CONFIG_H
50 #include <config.h>
51 #endif /*HAVE_CONFIG_H*/
52 #include <vips/intl.h>
53 
54 #include <stdio.h>
55 #include <stdlib.h>
56 #include <math.h>
57 
58 #include <vips/vips.h>
59 
60 #include "pmosaicing.h"
61 
62 /* vips__correl:
63  * @ref: reference image
64  * @sec: secondary image
65  * @xref: position in reference image
66  * @yref: position in reference image
67  * @xsec: position in secondary image
68  * @ysec: position in secondary image
69  * @hwindowsize: half window size
70  * @hsearchsize: half search size
71  * @correlation: return detected correlation
72  * @x: return found position
73  * @y: return found position
74  *
75  * This operation finds the position of @sec within @ref.
76  *
77  * The area around
78  * (@xsec, @ysec) is searched for the best match to the area around (@xref,
79  * @yref). It  searches an area of size @hsearchsize for a
80  * match of size @hwindowsize.  The position of the best match is
81  * returned, together with the correlation at that point.
82  *
83  * Only  the  first  band  of each image is correlated. @ref and @sec may be
84  * very large --- the function  extracts  and  generates  just  the
85  * parts needed.  Correlation is done with vips_spcor(); the position of
86  * the maximum is found with vips_max().
87  *
88  * See also: vips_match(), vips__lrmosaic().
89  *
90  * Returns: 0 on success, -1 on error
91  */
92 int
vips__correl(VipsImage * ref,VipsImage * sec,int xref,int yref,int xsec,int ysec,int hwindowsize,int hsearchsize,double * correlation,int * x,int * y)93 vips__correl( VipsImage *ref, VipsImage *sec,
94 	int xref, int yref, int xsec, int ysec,
95 	int hwindowsize, int hsearchsize,
96 	double *correlation, int *x, int *y )
97 {
98 	VipsImage *surface = vips_image_new();
99 	VipsImage **t = (VipsImage **)
100 		vips_object_local_array( VIPS_OBJECT( surface ), 5 );
101 
102 	VipsRect refr, secr;
103 	VipsRect winr, srhr;
104 	VipsRect wincr, srhcr;
105 
106 	/* Find position of window and search area, and clip against image
107 	 * size.
108 	 */
109 	refr.left = 0;
110 	refr.top = 0;
111 	refr.width = ref->Xsize;
112 	refr.height = ref->Ysize;
113 	winr.left = xref - hwindowsize;
114 	winr.top = yref - hwindowsize;
115 	winr.width = hwindowsize * 2 + 1;
116 	winr.height = hwindowsize * 2 + 1;
117 	vips_rect_intersectrect( &refr, &winr, &wincr );
118 
119 	secr.left = 0;
120 	secr.top = 0;
121 	secr.width = sec->Xsize;
122 	secr.height = sec->Ysize;
123 	srhr.left = xsec - hsearchsize;
124 	srhr.top = ysec - hsearchsize;
125 	srhr.width = hsearchsize * 2 + 1;
126 	srhr.height = hsearchsize * 2 + 1;
127 	vips_rect_intersectrect( &secr, &srhr, &srhcr );
128 
129 	/* Extract window and search area.
130 	 */
131 	if( vips_extract_area( ref, &t[0],
132 			wincr.left, wincr.top, wincr.width, wincr.height,
133 			NULL ) ||
134 		vips_extract_area( sec, &t[1],
135 			srhcr.left, srhcr.top, srhcr.width, srhcr.height,
136 			NULL ) ) {
137 		g_object_unref( surface );
138 		return( -1 );
139 	}
140 	ref = t[0];
141 	sec = t[1];
142 
143 	/* Make sure we have just one band. From vips_*mosaic() we will, but
144 	 * from vips_match() etc. we may not.
145 	 */
146 	if( ref->Bands != 1 ) {
147 		if( vips_extract_band( ref, &t[2], 0, NULL ) ) {
148 			g_object_unref( surface );
149 			return( -1 );
150 		}
151 		ref = t[2];
152 	}
153 	if( sec->Bands != 1 ) {
154 		if( vips_extract_band( sec, &t[3], 0, NULL ) ) {
155 			g_object_unref( surface );
156 			return( -1 );
157 		}
158 		sec = t[3];
159 	}
160 
161 	/* Search!
162 	 */
163 	if( vips_spcor( sec, ref, &t[4], NULL ) ) {
164 		g_object_unref( surface );
165 		return( -1 );
166 	}
167 
168 	/* Find maximum of correlation surface.
169 	 */
170 	if( vips_max( t[4], correlation, "x", x, "y", y, NULL ) ) {
171 		g_object_unref( surface );
172 		return( -1 );
173 	}
174 	g_object_unref( surface );
175 
176 	/* Translate back to position within sec.
177 	 */
178 	*x += srhcr.left;
179 	*y += srhcr.top;
180 
181 	return( 0 );
182 }
183 
184 int
vips__chkpair(VipsImage * ref,VipsImage * sec,TiePoints * points)185 vips__chkpair( VipsImage *ref, VipsImage *sec, TiePoints *points )
186 {
187 	int i;
188 	int x, y;
189 	double correlation;
190 
191 	const int hcor = points->halfcorsize;
192 	const int harea = points->halfareasize;
193 
194 	/* Check images.
195 	 */
196 	if( vips_image_wio_input( ref ) || vips_image_wio_input( sec ) )
197 		return( -1 );
198 	if( ref->Bands != sec->Bands || ref->BandFmt != sec->BandFmt ||
199 		ref->Coding != sec->Coding ) {
200 		vips_error( "vips_chkpair", "%s", _( "inputs incompatible" ) );
201 		return( -1 );
202 	}
203 	if( ref->Bands != 1 || ref->BandFmt != VIPS_FORMAT_UCHAR ) {
204 		vips_error( "vips_chkpair", "%s", _( "help!" ) );
205 		return( -1 );
206 	}
207 
208 	for( i = 0; i < points->nopoints; i++ ) {
209 		/* Find correlation point.
210 		 */
211 		if( vips__correl( ref, sec,
212 			points->x_reference[i], points->y_reference[i],
213 			points->x_reference[i], points->y_reference[i],
214 			hcor, harea,
215 			&correlation, &x, &y ) )
216 			return( -1 );
217 
218 		/* And note in x_secondary.
219 		 */
220 		points->x_secondary[i] = x;
221 		points->y_secondary[i] = y;
222 		points->correlation[i] = correlation;
223 
224 		/* Note each dx, dy too.
225 		 */
226 		points->dx[i] =
227 			points->x_secondary[i] - points->x_reference[i];
228 		points->dy[i] =
229 			points->y_secondary[i] - points->y_reference[i];
230 	}
231 
232 	return( 0 );
233 }
234