标签: great-circle

我有一个工作的PHP脚本,它获取经度和纬度值,然后将它们输入到MySQL查询中.我想把它做成MySQL.这是我目前的PHP代码:

if ($distance != "Any" && $customer_zip != "") { //get the great circle distance

    //get the origin zip code info
    $zip_sql = "SELECT * FROM zip_code WHERE zip_code = '$customer_zip'";
    $result = mysql_query($zip_sql);
    $row = mysql_fetch_array($result);
    $origin_lat = $row['lat'];
    $origin_lon = $row['lon'];

    //get the range
    $lat_range = $distance/69.172;
    $lon_range = abs($distance/(cos($details[0]) * 69.172));
    $min_lat = number_format($origin_lat - $lat_range, "4", ".", "");
    $max_lat = number_format($origin_lat + $lat_range, "4", ".", "");
    $min_lon = number_format($origin_lon - $lon_range, "4", ".", "");
    $max_lon = …

我正在寻找一个能够对lat/lon数据集进行空间计算的Java库.以下是我正在寻找的一些功能:

• 计算两点之间的大圆距离
• 确定点是否位于简单的闭合多边形内,其中多边形由有序的点列表定义
• 确定点"A"和点"B"之间的直线是否与简单的闭合多边形相交,其中多边形由有序的点列表定义
• 确定点"A"是否在点"B"的某个半径范围内

我不想要的:

• 我不想要一个依赖于数据库地理空间组件的库,例如Oracle Spatial,并且不能独立运行.
• 我不想要一个图书馆,其目的是生成图形/地图等.我正在为现有应用程序构建分析模块,最终目标不是创建漂亮的图片.
• 我不想要一个库来搜索大量的空间数据.如果它也碰巧这样做,那没关系,但我不会使用该功能.

该组织计划最终获得Oracle Spatial的许可证(因此,此时将覆盖空间搜索),但是现在我需要在不依赖数据库提供的空间的情况下在小数据集上实现我上面提到的分析功能.支持.

我正在尝试将日间/夜间线路添加到世界地图ggplot中,以指示日夜区域; 这样的事情:

计划是在24小时的周期内为我的地图制作动画,如下所示:

上面的动画是用正弦波实现的,我知道它是完全不准确的.我知道这geosphere::gcIntermediate可以让我绘制大圆线,就像这样:

library(ggplot2)
library(ggthemes)
library(geosphere)

sunPath1 <- data.frame(gcIntermediate(c(-179, -30), c(0, 30), n=100))
sunPath2 <- data.frame(gcIntermediate(c(0, 30), c(179, -30), n=100))
sunPath <- rbind(sunPath1, sunPath2)

ggplot(sunPath) +
  borders("world", colour = "gray95", fill = "gray90") +
  geom_ribbon(aes(lon, ymax = lat), ymin=-180, fill="black", alpha=0.2) +
  theme_map()

虽然我不确定是否有可能在一年中的不同点画出所需的线条,例如在三月看起来如此:

我没有找到解决方案的运气,但我猜我不知道要搜索的关键字,因为这超出了我的知识范围.我认为答案可能在日出方程中的某处,但我不知道如何应用这些来找到解决方案,我也不知道如何在一年中改变这些参数.这个网站(用于上面的情节)似乎也很有用,但我还不确定如何!

我熟悉计算两点之间的大圆距离的公式.

即

<?php
$theta = $lon1 - $lon2; 
$dist = sin(deg2rad($lat1)) * sin(deg2rad($lat2)) +  cos(deg2rad($lat1)) * cos(deg2rad($lat2)) * cos(deg2rad($theta)); 
$dist = acos($dist); 
$dist = rad2deg($dist); 
//convert degrees to distance depending on units desired
?>

我需要的是与此相反的.给定起点,距离和简单的基数NSEW方向,以计算目标点的位置.自从我上数学课以来已经很久了.;)

所以我有一个表中有一堆不同的地址.我需要一个proc,它将选择该表中与传入的lat/long值相距指定距离的地址.

我的表的例子如下:

- messageId
- lat (float)
- long (float)

Proc正在传入另一个纬度/长度对(同时float也是一对)以及一个int(英里)

