summaryrefslogtreecommitdiffstats
path: root/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX
diff options
context:
space:
mode:
Diffstat (limited to 'SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX')
-rwxr-xr-xSD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/-MacReadMe1
-rwxr-xr-xSD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/.AppleDouble/.Parentbin589 -> 0 bytes
-rwxr-xr-xSD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/.FBCIndexbin258048 -> 0 bytes
-rwxr-xr-xSD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/.FBCLockFolder/.FBCSemaphoreFilebin6 -> 0 bytes
-rwxr-xr-xSD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/Makefile-linux39
-rwxr-xr-xSD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/Makefile-ml6-linux39
-rwxr-xr-xSD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/Makefile-osx39
-rwxr-xr-xSD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/Makefile-solaris38
-rwxr-xr-xSD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/Makefile-sun439
-rwxr-xr-xSD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/convolve.c325
-rwxr-xr-xSD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/convolve.h55
-rwxr-xr-xSD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/corrDn.c145
-rwxr-xr-xSD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/edges-orig.c494
-rwxr-xr-xSD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/edges.c647
-rwxr-xr-xSD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/histo.c140
-rwxr-xr-xSD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/innerProd.c52
-rwxr-xr-xSD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/innerProd.dllbin41984 -> 0 bytes
-rwxr-xr-xSD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/innerProd.mexglxbin17704 -> 0 bytes
-rwxr-xr-xSD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/innerProd.mexlxbin4479 -> 0 bytes
-rwxr-xr-xSD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/innerProd.mexmacbin21600 -> 0 bytes
-rwxr-xr-xSD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/innerProd.mexsolbin3948 -> 0 bytes
-rwxr-xr-xSD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/pointOp.c126
-rwxr-xr-xSD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/range2.c56
-rwxr-xr-xSD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/upConv.c195
-rwxr-xr-xSD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/wrap.c281
25 files changed, 0 insertions, 2711 deletions
diff --git a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/-MacReadMe b/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/-MacReadMe
deleted file mode 100755
index 898dc0c..0000000
--- a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/-MacReadMe
+++ /dev/null
@@ -1 +0,0 @@
1MacReadMe How to compile a MEX file for Macintosh (Based on a similar note written by David Brainard and Denis Pelli to accompany the Psychophysics Toolbox.) A MEX file (short for "MATLAB extension") contains code that implements a new MATLAB function, allowing users to use the full power of the C language and the Macintosh hardware and operating system in implementing a new MATLAB function. This document explains how to produce a MEX file that is "fat" (i.e. will run as native code on both 68K and PowerPC Macs) and is compatible with both MATLAB 4 and 5. 1) To produce a MATLAB MEX file with PowerPC code you must have the Metrowerks CodeWarrior C compiler (version 10 or better, abbreviated as CW below). To produce 68K code we still use the Symantec THINK C compiler (version from Symantec C++ 8 CD-ROM release 5), but we will soon be switching to Metrowerks CodeWarrior. (See note A below.) 2) Place a copy of the MATLAB 4:Extern folder, supplied by Mathworks, on your compiler's search path. We suggest that you name the copy "MEX V4". (See notes B and C, below.) 3) Build any of the MEX files simply by opening its project file and asking your compiler to "Build Code Resource" (in THINK C) or to "Make" it (in CW). For each MEX file, e.g. histo.mex, there are two associated projects, e.g. histo.µ for CW, and histo.¹.4 for THINK C. To build a "fat" MEX, that runs native on both 68K and PowerPC, you should first compile in THINK C, and then in CW. (See note A, below.) Denis Pelli April 2, 1997 Notes A) The Mathworks support only the THINK C compiler to make 68K MEX code for MATLAB version 4 and only the CW compiler to make PPC MEX files for MATLAB 4 and both 68K and PPC for MATLAB 5. This archive includes THINK and CW projects. To build a fat MEX file for MATLAB 4, first "make" the THINK C version (e.g. histo.¹.4), producing a file with a .rsrc extension (e.g. histo.µ.rsrc). This is the 68K MEX file. When you then "make" histo.µ, the CW project incorporates the .rsrc file and generates a "fat" MEX file that will run native (i.e. fast) on both 68K and PowerPC. To make a 68K-only MEX file, simply rename, e.g., histo.µ.rsrc to histo.mex after you make the THINK project, and set the file type and creator to match the other MEX files. THINK C is slow and hard to work with. Symantec hasn't significantly upgraded in it many years. There is an error in the math.h header (version from Symantec C++ 8 CD-ROM release 5). We fix that error by some tricky preprocessor defines and undefines in the THINK C Prefix in each of the THINK projects. B) The easiest way to put a folder on your compilerÕs search path is simply to put the folder in the same folder as the compiler itself. If you want to use both CW and THINK C, then put the folder under CW, make an alias of it, and put the alias in THINK C's "Aliases" folder. C) Happily, MATLAB 5 is capable of running both V4 and V5 MEX files. Thus we are currently distributing sources that compile into V4 MEX files. The resulting MEX files run both under V4 and V5. In the future we will drop support for V4 and THINK C. (See note A above.) \ No newline at end of file
diff --git a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/.AppleDouble/.Parent b/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/.AppleDouble/.Parent
deleted file mode 100755
index f242a99..0000000
--- a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/.AppleDouble/.Parent
+++ /dev/null
Binary files differ
diff --git a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/.FBCIndex b/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/.FBCIndex
deleted file mode 100755
index 848736b..0000000
--- a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/.FBCIndex
+++ /dev/null
Binary files differ
diff --git a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/.FBCLockFolder/.FBCSemaphoreFile b/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/.FBCLockFolder/.FBCSemaphoreFile
deleted file mode 100755
index ab2c684..0000000
--- a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/.FBCLockFolder/.FBCSemaphoreFile
+++ /dev/null
Binary files differ
diff --git a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/Makefile-linux b/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/Makefile-linux
deleted file mode 100755
index 726dd31..0000000
--- a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/Makefile-linux
+++ /dev/null
@@ -1,39 +0,0 @@
1MLAB = /usr/local/matlab5.1
2
3MXSFX = mexlx
4MEX = ${MLAB}/bin/mex
5
6MFLAGS = -V4
7INC = -I ${MLAB}/extern/include
8LIB = -L ${MLAB}/extern/lib
9
10CC = gcc -Wall -pedantic
11C_OPTIMIZE_SWITCH = -O2 ## For GCC
12CFLAGS = ${C_OPTIMIZE_SWITCH} ${INC} ${LIB}
13
14all: corrDn.${MXSFX} upConv.${MXSFX} pointOp.${MXSFX} \
15 histo.${MXSFX} range2.${MXSFX}
16
17clean:
18 /bin/rm *.o
19
20corrDn.${MXSFX}: corrDn.o wrap.o convolve.o edges.o
21 ${MEX} ${MFLAGS} corrDn.o wrap.o convolve.o edges.o
22
23upConv.${MXSFX}: upConv.o wrap.o convolve.o edges.o
24 ${MEX} ${MFLAGS} upConv.o wrap.o convolve.o edges.o
25
26pointOp.${MXSFX}: pointOp.o
27 ${MEX} ${MFLAGS} pointOp.o
28
29histo.${MXSFX}: histo.o
30 ${MEX} ${MFLAGS} histo.o
31
32range2.${MXSFX}: range2.o
33 ${MEX} ${MFLAGS} range2.o
34
35convolve.o wrap.o edges.o: convolve.h
36
37%.o : %.c
38 ${CC} -c ${CFLAGS} $<
39
diff --git a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/Makefile-ml6-linux b/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/Makefile-ml6-linux
deleted file mode 100755
index f596ad7..0000000
--- a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/Makefile-ml6-linux
+++ /dev/null
@@ -1,39 +0,0 @@
1MLAB = /usr/local/matlab6.0
2
3MXSFX = mexglx
4MEX = ${MLAB}/bin/mex
5
6MFLAGS = -V4
7INC = -I ${MLAB}/extern/include
8LIB = -L ${MLAB}/extern/lib
9
10CC = gcc -Wall -pedantic
11C_OPTIMIZE_SWITCH = -O2 ## For GCC
12CFLAGS = ${C_OPTIMIZE_SWITCH} ${INC} ${LIB}
13
14all: corrDn.${MXSFX} upConv.${MXSFX} pointOp.${MXSFX} \
15 histo.${MXSFX} range2.${MXSFX}
16
17clean:
18 /bin/rm *.o
19
20corrDn.${MXSFX}: corrDn.o wrap.o convolve.o edges.o
21 ${MEX} ${MFLAGS} corrDn.o wrap.o convolve.o edges.o
22
23upConv.${MXSFX}: upConv.o wrap.o convolve.o edges.o
24 ${MEX} ${MFLAGS} upConv.o wrap.o convolve.o edges.o
25
26pointOp.${MXSFX}: pointOp.o
27 ${MEX} ${MFLAGS} pointOp.o
28
29histo.${MXSFX}: histo.o
30 ${MEX} ${MFLAGS} histo.o
31
32range2.${MXSFX}: range2.o
33 ${MEX} ${MFLAGS} range2.o
34
35convolve.o wrap.o edges.o: convolve.h
36
37%.o : %.c
38 ${CC} -c ${CFLAGS} $<
39
diff --git a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/Makefile-osx b/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/Makefile-osx
deleted file mode 100755
index 352d15b..0000000
--- a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/Makefile-osx
+++ /dev/null
@@ -1,39 +0,0 @@
1MLAB = /share/wotan/matlab13
2
3MXSFX = mexosx
4MEX = ${MLAB}/bin/mex
5
6MFLAGS =
7INC = -I ${MLAB}/extern/include
8LIB = -L ${MLAB}/extern/lib
9
10CC = cc -Wall -pedantic -no-cpp-precomp
11C_OPTIMIZE_SWITCH = -O2 ## For GCC
12CFLAGS = ${C_OPTIMIZE_SWITCH} ${INC} ${LIB}
13
14all: corrDn.${MXSFX} upConv.${MXSFX} pointOp.${MXSFX} \
15 histo.${MXSFX} range2.${MXSFX}
16
17clean:
18 /bin/rm *.o
19
20corrDn.${MXSFX}: corrDn.o wrap.o convolve.o edges.o
21 ${MEX} ${MFLAGS} corrDn.o wrap.o convolve.o edges.o
22
23upConv.${MXSFX}: upConv.o wrap.o convolve.o edges.o
24 ${MEX} ${MFLAGS} upConv.o wrap.o convolve.o edges.o
25
26pointOp.${MXSFX}: pointOp.o
27 ${MEX} ${MFLAGS} pointOp.o
28
29histo.${MXSFX}: histo.o
30 ${MEX} ${MFLAGS} histo.o
31
32range2.${MXSFX}: range2.o
33 ${MEX} ${MFLAGS} range2.o
34
35convolve.o wrap.o edges.o: convolve.h
36
37%.o : %.c
38 ${CC} -c ${CFLAGS} $<
39
diff --git a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/Makefile-solaris b/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/Makefile-solaris
deleted file mode 100755
index 2be2bdb..0000000
--- a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/Makefile-solaris
+++ /dev/null
@@ -1,38 +0,0 @@
1MLAB = /export/home/Solaris2/matlab5.1
2
3MXSFX = mexsol
4MEX = ${MLAB}/bin/mex
5
6MFLAGS = -V4
7INC = -I ${MLAB}/extern/include
8LIB = -L ${MLAB}/extern/lib
9
10CC = gcc -Wall -pedantic
11C_OPTIMIZE_SWITCH = -O2 ## For GCC
12CFLAGS = ${C_OPTIMIZE_SWITCH} ${INC} ${LIB}
13
14all: corrDn.${MXSFX} upConv.${MXSFX} pointOp.${MXSFX} \
15 histo.${MXSFX} range2.${MXSFX}
16
17clean:
18 /bin/rm *.o
19
20corrDn.${MXSFX}: corrDn.o wrap.o convolve.o edges.o
21 ${MEX} ${MFLAGS} corrDn.o wrap.o convolve.o edges.o
22
23upConv.${MXSFX}: upConv.o wrap.o convolve.o edges.o
24 ${MEX} ${MFLAGS} upConv.o wrap.o convolve.o edges.o
25
26pointOp.${MXSFX}: pointOp.o
27 ${MEX} ${MFLAGS} pointOp.o
28
29histo.${MXSFX}: histo.o
30 ${MEX} ${MFLAGS} histo.o
31
32range2.${MXSFX}: range2.o
33 ${MEX} ${MFLAGS} range2.o
34
35convolve.o wrap.o edges.o: convolve.h
36
37%.o : %.c
38 ${CC} -c ${CFLAGS} $<
diff --git a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/Makefile-sun4 b/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/Makefile-sun4
deleted file mode 100755
index 432b181..0000000
--- a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/Makefile-sun4
+++ /dev/null
@@ -1,39 +0,0 @@
1MLAB = /home/alberich/matlab4
2
3MXSFX = mex4
4MEX = ${MLAB}/bin/cmex
5
6MFLAGS =
7INC = -I ${MLAB}/extern/include
8LIB = -L ${MLAB}/extern/lib/sun4
9
10CC = gcc
11C_OPTIMIZE_SWITCH = -O2 ## For GCC
12CFLAGS = ${C_OPTIMIZE_SWITCH} ${INC} ${LIB}
13
14all: corrDn.${MXSFX} upConv.${MXSFX} pointOp.${MXSFX} \
15 histo.${MXSFX} range2.${MXSFX}
16
17clean:
18 /bin/rm *.o
19
20corrDn.${MXSFX}: corrDn.o wrap.o convolve.o edges.o
21 ${MEX} ${MFLAGS} corrDn.o wrap.o convolve.o edges.o
22
23upConv.${MXSFX}: upConv.o wrap.o convolve.o edges.o
24 ${MEX} ${MFLAGS} upConv.o wrap.o convolve.o edges.o
25
26pointOp.${MXSFX}: pointOp.o
27 ${MEX} ${MFLAGS} pointOp.o
28
29histo.${MXSFX}: histo.o
30 ${MEX} ${MFLAGS} histo.o
31
32range2.${MXSFX}: range2.o
33 ${MEX} ${MFLAGS} range2.o
34
35convolve.o wrap.o edges.o: convolve.h
36
37%.o : %.c
38 ${CC} -c ${CFLAGS} $<
39
diff --git a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/convolve.c b/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/convolve.c
deleted file mode 100755
index 60a11a4..0000000
--- a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/convolve.c
+++ /dev/null
@@ -1,325 +0,0 @@
1/*
2;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
3;;; File: convolve.c
4;;; Author: Eero Simoncelli
5;;; Description: General convolution code for 2D images
6;;; Creation Date: Spring, 1987.
7;;; MODIFICATIONS:
8;;; 10/89: approximately optimized the choice of register vars on SPARCS.
9;;; 6/96: Switched array types to double float.
10;;; 2/97: made more robust and readable. Added STOP arguments.
11;;; 8/97: Bug: when calling internal_reduce with edges in {reflect1,repeat,
12;;; extend} and an even filter dimension. Solution: embed the filter
13;;; in the upper-left corner of a filter with odd Y and X dimensions.
14;;; ----------------------------------------------------------------
15;;; Object-Based Vision and Image Understanding System (OBVIUS),
16;;; Copyright 1988, Vision Science Group, Media Laboratory,
17;;; Massachusetts Institute of Technology.
18;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
19*/
20
21#include <stdio.h>
22#include <math.h>
23#include "convolve.h"
24
25/*
26 --------------------------------------------------------------------
27 Correlate FILT with IMAGE, subsampling according to START, STEP, and
28 STOP parameters, with values placed into RESULT array. RESULT
29 dimensions should be ceil((stop-start)/step). TEMP should be a
30 pointer to a temporary double array the size of the filter.
31 EDGES is a string specifying how to handle boundaries -- see edges.c.
32 The convolution is done in 9 sections, where the border sections use
33 specially computed edge-handling filters (see edges.c). The origin
34 of the filter is assumed to be (floor(x_fdim/2), floor(y_fdim/2)).
35------------------------------------------------------------------------ */
36
37/* abstract out the inner product computation */
38#define INPROD(XCNR,YCNR) \
39 { \
40 sum=0.0; \
41 for (im_pos=YCNR*x_dim+XCNR, filt_pos=0, x_filt_stop=x_fdim; \
42 x_filt_stop<=filt_size; \
43 im_pos+=(x_dim-x_fdim), x_filt_stop+=x_fdim) \
44 for (; \
45 filt_pos<x_filt_stop; \
46 filt_pos++, im_pos++) \
47 sum+= image[im_pos]*temp[filt_pos]; \
48 result[res_pos] = sum; \
49 }
50
51int internal_reduce(image, x_dim, y_dim, filt, temp, x_fdim, y_fdim,
52 x_start, x_step, x_stop, y_start, y_step, y_stop,
53 result, edges)
54 register image_type *image, *temp;
55 register int x_fdim, x_dim;
56 register image_type *result;
57 register int x_step, y_step;
58 int x_start, y_start;
59 int x_stop, y_stop;
60 image_type *filt;
61 int y_dim, y_fdim;
62 char *edges;
63 {
64 register double sum;
65 register int filt_pos, im_pos, x_filt_stop;
66 register int x_pos, filt_size = x_fdim*y_fdim;
67 register int y_pos, res_pos;
68 register int y_ctr_stop = y_dim - ((y_fdim==1)?0:y_fdim);
69 register int x_ctr_stop = x_dim - ((x_fdim==1)?0:x_fdim);
70 register int x_res_dim = (x_stop-x_start+x_step-1)/x_step;
71 int x_ctr_start = ((x_fdim==1)?0:1);
72 int y_ctr_start = ((y_fdim==1)?0:1);
73 int x_fmid = x_fdim/2;
74 int y_fmid = y_fdim/2;
75 int base_res_pos;
76 fptr reflect = edge_function(edges); /* look up edge-handling function */
77
78 if (!reflect) return(-1);
79
80 /* shift start/stop coords to filter upper left hand corner */
81 x_start -= x_fmid; y_start -= y_fmid;
82 x_stop -= x_fmid; y_stop -= y_fmid;
83
84 if (x_stop < x_ctr_stop) x_ctr_stop = x_stop;
85 if (y_stop < y_ctr_stop) y_ctr_stop = y_stop;
86
87 for (res_pos=0, y_pos=y_start; /* TOP ROWS */
88 y_pos<y_ctr_start;
89 y_pos+=y_step)
90 {
91 for (x_pos=x_start; /* TOP-LEFT CORNER */
92 x_pos<x_ctr_start;
93 x_pos+=x_step, res_pos++)
94 {
95 (*reflect)(filt,x_fdim,y_fdim,x_pos-1,y_pos-1,temp,REDUCE);
96 INPROD(0,0)
97 }
98
99 (*reflect)(filt,x_fdim,y_fdim,0,y_pos-1,temp,REDUCE);
100 for (; /* TOP EDGE */
101 x_pos<x_ctr_stop;
102 x_pos+=x_step, res_pos++)
103 INPROD(x_pos,0)
104
105 for (; /* TOP-RIGHT CORNER */
106 x_pos<x_stop;
107 x_pos+=x_step, res_pos++)
108 {
109 (*reflect)(filt,x_fdim,y_fdim,x_pos-x_ctr_stop+1,y_pos-1,temp,REDUCE);
110 INPROD(x_ctr_stop,0)
111 }
112 } /* end TOP ROWS */
113
114 y_ctr_start = y_pos; /* hold location of top */
115 for (base_res_pos=res_pos, x_pos=x_start; /* LEFT EDGE */
116 x_pos<x_ctr_start;
117 x_pos+=x_step, base_res_pos++)
118 {
119 (*reflect)(filt,x_fdim,y_fdim,x_pos-1,0,temp,REDUCE);
120 for (y_pos=y_ctr_start, res_pos=base_res_pos;
121 y_pos<y_ctr_stop;
122 y_pos+=y_step, res_pos+=x_res_dim)
123 INPROD(0,y_pos)
124 }
125
126 (*reflect)(filt,x_fdim,y_fdim,0,0,temp,REDUCE);
127 for (; /* CENTER */
128 x_pos<x_ctr_stop;
129 x_pos+=x_step, base_res_pos++)
130 for (y_pos=y_ctr_start, res_pos=base_res_pos;
131 y_pos<y_ctr_stop;
132 y_pos+=y_step, res_pos+=x_res_dim)
133 INPROD(x_pos,y_pos)
134
135 for (; /* RIGHT EDGE */
136 x_pos<x_stop;
137 x_pos+=x_step, base_res_pos++)
138 {
139 (*reflect)(filt,x_fdim,y_fdim,x_pos-x_ctr_stop+1,0,temp,REDUCE);
140 for (y_pos=y_ctr_start, res_pos=base_res_pos;
141 y_pos<y_ctr_stop;
142 y_pos+=y_step, res_pos+=x_res_dim)
143 INPROD(x_ctr_stop,y_pos)
144 }
145
146 for (res_pos-=(x_res_dim-1);
147 y_pos<y_stop; /* BOTTOM ROWS */
148 y_pos+=y_step)
149 {
150 for (x_pos=x_start; /* BOTTOM-LEFT CORNER */
151 x_pos<x_ctr_start;
152 x_pos+=x_step, res_pos++)
153 {
154 (*reflect)(filt,x_fdim,y_fdim,x_pos-1,y_pos-y_ctr_stop+1,temp,REDUCE);
155 INPROD(0,y_ctr_stop)
156 }
157
158 (*reflect)(filt,x_fdim,y_fdim,0,y_pos-y_ctr_stop+1,temp,REDUCE);
159 for (; /* BOTTOM EDGE */
160 x_pos<x_ctr_stop;
161 x_pos+=x_step, res_pos++)
162 INPROD(x_pos,y_ctr_stop)
163
164 for (; /* BOTTOM-RIGHT CORNER */
165 x_pos<x_stop;
166 x_pos+=x_step, res_pos++)
167 {
168 (*reflect)(filt,x_fdim,y_fdim,x_pos-x_ctr_stop+1,y_pos-y_ctr_stop+1,temp,REDUCE);
169 INPROD(x_ctr_stop,y_ctr_stop)
170 }
171 } /* end BOTTOM */
172 return(0);
173 } /* end of internal_reduce */
174
175
176/*
177 --------------------------------------------------------------------
178 Upsample IMAGE according to START,STEP, and STOP parameters and then
179 convolve with FILT, adding values into RESULT array. IMAGE
180 dimensions should be ceil((stop-start)/step). See
181 description of internal_reduce (above).
182
183 WARNING: this subroutine destructively modifies the RESULT array!
184 ------------------------------------------------------------------------ */
185
186/* abstract out the inner product computation */
187#define INPROD2(XCNR,YCNR) \
188 { \
189 val = image[im_pos]; \
190 for (res_pos=YCNR*x_dim+XCNR, filt_pos=0, x_filt_stop=x_fdim; \
191 x_filt_stop<=filt_size; \
192 res_pos+=(x_dim-x_fdim), x_filt_stop+=x_fdim) \
193 for (; \
194 filt_pos<x_filt_stop; \
195 filt_pos++, res_pos++) \
196 result[res_pos] += val*temp[filt_pos]; \
197 }
198
199int internal_expand(image,filt,temp,x_fdim,y_fdim,
200 x_start,x_step,x_stop,y_start,y_step,y_stop,
201 result,x_dim,y_dim,edges)
202 register image_type *result, *temp;
203 register int x_fdim, x_dim;
204 register int x_step, y_step;
205 register image_type *image;
206 int x_start, y_start;
207 image_type *filt;
208 int y_fdim, y_dim;
209 char *edges;
210 {
211 register double val;
212 register int filt_pos, res_pos, x_filt_stop;
213 register int x_pos, filt_size = x_fdim*y_fdim;
214 register int y_pos, im_pos;
215 register int x_ctr_stop = x_dim - ((x_fdim==1)?0:x_fdim);
216 int y_ctr_stop = (y_dim - ((y_fdim==1)?0:y_fdim));
217 int x_ctr_start = ((x_fdim==1)?0:1);
218 int y_ctr_start = ((y_fdim==1)?0:1);
219 int x_fmid = x_fdim/2;
220 int y_fmid = y_fdim/2;
221 int base_im_pos, x_im_dim = (x_stop-x_start+x_step-1)/x_step;
222 fptr reflect = edge_function(edges); /* look up edge-handling function */
223
224 if (!reflect) return(-1);
225
226 /* shift start/stop coords to filter upper left hand corner */
227 x_start -= x_fmid; y_start -= y_fmid;
228 x_stop -= x_fmid; y_stop -= y_fmid;
229
230 if (x_stop < x_ctr_stop) x_ctr_stop = x_stop;
231 if (y_stop < y_ctr_stop) y_ctr_stop = y_stop;
232
233 for (im_pos=0, y_pos=y_start; /* TOP ROWS */
234 y_pos<y_ctr_start;
235 y_pos+=y_step)
236 {
237 for (x_pos=x_start; /* TOP-LEFT CORNER */
238 x_pos<x_ctr_start;
239 x_pos+=x_step, im_pos++)
240 {
241 (*reflect)(filt,x_fdim,y_fdim,x_pos-1,y_pos-1,temp,EXPAND);
242 INPROD2(0,0)
243 }
244
245 (*reflect)(filt,x_fdim,y_fdim,0,y_pos-1,temp,EXPAND);
246 for (; /* TOP EDGE */
247 x_pos<x_ctr_stop;
248 x_pos+=x_step, im_pos++)
249 INPROD2(x_pos,0)
250
251 for (; /* TOP-RIGHT CORNER */
252 x_pos<x_stop;
253 x_pos+=x_step, im_pos++)
254 {
255 (*reflect)(filt,x_fdim,y_fdim,x_pos-x_ctr_stop+1,y_pos-1,temp,EXPAND);
256 INPROD2(x_ctr_stop,0)
257 }
258 } /* end TOP ROWS */
259
260 y_ctr_start = y_pos; /* hold location of top */
261 for (base_im_pos=im_pos, x_pos=x_start; /* LEFT EDGE */
262 x_pos<x_ctr_start;
263 x_pos+=x_step, base_im_pos++)
264 {
265 (*reflect)(filt,x_fdim,y_fdim,x_pos-1,0,temp,EXPAND);
266 for (y_pos=y_ctr_start, im_pos=base_im_pos;
267 y_pos<y_ctr_stop;
268 y_pos+=y_step, im_pos+=x_im_dim)
269 INPROD2(0,y_pos)
270 }
271
272 (*reflect)(filt,x_fdim,y_fdim,0,0,temp,EXPAND);
273 for (; /* CENTER */
274 x_pos<x_ctr_stop;
275 x_pos+=x_step, base_im_pos++)
276 for (y_pos=y_ctr_start, im_pos=base_im_pos;
277 y_pos<y_ctr_stop;
278 y_pos+=y_step, im_pos+=x_im_dim)
279 INPROD2(x_pos,y_pos)
280
281 for (; /* RIGHT EDGE */
282 x_pos<x_stop;
283 x_pos+=x_step, base_im_pos++)
284 {
285 (*reflect)(filt,x_fdim,y_fdim,x_pos-x_ctr_stop+1,0,temp,EXPAND);
286 for (y_pos=y_ctr_start, im_pos=base_im_pos;
287 y_pos<y_ctr_stop;
288 y_pos+=y_step, im_pos+=x_im_dim)
289 INPROD2(x_ctr_stop,y_pos)
290 }
291
292 for (im_pos-=(x_im_dim-1);
293 y_pos<y_stop; /* BOTTOM ROWS */
294 y_pos+=y_step)
295 {
296 for (x_pos=x_start; /* BOTTOM-LEFT CORNER */
297 x_pos<x_ctr_start;
298 x_pos+=x_step, im_pos++)
299 {
300 (*reflect)(filt,x_fdim,y_fdim,x_pos-1,y_pos-y_ctr_stop+1,temp,EXPAND);
301 INPROD2(0,y_ctr_stop)
302 }
303
304 (*reflect)(filt,x_fdim,y_fdim,0,y_pos-y_ctr_stop+1,temp,EXPAND);
305 for (; /* BOTTOM EDGE */
306 x_pos<x_ctr_stop;
307 x_pos+=x_step, im_pos++)
308 INPROD2(x_pos,y_ctr_stop)
309
310 for (; /* BOTTOM-RIGHT CORNER */
311 x_pos<x_stop;
312 x_pos+=x_step, im_pos++)
313 {
314 (*reflect)(filt,x_fdim,y_fdim,x_pos-x_ctr_stop+1,y_pos-y_ctr_stop+1,temp,EXPAND);
315 INPROD2(x_ctr_stop,y_ctr_stop)
316 }
317 } /* end BOTTOM */
318 return(0);
319 } /* end of internal_expand */
320
321
322/* Local Variables: */
323/* buffer-read-only: t */
324/* End: */
325
diff --git a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/convolve.h b/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/convolve.h
deleted file mode 100755
index 48d55f7..0000000
--- a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/convolve.h
+++ /dev/null
@@ -1,55 +0,0 @@
1/*
2;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
3;;; File: convolve.h
4;;; Author: Simoncelli
5;;; Description: Header file for convolve.c
6;;; Creation Date:
7;;; ----------------------------------------------------------------
8;;; Object-Based Vision and Image Understanding System (OBVIUS),
9;;; Copyright 1988, Vision Science Group, Media Laboratory,
10;;; Massachusetts Institute of Technology.
11;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
12*/
13
14#include <stdio.h>
15#include <stdlib.h>
16
17#define ABS(x) (((x)>=0) ? (x) : (-(x)))
18#define ROOT2 1.4142135623730951
19#define REDUCE 0
20#define EXPAND 1
21#define IS ==
22#define ISNT !=
23#define AND &&
24#define OR ||
25
26typedef int (*fptr)();
27
28typedef struct
29 {
30 char *name;
31 fptr func;
32 } EDGE_HANDLER;
33
34typedef double image_type;
35
36fptr edge_function(char *edges);
37int internal_reduce(image_type *image, int x_idim, int y_idim,
38 image_type *filt, image_type *temp, int x_fdim, int y_fdim,
39 int x_start, int x_step, int x_stop,
40 int y_start, int y_step, int y_stop,
41 image_type *result, char *edges);
42int internal_expand(image_type *image,
43 image_type *filt, image_type *temp, int x_fdim, int y_fdim,
44 int x_start, int x_step, int x_stop,
45 int y_start, int y_step, int y_stop,
46 image_type *result, int x_rdim, int y_rdim, char *edges);
47int internal_wrap_reduce(image_type *image, int x_idim, int y_idim,
48 image_type *filt, int x_fdim, int y_fdim,
49 int x_start, int x_step, int x_stop,
50 int y_start, int y_step, int y_stop,
51 image_type *result);
52int internal_wrap_expand(image_type *image, image_type *filt, int x_fdim, int y_fdim,
53 int x_start, int x_step, int x_stop,
54 int y_start, int y_step, int y_stop,
55 image_type *result, int x_rdim, int y_rdim);
diff --git a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/corrDn.c b/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/corrDn.c
deleted file mode 100755
index c74df1f..0000000
--- a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/corrDn.c
+++ /dev/null
@@ -1,145 +0,0 @@
1/*
2RES = corrDn(IM, FILT, EDGES, STEP, START, STOP);
3 >>> See corrDn.m for documentation <<<
4 This is a matlab interface to the internal_reduce function.
5 EPS, 7/96.
6*/
7
8#define V4_COMPAT
9#include <matrix.h> /* Matlab matrices */
10#include <mex.h>
11
12#include "convolve.h"
13
14#define notDblMtx(it) (!mxIsNumeric(it) || !mxIsDouble(it) || mxIsSparse(it) || mxIsComplex(it))
15
16void mexFunction(int nlhs, /* Num return vals on lhs */
17 mxArray *plhs[], /* Matrices on lhs */
18 int nrhs, /* Num args on rhs */
19 const mxArray *prhs[] /* Matrices on rhs */
20 )
21 {
22 double *image,*filt, *temp, *result;
23 int x_fdim, y_fdim, x_idim, y_idim;
24 int x_rdim, y_rdim;
25 int x_start = 1;
26 int x_step = 1;
27 int y_start = 1;
28 int y_step = 1;
29 int x_stop, y_stop;
30 mxArray *arg;
31 double *mxMat;
32 char edges[15] = "reflect1";
33
34 if (nrhs<2) mexErrMsgTxt("requres at least 2 args.");
35
36 /* ARG 1: IMAGE */
37 arg = prhs[0];
38 if notDblMtx(arg) mexErrMsgTxt("IMAGE arg must be a non-sparse double float matrix.");
39 image = mxGetPr(arg);
40 x_idim = (int) mxGetM(arg); /* X is inner index! */
41 y_idim = (int) mxGetN(arg);
42
43 /* ARG 2: FILTER */
44 arg = prhs[1];
45 if notDblMtx(arg) mexErrMsgTxt("FILTER arg must be non-sparse double float matrix.");
46 filt = mxGetPr(arg);
47 x_fdim = (int) mxGetM(arg);
48 y_fdim = (int) mxGetN(arg);
49
50 if ((x_fdim > x_idim) || (y_fdim > y_idim))
51 {
52 mexPrintf("Filter: [%d %d], Image: [%d %d]\n",x_fdim,y_fdim,x_idim,y_idim);
53 mexErrMsgTxt("FILTER dimensions larger than IMAGE dimensions.");
54 }
55
56 /* ARG 3 (optional): EDGES */
57 if (nrhs>2)
58 {
59 if (!mxIsChar(prhs[2]))
60 mexErrMsgTxt("EDGES arg must be a string.");
61 mxGetString(prhs[2],edges,15);
62 }
63
64 /* ARG 4 (optional): STEP */
65 if (nrhs>3)
66 {
67 arg = prhs[3];
68 if notDblMtx(arg) mexErrMsgTxt("STEP arg must be a double float matrix.");
69 if (mxGetM(arg) * mxGetN(arg) != 2)
70 mexErrMsgTxt("STEP arg must contain two elements.");
71 mxMat = mxGetPr(arg);
72 x_step = (int) mxMat[0];
73 y_step = (int) mxMat[1];
74 if ((x_step<1) || (y_step<1))
75 mexErrMsgTxt("STEP values must be greater than zero.");
76 }
77
78 /* ARG 5 (optional): START */
79 if (nrhs>4)
80 {
81 arg = prhs[4];
82 if notDblMtx(arg) mexErrMsgTxt("START arg must be a double float matrix.");
83 if (mxGetM(arg) * mxGetN(arg) != 2)
84 mexErrMsgTxt("START arg must contain two elements.");
85 mxMat = mxGetPr(arg);
86 x_start = (int) mxMat[0];
87 y_start = (int) mxMat[1];
88 if ((x_start<1) || (x_start>x_idim) ||
89 (y_start<1) || (y_start>y_idim))
90 mexErrMsgTxt("START values must lie between 1 and the image dimensions.");
91 }
92 x_start--; /* convert from Matlab to standard C indexes */
93 y_start--;
94
95 /* ARG 6 (optional): STOP */
96 if (nrhs>5)
97 {
98 if notDblMtx(prhs[5]) mexErrMsgTxt("STOP arg must be double float matrix.");
99 if (mxGetM(prhs[5]) * mxGetN(prhs[5]) != 2)
100 mexErrMsgTxt("STOP arg must contain two elements.");
101 mxMat = mxGetPr(prhs[5]);
102 x_stop = (int) mxMat[0];
103 y_stop = (int) mxMat[1];
104 if ((x_stop<x_start) || (x_stop>x_idim) ||
105 (y_stop<y_start) || (y_stop>y_idim))
106 mexErrMsgTxt("STOP values must lie between START and the image dimensions.");
107 }
108 else
109 {
110 x_stop = x_idim;
111 y_stop = y_idim;
112 }
113
114 x_rdim = (x_stop-x_start+x_step-1) / x_step;
115 y_rdim = (y_stop-y_start+y_step-1) / y_step;
116
117 /* mxFreeMatrix(plhs[0]); */
118 plhs[0] = (mxArray *) mxCreateDoubleMatrix(x_rdim,y_rdim,mxREAL);
119 if (plhs[0] == NULL) mexErrMsgTxt("Cannot allocate result matrix");
120 result = mxGetPr(plhs[0]);
121
122 temp = mxCalloc(x_fdim*y_fdim, sizeof(double));
123 if (temp == NULL)
124 mexErrMsgTxt("Cannot allocate necessary temporary space");
125
126 /*
127 printf("i(%d, %d), f(%d, %d), r(%d, %d), X(%d, %d, %d), Y(%d, %d, %d), %s\n",
128 x_idim,y_idim,x_fdim,y_fdim,x_rdim,y_rdim,
129 x_start,x_step,x_stop,y_start,y_step,y_stop,edges);
130 */
131
132 if (strcmp(edges,"circular") == 0)
133 internal_wrap_reduce(image, x_idim, y_idim, filt, x_fdim, y_fdim,
134 x_start, x_step, x_stop, y_start, y_step, y_stop,
135 result);
136 else internal_reduce(image, x_idim, y_idim, filt, temp, x_fdim, y_fdim,
137 x_start, x_step, x_stop, y_start, y_step, y_stop,
138 result, edges);
139
140 mxFree((char *) temp);
141 return;
142 }
143
144
145
diff --git a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/edges-orig.c b/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/edges-orig.c
deleted file mode 100755
index 1f6a98b..0000000
--- a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/edges-orig.c
+++ /dev/null
@@ -1,494 +0,0 @@
1/*
2;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
3;;; File: edges.c
4;;; Author: Eero Simoncelli
5;;; Description: Boundary handling routines for use with convolve.c
6;;; Creation Date: Spring 1987.
7;;; MODIFIED, 6/96, to operate on double float arrays.
8;;; MODIFIED by dgp, 4/1/97, to support THINK C.
9;;; ----------------------------------------------------------------
10;;; Object-Based Vision and Image Understanding System (OBVIUS),
11;;; Copyright 1988, Vision Science Group, Media Laboratory,
12;;; Massachusetts Institute of Technology.
13;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
14*/
15
16/*
17This file contains functions which determine how edges are to be
18handled when performing convolutions of images with linear filters.
19Any edge handling function which is local and linear may be defined,
20except (unfortunately) constants cannot be added. So to treat the
21edges as if the image is surrounded by a gray field, you must paste it
22into a gray image, convolve, and crop it out...
23The main convolution function is called internal_filter and is defined
24in the file convolve.c. The idea is that the convolution function
25calls the edge handling function which computes a new filter based on
26the old filter and the distance to the edge of the image. For
27example, reflection is done by reflecting the filter through the
28appropriate axis and summing. Currently defined functions are listed
29below.
30*/
31
32/*
33#define DEBUG
34*/
35
36#include <stdio.h>
37#include <math.h>
38#include <string.h>
39#include "convolve.h"
40
41#define sgn(a) ( ((a)>0)?1:(((a)<0)?-1:0) )
42#define clip(a,mn,mx) ( ((a)<(mn))?(mn):(((a)>=(mx))?(mx-1):(a)) )
43
44int reflect1(), reflect2(), repeat(), zero(), Extend(), nocompute();
45int ereflect(), predict();
46
47/* Lookup table matching a descriptive string to the edge-handling function */
48#if !THINK_C
49 static EDGE_HANDLER edge_foos[] =
50 {
51 { "dont-compute", nocompute }, /* zero output for filter touching edge */
52 { "zero", zero }, /* zero outside of image */
53 { "repeat", repeat }, /* repeat edge pixel */
54 { "reflect1", reflect1 }, /* reflect about edge pixels */
55 { "reflect2", reflect2 }, /* reflect image, including edge pixels */
56 { "extend", Extend }, /* extend (reflect & invert) */
57 { "predict", predict }, /* predict based on portion covered by filt */
58 { "ereflect", ereflect }, /* orthogonal QMF reflection */
59 };
60#else
61 /*
62 This is really stupid, but THINK C won't allow initialization of static variables in
63 a code resource with string addresses. So we do it this way.
64 The 68K code for a MATLAB 4 MEX file can only be created by THINK C.
65 However, for MATLAB 5, we'll be able to use Metrowerks CodeWarrior for both 68K and PPC, so this
66 cludge can be dropped when we drop support for MATLAB 4.
67 Denis Pelli, 4/1/97.
68 */
69 static EDGE_HANDLER edge_foos[8];
70
71 void InitializeTable(EDGE_HANDLER edge_foos[])
72 {
73 static int i=0;
74
75 if(i>0) return;
76 edge_foos[i].name="dont-compute";
77 edge_foos[i++].func=nocompute;
78 edge_foos[i].name="zero";
79 edge_foos[i++].func=zero;
80 edge_foos[i].name="repeat";
81 edge_foos[i++].func=repeat;
82 edge_foos[i].name="reflect1";
83 edge_foos[i++].func=reflect1;
84 edge_foos[i].name="reflect2";
85 edge_foos[i++].func=reflect2;
86 edge_foos[i].name="extend";
87 edge_foos[i++].func=Extend;
88 edge_foos[i].name="predict";
89 edge_foos[i++].func=predict;
90 edge_foos[i].name="ereflect";
91 edge_foos[i++].func=ereflect;
92 }
93#endif
94
95/*
96Function looks up an edge handler id string in the structure above, and
97returns the associated function
98*/
99fptr edge_function(char *edges)
100 {
101 int i;
102
103#if THINK_C
104 InitializeTable(edge_foos);
105#endif
106 for (i = 0; i<sizeof(edge_foos)/sizeof(EDGE_HANDLER); i++)
107 if (strcmp(edges,edge_foos[i].name) == 0)
108 return(edge_foos[i].func);
109 printf("Error: '%s' is not the name of a valid edge-handler!\n",edges);
110 for (i=0; i<sizeof(edge_foos)/sizeof(EDGE_HANDLER); i++)
111 {
112 if (i==0) printf(" Options are: ");
113 else printf(", ");
114 printf("%s",edge_foos[i].name);
115 }
116 printf("\n");
117 return(0);
118 }
119
120/*
121---------------- EDGE HANDLER ARGUMENTS ------------------------
122filt - array of filter taps.
123x_dim, y_dim - x and y dimensions of filt.
124x_pos - position of filter relative to the horizontal image edges. Negative
125 values indicate left edge, positive indicate right edge. Zero
126 indicates that the filter is not touching either edge. An absolute
127 value of 1 indicates that the edge tap of the filter is over the
128 edge pixel of the image.
129y_pos - analogous to x_pos.
130result - array where the resulting filter will go. The edge
131 of this filter will be aligned with the image for application...
132f_or_e - equal to one of the two constants EXPAND or FILTER.
133--------------------------------------------------------------------
134*/
135
136
137/* --------------------------------------------------------------------
138nocompute() - Return zero for values where filter hangs over the edge.
139*/
140
141int nocompute(filt,x_dim,y_dim,x_pos,y_pos,result,f_or_e)
142 register double *filt, *result;
143 register int x_dim;
144 int y_dim, x_pos, y_pos, f_or_e;
145 {
146 register int i;
147 register int size = x_dim*y_dim;
148
149 if ( (x_pos>1) OR (x_pos<-1) OR (y_pos>1) OR (y_pos<-1) )
150 for (i=0; i<size; i++) result[i] = 0.0;
151 else
152 for (i=0; i<size; i++) result[i] = filt[i];
153 return(0);
154 }
155
156/* --------------------------------------------------------------------
157zero() - Zero outside of image. Discontinuous, but adds zero energy. */
158
159int zero(filt,x_dim,y_dim,x_pos,y_pos,result,f_or_e)
160 register double *filt, *result;
161 register int x_dim;
162 int y_dim, x_pos, y_pos, f_or_e;
163 {
164 register int y_filt,x_filt, y_res,x_res;
165 int filt_sz = x_dim*y_dim;
166 int x_start = ((x_pos>0)?(x_pos-1):((x_pos<0)?(x_pos+1):0));
167 int y_start = x_dim * ((y_pos>0)?(y_pos-1):((y_pos<0)?(y_pos+1):0));
168 int i;
169
170 for (i=0; i<filt_sz; i++) result[i] = 0.0;
171
172 for (y_filt=0, y_res=y_start;
173 y_filt<filt_sz;
174 y_filt+=x_dim, y_res+=x_dim)
175 if ((y_res >= 0) AND (y_res < filt_sz))
176 for (x_filt=y_filt, x_res=x_start;
177 x_filt<y_filt+x_dim;
178 x_filt++, x_res++)
179 if ((x_res >= 0) AND (x_res < x_dim))
180 result[y_res+x_res] = filt[x_filt];
181 return(0);
182 }
183
184/* --------------------------------------------------------------------
185repeat() - repeat edge pixel. Continuous, but content is usually
186different from image.
187*/
188
189int repeat(filt,x_dim,y_dim,x_pos,y_pos,result,f_or_e)
190 register double *filt, *result;
191 register int x_dim;
192 int y_dim, x_pos, y_pos, f_or_e;
193 {
194 register int y_filt,x_filt, y_res,x_res;
195 int filt_sz = x_dim*y_dim;
196 int x_start = ((x_pos>0)?(x_pos-1):((x_pos<0)?(x_pos+1):0));
197 int y_start = x_dim * ((y_pos>0)?(y_pos-1):((y_pos<0)?(y_pos+1):0));
198 int i;
199
200 for (i=0; i<filt_sz; i++) result[i] = 0.0;
201
202 for (y_filt=0, y_res=y_start;
203 y_filt<filt_sz;
204 y_filt+=x_dim, y_res+=x_dim)
205 for (x_filt=y_filt, x_res=x_start;
206 x_filt<y_filt+x_dim;
207 x_filt++, x_res++)
208 result[((y_res>=0)?((y_res<filt_sz)?y_res:(filt_sz-x_dim)):0)
209 + ((x_res>=0)?((x_res<x_dim)?x_res:(x_dim-1)):0)]
210 += filt[x_filt];
211 return(0);
212 }
213
214/* --------------------------------------------------------------------
215reflect2() - "Normal" image reflection. The edge pixel is repeated,
216then the next pixel, etc. Continuous, attempting to maintain
217"similar" content, but discontinuous first derivative.
218*/
219
220int reflect2(filt,x_dim,y_dim,x_pos,y_pos,result,f_or_e)
221 register double *filt, *result;
222 register int x_dim;
223 int y_dim, x_pos, y_pos, f_or_e;
224 {
225 register int y_filt,x_filt, y_edge,x_edge;
226 register int x_base = (x_pos>0)?(x_dim-1):0;
227 register int y_base = (y_pos>0)?(x_dim*(y_dim-1)):0;
228 int filt_sz = x_dim*y_dim;
229 int x_edge_dist = (x_pos>0)?(x_pos-x_dim-1):(x_pos+1);
230 int y_edge_dist = x_dim * ((y_pos>0)?(y_pos-y_dim-1):(y_pos+1));
231 int i;
232
233 #ifdef DEBUG
234 printf("(%d,%d) ",y_pos,x_pos);
235 if (x_pos==0) printf("\n");
236 #endif
237
238 for (i=0; i<filt_sz; i++) result[i] = 0.0;
239
240 for (y_filt=0, y_edge=y_edge_dist;
241 y_filt<filt_sz;
242 y_filt+=x_dim, y_edge+=x_dim)
243 {
244 if (y_edge IS 0) y_edge+=x_dim;
245 for (x_filt=y_filt, x_edge=x_edge_dist;
246 x_filt<y_filt+x_dim;
247 x_filt++, x_edge++)
248 {
249 if (x_edge IS 0) x_edge++;
250 result[ABS(y_base-ABS(y_edge)+x_dim) + ABS(x_base-ABS(x_edge)+1)]
251 += filt[x_filt];
252 }
253 }
254 return(0);
255 }
256
257/* --------------------------------------------------------------------
258reflect1() - Reflection through the edge pixels. This is the right thing
259to do if you are subsampling by 2, since it maintains parity (even
260pixels positions remain even, odd ones remain odd). (note: procedure differs
261depending on f_or_e parameter). */
262
263int reflect1(filt,x_dim,y_dim,x_pos,y_pos,result,f_or_e)
264 register double *filt, *result;
265 register int x_dim;
266 int y_dim, x_pos, y_pos, f_or_e;
267 {
268 int filt_sz = x_dim*y_dim;
269 register int x_start = 0, y_start = 0, x_stop = x_dim, y_stop = filt_sz;
270 register int y_filt,x_filt, y_edge,x_edge;
271 register int x_base = (x_pos>0)?(x_dim-1):0;
272 register int y_base = (y_pos>0)?(x_dim*(y_dim-1)):0;
273 int x_edge_dist = (x_pos>0)?(x_pos-x_dim):((x_pos<0)?(x_pos+1):0);
274 int y_edge_dist = x_dim * ((y_pos>0)?(y_pos-y_dim):((y_pos<0)?(y_pos+1):0));
275 int i;
276 int mx_pos = (x_dim/2)+1;
277 int my_pos = (y_dim/2)+1;
278
279 #ifdef DEBUG
280 printf("(%d,%d) ",y_pos,x_pos);
281 if (x_pos==0) printf("\n");
282 #endif
283
284 for (i=0; i<filt_sz; i++) result[i] = 0.0;
285
286 /* if EXPAND and filter is centered on image edge, do not reflect */
287 if (f_or_e IS EXPAND)
288 {
289 if (x_pos IS mx_pos) x_stop = (x_dim+1)/2;
290 else if (x_pos IS -mx_pos) { x_start = x_dim/2; x_edge_dist = 0; }
291
292 if (y_pos IS my_pos) y_stop = x_dim*((y_dim+1)/2);
293 else if (y_pos IS -my_pos) { y_start = x_dim*(y_dim/2); y_edge_dist = 0;}
294 }
295
296 /* reflect at boundary of image */
297 for (y_filt=y_start, y_edge=y_edge_dist;
298 y_filt<y_stop;
299 y_filt+=x_dim, y_edge+=x_dim)
300 for (x_filt=y_filt+x_start, x_edge=x_edge_dist;
301 x_filt<y_filt+x_stop;
302 x_filt++, x_edge++)
303 result[ABS(y_base-ABS(y_edge)) + ABS(x_base-ABS(x_edge))]
304 += filt[x_filt];
305
306 /* if EXPAND and filter is not centered on image edge, mult edge by 2 */
307 if (f_or_e IS EXPAND)
308 {
309 if ( (ABS(x_pos) ISNT mx_pos) AND (x_pos ISNT 0) )
310 for (y_filt=x_base; y_filt<filt_sz; y_filt+=x_dim)
311 result[y_filt] += result[y_filt];
312 if ( (ABS(y_pos) ISNT my_pos) AND (y_pos ISNT 0) )
313 for (x_filt=y_base; x_filt<y_base+x_dim; x_filt++)
314 result[x_filt] += result[x_filt];
315 }
316 return(0);
317 }
318
319/* --------------------------------------------------------------------
320Extend() - Extend image by reflecting and inverting about edge pixel
321value. Maintains continuity in intensity AND first derivative (but
322not higher derivs).
323*/
324
325int Extend(filt,x_dim,y_dim,x_pos,y_pos,result,f_or_e)
326 register double *filt, *result;
327 register int x_dim;
328 int y_dim, x_pos, y_pos, f_or_e;
329 {
330 int filt_sz = x_dim*y_dim;
331 register int x_start = 0, y_start = 0, x_stop = x_dim, y_stop = filt_sz;
332 register int y_filt,x_filt, y_edge,x_edge;
333 register int x_base = (x_pos>0)?(x_dim-1):0;
334 register int y_base = (y_pos>0)?(x_dim*(y_dim-1)):0;
335 int x_edge_dist = (x_pos>0)?(x_pos-x_dim):((x_pos<-1)?(x_pos+1):0);
336 int y_edge_dist = x_dim * ((y_pos>0)?(y_pos-y_dim):((y_pos<-1)?(y_pos+1):0));
337 int i;
338 int mx_pos = (x_dim/2)+1;
339 int my_pos = (y_dim/2)+1;
340
341 for (i=0; i<filt_sz; i++) result[i] = 0.0;
342
343 /* if EXPAND and filter is centered on image edge, do not reflect */
344 if (f_or_e IS EXPAND)
345 {
346 if (x_pos IS mx_pos) x_stop = (x_dim+1)/2;
347 else if (x_pos IS -mx_pos) { x_start = x_dim/2; x_edge_dist = 0; }
348
349 if (y_pos IS my_pos) y_stop = x_dim*((y_dim+1)/2);
350 else if (y_pos IS -my_pos) { y_start = x_dim*(y_dim/2); y_edge_dist = 0;}
351 }
352
353 /* reflect at boundary of image */
354 for (y_filt=y_start, y_edge=y_edge_dist;
355 y_filt<y_stop;
356 y_filt+=x_dim, y_edge+=x_dim)
357 for (x_filt=y_filt+x_start, x_edge=x_edge_dist;
358 x_filt<y_filt+x_stop;
359 x_filt++, x_edge++)
360 if (((!y_base AND (sgn(y_edge) IS -1)) /* y overhanging */
361 OR
362 (y_base AND (sgn(y_edge) IS 1)))
363 ISNT /* XOR */
364 ((!x_base AND (sgn(x_edge) IS -1)) /* x overhanging */
365 OR
366 (x_base AND (sgn(x_edge) IS 1))))
367 {
368 result[ABS(y_base-ABS(y_edge)) + ABS(x_base-ABS(x_edge))]
369 -= filt[x_filt];
370 result[clip(y_base+y_edge,0,y_dim) + clip(x_base+x_edge,0,x_dim)]
371 += filt[x_filt] + filt[x_filt];
372 }
373 else result[ABS(y_base-ABS(y_edge)) + ABS(x_base-ABS(x_edge))]
374 += filt[x_filt];
375 return(0);
376 }
377
378/* --------------------------------------------------------------------
379predict() - Simple prediction. Like zero, but multiplies the result
380by the reciprocal of the percentage of filter being used. (i.e. if
38150% of the filter is hanging over the edge of the image, multiply the
382taps being used by 2). */
383
384int predict(filt,x_dim,y_dim,x_pos,y_pos,result,f_or_e)
385 register double *filt, *result;
386 register int x_dim;
387 int y_dim, x_pos, y_pos, f_or_e;
388 {
389 register int y_filt,x_filt, y_res,x_res;
390 register double taps_used = 0.0; /* int *** */
391 register double fraction = 0.0;
392 int filt_sz = x_dim*y_dim;
393 int x_start = ((x_pos>0)?(x_pos-1):((x_pos<0)?(x_pos+1):0));
394 int y_start = x_dim * ((y_pos>0)?(y_pos-1):((y_pos<0)?(y_pos+1):0));
395 int i;
396
397 for (i=0; i<filt_sz; i++) result[i] = 0.0;
398
399 for (y_filt=0, y_res=y_start;
400 y_filt<filt_sz;
401 y_filt+=x_dim, y_res+=x_dim)
402 if ((y_res >= 0) AND (y_res < filt_sz))
403 for (x_filt=y_filt, x_res=x_start;
404 x_filt<y_filt+x_dim;
405 x_filt++, x_res++)
406 if ((x_res >= 0) AND (x_res < x_dim))
407 {
408 result[y_res+x_res] = filt[x_filt];
409 taps_used += ABS(filt[x_filt]);
410 }
411 printf("TU: %f\n",taps_used);
412 if (f_or_e IS FILTER)
413 {
414 /* fraction = ( (double) filt_sz ) / ( (double) taps_used ); */
415 for (i=0; i<filt_sz; i++) fraction += ABS(filt[i]);
416 fraction = ( fraction / taps_used );
417 for (i=0; i<filt_sz; i++) result[i] *= fraction;
418 }
419 return(0);
420 }
421
422
423/* --------------------------------------------------------------------
424Reflect, multiplying tap of filter which is at the edge of the image
425by root 2. This maintains orthogonality of odd-length linear-phase
426QMF filters, but it is not useful for most applications, since it
427alters the DC level. */
428
429int ereflect(filt,x_dim,y_dim,x_pos,y_pos,result,f_or_e)
430 register double *filt, *result;
431 register int x_dim;
432 int y_dim, x_pos, y_pos, f_or_e;
433 {
434 register int y_filt,x_filt, y_edge,x_edge;
435 register int x_base = (x_pos>0)?(x_dim-1):0;
436 register int y_base = x_dim * ( (y_pos>0)?(y_dim-1):0 );
437 int filt_sz = x_dim*y_dim;
438 int x_edge_dist = (x_pos>1)?(x_pos-x_dim):((x_pos<-1)?(x_pos+1):0);
439 int y_edge_dist = x_dim * ( (y_pos>1)?(y_pos-y_dim):((y_pos<-1)?(y_pos+1):0) );
440 int i;
441 double norm,onorm;
442
443 for (i=0; i<filt_sz; i++) result[i] = 0.0;
444
445 /* reflect at boundary */
446 for (y_filt=0, y_edge=y_edge_dist;
447 y_filt<filt_sz;
448 y_filt+=x_dim, y_edge+=x_dim)
449 for (x_filt=y_filt, x_edge=x_edge_dist;
450 x_filt<y_filt+x_dim;
451 x_filt++, x_edge++)
452 result[ABS(y_base-ABS(y_edge)) + ABS(x_base-ABS(x_edge))]
453 += filt[x_filt];
454
455 /* now multiply edge by root 2 */
456 if (x_pos ISNT 0)
457 for (y_filt=x_base; y_filt<filt_sz; y_filt+=x_dim)
458 result[y_filt] *= ROOT2;
459 if (y_pos ISNT 0)
460 for (x_filt=y_base; x_filt<y_base+x_dim; x_filt++)
461 result[x_filt] *= ROOT2;
462
463 /* now normalize to norm of original filter */
464 for (norm=0.0,i=0; i<filt_sz; i++)
465 norm += (result[i]*result[i]);
466 norm=sqrt(norm);
467
468 for (onorm=0.0,i=0; i<filt_sz; i++)
469 onorm += (filt[i]*filt[i]);
470 onorm = sqrt(onorm);
471
472 norm = norm/onorm;
473 for (i=0; i<filt_sz; i++)
474 result[i] /= norm;
475 return(0);
476 }
477
478
479/* ------- printout stuff for testing ------------------------------
480 printf("Xpos: %d, Ypos: %d", x_pos, y_pos);
481 for (y_filt=0; y_filt<y_dim; y_filt++)
482 {
483 printf("\n");
484 for (x_filt=0; x_filt<x_dim; x_filt++)
485 printf("%6.1f", result[y_filt*x_dim+x_filt]);
486 }
487 printf("\n");
488*/
489
490
491
492/* Local Variables: */
493/* buffer-read-only: t */
494/* End: */
diff --git a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/edges.c b/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/edges.c
deleted file mode 100755
index 98b377d..0000000
--- a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/edges.c
+++ /dev/null
@@ -1,647 +0,0 @@
1/*
2;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
3;;; File: edges.c
4;;; Author: Eero Simoncelli
5;;; Description: Boundary handling routines for use with convolve.c
6;;; Creation Date: Spring 1987.
7;;; MODIFIED, 6/96, to operate on double float arrays.
8;;; MODIFIED by dgp, 4/1/97, to support THINK C.
9;;; MODIFIED, 8/97: reflect1, reflect2, repeat, extend upgraded to
10;;; work properly for non-symmetric filters. Added qreflect2 to handle
11;;; even-length QMF's which broke under the reflect2 modification.
12;;; ----------------------------------------------------------------
13;;; Object-Based Vision and Image Understanding System (OBVIUS),
14;;; Copyright 1988, Vision Science Group, Media Laboratory,
15;;; Massachusetts Institute of Technology.
16;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
17*/
18
19/* This file contains functions which determine how edges are to be
20handled when performing convolutions of images with linear filters.
21Any edge handling function which is local and linear may be defined,
22except (unfortunately) constants cannot be added. So to treat the
23edges as if the image is surrounded by a gray field, you must paste it
24into a gray image, convolve, and crop it out... The main convolution
25functions are called internal_reduce and internal_expand and are
26defined in the file convolve.c. The idea is that the convolution
27function calls the edge handling function which computes a new filter
28based on the old filter and the distance to the edge of the image.
29For example, reflection is done by reflecting the filter through the
30appropriate axis and summing. Currently defined functions are listed
31below.
32*/
33
34
35#include <stdio.h>
36#include <math.h>
37#include <string.h>
38#include "convolve.h"
39
40#define sgn(a) ( ((a)>0)?1:(((a)<0)?-1:0) )
41#define clip(a,mn,mx) ( ((a)<(mn))?(mn):(((a)>=(mx))?(mx-1):(a)) )
42
43int reflect1(), reflect2(), qreflect2(), repeat(), zero(), Extend(), nocompute();
44int ereflect(), predict();
45
46/* Lookup table matching a descriptive string to the edge-handling function */
47#if !THINK_C
48 static EDGE_HANDLER edge_foos[] =
49 {
50 { "dont-compute", nocompute }, /* zero output for filter touching edge */
51 { "zero", zero }, /* zero outside of image */
52 { "repeat", repeat }, /* repeat edge pixel */
53 { "reflect1", reflect1 }, /* reflect about edge pixels */
54 { "reflect2", reflect2 }, /* reflect image, including edge pixels */
55 { "qreflect2", qreflect2 }, /* reflect image, including edge pixels
56 for even-length QMF decompositions */
57 { "extend", Extend }, /* extend (reflect & invert) */
58 { "ereflect", ereflect }, /* orthogonal QMF reflection */
59 };
60#else
61 /*
62 This is really stupid, but THINK C won't allow initialization of static variables in
63 a code resource with string addresses. So we do it this way.
64 The 68K code for a MATLAB 4 MEX file can only be created by THINK C.
65 However, for MATLAB 5, we'll be able to use Metrowerks CodeWarrior for both 68K and PPC, so this
66 cludge can be dropped when we drop support for MATLAB 4.
67 Denis Pelli, 4/1/97.
68 */
69 static EDGE_HANDLER edge_foos[8];
70
71 void InitializeTable(EDGE_HANDLER edge_foos[])
72 {
73 static int i=0;
74
75 if(i>0) return;
76 edge_foos[i].name="dont-compute";
77 edge_foos[i++].func=nocompute;
78 edge_foos[i].name="zero";
79 edge_foos[i++].func=zero;
80 edge_foos[i].name="repeat";
81 edge_foos[i++].func=repeat;
82 edge_foos[i].name="reflect1";
83 edge_foos[i++].func=reflect1;
84 edge_foos[i].name="reflect2";
85 edge_foos[i++].func=reflect2;
86 edge_foos[i].name="qreflect2";
87 edge_foos[i++].func=qreflect2;
88 edge_foos[i].name="extend";
89 edge_foos[i++].func=Extend;
90 edge_foos[i].name="ereflect";
91 edge_foos[i++].func=ereflect;
92 }
93#endif
94
95/*
96Function looks up an edge handler id string in the structure above, and
97returns the associated function
98*/
99fptr edge_function(char *edges)
100 {
101 int i;
102
103#if THINK_C
104 InitializeTable(edge_foos);
105#endif
106 for (i = 0; i<sizeof(edge_foos)/sizeof(EDGE_HANDLER); i++)
107 if (strcmp(edges,edge_foos[i].name) IS 0)
108 return(edge_foos[i].func);
109 printf("Error: '%s' is not the name of a valid edge-handler!\n",edges);
110 for (i=0; i<sizeof(edge_foos)/sizeof(EDGE_HANDLER); i++)
111 {
112 if (i IS 0) printf(" Options are: ");
113 else printf(", ");
114 printf("%s",edge_foos[i].name);
115 }
116 printf("\n");
117 return(0);
118 }
119
120/*
121---------------- EDGE HANDLER ARGUMENTS ------------------------
122filt - array of filter taps.
123x_dim, y_dim - x and y dimensions of filt.
124x_pos - position of filter relative to the horizontal image edges. Negative
125 values indicate left edge, positive indicate right edge. Zero
126 indicates that the filter is not touching either edge. An absolute
127 value of 1 indicates that the edge tap of the filter is over the
128 edge pixel of the image.
129y_pos - analogous to x_pos.
130result - array where the resulting filter will go. The edge
131 of this filter will be aligned with the image for application...
132r_or_e - equal to one of the two constants EXPAND or REDUCE.
133--------------------------------------------------------------------
134*/
135
136
137/* --------------------------------------------------------------------
138nocompute() - Return zero for values where filter hangs over the edge.
139*/
140
141int nocompute(filt,x_dim,y_dim,x_pos,y_pos,result,r_or_e)
142 register double *filt, *result;
143 register int x_dim;
144 int y_dim, x_pos, y_pos, r_or_e;
145 {
146 register int i;
147 register int size = x_dim*y_dim;
148
149 if ( (x_pos>1) OR (x_pos<-1) OR (y_pos>1) OR (y_pos<-1) )
150 for (i=0; i<size; i++) result[i] = 0.0;
151 else
152 for (i=0; i<size; i++) result[i] = filt[i];
153
154 return(0);
155 }
156
157/* --------------------------------------------------------------------
158zero() - Zero outside of image. Discontinuous, but adds zero energy. */
159
160int zero(filt,x_dim,y_dim,x_pos,y_pos,result,r_or_e)
161 register double *filt, *result;
162 register int x_dim;
163 int y_dim, x_pos, y_pos, r_or_e;
164 {
165 register int y_filt,x_filt, y_res,x_res;
166 int filt_sz = x_dim*y_dim;
167 int x_start = ((x_pos>0)?(x_pos-1):((x_pos<0)?(x_pos+1):0));
168 int y_start = x_dim * ((y_pos>0)?(y_pos-1):((y_pos<0)?(y_pos+1):0));
169 int i;
170
171 for (i=0; i<filt_sz; i++) result[i] = 0.0;
172
173 for (y_filt=0, y_res=y_start;
174 y_filt<filt_sz;
175 y_filt+=x_dim, y_res+=x_dim)
176 if ((y_res >= 0) AND (y_res < filt_sz))
177 for (x_filt=y_filt, x_res=x_start;
178 x_filt<y_filt+x_dim;
179 x_filt++, x_res++)
180 if ((x_res >= 0) AND (x_res < x_dim))
181 result[y_res+x_res] = filt[x_filt];
182 return(0);
183 }
184
185
186/* --------------------------------------------------------------------
187reflect1() - Reflection through the edge pixels. Continuous, but
188discontinuous first derivative. This is the right thing to do if you
189are subsampling by 2, since it maintains parity (even pixels positions
190remain even, odd ones remain odd).
191*/
192
193int reflect1(filt,x_dim,y_dim,x_pos,y_pos,result,r_or_e)
194 register double *filt, *result;
195 register int x_dim;
196 int y_dim, x_pos, y_pos, r_or_e;
197 {
198 int filt_sz = x_dim*y_dim;
199 register int y_filt,x_filt, y_res, x_res;
200 register int x_base = (x_pos>0)?(x_dim-1):0;
201 register int y_base = x_dim * ((y_pos>0)?(y_dim-1):0);
202 int x_overhang = (x_pos>0)?(x_pos-1):((x_pos<0)?(x_pos+1):0);
203 int y_overhang = x_dim * ((y_pos>0)?(y_pos-1):((y_pos<0)?(y_pos+1):0));
204 int i;
205 int mx_pos = (x_pos<0)?(x_dim/2):((x_dim-1)/2);
206 int my_pos = x_dim * ((y_pos<0)?(y_dim/2):((y_dim-1)/2));
207
208 for (i=0; i<filt_sz; i++) result[i] = 0.0;
209
210 if (r_or_e IS REDUCE)
211 for (y_filt=0, y_res=y_overhang-y_base;
212 y_filt<filt_sz;
213 y_filt+=x_dim, y_res+=x_dim)
214 for (x_filt=y_filt, x_res=x_overhang-x_base;
215 x_filt<y_filt+x_dim;
216 x_filt++, x_res++)
217 result[ABS(y_base-ABS(y_res)) + ABS(x_base-ABS(x_res))]
218 += filt[x_filt];
219 else {
220 y_overhang = ABS(y_overhang);
221 x_overhang = ABS(x_overhang);
222 for (y_res=y_base, y_filt = y_base-y_overhang;
223 y_filt > y_base-filt_sz;
224 y_filt-=x_dim, y_res-=x_dim)
225 {
226 for (x_res=x_base, x_filt=x_base-x_overhang;
227 x_filt > x_base-x_dim;
228 x_res--, x_filt--)
229 result[ABS(y_res)+ABS(x_res)] += filt[ABS(y_filt)+ABS(x_filt)];
230 if ((x_overhang ISNT mx_pos) AND (x_pos ISNT 0))
231 for (x_res=x_base, x_filt=x_base-2*mx_pos+x_overhang;
232 x_filt > x_base-x_dim;
233 x_res--, x_filt--)
234 result[ABS(y_res)+ABS(x_res)] += filt[ABS(y_filt)+ABS(x_filt)];
235 }
236 if ((y_overhang ISNT my_pos) AND (y_pos ISNT 0))
237 for (y_res=y_base, y_filt = y_base-2*my_pos+y_overhang;
238 y_filt > y_base-filt_sz;
239 y_filt-=x_dim, y_res-=x_dim)
240 {
241 for (x_res=x_base, x_filt=x_base-x_overhang;
242 x_filt > x_base-x_dim;
243 x_res--, x_filt--)
244 result[ABS(y_res)+ABS(x_res)] += filt[ABS(y_filt)+ABS(x_filt)];
245 if ((x_overhang ISNT mx_pos) AND (x_pos ISNT 0))
246 for (x_res=x_base, x_filt=x_base-2*mx_pos+x_overhang;
247 x_filt > x_base-x_dim;
248 x_res--, x_filt--)
249 result[ABS(y_res)+ABS(x_res)] += filt[ABS(y_filt)+ABS(x_filt)];
250 }
251 }
252
253 return(0);
254 }
255
256/* --------------------------------------------------------------------
257reflect2() - Reflect image at boundary. The edge pixel is repeated,
258then the next pixel, etc. Continuous, but discontinuous first
259derivative.
260*/
261
262int reflect2(filt,x_dim,y_dim,x_pos,y_pos,result,r_or_e)
263 register double *filt, *result;
264 register int x_dim;
265 int y_dim, x_pos, y_pos, r_or_e;
266 {
267 int filt_sz = x_dim*y_dim;
268 register int y_filt,x_filt, y_res, x_res;
269 register int x_base = (x_pos>0)?(x_dim-1):0;
270 register int y_base = x_dim * ((y_pos>0)?(y_dim-1):0);
271 int x_overhang = (x_pos>0)?(x_pos-1):((x_pos<0)?(x_pos+1):0);
272 int y_overhang = x_dim * ((y_pos>0)?(y_pos-1):((y_pos<0)?(y_pos+1):0));
273 int i;
274 int mx_pos = (x_pos<0)?(x_dim/2):((x_dim-1)/2);
275 int my_pos = x_dim * ((y_pos<0)?(y_dim/2):((y_dim-1)/2));
276
277 for (i=0; i<filt_sz; i++) result[i] = 0.0;
278
279 if (r_or_e IS REDUCE)
280 for (y_filt=0, y_res = (y_overhang-y_base) - ((y_pos>0)?x_dim:0);
281 y_filt<filt_sz;
282 y_filt+=x_dim, y_res+=x_dim)
283 {
284 if (y_res IS 0) y_res+=x_dim;
285 for (x_filt=y_filt, x_res = (x_overhang-x_base) - ((x_pos>0)?1:0);
286 x_filt<y_filt+x_dim;
287 x_filt++, x_res++)
288 {
289 if (x_res IS 0) x_res++;
290 result[ABS(y_base-ABS(y_res)+x_dim) + ABS(x_base-ABS(x_res)+1)]
291 += filt[x_filt];
292 }
293 }
294 else {
295 y_overhang = ABS(y_overhang);
296 x_overhang = ABS(x_overhang);
297 for (y_res=y_base, y_filt = y_base-y_overhang;
298 y_filt > y_base-filt_sz;
299 y_filt-=x_dim, y_res-=x_dim)
300 {
301 for (x_res=x_base, x_filt=x_base-x_overhang;
302 x_filt > x_base-x_dim;
303 x_res--, x_filt--)
304 result[ABS(y_res)+ABS(x_res)] += filt[ABS(y_filt)+ABS(x_filt)];
305 if (x_pos ISNT 0)
306 for (x_res=x_base, x_filt=x_base-2*mx_pos+x_overhang-1;
307 x_filt > x_base-x_dim;
308 x_res--, x_filt--)
309 result[ABS(y_res)+ABS(x_res)] += filt[ABS(y_filt)+ABS(x_filt)];
310 }
311 if (y_pos ISNT 0)
312 for (y_res=y_base, y_filt = y_base-2*my_pos+y_overhang-x_dim;
313 y_filt > y_base-filt_sz;
314 y_filt-=x_dim, y_res-=x_dim)
315 {
316 for (x_res=x_base, x_filt=x_base-x_overhang;
317 x_filt > x_base-x_dim;
318 x_res--, x_filt--)
319 result[ABS(y_res)+ABS(x_res)] += filt[ABS(y_filt)+ABS(x_filt)];
320 if (x_pos ISNT 0)
321 for (x_res=x_base, x_filt=x_base-2*mx_pos+x_overhang-1;
322 x_filt > x_base-x_dim;
323 x_res--, x_filt--)
324 result[ABS(y_res)+ABS(x_res)] += filt[ABS(y_filt)+ABS(x_filt)];
325 }
326 }
327
328 return(0);
329 }
330
331
332/* --------------------------------------------------------------------
333qreflect2() - Modified version of reflect2 that works properly for
334even-length QMF filters.
335*/
336
337int qreflect2(filt,x_dim,y_dim,x_pos,y_pos,result,r_or_e)
338 double *filt, *result;
339 int x_dim, y_dim, x_pos, y_pos, r_or_e;
340 {
341 reflect2(filt,x_dim,y_dim,x_pos,y_pos,result,0);
342 return(0);
343 }
344
345/* --------------------------------------------------------------------
346repeat() - repeat edge pixel. Continuous, with discontinuous first
347derivative.
348*/
349
350int repeat(filt,x_dim,y_dim,x_pos,y_pos,result,r_or_e)
351 register double *filt, *result;
352 register int x_dim;
353 int y_dim, x_pos, y_pos, r_or_e;
354 {
355 register int y_filt,x_filt, y_res,x_res, y_tmp, x_tmp;
356 register int x_base = (x_pos>0)?(x_dim-1):0;
357 register int y_base = x_dim * ((y_pos>0)?(y_dim-1):0);
358 int x_overhang = ((x_pos>0)?(x_pos-1):((x_pos<0)?(x_pos+1):0));
359 int y_overhang = x_dim * ((y_pos>0)?(y_pos-1):((y_pos<0)?(y_pos+1):0));
360 int filt_sz = x_dim*y_dim;
361 int mx_pos = (x_dim/2);
362 int my_pos = x_dim * (y_dim/2);
363 int i;
364
365 for (i=0; i<filt_sz; i++) result[i] = 0.0;
366
367 if (r_or_e IS REDUCE)
368 for (y_filt=0, y_res=y_overhang;
369 y_filt<filt_sz;
370 y_filt+=x_dim, y_res+=x_dim)
371 for (x_filt=y_filt, x_res=x_overhang;
372 x_filt<y_filt+x_dim;
373 x_filt++, x_res++)
374 /* Clip the index: */
375 result[((y_res>=0)?((y_res<filt_sz)?y_res:(filt_sz-x_dim)):0)
376 + ((x_res>=0)?((x_res<x_dim)?x_res:(x_dim-1)):0)]
377 += filt[x_filt];
378 else {
379 if ((y_base-y_overhang) ISNT my_pos)
380 for (y_res=y_base, y_filt=y_base-ABS(y_overhang);
381 y_filt > y_base-filt_sz;
382 y_filt-=x_dim, y_res-=x_dim)
383 if ((x_base-x_overhang) ISNT mx_pos)
384 for (x_res=x_base, x_filt=x_base-ABS(x_overhang);
385 x_filt > x_base-x_dim;
386 x_res--, x_filt--)
387 result[ABS(y_res)+ABS(x_res)] += filt[ABS(y_filt)+ABS(x_filt)];
388 else /* ((x_base-x_overhang) IS mx_pos) */
389 for (x_res=x_base, x_filt=x_base-ABS(x_overhang);
390 x_filt > x_base-x_dim;
391 x_filt--, x_res--)
392 for(x_tmp=x_filt; x_tmp > x_base-x_dim; x_tmp--)
393 result[ABS(y_res)+ABS(x_res)] += filt[ABS(y_filt)+ABS(x_tmp)];
394 else /* ((y_base-y_overhang) IS my_pos) */
395 for (y_res=y_base, y_filt=y_base-ABS(y_overhang);
396 y_filt > y_base-filt_sz;
397 y_filt-=x_dim, y_res-=x_dim)
398 for (y_tmp=y_filt; y_tmp > y_base-filt_sz; y_tmp-=x_dim)
399 if ((x_base-x_overhang) ISNT mx_pos)
400 for (x_res=x_base, x_filt=x_base-ABS(x_overhang);
401 x_filt > x_base-x_dim;
402 x_filt--, x_res--)
403 result[ABS(y_res)+ABS(x_res)] += filt[ABS(y_tmp)+ABS(x_filt)];
404 else /* ((x_base-x_overhang) IS mx_pos) */
405 for (x_res=x_base, x_filt=x_base-ABS(x_overhang);
406 x_filt > x_base-x_dim;
407 x_filt--, x_res--)
408 for (x_tmp=x_filt; x_tmp > x_base-x_dim; x_tmp--)
409 result[ABS(y_res)+ABS(x_res)] += filt[ABS(y_tmp)+ABS(x_tmp)];
410 } /* End, if r_or_e IS EXPAND */
411
412 return(0);
413 }
414
415/* --------------------------------------------------------------------
416extend() - Extend image by reflecting and inverting about edge pixel
417value. Maintains continuity in intensity AND first derivative (but
418not higher derivs).
419*/
420
421int Extend(filt,x_dim,y_dim,x_pos,y_pos,result,r_or_e)
422 register double *filt, *result;
423 register int x_dim;
424 int y_dim, x_pos, y_pos, r_or_e;
425 {
426 int filt_sz = x_dim*y_dim;
427 register int y_filt,x_filt, y_res,x_res, y_tmp, x_tmp;
428 register int x_base = (x_pos>0)?(x_dim-1):0;
429 register int y_base = x_dim * ((y_pos>0)?(y_dim-1):0);
430 int x_overhang = (x_pos>0)?(x_pos-1):((x_pos<0)?(x_pos+1):0);
431 int y_overhang = x_dim * ((y_pos>0)?(y_pos-1):((y_pos<0)?(y_pos+1):0));
432 int mx_pos = (x_pos<0)?(x_dim/2):((x_dim-1)/2);
433 int my_pos = x_dim * ((y_pos<0)?(y_dim/2):((y_dim-1)/2));
434 int i;
435
436 for (i=0; i<filt_sz; i++) result[i] = 0.0;
437
438 /* Modeled on reflect1() */
439 if (r_or_e IS REDUCE)
440 for (y_filt=0, y_res=y_overhang;
441 y_filt<filt_sz;
442 y_filt+=x_dim, y_res+=x_dim)
443 if ((y_res>=0) AND (y_res<filt_sz))
444 for (x_filt=y_filt, x_res=x_overhang;
445 x_filt<y_filt+x_dim;
446 x_filt++, x_res++)
447 if ((x_res>=0) AND (x_res<x_dim))
448 result[y_res+x_res] += filt[x_filt];
449 else
450 {
451 result[y_res+ABS(x_base-ABS(x_res-x_base))] -= filt[x_filt];
452 result[y_res+x_base] += 2*filt[x_filt];
453 }
454 else
455 for (x_filt=y_filt, x_res=x_overhang;
456 x_filt<y_filt+x_dim;
457 x_filt++, x_res++)
458 if ((x_res>=0) AND (x_res<x_dim))
459 {
460 result[ABS(y_base-ABS(y_res-y_base))+x_res] -= filt[x_filt];
461 result[y_base+x_res] += 2*filt[x_filt];
462 }
463 else
464 {
465 result[ABS(y_base-ABS(y_res-y_base))+ABS(x_base-ABS(x_res-x_base))]
466 -= filt[x_filt];
467 result[y_base+x_base] += 2*filt[x_filt];
468 }
469 else { /* r_or_e ISNT REDUCE */
470 y_overhang = ABS(y_overhang);
471 x_overhang = ABS(x_overhang);
472 for (y_res=y_base, y_filt = y_base-y_overhang;
473 y_filt > y_base-filt_sz;
474 y_filt-=x_dim, y_res-=x_dim)
475 {
476 for (x_res=x_base, x_filt=x_base-x_overhang;
477 x_filt > x_base-x_dim;
478 x_res--, x_filt--)
479 result[ABS(y_res)+ABS(x_res)] += filt[ABS(y_filt)+ABS(x_filt)];
480 if (x_pos ISNT 0)
481 if (x_overhang ISNT mx_pos)
482 for (x_res=x_base, x_filt=x_base-2*mx_pos+x_overhang;
483 x_filt > x_base-x_dim;
484 x_res--, x_filt--)
485 result[ABS(y_res)+ABS(x_res)] -= filt[ABS(y_filt)+ABS(x_filt)];
486 else /* x_overhang IS mx_pos */
487 for (x_res=x_base, x_filt=x_base-2*mx_pos+x_overhang-1;
488 x_filt > x_base-x_dim;
489 x_res--, x_filt--)
490 for (x_tmp=x_filt; x_tmp > x_base-x_dim; x_tmp--)
491 result[ABS(y_res)+ABS(x_res)] += 2*filt[ABS(y_filt)+ABS(x_tmp)];
492 }
493 if (y_pos ISNT 0)
494 if (y_overhang ISNT my_pos)
495 for (y_res=y_base, y_filt = y_base-2*my_pos+y_overhang;
496 y_filt > y_base-filt_sz;
497 y_filt-=x_dim, y_res-=x_dim)
498 {
499 for (x_res=x_base, x_filt=x_base-x_overhang;
500 x_filt > x_base-x_dim;
501 x_res--, x_filt--)
502 result[ABS(y_res)+ABS(x_res)] -= filt[ABS(y_filt)+ABS(x_filt)];
503 if ((x_pos ISNT 0) AND (x_overhang ISNT mx_pos))
504 for (x_res=x_base, x_filt=x_base-2*mx_pos+x_overhang;
505 x_filt > x_base-x_dim;
506 x_res--, x_filt--)
507 result[ABS(y_res)+ABS(x_res)] -= filt[ABS(y_filt)+ABS(x_filt)];
508 }
509 else /* y_overhang IS my_pos */
510 for (y_res=y_base, y_filt = y_base-2*my_pos+y_overhang-x_dim;
511 y_filt > y_base-filt_sz;
512 y_res-=x_dim, y_filt-=x_dim)
513 for (y_tmp=y_filt; y_tmp > y_base-filt_sz; y_tmp-=x_dim)
514 {
515 for (x_res=x_base, x_filt=x_base-x_overhang;
516 x_filt > x_base-x_dim;
517 x_res--, x_filt--)
518 result[ABS(y_res)+ABS(x_res)] += 2*filt[ABS(y_tmp)+ABS(x_filt)];
519 if ((x_pos ISNT 0) AND (x_overhang IS mx_pos))
520 for (x_res=x_base, x_filt=x_base-2*mx_pos+x_overhang-1;
521 x_filt > x_base-x_dim;
522 x_res--, x_filt--)
523 for (x_tmp=x_filt; x_tmp > x_base-x_dim; x_tmp--)
524 result[ABS(y_res)+ABS(x_res)] += 2*filt[ABS(y_tmp)+ABS(x_tmp)];
525 }
526 } /* r_or_e ISNT REDUCE */
527
528 return(0);
529 }
530
531/* --------------------------------------------------------------------
532predict() - Simple prediction. Like zero, but multiplies the result
533by the reciprocal of the percentage of filter being used. (i.e. if
53450% of the filter is hanging over the edge of the image, multiply the
535taps being used by 2). */
536
537int predict(filt,x_dim,y_dim,x_pos,y_pos,result,r_or_e)
538 register double *filt, *result;
539 register int x_dim;
540 int y_dim, x_pos, y_pos, r_or_e;
541 {
542 register int y_filt,x_filt, y_res,x_res;
543 register double taps_used = 0.0; /* int *** */
544 register double fraction = 0.0;
545 int filt_sz = x_dim*y_dim;
546 int x_start = ((x_pos>0)?(x_pos-1):((x_pos<0)?(x_pos+1):0));
547 int y_start = x_dim * ((y_pos>0)?(y_pos-1):((y_pos<0)?(y_pos+1):0));
548 int i;
549
550 for (i=0; i<filt_sz; i++) result[i] = 0.0;
551
552 for (y_filt=0, y_res=y_start;
553 y_filt<filt_sz;
554 y_filt+=x_dim, y_res+=x_dim)
555 if ((y_res >= 0) AND (y_res < filt_sz))
556 for (x_filt=y_filt, x_res=x_start;
557 x_filt<y_filt+x_dim;
558 x_filt++, x_res++)
559 if ((x_res >= 0) AND (x_res < x_dim))
560 {
561 result[y_res+x_res] = filt[x_filt];
562 taps_used += ABS(filt[x_filt]);
563 }
564
565 if (r_or_e IS REDUCE)
566 {
567 /* fraction = ( (double) filt_sz ) / ( (double) taps_used ); */
568 for (i=0; i<filt_sz; i++) fraction += ABS(filt[i]);
569 fraction = ( fraction / taps_used );
570 for (i=0; i<filt_sz; i++) result[i] *= fraction;
571 }
572 return(0);
573 }
574
575
576/* --------------------------------------------------------------------
577Reflect, multiplying tap of filter which is at the edge of the image
578by root 2. This maintains orthogonality of odd-length linear-phase
579QMF filters, but it is not useful for most applications, since it
580alters the DC level. */
581
582int ereflect(filt,x_dim,y_dim,x_pos,y_pos,result,r_or_e)
583 register double *filt, *result;
584 register int x_dim;
585 int y_dim, x_pos, y_pos, r_or_e;
586 {
587 register int y_filt,x_filt, y_res,x_res;
588 register int x_base = (x_pos>0)?(x_dim-1):0;
589 register int y_base = x_dim * ((y_pos>0)?(y_dim-1):0);
590 int filt_sz = x_dim*y_dim;
591 int x_overhang = (x_pos>1)?(x_pos-x_dim):((x_pos<-1)?(x_pos+1):0);
592 int y_overhang = x_dim * ( (y_pos>1)?(y_pos-y_dim):((y_pos<-1)?(y_pos+1):0) );
593 int i;
594 double norm,onorm;
595
596 for (i=0; i<filt_sz; i++) result[i] = 0.0;
597
598 /* reflect at boundary */
599 for (y_filt=0, y_res=y_overhang;
600 y_filt<filt_sz;
601 y_filt+=x_dim, y_res+=x_dim)
602 for (x_filt=y_filt, x_res=x_overhang;
603 x_filt<y_filt+x_dim;
604 x_filt++, x_res++)
605 result[ABS(y_base-ABS(y_res)) + ABS(x_base-ABS(x_res))]
606 += filt[x_filt];
607
608 /* now multiply edge by root 2 */
609 if (x_pos ISNT 0)
610 for (y_filt=x_base; y_filt<filt_sz; y_filt+=x_dim)
611 result[y_filt] *= ROOT2;
612 if (y_pos ISNT 0)
613 for (x_filt=y_base; x_filt<y_base+x_dim; x_filt++)
614 result[x_filt] *= ROOT2;
615
616 /* now normalize to norm of original filter */
617 for (norm=0.0,i=0; i<filt_sz; i++)
618 norm += (result[i]*result[i]);
619 norm=sqrt(norm);
620
621 for (onorm=0.0,i=0; i<filt_sz; i++)
622 onorm += (filt[i]*filt[i]);
623 onorm = sqrt(onorm);
624
625 norm = norm/onorm;
626 for (i=0; i<filt_sz; i++)
627 result[i] /= norm;
628 return(0);
629 }
630
631
632/* ------- printout stuff for testing ------------------------------
633 printf("Xpos: %d, Ypos: %d", x_pos, y_pos);
634 for (y_filt=0; y_filt<y_dim; y_filt++)
635 {
636 printf("\n");
637 for (x_filt=0; x_filt<x_dim; x_filt++)
638 printf("%6.1f", result[y_filt*x_dim+x_filt]);
639 }
640 printf("\n");
641*/
642
643
644
645/* Local Variables: */
646/* buffer-read-only: t */
647/* End: */
diff --git a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/histo.c b/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/histo.c
deleted file mode 100755
index 43a99c7..0000000
--- a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/histo.c
+++ /dev/null
@@ -1,140 +0,0 @@
1/*
2[N, X] = histo(MTX, NBINS_OR_BINSIZE, BIN_CENTER)
3 >>> See histo.m for documentation <<<
4 EPS, ported from OBVIUS, 3/97.
5*/
6
7#define V4_COMPAT
8#include <matrix.h> /* Matlab matrices */
9#include <mex.h>
10
11#include <stddef.h> /* NULL */
12#include <math.h> /* ceil */
13
14#define notDblMtx(it) (!mxIsNumeric(it) || !mxIsDouble(it) || mxIsSparse(it) || mxIsComplex(it))
15
16#define PAD 0.49999 /* A hair below 1/2, to avoid roundoff errors */
17#define MAXBINS 20000
18
19void mexFunction(int nlhs, /* Num return vals on lhs */
20 mxArray *plhs[], /* Matrices on lhs */
21 int nrhs, /* Num args on rhs */
22 const mxArray *prhs[] /* Matrices on rhs */
23 )
24 {
25 register double temp;
26 register int binnum, i, size;
27 register double *im, binsize;
28 register double origin, *hist, mn, mx, mean;
29 register int nbins;
30 double *bincenters;
31 mxArray *arg;
32 double *mxMat;
33
34 if (nrhs < 1 ) mexErrMsgTxt("requires at least 1 argument.");
35
36 /* ARG 1: MATRIX */
37 arg = prhs[0];
38 if notDblMtx(arg) mexErrMsgTxt("MTX arg must be a real non-sparse matrix.");
39 im = mxGetPr(arg);
40 size = (int) mxGetM(arg) * mxGetN(arg);
41
42 /* FIND min, max, mean values of MTX */
43 mn = *im; mx = *im; binsize = 0;
44 for (i=1; i<size; i++)
45 {
46 temp = im[i];
47 if (temp < mn)
48 mn = temp;
49 else if (temp > mx)
50 mx = temp;
51 binsize += temp;
52 }
53 mean = binsize / size;
54
55 /* ARG 3: BIN_CENTER */
56 if (nrhs > 2)
57 {
58 arg = prhs[2];
59 if notDblMtx(arg) mexErrMsgTxt("BIN_CENTER arg must be a real scalar.");
60 if (mxGetM(arg) * mxGetN(arg) != 1)
61 mexErrMsgTxt("BIN_CENTER must be a real scalar.");
62 mxMat= mxGetPr(arg);
63 origin = *mxMat;
64 }
65 else
66 origin = mean;
67
68 /* ARG 2: If positive, NBINS. If negative, -BINSIZE. */
69 if (nrhs > 1)
70 {
71 arg = prhs[1];
72 if notDblMtx(arg) mexErrMsgTxt("NBINS_OR_BINSIZE arg must be a real scalar.");
73 if (mxGetM(arg) * mxGetN(arg) != 1)
74 mexErrMsgTxt("NBINS_OR_BINSIZE must be a real scalar.");
75 mxMat= mxGetPr(arg);
76 binsize = *mxMat;
77 }
78 else
79 {
80 binsize = 101; /* DEFAULT: 101 bins */
81 }
82
83 /* --------------------------------------------------
84 Adjust origin, binsize, nbins such that
85 mx <= origin + (nbins-1)*binsize + PAD*binsize
86 mn >= origin - PAD*binsize
87 -------------------------------------------------- */
88 if (binsize < 0) /* user specified BINSIZE */
89 {
90 binsize = -binsize;
91 origin -= binsize * ceil((origin-mn-PAD*binsize)/binsize);
92 nbins = (int) ceil((mx-origin-PAD*binsize)/binsize) + 1;
93 }
94 else /* user specified NBINS */
95 {
96 nbins = (int) (binsize + 0.5); /* round to int */
97 if (nbins == 0)
98 mexErrMsgTxt("NBINS must be greater than zero.");
99 binsize = (mx-mn)/(nbins-1+2*PAD); /* start with lower bound */
100 i = ceil((origin-mn-binsize/2)/binsize);
101 if ( mn < (origin-i*binsize-PAD*binsize) )
102 binsize = (origin-mn)/(i+PAD);
103 else if ( mx > (origin+(nbins-1-i)*binsize+PAD*binsize) )
104 binsize = (mx-origin)/((nbins-1-i)+PAD);
105 origin -= binsize * ceil((origin-mn-PAD*binsize)/binsize);
106 }
107
108 if (nbins > MAXBINS)
109 {
110 mexPrintf("nbins: %d, MAXBINS: %d\n",nbins,MAXBINS);
111 mexErrMsgTxt("Number of histo bins has exceeded maximum");
112 }
113
114 /* Allocate hist and xvals */
115 plhs[0] = (mxArray *) mxCreateDoubleMatrix(1,nbins,mxREAL);
116 if (plhs[0] == NULL) mexErrMsgTxt("Error allocating result matrix");
117 hist = mxGetPr(plhs[0]);
118
119 if (nlhs > 1)
120 {
121 plhs[1] = (mxArray *) mxCreateDoubleMatrix(1,nbins,mxREAL);
122 if (plhs[1] == NULL) mexErrMsgTxt("Error allocating result matrix");
123 bincenters = mxGetPr(plhs[1]);
124 for (i=0, temp=origin; i<nbins; i++, temp+=binsize)
125 bincenters[i] = temp;
126 }
127
128 for (i=0; i<size; i++)
129 {
130 binnum = (int) ((im[i] - origin)/binsize + 0.5);
131 if ((binnum < nbins) && (binnum >= 0))
132 (hist[binnum]) += 1.0;
133 else
134 printf("HISTO warning: value %f outside of range [%f,%f]\n",
135 im[i], origin-0.5*binsize, origin+(nbins-0.5)*binsize);
136 }
137
138 return;
139 }
140
diff --git a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/innerProd.c b/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/innerProd.c
deleted file mode 100755
index 8fa1224..0000000
--- a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/innerProd.c
+++ /dev/null
@@ -1,52 +0,0 @@
1/*
2RES = innerProd(MAT);
3 Computes mat'*mat
4 Odelia Schwartz, 8/97.
5*/
6
7#define V4_COMPAT
8#include <matrix.h>
9
10#include <stdio.h>
11#include <ctype.h>
12#include <math.h>
13#include <strings.h>
14#include <stdlib.h>
15
16void mexFunction(int nlhs, /* Num return vals on lhs */
17 mxArray *plhs[], /* Matrices on lhs */
18 int nrhs, /* Num args on rhs */
19 const mxArray *prhs[] /* Matrices on rhs */
20 )
21{
22 register double *res, *mat, tmp;
23 register int len, wid, i, k, j, jlen, ilen, imat, jmat;
24 mxArray *arg;
25
26 /* get matrix input argument */
27 /* should be matrix in which num rows >= num columns */
28 arg=prhs[0];
29 mat= mxGetPr(arg);
30 len = (int) mxGetM(arg);
31 wid = (int) mxGetN(arg);
32 if ( wid > len )
33 printf("innerProd: Warning: width %d is greater than length %d.\n",wid,len);
34 plhs[0] = (mxArray *) mxCreateDoubleMatrix(wid,wid,mxREAL);
35 if (plhs[0] == NULL)
36 mexErrMsgTxt(sprintf("Error allocating %dx%d result matrix",wid,wid));
37 res = mxGetPr(plhs[0]);
38
39 for(i=0, ilen=0; i<wid; i++, ilen+=len)
40 {
41 for(j=i, jlen=ilen; j<wid; j++, jlen+=len)
42 {
43 tmp = 0.0;
44 for(k=0, imat=ilen, jmat=jlen; k<len; k++, imat++, jmat++)
45 tmp += mat[imat]*mat[jmat];
46 res[i*wid+j] = tmp;
47 res[j*wid+i] = tmp;
48 }
49 }
50 return;
51
52}
diff --git a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/innerProd.dll b/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/innerProd.dll
deleted file mode 100755
index 40ac896..0000000
--- a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/innerProd.dll
+++ /dev/null
Binary files differ
diff --git a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/innerProd.mexglx b/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/innerProd.mexglx
deleted file mode 100755
index 749a5b5..0000000
--- a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/innerProd.mexglx
+++ /dev/null
Binary files differ
diff --git a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/innerProd.mexlx b/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/innerProd.mexlx
deleted file mode 100755
index 151b08f..0000000
--- a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/innerProd.mexlx
+++ /dev/null
Binary files differ
diff --git a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/innerProd.mexmac b/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/innerProd.mexmac
deleted file mode 100755
index 6e11fcd..0000000
--- a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/innerProd.mexmac
+++ /dev/null
Binary files differ
diff --git a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/innerProd.mexsol b/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/innerProd.mexsol
deleted file mode 100755
index d761ed8..0000000
--- a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/innerProd.mexsol
+++ /dev/null
Binary files differ
diff --git a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/pointOp.c b/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/pointOp.c
deleted file mode 100755
index 6ffcb45..0000000
--- a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/pointOp.c
+++ /dev/null
@@ -1,126 +0,0 @@
1/*
2RES = pointOp(IM, LUT, ORIGIN, INCREMENT, WARNINGS)
3 >>> See pointOp.m for documentation <<<
4 EPS, ported from OBVIUS, 7/96.
5*/
6
7#define V4_COMPAT
8#include <matrix.h> /* Matlab matrices */
9#include <mex.h>
10
11#include <stddef.h> /* NULL */
12
13#define notDblMtx(it) (!mxIsNumeric(it) || !mxIsDouble(it) || mxIsSparse(it) || mxIsComplex(it))
14
15void internal_pointop();
16
17void mexFunction(int nlhs, /* Num return vals on lhs */
18 mxArray *plhs[], /* Matrices on lhs */
19 int nrhs, /* Num args on rhs */
20 const mxArray *prhs[] /* Matrices on rhs */
21 )
22 {
23 double *image, *lut, *res;
24 double origin, increment;
25 int x_dim, y_dim, lx_dim, ly_dim;
26 int warnings = 1;
27 mxArray *arg;
28 double *mxMat;
29
30 if (nrhs < 4 ) mexErrMsgTxt("requres at least 4 args.");
31
32 /* ARG 1: IMAGE */
33 arg = prhs[0];
34 if notDblMtx(arg) mexErrMsgTxt("IMAGE arg must be a real non-sparse matrix.");
35 image = mxGetPr(arg);
36 x_dim = (int) mxGetM(arg); /* X is inner index! */
37 y_dim = (int) mxGetN(arg);
38
39 /* ARG 2: Lookup table */
40 arg = prhs[1];
41 if notDblMtx(arg) mexErrMsgTxt("LUT arg must be a real non-sparse matrix.");
42 lut = mxGetPr(arg);
43 lx_dim = (int) mxGetM(arg); /* X is inner index! */
44 ly_dim = (int) mxGetN(arg);
45 if ( (lx_dim != 1) && (ly_dim != 1) )
46 mexErrMsgTxt("Lookup table must be a row or column vector.");
47
48 /* ARG 3: ORIGIN */
49 arg = prhs[2];
50 if notDblMtx(arg) mexErrMsgTxt("ORIGIN arg must be a real scalar.");
51 if (mxGetM(arg) * mxGetN(arg) != 1)
52 mexErrMsgTxt("ORIGIN arg must be a real scalar.");
53 mxMat = mxGetPr(arg);
54 origin = *mxMat;
55
56 /* ARG 4: INCREMENT */
57 arg = prhs[3];
58 if notDblMtx(arg) mexErrMsgTxt("INCREMENT arg must be a real scalar.");
59 if (mxGetM(arg) * mxGetN(arg) != 1)
60 mexErrMsgTxt("INCREMENT arg must be a real scalar.");
61 mxMat = mxGetPr(arg);
62 increment = *mxMat;
63
64 /* ARG 5: WARNINGS */
65 if (nrhs>4)
66 {
67 arg = prhs[4];
68 if notDblMtx(arg) mexErrMsgTxt("WARINGS arg must be a real scalar.");
69 if (mxGetM(arg) * mxGetN(arg) != 1)
70 mexErrMsgTxt("WARNINGS arg must be a real scalar.");
71 mxMat = mxGetPr(arg);
72 warnings = (int) *mxMat;
73 }
74
75 plhs[0] = (mxArray *) mxCreateDoubleMatrix(x_dim,y_dim,mxREAL);
76 if (plhs[0] == NULL) mexErrMsgTxt("Cannot allocate result matrix");
77 res = mxGetPr(plhs[0]);
78
79 internal_pointop(image, res, x_dim*y_dim, lut, lx_dim*ly_dim,
80 origin, increment, warnings);
81 return;
82 }
83
84
85/* Use linear interpolation on a lookup table.
86 Taken from OBVIUS. EPS, Spring, 1987.
87 */
88void internal_pointop (im, res, size, lut, lutsize, origin, increment, warnings)
89 register double *im, *res, *lut;
90 register double origin, increment;
91 register int size, lutsize, warnings;
92 {
93 register int i, index;
94 register double pos;
95 register int l_unwarned = warnings;
96 register int r_unwarned = warnings;
97
98 lutsize = lutsize - 2; /* Maximum index value */
99 if (increment > 0)
100 for (i=0; i<size; i++)
101 {
102 pos = (im[i] - origin) / increment;
103 index = (int) pos; /* Floor */
104 if (index < 0)
105 {
106 index = 0;
107 if (l_unwarned)
108 {
109 mexPrintf("Warning: Extrapolating to left of lookup table...\n");
110 l_unwarned = 0;
111 }
112 }
113 else if (index > lutsize)
114 {
115 index = lutsize;
116 if (r_unwarned)
117 {
118 mexPrintf("Warning: Extrapolating to right of lookup table...\n");
119 r_unwarned = 0;
120 }
121 }
122 res[i] = lut[index] + (lut[index+1] - lut[index]) * (pos - index);
123 }
124 else
125 for (i=0; i<size; i++) res[i] = *lut;
126 }
diff --git a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/range2.c b/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/range2.c
deleted file mode 100755
index b84f4e1..0000000
--- a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/range2.c
+++ /dev/null
@@ -1,56 +0,0 @@
1/*
2[MIN, MAX] = range2(MTX)
3 >>> See range2.m for documentation <<<
4 EPS, 3/97.
5*/
6
7#define V4_COMPAT
8#include <matrix.h> /* Matlab matrices */
9#include <mex.h>
10
11#include <stddef.h> /* NULL */
12
13#define notDblMtx(it) (!mxIsNumeric(it) || !mxIsDouble(it) || mxIsSparse(it) || mxIsComplex(it))
14
15void mexFunction(int nlhs, /* Num return vals on lhs */
16 mxArray *plhs[], /* Matrices on lhs */
17 int nrhs, /* Num args on rhs */
18 const mxArray *prhs[] /* Matrices on rhs */
19 )
20 {
21 register double temp, mn, mx;
22 register double *mtx;
23 register int i, size;
24 mxArray *arg;
25
26 if (nrhs != 1) mexErrMsgTxt("requires 1 argument.");
27
28 /* ARG 1: MATRIX */
29 arg = prhs[0];
30 if notDblMtx(arg) mexErrMsgTxt("MTX arg must be a real non-sparse matrix.");
31 mtx = mxGetPr(arg);
32 size = (int) mxGetM(arg) * mxGetN(arg);
33
34 /* FIND min, max values of MTX */
35 mn = *mtx; mx = *mtx;
36 for (i=1; i<size; i++)
37 {
38 temp = mtx[i];
39 if (temp < mn)
40 mn = temp;
41 else if (temp > mx)
42 mx = temp;
43 }
44
45 plhs[0] = (mxArray *) mxCreateDoubleMatrix(1,1,mxREAL);
46 if (plhs[0] == NULL) mexErrMsgTxt("Error allocating result matrix");
47 plhs[1] = (mxArray *) mxCreateDoubleMatrix(1,1,mxREAL);
48 if (plhs[1] == NULL) mexErrMsgTxt("Error allocating result matrix");
49 mtx = mxGetPr(plhs[0]);
50 mtx[0] = mn;
51 mtx = mxGetPr(plhs[1]);
52 mtx[0] = mx;
53
54 return;
55 }
56
diff --git a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/upConv.c b/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/upConv.c
deleted file mode 100755
index 3708f8a..0000000
--- a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/upConv.c
+++ /dev/null
@@ -1,195 +0,0 @@
1/*
2RES = upConv(IM, FILT, EDGES, STEP, START, STOP, RES);
3 >>> See upConv.m for documentation <<<
4 This is a matlab interface to the internal_expand function.
5 EPS, 7/96.
6*/
7
8#define V4_COMPAT
9#include <matrix.h> /* Matlab matrices */
10#include <mex.h>
11
12#include "convolve.h"
13
14#define notDblMtx(it) (!mxIsNumeric(it) || !mxIsDouble(it) || mxIsSparse(it) || mxIsComplex(it))
15
16void mexFunction(int nlhs, /* Num return vals on lhs */
17 mxArray *plhs[], /* Matrices on lhs */
18 int nrhs, /* Num args on rhs */
19 const mxArray *prhs[] /* Matrices on rhs */
20 )
21 {
22 double *image,*filt, *temp, *result, *orig_filt;
23 int x_fdim, y_fdim, x_idim, y_idim;
24 int orig_x = 0, orig_y, x, y;
25 int x_rdim, y_rdim;
26 int x_start = 1;
27 int x_step = 1;
28 int y_start = 1;
29 int y_step = 1;
30 int x_stop, y_stop;
31 mxArray *arg;
32 double *mxMat;
33 char edges[15] = "reflect1";
34
35 if (nrhs<2) mexErrMsgTxt("requres at least 2 args.");
36
37 /* ARG 1: IMAGE */
38 arg = prhs[0];
39 if notDblMtx(arg) mexErrMsgTxt("IMAGE arg must be a non-sparse double float matrix.");
40 image = mxGetPr(arg);
41 x_idim = (int) mxGetM(arg); /* X is inner index! */
42 y_idim = (int) mxGetN(arg);
43
44 /* ARG 2: FILTER */
45 arg = prhs[1];
46 if notDblMtx(arg) mexErrMsgTxt("FILTER arg must be non-sparse double float matrix."); filt = mxGetPr(arg);
47 x_fdim = (int) mxGetM(arg);
48 y_fdim = (int) mxGetN(arg);
49
50 /* ARG 3 (optional): EDGES */
51 if (nrhs>2)
52 {
53 if (!mxIsChar(prhs[2]))
54 mexErrMsgTxt("EDGES arg must be a string.");
55 mxGetString(prhs[2],edges,15);
56 }
57
58 /* ARG 4 (optional): STEP */
59 if (nrhs>3)
60 {
61 arg = prhs[3];
62 if notDblMtx(arg) mexErrMsgTxt("STEP arg must be double float matrix.");
63 if (mxGetM(arg) * mxGetN(arg) != 2)
64 mexErrMsgTxt("STEP arg must contain two elements.");
65 mxMat = mxGetPr(arg);
66 x_step = (int) mxMat[0];
67 y_step = (int) mxMat[1];
68 if ((x_step<1) || (y_step<1))
69 mexErrMsgTxt("STEP values must be greater than zero.");
70 }
71
72 /* ARG 5 (optional): START */
73 if (nrhs>4)
74 {
75 arg = prhs[4];
76 if notDblMtx(arg) mexErrMsgTxt("START arg must be double float matrix.");
77 if (mxGetM(arg) * mxGetN(arg) != 2)
78 mexErrMsgTxt("START arg must contain two elements.");
79 mxMat = mxGetPr(arg);
80 x_start = (int) mxMat[0];
81 y_start = (int) mxMat[1];
82 if ((x_start<1) || (y_start<1))
83 mexErrMsgTxt("START values must be greater than zero.");
84 }
85 x_start--; /* convert to standard C indexes */
86 y_start--;
87
88 /* ARG 6 (optional): STOP */
89 if (nrhs>5)
90 {
91 if notDblMtx(prhs[5]) mexErrMsgTxt("STOP arg must be double float matrix.");
92 if (mxGetM(prhs[5]) * mxGetN(prhs[5]) != 2)
93 mexErrMsgTxt("STOP arg must contain two elements.");
94 mxMat = mxGetPr(prhs[5]);
95 x_stop = (int) mxMat[0];
96 y_stop = (int) mxMat[1];
97 if ((x_stop<x_start) || (y_stop<y_start))
98 mexErrMsgTxt("STOP values must be greater than START values.");
99 }
100 else
101 { /* default: make res dims a multiple of STEP size */
102 x_stop = x_step * ((x_start/x_step) + x_idim);
103 y_stop = y_step * ((y_start/y_step) + y_idim);
104 }
105
106 /* ARG 6 (optional): RESULT image */
107 if (nrhs>6)
108 {
109 arg = prhs[6];
110 if notDblMtx(arg) mexErrMsgTxt("RES arg must be double float matrix.");
111
112 /* 7/10/97: Returning one of the args causes problems with Matlab's memory
113 manager, so we don't return anything if the result image is passed */
114 /* plhs[0] = arg; */
115 result = mxGetPr(arg);
116 x_rdim = (int) mxGetM(arg); /* X is inner index! */
117 y_rdim = (int) mxGetN(arg);
118 if ((x_stop>x_rdim) || (y_stop>y_rdim))
119 mexErrMsgTxt("STOP values must within image dimensions.");
120 }
121 else
122 {
123 x_rdim = x_stop;
124 y_rdim = y_stop;
125 /* x_rdim = x_step * ((x_stop+x_step-1)/x_step);
126 y_rdim = y_step * ((y_stop+y_step-1)/y_step); */
127
128 plhs[0] = (mxArray *) mxCreateDoubleMatrix(x_rdim,y_rdim,mxREAL);
129 if (plhs[0] == NULL) mexErrMsgTxt("Cannot allocate result matrix");
130 result = mxGetPr(plhs[0]);
131 }
132
133 if ( (((x_stop-x_start+x_step-1) / x_step) != x_idim) ||
134 (((y_stop-y_start+y_step-1) / y_step) != y_idim) )
135 {
136 mexPrintf("Im dims: [%d %d]\n",x_idim,y_idim);
137 mexPrintf("Start: [%d %d]\n",x_start,y_start);
138 mexPrintf("Step: [%d %d]\n",x_step,y_step);
139 mexPrintf("Stop: [%d %d]\n",x_stop,y_stop);
140 mexPrintf("Res dims: [%d %d]\n",x_rdim,y_rdim);
141 mexErrMsgTxt("Image sizes and upsampling args are incompatible!");
142 }
143
144 /* upConv has a bug for even-length kernels when using the
145 reflect1, extend, or repeat edge-handlers */
146 if ((!strcmp(edges,"reflect1") || !strcmp(edges,"extend") || !strcmp(edges,"repeat"))
147 &&
148 ((x_fdim%2 == 0) || (y_fdim%2 == 0)))
149 {
150 orig_filt = filt;
151 orig_x = x_fdim;
152 orig_y = y_fdim;
153 x_fdim = 2*(orig_x/2)+1;
154 y_fdim = 2*(orig_y/2)+1;
155 filt = mxCalloc(x_fdim*y_fdim, sizeof(double));
156 if (filt == NULL)
157 mexErrMsgTxt("Cannot allocate necessary temporary space");
158 for (y=0; y<orig_y; y++)
159 for (x=0; x<orig_x; x++)
160 filt[y*x_fdim + x] = orig_filt[y*orig_x + x];
161 }
162
163 if ((x_fdim > x_rdim) || (y_fdim > y_rdim))
164 {
165 mexPrintf("Filter: [%d %d], ",x_fdim,y_fdim);
166 mexPrintf("Result: [%d %d]\n",x_rdim,y_rdim);
167 mexErrMsgTxt("FILTER dimensions larger than RESULT dimensions.");
168 }
169
170 temp = mxCalloc(x_fdim*y_fdim, sizeof(double));
171 if (temp == NULL)
172 mexErrMsgTxt("Cannot allocate necessary temporary space");
173
174 /*
175 printf("(%d, %d), (%d, %d), (%d, %d), (%d, %d), (%d, %d), %s\n",
176 x_idim,y_idim,x_fdim,y_fdim,x_rdim,y_rdim,
177 x_start,x_step,y_start,y_step,edges);
178 */
179
180 if (strcmp(edges,"circular") == 0)
181 internal_wrap_expand(image, filt, x_fdim, y_fdim,
182 x_start, x_step, x_stop, y_start, y_step, y_stop,
183 result, x_rdim, y_rdim);
184 else internal_expand(image, filt, temp, x_fdim, y_fdim,
185 x_start, x_step, x_stop, y_start, y_step, y_stop,
186 result, x_rdim, y_rdim, edges);
187
188 if (orig_x) mxFree((char *) filt);
189 mxFree((char *) temp);
190
191 return;
192 }
193
194
195
diff --git a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/wrap.c b/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/wrap.c
deleted file mode 100755
index a081123..0000000
--- a/SD-VBS/benchmarks/texture_synthesis/src/matlabPyrTools/MEX/wrap.c
+++ /dev/null
@@ -1,281 +0,0 @@
1/*
2;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
3;;; File: wrap.c
4;;; Author: Eero Simoncelli
5;;; Description: Circular convolution on 2D images.
6;;; Creation Date: Spring, 1987.
7;;; MODIFICATIONS:
8;;; 6/96: Switched array types to double float.
9;;; 2/97: made more robust and readable. Added STOP arguments.
10;;; ----------------------------------------------------------------
11;;; Object-Based Vision and Image Understanding System (OBVIUS),
12;;; Copyright 1988, Vision Science Group, Media Laboratory,
13;;; Massachusetts Institute of Technology.
14;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
15*/
16
17#include <stdlib.h>
18
19#include "convolve.h"
20
21/*
22 --------------------------------------------------------------------
23 Performs correlation (i.e., convolution with filt(-x,-y)) of FILT
24 with IMAGE followed by subsampling (a.k.a. REDUCE in Burt&Adelson81).
25 The operations are combined to avoid unnecessary computation of the
26 convolution samples that are to be discarded in the subsampling
27 operation. The convolution is done in 9 sections so that mod
28 operations are not performed unnecessarily. The subsampling lattice
29 is specified by the START, STEP and STOP parameters.
30 -------------------------------------------------------------------- */
31
32/* abstract out the inner product computation */
33#define INPROD(YSTART,YIND,XSTART,XIND) \
34 { \
35 sum=0.0; \
36 for (y_im=YSTART, filt_pos=0, x_filt_stop=x_fdim; \
37 x_filt_stop<=filt_size; \
38 y_im++, x_filt_stop+=x_fdim) \
39 for (x_im=XSTART ; \
40 filt_pos<x_filt_stop; \
41 filt_pos++, x_im++) \
42 sum += imval[YIND][XIND] * filt[filt_pos]; \
43 result[res_pos] = sum; \
44 }
45
46int internal_wrap_reduce(image, x_dim, y_dim, filt, x_fdim, y_fdim,
47 x_start, x_step, x_stop, y_start, y_step, y_stop,
48 result)
49 register image_type *filt, *result;
50 register int x_dim, y_dim, x_fdim, y_fdim;
51 image_type *image;
52 int x_start, x_step, x_stop, y_start, y_step, y_stop;
53 {
54 register double sum;
55 register int filt_size = x_fdim*y_fdim;
56 image_type **imval;
57 register int filt_pos, x_im, y_im, x_filt_stop;
58 register int x_pos, y_pos, res_pos;
59 int x_ctr_stop = x_dim - x_fdim + 1;
60 int y_ctr_stop = y_dim - y_fdim + 1;
61 int x_ctr_start = 0;
62 int y_ctr_start = 0;
63 int x_fmid = x_fdim/2;
64 int y_fmid = y_fdim/2;
65
66 /* shift start/stop coords to filter upper left hand corner */
67 x_start -= x_fmid; y_start -= y_fmid;
68 x_stop -= x_fmid; y_stop -= y_fmid;
69
70 if (x_stop < x_ctr_stop) x_ctr_stop = x_stop;
71 if (y_stop < y_ctr_stop) y_ctr_stop = y_stop;
72
73 /* Set up pointer array for rows */
74 imval = (image_type **) malloc(y_dim*sizeof(image_type *));
75 if (imval IS NULL)
76 {
77 printf("INTERNAL_WRAP: Failed to allocate temp array!");
78 return(-1);
79 }
80 for (y_pos=y_im=0;y_pos<y_dim;y_pos++,y_im+=x_dim)
81 imval[y_pos] = (image+y_im);
82
83 for (res_pos=0, y_pos=y_start; /* TOP ROWS */
84 y_pos<y_ctr_start;
85 y_pos+=y_step)
86 {
87 for (x_pos=x_start;
88 x_pos<x_ctr_start;
89 x_pos+=x_step, res_pos++)
90 INPROD(y_pos+y_dim, y_im%y_dim, x_pos+x_dim, x_im%x_dim)
91
92 for (;
93 x_pos<x_ctr_stop;
94 x_pos+=x_step, res_pos++)
95 INPROD(y_pos+y_dim, y_im%y_dim, x_pos, x_im)
96
97 for (;
98 x_pos<x_stop;
99 x_pos+=x_step, res_pos++)
100 INPROD(y_pos+y_dim, y_im%y_dim, x_pos, x_im%x_dim)
101 } /* end TOP ROWS */
102
103 for (; /* MID ROWS */
104 y_pos<y_ctr_stop;
105 y_pos+=y_step)
106 {
107 for (x_pos=x_start;
108 x_pos<x_ctr_start;
109 x_pos+=x_step, res_pos++)
110 INPROD(y_pos, y_im, x_pos+x_dim, x_im%x_dim)
111
112 for (; /* CENTER SECTION */
113 x_pos<x_ctr_stop;
114 x_pos+=x_step, res_pos++)
115 INPROD(y_pos, y_im, x_pos, x_im)
116
117 for (;
118 x_pos<x_stop;
119 x_pos+=x_step, res_pos++)
120 INPROD(y_pos, y_im, x_pos, x_im%x_dim)
121 } /* end MID ROWS */
122
123 for (; /* BOTTOM ROWS */
124 y_pos<y_stop;
125 y_pos+=y_step)
126 {
127 for (x_pos=x_start;
128 x_pos<x_ctr_start;
129 x_pos+=x_step, res_pos++)
130 INPROD(y_pos, y_im%y_dim, x_pos+x_dim, x_im%x_dim)
131
132 for (;
133 x_pos<x_ctr_stop;
134 x_pos+=x_step, res_pos++)
135 INPROD(y_pos, y_im%y_dim, x_pos, x_im)
136
137 for (;
138 x_pos<x_stop;
139 x_pos+=x_step, res_pos++)
140 INPROD(y_pos, y_im%y_dim, x_pos, x_im%x_dim)
141 } /* end BOTTOM ROWS */
142
143 free ((image_type **) imval);
144
145 return(0);
146 } /* end of internal_wrap_reduce */
147
148
149
150/*
151 --------------------------------------------------------------------
152 Performs upsampling (padding with zeros) followed by convolution of
153 FILT with IMAGE (a.k.a. EXPAND in Burt&Adelson81). The operations
154 are combined to avoid unnecessary multiplication of filter samples
155 with zeros in the upsampled image. The convolution is done in 9
156 sections so that mod operation is not performed unnecessarily.
157 Arguments are described in the comment above internal_wrap_reduce.
158
159 WARNING: this subroutine destructively modifes the RESULT image, so
160 the user must zero the result before invocation!
161 -------------------------------------------------------------------- */
162
163/* abstract out the inner product computation */
164#define INPROD2(YSTART,YIND,XSTART,XIND) \
165 { \
166 val = image[im_pos]; \
167 for (y_res=YSTART, filt_pos=0, x_filt_stop=x_fdim; \
168 x_filt_stop<=filt_size; \
169 y_res++, x_filt_stop+=x_fdim) \
170 for (x_res=XSTART; \
171 filt_pos<x_filt_stop; \
172 filt_pos++, x_res++) \
173 imval[YIND][XIND] += val * filt[filt_pos]; \
174 }
175
176int internal_wrap_expand(image, filt, x_fdim, y_fdim,
177 x_start, x_step, x_stop, y_start, y_step, y_stop,
178 result, x_dim, y_dim)
179 register image_type *filt, *result;
180 register int x_fdim, y_fdim, x_dim, y_dim;
181 image_type *image;
182 int x_start, x_step, x_stop, y_start, y_step, y_stop;
183 {
184 register double val;
185 register int filt_size = x_fdim*y_fdim;
186 image_type **imval;
187 register int filt_pos, x_res, y_res, x_filt_stop;
188 register int x_pos, y_pos, im_pos;
189 int x_ctr_stop = x_dim - x_fdim + 1;
190 int y_ctr_stop = y_dim - y_fdim + 1;
191 int x_ctr_start = 0;
192 int y_ctr_start = 0;
193 int x_fmid = x_fdim/2;
194 int y_fmid = y_fdim/2;
195
196 /* shift start/stop coords to filter upper left hand corner */
197 x_start -= x_fmid; y_start -= y_fmid;
198 x_stop -= x_fmid; y_stop -= y_fmid;
199
200 if (x_stop < x_ctr_stop) x_ctr_stop = x_stop;
201 if (y_stop < y_ctr_stop) y_ctr_stop = y_stop;
202
203 /* Set up pointer array for rows */
204 imval = (image_type **) malloc(y_dim*sizeof(image_type *));
205 if (imval IS NULL)
206 {
207 printf("INTERNAL_WRAP: Failed to allocate temp array!");
208 return(-1);
209 }
210 for (y_pos=y_res=0;y_pos<y_dim;y_pos++,y_res+=x_dim)
211 imval[y_pos] = (result+y_res);
212
213 for (im_pos=0, y_pos=y_start; /* TOP ROWS */
214 y_pos<y_ctr_start;
215 y_pos+=y_step)
216 {
217 for (x_pos=x_start;
218 x_pos<x_ctr_start;
219 x_pos+=x_step, im_pos++)
220 INPROD2(y_pos+y_dim, y_res%y_dim, x_pos+x_dim, x_res%x_dim)
221
222 for (;
223 x_pos<x_ctr_stop;
224 x_pos+=x_step, im_pos++)
225 INPROD2(y_pos+y_dim, y_res%y_dim, x_pos, x_res)
226
227 for (;
228 x_pos<x_stop;
229 x_pos+=x_step, im_pos++)
230 INPROD2(y_pos+y_dim, y_res%y_dim, x_pos, x_res%x_dim)
231 } /* end TOP ROWS */
232
233 for (; /* MID ROWS */
234 y_pos<y_ctr_stop;
235 y_pos+=y_step)
236 {
237 for (x_pos=x_start;
238 x_pos<x_ctr_start;
239 x_pos+=x_step, im_pos++)
240 INPROD2(y_pos, y_res, x_pos+x_dim, x_res%x_dim)
241
242 for (; /* CENTER SECTION */
243 x_pos<x_ctr_stop;
244 x_pos+=x_step, im_pos++)
245 INPROD2(y_pos, y_res, x_pos, x_res)
246
247 for (;
248 x_pos<x_stop;
249 x_pos+=x_step, im_pos++)
250 INPROD2(y_pos, y_res, x_pos, x_res%x_dim)
251 } /* end MID ROWS */
252
253 for (; /* BOTTOM ROWS */
254 y_pos<y_stop;
255 y_pos+=y_step)
256 {
257 for (x_pos=x_start;
258 x_pos<x_ctr_start;
259 x_pos+=x_step, im_pos++)
260 INPROD2(y_pos, y_res%y_dim, x_pos+x_dim, x_res%x_dim)
261
262 for (;
263 x_pos<x_ctr_stop;
264 x_pos+=x_step, im_pos++)
265 INPROD2(y_pos, y_res%y_dim, x_pos, x_res)
266
267 for (;
268 x_pos<x_stop;
269 x_pos+=x_step, im_pos++)
270 INPROD2(y_pos, y_res%y_dim, x_pos, x_res%x_dim)
271 } /* end BOTTOM ROWS */
272
273 free ((image_type **) imval);
274 return(0);
275 } /* end of internal_wrap_expand */
276
277
278
279/* Local Variables: */
280/* buffer-read-only: t */
281/* End: */
generated by cgit v1.2.2 (git 2.25.0) at 2025-04-30 19:29:44 -0400