3D convolution

Jayashree_N_ — Sun, 02 Mar 2014 22:17:42 GMT

I am using MKL's VSL multi-dimensional convolution routines to perform a 3D convolution. As a test, I am using the input function as a 3D gaussian as follows:

f(x,y,z) = exp(- (x-2)**2) exp(-(y+2)**2) exp(- (z**2))

and the convolution kernel as

g(x,y,z) = exp(-(x**2 + y**2 + z**2))

The resulting function is formally C(x,y,z), also a gaussian. I have some trouble understanding the access of the exact indices of the ouput 3D array. Given below is an output for the non-zero values of the convoluted function C, in the last column, calculated with the following code (forgive the lengthy information).

26 26 26   26    3    1     0.00000018
27 27 27   27    3    1     0.00000496
28 28 28   28    3    1     0.00008922
29 29 29   29    3    1     0.00102816
30 30 30   30    3    1     0.00759755
31 31 31   31    3    1     0.03599289
32 32 32   32    3    1     0.10934334
33 33 33   33    3    1     0.21295892
34 34 34   34    3    1     0.26596953
35 35 35   35    3    1     0.21295892
36 36 36   36    3    1     0.10934334
37 37 37   37    3    1     0.03599289
38 38 38   38    3    1     0.00759755
39 39 39   39    3    1     0.00102816
40 40 40   40    3    1     0.00008922
41 41 41   41    3    1     0.00000496
42 42 42   42    3    1     0.00000018

      !n is number of points on the 3D grid.
      n=30
      rank=3
      cdim=n+n-1
      r1shape(1)=n;r1shape(2)=n;r1shape(3)=n;r2shape=r1shape;
      convshape(1)=cdim;convshape(2)=cdim;convshape(3)=cdim
      r1stride=1; r2stride=1; convstride=1
      allocate(conv_mkl(cdim,cdim,cdim));conv_mkl=0.d0
      ok=.true.
      precis=VSL_CONV_PRECISION_SINGLE
      stat=vslconvsetmode(task,mode)
      stat=vslconvsetinternalprecision(task,precis)
      stat=vsldconvnewtask(task,mode,rank,r1shape,r2shape,convshape)
      stat=vsldconvexec(task,r1,r1stride,r2,r2stride,conv_mkl,convstride)
      ier=vslconvdeletetask(task)
      !normalise this
      conv_mkl=conv_mkl/sum(conv_mkl)
      if(stat/=VSL_STATUS_OK) then
       print*,'Error: Bad status : ',stat
       ok=.false.
      end if
      if(ier/=0) then
       print*,'Failed to destroy task'
      end if
      write(*,*)'Successful? ',ok
      do i=1,cdim
       do j=1,cdim
        do k=1,cdim
         write(23,'(3I4,3I5,F15.8)')i,i-j+3,i-j-k+4,i,j,k,conv_mkl(i,j,k)
        end do
       end do
      end do

I find that the value of j and k are always 3 and 1 respectively, even when I change parameters in function f. What is so special about 3 and 1?

Thanks.

Please take a look at the

Zhang_Z_Intel — Mon, 03 Mar 2014 20:21:07 GMT

Please take a look at the "Data Allocation" section in the MKL reference manual that describes in great details about multi-dimensional input and output data vectors in convolution routines: http://software.intel.com/en-us/node/470762

topic 3D convolution in Intel® oneAPI Math Kernel Library

3D convolution

Please take a look at the