2011-10-24 2 views
12

나는 easting/northing 형식으로 위치 좌표를 가졌지 만 빙 표로 중심을 맞추려면 적절한 위도로 변환해야합니다. 모든 공식 또는 세부 정보를 east/northing을 lat/lon으로 변환하는 방법은 무엇입니까?위도 경도로 북쪽으로 이동

편집 : 구체적으로는, 나는이

답변

13

Eastings 같은 것을 사용하여 코드의 변환을 수행 할 수 있습니다 북쪽은 각각 기점의 동쪽과 북쪽 거리입니다. 기준점은 일반적으로 위도와 경도이며, eastings와 northings는 일반적으로 미터 또는 피트로 표시됩니다. 그러나 easting과 northing은 보통 그들을 긍정적으로 만들고 그들이베이스 포인트의 서쪽과 남쪽을 표현할 수 있도록 특별한 가치를 상쇄합니다.

일반적으로 하나의 좌표계에서 다른 좌표계로의 변환은 서로 다른 타원체 (지구 모델) 및 데이텀을 가질 수 있으므로 간단하지 않습니다. 필자가 이해하는 것처럼 한 좌표계에서 다른 좌표계로 변환하는 공식은 다소 복잡합니다.

SVY21 그러나 WGS84와 동일한 데이텀 및 타원체를 사용하므로 작업이 간단 해집니다. SVY21에서 eastings 및 northings의 기점은 Base 7 at Pierce Reservoir, 1 deg입니다. 22 분 02.9154 초. 북쪽과 103도 49 분 31.9752 초. (즉, 약 1.3674765 도의 위도와 약 103.8255487 도의 경도, 잘 알려진 텍스트는 각각 1.3666도 및 103.8333도를 사용합니다). easting의 오프셋은 28001.642 미터이고 northing의 오프셋은 38744.572 미터입니다. EPSG 코드는 3414입니다. 귀하의 eastings 및 northings는 미터로 표시됩니다. SVY21는 WGS84와 같은 시스템을 사용하기 때문에

