Подсчет точек внутри эллипса

Question

Я пытаюсь подсчитать заданные точки данных внутри каждого кольца эллипса:

Проблема в том, что у меня есть функция, которая проверяет, что: так для каждого эллипса,чтобы убедиться, что в ней есть точка, нужно вычислить три входа:

def get_focal_point(r1,r2,center_x):
    # f = square root of r1-squared - r2-squared
    focal_dist = sqrt((r1**2) - (r2**2))
    f1_x = center_x - focal_dist
    f2_x = center_x + focal_dist
    return f1_x, f2_x

def get_distance(f1,f2,center_y,t_x,t_y):
    d1 = sqrt(((f1-t_x)**2) + ((center_y - t_y)**2)) 
    d2 = sqrt(((f2-t_x)**2) + ((center_y - t_y)**2))
    return d1,d2

def in_ellipse(major_ax,d1,d2):
    if (d1+d2) <= 2*major_ax:
        return True
    else:
        return False

Сейчас я проверяю, находится ли он в эллипсе, по:

for i in range(len(data.latitude)):
    t_x = data.latitude[i] 
    t_y = data.longitude[i] 
    d1,d2 = get_distance(f1,f2,center_y,t_x,t_y)
    d1_array.append(d1)
    d2_array.append(d2)
    if in_ellipse(major_ax,d1,d2) == True:
        core_count += 1
        # if the point is not in core ellipse 
        # check the next ring up
    else:
        for i in range(loop):
            .....

Но тогда мне нужно было бы рассчитать каждую пару фокусных точек внешних петель. Есть ли более эффективный или умный способ сделать это?

Answer 1

Это может быть что-то похожее на то, что вы делаете.Я просто смотрю, если f (x, y) = x ^ 2 / r1 ^ 2 + y ^ 2 / r2 ^ 2 = 1.

Когда f (x, y) больше 1точка x, y находится вне эллипса.Когда оно меньше, оно находится внутри эллипса.Я перебираю каждый эллипс, чтобы найти тот, когда f (x, y) меньше 1.

Код также не учитывает эллипс, центрированный относительно начала координат.Это небольшое изменение, чтобы включить эту функцию.

import matplotlib.pyplot as plt
import matplotlib.patches as patches
import numpy as np

def inWhichEllipse(x,y,rads):
    '''
    With a list of (r1,r2) pairs, rads, return the index of the pair in which
    the point x,y resides. Return None as the index if it is outside all 
    Ellipses.
    '''
    xx = x*x
    yy = y*y

    count = 0
    ithEllipse =0
    while True:
        rx,ry = rads[count]
        ellips = xx/(rx*rx)+yy/(ry*ry)
        if ellips < 1:
            ithEllipse = count
            break
        count+=1
        if count >= len(rads):
            ithEllipse = None
            break

    return ithEllipse

rads = zip(np.arange(.5,10,.5),np.arange(.125,2.5,.25))

fig = plt.figure()
ax = fig.add_subplot(111)
ax.set_xlim(-15,15)
ax.set_ylim(-15,15)

# plot Ellipses
for rx,ry in rads:
    ellipse = patches.Ellipse((0,0),rx*2,ry*2,fc='none',ec='red')    
    ax.add_patch(ellipse)

x=3.0
y=1.0
idx = inWhichEllipse(x,y,rads)
rx,ry = rads[idx]
ellipse = patches.Ellipse((0,0),rx*2,ry*2,fc='none',ec='blue')    
ax.add_patch(ellipse)

if idx != None:
    circle = patches.Circle((x,y),.1)
    ax.add_patch(circle)

plt.show()

Этот код производит следующий рисунок:

Имейте в виду, это только отправная точка.Во-первых, вы можете изменить inWhichEllipse так, чтобы он принимал список квадратов пар r1 и r2, то есть (r1 * r1, r2 * r2), и это еще больше сократило бы вычисления.

Answer 2

Вы все усложняете. Нет необходимости вычислять фокусные точки и расстояния до фокусных точек и т. Д. В соответствии с геометрическим определением эллипса. Если вы знаете большую и меньшую ось (вы знаете), просто немного сожмите весь вопрос (например, чтобы оба были 1.0, путем деления x-centerx и y-centery на xaxis и yaxis), а затем задайте вопрос, является ли точка внутри эллипса просто

xnormalized**2 + ynormalized**2 <= 1

P.S .: В общем, хороший совет в этой области: нет sqrt, если вы можете сделать то же самое, фактически не вычисляя расстояние, но оставаясь комфортно в области его квадрата.

Answer 3

Вот несколько идей для вас:

• У вас есть правильная идея, перемещая код для вычисления фокусов за пределы цикла.
• Вычисления расстояния могут бытьускорение путем удаления квадратных корней.Другими словами, мы знаем, a < b подразумевает sqrt(a) < sqrt(b), поэтому нет необходимости вычислять квадратный корень.
• Если эллипсы концентрически, а главная ось параллельна оси X, вы можете упроститьот проблемы эллипса к проблеме окружности путем перерасчета значения x.

Кроме того, вот одна небольшая нота кодирования.Нет необходимости в выражении if для возврата True или False .Вместо этого вы можете вернуть само условное выражение:

def in_ellipse(major_ax,d1,d2):
    return (d1+d2) <= 2*major_ax: