Take the 2-minute tour ×
Mathematics Stack Exchange is a question and answer site for people studying math at any level and professionals in related fields. It's 100% free, no registration required.

I am not completely sure if this is where a MatLab question belongs, so if not, please direct me where I should ask.

But onto my question. I am working on trying to deconvolution a signal with noise. So I have $h(x)=f(x) \ast g(x) +n(x) $. I want to find what the function g is. n(x) is white gaussian noise. I am using the Wiener Deconvolution Method. The function I am using for h is a box:

x=-20:.01:20;
h=zeros(size(x));
h(1900:2100)=1; 

I know that the Fourier transform of the box should be a Sinc function, but when I use fft(h), I get nothing that resembles this. Additionally, this contains zeros, so when I use the linked method, I end up dividing by zero. So I am wondering how I can go about doing this correctly. (Note: I do not have the symbolic toolbox which was recommended me by a colleague).

share|improve this question
2  
How do you visualize this? Doing plot(x,fftshift(fft(h))) seems to give something that "looks like" sinc, but maybe not enough? –  Cocopuffs Jun 17 '13 at 15:39
    
The use of fftshift makes it look alot nicer. I still have to deal with the zeros... –  yankeefan11 Jun 17 '13 at 15:42
    
And it also contains an imaginary part. –  yankeefan11 Jun 17 '13 at 15:45
    
If a question is mathematical, this is a fine place for MatLab questions as far as I know. This particular question might also fit well (better?) on dsp.stackexchange.com. –  Bjorn Roche Jun 17 '13 at 18:09
    
@JamesMaslek The imaginary part is from a linear phase term that is related to the location of your signal in the time domain and shouldn't really affect what you are doing. Also, make sure you plot the magnitude if the sinc pattern is what you want to see: plot(abs(fftshift(fft(h)))). If you plot the unwrapped angle, you should see a nice line: plot(unwrap(angle(fftshift(fft(h))))) –  AnonSubmitter85 Jun 23 '13 at 8:55

1 Answer 1

up vote 1 down vote accepted

If you want to keep the length of $h$ unchanged (4001 samples) you can define it like

h=[ones(101,1);zeros(3800,1);ones(100,1)]

This will guarantee that its FFT is real-valued (up to numerical noise). The reason for this is that Matlab assumes that the first value of your signal, i.e. $h(1)$, corresponds to the actual index $0$ (i.e. the center of symmetry for a signal with a real-valued FFT). So if you want a signal which is symmetric, you need to shift the part corresponding to negative indices to the right (the 'end' of the signal). This is because the FFT implies a periodic continuation of the signal. Note that you still have to take the real part of the FFT of $h$ because due to small numerical errors, the imaginary part of the FFT of $h$ is in the order of $10^{-14}$, which corresponds to Matlab's floating point accuracy.

share|improve this answer
    
fftshift(fft(fftshift(x))) is often used to make the center pixel in the time/image domain represent zero linear phase in the frequency domain. –  AnonSubmitter85 Jun 23 '13 at 8:50

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.