Module:Ordnance Survey coordinates/sandbox: Difference between revisions

Revision as of 17:32, 16 December 2023

Revision as of 17:38, 16 December 2023

Line 157:

phi = atan ( z2 / p2*(1.-e2) )

local limitn0 = 1. / b2101

local n0 = 11

local phi_alt

Line 166 ⟶ 165:

phi = atan ( (z2 + e2 * ny * sin(phi)) / p2 )

n0 = n0 - 1

until abs ( phi - phi_alt ) <= limit0.000000000001 or n0 <= 0

return {lat=rad2deg(phi),lon=rad2deg(lambda)}

Line 220 ⟶ 219:

local N0uk=northingGB

-- semi-major axis (in line to equator) 0.996012717 is yer scale @ central meridian

local auk=Airyradius OSGBglobe.semimajor*F0uk

--semi-minor axis (in line to poles)

local buk=6356256OSGBglobe.91semiminor*F0uk

-- flatness=a1-b1/(a1+b1)

local n1uk=0.00167322025032508731869331280635710896296

-- first eccentricity squared=2*f-f^2where f=(a1-b1)/a1

local e2uk=AiryeccOSGBglobe.ecc

local k=(n-N0uk)/auk+lambda0uk/dr

local nextcounter=0

Line 256 ⟶ 255:

local long=phi0uk+dr*(y1*j6-y1*y1*y1*j7+pow(y1,5.0)*j8-pow(y1,7.0)*j9)

local lat=k9*dr

local v=AiryradiusOSGBglobe.semimajor/sqrt(1.0-e2uk*pow(sin(k),2.0))

local cartxa=v*cos(k9)*cos(long/dr)

local cartya=v*cos(k9)*sin(long/dr)

Line 269 ⟶ 268:

local cartzb=542.06+roty*cartxa-rotx*cartya+(1.0+scale)*cartza

-- convert Cartesian to long/lat

local awgs84=wgs84radiusWGSglobe.semimajor

local bwgs84=6356752WGSglobe.3141semiminor

local e2wgs84=wgs84eccWGSglobe.ecc

local lambdaradwgs84=atan(cartyb/cartxb)

long=lambdaradwgs84*dr

Line 367 ⟶ 366:

local cartzb=564.557+roty*cartxa-rotx*cartya+(1.0+scale)*cartza

-- convert Cartesian to long/lat

local awgs84=wgs84radiusWGSglobe.semimajor

local bwgs84=6356752WGSglobe.3141semiminor

local e2wgs84=wgs84eccWGSglobe.ecc

local lambdaradwgs84=atan(cartyb/cartxb)

long=lambdaradwgs84*dr

Line 583 ⟶ 582:

local function find_M ( lat_rad )

local e = AiryeccOSGBglobe.ecc

return AiryradiusOSGBglobe.semimajor * ( ( 1 - e/4 - 3 * e * e/64

- 5 * e * e * e/256

) * lat_rad

Line 613 ⟶ 612:

local lat_rad = deg2rad( latitude )

local e = AiryeccOSGBglobe.ecc

local e_prime_sq = e / (1-e)

local v = AiryradiusOSGBglobe.semimajor / sqrt(1 - e * sin(lat_rad)*sin(lat_rad))

local T = pow( tan(lat_rad), 2)

local C = e_prime_sq * pow( cos(lat_rad), 2)