Thread Subject:

2-D Convolution

Hi all,

I have a question about convolving two matrices together using conv2.
I have a matrix of coded bits say 111 X 111 and I'm convolving it with a 2-D sinc function i.e. sinc(X) .* sinc(Y).
If sampling time T = 1 ( e.g. [X,Y] =meshgrid(-50:50) ) then simply conv2(... , ...) does that. The question is how to perform the same thing with higher sampling rates e.g. T = 0.5 (i.e. the step size in X and Y is 0.5)
So far I am zero padding the coded matrix by adding 2 columns and rows of zeros in between each row and column before convolution. However, based on the final results I'm suspicious that this may not be correct. What do you think?

Elnaz

"Elnaz " <ebsadeghian@gmail.com> wrote in message <jsd72f$ee9$1@newscl01ah.mathworks.com>...
> Hi all,
>
> I have a question about convolving two matrices together using conv2.
> I have a matrix of coded bits say 111 X 111 and I'm convolving it with a 2-D sinc function i.e. sinc(X) .* sinc(Y).
> If sampling time T = 1 ( e.g. [X,Y] =meshgrid(-50:50) ) then simply conv2(... , ...) does that. The question is how to perform the same thing with higher sampling rates e.g. T = 0.5 (i.e. the step size in X and Y is 0.5)
> So far I am zero padding the coded matrix by adding 2 columns and rows of zeros in between each row and column before convolution. However, based on the final results I'm suspicious that this may not be correct. What do you think?
>
> Elnaz

It is quite confusing how convolution is related with sampling time T in your case. Sampling time has nothing to do with convolution. Convolution can be performed between any two arrays of data, whose size may or may not be equal. So, your problem might be somewhere else not in the convolution.

Yumnam Kirani Singh
Tronglaobi Awang Leikai

> It is quite confusing how convolution is related with sampling time T in your case. Sampling time has nothing to do with convolution. Convolution can be performed between any two arrays of data, whose size may or may not be equal. So, your problem might be somewhere else not in the convolution.

Why confusing?
"any array of data" is a sampled waveform (signal). If you think about what we do on paper to convolve two signals in time and then compare it with the spacial convolution in MATLAB you see how these two are different. Since MATLAB does conv on samples regardless of the timing, everything must be translated to sample unit. Here, with higher sampling rates the sinc function will be denser i.e. more samples while zerocrossing times (symbol period) stays the same. In other words, number of samples increase while time stays the same.

"Elnaz " <ebsadeghian@gmail.com> wrote in message <jsd72f$ee9$1@newscl01ah.mathworks.com>...
> Hi all,
>
> I have a question about convolving two matrices together using conv2.
> I have a matrix of coded bits say 111 X 111 and I'm convolving it with a 2-D sinc function i.e. sinc(X) .* sinc(Y).
> If sampling time T = 1 ( e.g. [X,Y] =meshgrid(-50:50) ) then simply conv2(... , ...) does that. The question is how to perform the same thing with higher sampling rates e.g. T = 0.5 (i.e. the step size in X and Y is 0.5)
> So far I am zero padding the coded matrix by adding 2 columns and rows of zeros in between each row and column before convolution. However, based on the final results I'm suspicious that this may not be correct. What do you think?
================

I don't really understand the rationale behind inserting zeros. When you gather samples of bits at a higher rate, why are the samples at the new sampling times all zero? Also, you have to upsample both signals that participate in the convolution. Both the bit matrix and the sinc need to have more samples.

Incidentally, I hope you are doing the convolution separably in X and Y using conv2's 3 argument form. Otherwise, it will be painfully slow. In fact, it would make more sense to perform convolution using ffts, since sincs are more sparse in the Fourier domain. Better still, many people approximate sinc filtering with cubic B-splines filtering, because cubic B-splines are sparser than sincs.



 

> I don't really understand the rationale behind inserting zeros. When you gather samples of bits at a higher rate, why are the samples at the new sampling times all zero? Also, you have to upsample both signals that participate in the convolution. Both the bit matrix and the sinc need to have more samples.
>
The number of bits is constant. There is no signal to upsample there. Upsampling is with the sinc function. I'm convolving a matrix of bits with a 2-D sinc function i.e. each bit is modulated with a sinc and then add those sincs all together:
r(t1, t2) = ?i?jai,j g(t1 – iT , t2 – jT)

ai,j are the bits, and g(t1,t2) is the 2-D sinc (sinc(t1) .* sinc(t2)).

The reason for zeropadding is that when each sinc is multiplied by a bit the next sinc should be at its first zerocrossing (T sec apart from each other), if you know what I mean, in order to correctly implement the convolution. Here T refers to symbol period. If symbol period is the same as sample period then there is no problem there; we simply can write conv2(matrix of bits, g) but if we have e.g. 2 samples at each T (upsampling) then I think we should add zeros in between bits to accommodate for the extra samples of the upsampled sinc.
Did I succeed in explaining the problem?

Thanks,
Elnaz