我发现这个http://www.sqlteam.com/forums/topic.asp?TOPIC_ID=81360来计算实际公式,但我无法弄清楚在proc表单中修改它以便能够通过整个列表地址,只给我一个Id<=里程(我传入的)的地址,来自我传入的纬度/经度.

我可以在这里得到任何帮助吗？

谢谢!

我在寻找python的地理库.我需要能够做到以下几点:

1. 使用大圆距离(不是线性距离计算)获取2点(以米为单位)之间的距离
2. 检查点是否在多边形内
3. 每秒执行1到2次,每次数千次

一开始我看了这篇文章: Python模块用于存储和查询地理坐标并开始使用geopy.我遇到了两个问题:

1. Geopy不支持多边形
2. geoPy的高CPU使用率(计算一个点与相对5000点之间的距离需要大约140ms的CPU)

我一直在寻找并发现最佳Python GIS库？和https://gis.stackexchange.com/.它看起来很有希望因为地理位置使用的是编译的C代码,它应该更快,并且可以更好地支持多边形.问题是geos/OGR执行线性距离计算而不是球体.这消除了所有其他基于地理的模块(如GEODjango和形状).我在这里错过了什么吗？我不认为我是第一个使用python进行GIS计算并希望获得准确结果的人.

我在流动数据大圆映射教程中学到了一些我学到的东西,并将它们与评论中链接的代码结合起来,以防止在R绘制跨赤道大圆时发生奇怪的事情.这给了我这个:

airports <- read.csv("/home/geoff/Desktop/DissertationData/airports.csv", header=TRUE) 
flights <- read.csv("/home/geoff/Desktop/DissertationData/ATL.csv", header=TRUE, as.is=TRUE)

library(maps)
library(geosphere)

