-
-
Save esisa/573b2f22e437a831cf40 to your computer and use it in GitHub Desktop.
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| /** | |
| Simple script to transform coordinates from WGS84 Geographic to UTM Zone 33N. | |
| Based on Holsen.js, which in turn is based on "Holsens småprogrammer" | |
| License: MIT | |
| Written by: Atle F. Sveen | |
| **/ | |
| var toUtm33 = (function () { | |
| 'use strict'; | |
| //make these easier to type/read | |
| var cos = Math.cos; | |
| var sin = Math.sin; | |
| var tan = Math.tan; | |
| var pow = Math.pow; | |
| var sqrt = Math.sqrt; | |
| //Ellipsoid constants for WGS84-ellipsoid | |
| var ellipsoid = { | |
| a: 6378137.000, | |
| f: 1 / 298.257223563 | |
| }; | |
| ellipsoid.b = -1 * ((ellipsoid.f * ellipsoid.a) - ellipsoid.a); | |
| //Values specific for Zone 33N | |
| var coordsystem = { | |
| scaleFactor: 0.9996, | |
| falseEasting: 500000, | |
| falseNorthing: 0, | |
| centralMeridian: 15, | |
| latitudeOfOrigin: 0 | |
| }; | |
| function toRad(deg) { | |
| return deg / (180 / Math.PI); | |
| } | |
| function ellipsoidParams(ellipsoid) { | |
| var n = (ellipsoid.a - ellipsoid.b) / (ellipsoid.a + ellipsoid.b); | |
| return { | |
| n: n, | |
| k: ellipsoid.a / (n + 1), | |
| k1: 1 + n * (n / 4) + pow(n, 4) / 64, | |
| k2: (n - n * n * n / 8) * 3, | |
| k3: (n * n - pow(n, 4) / 4) * 15 / 8 | |
| }; | |
| } | |
| function getGeodeticArch(lon1, lon2, ell) { | |
| var db = lon1 - lon2; | |
| var bm = lon1 - db / 2; | |
| return (ell.k * (ell.k1 * db - ell.k2 * cos(2 * bm) * sin(db) + | |
| ell.k3 * cos(4 * bm) * sin(2 * db)) - | |
| ell.k * (pow(ell.n, 3) * cos(6 * bm) * sin(3 * db) * 35 / 24) + | |
| ell.k * (pow(ell.n, 4) * cos(8 * bm) * sin(4 * db) * 315) / 256); | |
| } | |
| return function (lon, lat) { | |
| lon = toRad(lon) - toRad(coordsystem.centralMeridian); | |
| lat = toRad(lat); | |
| var a = ellipsoid.a; | |
| var b = ellipsoid.b; | |
| var et = (pow(a, 2) - pow(b, 2)) / pow(b, 2) * pow(cos(lat), 2); | |
| var n1 = pow(a, 2) / (sqrt(1 + et) * b); | |
| var t = tan(lat); | |
| var a1 = n1 * cos(lat); | |
| var a2 = -(n1 * t * pow(cos(lat), 2)) / 2.0; | |
| var a3 = -(n1 * pow(cos(lat), 3)) * (1 - pow(t, 2) + et) / 6; | |
| var a4 = n1 * t * pow(cos(lat), 4) * (5 - pow(t, 2) + 9 * et + 4 * pow(et, 2)) / 24; | |
| var a5 = n1 * pow(cos(lat), 5) * (5 - 18 * pow(t, 2) + pow(t, 4) + 14 * et - 58 * et * pow(t, 2)) / 120; | |
| var a6 = n1 * t * pow(cos(lat), 6) * (61 - 58 * pow(t, 2) + pow(t, 4) + 270 * et - 330 * et * pow(t, 2)) / 720; | |
| var north = -a2 * pow(lon, 2) + a4 * pow(lon, 4) + a6 * pow(lon, 6); | |
| var east = a1 * lon - a3 * pow(lon, 3) + a5 * pow(lon, 5); | |
| north = north + getGeodeticArch( | |
| lat, | |
| toRad(coordsystem.latitudeOfOrigin), | |
| ellipsoidParams(ellipsoid) | |
| ); | |
| //apply UTM "tricks" | |
| return { | |
| y: north * coordsystem.scaleFactor + coordsystem.falseNorthing, | |
| x: east * coordsystem.scaleFactor + coordsystem.falseEasting | |
| }; | |
| }; | |
| }()); |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment