De onde vem o raio padrão da Terra no ST_Distance_Sphere?

15

MySQL diz nos documentos para ST_Distance_Sphere

Os cálculos usam uma terra esférica e um raio configurável. O argumento do raio opcional deve ser dado em metros. Se omitido, o raio padrão é 6.370.986 metros. Se o argumento do raio estiver presente, mas não positivo, ER_WRONG_ARGUMENTSocorrerá um erro.

O PostGIS diz nos documentos de ST_Distance_Sphere(embora os documentos não sejam mais precisos )

Usa uma terra esférica e raio de 6370986 metros.

De onde eles tiraram os 6.370.986 metros padrão? WGS84 diz que o raio do eixo principal é 6.378.137,0 m. O PostGIS que agora usa um raio médio usa essencialmente 6371008.

Olhando o código

#define WGS84_MAJOR_AXIS 6378137.0
#define WGS84_INVERSE_FLATTENING 298.257223563
#define WGS84_MINOR_AXIS (WGS84_MAJOR_AXIS - WGS84_MAJOR_AXIS / WGS84_INVERSE_FLATTENING)
#define WGS84_RADIUS ((2.0 * WGS84_MAJOR_AXIS + WGS84_MINOR_AXIS ) / 3.0)

que significa

-- SELECT 6378137.0 - 6378137.0 / 298.257223563;
WGS84_MINOR_AXIS = 6356752.314245179498
-- SELECT ( 2.0 * 6378137.0 + ( 6378137.0 - 6378137.0 / 298.257223563) ) / 3.0;
WGS84_RADIUS = 6371008.771415059833

As versões mais recentes são muito menos eficientes, mais complexas e usam o Pro4j, mas parecem fazer a mesma coisa.

Ainda de onde vem o 6370986?

Evan Carroll
fonte
1
Ele representa o raio da Terra média, que deve ser (2*minorAxis+majorAxis)/3 ... embora esse valor para WGS84 ainda é a poucos metros maior (6,371,008.771)
JGH
Sim, essa é a pergunta por que a discrepância.
Evan Carroll
2
Algum desenvolvedor procurou na net? A fonte do postgis pode lançar alguma luz sobre ele
Ian Turton
2
@IanTurton A maioria dos erros pode ser reduzida para "algum desenvolvedor fez alguma coisa e a fonte pode esclarecer isso". Eu pretendia fazer o trabalho, imaginando que seria necessário se ninguém conhecesse a história. Veja a resposta abaixo.
Evan Carroll
1
Talvez tenha havido um erro de digitação e eles significaram 6370996 ... isso é muito próximo do raio autálico de Clarke 1866.
mkennedy

Respostas:

21

Ok, isso é hilarriuusss . Eu localizei isso. Em uma cópia antiga do lwgeom/lwgeom_spheroid.cPostGIS 1.0.0rc4, você pode ver isso,

/*
 * This algorithm was taken from the geo_distance function of the 
 * earthdistance package contributed by Bruno Wolff III.
 * It was altered to accept GEOMETRY objects and return results in
 * meters.
 */
PG_FUNCTION_INFO_V1(LWGEOM_distance_sphere);
Datum LWGEOM_distance_sphere(PG_FUNCTION_ARGS)
{
        const double EARTH_RADIUS = 6370986.884258304;

Passando para os documentos de earthdistance, você encontrará o seguinte:

Observe que, ao contrário da parte do módulo baseada em cubo, as unidades são conectadas aqui: alterar a earth()função não afetará os resultados desse operador.

E esse número conectado: EARTH_RADIUSpode ser visto aqui

/* Earth's radius is in statute miles. */
static const double EARTH_RADIUS = 3958.747716;

Então você pode fazer um simples.

EARTH_RADIUS * MILES_TO_METERS = EARTH_RADIUS_IN_METERS
 3958.747716 * 1609.344        = 6370986.884258304

E você tem o seu 6370986.884258304. Claro, apenas trunque isso e guarde-o em um longporque porque não.

Então, em essência, o raio no MySQL foi aumentado de um trabalho de cópia preguiçoso do PostGIS que converteu um raio em milhas em metros de uma constante obscura de um módulo aleatório PostgreSQL de 20 anos .

earth_distanceé um módulo pré-PostGIS de Bruce Momjian. Proclamo por este meio 6370986 a constante Bmomjian: uma boa aproximação da Terra em metros para satisfazer o MySQL. Embora talvez não por muito tempo.

Evan Carroll
fonte
2
Mas de onde veio esse número muito preciso 3958.747716? O mais próximo que posso encontrar é 3958.74795, que é o número de milhas de pesquisa nos EUA em 6371 quilômetros, mas que ainda resta cerca de 37 cm de
distância
1
@HenningMakholm continua lutando contra a boa luta, sem idéia. ;)
Evan Carroll
2
Very nice find!
Paul Ramsey