将一维高斯插值到二维高斯

发布于 2024-08-25 11:09:03 字数 561 浏览 2 评论 0原文

假设我有一个一维高斯函数。 它的长度是 600。

我想将其插值成大小为 600 X 600 的 2D 高斯。

这是我编写的代码(OTFx 是高斯函数,OTF - 2d 插值函数):

[x, y] = meshgrid([-300:299], [-300:299]);
r = sqrt((x .^ 2) + (y .^ 2));

OTF = interp1([-300:299], OTFx, r(:), 'spline');
OTF = reshape(OTF, [600, 600]);

问题是我最后得到了 Overshoot : 替代文本

如何防止这种超调? 对于单调递减函数是否有更好的插值算法?

注意:我正在寻找一种将一维函数插值到二维径向对称函数的通用解决方案,高斯只是一个例子。

Let's say I have a 1D Gaussian function.
Its length is 600.

I want to interpolate it into a 2D Gaussian of the size 600 X 600.

This is the code I wrote (OTFx is the Gaussian Function, OTF - 2d Interpolated Function):

[x, y] = meshgrid([-300:299], [-300:299]);
r = sqrt((x .^ 2) + (y .^ 2));

OTF = interp1([-300:299], OTFx, r(:), 'spline');
OTF = reshape(OTF, [600, 600]);

The problem is I get Overshoot at the end:
alt text

How can I prevent this overshoot?
Is there better interpolating algorithm for monotonic descending functions?

NOTE: I am looking for a generic solution for interpolating a 1D function into a 2D radially symmetric function, the Gaussian is just an example.

如果你对这篇内容有疑问,欢迎到本站社区发帖提问 参与讨论,获取更多帮助,或者扫码二维码加入 Web 技术交流群。

扫码二维码加入Web技术交流群

发布评论

需要 登录 才能够评论, 你可以免费 注册 一个本站的账号。

