Восток Север Север широта долгота

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

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

5 ответов

Решение

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

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

Однако 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 градусов, если отрицательное. Эта ссылка содержит соответствующие формулы. (Для этого расчета вы можете использовать либо сферический закон косинусов, как указано в разделе "Точка назначения с учетом расстояния и пеленга от начальной точки", либо более точную прямую формулу Винсенти. Однако первая связанная страница не используйте формулу Haversine для этого расчета.)
  • Найдите широту данной точки, найдя точку назначения с учетом базовой точки, абсолютного значения направления на север и направления 0 градусов, если направление на север положительное, или 180 градусов, если оно отрицательное.

Я преобразовал реализацию 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

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

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

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

В 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, чтобы проанализировать данные, извлечь восточные / северные ГБ и преобразовать их в широту / долготу.

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

пожалуйста, запустите это на bigquery

Другие вопросы по тегам