2013年5月12日日曜日

MySQLの&GIS&二点間の半正矢式と距離。

Original post: http://anothermysqldba.blogspot.com/2013/05/mysql-gis-haversine-formula-and.html


MySQLは、人々がGISを考えるときに最初に頭に浮かぶのデータベースではありません。 以下に記載されているデータベースは、以下のとおりです。

  • オラクル
  • Microsoft SQL Serverの
  • IBM DB2
  • IBM Informixの
  • PostgreSQLの

http://webhelp.esri.com/arcgisserver/9.3/java/index.htm#ジオ/ types_of_geodatabases.htm
#geodatabases/determ-1479045992.htm http://webhelp.esri.com/arcgisserver/9.3/java/index.htm


MySQLは今いつかのためにGISで取り組んできた、それは、MySQL 4.1に戻ります:
上記マークマウンダーの記事は、GISで素晴らしい外観です。 このブログの記事は、式に集中しますが、私はMarkの素敵なポストを指摘したかったのです。

空間インデックスの使用は、MySQL 5.5上にあるか、上記念頭にInnoDBがデフォルトのストレージエンジンであるとことを続ければそう、MyISAMストレージエンジンで動作します。

> CREATE TABLE geom (
-> lat float(10,7) NOT NULL,
-> lon float(10,7) NOT NULL,
-> g GEOMETRY NOT NULL,
-> SPATIAL INDEX(g)
-> ) ENGINE=Innodb;
ERROR 1464 (HY000): The used table type doesn't support SPATIAL indexes
> CREATE TABLE geom (
-> lat float(10,7) NOT NULL,
-> lon float(10,7) NOT NULL,
-> g GEOMETRY NOT NULL,
-> SPATIAL INDEX(g)
-> ) ENGINE=MyISAM; 

だからマルコスキーマ設計で私は中国の一部の都市の緯度経度を移植。
このデータは、ここを経由して収集されたhttp://www.infoplease.com/ipa/A0001769.html

CREATE TABLE china (
cityname varchar(50) DEFAULT NULL,
lat float(10,7) NOT NULL,
lon float(10,7) NOT NULL,
g GEOMETRY NOT NULL,
SPATIAL INDEX(g)
) ENGINE=MyISAM; 

INSERT INTO china VALUES ('Beijing', 39.55, 116.25, GeomFromText('POINT(39.55 116.25)'));
INSERT INTO china VALUES ('Canton', 23.7, 113.15, GeomFromText('POINT(23.7 113.15)'));
INSERT INTO china VALUES ('Chongqing', 29.46, 106.34, GeomFromText('POINT(29.46, 106.34)'));
INSERT INTO china VALUES ('Hong Kong', 22.20, 114.11, GeomFromText('POINT( 22.20 114.11)'));
INSERT INTO china VALUES ('Shanghai', 31.10, 121.28, GeomFromText('POINT(31.10 121.28)'));

ただ、それはすべての仕事を確認するため....

select lat, lon from china;
+------------+-------------+
| lat | lon |
+------------+-------------+
| 39.5499992 | 116.2500000 |
| 23.7000008 | 113.1500015 |
| 22.2000008 | 114.1100006 |
| 31.1000004 | 121.2799988 |
+------------+-------------+


だから私たちは北京から香港までの距離を確認することができます。
クイックルックオンラインは 、それが1963キロですが、私たちのデータを確認してみましょうを教えてくれる。


今では距離式の世界に掘るための時間だった.. 私はGISを当てDBAではない、私はそれは問題を認めることはありません。 私が遭遇した数学の数式の答えの範囲をもっと楽しんで、そうでなければ私はもっと上のテーブル情報を使用していただろう。 代わりに、あなたのためだけのリファレンスです。 私には2つの緯度/長い点間の距離を計算するために広く議論の方法を見に行った。 あなたは、異なる関数、プロシージャ、これを計算するオンラインなどの多くを見つけるでしょう。

私は式をテストしたかったので、最初に私はいくつかの変数を設定する

SET @lat1 =39.55;
SET @long1 =116.25;
SET @lat2 =22.20;
SET @long2 =114.11; 

たとえば、次のように

このウェブサイトでは、距離のための機能を持っています。 私はあなたがあなたのための場所で区切りオプションでそれをカット&ペーストできるようにするには、以下のことを追加した。 しかし、それは動作しますか? ところで私も3959に地球の半径値を更新した。
http://www.sqlexamples.info/SPAT/mysql_distance.htm

delimiter //
CREATE FUNCTION fn_distance
(p_x1 FLOAT, p_y1 FLOAT, p_x2 FLOAT, p_y2 FLOAT)
RETURNS FLOAT
DETERMINISTIC
BEGIN
DECLARE v_dist FLOAT;
DECLARE A FLOAT; DECLARE B FLOAT;
DECLARE C FLOAT; DECLARE D FLOAT;
/*
returns distance calculation between two points in LAT-LONG coordinates
*/

SET v_dist = 0;

-- convert to radians
SET A = p_x1 / 57.29577951;
SET B = p_y1 / 57.29577951;
SET C = p_x2 / 57.29577951;
SET D = p_y2 / 57.29577951;

IF (A = C && B = D) THEN
SET v_dist = 0;
ELSEIF ((sin(A)*sin(C)+cos(A)*cos(C)*cos(B - D)) > 1) THEN
SET v_dist = 3959 * acos(1);
ELSE
SET v_dist = 3959 *acos(sin(A)*sin(C) + cos(A)*cos(C)*cos(B - D));
END IF;