评论(3

落在眉间の轻吻 2024-09-01 11:09:03

编辑:根据您的澄清,很清楚发生了什么。您正在尝试对超出可用数据范围的函数进行插值 - 即您将从插值转向外推。样条线将导致您观察到的超调。解决方案只是确保您的一维函数具有区间 [min(r), max(r)] 内的值。请注意,在原始数据中, max(r) 约为 424,而您要插值的函数定义在范围 [-300,299]

% Simulated overshoot, see left figure:
x1d = [-300:299];
[x,y]=meshgrid(x1d,x1d);
r = sqrt(x.^2+y.^2);
gsn1d = exp(-x1d.^2/500);
lowpass = @(x)(x1d > -x & x1d < x);
gsn1dcutoff = ifft(fftshift(lowpass(10).*fftshift(fft(gsn1d))));
plot(gsn1dcutoff)
OTF2d = reshape(interp1(x1d,gsn1dcutoff,r(:),'spline'),[length(x1d),length(x1d)]);
mesh(OTF2d)

% Quick and dirty fix, see right figure:
x1dExtended = linspace(min(x1d*sqrt(2)),max(x1d*sqrt(2)),ceil(length(x1d)*sqrt(2)));
gsn1dE = exp(-x1dExtended.^2/500);
% ^^^ note that this has 600*sqrt(2) points and is defined on the diagonal of your square.   Now we can low-pass filter in the freq. domain to add ripple in space domain:
lowpass = @(x)(x1dExtended > -x & x1dExtended < x);
gsn1dcutoff = -real(ifft(fftshift(lowpass(10).*fftshift(fft(gsn1dE)))));
plot(gsn1dcutoff)
OTF2d = reshape(interp1(x1dExtended,gsn1dcutoff,r(:),'spline'),[length(x1d),length(x1d)]);
mesh(OTF2d)

替代文本 http://img54.imageshack.us/img54/8255/clipboard01vz.png

EDIT: based on your clarification, it's clear what's going on. You are trying to interpolate a function beyond the range of available data -- i.e. you are going from interpolation to extrapolation. Splines are going to result in the overshoot that you are observing. The solution is simply to make sure that your 1D function has values in the interval [min(r), max(r)]. Note that in the original data, max(r) is about 424, while the function you are interpolating is defined on the range [-300,299]

% Simulated overshoot, see left figure:
x1d = [-300:299];
[x,y]=meshgrid(x1d,x1d);
r = sqrt(x.^2+y.^2);
gsn1d = exp(-x1d.^2/500);
lowpass = @(x)(x1d > -x & x1d < x);
gsn1dcutoff = ifft(fftshift(lowpass(10).*fftshift(fft(gsn1d))));
plot(gsn1dcutoff)
OTF2d = reshape(interp1(x1d,gsn1dcutoff,r(:),'spline'),[length(x1d),length(x1d)]);
mesh(OTF2d)

% Quick and dirty fix, see right figure:
x1dExtended = linspace(min(x1d*sqrt(2)),max(x1d*sqrt(2)),ceil(length(x1d)*sqrt(2)));
gsn1dE = exp(-x1dExtended.^2/500);
% ^^^ note that this has 600*sqrt(2) points and is defined on the diagonal of your square.   Now we can low-pass filter in the freq. domain to add ripple in space domain:
lowpass = @(x)(x1dExtended > -x & x1dExtended < x);
gsn1dcutoff = -real(ifft(fftshift(lowpass(10).*fftshift(fft(gsn1dE)))));
plot(gsn1dcutoff)
OTF2d = reshape(interp1(x1dExtended,gsn1dcutoff,r(:),'spline'),[length(x1d),length(x1d)]);
mesh(OTF2d)

alt text http://img54.imageshack.us/img54/8255/clipboard01vz.png

沙沙粒小 2024-09-01 11:09:03

利奥的诊断是正确的。我想建议一个更简单(我希望)的补救措施:做你想做的事(基本上是围绕其对称轴旋转高斯)并在 600x600 正方形中获得合理的答案,你需要一个高斯 600*sqrt (2)=849 像素长。如果你能做到这一点,那么所有进一步的http://stackoverflow.com/questions/2443046/interpolating-1d-gaussian-into-2d-gaussianrickery都是不必要的。

编辑:换句话说:如果围绕其中心旋转 600 像素宽的物体,就会得到一个直径 600 像素的圆。您想要覆盖 600x600正方形。为此,您需要一个直径为 849 像素的圆,因为这是正方形的对角线。

Leo is right in his diagnosis. I'd like to suggest a simpler (I hope) remedy: to do what you want to do (which is basically rotate a Gaussian around its symmetry axis) and to get a reasonable answer in a 600x600 square you need a Gaussian 600*sqrt(2)=849 pixels long. If you can do that, then all further thttp://stackoverflow.com/questions/2443046/interpolating-1d-gaussian-into-2d-gaussianrickery is not necessary.

Edit: In other words: if you rotate something 600 pixels wide around its center, you get a circle 600 pixels in diameter. You want to cover a 600x600 square. For that you need a circle 849 pixels in diameter, since this is the diagonal of the square.

浅浅淡淡 2024-09-01 11:09:03

在高斯的特殊情况下,您可以利用它是可分离的这一事实来计算高斯:

OTF2(x,y) = exp( - x^2 - y^2) = exp( - x^2) * exp( - y^2) = OTFx(x) * OTFx(y)

因此您仍然只需要在内存中存储 OTFx 即可。

In the particular case of the gaussian, you may compute the gaussian by using the fact that it is separable:

OTF2(x,y) = exp( - x^2 - y^2) = exp( - x^2) * exp( - y^2) = OTFx(x) * OTFx(y)

So you still just need to store just OTFx in memory.

~没有更多了~
我们使用 Cookies 和其他技术来定制您的体验包括您的登录状态等。通过阅读我们的 隐私政策 了解更多相关信息。 单击 接受 或继续使用网站,即表示您同意使用 Cookies 和您的相关数据。
原文