123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281 |
- /*
- ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
- ;;; File: wrap.c
- ;;; Author: Eero Simoncelli
- ;;; Description: Circular convolution on 2D images.
- ;;; Creation Date: Spring, 1987.
- ;;; MODIFICATIONS:
- ;;; 6/96: Switched array types to double float.
- ;;; 2/97: made more robust and readable. Added STOP arguments.
- ;;; ----------------------------------------------------------------
- ;;; Object-Based Vision and Image Understanding System (OBVIUS),
- ;;; Copyright 1988, Vision Science Group, Media Laboratory,
- ;;; Massachusetts Institute of Technology.
- ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
- */
- #include <stdlib.h>
- #include "convolve.h"
- /*
- --------------------------------------------------------------------
- Performs correlation (i.e., convolution with filt(-x,-y)) of FILT
- with IMAGE followed by subsampling (a.k.a. REDUCE in Burt&Adelson81).
- The operations are combined to avoid unnecessary computation of the
- convolution samples that are to be discarded in the subsampling
- operation. The convolution is done in 9 sections so that mod
- operations are not performed unnecessarily. The subsampling lattice
- is specified by the START, STEP and STOP parameters.
- -------------------------------------------------------------------- */
- /* abstract out the inner product computation */
- #define INPROD(YSTART,YIND,XSTART,XIND) \
- { \
- sum=0.0; \
- for (y_im=YSTART, filt_pos=0, x_filt_stop=x_fdim; \
- x_filt_stop<=filt_size; \
- y_im++, x_filt_stop+=x_fdim) \
- for (x_im=XSTART ; \
- filt_pos<x_filt_stop; \
- filt_pos++, x_im++) \
- sum += imval[YIND][XIND] * filt[filt_pos]; \
- result[res_pos] = sum; \
- }
- int internal_wrap_reduce(image, x_dim, y_dim, filt, x_fdim, y_fdim,
- x_start, x_step, x_stop, y_start, y_step, y_stop,
- result)
- register image_type *filt, *result;
- register int x_dim, y_dim, x_fdim, y_fdim;
- image_type *image;
- int x_start, x_step, x_stop, y_start, y_step, y_stop;
- {
- register double sum;
- register int filt_size = x_fdim*y_fdim;
- image_type **imval;
- register int filt_pos, x_im, y_im, x_filt_stop;
- register int x_pos, y_pos, res_pos;
- int x_ctr_stop = x_dim - x_fdim + 1;
- int y_ctr_stop = y_dim - y_fdim + 1;
- int x_ctr_start = 0;
- int y_ctr_start = 0;
- int x_fmid = x_fdim/2;
- int y_fmid = y_fdim/2;
-
- /* shift start/stop coords to filter upper left hand corner */
- x_start -= x_fmid; y_start -= y_fmid;
- x_stop -= x_fmid; y_stop -= y_fmid;
- if (x_stop < x_ctr_stop) x_ctr_stop = x_stop;
- if (y_stop < y_ctr_stop) y_ctr_stop = y_stop;
- /* Set up pointer array for rows */
- imval = (image_type **) malloc(y_dim*sizeof(image_type *));
- if (imval IS NULL)
- {
- printf("INTERNAL_WRAP: Failed to allocate temp array!");
- return(-1);
- }
- for (y_pos=y_im=0;y_pos<y_dim;y_pos++,y_im+=x_dim)
- imval[y_pos] = (image+y_im);
-
- for (res_pos=0, y_pos=y_start; /* TOP ROWS */
- y_pos<y_ctr_start;
- y_pos+=y_step)
- {
- for (x_pos=x_start;
- x_pos<x_ctr_start;
- x_pos+=x_step, res_pos++)
- INPROD(y_pos+y_dim, y_im%y_dim, x_pos+x_dim, x_im%x_dim)
- for (;
- x_pos<x_ctr_stop;
- x_pos+=x_step, res_pos++)
- INPROD(y_pos+y_dim, y_im%y_dim, x_pos, x_im)
- for (;
- x_pos<x_stop;
- x_pos+=x_step, res_pos++)
- INPROD(y_pos+y_dim, y_im%y_dim, x_pos, x_im%x_dim)
- } /* end TOP ROWS */
-
- for (; /* MID ROWS */
- y_pos<y_ctr_stop;
- y_pos+=y_step)
- {
- for (x_pos=x_start;
- x_pos<x_ctr_start;
- x_pos+=x_step, res_pos++)
- INPROD(y_pos, y_im, x_pos+x_dim, x_im%x_dim)
- for (; /* CENTER SECTION */
- x_pos<x_ctr_stop;
- x_pos+=x_step, res_pos++)
- INPROD(y_pos, y_im, x_pos, x_im)
- for (;
- x_pos<x_stop;
- x_pos+=x_step, res_pos++)
- INPROD(y_pos, y_im, x_pos, x_im%x_dim)
- } /* end MID ROWS */
-
- for (; /* BOTTOM ROWS */
- y_pos<y_stop;
- y_pos+=y_step)
- {
- for (x_pos=x_start;
- x_pos<x_ctr_start;
- x_pos+=x_step, res_pos++)
- INPROD(y_pos, y_im%y_dim, x_pos+x_dim, x_im%x_dim)
- for (;
- x_pos<x_ctr_stop;
- x_pos+=x_step, res_pos++)
- INPROD(y_pos, y_im%y_dim, x_pos, x_im)
- for (;
- x_pos<x_stop;
- x_pos+=x_step, res_pos++)
- INPROD(y_pos, y_im%y_dim, x_pos, x_im%x_dim)
- } /* end BOTTOM ROWS */
- free ((image_type **) imval);
- return(0);
- } /* end of internal_wrap_reduce */
- /*
- --------------------------------------------------------------------
- Performs upsampling (padding with zeros) followed by convolution of
- FILT with IMAGE (a.k.a. EXPAND in Burt&Adelson81). The operations
- are combined to avoid unnecessary multiplication of filter samples
- with zeros in the upsampled image. The convolution is done in 9
- sections so that mod operation is not performed unnecessarily.
- Arguments are described in the comment above internal_wrap_reduce.
- WARNING: this subroutine destructively modifes the RESULT image, so
- the user must zero the result before invocation!
- -------------------------------------------------------------------- */
- /* abstract out the inner product computation */
- #define INPROD2(YSTART,YIND,XSTART,XIND) \
- { \
- val = image[im_pos]; \
- for (y_res=YSTART, filt_pos=0, x_filt_stop=x_fdim; \
- x_filt_stop<=filt_size; \
- y_res++, x_filt_stop+=x_fdim) \
- for (x_res=XSTART; \
- filt_pos<x_filt_stop; \
- filt_pos++, x_res++) \
- imval[YIND][XIND] += val * filt[filt_pos]; \
- }
- int internal_wrap_expand(image, filt, x_fdim, y_fdim,
- x_start, x_step, x_stop, y_start, y_step, y_stop,
- result, x_dim, y_dim)
- register image_type *filt, *result;
- register int x_fdim, y_fdim, x_dim, y_dim;
- image_type *image;
- int x_start, x_step, x_stop, y_start, y_step, y_stop;
- {
- register double val;
- register int filt_size = x_fdim*y_fdim;
- image_type **imval;
- register int filt_pos, x_res, y_res, x_filt_stop;
- register int x_pos, y_pos, im_pos;
- int x_ctr_stop = x_dim - x_fdim + 1;
- int y_ctr_stop = y_dim - y_fdim + 1;
- int x_ctr_start = 0;
- int y_ctr_start = 0;
- int x_fmid = x_fdim/2;
- int y_fmid = y_fdim/2;
-
- /* shift start/stop coords to filter upper left hand corner */
- x_start -= x_fmid; y_start -= y_fmid;
- x_stop -= x_fmid; y_stop -= y_fmid;
- if (x_stop < x_ctr_stop) x_ctr_stop = x_stop;
- if (y_stop < y_ctr_stop) y_ctr_stop = y_stop;
- /* Set up pointer array for rows */
- imval = (image_type **) malloc(y_dim*sizeof(image_type *));
- if (imval IS NULL)
- {
- printf("INTERNAL_WRAP: Failed to allocate temp array!");
- return(-1);
- }
- for (y_pos=y_res=0;y_pos<y_dim;y_pos++,y_res+=x_dim)
- imval[y_pos] = (result+y_res);
-
- for (im_pos=0, y_pos=y_start; /* TOP ROWS */
- y_pos<y_ctr_start;
- y_pos+=y_step)
- {
- for (x_pos=x_start;
- x_pos<x_ctr_start;
- x_pos+=x_step, im_pos++)
- INPROD2(y_pos+y_dim, y_res%y_dim, x_pos+x_dim, x_res%x_dim)
-
- for (;
- x_pos<x_ctr_stop;
- x_pos+=x_step, im_pos++)
- INPROD2(y_pos+y_dim, y_res%y_dim, x_pos, x_res)
- for (;
- x_pos<x_stop;
- x_pos+=x_step, im_pos++)
- INPROD2(y_pos+y_dim, y_res%y_dim, x_pos, x_res%x_dim)
- } /* end TOP ROWS */
-
- for (; /* MID ROWS */
- y_pos<y_ctr_stop;
- y_pos+=y_step)
- {
- for (x_pos=x_start;
- x_pos<x_ctr_start;
- x_pos+=x_step, im_pos++)
- INPROD2(y_pos, y_res, x_pos+x_dim, x_res%x_dim)
-
- for (; /* CENTER SECTION */
- x_pos<x_ctr_stop;
- x_pos+=x_step, im_pos++)
- INPROD2(y_pos, y_res, x_pos, x_res)
- for (;
- x_pos<x_stop;
- x_pos+=x_step, im_pos++)
- INPROD2(y_pos, y_res, x_pos, x_res%x_dim)
- } /* end MID ROWS */
-
- for (; /* BOTTOM ROWS */
- y_pos<y_stop;
- y_pos+=y_step)
- {
- for (x_pos=x_start;
- x_pos<x_ctr_start;
- x_pos+=x_step, im_pos++)
- INPROD2(y_pos, y_res%y_dim, x_pos+x_dim, x_res%x_dim)
-
- for (;
- x_pos<x_ctr_stop;
- x_pos+=x_step, im_pos++)
- INPROD2(y_pos, y_res%y_dim, x_pos, x_res)
- for (;
- x_pos<x_stop;
- x_pos+=x_step, im_pos++)
- INPROD2(y_pos, y_res%y_dim, x_pos, x_res%x_dim)
- } /* end BOTTOM ROWS */
- free ((image_type **) imval);
- return(0);
- } /* end of internal_wrap_expand */
- /* Local Variables: */
- /* buffer-read-only: t */
- /* End: */
|