, 당신이해야 할 모든은 다음과 같습니다

  • 는 각각의 오프셋 값에 의해 Y 좌표와 X 좌표를 뺍니다. 값은 미터로 표시됩니다.
  • 기준점, easting의 절대 값 및 easting이 양수인 경우 90 도의 방위 또는 270을 사용하여 대상 점을 찾아 주어진 점의 경도를 찾습니다. 네가 음수이면 This link에는 관련 수식이 들어 있습니다. (이 계산에서는 "주어진 거리와 시작점으로부터의 거리"섹션에서 주어진 것처럼 코사인의 구형 법칙을 사용하거나 더 정확하게 Vincenty's direct formula을 사용할 수 있습니다. 그러나 첫 번째 연결된 페이지는 이 계산에 대한 Haversine 공식)
  • 기준점, 북쪽 좌표의 절대 값, 북쪽 좌표가 양수이면 0도, 북쪽 좌표가 180 도인 경우의 대상 점을 찾아 주어진 점의 위도를 찾습니다. 그건 부정적이야.
+0

당신은 Haversine의 공식을 언급하고 있습니까? 동/북에서 위도/경도로가는 것은 어떻습니까? – Bahamut

+0

아니요, 저는 코사인의 구형 법칙을 언급하고 있습니다 만, Vincenty의 직접 공식에서도 사용할 수 있습니다. –

+0

나는 그것을 시험해 볼 것이다. – Bahamut

3

다른 수백 개의 좌표계가 있습니다 WGS84하기에 SVY21 좌표 변환해야 - Y 좌표/X 좌표와 위도/긴 유형의 좌표이 있지만, 그들은 ' 이러한 좌표를 얻는 시스템을 고유하게 식별하기에는 충분하지 않습니다.

소스 코드와 소스 코드 모두에 대해 EPSG 코드 (예 : 4326, 4269, 27700, 32701) 또는 공간 참조 시스템의 세부 정보 (데이텀, 프로젝션, 본초 자오선 및 측정 단위)가 필요합니다 선택한 대상 형식. 질문 제목에 "GPS"라고 말하면 원하는 위도/경도가 전역 위치 결정 시스템에서 사용하는 WGS84 데이터에 상대적으로 정의되어 있다고 가정하고 있지만 다른 데이터로 이어질 수있는 데이터의 예상이 여전히 많이 있습니다 동쪽/북쪽 가치.

사용한 투사의 세부 사항을 가지고 후에는 PROJ.4 라이브러리 (http://trac.osgeo.org/proj/)

+0

내가 얻는 결과는 SVY21 형식이며 WGS84로 변환하고 빙지도에 그려야합니다. 라이브러리 사용이 가능합니까? – Bahamut

2

펄에서 상대적으로 간단한 해결책이있다 :

것은 그래서, 우선, 펄이 설치되어 있는지 확인하십시오가. 그런 다음, 다음과 같은 4 개 개의 모듈을 설치합니다

지오 :: HelmertTransform 지리 :: NationalGrid CAM :: DBF mySociety :: GeoUtil

당신은 여러 가지 방법으로이 작업을 수행 할 수 있습니다. 여기에 어떻게 내가 해냈어 :

# Geo::HelmertTransform 
wget http://search.cpan.org/CPAN/authors/id/M/MY/MYSOCIETY/Geo-HelmertTransform-1.13.tar.gz 
tar xzf Geo-HelmertTransform-1.13.tar.gz 
perl Makefile.PL 
make 
make install 

# Geography::NationalGrid 
http://search.cpan.org/CPAN/authors/id/P/PK/PKENT/Geography-NationalGrid-1.6.tar.gz 
tar xzf Geography-NationalGrid-1.6.tar.gz 
perl Makefile.PL 
make 
make install 

# CAM::DBF 
wget http://search.cpan.org/CPAN/authors/id/C/CL/CLOTHO/CAM-DBF-1.02.tgz 
tar xzf CAM-DBF-1.02.tgz 
perl Makefile.PL 
make 
make install 

# mySociety::GeoUtil 
# See: http://parlvid.mysociety.org:81/os/ -> https://github.com/mysociety/commonlib/blob/master/perllib/mySociety/GeoUtil.pm 
mkdir -p mySociety 
wget -O mySociety/GeoUtil.pm 'https://raw.githubusercontent.com/mysociety/commonlib/master/perllib/mySociety/GeoUtil.pm' 
  1. 는 GB의 데이터를 가져옵니다.

을 클릭하고 지침에 따라 대 브리튼 "Code-Point® Open"데이터 세트를 다운로드하십시오. 당신이 codepo_gb.zip 다운로드하면 다음과 같이 당신이 그것을 추출 할 수 있습니다 :

압축 해제 codepo_gb.zip

을 압축 해제 된 파일이 현재 디렉토리에있다 염치, 당신은 다음 순서에 따라 펄 스크립트를 실행할 수 있습니다 데이터를 구문 분석하려면 GB eastings/northings를 추출하여 위도/경도로 변환하십시오.

use strict; 
use mySociety::GeoUtil qw/national_grid_to_wgs84/; 

while (<>) { 
    my @x=split(/,/); # split csv 
    my ($pc, $east, $north) = ($x[0], $x[10], $x[11]); 
    $pc=~s/\"//g; # remove quotes around postcode 
    my ($lat, $lng) = national_grid_to_wgs84($east, $north, "G"); # "G" means Great Britain 
    print "$pc,$lat,$lng\n"; 
} 

가 (호출하려면, 또한 기억 ... $를이 .pl 파일에 마지막 코드 블록을 저장 한 다음 perl script.pl your.csv 전화 X [0], $의 X [10]와 $ X [(11) ] 우편 번호, Y 좌표 각각 X 좌표의 열 번호이어야한다. 나는/경도 값을 위도 할 수있는 WGS84에 대한 T-SQL의 기능에 자바 스크립트 구현을 변환 한

전체 신용

+0

은 좋을 것이지만 나는 이식을위한 변환을위한 계산 및/또는 공식을 선호한다. 둘째, svy21은 싱가포르이므로 영국의 SVY에서 WGS84 로의 변환은 적용 할 수 없습니다. – Bahamut

+0

내 문제는, 바라건대 누군가를 도울 수 있지만, 해결책을 찾기 위해 서두르는 동안이 질문에 비틀 거리며 이것이 영국식 전환을위한 가장 쉬운 방법이었습니다. – rickyduck

+1

사실. 그래서 누군가가 이것을 유용하다고 생각하기 때문에 나는 그것을 downvote하지 않았다. – Bahamut

2

http://baroque.posterous.com/uk-postcode-latitudelongitude에.에 자유롭게 다른 좌표계가 필요한 경우 위스콘신 대학교 (University of Wisconsin)를 확인하십시오. 내가 소스로 사용하고 업데이트 된 상수를 얻는 그린 베이 웹 페이지.

drop function UF_utm_to_lat 
go 
create function UF_utm_to_lat(@utmz float, @x float, @y float) returns float 
as 
begin 
    --Based on code from this page: http://www.uwgb.edu/dutchs/usefuldata/ConvertUTMNoOZ.HTM 
    declare @latitude float; 
    declare @longitude float; 
    set @latitude = 0.00; 
    set @longitude = 0.00; 

    --Declarations 
    declare @a float; 
    declare @f float; 
    declare @drad float; 
    declare @k0 float; 
    declare @b float; 
    declare @e float; 
    declare @e0 float; 
    declare @esq float; 
    declare @e0sq float; 
    declare @zcm float; 
    declare @e1 float; 
    declare @M float; 
    declare @mu float; 
    declare @phi1 float; 
    declare @C1 float; 
    declare @T1 float; 
    declare @N1 float; 
    declare @R1 float; 
    declare @D float; 
    declare @phi float; 
    declare @lng float; 
    declare @lngd float; 

    --Datum Info here: Name, a, b, f, 1/f 
    --WGS 84 6,378,137.0 6356752.314 0.003352811 298.2572236 

    set @a = 6378137.0; 
    set @b = 6356752.314; 
    set @f = 0.003352811; 
    set @drad = PI()/180.0; 
    set @k0 = 0.9996; --scale on central meridian 

    set @e = SQRT(1.0 - (@b/@a)*(@b/@a)); --Eccentricity 
    --e = Math.sqrt(1 - (b/a)*(b/a));//eccentricity 
    set @e0 = @e/SQRT(1.0 - @e*@e); --Called e prime in reference 
    --e0 = e/Math.sqrt(1 - e*e);//Called e prime in reference 
    set @esq = (1.0 - (@b/@a)*(@b/@a)); --e squared for use in expansions 
    --esq = (1 - (b/a)*(b/a));//e squared for use in expansions 
    set @e0sq = @e*@e/([email protected]*@e); --e0 squared - always even powers 
    --e0sq = e*e/(1-e*e);// e0 squared - always even powers 
    set @zcm = 3.0 + 6.0*(@utmz-1.0) - 180.0; --Central meridian of zone 
    --zcm = 3 + 6*(utmz-1) - 180;//Central meridian of zone 
    set @e1 = (1.0 - SQRT(1.0 - @e*@e))/(1.0 + SQRT(1.0 - @e*@e)); --Called e1 in USGS PP 1395 also 
    --e1 = (1 - Math.sqrt(1 - e*e))/(1 + Math.sqrt(1 - e*e));//Called e1 in USGS PP 1395 also 
    set @M = 0.0 + @y/@k0; --Arc length along standard meridian 
    --M = M0 + y/k0;//Arc length along standard meridian. 
    set @mu = @M/(@a*(1.0 - @esq*(1.0/4.0 + @esq*(3.0/64.0 + 5.0*@esq/256.0)))); 
    --mu = M/(a*(1 - esq*(1/4 + esq*(3/64 + 5*esq/256)))); 
    set @phi1 = @mu + @e1*(3.0/2.0 - 27.0*@e1*@e1/32.0)*SIN(2.0*@mu) + @e1*@e1*(21.0/16.0 - 55.0*@e1*@e1/32.0)*SIN(4.0*@mu); --Footprint Latitude 
    --phi1 = mu + e1*(3/2 - 27*e1*e1/32)*Math.sin(2*mu) + e1*e1*(21/16 -55*e1*e1/32)*Math.sin(4*mu);//Footprint Latitude 
    set @phi1 = @phi1 + @e1*@e1*@e1*(SIN(6.0*@mu)*151.0/96.0 + @e1*SIN(8.0*@mu)*1097.0/512.0); 
    --phi1 = phi1 + e1*e1*e1*(Math.sin(6*mu)*151/96 + e1*Math.sin(8*mu)*1097/512); 
    set @C1 = @e0sq*POWER(COS(@phi1),2.0); 
    --C1 = e0sq*Math.pow(Math.cos(phi1),2); 
    set @T1 = POWER(TAN(@phi1),2.0); 
    --T1 = Math.pow(Math.tan(phi1),2); 
    set @N1 = @a/SQRT(1.0-POWER(@e*SIN(@phi1),2.0)); 
    --N1 = a/Math.sqrt(1-Math.pow(e*Math.sin(phi1),2)); 
    set @R1 = @N1*([email protected]*@e)/(1.0-POWER(@e*SIN(@phi1),2.0)); 
    --R1 = N1*(1-e*e)/(1-Math.pow(e*Math.sin(phi1),2)); 
    set @D = (@x-500000.0)/(@N1*@k0); 
    --D = (x-500000)/(N1*k0); 
    set @phi = (@D*@D)*(1.0/2.0 - @D*@D*(5.0 + 3.0*@T1 + 10.0*@C1 - 4.0*@C1*@C1 - 9.0*@e0sq)/24.0); 
    --phi = (D*D)*(1/2 - D*D*(5 + 3*T1 + 10*C1 - 4*C1*C1 - 9*e0sq)/24); 
    set @phi = @phi + POWER(@D,6.0)*(61.0 + 90.0*@T1 + 298.0*@C1 + 45.0*@T1*@T1 - 252.0*@e0sq - 3.0*@C1*@C1)/720.0; 
    --phi = phi + Math.pow(D,6)*(61 + 90*T1 + 298*C1 + 45*T1*T1 -252*e0sq - 3*C1*C1)/720; 
    set @phi = @phi1 - (@N1*TAN(@phi1)/@R1)*@phi; 
    --phi = phi1 - (N1*Math.tan(phi1)/R1)*phi; 


    set @latitude = FLOOR(1000000.0*@phi/@drad)/1000000.0; 

    set @lng = @D*(1.0 + @D*@D*((-1.0 - 2.0*@T1 - @C1)/6.0 + @D*@D*(5.0 - 2.0*@C1 + 28.0*@T1 - 3.0*@C1*@C1 + 8.0*@e0sq + 24.0*@T1*@T1)/120))/COS(@phi1); 
    set @lngd = @[email protected]/@drad; 
    set @longitude = FLOOR(1000000.0*@lngd)/1000000.0; 


    return @latitude; 
end 
go 
drop function UF_utm_to_long 
go 
create function UF_utm_to_long(@utmz float, @x float, @y float) returns float 
as 
begin 
    --Based on code from this page: http://www.uwgb.edu/dutchs/usefuldata/ConvertUTMNoOZ.HTM 
    declare @latitude float; 
    declare @longitude float; 
    set @latitude = 0.00; 
    set @longitude = 0.00; 

    --Declarations 
    declare @a float; 
    declare @f float; 
    declare @drad float; 
    declare @k0 float; 
    declare @b float; 
    declare @e float; 
    declare @e0 float; 
    declare @esq float; 
    declare @e0sq float; 
    declare @zcm float; 
    declare @e1 float; 
    declare @M float; 
    declare @mu float; 
    declare @phi1 float; 
    declare @C1 float; 
    declare @T1 float; 
    declare @N1 float; 
    declare @R1 float; 
    declare @D float; 
    declare @phi float; 
    declare @lng float; 
    declare @lngd float; 

    --Datum Info here: Name, a, b, f, 1/f 
    --WGS 84 6,378,137.0 6356752.314 0.003352811 298.2572236 

    set @a = 6378137.0; 
    set @b = 6356752.314; 
    set @f = 0.003352811; 
    set @drad = PI()/180.0; 
    set @k0 = 0.9996; --scale on central meridian 

    set @e = SQRT(1.0 - (@b/@a)*(@b/@a)); --Eccentricity 
    --e = Math.sqrt(1 - (b/a)*(b/a));//eccentricity 
    set @e0 = @e/SQRT(1.0 - @e*@e); --Called e prime in reference 
    --e0 = e/Math.sqrt(1 - e*e);//Called e prime in reference 
    set @esq = (1.0 - (@b/@a)*(@b/@a)); --e squared for use in expansions 
    --esq = (1 - (b/a)*(b/a));//e squared for use in expansions 
    set @e0sq = @e*@e/([email protected]*@e); --e0 squared - always even powers 
    --e0sq = e*e/(1-e*e);// e0 squared - always even powers 
    set @zcm = 3.0 + 6.0*(@utmz-1.0) - 180.0; --Central meridian of zone 
    --zcm = 3 + 6*(utmz-1) - 180;//Central meridian of zone 
    set @e1 = (1.0 - SQRT(1.0 - @e*@e))/(1.0 + SQRT(1.0 - @e*@e)); --Called e1 in USGS PP 1395 also 
    --e1 = (1 - Math.sqrt(1 - e*e))/(1 + Math.sqrt(1 - e*e));//Called e1 in USGS PP 1395 also 
    set @M = 0.0 + @y/@k0; --Arc length along standard meridian 
    --M = M0 + y/k0;//Arc length along standard meridian. 
    set @mu = @M/(@a*(1.0 - @esq*(1.0/4.0 + @esq*(3.0/64.0 + 5.0*@esq/256.0)))); 
    --mu = M/(a*(1 - esq*(1/4 + esq*(3/64 + 5*esq/256)))); 
    set @phi1 = @mu + @e1*(3.0/2.0 - 27.0*@e1*@e1/32.0)*SIN(2.0*@mu) + @e1*@e1*(21.0/16.0 - 55.0*@e1*@e1/32.0)*SIN(4.0*@mu); --Footprint Latitude 
    --phi1 = mu + e1*(3/2 - 27*e1*e1/32)*Math.sin(2*mu) + e1*e1*(21/16 -55*e1*e1/32)*Math.sin(4*mu);//Footprint Latitude 
    set @phi1 = @phi1 + @e1*@e1*@e1*(SIN(6.0*@mu)*151.0/96.0 + @e1*SIN(8.0*@mu)*1097.0/512.0); 
    --phi1 = phi1 + e1*e1*e1*(Math.sin(6*mu)*151/96 + e1*Math.sin(8*mu)*1097/512); 
    set @C1 = @e0sq*POWER(COS(@phi1),2.0); 
    --C1 = e0sq*Math.pow(Math.cos(phi1),2); 
    set @T1 = POWER(TAN(@phi1),2.0); 
    --T1 = Math.pow(Math.tan(phi1),2); 
    set @N1 = @a/SQRT(1.0-POWER(@e*SIN(@phi1),2.0)); 
    --N1 = a/Math.sqrt(1-Math.pow(e*Math.sin(phi1),2)); 
    set @R1 = @N1*([email protected]*@e)/(1.0-POWER(@e*SIN(@phi1),2.0)); 
    --R1 = N1*(1-e*e)/(1-Math.pow(e*Math.sin(phi1),2)); 
    set @D = (@x-500000.0)/(@N1*@k0); 
    --D = (x-500000)/(N1*k0); 
    set @phi = (@D*@D)*(1.0/2.0 - @D*@D*(5.0 + 3.0*@T1 + 10.0*@C1 - 4.0*@C1*@C1 - 9.0*@e0sq)/24.0); 
    --phi = (D*D)*(1/2 - D*D*(5 + 3*T1 + 10*C1 - 4*C1*C1 - 9*e0sq)/24); 
    set @phi = @phi + POWER(@D,6.0)*(61.0 + 90.0*@T1 + 298.0*@C1 + 45.0*@T1*@T1 - 252.0*@e0sq - 3.0*@C1*@C1)/720.0; 
    --phi = phi + Math.pow(D,6)*(61 + 90*T1 + 298*C1 + 45*T1*T1 -252*e0sq - 3*C1*C1)/720; 
    set @phi = @phi1 - (@N1*TAN(@phi1)/@R1)*@phi; 
    --phi = phi1 - (N1*Math.tan(phi1)/R1)*phi; 

    set @latitude = FLOOR(1000000.0*@phi/@drad)/1000000.0; 

    set @lng = @D*(1.0 + @D*@D*((-1.0 - 2.0*@T1 - @C1)/6.0 + @D*@D*(5.0 - 2.0*@C1 + 28.0*@T1 - 3.0*@C1*@C1 + 8.0*@e0sq + 24.0*@T1*@T1)/120))/COS(@phi1); 
    set @lngd = @[email protected]/@drad; 
    set @longitude = FLOOR(1000000.0*@lngd)/1000000.0; 


    return @longitude; 
end 
관련 문제