checkDateLine <- function(l){
  n<-0
  k<-length(l)
  k<-k-1
  for (j in 1:k){
    n[j] <- l[j+1] - l[j]
  }
  n <- abs(n)
  m<-max(n, rm.na=TRUE)
  ifelse(m > 30, TRUE, FALSE)
}
clean.Inter <- function(p1, p2, n, addStartEnd){
  inter <- gcIntermediate(p1, p2, n=n, addStartEnd=addStartEnd)
  if (checkDateLine(inter[,1])){
    m1 <- midPoint(p1, p2)
    m1[,1] <- (m1[,1]+180)%%360 - 180
    a1 <- antipode(m1)
    l1 <- gcIntermediate(p1, a1, n=n, addStartEnd=addStartEnd)
    l2 <- gcIntermediate(a1, p2, n=n, addStartEnd=addStartEnd)
    l3 …

我正在建立一个小程序,从用户那里获取2个地理坐标,然后计算它们之间的距离(考虑到地球的曲率).所以我查看维基百科,了解这里的公式.

基本上我基本上设置了我的python函数,这就是我想出的:

def geocalc(start_lat, start_long, end_lat, end_long):
    start_lat = math.radians(start_lat)
    start_long = math.radians(start_long)
    end_lat = math.radians(end_long)
    end_long = math.radians(end_long)

    d_lat = start_lat - end_lat
    d_long = start_long - end_long

    EARTH_R = 6372.8

    c = math.atan((math.sqrt( (math.cos(end_lat)*d_long)**2 +( (math.cos(start_lat)*math.sin(end_lat)) - (math.sin(start_lat)*math.cos(end_lat)*math.cos(d_long)))**2)) / ((math.sin(start_lat)*math.sin(end_lat)) + (math.cos(start_lat)*math.cos(end_lat)*math.cos(d_long))) )

    return EARTH_R*c

问题是结果真的不准确.我是python的新手,所以一些帮助或建议将不胜感激!

问题,简而言之:如何使用Google Maps API之外的其他东西在Mercator中准确地投影包裹大圆半径？

问题,长期:

所以我有一个难题.我运行了一个地图应用程序,它使用谷歌地图API将巨大的圆圈投射到墨卡托地图上 - 它试图显示非常大的,准确的半径,比如大约13,000公里.但我不想再使用谷歌地图API,因为谷歌的新定价方案是疯了.所以我试图将代码转换为Leaflet,Mapbox或任何非谷歌,并且没有任何东西可以正确处理这些圈子.

以下是Google Maps API处理非洲北部中心13,000公里半径的测地圆的方法:

这看起来很直观,但是是正确的.波浪形图案是由环绕地球的圆圈引起的.

D3.js可以在正交投影中正确渲染.所以这是在地球上用d3.geo.circle()在D3.js中渲染的相同圆圈,两次旋转:

这使2D-"波浪"图案更有意义,对吧？对.我喜欢它.完全符合我的科学传播目的和所有这些.

但是当我将我的代码转换为Leaflet时,它根本不起作用.为什么？因为Leaflet的圆圈类不是一个很棒的圆圈.相反,它似乎只是一个椭圆,它与纬度有点扭曲,但不是真正的测地线.相同的圆,相同的半径,相同的原点,我们得到这个:

太错了,错了!除了看起来完全不现实之外,它只是不正确 - 澳大利亚不会在这样的圆形范围内.这对我的申请很重要!这不可能.

好吧,我想,也许诀窍就是尝试实现我自己的大圆圈课程.我采用的方法是迭代圆点作为距离原点的距离,但是使用这个非常有用的网站上的"目标点给定距离和从起点承载"计算距离,然后将它们投影为多边形.传单.这就是我得到的结果:

这看起来很糟糕但实际上更接近于准确!我们得到了波浪效应,这是正确的.像我一样,你可能会问,"这里到底发生了什么？" 所以我做了一个允许我突出显示每个迭代点的版本:

你可以非常清楚地看到它正确地渲染了圆圈,但是多边形错误地连接了它.它应该做什么(人们可能天真地想到)是围绕墨卡托地图投影的多个实例包裹该波形图,而不是天真地将它们连接在顶部,而是将它们球形地连接起来.像这样粗糙的Photoshop渲染:

然后多边形将以一种方式"关闭",表明多边形上方的所有东西都包含在其中.

我不知道如何在Leaflet中实现这样的东西.或其他任何事情.也许我必须以某种方式自己处理原始SVG,考虑到缩放状态？或者其他的东西？在我走进那些奸诈的杂草之前,我想我会要求任何建议/想法/等等.也许有一些更明显的方法.

哦,我尝试了另外一件事:使用相同的d3.geo.circle构造函数,该构造函数在墨卡托/ Leaflet投影的正交投影中运行良好.它产生或多或少与我的"天真"Leaflet大圆实现相同的结果:

我想这是有希望的.但是,如果移动原点的经度,则D3.js版本以更奇怪的方式包装(D3.js为红色,我的Leaflet类为绿松石):

如果在D3.js中有某种方式可以改变其工作方式,我不会感到惊讶,但我没有完全放下D3.js兔子洞.我希望D3.js会让这个"更容易"(因为它是比Leaflet更完整的制图工具),所以我会继续研究这个.

我还没有尝试在Mapbox-gl中这样做(我想这是"尝试"列表中的下一个).

无论如何.谢谢阅读.重申一个问题:如何使用Google Maps API之外的其他内容准确地在Mercator中投影包裹大圆半径？

我正在研究R中的通勤旅行模式(起源 - 目的地)流程图.我所拥有的数据是通勤者的日常交易(Date,Card,Entry_lat,Entry_Long,Exit_Lat,Exit_Long).旅行路径可能类似(因为他们通勤上班).

我需要绘制这个map (great circles).如果原点和目的地相同 - 连接线的不间断应指示相同的原点 - 目的地.

structure(list(business_date = structure(c(17245, 17245, 17245, 
17245, 17245, 17245, 17245, 17245, 17245, 17245, 17245, 17245, 
17245, 17245, 17245, 17245, 17245, 17245, 17245, 17245, 17245, 
17245, 17245, 17245, 17245, 17245, 17245, 17245, 17245, 17245, 
17245, 17245, 17245, 17245, 17245, 17245, 17245, 17245, 17245, 
17245, 17245, 17245, 17245, 17245, 17245, 17245, 17245, 17245, 
17245, 17245, 17245, 17245, 17245, 17245, 17245, 17245, 17245, 
17245, 17245, 17245, 17245, 17245, 17245, 17245, …