Восток север север широта долгота - PullRequest
14 голосов
/ 24 октября 2011

У меня есть координаты местоположения в формате восток / север, но мне нужно преобразовать его в правильный широту, чтобы центрировать его в картах Bing. Любая формула или детали, как конвертировать восток / север в широту / долготу?

РЕДАКТИРОВАТЬ: Чтобы быть более точным, мне нужно преобразовать координаты SVY21 в WGS84

Ответы [ 4 ]

14 голосов
/ 29 октября 2011

Восток и север - это расстояния на восток и север, соответственно, от базовой точки.Базовая точка обычно представляет собой широту и долготу, а восточные и северные направления обычно выражаются в метрах или футах.Восток и север, однако, обычно смещают определенное значение, чтобы сделать их положительными и позволить им выражать места к западу и югу от базовой точки.

В общем, преобразование из одной системы координат в другую не является простым, поскольку оба могут иметь разные эллипсоиды (модели Земли) и датумы.Как я понимаю, формулы для преобразования из одной системы координат в другую довольно сложны.

SVY21 , однако, использует те же данные и эллипсоид, что и WGS84, что упрощает задачу.В SVY21 базовая точка для восточных и северных направлений составляет База 7 в Пирсовом водохранилище , 1 град.22 мин.02,9154 сек.север и 103 град.49 мин. 31,9752 сек.восток (то есть, широта около 1,367765 градусов и долгота около 103,8255487 градусов; однако в хорошо известном тексте используются 1,3666 ... градусов и 103,8333 ... градусов соответственно).Смещение на восток составляет 28001,642 метра, а смещение на север - 38744,572 метра.Код EPSG - 3414. Я предполагаю, что ваши восточные и северные значения выражены в метрах.

Поскольку SVY21 использует ту же систему, что и WGS84, все, что вам нужно сделать, это:

  • Вычестьвосток и север на их соответствующие значения смещения.(Значения будут в метрах.)
  • Найдите долготу данной точки, найдя точку назначения с учетом базовой точки, абсолютное значение восточного направления и направление 90 градусов, если восточное направление положительноеили 270 градусов, если это отрицательно. Эта ссылка содержит соответствующие формулы.(Для этого расчета вы можете использовать либо сферический закон косинусов, как указано в разделе «Точка назначения с учетом расстояния и направления от начальной точки», либо более точную прямую формулу Винсенти . Первая связанная страницаоднако для этого расчета формула Хаверсайна не используется.)
  • Найдите широту данной точки, найдя точку назначения с учетом базовой точки, абсолютное значение северного направления и азимут 0градусов, если направление на север положительное, или 180 градусов, если оно отрицательное.
3 голосов
/ 24 октября 2011

Существуют сотни различных систем координат - Easting / Northing и Lat / Long имеют типов координат, но их недостаточно для однозначной идентификации системы, из которой получены эти координаты.

Вам нужен либо код EPSG (например, 4326, 4269, 27700, 32701), либо, альтернативно, детали системы пространственной привязки (данные, проекция, начальный меридиан и единица измерения) как для источника, так и для выбранный формат назначения. Вы упомянули «GPS» в заголовке своего вопроса, поэтому я предполагаю, что требуемая широта / долгота определяется относительно базовой точки WGS84, используемой системами глобального позиционирования, но все еще есть много проекций этой базовой точки, которые могут привести к разным восточным значениям. / Север значения.

Получив сведения об используемой проекции, вы можете выполнить преобразование в коде, используя что-то вроде библиотеки Proj.4 (http://trac.osgeo.org/proj/)

2 голосов
/ 08 марта 2015

Я преобразовал реализацию Javascript в функции T-SQL для WGS84 в значения широты / долготы. Не стесняйтесь использовать по своему усмотрению. Если вам нужна другая система координат, посетите веб-страницу Университета Висконсин - Грин-Бей, которую я использовал в качестве источника, и получите обновленные константы.

    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/(1.0-@e*@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*(1.0-@e*@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 = @zcm+@lng/@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/(1.0-@e*@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*(1.0-@e*@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 = @zcm+@lng/@drad;
    set @longitude = FLOOR(1000000.0*@lngd)/1000000.0;


    return @longitude;
end
2 голосов
/ 19 января 2012

В Perl есть относительно простое решение:

Итак, прежде всего, убедитесь, что у вас установлен Perl. Затем установите следующие четыре модуля:

Geo :: HelmertTransform Geography :: 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. Получить данные в ГБ.

Загрузите набор данных "Code-Point® Open" из Великобритании, нажав здесь и следуя инструкциям. После того, как вы скачали codepo_gb.zip Вы можете извлечь его следующим образом:

распаковать codepo_gb.zip

Предполагая, что разархивированные файлы теперь находятся в текущем каталоге, Затем вы можете запустить следующий perlscript для анализа данных, извлечь GB восточных / северных и преобразовать их в широта / долгота.

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] должны быть номера столбцов почтовых индексов, восток и север соответственно.

Полный кредит до http://baroque.posterous.com/uk-postcode-latitudelongitude

...