"Elnaz " <ebsadeghian@gmail.com> wrote in message <jsff10$7lj$1@newscl01ah.mathworks.com>...
>
> The reason for zeropadding is that when each sinc is multiplied by a bit the next sinc should be at its first zerocrossing (T sec apart from each other), if you know what I mean, in order to correctly implement the convolution. Here T refers to symbol period. If symbol period is the same as sample period then there is no problem there; we simply can write conv2(matrix of bits, g) but if we have e.g. 2 samples at each T (upsampling) then I think we should add zeros in between bits to accommodate for the extra samples of the upsampled sinc.
> Did I succeed in explaining the problem?
================

Yes. So basically you're trying to sinc-interpolate between the bits... You can implement this as a convolution by interspersing zeros, but then it's important that you sample the sinc more finely so that the zero crossings are still at the locations of the bits. If that's what you've been doing, I don't know why you wouldn't be getting the results you want.

My issue with implementing as a convolution is that it can be a bit computationally inefficient. If you have N symbols per period, the convolution is doing N times more computation than as really needed because the inserted zeros contribute 0 to the convolution. The following FEX file will let you do this more efficiently using matrix methods

http://www.mathworks.com/matlabcentral/fileexchange/26292-regular-control-point-interpolation-matrix-with-boundary-conditions


Below is an example. The FEX download also contains examples of how to do something similar with cubic B-splines, which as I said is often used as a cheaper approximation to sinc interp.



%%Fake data
SamplesperPeriod=2;
matrixdim=100;
BitMatrix=randi(2,matrixdim)-1;


%%Engine
sinc1D=sinc(pi/SamplesperPeriod*(-matrixdim:matrixdim));

T=interpMatrix(sinc1D, 'max', matrixdim, SamplesperPeriod);
 
result=T*BitMatrix*T.';

"Matt J" wrote in message <jsfgnb$f9t$1@newscl01ah.mathworks.com>...
>
> T=interpMatrix(sinc1D, 'max', matrixdim, SamplesperPeriod);
============

In this case, because your kernels are not short, it will be more efficient to express T as a full matrix, instead of sparse:

 T=full(interpMatrix(sinc1D, 'max', matrixdim, SamplesperPeriod));
  

"Elnaz " <ebsadeghian@gmail.com> wrote in message <jsf91n$a9d$1@newscl01ah.mathworks.com>...
> > It is quite confusing how convolution is related with sampling time T in your case. Sampling time has nothing to do with convolution. Convolution can be performed between any two arrays of data, whose size may or may not be equal. So, your problem might be somewhere else not in the convolution.
>
> Why confusing?
> "any array of data" is a sampled waveform (signal). If you think about what we do on paper to convolve two signals in time and then compare it with the spacial convolution in MATLAB you see how these two are different. Since MATLAB does conv on samples regardless of the timing, everything must be translated to sample unit. Here, with higher sampling rates the sinc function will be denser i.e. more samples while zerocrossing times (symbol period) stays the same. In other words, number of samples increase while time stays the same.

Sampling is altogether a different process from convolution. Sampling has something to do with the resolution. Yes, the number of samples increases or decreases with the decrease or increase in sampling period. That could be a part of data acquisition or digitzation process. It is not always mendatory that the sampling should be followed by convolution. We may or may not apply convolution after sampling the signals at different various sampling rate.
I do not know how you do convolution using pen and paper that gives different results with the matlab conv function. I get the almost the same result.
For example:
>>x=[2 3 5 1 7 2 5];
>>y=[1 2 1];
The convolution of x and y gives [2 7 13 14 14 17 16 12 5]
If you get different result, I think you need to recheck your way of doing convolution using pen and paper.

Yumnam Kirani Singh
Tronglaobi Awang Leikai

 

> Yes. So basically you're trying to sinc-interpolate between the bits... You can implement this as a convolution by interspersing zeros, but then it's important that you sample the sinc more finely so that the zero crossings are still at the locations of the bits. If that's what you've been doing, I don't know why you wouldn't be getting the results you want.

Exactly. However, I found a clue why I'm not getting the result I want:
So, I try two convolutions one with without oversampling i.e. Tsymbol=1 and Tsample=1, and the second with upsampling e.g. Tsymbol=1 and Tsample=0.5. Therefore, I expect that the second convolution to be the interpolated version of the first one. In other words, I expect to have the exact same elements of the 1st convolution in the 2nd convolution but with 1 row and column more data in between each.
The thing is that, those supposed-to-be equal elements are not EXACTLY equal. They differ in very small fractional digits and this happens with larger number of bits in the coded matrix.
Can it be a computational error/approximation in MATLAB?

I include the example code here:
tbyT = (-50:50).';
[t2byT,t1byT] = meshgrid(tbyT);
a = (-1).^(rand(101,101) > 0.5);
 s = conv2(a, sinc(t1byT-0.1).*sinc(t2byT-0.2) );

tbyT1 = (-50:0.5:50).'; %over-sampling by 2
[t2byT1,t1byT1] = meshgrid(tbyT1);
 bits(1:2:201,1:2:201) = a; %zero-pad coded matrix "a"
s1 = conv2(bits, sinc(t1byT1 - 0.1).*sinc(t2byT1 - 0.2) );

