--
local oscoord = local getArgs = require('Module:Arguments').getArgslocal yesno = require('Module:Yesno')local namespace = mw.title.getCurrentTitle.namespace
local pow = math.powlocal sqrt = math.sqrtlocal abs = math.abslocal floor = math.floorlocal ceil = math.ceillocal sin = math.sinlocal cos = math.coslocal tan = math.tanlocal atan = math.atanlocal deg2rad = math.radlocal rad2deg = math.deg
--
-- Datum parameters
local OSGBglobe =
local IEglobe =
local WGSglobe =
local WGS2OSGBparam =
local OSGB2WGSparam =
local IE2WGSparam =
local function HelmertDatumShift (latitude, longitude, param) -- latitude and longitude are in degrees -- return if latitude or longitude out of range
if abs(latitude) > 89. or abs(longitude) > 89. then return end
param = param or WGS2OSGBparam
-- parameters for ellipsoids -- Annex A.1, Ordnance Survey document --
local a1 = param.semimajor_src local b1 = param.semiminor_src local e1 = param.ecc_src
local a2 = param.semimajor_dst local b2 = param.semiminor_dst local e2 = param.ecc_dst
-- convert latitude and longitude to cartesian coordinates -- math in Annex B.1, Ordnance Survey document local phi = deg2rad (latitude) local lambda = deg2rad (longitude) local cosphi = cos (phi) local sinphi = sin (phi) local coslambda = cos (lambda) local sinlambda = sin (lambda)
local ny = a1 / sqrt (1. - e1*sinphi*sinphi)
local x1 = ny * cosphi * coslambda local y1 = ny * cosphi * sinlambda local z1 = ny * (1. - e1) * sinphi -- the helmert transformation proper -- Equation 3, Section 6.2, Ordnance Survey document
local x2 = x1 + param.tx + (param.s0 * x1 - param.rz * y1 + param.ry * z1) local y2 = y1 + param.ty + (param.rz * x1 + param.s0 * y1 - param.rx * z1) local z2 = z1 + param.tz + (-param.ry *x1 + param.rx * y1 + param.s0 * z1) -- convert cartesian coordinates to latitude and longitude -- math in Annex B.2, of Ordnance Survey document
lambda = atan (y2 / x2)
local p2 = sqrt (x2*x2 + y2*y2) phi = atan (z2 / p2*(1.-e2))
local n0 = 101
local phi_alt repeat phi_alt = phi ny = a2 / sqrt (1. - e2 * sin(phi) * sin(phi)) phi = atan ((z2 + e2 * ny * sin(phi)) / p2) n0 = n0 - 1 until abs (phi - phi_alt) <= 0.000000000001 or n0 <= 0
return
end
-- LAT LONG TO OS GRID LIBRARY RESUMES HERE
local function northeast(lett,num,shift) -- split into northings and eastings local le=mw.ustring.len(num) if le%2
and 0 or n e = e
local function EN2LL(e,n,datum) local a = datum.semimajor*datum.scale local b = datum.semiminor*datum.scale local lat0rad = deg2rad(datum.lat0) local n1 = datum.n1 local n12 = n1*n1 local n13 = n12*n1 local k=(n-datum.n0)/a+lat0rad local nextcounter=0 local j3, j4, j5, j6, m repeat nextcounter=nextcounter+1 local k3=k-lat0rad local k4=k+lat0rad j3=(1.0+n1+1.25*n12+1.25*n13)*k3 j4=(3.0*n1+3.0*n12+2.625*n13)*sin(k3)*cos(k4) j5=(1.875*n12+1.875*n13)*sin(2.0*k3)*cos(2.0*k4) j6=35.0/24.0*n13*sin(3.0*k3)*cos(3.0*k4) m=b*(j3-j4+j5-j6) k=k+(n-datum.n0-m)/a until abs(n-datum.n0-m)<=0.000000001 or nextcounter>=100 local x = 1.0-datum.ecc*pow(sin(k),2.0) local v=a/sqrt(x) local r=v*(1.0-datum.ecc)/x local h2=v/r-1.0 local y1=e-datum.e0 local tank = tan(k) local tank2 = tank*tank local tank4 = tank2*tank2 local tank6 = tank2*tank4 local cosk = cos(k) local yv = y1/v -- dimensionless quantity in series expansion local yv2 = yv*yv local yv3 = yv*yv2 local yv5 = yv3*yv2 local yv7 = yv5*yv2 j3=tank/(2.0*r) j4=tank/(24.0*r)*(5.0+3.0*tank2+h2-9.0*tank2*h2) j5=tank/(720.0*r)*(61.0+90.0*tank2+45.0*tank4) local k9=k-y1*(yv*j3+yv3*j4-yv5*j5) j6=1.0/(cosk) local j7=1.0/(cosk*6.0)*(v/r+2.0*tank2) local j8=1.0/(cosk*120.0)*(5.0+28.0*tank2+24.0*tank4) local j9=1.0/(cosk*5040.0) j9=j9*(61.0+662.0*tank2+1320.0*tank4+720.0*tank6) local long=datum.lon0+rad2deg(yv*j6-yv3*j7+yv5*j8-yv7*j9) local lat=rad2deg(k9) return end
local function GBEN2LL(e,n) local latlong = EN2LL(e,n,OSGBglobe) local helmert = HelmertDatumShift (latlong.lat, latlong.lon, OSGB2WGSparam) return end
local function GB2LL(lett,num) -- British OS to Lat+Long -- first caclulate e,n -- computing e and n exactly, to get SW corner of box local ne = northeast(lett,num) if ne.err then return end -- use British definition of e and n local e=500000.0*(ne.T1%5)+100000.0*(ne.T2%5)-1000000.0+ne.e*ne.pr local n=1900000.0-500000.0*floor(ne.T1/5)-100000.0*floor(ne.T2/5)+ne.n*ne.pr local result = GBEN2LL(e,n) result.prec = 0.8165*ne.pr return resultend local function IrishEN2LL(e,n) local latlong = EN2LL(e,n,IEglobe) local helmert = HelmertDatumShift (latlong.lat, latlong.lon, IE2WGSparam) return end
local function Irish2LL(lett,num) -- Irish OS to Lat+Long -- first caclulate e,n -- computing e and n exactly, to get SW corner of box local ne = northeast(lett,num) if ne.err then return end -- use Irish definition of northing and easting local e=100000.0*(ne.T1%5.0)+ne.e*ne.pr local n=ne.n*ne.pr+100000.0*(4.0-floor(ne.T1/5.0)) local result = IrishEN2LL(e,n) result.prec = 0.8165*ne.pr -- useful @ Commons return resultend
local function empty(s) return not s or s
local function NGR2LL(ngr) local result = local _ = 0 ngr, _ = mw.ustring.gsub(mw.ustring.upper(ngr),"[%s%p]","") local first, last, lett, num = mw.ustring.find(ngr,"^([A-Z]+)(%d+)$") if not first or empty(lett) or empty(num) or mw.ustring.len(lett) > 2 then return end if mw.ustring.len(lett)
local function split(s,sep)-- split a string s into chunks, separated by sep sep = sep or "%s" local t = for chunk in mw.ustring.gmatch(s,"([^"..sep.."]+)") do table.insert(t,chunk) end return tend
local function trim(s) local _ = 0 s, _ = mw.ustring.gsub(s,"^%s+","") s, _ = mw.ustring.gsub(s,"%s+$","") return send
local function alldigits(s) return not mw.ustring.find(s,"%D")end
local function warning(errmsg) local preview = require('Module:If preview') local msg = errmsg or 'Empty OS grid ref'
local html = preview._warning
if namespace
function oscoord._main(args) local input = args[1] if empty(input) then return warning(nil) end local linktitle = args[2] local namearg = args["name"] local rawurl = yesno(args["rawurl"]) local args = split(input,'_') local LL local restargs = 1 local current_page = mw.title.getCurrentTitle local pagename = mw.uri.encode(current_page.prefixedText, 'WIKI') if #args >= 2 and alldigits(args[2]) then if mw.ustring.sub(args[1],1,1)
function oscoord._oscoord(args) local output = '
' .. oscoord._main(args) .. '' if namespacefunction oscoord.main(frame) local args = getArgs(frame) return oscoord._main(args)end
function oscoord.oscoord(frame) local args = getArgs(frame) return oscoord._oscoord(args)end
--
local function find_M (lat_rad) local e = OSGBglobe.ecc local e2 = e*e local e3 = e*e*e
return OSGBglobe.semimajor * ((1 - e/4 - 3 * e2/64 - 5 * e3/256 ) * lat_rad - (3 * e/8 + 3 * e2/32 + 45 * e3/1024 ) * sin(2 * lat_rad) + (15 * e2/256 + 45 * e3/1024 ) * sin(4 * lat_rad) - (35 * e3/3072 ) * sin(6 * lat_rad))end
--
local longitude2 = longitude - floor((longitude + 180)/360) * 360
local lat_rad = deg2rad(latitude)
local e = OSGBglobe.ecc local e_prime_sq = e / (1-e)
local v = OSGBglobe.semimajor / sqrt(1 - e * sin(lat_rad)*sin(lat_rad)) local tank = tan(lat_rad) local T = tank*tank local T2 = T*T local C = e_prime_sq * pow(cos(lat_rad), 2) local A = deg2rad(longitude2 -OSGBglobe.lon0) * cos(lat_rad) local A2 = A*A local A3 = A2*A local A4 = A2*A2 local A5 = A3*A2 local A6 = A3*A3 local M = find_M(lat_rad) local M0 = 0.0 if latitude_origin ~= 0 then M0 = find_M(deg2rad(OSGBglobe.lat0)) end
local northing = OSGBglobe.n0 + OSGBglobe.scale * ((M - M0) + v*tan(lat_rad) * (A2/2 + (5 - T + 9*C + 4*C*C) * A4/24 + (61 - 58*T + T2 + 600*C - 330*e_prime_sq) * A6/720))
local easting = OSGBglobe.e0 + OSGBglobe.scale * v * (A + (1-T+C)*A3/6 + (5 - 18*T + T2 + 72*C - 58 * e_prime_sq)*A5/120)
return end
--
-- fix by Roger W Haworth local grid_x = floor(tm.easting / 100000) local grid_y = floor(tm.northing / 100000) if grid_x < 0 or grid_x > 6 or grid_y < 0 or grid_y > 12 then return end -- 0000000001111111111222222 -- 1234567890123456789012345 local letters = "ABCDEFGHJKLMNOPQRSTUVWXYZ" local indx1 = 18-floor(grid_y/5)*5+floor(grid_x/5) local indx2 = 21-(grid_y%5)*5+grid_x%5 local c1 = mw.ustring.sub(letters,indx1,indx1) local c2 = mw.ustring.sub(letters,indx2,indx2)
local easting = tm.easting%100000 local northing = tm.northing%100000 local grid = pow(10.0,5.0-prec) easting = floor(easting/grid) northing = floor(northing/grid) local fmt = '%0'..prec..'d' local e = mw.ustring.format(fmt, easting) local n = mw.ustring.format(fmt, northing)
return c1..c2..e..n
end
function oscoord._WGS2OSGB(lat,lon,prec) return LatLon2OSGB36(HelmertDatumShift(lat,lon,WGS2OSGBparam),prec)end
function oscoord.WGS2OSGB(frame) local args = getArgs(frame) return args[1] and args[2] and oscoord._WGS2OSGB(args[1],args[2],args[3]) or end
function oscoord.LL2OS(frame) local args = getArgs(frame) if not args[1] or not args[2] then return end local gridRef = oscoord._WGS2OSGB(args[1],args[2],args.prec) if not gridRef or #gridRef
return oscoord