Xc@`sddlmZmZmZddlZddlmZddlZddl m Z ddl m Z mZddlmZddlmZddlmZd ejfd YZeZd Zd ZdS( i(tabsolute_importtprint_functiontdivisionN(txrange(tbasic(tblas_header_texttblas_header_version(tldflags(tstrutil(tgrad_undefinedtConv3DcB`sweZdZd ZdZdZdZdZdZdZ dZ dZ d Z d Z d ZRS( s 3D `convolution` of multiple filters on a minibatch. Notes ----- Does not flip the kernel, moves kernel with a user specified stride. cC`s dtfS(Ni(R(tself((s9/tmp/pip-build-X4mzal/theano/theano/tensor/nnet/Conv3D.pytc_code_cache_version9sc C`stj|}tj|}tj|}tj|}|jdttt|jdf} tj|d||||gdtj|j| g} | S(sK Parameters ---------- V Visible unit, input(batch,row,column,time,in channel) W Weights, filter(out channel,row,column,time,in channel) b bias, shape == (W.shape[0],) d strides when moving the filter over the input(dx,dy,dt) itinputstoutputs(tTtas_tensor_variablet broadcastabletFalsettheanotApplyt TensorTypetdtype( R tVtWtbtdtV_tW_tb_td_tbcasttnode((s9/tmp/pip-build-X4mzal/theano/theano/tensor/nnet/Conv3D.pyt make_node<s#c C`s;|\}}}}|\}tjjj|tj|ddddddf|||jdd!}tj||j}|j} tjjj ||| |} tj| |j} tj |dd} tj| |j} t |d|dd} dt |kr)|j dk r)|j } nd } dt |kr\|j dk r\|j }nd }dt |kr|j dk r|j }nd }dt |kr|j dk r|j }nd }d | d|d|_ d| d|d|d| _ d| d|d|d|d| _ || | | gS(Niiitaxisiis~The gradient of Conv3D with respect to the convolution stride is undefined because Conv3D is only defined for integer strides.tnamet anon_dCdHtanon_Vtanon_Wtanon_bsConv3D_dCdV(dCdH=s,V=t)sConv3D_dCdW(dCdH=s,W=sConv3D_dCdb(dCdH=s,b=(iiii(Rttensortnnett convTransp3DRt zeros_liketshapetpatternbroadcastRt convGrad3DtsumR tdirR#tNone(R R toutput_gradientsRRRRtdCdHtdCdVtWShapetdCdWtdCdbtdCddt dCdH_nametV_nametW_nametb_name((s9/tmp/pip-build-X4mzal/theano/theano/tensor/nnet/Conv3D.pytgradWs: A   ! ! ! ! !)cC`s3|\}}}}t|||||dd //printf(" Conv3D c code\n"); //Check dimensionality of inputs if (PyArray_NDIM(%(W)s) != 5) { PyErr_Format(PyExc_ValueError, "Conv3D: W must be a 5 dimensional tensor"); %(fail)s } if (PyArray_NDIM(%(V)s) != 5) { PyErr_Format(PyExc_ValueError, "Conv3D: V must be a 5 dimensional tensor"); %(fail)s } if (PyArray_NDIM(%(b)s) != 1) { PyErr_Format(PyExc_ValueError,"Conv3D: b must be a vector."); %(fail)s } if (PyArray_NDIM(%(d)s) != 1) { PyErr_Format(PyExc_ValueError,"Conv3D: d must be a vector."); %(fail)s } if (PyArray_DIMS(%(d)s)[0] != 3) { PyErr_Format(PyExc_ValueError,"Conv3D: 3 stride length arguments expected (row, col, time) but %%li were given", (long)PyArray_DIMS(%(d)s)[0]); %(fail)s } //Read and check sizes of inputs { // exta scope so error handler jumps don't cause errors const int batchSize = PyArray_DIMS(%(V)s)[0]; const int outputChannels = PyArray_DIMS(%(W)s)[0]; const int inputChannels = PyArray_DIMS(%(V)s)[4]; if (PyArray_DIMS(%(W)s)[4] != inputChannels) { PyErr_Format(PyExc_ValueError, "Conv3D: W operates on a %%ld channel image but the image has %%d channels. Overall shape of input: (%%ld,%%ld,%%ld,%%ld,%%ld)", (long)PyArray_DIMS(%(W)s)[4], inputChannels, (long)PyArray_DIMS(%(V)s)[0], (long)PyArray_DIMS(%(V)s)[1], (long)PyArray_DIMS(%(V)s)[2], (long)PyArray_DIMS(%(V)s)[3], (long)PyArray_DIMS(%(V)s)[4]); %(fail)s } if (PyArray_DIMS(%(b)s)[0] != outputChannels) { PyErr_Format(PyExc_ValueError, "Conv3D: b adds to a(n) %%ld channel output image but the output has %%d channels", (long)PyArray_DIMS(%(b)s)[0], outputChannels); %(fail)s } { //extra scope so error handler jumps don't cause errors const int filterHeight = PyArray_DIMS(%(W)s)[1]; const int filterWidth = PyArray_DIMS(%(W)s)[2]; const int filterDur = PyArray_DIMS(%(W)s)[3]; const int vidHeight = PyArray_DIMS(%(V)s)[1]; const int vidWidth = PyArray_DIMS(%(V)s)[2]; const int vidDur = PyArray_DIMS(%(V)s)[3]; if (vidHeight < filterHeight) { PyErr_Format(PyExc_ValueError, "W has a height of %%i but V is only %%i pixels tall",filterHeight,vidHeight); %(fail)s } { // extra scope so fail works if (vidWidth < filterWidth) { PyErr_Format(PyExc_ValueError, "W has a width of %%i but V is only %%i pixels wide",filterWidth,vidWidth); %(fail)s } { // extra scope so fail works if (vidDur < filterDur) { PyErr_Format(PyExc_ValueError, "W has a duration of %%i but V is only %%i pixels long",filterDur,vidDur); %(fail)s } { // extra scope so fail works //Read and check stride arguments const int dr = *(dtype_%(d)s*) PyArray_GETPTR1(%(d)s,0); const int dc = *(dtype_%(d)s*) PyArray_GETPTR1(%(d)s,1); const int dt = *(dtype_%(d)s*) PyArray_GETPTR1(%(d)s,2); if (dr <= 0 || dc <= 0 || dt <= 0) { PyErr_Format(PyExc_ValueError,"Conv3D: Strides must all be positive but are %%i, %%i, %%i",dr,dc,dt); %(fail)s } { // extra scope so fail works //Make correctly sized output const long long outputHeight = int( (vidHeight - filterHeight) / dr )+1; const long long outputWidth = int( (vidWidth - filterWidth) / dc )+1; const long long outputDur = int( (vidDur - filterDur) / dt ) +1; npy_intp dims[5]; dims[0] = batchSize; dims[4] = outputChannels; dims[1] = outputHeight; dims[2] = outputWidth; dims[3] = outputDur; if(!(%(H)s) || PyArray_DIMS(%(H)s)[0]!=dims[0] || PyArray_DIMS(%(H)s)[1]!=dims[1] || PyArray_DIMS(%(H)s)[2]!=dims[2] || PyArray_DIMS(%(H)s)[3]!=dims[3] || PyArray_DIMS(%(H)s)[4]!=dims[4]){ Py_XDECREF(%(H)s); %(H)s = (PyArrayObject *) PyArray_SimpleNew(5, dims, PyArray_DESCR(%(V)s)->type_num); if (!(%(H)s)) { PyErr_Format(PyExc_MemoryError,"Conv3D: Could not allocate output."); %(fail)s } } { // extra scope so fail works #define ELEM_AT(x, i) * ( dtype_ ## x *) ( PyArray_BYTES(x) + (i) ) const int ws0 = PyArray_STRIDES(%(W)s)[0]; const int ws1 = PyArray_STRIDES(%(W)s)[1]; const int ws2 = PyArray_STRIDES(%(W)s)[2]; const int vs1 = PyArray_STRIDES(%(V)s)[1]; const int ws4 = PyArray_STRIDES(%(W)s)[4]; const int vs4 = PyArray_STRIDES(%(V)s)[4]; const int ws3 = PyArray_STRIDES(%(W)s)[3]; const int vs3 = PyArray_STRIDES(%(V)s)[3]; const int vs2 = PyArray_STRIDES(%(V)s)[2]; const int bs = PyArray_STRIDES(%(b)s)[0]; const int hs4 = PyArray_STRIDES(%(H)s)[4]; // Compute H //H[i,j,x,y,t] = b_j + sum_k sum_l sum_m sum_z W[j,z,k,l,m] V[i,z, dr*r+k,dc*c+l,dt*t+m] //TODO: add special cases // ex: filterDur == 1 && batchSize == 1 && dt = 1 (for SFA) // ex: inputChannels == 1 tfloat64tdgemv_tfloat32tsgemv_s#Unrecognized dtype for convolution s if (inputChannels > 20 && outputChannels > 20 && ws4 == sizeof(ELEM_AT(%(W)s,0))) { //std::cout << "lots of channels special case code" << std::endl; #define blas_type dtype_ ## %(V)s const blas_type constant_one = 1.0; char N = 'T'; int ws0e = ws0 / sizeof(ELEM_AT(%(W)s,0)); int vs4e = vs4 / sizeof(ELEM_AT(%(V)s,4)); int hs4e = hs4 / sizeof(ELEM_AT(%(H)s,4)); //special case code for the "lots of channels" case //uses a BLAS matrix vector multiply to compute the contribute for //all channels of an input pixel to all channels of an output pixel //simultaneously long long Hpos = 0; long long Vpos = 0; for (int i = 0; i < batchSize; i++) { long long Hposi = Hpos; long long Vposi = Vpos; for (int r = 0; r < outputHeight; r++) { long long Hposr = Hpos; long long Vposr = Vpos; for (int c = 0; c < outputWidth; c++) { long long Hposc = Hpos; long long Vposc = Vpos; for (int t = 0; t < outputDur; t++) { long long Hpost = Hpos; long long Vpost = Vpos; //of the loops so far, j should be the innermost, because //each loop through j visits the same elements of V //this implies that the last index of H should be the j index //since V and H should have the same format, this means //z should be the last index in v, and therefore the innermost //of the next set of for loops int Wpos = 0; int bPos = 0; long long Hposj = Hpos; for (int j = 0; j < outputChannels; j++) { // H[i,r,c,t,j] = b[j] ELEM_AT(%(H)s,Hposj) = ELEM_AT(%(b)s,bPos); Hposj += hs4; bPos += bs; } dtype_%(H)s * writePos = & ELEM_AT(%(H)s,Hpos); for (int k =0; k < filterHeight; k++) { int Wposk = Wpos; long long Vposk = Vpos; for (int l = 0; l < filterWidth; l++) { int Wposl = Wpos; long long Vposl = Vpos; for (int m = 0; m < filterDur; m++) { //H[i,r,c,t,:] += N.dot(W[:,k,l,m,:],V[i,dr*r+k,dc*c+l,dt*t+m,:]) //note: changing the weights so that outputChannels and inputChannels were the last two rather than //the first and last elements did not speed this up, even for extremely large input sizes %(gemv)s(&N, & inputChannels, & outputChannels, &constant_one, & ELEM_AT( %(W)s , Wpos),& ws0e, & ELEM_AT(%(V)s, Vpos),& vs4e, &constant_one, writePos,& hs4e); Wpos += ws3; Vpos += vs3; } // close m Wpos = Wposl + ws2; Vpos = Vposl + vs2; } //close l Wpos = Wposk + PyArray_STRIDES(%(W)s)[1]; Vpos = Vposk + PyArray_STRIDES(%(V)s)[1]; } //close k Hpos = Hpost + PyArray_STRIDES(%(H)s)[3]; Vpos = Vpost + vs3 * dt; } //close t Hpos = Hposc + PyArray_STRIDES(%(H)s)[2]; Vpos = Vposc + vs2 * dc; } //close c Hpos = Hposr + PyArray_STRIDES(%(H)s)[1]; Vpos = Vposr + PyArray_STRIDES(%(V)s)[1] * dr; } //closes r Hpos = Hposi + PyArray_STRIDES(%(H)s)[0]; Vpos = Vposi + PyArray_STRIDES(%(V)s)[0]; } //closes i } //closes "lots of channels" special case code else s { //General case code //std::cout << "general case code" << std::endl; long long Hpos = 0; long long Vpos = 0; for (int i = 0; i < batchSize; i++) { long long Hposi = Hpos; long long Vposi = Vpos; for (int r = 0; r < outputHeight; r++) { long long Hposr = Hpos; long long Vposr = Vpos; for (int c = 0; c < outputWidth; c++) { long long Hposc = Hpos; long long Vposc = Vpos; for (int t = 0; t < outputDur; t++) { long long Hpost = Hpos; long long Vpost = Vpos; //of the loops so far, j should be the innermost, because //each loop through j visits the same elements of V //this implies that the last index of H should be the j index //since V and H should have the same format, this means //z should be the last index in v, and therefore the innermost //of the next set of for loops int Wpos = 0; int bPos = 0; for (int j = 0; j < outputChannels; j++) { long long Hposj = Hpos; long long Vposj = Vpos; int Wposj = Wpos; // H[i,r,c,t,j] = b[j] dtype_%(H)s & writePos = ELEM_AT(%(H)s,Hpos); writePos = ELEM_AT(%(b)s,bPos); for (int k =0; k < filterHeight; k++) { int Wposk = Wpos; long long Vposk = Vpos; for (int l = 0; l < filterWidth; l++) { int Wposl = Wpos; long long Vposl = Vpos; for (int m = 0; m < filterDur; m++) { int Wposm = Wpos; long long Vposm = Vpos; for (int z = 0; z < inputChannels; z++) { //H[i,r,c,t,j] += W[j,z,k,l,m] * V[i,dr*r+k, dc*c+l, dt*t+m,z] writePos += ELEM_AT(%(W)s,Wpos) * ELEM_AT(%(V)s,Vpos); Wpos += ws4; Vpos += vs4; } // close z Wpos = Wposm + ws3; Vpos = Vposm + vs3; } // close m Wpos = Wposl + ws2; Vpos = Vposl + vs2; } //close l Wpos = Wposk + PyArray_STRIDES(%(W)s)[1]; Vpos = Vposk + PyArray_STRIDES(%(V)s)[1]; } //close k bPos += bs; Wpos = Wposj + ws0; Hpos = Hposj + hs4; Vpos = Vposj; //std::cout << "incremented Wpos by " << ws0 << std::endl; //std::cout << "incremented Hpos by " << hs4 << std::endl; } //close j Hpos = Hpost + PyArray_STRIDES(%(H)s)[3]; Vpos = Vpost + vs3 * dt; } //close t Hpos = Hposc + PyArray_STRIDES(%(H)s)[2]; Vpos = Vposc + vs2 * dc; } //close c Hpos = Hposr + PyArray_STRIDES(%(H)s)[1]; Vpos = Vposr + PyArray_STRIDES(%(V)s)[1] * dr; } //closes r Hpos = Hposi + PyArray_STRIDES(%(H)s)[0]; Vpos = Vposi + PyArray_STRIDES(%(V)s)[0]; } //closes i } //closes general case code }}}}}}} //extra scope so error handler jumps don't cross declarations ///////////// < /code generated by Conv3D > ( R RRtconfigtblasRRt ExceptiontvalueRt render_stringtlocals(R R tnodenameR RtsubRRRRRatHt codeSourcetVVtWVtbvtdvtHVtgemv((s9/tmp/pip-build-X4mzal/theano/theano/tensor/nnet/Conv3D.pytc_codes$   $  a a((t__name__t __module__t__doc__t __props__R R!R>RARVRWRXR\R^R`Rv(((s9/tmp/pip-build-X4mzal/theano/theano/tensor/nnet/Conv3D.pyR .s   5       cC`st||||S(s 3D "convolution" of multiple filters on a minibatch. (does not flip the kernel, moves kernel with a user specified stride) Parameters ---------- V Visible unit, input. Dimensions: (batch, row, column, time, in channel). W Weights, filter. Dimensions: (out channel, row, column, time ,in channel). b Bias, shape == (W.shape[0],). d Strides when moving the filter over the input(dx, dy, dt). Notes ----- The order of dimensions does not correspond to the one in `conv2d`. This is for optimization. The GPU implementation is very slow. You should use :func:`conv3d2d ` or :func:`conv3d_fft ` for a GPU graph instead. See Also -------- Someone made a script that shows how to swap the axes between both 3d convolution implementations in Theano. See the last `attachment `_ (t_conv3D(RRRR((s9/tmp/pip-build-X4mzal/theano/theano/tensor/nnet/Conv3D.pytconv3D(s$cC`st|jdkstt|jdks6tt|jdkrgt|jtsgtnt|dkst|jd}|jd}|jd}|jd|krtdt|jddt|n|jd}|jd}|jd} |jd} |jd} |jd} | |ksGt| |ksYt| | kskt|\} }}| dkst|dkst|dkstt| || d}t| ||d}t| | |d}tj |||||fd |j }xt d|jdD]}xt d|jdD]g}x^t d|jdD]F}x=t d|jdD]%}xt d|jdD]}||||||||fs  "  '