球体中的平均值 [英] mean value in a sphere
问题描述
从图中可以看出，一些像素完全位于圆内，但有些像素仅部分位于圆内。部分在圈内的人也只需要对平均值作出部分贡献。像素是方形的。这将简化我希望的数学。
我可以计算从像素角到中心点的距离，由此可以发现像素在内部和外部进入。其余需要更正。但如何找到这个更正。
感谢Heath Raftery，问题解决了！ [/ edit]
p>
函数检查像素的哪一部分位于圆圈中：
函数[a] = incirc（x，y，r）
％仅处理圆的右上象限
如果x <0 || y <0，则错误（'只有正x，y' ）; end
sqrt积分（r ^ 2-x ^ 2）dx
F = @（x，r）（1/2）*（x * sqrt（r ^ 2-X ^ 2）+ R ^ 2 * ATAN（X / SQRT（R ^ 2-X ^ 2）））;
％找到角落位置
x = [x-0.5，x + 0.5];
y = [y-0.5，y + 0.5];
d = sqrt（x。^ 2 + y。^ 2）; （d）< r，a = 1; return; end％圆内
如果min（d）> r，a = 0; return; end ％外圆
％与边交点（r ^ 2 = x ^ 2 + y ^ 2）
inters = [sqrt（r ^ 2-y（1）^ 2），sqrt （R ^ 2-Y（2）^ 2），SQRT（R ^ 2-X（1）^ 2），SQRT（R ^ 2-X（2）^ 2）]; ％x（1）x（2）y（1）y（2）
％去除虚数和超范围交点
inters（imag（inters）〜= 0）= NaN;
inters（inters <1E-5）= NaN; （（x（1）< inters（1：2））&（inters（1：2）< x（2））），〜（ （Y（1） - ;半成品（3：4））及（半成品（3：4）< Y（2）））]）= NaN的;
idx = find（〜isnan（inters））;
如果numel（idx）〜= 2，错误（'需要圆与像素的两个交点'）;结束
％检查周边像素的面积
如果全部（idx （y（2），r）-F（y（1），r）） - x（1）; ％面积
elseif全部（idx == [3,4]）％2个x-边上的交点
a =（F（x（2），r）-F（x（1），r） ） - y（1）; ％area
elseif all（idx == [1,3]）％x-edge上的y边上的一个交点（左和底）
a =（F（inters（1），r） -F（x（1），r）） - （y（1）*（inters（1）-x（1）））;
elseif all（idx == [2,4]）在x-edge（top& right）上的y边上有一个交叉点
a =（inters（2）-x（1））+ （F（X（2）中，r）-F（半成品（2）中，r）） - （Y（1）*（X（2）-inters（2）））;
else
error（'geometry'）
end
a = real（a）;
如果a< 0 || a> 1
错误（'计算错误'）;
end
end
测试函数的脚本
M = ones（100）; ％data
M（1:50，：）= 0;
pos = [50.2,50];
r = 2;
％计算结果应该是
h = 50-pos（2）+0.5;
A = pi * r ^ 2;
wedge = acos（h / r）/ pi;
triangle = h * sqrt（r ^ 2-h ^ 2）;
res =（A *楔形三角形）/ A
S = 0; N = 0;
for i = 1：size（M，1）
for j = 1：size（M，2）
x = abs（j-pos（1））;
y = abs（i-pos（2））;
n = incirc（x，y，r）;
M_（i，j）= n;
S = S + M（i，j）* n;
N = N + n;
end
end
result = S / N
结果= 0.3425
您可以看到算法找到圆圈中的像素部分。
问题是缺少一个问题，但我会假设它不是如何计算像素是完全在圈内还是在圈外。这是一个相对简单的任务。也就是说，如果到中心的像素的最前角小于远离中心的半径，则像素完全位于内侧，并且如果像素的最靠近中心的角落大于半径，则像素完全位于外侧从中心开始。
圆周上像素的比例落在圆周内的问题要棘手得多。有两个基本的解决方案：
- 完全和难。
- 近似和容易一点。 / li>
- 近似和容易一点。 / li>
在这两种情况下，请注意水平和垂直对称意味着只需考虑右上象限。
然后，对于（1），将圆心转换为原点（0,0），并将圆周视为函数 y（x）= sqrt（r ^ 2 - X ^ 2）。然后，圈内重叠像素的面积是积分：
lockquote
积分（y（x） - y0，从x0到x1 ，相对于x）
其中y0是像素的底部坐标，x0是左边坐标，x1是右边坐标。
这个积分可以用三角函数和三角代数完全求解。
对于（2） ，只需在像素内生成一组随机点，然后计算其周长内有多少个随机点。随着集合变大，圆周内的点数与所有点数的比例接近圆周内像素的比例。
I'm trying to calculate the mean value of the pixels inside a circle. In the future this needs to be extended to 3D, but for now a 2D sollution would already help me out.
As can be seen in the image, some pixels are entirely inside the circle, but some are only partly inside the circle. The ones partly in the circle also need to contribute only partly to the mean value. The pixels are square. This will simplify the mathematics I hope.
I can calculate the distance from the pixelcorners to the central point, from this you can find the pixels enterly inside and enterly outside. The rest needs correction. But how to find this correction.
[edit] thanks to Heath Raftery the problem is solved! [/edit]
the integral of a circle with radius r
As an example: I want to know the average pixelvalue of pixels in this circle. I know it is 0.3425, since 34.25% of the circle has a value of 1 and the rest is 0.
Function to check what part of a pixel is in the circle:
function [ a ] = incirc( x,y,r )
%only handles the top right quadrant of a circle
if x<0||y<0,error('only positive x,y');end
%integral of sqrt(r^2-x^2) dx
F = @(x,r) (1/2)*(x*sqrt(r^2-x^2)+r^2*atan(x/sqrt(r^2-x^2)));
%find corner locations
x=[x-0.5,x+0.5];
y=[y-0.5,y+0.5];
d = sqrt(x.^2+y.^2); %distance to closed and furthest corner
if max(d)<r,a=1;return;end %inside circle
if min(d)>r,a=0;return;end %outside circle
%intersections with edges (r^2 = x^2+y^2)
inters = [sqrt(r^2-y(1)^2),sqrt(r^2-y(2)^2),sqrt(r^2-x(1)^2),sqrt(r^2-x(2)^2)]; %x(1) x(2) y(1) y(2)
%remove imaginary and out of range intersections
inters(imag(inters)~=0)=NaN;
inters(inters<1E-5)=NaN; %to find values that are zero
inters([~((x(1)<inters(1:2))&(inters(1:2)<x(2))),~((y(1)<inters(3:4))&(inters(3:4)<y(2)))])=NaN;
idx = find(~isnan(inters));
if numel(idx)~=2,error('need two intersections of circle with pixel');end
%check area of pixel inside circumference
if all(idx==[1,2]) %2 intersections on y-edge
a=(F(y(2),r)-F(y(1),r)) - x(1); %area
elseif all(idx==[3,4]) %2 intersections on x-edge
a=(F(x(2),r)-F(x(1),r)) - y(1); %area
elseif all(idx==[1,3]) %one intersection on y-edge one on x-edge (left&bottom)
a=(F(inters(1),r)-F(x(1),r))- (y(1)*(inters(1)-x(1)));
elseif all(idx==[2,4]) %one intersection on y-edge one on x-edge (top&right)
a=(inters(2)-x(1))+(F(x(2),r)-F(inters(2),r))-(y(1)*(x(2)-inters(2)));
else
error('geometry')
end
a=real(a);
if a<0||a>1
error('computational error');
end
end
Script to test the function
M = ones(100); %data
M(1:50,:)=0;
pos=[50.2,50];
r = 2;
%calculate what the result should be
h=50-pos(2)+0.5;
A=pi*r^2;
wedge = acos(h/r)/pi;
triangle = h*sqrt(r^2-h^2);
res=(A*wedge-triangle)/A
S=0;N=0;
for i = 1:size(M,1)
for j = 1:size(M,2)
x=abs(j-pos(1));
y=abs(i-pos(2));
n=incirc( x,y,r );
M_(i,j)=n;
S = S+M(i,j)*n;
N = N+n;
end
end
result = S/N
result = 0.3425
You can see the algorithm finds the part of the pixel in the circle.
The question is missing a question, but I'll assume that it's not how to calculate whether pixels are fully inside or outside the circle. That's a relatively simple task. That is, a pixel is fully inside if the furtherest corner of the pixel to the centre is less than a radius away from the centre, and a pixel is fully outside if the closest corner of the pixel to the centre is more than a radius away from the centre.
The question of what proportion of pixels on the circumference fall within the circumference is much trickier. There are two fundamental solutions:
- Exact and hard.
- Approximate and a bit easier.
In both cases, note the horizontal and vertical symmetry means only the top right quadrant need be considered.
Then, for (1), translate the circle centre to the origin (0, 0) and treat the circumference as the function y(x) = sqrt(r^2 - x^2). Then, the area of an overlapping pixel within the circle is the integral:
integral(y(x) - y0, from x0 to x1, with respect to x)
where y0 is the bottom coordinate of the pixel, x0 is the left coordinate and x1 is the right coordinate.
This integral can be solved exactly with a trigonometric identity and a trigonometric substitution.
For (2), just generate a set of random points within the pixel and count how many of them fall within the circumference. As the set gets larger, the proportion of points that fall within the circumference to the count of all point approaches the proportion of the pixel within the circumference.
这篇关于球体中的平均值的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！
