From f618466c25d43f3bae9e40920273bf77de1e1149 Mon Sep 17 00:00:00 2001 From: leochanj105 Date: Mon, 19 Oct 2020 23:09:30 -0400 Subject: initial sd-vbs initial sd-vbs add sd-vbs sd-vbs --- SD-VBS/benchmarks/sift/src/matlab/siftrefinemx.c | 342 +++++++++++++++++++++++ 1 file changed, 342 insertions(+) create mode 100644 SD-VBS/benchmarks/sift/src/matlab/siftrefinemx.c (limited to 'SD-VBS/benchmarks/sift/src/matlab/siftrefinemx.c') diff --git a/SD-VBS/benchmarks/sift/src/matlab/siftrefinemx.c b/SD-VBS/benchmarks/sift/src/matlab/siftrefinemx.c new file mode 100644 index 0000000..a1ef8c7 --- /dev/null +++ b/SD-VBS/benchmarks/sift/src/matlab/siftrefinemx.c @@ -0,0 +1,342 @@ +/* file: siftrefinemx.c + +Benchmark - sift +Data set - sqcif + + < M A T L A B > + Copyright 1984-2006 The MathWorks, Inc. + Version 7.3.0.298 (R2006b) + August 03, 2006 + + + To get started, type one of these: helpwin, helpdesk, or demo. + For product information, visit www.mathworks.com. + +Warning: Function /h/g2/kvs/checkParallel/sdvbs-svn/common/matlab/randn.m has the same name as a MATLAB builtin. We suggest you rename the function to avoid a potential name conflict. +> In path at 113 + In script_run_profile at 3 +Warning: You are using gcc version "4.1.1". The earliest gcc version supported +with mex is "3.4.0". The latest version tested for use with mex is "3.4.5". +To download a different version of gcc, visit http://gcc.gnu.org +Warning: You are using gcc version "4.1.1". The earliest gcc version supported +with mex is "3.4.0". The latest version tested for use with mex is "3.4.5". +To download a different version of gcc, visit http://gcc.gnu.org +Warning: You are using gcc version "4.1.1". The earliest gcc version supported +with mex is "3.4.0". The latest version tested for use with mex is "3.4.5". +To download a different version of gcc, visit http://gcc.gnu.org +Warning: You are using gcc version "4.1.1". The earliest gcc version supported +with mex is "3.4.0". The latest version tested for use with mex is "3.4.5". +To download a different version of gcc, visit http://gcc.gnu.org +Warning: You are using gcc version "4.1.1". The earliest gcc version supported +with mex is "3.4.0". The latest version tested for use with mex is "3.4.5". +To download a different version of gcc, visit http://gcc.gnu.org +Warning: You are using gcc version "4.1.1". The earliest gcc version supported +with mex is "3.4.0". The latest version tested for use with mex is "3.4.5". +To download a different version of gcc, visit http://gcc.gnu.org +Input size - (96x96) +** author: Andrea Vedaldi +** description: Subpixel localization, thresholding and off-edge test. +**/ + +/* AUTORIGHTS +Copyright (c) 2006 The Regents of the University of California. +All Rights Reserved. + +Created by Andrea Vedaldi +UCLA Vision Lab - Department of Computer Science + +Permission to use, copy, modify, and distribute this software and its +documentation for educational, research and non-profit purposes, +without fee, and without a written agreement is hereby granted, +provided that the above copyright notice, this paragraph and the +following three paragraphs appear in all copies. + +This software program and documentation are copyrighted by The Regents +of the University of California. The software program and +documentation are supplied "as is", without any accompanying services +from The Regents. The Regents does not warrant that the operation of +the program will be uninterrupted or error-free. The end-user +understands that the program was developed for research purposes and +is advised not to rely exclusively on the program for any reason. + +This software embodies a method for which the following patent has +been issued: "Method and apparatus for identifying scale invariant +features in an image and use of same for locating an object in an +image," David G. Lowe, US Patent 6,711,293 (March 23, +2004). Provisional application filed March 8, 1999. Asignee: The +University of British Columbia. + +IN NO EVENT SHALL THE UNIVERSITY OF CALIFORNIA BE LIABLE TO ANY PARTY +FOR DIRECT, INDIRECT, SPECIAL, INCIDENTAL, OR CONSEQUENTIAL DAMAGES, +INCLUDING LOST PROFITS, ARISING OUT OF THE USE OF THIS SOFTWARE AND +ITS DOCUMENTATION, EVEN IF THE UNIVERSITY OF CALIFORNIA HAS BEEN +ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. THE UNIVERSITY OF +CALIFORNIA SPECIFICALLY DISCLAIMS ANY WARRANTIES, INCLUDING, BUT NOT +LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR +A PARTICULAR PURPOSE. THE SOFTWARE PROVIDED HEREUNDER IS ON AN "AS IS" +BASIS, AND THE UNIVERSITY OF CALIFORNIA HAS NO OBLIGATIONS TO PROVIDE +MAINTENANCE, SUPPORT, UPDATES, ENHANCEMENTS, OR MODIFICATIONS. +*/ + +#include"mex.h" + +#include + +#include +#include + +/* Prototype of DGESV LAPACK function for the solution of a linear system. */ +#ifdef WINDOWS +#define DGESV dgesv +#undef min +#undef max +#else +#define DGESV dgesv_ +#endif + +#ifdef WINDOWS +#ifdef __cplusplus__ +extern "C" { + extern int DGESV(int *n, int *nrhs, double *a, int *lda, + int *ipiv, double *b, int *ldb, int *info) ; +} +#else + extern int DGESV(int *n, int *nrhs, double *a, int *lda, + int *ipiv, double *b, int *ldb, int *info) ; +#define sqrtf(x) ((float)sqrt((double)x) +#define powf(x) ((float)pow((double)x) +#define fabsf(x) ((float)fabs((double)x) +#endif +#else +extern int DGESV(int *n, int *nrhs, double *a, int *lda, + int *ipiv, double *b, int *ldb, int *info) ; +#endif + +#define greater(a,b) ((a) > (b)) +#define min(a,b) (((a)<(b))?(a):(b)) +#define max(a,b) (((a)>(b))?(a):(b)) + +const int max_iter = 5 ; + +void +mexFunction(int nout, mxArray *out[], + int nin, const mxArray *in[]) +{ + int M,N,S,smin,K ; + const int* dimensions ; + const double* P_pt ; + const double* D_pt ; + double threshold = 0.01 ; /*0.02 ;*/ + double r = 10.0 ; + double* result ; + enum {IN_P=0,IN_D,IN_SMIN,IN_THRESHOLD,IN_R} ; + enum {OUT_Q=0} ; + + /* ----------------------------------------------------------------- + ** Check the arguments + ** -------------------------------------------------------------- */ + if (nin < 3) { + mexErrMsgTxt("At least three input arguments required."); + } else if (nout > 1) { + mexErrMsgTxt("Too many output arguments."); + } + + if( !uIsRealMatrix(in[IN_P],3,-1) ) { + mexErrMsgTxt("P must be a 3xK real matrix") ; + } + + if( !mxIsDouble(in[IN_D]) || mxGetNumberOfDimensions(in[IN_D]) != 3) { + mexErrMsgTxt("G must be a three dimensional real array.") ; + } + + if( !uIsRealScalar(in[IN_SMIN]) ) { + mexErrMsgTxt("SMIN must be a real scalar.") ; + } + + if(nin >= 4) { + if(!uIsRealScalar(in[IN_THRESHOLD])) { + mexErrMsgTxt("THRESHOLD must be a real scalar.") ; + } + threshold = *mxGetPr(in[IN_THRESHOLD]) ; + } + + if(nin >= 5) { + if(!uIsRealScalar(in[IN_R])) { + mexErrMsgTxt("R must be a real scalar.") ; + } + r = *mxGetPr(in[IN_R]) ; + } + + dimensions = mxGetDimensions(in[IN_D]) ; + M = dimensions[0] ; + N = dimensions[1] ; + S = dimensions[2] ; + smin = (int)(*mxGetPr(in[IN_SMIN])) ; + + if(S < 3 || M < 3 || N < 3) { + mexErrMsgTxt("All dimensions of DOG must be not less than 3.") ; + } + + K = mxGetN(in[IN_P]) ; + P_pt = mxGetPr(in[IN_P]) ; + D_pt = mxGetPr(in[IN_D]) ; + + /* If the input array is empty, then output an empty array as well. */ + if( K == 0) { + out[OUT_Q] = mxDuplicateArray(in[IN_P]) ; + return ; + } + + /* ----------------------------------------------------------------- + * Do the job + * -------------------------------------------------------------- */ + { + double* buffer = (double*) mxMalloc(K*3*sizeof(double)) ; + double* buffer_iterator = buffer ; + int p ; + const int yo = 1 ; + const int xo = M ; + const int so = M*N ; + +/* printf("Actual values = %d\n\n", K); +*/ + for(p = 0 ; p < K ; ++p) { + int x = ((int)*P_pt++) ; + int y = ((int)*P_pt++) ; +/* printf("%d\t%d\n", ((int)*P_pt), smin); +*/ + int s = ((int)*P_pt++) - smin ; + int iter ; + double b[3] ; + + /* Local maxima extracted from the DOG + * have coorrinates 1<=x<=N-2, 1<=y<=M-2 + * and 1<=s-mins<=S-2. This is also the range of the points + * that we can refine. + */ +/* printf("%d\t%d\t%d\t%d\t%d\t%d\n", x, N-2,y,M-2, s, S-2); +*/ + if(x < 1 || x > N-2 || + y < 1 || y > M-2 || + s < 1 || s > S-2) { + continue ; + } + +#define at(dx,dy,ds) (*(pt + (dx)*xo + (dy)*yo + (ds)*so)) + + { + const double* pt = D_pt + y*yo + x*xo + s*so ; + + double Dx=0,Dy=0,Ds=0,Dxx=0,Dyy=0,Dss=0,Dxy=0,Dxs=0,Dys=0 ; + int dx = 0 ; + int dy = 0 ; +/* printf("%d\t%d\t%d\t%d\t%d\t%d\t%d\t%f\t%f\n",S, y, yo, x, xo, s, so, *D_pt, *pt); +*/ + for(iter = 0 ; iter < max_iter ; ++iter) { + + double A[3*3] ; + int ipiv[3] ; + int n = 3 ; + int one = 1 ; + int info = 0 ; + +#define Aat(i,j) (A[(i)+(j)*3]) + + x += dx ; + y += dy ; + pt = D_pt + y*yo + x*xo + s*so ; + + /* Compute the gradient. */ + Dx = 0.5 * (at(+1,0,0) - at(-1,0,0)) ; + Dy = 0.5 * (at(0,+1,0) - at(0,-1,0)); + Ds = 0.5 * (at(0,0,+1) - at(0,0,-1)) ; + + /* Compute the Hessian. */ + Dxx = (at(+1,0,0) + at(-1,0,0) - 2.0 * at(0,0,0)) ; + Dyy = (at(0,+1,0) + at(0,-1,0) - 2.0 * at(0,0,0)) ; + Dss = (at(0,0,+1) + at(0,0,-1) - 2.0 * at(0,0,0)) ; + + Dxy = 0.25 * ( at(+1,+1,0) + at(-1,-1,0) - at(-1,+1,0) - at(+1,-1,0) ) ; + Dxs = 0.25 * ( at(+1,0,+1) + at(-1,0,-1) - at(-1,0,+1) - at(+1,0,-1) ) ; + Dys = 0.25 * ( at(0,+1,+1) + at(0,-1,-1) - at(0,-1,+1) - at(0,+1,-1) ) ; + + /* Solve linear system. */ + Aat(0,0) = Dxx ; + Aat(1,1) = Dyy ; + Aat(2,2) = Dss ; + Aat(0,1) = Aat(1,0) = Dxy ; + Aat(0,2) = Aat(2,0) = Dxs ; + Aat(1,2) = Aat(2,1) = Dys ; + + b[0] = - Dx ; + b[1] = - Dy ; + b[2] = - Ds ; + + /* DGESV (&n, &one, A, &n, ipiv, b, &n, &info) ; + */ + /* If the translation of the keypoint is big, move the keypoint + * and re-iterate the computation. Otherwise we are all set. + */ + dx= ((b[0] > 0.6 && x < N-2) ? 1 : 0 ) + + ((b[0] < -0.6 && x > 1 ) ? -1 : 0 ) ; + + dy= ((b[1] > 0.6 && y < M-2) ? 1 : 0 ) + + ((b[1] < -0.6 && y > 1 ) ? -1 : 0 ) ; + + if( dx == 0 && dy == 0 ) break ; + + } + + { + double val = at(0,0,0) + 0.5 * (Dx * b[0] + Dy * b[1] + Ds * b[2]) ; + double score = (Dxx+Dyy)*(Dxx+Dyy) / (Dxx*Dyy - Dxy*Dxy) ; + double xn = x + b[0] ; + double yn = y + b[1] ; + double sn = s + b[2] ; + +/* printf("%f\t%f\t%f\t%f\t%f\t%f\t%f\t%f\t%f\t%f\t%f\n", fabs(val),threshold,score,(r+1)*(r+1)/r,fabs(b[0]), fabs(b[1]), fabs(b[2]),xn,yn,sn,r); +*/ + if(fabs(val) > threshold && + score < (r+1)*(r+1)/r && + score >= 0 && + fabs(b[0]) < 1.5 && + fabs(b[1]) < 1.5 && + fabs(b[2]) < 1.5 && + xn >= 0 && + xn <= N-1 && + yn >= 0 && + yn <= M-1 && + sn >= 0 && + sn <= S-1) { + *buffer_iterator++ = xn ; + *buffer_iterator++ = yn ; + *buffer_iterator++ = sn+smin ; + } + } + } + } + + /* Copy the result into an array. */ + { + int i; + int NL = (buffer_iterator - buffer)/3 ; +/* printf("%NL VALUE = %d\t%d\t%d\n", NL, buffer_iterator, buffer); +*/ + out[OUT_Q] = mxCreateDoubleMatrix(3, NL, mxREAL) ; + result = mxGetPr(out[OUT_Q]); + for(i=0; i<(3*NL); i++) + { + result[i] = buffer[i]; +/* printf("%f\t", buffer[i]); +*/ + } +/* printf("\n"); + memcpy(result, buffer, sizeof(double) * 3 * NL) ; +*/ + } + mxFree(buffer) ; + } + +} + + -- cgit v1.2.2