| Re: Convolution Kernel [message #73899 is a reply to message #73813] |
Fri, 03 December 2010 10:06   |
Jeremy Bailin
Messages: 618
Registered: April 2008
|
Senior Member |
|
|
On Dec 3, 5:23 am, MC <morefl...@gmail.com> wrote:
> On Dec 3, 6:35 am, Gray <grayliketheco...@gmail.com> wrote:
>
>
>
>
>
>
>
>
>
>> Hi all,
>
>> Maybe my calculus is screwy, but this doesn't make sense to me.
>> Here's my issue:
>
>> I have two astronomical images (of stars). I've fit an average PSF as
>> a Moffat profile for each of the two images. I want to find the
>> optimal convolution kernel to match the two psfs, so I call on my old
>> friend Mr. Fourier. If MA is the Moffat profile for image A and MB is
>> the Moffat profile for image B (both 2d), and K is my optimal kernel,
>> then I can do this:
>
>> MA ** K = MB --> ** is convolution in this scenario
>> F(MA**K) = F(MB) --> F() is the Fourier transform
>> F(MA) * F(K) = F(MB)
>> K = F^-1(F(MB)/F(MA))
>
>> With me so far? So I do this in IDL.
>> IDL> ma = moffat(params_a)
>> IDL> mb = moffat(params_b)
>> IDL> fma = fft(ma) & fmb = fft(mb)
>> IDL> k = fft(fma/fmb,/inverse)
>> IDL> mc = convol(ma,k)
>
>> What I get, however, is that MC is a 2d delta function. Why...? It
>> happens with 2d Gaussians, as well. Thanks for your help!
>
> Hmm, not sure what Moffat does, but maybe the division by fmb terms
> with small amplitudes could be the problem? What happens if you add a
> small constant to fmb? If you want to match the OTF of the two images
> could you just convolve a with the psf of b and convolve b with the
> psf of a?
>
> Hope this helps
> MC
(I thought I posted this this morning, but it never showed up)
This is certainly what's going on. I tried this with some Gaussians,
but used
fma/(fmb > 1e-10)
instead inside the k=fft(...) line, and it changed from nonsense to
working beautifully. Outside of the support of the Moffat fit to b,
you're just getting amplified numerical noise. Try looking at the
difference between those ratios (with TVIMAGE or SURFACE or something)
and you'll see what's going on.
-Jeremy.
|
|
|
|
comp.lang.idl-pvwave archive
Members
Search
Help
Login
Home
