地圖上繪制任意角度的橢圓
或者,如何選擇下班后去海灘的最佳方式 (Or, how to choose the best way to walk to the beach after work)
It was a cool autumn evening when Hila Kloper and I were thinking of going to the beach after work. The beach is about 2.5Km away from the office.
這是一個涼爽的秋天晚上,希拉·克洛珀和我想下班后去海灘。 海灘距離辦公室約2.5公里。
We were even considering strolling down the streets of Tel Aviv, willing to stretch our path to 3Km, and thinking to ourselves “mmm we wonder how far this stretch can take us?”
我們甚至正在考慮漫步在特拉維夫的街道上,愿意將我們的道路延伸到3Km,然后對自己進行思考:“嗯,我們想知道這條延伸能帶我們走多遠?”
Well, long story short, we didn’t go to the beach. Instead, we wrote a script that draws an ellipse around the office and the beach. The ellipse covers the city area we may go through if we ever decide to go to the beach after work.
好吧,長話短說,我們沒有去海灘。 相反,我們編寫了一個腳本,在辦公室和海灘周圍繪制了一個橢圓形。 如果我們決定下班后去海灘,那橢圓形可能覆蓋我們可能經歷的城市地區。
這是生命的橢圓 (It’s the Ellipse of Life)
Or - why should we care about ellipses?
或者-為什么我們要關心橢圓?
A circle is in some way the “natural” area around one point. An ellipse is the “natural” area around two points or a line. To name a few examples, bodies of mass move in elliptic orbits, ellipses represent the distortion caused by projecting a 3D map on 2D, and ellipses are also an accurate way to plot confidence of noisy GPS data (and confidence areas in 2D data in general).
在某種程度上,圓是圍繞一個點的“自然”區域。 橢圓是圍繞兩個點或一條線的“自然”區域。 僅舉幾個例子, 質量塊在橢圓軌道上移動 ,橢圓表示在2D上投影3D地圖所引起的變形 ,橢圓也是繪制嘈雜GPS數據 置信度 (以及通常在2D數據中置信度區域)的準確方法)。
In our case, we wanted to draw an area around the line starting at our office and ending at the beach. The easiest solution we found for drawing ellipses involved shapely and pyplot. It still required some modifications due to our GPS and map constraints.
在我們的案例中,我們想在直線周圍繪制一個區域,該區域從我們的辦公室開始,一直到海灘。 我們找到的繪制橢圓的最簡單解決方案涉及形狀和pyplot。 由于我們的GPS和地圖限制,仍然需要進行一些修改。
So, if you are here because you are looking for an easy copy-pastable code that draws an ellipse on a map — you can go to this repository we made. If you are also interested to learn how we found the complete solution to our problem, you are welcome to join us for the ride. We rediscover elementary geometry, learn about coordinates systems, and play around with some math code.
因此,如果您在這里是因為您正在尋找一個易于復制的可復制代碼,該代碼可在地圖上繪制橢圓形,那么您可以轉到我們創建的此存儲庫 。 如果您還想了解我們如何找到解決問題的完整解決方案,歡迎您加入我們。 我們重新發現基本幾何,了解坐標系,并試用一些數學代碼。
女孩只是想橢圓 (Girls Just Wanna Have Ellipses)
Us: “Oh there must be a package somewhere that can draw an ellipse on a map!”
我們: “哦,必須在某個地方可以在地圖上繪制橢圓的包裝!”
Us: “Oh there must be a package somewhere that can draw an ellipse on a map!”The Internet: “No there isn’t.”
我們: “哦,必須在某個地方可以在地圖上繪制橢圓的包裝!” 互聯網: “不,沒有。”
Us: “Oh there must be a package somewhere that can draw an ellipse on a map!”The Internet: “No there isn’t.”Us: “But there must be a simple copy-pastable code somewhere!”
我們: “哦,在某處必須有可以在地圖上繪制橢圓的包裝!” 互聯網: “不,沒有。” 我們: “但是某個地方必須有一個簡單的,可復制的代碼!”
Us: “Oh there must be a package somewhere that can draw an ellipse on a map!”The Internet: “No there isn’t.”Us: “But there must be a simple copy-pastable code somewhere!”Stackoverflow: “I have incomprehensible ones if you want.”
我們: “哦,在某處必須有可以在地圖上繪制橢圓的包裝!” 互聯網: “不,沒有。” 我們: “但是某個地方必須有一個簡單的,可復制的代碼!” Stackoverflow: “如果您愿意,我還有一些難以理解的內容。”
Us: “Oh there must be a package somewhere that can draw an ellipse on a map!”The Internet: “No there isn’t.”Us: “But there must be a simple copy-pastable code somewhere!”Stackoverflow: “I have incomprehensible ones if you want.”Us: “Well, OK, we’ll take an hour and make one ourselves!”
我們: “哦,在某處必須有可以在地圖上繪制橢圓的包裝!” 互聯網: “不,沒有。” 我們: “但是某個地方必須有一個簡單的,可復制的代碼!” Stackoverflow: “如果您愿意,我還有一些難以理解的內容。” 我們: “好吧,我們要花一個小時自己做一個!”
Us: “Oh there must be a package somewhere that can draw an ellipse on a map!”The Internet: “No there isn’t.”Us: “But there must be a simple copy-pastable code somewhere!”Stackoverflow: “I have incomprehensible ones if you want.”Us: “Well, OK, we’ll take an hour and make one ourselves!”Reality: “…”
我們: “哦,必須在某個地方可以在地圖上繪制橢圓的包裝!” 互聯網: “不,沒有。” 我們: “但是某個地方必須有一個簡單的,可復制的代碼!” Stackoverflow: “如果您愿意,我還有一些難以理解的內容。” 我們: “好吧,我們要花一個小時自己做一個!” 現實: “…”
Figuring out what an ellipse actually is was the first challenge.Wolfram Alpha told us that an ellipse is the set of points that have the same sum-of-distances from two mutual centers. Or something like that. Wolfram Alpha can be rather cryptic sometimes. But they have a gif so that’s nice.
弄清楚橢圓實際上是什么是第一個挑戰。 沃爾夫勒姆·阿爾法 ( Wolfram Alpha)告訴我們,橢圓是與兩個相互中心的距離之和相同的點集。 或類似的東西。 沃爾夫勒姆·阿爾法有時可能很神秘。 但是他們有一個gif,所以很好。
So all we have to do is make sure we have these inputs:
因此,我們要做的就是確保獲得以下輸入:
p1, p2
- GPS coordinates.p1, p2
-GPS坐標。r
- the radius which is actually the sum of distances from a point on the ellipse to the two centers.r
半徑,實際上是從橢圓上的點到兩個中心的距離之和。
Then follow this plan:
然后執行以下計劃:
Find
a
andb
, the axes of the ellipse.找到
a
和b
,橢圓的軸。Draw an ellipse around the origin
(0,0)
measured in meters.在以米為單位的原點
(0,0)
周圍繪制一個橢圓。- Move the ellipse to the center between the input GPS locations. 將橢圓移動到輸入GPS位置之間的中心。
- Rotate according to the angle between the input GPS locations. 根據輸入的GPS位置之間的角度旋轉。
And that’s it.
就是這樣。
Sounds simple enough, right?
聽起來很簡單,對吧?
漫長而曲折的道路(通往橢圓形) (The Long and Winding Road (That Leads to Ellipse))
步驟1.找到軸。 (Step 1. Find the axes.)
Here we needed to do some basic Pythagorean algebra. This image from Wikipedia was somewhat helpful:
在這里,我們需要做一些基本的畢達哥拉斯代數。 維基百科的這張圖片有些幫助:
Computing c
from GPS coordinates was easy thanks to haversine package.
多虧了hasersine軟件包,從GPS坐標計算c
很容易。
def GetEllipseAxisLengths(p1_lat, p1_lng, p2_lat, p2_lng,radius_in_meters):c2 = haversine((p1_lat, p1_lng), (p2_lat, p2_lng)) * 1000.0if radius_in_meters < c2:raise ValueError("Please specify radius larger than the distance between the two input points.")a = radius_in_meters / 2.0b = sqrt(pow(a, 2) - pow(c2 / 2.0, 2))return a, b
步驟2.在原點周圍繪制一個橢圓。 (Step 2. Draw an Ellipse Around the Origin.)
What we did here is that we took evenly spaced points on the x
axis and for each one found the two points on the ellipse that project to it:
我們在這里所做的是,我們在x
軸上獲取了均勻間隔的點,并且每個點在橢圓上找到了投射到該點的兩個點:
def GetEllipsePointInMeters(a, b, num_points):""":param a: length of "horizontal" axis in meters:param b: length of "vertical" axis in meters:param num_points: (half the) number of points to draw:return: List of tuples of perimeter points on the ellipse, centered around (0,0), in m."""x_points = list(np.linspace(-a, a, num_points))[1:-1]y_points_pos = [sqrt(pow(a, 2) - pow(x, 2)) * (float(b) / float(a))for x in x_points]y_points_neg = [-y for y in y_points_pos]perimeter_points_in_meters = [tuple([-a, 0])] + \[tuple([x, y]) for x, y in zip(x_points, y_points_pos)] + \[tuple([a, 0])] + \list(reversed([tuple([x, y]) for x, y in zip(x_points, y_points_neg)]))return perimeter_points_in_meters
第3步。如何將儀表添加到GPS? (Step 3. How Do You Even Add Meters to GPS?)
Well this was a tricky one, and the answer lies in understanding that
嗯,這是一個棘手的問題,答案在于了解
Fun fact: the Earth’s radius is around 6371000 meters on average!
有趣的事實:地球的半徑平均約為6371000米!
def AddMetersToPoint(center_lng, center_lat, dx, dy):""":param center_lng, center_lat: GPS coordinates of the center between the two input points.:param dx: distance to add to x-axis (lng) in meters:param dy: distance to add to y-axis (lat) in meters"""new_x = (center_lng + (dx / R_EARTH) * (180 / pi) / np.cos(center_lat * pi/180))new_y = center_lat + (dy / R_EARTH) * (180 / pi)return tuple([new_x, new_y])
步驟4.旋轉。 (Step 4. Rotate.)
This time, the wonders of the internet did not fail us (as they did on our major ellipse-drawing task). We found shapely package to do the rotation for us. The one trick to remember here is that you can’t rotate the points one by one. Rather you should form a shape first, and then rotate the entire shape.
這次,互聯網的奇跡并沒有使我們失望(就像他們在完成主要橢圓繪制任務時所做的那樣)。 我們發現勻稱的包裝可以為我們進行輪換。 這里要記住的一個技巧是不能一一旋轉點。 而是應該先形成一個形狀 ,然后再旋轉整個形狀 。
def GetEllipsePoints(p1_lat, p1_lng, p2_lat, p2_lng, perimeter_points_in_meters):"""Enter ellipse centers in lat-lng and ellipse perimeter points around the origin (0,0), and get points on the perimeter of the ellipse around the centers in lat-lng.:param p1_lat: lat coordinates of center point 1:param p1_lng: lng coordinates of center point 1:param p2_lat: lat coordinates of center point 2:param p2_lng: lng coordinates of center point 2:param perimeter_points_in_meters: List of tuples of perimeter points on the ellipse, centered around (0,0), in m.:return: List of the points we really want, tuples of (lat,lng)"""center_lng = (p1_lng + p2_lng) / 2.0center_lat = (p1_lat + p2_lat) / 2.0perimeter_points_in_lng_lat = \[AddMetersToPoint(center_lng, center_lat, p[0], p[1])for p in perimeter_points_in_meters]ellipse = LineString(perimeter_points_in_lng_lat)angle = degrees(atan2(p2_lat - p1_lat, p2_lng - p1_lng))ellipse_rotated = affinity.rotate(ellipse, angle)ellipse_points_lng_lat = list(ellipse_rotated.coords)ellipse_points = [tuple([p[1], p[0]]) for p in ellipse_points_lng_lat]return ellipse_points
驚喜! 步驟5.在s2 Map上繪制! (Surprise! Step 5. Draw on s2 Map!)
We wanted to present the ellipse nicely on an s2map. Apparently you can do that by opening the URL from inside your script. We used subprocess to do that.
我們想在s2map上很好地呈現橢圓。 顯然,您可以通過從腳本內部打開URL來實現。 我們使用子過程來做到這一點。
def OpenS2Map(points):url = \"http://s2map.com/#order=latlng&mode=polygon&s2=false" \"&points={}".format(str(points).replace(" ", ","))cmd = ["python", "-m", "webbrowser", "-t", url]subprocess.Popen(cmd, stdout=subprocess.PIPE, stderr=subprocess.STDOUT).communicate()
鄰居的橢圓更圓 (The Neighbor’s Ellipse Is Rounder)
You might notice our ellipse is not perfect. The points are evenly spaced on the axis between the centers, but they are not evenly spaced on the perimeter of the ellipse. The GPS->meters->GPS
transformation might result in loss of meters here and there. But hey, done is better than perfect, and we have to leave something to do for the next time we want to go to the beach, right?
您可能會注意到我們的橢圓不是完美的。 這些點在中心之間的軸上均勻分布,但在橢圓的周長上卻不均勻分布。 GPS->meters->GPS
轉換可能會導致此處和那里的電表丟失。 但是,嘿,做起來比完美還好,下次我們要去海灘時,我們必須做點事情,對嗎?
翻譯自: https://www.freecodecamp.org/news/a-total-ellipse-on-the-map-9e30d5235078/
地圖上繪制任意角度的橢圓