sepConv2PreserveEnergy - energy preserving 2d convolution with separable filter result = sepConv2PreserveEnergy(yFilter,xFilter,data) Convolves data with the separable filters xFilter and yFilter. For border elements, the filters are truncated and the filter energy is renormalized to avoid bleeding over the edge. The result has the same size as data. See also mexConv2PreserveEnergy, maxNormalizeIterative.
0001 % sepConv2PreserveEnergy - energy preserving 2d convolution with separable filter 0002 % 0003 % result = sepConv2PreserveEnergy(yFilter,xFilter,data) 0004 % 0005 % Convolves data with the separable filters xFilter and 0006 % yFilter. For border elements, the filters are truncated 0007 % and the filter energy is renormalized to avoid bleeding 0008 % over the edge. The result has the same size as data. 0009 % 0010 % See also mexConv2PreserveEnergy, maxNormalizeIterative. 0011 0012 % This file is part of the SaliencyToolbox - Copyright (C) 2006-2008 0013 % by Dirk B. Walther and the California Institute of Technology. 0014 % See the enclosed LICENSE.TXT document for the license agreement. 0015 % More information about this project is available at: 0016 % http://www.saliencytoolbox.net 0017 0018 function result = sepConv2PreserveEnergy(filter1,filter2,data) 0019 0020 [sd1,sd2] = size(data); 0021 0022 % 2d convolution with zero padding from Matlab 0023 result = conv2(filter1,filter2,data,'same'); 0024 0025 % filter length and half-length for the vertical filter 0026 fl1 = length(filter1); 0027 fl1b = floor((fl1-1)/2); 0028 0029 % cumulative sums forward and backward 0030 fc1a = cumsum(filter1); 0031 fc1b = cumsum(filter1(end:-1:1)); 0032 fs1 = sum(filter1); 0033 0034 % Is data field very small? 0035 if (fl1 > sd1) 0036 % very small -> need to find the various sums of filter 0037 % elements used for each convolution 0038 tmp = [zeros(1,fl1) fc1a(:)' fs1*ones(1,fl1)]; 0039 range = fl1+fl1b+(1:sd1); 0040 ff1 = repmat((tmp(range)-tmp(range-sd1))',1,sd2); 0041 0042 % normalize result by (actual filter size / full filter size) 0043 result = result * fs1 ./ ff1; 0044 else 0045 % border with incomplete filter coverage at top 0046 ff1a = repmat(fc1a(fl1b+1:end-1)',1,sd2); 0047 result(1:fl1b,:) = result(1:fl1b,:) * fs1 ./ ff1a; 0048 0049 % border with incomplete filter coverage at bottom 0050 ff1b = repmat(fc1b(fl1b+1:end-1)',1,sd2); 0051 result(end:-1:end-fl1b+1,:) = result(end:-1:end-fl1b+1,:) * fs1 ./ ff1b; 0052 end 0053 0054 % now the same again for the horizontal filter 0055 0056 % filter length and half-length 0057 fl2 = length(filter2); 0058 fl2b = floor((fl2-1)/2); 0059 0060 % cumulative sums forward and backward 0061 fc2a = cumsum(filter2); 0062 fc2b = cumsum(filter2(end:-1:1)); 0063 fs2 = sum(filter2); 0064 0065 % Is data field very small? 0066 if (fl2 > sd2) 0067 % very small -> need to find the various sums of filter 0068 % elements used for each convolution 0069 tmp = [zeros(1,fl2) fc2a(:)' fs2*ones(1,fl2)]; 0070 range = fl2+fl2b+(1:sd2); 0071 ff2 = repmat((tmp(range)-tmp(range-sd2)),sd1,1); 0072 0073 % normalize result by (actual filter size / full filter size) 0074 result = result * fs2 ./ ff2; 0075 else 0076 % border with incomplete filter coverage at the left 0077 ff2a = repmat(fc2a(fl2b+1:end-1),sd1,1); 0078 result(:,1:fl2b) = result(:,1:fl2b) * fs2 ./ ff2a; 0079 0080 % border with incomplete filter coverage at the right 0081 ff2b = repmat(fc2b(fl2b+1:end-1),sd1,1); 0082 result(:,end:-1:end-fl2b+1) = result(:,end:-1:end-fl2b+1) * fs2 ./ ff2b; 0083 end