Here we expect that: s(51,51) ==s1(101,101) but this returns ans=0 and we see that the two differ in the last few fractional digits.

"Elnaz " <ebsadeghian@gmail.com> wrote in message <jssecr$9s6$1@newscl01ah.mathworks.com>...
>
> The thing is that, those supposed-to-be equal elements are not EXACTLY equal. They differ in very small fractional digits and this happens with larger number of bits in the coded matrix.
> Can it be a computational error/approximation in MATLAB?
>
> I include the example code here:
> tbyT = (-50:50).';
> [t2byT,t1byT] = meshgrid(tbyT);
> a = (-1).^(rand(101,101) > 0.5);
> s = conv2(a, sinc(t1byT-0.1).*sinc(t2byT-0.2) );
>
> tbyT1 = (-50:0.5:50).'; %over-sampling by 2
> [t2byT1,t1byT1] = meshgrid(tbyT1);
> bits(1:2:201,1:2:201) = a; %zero-pad coded matrix "a"
> s1 = conv2(bits, sinc(t1byT1 - 0.1).*sinc(t2byT1 - 0.2) );
>
> Here we expect that: s(51,51) ==s1(101,101) but this returns ans=0 and we see that the two differ in the last few fractional digits.
============

I see many errors or potential errors in the code:

(1) You haven't shown the implementation of your sinc() function. Is it time-scaled so that the zero crossings of the sinc occur at integer locations?

(2) Your sincs are shifted off center wrt the sampling locations by 0.1 and 0.2. This should make it so that the zero-crossings don't fall at the sample locations.

(3) You are not using CONV2 with the 'same' option. That means the sample locations of "s1" won't coincide with the sample locations of "bits".

(4) Expecting EXACT equality between s(51,51) and s1(101,101) may not be so realistic. There are always floating point approximations in the evaluation of sincs because a true sinc is only exactly zero at multiples of pi, which cannot be exactly evaluated in finite precision. Possibly, you could base your sinc implementation on SIND instead of SIN and that could force them to have exact zero crossings. But have you done that?

> I see many errors or potential errors in the code:
>
> (1) You haven't shown the implementation of your sinc() function. Is it time-scaled so that the zero crossings of the sinc occur at integer locations?

This is MATLAB's sinc. I don't have another implementation. So, to answer your question, yes, zero crossings occur at integer places i.e. multiples of 1.
Also, zero-crossings as you said are not exactly zero but 1e-17.
 
> (2) Your sincs are shifted off center wrt the sampling locations by 0.1 and 0.2. This should make it so that the zero-crossings don't fall at the sample locations.

Yes, this is due to the physical aspect of the problem I am modeling in MATLAB. And, this is fine. The only effect is that here instead of +1/-1, we'll have e.g. .95 or 1.12 but, again, whatever the number is it should appear exactly the same in both s and s1.

> (3) You are not using CONV2 with the 'same' option. That means the sample locations of "s1" won't coincide with the sample locations of "bits".

I don't understand what you mean here. Could you explain a bit more?
 
> (4) Expecting EXACT equality between s(51,51) and s1(101,101) may not be so realistic. There are always floating point approximations in the evaluation of sincs because a true sinc is only exactly zero at multiples of pi, which cannot be exactly evaluated in finite precision. Possibly, you could base your sinc implementation on SIND instead of SIN and that could force them to have exact zero crossings. But have you done that?

No, I did not implement sinc as sin(t)/t; I'm using MATLAB's sinc.
One thing that may shed some light on this problem is that when I said earlier that this problem happens with larger number of bits in the coded matrix, it exaclty happens with number of bits bigger than 100 (half the length of the sinc). In the above code, we have 101 X 101 bits in the coded matrix. But I can't figure out why?

"Elnaz " <ebsadeghian@gmail.com> wrote in message <jssinl$r5s$1@newscl01ah.mathworks.com>...
>
> > (2) Your sincs are shifted off center wrt the sampling locations by 0.1 and 0.2. This should make it so that the zero-crossings don't fall at the sample locations.
>
> Yes, this is due to the physical aspect of the problem I am modeling in MATLAB. And, this is fine. The only effect is that here instead of +1/-1, we'll have e.g. .95 or 1.12 but, again, whatever the number is it should appear exactly the same in both s and s1.
===========

OK, I've run your code. I see a very trivial difference between s and s1, a difference easily attributable to machine precision limits:

>> percentError=abs(s(51,51) - s1(101,101))/s(51,51)*100

percentError =

 -8.0205e-014


Not sure why you expect to do better than this. I'm using my own implementation of sinc instead of the Signal Processing Toolboxes (I don't have that toolbox), but I doubt it matters. I too get a zero crossing error on the order of 1e-17.





> > (3) You are not using CONV2 with the 'same' option. That means the sample locations of "s1" won't coincide with the sample locations of "bits".
>
> I don't understand what you mean here. Could you explain a bit more?
===========

In hindsight, I don't think it matter. However, the way you're running CONV2 right now, you are extrapolating as well as interpolating. If you only want the convolution results in between the bits sample locations, you would use the syntax CONV2(A,B,'same')