summaryrefslogtreecommitdiffstats
path: root/SD-VBS/benchmarks/sift/src/matlab/siftlocalmax.c
diff options
context:
space:
mode:
Diffstat (limited to 'SD-VBS/benchmarks/sift/src/matlab/siftlocalmax.c')
-rw-r--r--SD-VBS/benchmarks/sift/src/matlab/siftlocalmax.c291
1 files changed, 0 insertions, 291 deletions
diff --git a/SD-VBS/benchmarks/sift/src/matlab/siftlocalmax.c b/SD-VBS/benchmarks/sift/src/matlab/siftlocalmax.c
deleted file mode 100644
index 3646692..0000000
--- a/SD-VBS/benchmarks/sift/src/matlab/siftlocalmax.c
+++ /dev/null
@@ -1,291 +0,0 @@
1/* file: siftlocalmax.c
2** author: Andrea Vedaldi
3** description: Find local maximizer of multi-dimensional array.
4**/
5
6/* AUTORIGHTS
7Copyright (c) 2006 The Regents of the University of California.
8All Rights Reserved.
9
10Created by Andrea Vedaldi
11UCLA Vision Lab - Department of Computer Science
12
13Permission to use, copy, modify, and distribute this software and its
14documentation for educational, research and non-profit purposes,
15without fee, and without a written agreement is hereby granted,
16provided that the above copyright notice, this paragraph and the
17following three paragraphs appear in all copies.
18
19This software program and documentation are copyrighted by The Regents
20of the University of California. The software program and
21documentation are supplied "as is", without any accompanying services
22from The Regents. The Regents does not warrant that the operation of
23the program will be uninterrupted or error-free. The end-user
24understands that the program was developed for research purposes and
25is advised not to rely exclusively on the program for any reason.
26
27This software embodies a method for which the following patent has
28been issued: "Method and apparatus for identifying scale invariant
29features in an image and use of same for locating an object in an
30image," David G. Lowe, US Patent 6,711,293 (March 23,
312004). Provisional application filed March 8, 1999. Asignee: The
32University of British Columbia.
33
34IN NO EVENT SHALL THE UNIVERSITY OF CALIFORNIA BE LIABLE TO ANY PARTY
35FOR DIRECT, INDIRECT, SPECIAL, INCIDENTAL, OR CONSEQUENTIAL DAMAGES,
36INCLUDING LOST PROFITS, ARISING OUT OF THE USE OF THIS SOFTWARE AND
37ITS DOCUMENTATION, EVEN IF THE UNIVERSITY OF CALIFORNIA HAS BEEN
38ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. THE UNIVERSITY OF
39CALIFORNIA SPECIFICALLY DISCLAIMS ANY WARRANTIES, INCLUDING, BUT NOT
40LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
41A PARTICULAR PURPOSE. THE SOFTWARE PROVIDED HEREUNDER IS ON AN "AS IS"
42BASIS, AND THE UNIVERSITY OF CALIFORNIA HAS NO OBLIGATIONS TO PROVIDE
43MAINTENANCE, SUPPORT, UPDATES, ENHANCEMENTS, OR MODIFICATIONS.
44*/
45
46#include"mex.h"
47#include<mexutils.c>
48#include<stdlib.h>
49
50/** Matlab driver.
51 **/
52#define greater(a,b) ((a) > (b)+threshold)
53
54void
55mexFunction(int nout, mxArray *out[],
56 int nin, const mxArray *in[])
57{
58 int M, N ;
59 const double* F_pt ;
60 int ndims ;
61 int pdims = -1 ;
62 int* offsets ;
63 int* midx ;
64 int* neighbors ;
65 int nneighbors ;
66 int* dims ;
67 enum {F=0,THRESHOLD,P} ;
68 enum {MAXIMA=0} ;
69 double threshold = - mxGetInf() ;
70
71
72 /* ------------------------------------------------------------------
73 * Check the arguments
74 * --------------------------------------------------------------- */
75 if(nin > 1) {
76 if(!uIsRealScalar(in[THRESHOLD])) {
77 mexErrMsgTxt("THRESHOLD must be a real scalar.") ;
78 }
79 threshold = *mxGetPr(in[THRESHOLD]) ;
80 }
81
82 ndims = mxGetNumberOfDimensions(in[F]) ;
83
84if(ndims !=3)
85 printf("Error Error: NDIMS not equal to 3. Not handled by C version\n");
86
87 {
88 /* We need to make a copy because in one special case (see below)
89 we need to adjust dims[].
90 */
91 int d ;
92 const int* const_dims = (int*) mxGetDimensions(in[F]) ;
93 dims = mxMalloc(sizeof(int)*ndims) ;
94 for(d=0 ; d < ndims ; ++d)
95 {
96/* printf("Const Dimensions = %d\n", const_dims[d]);
97*/
98 dims[d] = const_dims[d] ;
99 }
100 }
101 M = dims[0] ;
102 N = dims[1] ;
103 F_pt = mxGetPr(in[F]) ;
104
105 /*
106 If there are only two dimensions and if one is singleton, then
107 assume that a vector has been provided as input (and treat this
108 as a COLUMN matrix with p=1). We do this because Matlab does not
109 distinguish between vectors and 1xN or Mx1 matrices and because
110 the cases 1xN and Mx1 are trivial (the result is alway empty).
111 */
112 if((ndims == 2) && (pdims < 0) && (M == 1 || N == 1)) {
113
114 printf("ERROR ERROR: Entered a different loop here. Not handled by C version");
115
116 pdims = 1 ;
117 M = (M>N)?M:N ;
118 N = 1 ;
119 dims[0]=M ;
120 dims[1]=N ;
121 }
122
123 /* search the local maxima along the first p dimensions only */
124 if(pdims < 0)
125 {
126 pdims = ndims ;
127 }
128
129 if(pdims > ndims) {
130 mxFree(dims) ;
131 mexErrMsgTxt("P must not be greater than the number of dimensions") ;
132 }
133
134 /* ------------------------------------------------------------------
135 * Do the job
136 * --------------------------------------------------------------- */
137 {
138 int maxima_size = M*N ;
139 int* maxima_start = (int*) mxMalloc(sizeof(int) * maxima_size) ;
140 int* maxima_iterator = maxima_start ;
141 int* maxima_end = maxima_start + maxima_size ;
142 int i,j,h,o ;
143 const double* pt = F_pt ;
144
145 /* Compute the offsets between dimensions. */
146 offsets = (int*) mxMalloc(sizeof(int) * ndims) ;
147 offsets[0] = 1 ;
148
149 for(h = 1 ; h < ndims ; ++h)
150 {
151 offsets[h] = offsets[h-1]*dims[h-1] ;
152/* printf("%d:%d\t%d\n", h, offsets[h], dims[h-1]);
153*/
154 }
155
156 /* Multi-index. */
157 midx = (int*) mxMalloc(sizeof(int) * ndims) ;
158 for(h = 0 ; h < ndims ; ++h)
159 midx[h] = 1 ;
160
161 /* Neighbors. */
162 nneighbors = 1 ;
163 o=0 ;
164 for(h = 0 ; h < pdims ; ++h) {
165 nneighbors *= 3 ;
166 midx[h] = -1 ;
167 o -= offsets[h] ;
168 }
169 nneighbors -= 1 ;
170 neighbors = (int*) mxMalloc(sizeof(int) * nneighbors) ;
171
172 /* Precompute offsets from offset(-1,...,-1,0,...0) to
173 * offset(+1,...,+1,0,...,0). */
174 i = 0 ;
175
176 while(true) {
177 if(o != 0)
178 neighbors[i++] = o ;
179 h = 0 ;
180 while( o += offsets[h], (++midx[h]) > 1 ) {
181 o -= 3*offsets[h] ;
182 midx[h] = -1 ;
183 if(++h >= pdims)
184 goto stop ;
185 }
186 }
187 stop: ;
188
189 /* Starts at the corner (1,1,...,1,0,0,...0) */
190 for(h = 0 ; h < pdims ; ++h) {
191 midx[h] = 1 ;
192 pt += offsets[h] ;
193/* printf("%d:%x\t%d\n", h, pt, offsets[h]);
194*/
195 }
196
197 for(h = pdims ; h < ndims ; ++h) {
198 midx[h] = 0 ;
199 }
200
201 /* ---------------------------------------------------------------
202 * Loop
203 * ------------------------------------------------------------ */
204
205 /*
206 If any dimension in the first P is less than 3 elements wide
207 then just return the empty matrix (if we proceed without doing
208 anything we break the carry reporting algorithm below).
209 */
210 for(h=0 ; h < pdims ; ++h)
211 if(dims[h] < 3) goto end ;
212
213 while(true) {
214
215 /* Propagate carry along multi index midx */
216 h = 0 ;
217 while((midx[h]) >= dims[h] - 1) {
218/* pt += 2*offsets[h] ; skip first and last el. */
219 midx[h] = 1 ;
220 if(++h >= pdims)
221 goto next_layer ;
222 ++midx[h] ;
223 }
224
225 /* Scan neighbors */
226 {
227 double v = *pt ;
228 bool is_greater = (v >= threshold) ;
229 /* printf("%f\n", v);
230 */
231 i = 0 ;
232 while(is_greater && i < nneighbors)
233 {
234 is_greater &= v > *(pt + neighbors[i++]) ;
235 }
236
237 /* Add the local maximum */
238 if(is_greater) {
239 /* Need more space? */
240 if(maxima_iterator == maxima_end) {
241 maxima_size += M*N ;
242 maxima_start = (int*) mxRealloc(maxima_start,
243 maxima_size*sizeof(int)) ;
244 maxima_end = maxima_start + maxima_size ;
245 maxima_iterator = maxima_end - M*N ;
246 }
247
248 *maxima_iterator++ = pt - F_pt + 1 ;
249 }
250
251 /* Go to next element */
252 pt += 1 ;
253 ++midx[0] ;
254 continue ;
255
256 next_layer: ;
257 if( h >= ndims )
258 goto end ;
259
260 while((++midx[h]) >= dims[h]) {
261 midx[h] = 0 ;
262 if(++h >= ndims)
263 goto end ;
264 }
265 }
266 }
267 end:;
268 /* Return. */
269 {
270 double* M_pt ;
271/* printf("%x\t%x\t%d\n", maxima_iterator, maxima_start, maxima_iterator-maxima_start);
272*/
273 out[MAXIMA] = mxCreateDoubleMatrix
274 (1, maxima_iterator-maxima_start, mxREAL) ;
275 maxima_end = maxima_iterator ;
276 maxima_iterator = maxima_start ;
277 M_pt = mxGetPr(out[MAXIMA]) ;
278 while(maxima_iterator != maxima_end)
279 {
280 *M_pt++ = *maxima_iterator++ ;
281 }
282 }
283
284 /* Release space. */
285 mxFree(offsets) ;
286 mxFree(neighbors) ;
287 mxFree(midx) ;
288 mxFree(maxima_start) ;
289 }
290 mxFree(dims) ;
291}
generated by cgit v1.2.2 (git 2.25.0) at 2025-04-30 17:25:49 -0400