SET v_dist = v_dist * 1.609;

/* return distance in km. */
RETURN v_dist;

END;
//
delimiter ;

> SELECT fn_distance (@lat1, @long1, @lat2 , @long2) AS dist_km;
+-------------------+
| dist_km |
+-------------------+
| 1939.5457763671875 |
+-------------------+

別のクエリのテストを示してい 

> SELECT ( GLength( LineString(( PointFromWKB( POINT( @lat1, @long1 ))), ( PointFromWKB( POINT( @lat2, @long2 ) ))))) * 100 AS distance;
+--------------------+
| distance |
+--------------------+
| 1748.1478770401545 |
+--------------------+ 

しかし、もう一つはここで見つけるhttp://www.posteet.com/view/1555を : 

set log_bin_trust_function_creators=TRUE;
DELIMITER |
CREATE FUNCTION GeoDistKM( lat1 FLOAT, lon1 FLOAT, lat2 FLOAT, lon2 FLOAT ) RETURNS float
BEGIN
DECLARE pi, q1, q2, q3 FLOAT;
DECLARE rads FLOAT DEFAULT 0;
SET pi = PI();
SET lat1 = lat1 * pi / 180;
SET lon1 = lon1 * pi / 180;
SET lat2 = lat2 * pi / 180;
SET lon2 = lon2 * pi / 180;
SET q1 = COS(lon1-lon2);
SET q2 = COS(lat1-lat2);
SET q3 = COS(lat1+lat2);
SET rads = ACOS( 0.5*((1.0+q1)*q2 - (1.0-q1)*q3) );
RETURN 6378.388 * rads;
END;
|
DELIMITER ;
select geodistkm(
 set log_bin_trust_function_creators=TRUE;
DELIMITER |
CREATE FUNCTION GeoDistKM( lat1 FLOAT, lon1 FLOAT, lat2 FLOAT, lon2 FLOAT ) RETURNS float
BEGIN
DECLARE pi, q1, q2, q3 FLOAT;
DECLARE rads FLOAT DEFAULT 0;
SET pi = PI();
SET lat1 = lat1 * pi / 180;
SET lon1 = lon1 * pi / 180;
SET lat2 = lat2 * pi / 180;
SET lon2 = lon2 * pi / 180;
SET q1 = COS(lon1-lon2);
SET q2 = COS(lat1-lat2);
SET q3 = COS(lat1+lat2);
SET rads = ACOS( 0.5*((1.0+q1)*q2 - (1.0-q1)*q3) );
RETURN 6378.388 * rads;
END;
|
DELIMITER ;
select geodistkm(
 @ LAT1、@ long1、@ LAT2、@ long2 ) as distance;
+----------------------------------------+
| distance |
+----------------------------------------+
| 1942.0909423828125 |
+----------------------------------------+ 
) as distance;
+----------------------------------------+
| distance |
+----------------------------------------+
| 1942.0909423828125 |
+----------------------------------------+ 


しかし、これはまだ間違っている。 北京から香港王までの距離は、このウェブサイトによると、彼の1963キロと1224.9マイルですので:http://www.timeanddate.com/worldclock/distances.html?n=102 。 だから、最初と最後の結果は非常に近くにあります。 

だから、多くの時間を無駄にした後、はい、私は、私はまだ私は簡単と比較できることを数式を使用して結果を望んでいたことを認めない。 

座標{LAT1、lon1}と{LAT2、lon2}を持つ2つのポイント間の大圏距離dは次式で与えられます。 
D = ACOS(SIN(LAT1)*罪(LAT2)+ COS(LAT1)* COS(LAT2)* COS(lon1-lon2)) 
> SELECT ACOS(
-> SIN(@lat1) * SIN(@lat2) + COS(@lat1) * COS(@lat2) * COS( @long1 - @long2 )
-> ) *1000 as Aviation_forumula_DISTANCE;
+----------------------------+
| Aviation_forumula_DISTANCE |
+----------------------------+
| 1923.0473470093848 |
+----------------------------+


別の式が半正矢式です 
> SELECT 3956* 2 * ASIN ( SQRT (POWER(SIN((@lat1 - @lat2)*pi()/180 / 2),2) + COS(@lat1 * pi()/180) * COS(@lat2 *pi()/180) * POWER(SIN((@long1 - @long2) *pi()/180 / 2), 2) ) ) as Haversine_Formula_distance;
+----------------------------+
| Haversine_Formula_distance |
+----------------------------+
| 1204.5222518763514 |
+----------------------------+ 

もう一度、あなたは非常に異なる結果を参照してください。 私は、 半正矢 でより良い結果を期待した 
だから、それは私が正しく実装するときに使用する正しい式であると信じる 半正矢 式に最も近いように私には表示されますので、私は GeoDistKM 機能と一緒に行ったこれらの式で遊んで後。 これは、ウィリアムズの情報に基づいて私が書いた航空式に僅差で続いている。 

1942年は、私は、彼らは同様にそれを正しく計算したと言うことですオンライン検索を介して集まった1963年の結果ではありませんが。 極で近づい一緒に来て、地球と緯度の曲、LONは異なる式に多少の誤差が可能になります。 だから私は、現在、このに固執します: 

select geodistkm( @ LAT1、@ long1、@ LAT2、@ long2 ) as distance;
+----------------------------------------+
| distance |
+----------------------------------------+
| 1942.0909423828125 |
+----------------------------------------+
 ) as distance;
+----------------------------------------+
| distance |
+----------------------------------------+
| 1942.0909423828125 |
+----------------------------------------+