You can not select more than 25 topics Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.

7356 lines
220 KiB

This file contains ambiguous Unicode characters!

This file contains ambiguous Unicode characters that may be confused with others in your current locale. If your use case is intentional and legitimate, you can safely ignore this warning. Use the Escape button to highlight these characters.

(function (global, factory) {
typeof exports === 'object' && typeof module !== 'undefined' ? module.exports = factory() :
typeof define === 'function' && define.amd ? define(factory) :
(global.proj4 = factory());
}(this, (function () { 'use strict';
var globals = function(defs) {
defs('EPSG:4326', "+title=WGS 84 (long/lat) +proj=longlat +ellps=WGS84 +datum=WGS84 +units=degrees");
defs('EPSG:4269', "+title=NAD83 (long/lat) +proj=longlat +a=6378137.0 +b=6356752.31414036 +ellps=GRS80 +datum=NAD83 +units=degrees");
defs('EPSG:3857', "+title=WGS 84 / Pseudo-Mercator +proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null +no_defs");
defs.WGS84 = defs['EPSG:4326'];
defs['EPSG:3785'] = defs['EPSG:3857']; // maintain backward compat, official code is 3857
defs.GOOGLE = defs['EPSG:3857'];
defs['EPSG:900913'] = defs['EPSG:3857'];
defs['EPSG:102113'] = defs['EPSG:3857'];
};
var PJD_3PARAM = 1;
var PJD_7PARAM = 2;
var PJD_GRIDSHIFT = 3;
var PJD_WGS84 = 4; // WGS84 or equivalent
var PJD_NODATUM = 5; // WGS84 or equivalent
var SRS_WGS84_SEMIMAJOR = 6378137.0; // only used in grid shift transforms
var SRS_WGS84_SEMIMINOR = 6356752.314; // only used in grid shift transforms
var SRS_WGS84_ESQUARED = 0.0066943799901413165; // only used in grid shift transforms
var SEC_TO_RAD = 4.84813681109535993589914102357e-6;
var HALF_PI = Math.PI/2;
// ellipoid pj_set_ell.c
var SIXTH = 0.1666666666666666667;
/* 1/6 */
var RA4 = 0.04722222222222222222;
/* 17/360 */
var RA6 = 0.02215608465608465608;
var EPSLN = 1.0e-10;
// you'd think you could use Number.EPSILON above but that makes
// Mollweide get into an infinate loop.
var D2R = 0.01745329251994329577;
var R2D = 57.29577951308232088;
var FORTPI = Math.PI/4;
var TWO_PI = Math.PI * 2;
// SPI is slightly greater than Math.PI, so values that exceed the -180..180
// degree range by a tiny amount don't get wrapped. This prevents points that
// have drifted from their original location along the 180th meridian (due to
// floating point error) from changing their sign.
var SPI = 3.14159265359;
var exports$1 = {};
exports$1.greenwich = 0.0; //"0dE",
exports$1.lisbon = -9.131906111111; //"9d07'54.862\"W",
exports$1.paris = 2.337229166667; //"2d20'14.025\"E",
exports$1.bogota = -74.080916666667; //"74d04'51.3\"W",
exports$1.madrid = -3.687938888889; //"3d41'16.58\"W",
exports$1.rome = 12.452333333333; //"12d27'8.4\"E",
exports$1.bern = 7.439583333333; //"7d26'22.5\"E",
exports$1.jakarta = 106.807719444444; //"106d48'27.79\"E",
exports$1.ferro = -17.666666666667; //"17d40'W",
exports$1.brussels = 4.367975; //"4d22'4.71\"E",
exports$1.stockholm = 18.058277777778; //"18d3'29.8\"E",
exports$1.athens = 23.7163375; //"23d42'58.815\"E",
exports$1.oslo = 10.722916666667; //"10d43'22.5\"E"
var units = {
ft: {to_meter: 0.3048},
'us-ft': {to_meter: 1200 / 3937}
};
var ignoredChar = /[\s_\-\/\(\)]/g;
function match(obj, key) {
if (obj[key]) {
return obj[key];
}
var keys = Object.keys(obj);
var lkey = key.toLowerCase().replace(ignoredChar, '');
var i = -1;
var testkey, processedKey;
while (++i < keys.length) {
testkey = keys[i];
processedKey = testkey.toLowerCase().replace(ignoredChar, '');
if (processedKey === lkey) {
return obj[testkey];
}
}
}
var parseProj = function(defData) {
var self = {};
var paramObj = defData.split('+').map(function(v) {
return v.trim();
}).filter(function(a) {
return a;
}).reduce(function(p, a) {
var split = a.split('=');
split.push(true);
p[split[0].toLowerCase()] = split[1];
return p;
}, {});
var paramName, paramVal, paramOutname;
var params = {
proj: 'projName',
datum: 'datumCode',
rf: function(v) {
self.rf = parseFloat(v);
},
lat_0: function(v) {
self.lat0 = v * D2R;
},
lat_1: function(v) {
self.lat1 = v * D2R;
},
lat_2: function(v) {
self.lat2 = v * D2R;
},
lat_ts: function(v) {
self.lat_ts = v * D2R;
},
lon_0: function(v) {
self.long0 = v * D2R;
},
lon_1: function(v) {
self.long1 = v * D2R;
},
lon_2: function(v) {
self.long2 = v * D2R;
},
alpha: function(v) {
self.alpha = parseFloat(v) * D2R;
},
gamma: function(v) {
self.rectified_grid_angle = parseFloat(v);
},
lonc: function(v) {
self.longc = v * D2R;
},
x_0: function(v) {
self.x0 = parseFloat(v);
},
y_0: function(v) {
self.y0 = parseFloat(v);
},
k_0: function(v) {
self.k0 = parseFloat(v);
},
k: function(v) {
self.k0 = parseFloat(v);
},
a: function(v) {
self.a = parseFloat(v);
},
b: function(v) {
self.b = parseFloat(v);
},
r_a: function() {
self.R_A = true;
},
zone: function(v) {
self.zone = parseInt(v, 10);
},
south: function() {
self.utmSouth = true;
},
towgs84: function(v) {
self.datum_params = v.split(",").map(function(a) {
return parseFloat(a);
});
},
to_meter: function(v) {
self.to_meter = parseFloat(v);
},
units: function(v) {
self.units = v;
var unit = match(units, v);
if (unit) {
self.to_meter = unit.to_meter;
}
},
from_greenwich: function(v) {
self.from_greenwich = v * D2R;
},
pm: function(v) {
var pm = match(exports$1, v);
self.from_greenwich = (pm ? pm : parseFloat(v)) * D2R;
},
nadgrids: function(v) {
if (v === '@null') {
self.datumCode = 'none';
}
else {
self.nadgrids = v;
}
},
axis: function(v) {
var legalAxis = "ewnsud";
if (v.length === 3 && legalAxis.indexOf(v.substr(0, 1)) !== -1 && legalAxis.indexOf(v.substr(1, 1)) !== -1 && legalAxis.indexOf(v.substr(2, 1)) !== -1) {
self.axis = v;
}
},
approx: function() {
self.approx = true;
}
};
for (paramName in paramObj) {
paramVal = paramObj[paramName];
if (paramName in params) {
paramOutname = params[paramName];
if (typeof paramOutname === 'function') {
paramOutname(paramVal);
}
else {
self[paramOutname] = paramVal;
}
}
else {
self[paramName] = paramVal;
}
}
if(typeof self.datumCode === 'string' && self.datumCode !== "WGS84"){
self.datumCode = self.datumCode.toLowerCase();
}
return self;
};
var NEUTRAL = 1;
var KEYWORD = 2;
var NUMBER = 3;
var QUOTED = 4;
var AFTERQUOTE = 5;
var ENDED = -1;
var whitespace = /\s/;
var latin = /[A-Za-z]/;
var keyword = /[A-Za-z84_]/;
var endThings = /[,\]]/;
var digets = /[\d\.E\-\+]/;
// const ignoredChar = /[\s_\-\/\(\)]/g;
function Parser(text) {
if (typeof text !== 'string') {
throw new Error('not a string');
}
this.text = text.trim();
this.level = 0;
this.place = 0;
this.root = null;
this.stack = [];
this.currentObject = null;
this.state = NEUTRAL;
}
Parser.prototype.readCharicter = function() {
var char = this.text[this.place++];
if (this.state !== QUOTED) {
while (whitespace.test(char)) {
if (this.place >= this.text.length) {
return;
}
char = this.text[this.place++];
}
}
switch (this.state) {
case NEUTRAL:
return this.neutral(char);
case KEYWORD:
return this.keyword(char)
case QUOTED:
return this.quoted(char);
case AFTERQUOTE:
return this.afterquote(char);
case NUMBER:
return this.number(char);
case ENDED:
return;
}
};
Parser.prototype.afterquote = function(char) {
if (char === '"') {
this.word += '"';
this.state = QUOTED;
return;
}
if (endThings.test(char)) {
this.word = this.word.trim();
this.afterItem(char);
return;
}
throw new Error('havn\'t handled "' +char + '" in afterquote yet, index ' + this.place);
};
Parser.prototype.afterItem = function(char) {
if (char === ',') {
if (this.word !== null) {
this.currentObject.push(this.word);
}
this.word = null;
this.state = NEUTRAL;
return;
}
if (char === ']') {
this.level--;
if (this.word !== null) {
this.currentObject.push(this.word);
this.word = null;
}
this.state = NEUTRAL;
this.currentObject = this.stack.pop();
if (!this.currentObject) {
this.state = ENDED;
}
return;
}
};
Parser.prototype.number = function(char) {
if (digets.test(char)) {
this.word += char;
return;
}
if (endThings.test(char)) {
this.word = parseFloat(this.word);
this.afterItem(char);
return;
}
throw new Error('havn\'t handled "' +char + '" in number yet, index ' + this.place);
};
Parser.prototype.quoted = function(char) {
if (char === '"') {
this.state = AFTERQUOTE;
return;
}
this.word += char;
return;
};
Parser.prototype.keyword = function(char) {
if (keyword.test(char)) {
this.word += char;
return;
}
if (char === '[') {
var newObjects = [];
newObjects.push(this.word);
this.level++;
if (this.root === null) {
this.root = newObjects;
} else {
this.currentObject.push(newObjects);
}
this.stack.push(this.currentObject);
this.currentObject = newObjects;
this.state = NEUTRAL;
return;
}
if (endThings.test(char)) {
this.afterItem(char);
return;
}
throw new Error('havn\'t handled "' +char + '" in keyword yet, index ' + this.place);
};
Parser.prototype.neutral = function(char) {
if (latin.test(char)) {
this.word = char;
this.state = KEYWORD;
return;
}
if (char === '"') {
this.word = '';
this.state = QUOTED;
return;
}
if (digets.test(char)) {
this.word = char;
this.state = NUMBER;
return;
}
if (endThings.test(char)) {
this.afterItem(char);
return;
}
throw new Error('havn\'t handled "' +char + '" in neutral yet, index ' + this.place);
};
Parser.prototype.output = function() {
while (this.place < this.text.length) {
this.readCharicter();
}
if (this.state === ENDED) {
return this.root;
}
throw new Error('unable to parse string "' +this.text + '". State is ' + this.state);
};
function parseString(txt) {
var parser = new Parser(txt);
return parser.output();
}
function mapit(obj, key, value) {
if (Array.isArray(key)) {
value.unshift(key);
key = null;
}
var thing = key ? {} : obj;
var out = value.reduce(function(newObj, item) {
sExpr(item, newObj);
return newObj
}, thing);
if (key) {
obj[key] = out;
}
}
function sExpr(v, obj) {
if (!Array.isArray(v)) {
obj[v] = true;
return;
}
var key = v.shift();
if (key === 'PARAMETER') {
key = v.shift();
}
if (v.length === 1) {
if (Array.isArray(v[0])) {
obj[key] = {};
sExpr(v[0], obj[key]);
return;
}
obj[key] = v[0];
return;
}
if (!v.length) {
obj[key] = true;
return;
}
if (key === 'TOWGS84') {
obj[key] = v;
return;
}
if (key === 'AXIS') {
if (!(key in obj)) {
obj[key] = [];
}
obj[key].push(v);
return;
}
if (!Array.isArray(key)) {
obj[key] = {};
}
var i;
switch (key) {
case 'UNIT':
case 'PRIMEM':
case 'VERT_DATUM':
obj[key] = {
name: v[0].toLowerCase(),
convert: v[1]
};
if (v.length === 3) {
sExpr(v[2], obj[key]);
}
return;
case 'SPHEROID':
case 'ELLIPSOID':
obj[key] = {
name: v[0],
a: v[1],
rf: v[2]
};
if (v.length === 4) {
sExpr(v[3], obj[key]);
}
return;
case 'PROJECTEDCRS':
case 'PROJCRS':
case 'GEOGCS':
case 'GEOCCS':
case 'PROJCS':
case 'LOCAL_CS':
case 'GEODCRS':
case 'GEODETICCRS':
case 'GEODETICDATUM':
case 'EDATUM':
case 'ENGINEERINGDATUM':
case 'VERT_CS':
case 'VERTCRS':
case 'VERTICALCRS':
case 'COMPD_CS':
case 'COMPOUNDCRS':
case 'ENGINEERINGCRS':
case 'ENGCRS':
case 'FITTED_CS':
case 'LOCAL_DATUM':
case 'DATUM':
v[0] = ['name', v[0]];
mapit(obj, key, v);
return;
default:
i = -1;
while (++i < v.length) {
if (!Array.isArray(v[i])) {
return sExpr(v, obj[key]);
}
}
return mapit(obj, key, v);
}
}
var D2R$1 = 0.01745329251994329577;
function rename(obj, params) {
var outName = params[0];
var inName = params[1];
if (!(outName in obj) && (inName in obj)) {
obj[outName] = obj[inName];
if (params.length === 3) {
obj[outName] = params[2](obj[outName]);
}
}
}
function d2r(input) {
return input * D2R$1;
}
function cleanWKT(wkt) {
if (wkt.type === 'GEOGCS') {
wkt.projName = 'longlat';
} else if (wkt.type === 'LOCAL_CS') {
wkt.projName = 'identity';
wkt.local = true;
} else {
if (typeof wkt.PROJECTION === 'object') {
wkt.projName = Object.keys(wkt.PROJECTION)[0];
} else {
wkt.projName = wkt.PROJECTION;
}
}
if (wkt.AXIS) {
var axisOrder = '';
for (var i = 0, ii = wkt.AXIS.length; i < ii; ++i) {
var axis = [wkt.AXIS[i][0].toLowerCase(), wkt.AXIS[i][1].toLowerCase()];
if (axis[0].indexOf('north') !== -1 || ((axis[0] === 'y' || axis[0] === 'lat') && axis[1] === 'north')) {
axisOrder += 'n';
} else if (axis[0].indexOf('south') !== -1 || ((axis[0] === 'y' || axis[0] === 'lat') && axis[1] === 'south')) {
axisOrder += 's';
} else if (axis[0].indexOf('east') !== -1 || ((axis[0] === 'x' || axis[0] === 'lon') && axis[1] === 'east')) {
axisOrder += 'e';
} else if (axis[0].indexOf('west') !== -1 || ((axis[0] === 'x' || axis[0] === 'lon') && axis[1] === 'west')) {
axisOrder += 'w';
}
}
if (axisOrder.length === 2) {
axisOrder += 'u';
}
if (axisOrder.length === 3) {
wkt.axis = axisOrder;
}
}
if (wkt.UNIT) {
wkt.units = wkt.UNIT.name.toLowerCase();
if (wkt.units === 'metre') {
wkt.units = 'meter';
}
if (wkt.UNIT.convert) {
if (wkt.type === 'GEOGCS') {
if (wkt.DATUM && wkt.DATUM.SPHEROID) {
wkt.to_meter = wkt.UNIT.convert*wkt.DATUM.SPHEROID.a;
}
} else {
wkt.to_meter = wkt.UNIT.convert;
}
}
}
var geogcs = wkt.GEOGCS;
if (wkt.type === 'GEOGCS') {
geogcs = wkt;
}
if (geogcs) {
//if(wkt.GEOGCS.PRIMEM&&wkt.GEOGCS.PRIMEM.convert){
// wkt.from_greenwich=wkt.GEOGCS.PRIMEM.convert*D2R;
//}
if (geogcs.DATUM) {
wkt.datumCode = geogcs.DATUM.name.toLowerCase();
} else {
wkt.datumCode = geogcs.name.toLowerCase();
}
if (wkt.datumCode.slice(0, 2) === 'd_') {
wkt.datumCode = wkt.datumCode.slice(2);
}
if (wkt.datumCode === 'new_zealand_geodetic_datum_1949' || wkt.datumCode === 'new_zealand_1949') {
wkt.datumCode = 'nzgd49';
}
if (wkt.datumCode === 'wgs_1984' || wkt.datumCode === 'world_geodetic_system_1984') {
if (wkt.PROJECTION === 'Mercator_Auxiliary_Sphere') {
wkt.sphere = true;
}
wkt.datumCode = 'wgs84';
}
if (wkt.datumCode.slice(-6) === '_ferro') {
wkt.datumCode = wkt.datumCode.slice(0, - 6);
}
if (wkt.datumCode.slice(-8) === '_jakarta') {
wkt.datumCode = wkt.datumCode.slice(0, - 8);
}
if (~wkt.datumCode.indexOf('belge')) {
wkt.datumCode = 'rnb72';
}
if (geogcs.DATUM && geogcs.DATUM.SPHEROID) {
wkt.ellps = geogcs.DATUM.SPHEROID.name.replace('_19', '').replace(/[Cc]larke\_18/, 'clrk');
if (wkt.ellps.toLowerCase().slice(0, 13) === 'international') {
wkt.ellps = 'intl';
}
wkt.a = geogcs.DATUM.SPHEROID.a;
wkt.rf = parseFloat(geogcs.DATUM.SPHEROID.rf, 10);
}
if (geogcs.DATUM && geogcs.DATUM.TOWGS84) {
wkt.datum_params = geogcs.DATUM.TOWGS84;
}
if (~wkt.datumCode.indexOf('osgb_1936')) {
wkt.datumCode = 'osgb36';
}
if (~wkt.datumCode.indexOf('osni_1952')) {
wkt.datumCode = 'osni52';
}
if (~wkt.datumCode.indexOf('tm65')
|| ~wkt.datumCode.indexOf('geodetic_datum_of_1965')) {
wkt.datumCode = 'ire65';
}
if (wkt.datumCode === 'ch1903+') {
wkt.datumCode = 'ch1903';
}
if (~wkt.datumCode.indexOf('israel')) {
wkt.datumCode = 'isr93';
}
}
if (wkt.b && !isFinite(wkt.b)) {
wkt.b = wkt.a;
}
function toMeter(input) {
var ratio = wkt.to_meter || 1;
return input * ratio;
}
var renamer = function(a) {
return rename(wkt, a);
};
var list = [
['standard_parallel_1', 'Standard_Parallel_1'],
['standard_parallel_1', 'Latitude of 1st standard parallel'],
['standard_parallel_2', 'Standard_Parallel_2'],
['standard_parallel_2', 'Latitude of 2nd standard parallel'],
['false_easting', 'False_Easting'],
['false_easting', 'False easting'],
['false-easting', 'Easting at false origin'],
['false_northing', 'False_Northing'],
['false_northing', 'False northing'],
['false_northing', 'Northing at false origin'],
['central_meridian', 'Central_Meridian'],
['central_meridian', 'Longitude of natural origin'],
['central_meridian', 'Longitude of false origin'],
['latitude_of_origin', 'Latitude_Of_Origin'],
['latitude_of_origin', 'Central_Parallel'],
['latitude_of_origin', 'Latitude of natural origin'],
['latitude_of_origin', 'Latitude of false origin'],
['scale_factor', 'Scale_Factor'],
['k0', 'scale_factor'],
['latitude_of_center', 'Latitude_Of_Center'],
['latitude_of_center', 'Latitude_of_center'],
['lat0', 'latitude_of_center', d2r],
['longitude_of_center', 'Longitude_Of_Center'],
['longitude_of_center', 'Longitude_of_center'],
['longc', 'longitude_of_center', d2r],
['x0', 'false_easting', toMeter],
['y0', 'false_northing', toMeter],
['long0', 'central_meridian', d2r],
['lat0', 'latitude_of_origin', d2r],
['lat0', 'standard_parallel_1', d2r],
['lat1', 'standard_parallel_1', d2r],
['lat2', 'standard_parallel_2', d2r],
['azimuth', 'Azimuth'],
['alpha', 'azimuth', d2r],
['srsCode', 'name']
];
list.forEach(renamer);
if (!wkt.long0 && wkt.longc && (wkt.projName === 'Albers_Conic_Equal_Area' || wkt.projName === 'Lambert_Azimuthal_Equal_Area')) {
wkt.long0 = wkt.longc;
}
if (!wkt.lat_ts && wkt.lat1 && (wkt.projName === 'Stereographic_South_Pole' || wkt.projName === 'Polar Stereographic (variant B)')) {
wkt.lat0 = d2r(wkt.lat1 > 0 ? 90 : -90);
wkt.lat_ts = wkt.lat1;
}
}
var wkt = function(wkt) {
var lisp = parseString(wkt);
var type = lisp.shift();
var name = lisp.shift();
lisp.unshift(['name', name]);
lisp.unshift(['type', type]);
var obj = {};
sExpr(lisp, obj);
cleanWKT(obj);
return obj;
};
function defs(name) {
/*global console*/
var that = this;
if (arguments.length === 2) {
var def = arguments[1];
if (typeof def === 'string') {
if (def.charAt(0) === '+') {
defs[name] = parseProj(arguments[1]);
}
else {
defs[name] = wkt(arguments[1]);
}
} else {
defs[name] = def;
}
}
else if (arguments.length === 1) {
if (Array.isArray(name)) {
return name.map(function(v) {
if (Array.isArray(v)) {
defs.apply(that, v);
}
else {
defs(v);
}
});
}
else if (typeof name === 'string') {
if (name in defs) {
return defs[name];
}
}
else if ('EPSG' in name) {
defs['EPSG:' + name.EPSG] = name;
}
else if ('ESRI' in name) {
defs['ESRI:' + name.ESRI] = name;
}
else if ('IAU2000' in name) {
defs['IAU2000:' + name.IAU2000] = name;
}
else {
console.log(name);
}
return;
}
}
globals(defs);
function testObj(code){
return typeof code === 'string';
}
function testDef(code){
return code in defs;
}
var codeWords = ['PROJECTEDCRS', 'PROJCRS', 'GEOGCS','GEOCCS','PROJCS','LOCAL_CS', 'GEODCRS', 'GEODETICCRS', 'GEODETICDATUM', 'ENGCRS', 'ENGINEERINGCRS'];
function testWKT(code){
return codeWords.some(function (word) {
return code.indexOf(word) > -1;
});
}
var codes = ['3857', '900913', '3785', '102113'];
function checkMercator(item) {
var auth = match(item, 'authority');
if (!auth) {
return;
}
var code = match(auth, 'epsg');
return code && codes.indexOf(code) > -1;
}
function checkProjStr(item) {
var ext = match(item, 'extension');
if (!ext) {
return;
}
return match(ext, 'proj4');
}
function testProj(code){
return code[0] === '+';
}
function parse(code){
if (testObj(code)) {
//check to see if this is a WKT string
if (testDef(code)) {
return defs[code];
}
if (testWKT(code)) {
var out = wkt(code);
// test of spetial case, due to this being a very common and often malformed
if (checkMercator(out)) {
return defs['EPSG:3857'];
}
var maybeProjStr = checkProjStr(out);
if (maybeProjStr) {
return parseProj(maybeProjStr);
}
return out;
}
if (testProj(code)) {
return parseProj(code);
}
}else{
return code;
}
}
var extend = function(destination, source) {
destination = destination || {};
var value, property;
if (!source) {
return destination;
}
for (property in source) {
value = source[property];
if (value !== undefined) {
destination[property] = value;
}
}
return destination;
};
var msfnz = function(eccent, sinphi, cosphi) {
var con = eccent * sinphi;
return cosphi / (Math.sqrt(1 - con * con));
};
var sign = function(x) {
return x<0 ? -1 : 1;
};
var adjust_lon = function(x) {
return (Math.abs(x) <= SPI) ? x : (x - (sign(x) * TWO_PI));
};
var tsfnz = function(eccent, phi, sinphi) {
var con = eccent * sinphi;
var com = 0.5 * eccent;
con = Math.pow(((1 - con) / (1 + con)), com);
return (Math.tan(0.5 * (HALF_PI - phi)) / con);
};
var phi2z = function(eccent, ts) {
var eccnth = 0.5 * eccent;
var con, dphi;
var phi = HALF_PI - 2 * Math.atan(ts);
for (var i = 0; i <= 15; i++) {
con = eccent * Math.sin(phi);
dphi = HALF_PI - 2 * Math.atan(ts * (Math.pow(((1 - con) / (1 + con)), eccnth))) - phi;
phi += dphi;
if (Math.abs(dphi) <= 0.0000000001) {
return phi;
}
}
//console.log("phi2z has NoConvergence");
return -9999;
};
function init() {
var con = this.b / this.a;
this.es = 1 - con * con;
if(!('x0' in this)){
this.x0 = 0;
}
if(!('y0' in this)){
this.y0 = 0;
}
this.e = Math.sqrt(this.es);
if (this.lat_ts) {
if (this.sphere) {
this.k0 = Math.cos(this.lat_ts);
}
else {
this.k0 = msfnz(this.e, Math.sin(this.lat_ts), Math.cos(this.lat_ts));
}
}
else {
if (!this.k0) {
if (this.k) {
this.k0 = this.k;
}
else {
this.k0 = 1;
}
}
}
}
/* Mercator forward equations--mapping lat,long to x,y
--------------------------------------------------*/
function forward(p) {
var lon = p.x;
var lat = p.y;
// convert to radians
if (lat * R2D > 90 && lat * R2D < -90 && lon * R2D > 180 && lon * R2D < -180) {
return null;
}
var x, y;
if (Math.abs(Math.abs(lat) - HALF_PI) <= EPSLN) {
return null;
}
else {
if (this.sphere) {
x = this.x0 + this.a * this.k0 * adjust_lon(lon - this.long0);
y = this.y0 + this.a * this.k0 * Math.log(Math.tan(FORTPI + 0.5 * lat));
}
else {
var sinphi = Math.sin(lat);
var ts = tsfnz(this.e, lat, sinphi);
x = this.x0 + this.a * this.k0 * adjust_lon(lon - this.long0);
y = this.y0 - this.a * this.k0 * Math.log(ts);
}
p.x = x;
p.y = y;
return p;
}
}
/* Mercator inverse equations--mapping x,y to lat/long
--------------------------------------------------*/
function inverse(p) {
var x = p.x - this.x0;
var y = p.y - this.y0;
var lon, lat;
if (this.sphere) {
lat = HALF_PI - 2 * Math.atan(Math.exp(-y / (this.a * this.k0)));
}
else {
var ts = Math.exp(-y / (this.a * this.k0));
lat = phi2z(this.e, ts);
if (lat === -9999) {
return null;
}
}
lon = adjust_lon(this.long0 + x / (this.a * this.k0));
p.x = lon;
p.y = lat;
return p;
}
var names$1 = ["Mercator", "Popular Visualisation Pseudo Mercator", "Mercator_1SP", "Mercator_Auxiliary_Sphere", "merc"];
var merc = {
init: init,
forward: forward,
inverse: inverse,
names: names$1
};
function init$1() {
//no-op for longlat
}
function identity(pt) {
return pt;
}
var names$2 = ["longlat", "identity"];
var longlat = {
init: init$1,
forward: identity,
inverse: identity,
names: names$2
};
var projs = [merc, longlat];
var names = {};
var projStore = [];
function add(proj, i) {
var len = projStore.length;
if (!proj.names) {
console.log(i);
return true;
}
projStore[len] = proj;
proj.names.forEach(function(n) {
names[n.toLowerCase()] = len;
});
return this;
}
function get(name) {
if (!name) {
return false;
}
var n = name.toLowerCase();
if (typeof names[n] !== 'undefined' && projStore[names[n]]) {
return projStore[names[n]];
}
}
function start() {
projs.forEach(add);
}
var projections = {
start: start,
add: add,
get: get
};
var exports$2 = {};
exports$2.MERIT = {
a: 6378137.0,
rf: 298.257,
ellipseName: "MERIT 1983"
};
exports$2.SGS85 = {
a: 6378136.0,
rf: 298.257,
ellipseName: "Soviet Geodetic System 85"
};
exports$2.GRS80 = {
a: 6378137.0,
rf: 298.257222101,
ellipseName: "GRS 1980(IUGG, 1980)"
};
exports$2.IAU76 = {
a: 6378140.0,
rf: 298.257,
ellipseName: "IAU 1976"
};
exports$2.airy = {
a: 6377563.396,
b: 6356256.910,
ellipseName: "Airy 1830"
};
exports$2.APL4 = {
a: 6378137,
rf: 298.25,
ellipseName: "Appl. Physics. 1965"
};
exports$2.NWL9D = {
a: 6378145.0,
rf: 298.25,
ellipseName: "Naval Weapons Lab., 1965"
};
exports$2.mod_airy = {
a: 6377340.189,
b: 6356034.446,
ellipseName: "Modified Airy"
};
exports$2.andrae = {
a: 6377104.43,
rf: 300.0,
ellipseName: "Andrae 1876 (Den., Iclnd.)"
};
exports$2.aust_SA = {
a: 6378160.0,
rf: 298.25,
ellipseName: "Australian Natl & S. Amer. 1969"
};
exports$2.GRS67 = {
a: 6378160.0,
rf: 298.2471674270,
ellipseName: "GRS 67(IUGG 1967)"
};
exports$2.bessel = {
a: 6377397.155,
rf: 299.1528128,
ellipseName: "Bessel 1841"
};
exports$2.bess_nam = {
a: 6377483.865,
rf: 299.1528128,
ellipseName: "Bessel 1841 (Namibia)"
};
exports$2.clrk66 = {
a: 6378206.4,
b: 6356583.8,
ellipseName: "Clarke 1866"
};
exports$2.clrk80 = {
a: 6378249.145,
rf: 293.4663,
ellipseName: "Clarke 1880 mod."
};
exports$2.clrk80ign = {
a: 6378249.2,
b: 6356515,
rf: 293.4660213,
ellipseName: "Clarke 1880 (IGN)"
};
exports$2.clrk58 = {
a: 6378293.645208759,
rf: 294.2606763692654,
ellipseName: "Clarke 1858"
};
exports$2.CPM = {
a: 6375738.7,
rf: 334.29,
ellipseName: "Comm. des Poids et Mesures 1799"
};
exports$2.delmbr = {
a: 6376428.0,
rf: 311.5,
ellipseName: "Delambre 1810 (Belgium)"
};
exports$2.engelis = {
a: 6378136.05,
rf: 298.2566,
ellipseName: "Engelis 1985"
};
exports$2.evrst30 = {
a: 6377276.345,
rf: 300.8017,
ellipseName: "Everest 1830"
};
exports$2.evrst48 = {
a: 6377304.063,
rf: 300.8017,
ellipseName: "Everest 1948"
};
exports$2.evrst56 = {
a: 6377301.243,
rf: 300.8017,
ellipseName: "Everest 1956"
};
exports$2.evrst69 = {
a: 6377295.664,
rf: 300.8017,
ellipseName: "Everest 1969"
};
exports$2.evrstSS = {
a: 6377298.556,
rf: 300.8017,
ellipseName: "Everest (Sabah & Sarawak)"
};
exports$2.fschr60 = {
a: 6378166.0,
rf: 298.3,
ellipseName: "Fischer (Mercury Datum) 1960"
};
exports$2.fschr60m = {
a: 6378155.0,
rf: 298.3,
ellipseName: "Fischer 1960"
};
exports$2.fschr68 = {
a: 6378150.0,
rf: 298.3,
ellipseName: "Fischer 1968"
};
exports$2.helmert = {
a: 6378200.0,
rf: 298.3,
ellipseName: "Helmert 1906"
};
exports$2.hough = {
a: 6378270.0,
rf: 297.0,
ellipseName: "Hough"
};
exports$2.intl = {
a: 6378388.0,
rf: 297.0,
ellipseName: "International 1909 (Hayford)"
};
exports$2.kaula = {
a: 6378163.0,
rf: 298.24,
ellipseName: "Kaula 1961"
};
exports$2.lerch = {
a: 6378139.0,
rf: 298.257,
ellipseName: "Lerch 1979"
};
exports$2.mprts = {
a: 6397300.0,
rf: 191.0,
ellipseName: "Maupertius 1738"
};
exports$2.new_intl = {
a: 6378157.5,
b: 6356772.2,
ellipseName: "New International 1967"
};
exports$2.plessis = {
a: 6376523.0,
rf: 6355863.0,
ellipseName: "Plessis 1817 (France)"
};
exports$2.krass = {
a: 6378245.0,
rf: 298.3,
ellipseName: "Krassovsky, 1942"
};
exports$2.SEasia = {
a: 6378155.0,
b: 6356773.3205,
ellipseName: "Southeast Asia"
};
exports$2.walbeck = {
a: 6376896.0,
b: 6355834.8467,
ellipseName: "Walbeck"
};
exports$2.WGS60 = {
a: 6378165.0,
rf: 298.3,
ellipseName: "WGS 60"
};
exports$2.WGS66 = {
a: 6378145.0,
rf: 298.25,
ellipseName: "WGS 66"
};
exports$2.WGS7 = {
a: 6378135.0,
rf: 298.26,
ellipseName: "WGS 72"
};
var WGS84 = exports$2.WGS84 = {
a: 6378137.0,
rf: 298.257223563,
ellipseName: "WGS 84"
};
exports$2.sphere = {
a: 6370997.0,
b: 6370997.0,
ellipseName: "Normal Sphere (r=6370997)"
};
function eccentricity(a, b, rf, R_A) {
var a2 = a * a; // used in geocentric
var b2 = b * b; // used in geocentric
var es = (a2 - b2) / a2; // e ^ 2
var e = 0;
if (R_A) {
a *= 1 - es * (SIXTH + es * (RA4 + es * RA6));
a2 = a * a;
es = 0;
} else {
e = Math.sqrt(es); // eccentricity
}
var ep2 = (a2 - b2) / b2; // used in geocentric
return {
es: es,
e: e,
ep2: ep2
};
}
function sphere(a, b, rf, ellps, sphere) {
if (!a) { // do we have an ellipsoid?
var ellipse = match(exports$2, ellps);
if (!ellipse) {
ellipse = WGS84;
}
a = ellipse.a;
b = ellipse.b;
rf = ellipse.rf;
}
if (rf && !b) {
b = (1.0 - 1.0 / rf) * a;
}
if (rf === 0 || Math.abs(a - b) < EPSLN) {
sphere = true;
b = a;
}
return {
a: a,
b: b,
rf: rf,
sphere: sphere
};
}
var exports$3 = {};
exports$3.wgs84 = {
towgs84: "0,0,0",
ellipse: "WGS84",
datumName: "WGS84"
};
exports$3.ch1903 = {
towgs84: "674.374,15.056,405.346",
ellipse: "bessel",
datumName: "swiss"
};
exports$3.ggrs87 = {
towgs84: "-199.87,74.79,246.62",
ellipse: "GRS80",
datumName: "Greek_Geodetic_Reference_System_1987"
};
exports$3.nad83 = {
towgs84: "0,0,0",
ellipse: "GRS80",
datumName: "North_American_Datum_1983"
};
exports$3.nad27 = {
nadgrids: "@conus,@alaska,@ntv2_0.gsb,@ntv1_can.dat",
ellipse: "clrk66",
datumName: "North_American_Datum_1927"
};
exports$3.potsdam = {
towgs84: "598.1,73.7,418.2,0.202,0.045,-2.455,6.7",
ellipse: "bessel",
datumName: "Potsdam Rauenberg 1950 DHDN"
};
exports$3.carthage = {
towgs84: "-263.0,6.0,431.0",
ellipse: "clark80",
datumName: "Carthage 1934 Tunisia"
};
exports$3.hermannskogel = {
towgs84: "577.326,90.129,463.919,5.137,1.474,5.297,2.4232",
ellipse: "bessel",
datumName: "Hermannskogel"
};
exports$3.osni52 = {
towgs84: "482.530,-130.596,564.557,-1.042,-0.214,-0.631,8.15",
ellipse: "airy",
datumName: "Irish National"
};
exports$3.ire65 = {
towgs84: "482.530,-130.596,564.557,-1.042,-0.214,-0.631,8.15",
ellipse: "mod_airy",
datumName: "Ireland 1965"
};
exports$3.rassadiran = {
towgs84: "-133.63,-157.5,-158.62",
ellipse: "intl",
datumName: "Rassadiran"
};
exports$3.nzgd49 = {
towgs84: "59.47,-5.04,187.44,0.47,-0.1,1.024,-4.5993",
ellipse: "intl",
datumName: "New Zealand Geodetic Datum 1949"
};
exports$3.osgb36 = {
towgs84: "446.448,-125.157,542.060,0.1502,0.2470,0.8421,-20.4894",
ellipse: "airy",
datumName: "Airy 1830"
};
exports$3.s_jtsk = {
towgs84: "589,76,480",
ellipse: 'bessel',
datumName: 'S-JTSK (Ferro)'
};
exports$3.beduaram = {
towgs84: '-106,-87,188',
ellipse: 'clrk80',
datumName: 'Beduaram'
};
exports$3.gunung_segara = {
towgs84: '-403,684,41',
ellipse: 'bessel',
datumName: 'Gunung Segara Jakarta'
};
exports$3.rnb72 = {
towgs84: "106.869,-52.2978,103.724,-0.33657,0.456955,-1.84218,1",
ellipse: "intl",
datumName: "Reseau National Belge 1972"
};
function datum(datumCode, datum_params, a, b, es, ep2, nadgrids) {
var out = {};
if (datumCode === undefined || datumCode === 'none') {
out.datum_type = PJD_NODATUM;
} else {
out.datum_type = PJD_WGS84;
}
if (datum_params) {
out.datum_params = datum_params.map(parseFloat);
if (out.datum_params[0] !== 0 || out.datum_params[1] !== 0 || out.datum_params[2] !== 0) {
out.datum_type = PJD_3PARAM;
}
if (out.datum_params.length > 3) {
if (out.datum_params[3] !== 0 || out.datum_params[4] !== 0 || out.datum_params[5] !== 0 || out.datum_params[6] !== 0) {
out.datum_type = PJD_7PARAM;
out.datum_params[3] *= SEC_TO_RAD;
out.datum_params[4] *= SEC_TO_RAD;
out.datum_params[5] *= SEC_TO_RAD;
out.datum_params[6] = (out.datum_params[6] / 1000000.0) + 1.0;
}
}
}
if (nadgrids) {
out.datum_type = PJD_GRIDSHIFT;
out.grids = nadgrids;
}
out.a = a; //datum object also uses these values
out.b = b;
out.es = es;
out.ep2 = ep2;
return out;
}
/**
* Resources for details of NTv2 file formats:
* - https://web.archive.org/web/20140127204822if_/http://www.mgs.gov.on.ca:80/stdprodconsume/groups/content/@mgs/@iandit/documents/resourcelist/stel02_047447.pdf
* - http://mimaka.com/help/gs/html/004_NTV2%20Data%20Format.htm
*/
var loadedNadgrids = {};
/**
* Load a binary NTv2 file (.gsb) to a key that can be used in a proj string like +nadgrids=<key>. Pass the NTv2 file
* as an ArrayBuffer.
*/
function nadgrid(key, data) {
var view = new DataView(data);
var isLittleEndian = detectLittleEndian(view);
var header = readHeader(view, isLittleEndian);
if (header.nSubgrids > 1) {
console.log('Only single NTv2 subgrids are currently supported, subsequent sub grids are ignored');
}
var subgrids = readSubgrids(view, header, isLittleEndian);
var nadgrid = {header: header, subgrids: subgrids};
loadedNadgrids[key] = nadgrid;
return nadgrid;
}
/**
* Given a proj4 value for nadgrids, return an array of loaded grids
*/
function getNadgrids(nadgrids) {
// Format details: http://proj.maptools.org/gen_parms.html
if (nadgrids === undefined) { return null; }
var grids = nadgrids.split(',');
return grids.map(parseNadgridString);
}
function parseNadgridString(value) {
if (value.length === 0) {
return null;
}
var optional = value[0] === '@';
if (optional) {
value = value.slice(1);
}
if (value === 'null') {
return {name: 'null', mandatory: !optional, grid: null, isNull: true};
}
return {
name: value,
mandatory: !optional,
grid: loadedNadgrids[value] || null,
isNull: false
};
}
function secondsToRadians(seconds) {
return (seconds / 3600) * Math.PI / 180;
}
function detectLittleEndian(view) {
var nFields = view.getInt32(8, false);
if (nFields === 11) {
return false;
}
nFields = view.getInt32(8, true);
if (nFields !== 11) {
console.warn('Failed to detect nadgrid endian-ness, defaulting to little-endian');
}
return true;
}
function readHeader(view, isLittleEndian) {
return {
nFields: view.getInt32(8, isLittleEndian),
nSubgridFields: view.getInt32(24, isLittleEndian),
nSubgrids: view.getInt32(40, isLittleEndian),
shiftType: decodeString(view, 56, 56 + 8).trim(),
fromSemiMajorAxis: view.getFloat64(120, isLittleEndian),
fromSemiMinorAxis: view.getFloat64(136, isLittleEndian),
toSemiMajorAxis: view.getFloat64(152, isLittleEndian),
toSemiMinorAxis: view.getFloat64(168, isLittleEndian),
};
}
function decodeString(view, start, end) {
return String.fromCharCode.apply(null, new Uint8Array(view.buffer.slice(start, end)));
}
function readSubgrids(view, header, isLittleEndian) {
var gridOffset = 176;
var grids = [];
for (var i = 0; i < header.nSubgrids; i++) {
var subHeader = readGridHeader(view, gridOffset, isLittleEndian);
var nodes = readGridNodes(view, gridOffset, subHeader, isLittleEndian);
var lngColumnCount = Math.round(
1 + (subHeader.upperLongitude - subHeader.lowerLongitude) / subHeader.longitudeInterval);
var latColumnCount = Math.round(
1 + (subHeader.upperLatitude - subHeader.lowerLatitude) / subHeader.latitudeInterval);
// Proj4 operates on radians whereas the coordinates are in seconds in the grid
grids.push({
ll: [secondsToRadians(subHeader.lowerLongitude), secondsToRadians(subHeader.lowerLatitude)],
del: [secondsToRadians(subHeader.longitudeInterval), secondsToRadians(subHeader.latitudeInterval)],
lim: [lngColumnCount, latColumnCount],
count: subHeader.gridNodeCount,
cvs: mapNodes(nodes)
});
}
return grids;
}
function mapNodes(nodes) {
return nodes.map(function (r) {return [secondsToRadians(r.longitudeShift), secondsToRadians(r.latitudeShift)];});
}
function readGridHeader(view, offset, isLittleEndian) {
return {
name: decodeString(view, offset + 8, offset + 16).trim(),
parent: decodeString(view, offset + 24, offset + 24 + 8).trim(),
lowerLatitude: view.getFloat64(offset + 72, isLittleEndian),
upperLatitude: view.getFloat64(offset + 88, isLittleEndian),
lowerLongitude: view.getFloat64(offset + 104, isLittleEndian),
upperLongitude: view.getFloat64(offset + 120, isLittleEndian),
latitudeInterval: view.getFloat64(offset + 136, isLittleEndian),
longitudeInterval: view.getFloat64(offset + 152, isLittleEndian),
gridNodeCount: view.getInt32(offset + 168, isLittleEndian)
};
}
function readGridNodes(view, offset, gridHeader, isLittleEndian) {
var nodesOffset = offset + 176;
var gridRecordLength = 16;
var gridShiftRecords = [];
for (var i = 0; i < gridHeader.gridNodeCount; i++) {
var record = {
latitudeShift: view.getFloat32(nodesOffset + i * gridRecordLength, isLittleEndian),
longitudeShift: view.getFloat32(nodesOffset + i * gridRecordLength + 4, isLittleEndian),
latitudeAccuracy: view.getFloat32(nodesOffset + i * gridRecordLength + 8, isLittleEndian),
longitudeAccuracy: view.getFloat32(nodesOffset + i * gridRecordLength + 12, isLittleEndian),
};
gridShiftRecords.push(record);
}
return gridShiftRecords;
}
function Projection(srsCode,callback) {
if (!(this instanceof Projection)) {
return new Projection(srsCode);
}
callback = callback || function(error){
if(error){
throw error;
}
};
var json = parse(srsCode);
if(typeof json !== 'object'){
callback(srsCode);
return;
}
var ourProj = Projection.projections.get(json.projName);
if(!ourProj){
callback(srsCode);
return;
}
if (json.datumCode && json.datumCode !== 'none') {
var datumDef = match(exports$3, json.datumCode);
if (datumDef) {
json.datum_params = json.datum_params || (datumDef.towgs84 ? datumDef.towgs84.split(',') : null);
json.ellps = datumDef.ellipse;
json.datumName = datumDef.datumName ? datumDef.datumName : json.datumCode;
}
}
json.k0 = json.k0 || 1.0;
json.axis = json.axis || 'enu';
json.ellps = json.ellps || 'wgs84';
json.lat1 = json.lat1 || json.lat0; // Lambert_Conformal_Conic_1SP, for example, needs this
var sphere_ = sphere(json.a, json.b, json.rf, json.ellps, json.sphere);
var ecc = eccentricity(sphere_.a, sphere_.b, sphere_.rf, json.R_A);
var nadgrids = getNadgrids(json.nadgrids);
var datumObj = json.datum || datum(json.datumCode, json.datum_params, sphere_.a, sphere_.b, ecc.es, ecc.ep2,
nadgrids);
extend(this, json); // transfer everything over from the projection because we don't know what we'll need
extend(this, ourProj); // transfer all the methods from the projection
// copy the 4 things over we calculated in deriveConstants.sphere
this.a = sphere_.a;
this.b = sphere_.b;
this.rf = sphere_.rf;
this.sphere = sphere_.sphere;
// copy the 3 things we calculated in deriveConstants.eccentricity
this.es = ecc.es;
this.e = ecc.e;
this.ep2 = ecc.ep2;
// add in the datum object
this.datum = datumObj;
// init the projection
this.init();
// legecy callback from back in the day when it went to spatialreference.org
callback(null, this);
}
Projection.projections = projections;
Projection.projections.start();
'use strict';
function compareDatums(source, dest) {
if (source.datum_type !== dest.datum_type) {
return false; // false, datums are not equal
} else if (source.a !== dest.a || Math.abs(source.es - dest.es) > 0.000000000050) {
// the tolerance for es is to ensure that GRS80 and WGS84
// are considered identical
return false;
} else if (source.datum_type === PJD_3PARAM) {
return (source.datum_params[0] === dest.datum_params[0] && source.datum_params[1] === dest.datum_params[1] && source.datum_params[2] === dest.datum_params[2]);
} else if (source.datum_type === PJD_7PARAM) {
return (source.datum_params[0] === dest.datum_params[0] && source.datum_params[1] === dest.datum_params[1] && source.datum_params[2] === dest.datum_params[2] && source.datum_params[3] === dest.datum_params[3] && source.datum_params[4] === dest.datum_params[4] && source.datum_params[5] === dest.datum_params[5] && source.datum_params[6] === dest.datum_params[6]);
} else {
return true; // datums are equal
}
} // cs_compare_datums()
/*
* The function Convert_Geodetic_To_Geocentric converts geodetic coordinates
* (latitude, longitude, and height) to geocentric coordinates (X, Y, Z),
* according to the current ellipsoid parameters.
*
* Latitude : Geodetic latitude in radians (input)
* Longitude : Geodetic longitude in radians (input)
* Height : Geodetic height, in meters (input)
* X : Calculated Geocentric X coordinate, in meters (output)
* Y : Calculated Geocentric Y coordinate, in meters (output)
* Z : Calculated Geocentric Z coordinate, in meters (output)
*
*/
function geodeticToGeocentric(p, es, a) {
var Longitude = p.x;
var Latitude = p.y;
var Height = p.z ? p.z : 0; //Z value not always supplied
var Rn; /* Earth radius at location */
var Sin_Lat; /* Math.sin(Latitude) */
var Sin2_Lat; /* Square of Math.sin(Latitude) */
var Cos_Lat; /* Math.cos(Latitude) */
/*
** Don't blow up if Latitude is just a little out of the value
** range as it may just be a rounding issue. Also removed longitude
** test, it should be wrapped by Math.cos() and Math.sin(). NFW for PROJ.4, Sep/2001.
*/
if (Latitude < -HALF_PI && Latitude > -1.001 * HALF_PI) {
Latitude = -HALF_PI;
} else if (Latitude > HALF_PI && Latitude < 1.001 * HALF_PI) {
Latitude = HALF_PI;
} else if (Latitude < -HALF_PI) {
/* Latitude out of range */
//..reportError('geocent:lat out of range:' + Latitude);
return { x: -Infinity, y: -Infinity, z: p.z };
} else if (Latitude > HALF_PI) {
/* Latitude out of range */
return { x: Infinity, y: Infinity, z: p.z };
}
if (Longitude > Math.PI) {
Longitude -= (2 * Math.PI);
}
Sin_Lat = Math.sin(Latitude);
Cos_Lat = Math.cos(Latitude);
Sin2_Lat = Sin_Lat * Sin_Lat;
Rn = a / (Math.sqrt(1.0e0 - es * Sin2_Lat));
return {
x: (Rn + Height) * Cos_Lat * Math.cos(Longitude),
y: (Rn + Height) * Cos_Lat * Math.sin(Longitude),
z: ((Rn * (1 - es)) + Height) * Sin_Lat
};
} // cs_geodetic_to_geocentric()
function geocentricToGeodetic(p, es, a, b) {
/* local defintions and variables */
/* end-criterium of loop, accuracy of sin(Latitude) */
var genau = 1e-12;
var genau2 = (genau * genau);
var maxiter = 30;
var P; /* distance between semi-minor axis and location */
var RR; /* distance between center and location */
var CT; /* sin of geocentric latitude */
var ST; /* cos of geocentric latitude */
var RX;
var RK;
var RN; /* Earth radius at location */
var CPHI0; /* cos of start or old geodetic latitude in iterations */
var SPHI0; /* sin of start or old geodetic latitude in iterations */
var CPHI; /* cos of searched geodetic latitude */
var SPHI; /* sin of searched geodetic latitude */
var SDPHI; /* end-criterium: addition-theorem of sin(Latitude(iter)-Latitude(iter-1)) */
var iter; /* # of continous iteration, max. 30 is always enough (s.a.) */
var X = p.x;
var Y = p.y;
var Z = p.z ? p.z : 0.0; //Z value not always supplied
var Longitude;
var Latitude;
var Height;
P = Math.sqrt(X * X + Y * Y);
RR = Math.sqrt(X * X + Y * Y + Z * Z);
/* special cases for latitude and longitude */
if (P / a < genau) {
/* special case, if P=0. (X=0., Y=0.) */
Longitude = 0.0;
/* if (X,Y,Z)=(0.,0.,0.) then Height becomes semi-minor axis
* of ellipsoid (=center of mass), Latitude becomes PI/2 */
if (RR / a < genau) {
Latitude = HALF_PI;
Height = -b;
return {
x: p.x,
y: p.y,
z: p.z
};
}
} else {
/* ellipsoidal (geodetic) longitude
* interval: -PI < Longitude <= +PI */
Longitude = Math.atan2(Y, X);
}
/* --------------------------------------------------------------
* Following iterative algorithm was developped by
* "Institut for Erdmessung", University of Hannover, July 1988.
* Internet: www.ife.uni-hannover.de
* Iterative computation of CPHI,SPHI and Height.
* Iteration of CPHI and SPHI to 10**-12 radian resp.
* 2*10**-7 arcsec.
* --------------------------------------------------------------
*/
CT = Z / RR;
ST = P / RR;
RX = 1.0 / Math.sqrt(1.0 - es * (2.0 - es) * ST * ST);
CPHI0 = ST * (1.0 - es) * RX;
SPHI0 = CT * RX;
iter = 0;
/* loop to find sin(Latitude) resp. Latitude
* until |sin(Latitude(iter)-Latitude(iter-1))| < genau */
do {
iter++;
RN = a / Math.sqrt(1.0 - es * SPHI0 * SPHI0);
/* ellipsoidal (geodetic) height */
Height = P * CPHI0 + Z * SPHI0 - RN * (1.0 - es * SPHI0 * SPHI0);
RK = es * RN / (RN + Height);
RX = 1.0 / Math.sqrt(1.0 - RK * (2.0 - RK) * ST * ST);
CPHI = ST * (1.0 - RK) * RX;
SPHI = CT * RX;
SDPHI = SPHI * CPHI0 - CPHI * SPHI0;
CPHI0 = CPHI;
SPHI0 = SPHI;
}
while (SDPHI * SDPHI > genau2 && iter < maxiter);
/* ellipsoidal (geodetic) latitude */
Latitude = Math.atan(SPHI / Math.abs(CPHI));
return {
x: Longitude,
y: Latitude,
z: Height
};
} // cs_geocentric_to_geodetic()
/****************************************************************/
// pj_geocentic_to_wgs84( p )
// p = point to transform in geocentric coordinates (x,y,z)
/** point object, nothing fancy, just allows values to be
passed back and forth by reference rather than by value.
Other point classes may be used as long as they have
x and y properties, which will get modified in the transform method.
*/
function geocentricToWgs84(p, datum_type, datum_params) {
if (datum_type === PJD_3PARAM) {
// if( x[io] === HUGE_VAL )
// continue;
return {
x: p.x + datum_params[0],
y: p.y + datum_params[1],
z: p.z + datum_params[2],
};
} else if (datum_type === PJD_7PARAM) {
var Dx_BF = datum_params[0];
var Dy_BF = datum_params[1];
var Dz_BF = datum_params[2];
var Rx_BF = datum_params[3];
var Ry_BF = datum_params[4];
var Rz_BF = datum_params[5];
var M_BF = datum_params[6];
// if( x[io] === HUGE_VAL )
// continue;
return {
x: M_BF * (p.x - Rz_BF * p.y + Ry_BF * p.z) + Dx_BF,
y: M_BF * (Rz_BF * p.x + p.y - Rx_BF * p.z) + Dy_BF,
z: M_BF * (-Ry_BF * p.x + Rx_BF * p.y + p.z) + Dz_BF
};
}
} // cs_geocentric_to_wgs84
/****************************************************************/
// pj_geocentic_from_wgs84()
// coordinate system definition,
// point to transform in geocentric coordinates (x,y,z)
function geocentricFromWgs84(p, datum_type, datum_params) {
if (datum_type === PJD_3PARAM) {
//if( x[io] === HUGE_VAL )
// continue;
return {
x: p.x - datum_params[0],
y: p.y - datum_params[1],
z: p.z - datum_params[2],
};
} else if (datum_type === PJD_7PARAM) {
var Dx_BF = datum_params[0];
var Dy_BF = datum_params[1];
var Dz_BF = datum_params[2];
var Rx_BF = datum_params[3];
var Ry_BF = datum_params[4];
var Rz_BF = datum_params[5];
var M_BF = datum_params[6];
var x_tmp = (p.x - Dx_BF) / M_BF;
var y_tmp = (p.y - Dy_BF) / M_BF;
var z_tmp = (p.z - Dz_BF) / M_BF;
//if( x[io] === HUGE_VAL )
// continue;
return {
x: x_tmp + Rz_BF * y_tmp - Ry_BF * z_tmp,
y: -Rz_BF * x_tmp + y_tmp + Rx_BF * z_tmp,
z: Ry_BF * x_tmp - Rx_BF * y_tmp + z_tmp
};
} //cs_geocentric_from_wgs84()
}
function checkParams(type) {
return (type === PJD_3PARAM || type === PJD_7PARAM);
}
var datum_transform = function(source, dest, point) {
// Short cut if the datums are identical.
if (compareDatums(source, dest)) {
return point; // in this case, zero is sucess,
// whereas cs_compare_datums returns 1 to indicate TRUE
// confusing, should fix this
}
// Explicitly skip datum transform by setting 'datum=none' as parameter for either source or dest
if (source.datum_type === PJD_NODATUM || dest.datum_type === PJD_NODATUM) {
return point;
}
// If this datum requires grid shifts, then apply it to geodetic coordinates.
var source_a = source.a;
var source_es = source.es;
if (source.datum_type === PJD_GRIDSHIFT) {
var gridShiftCode = applyGridShift(source, false, point);
if (gridShiftCode !== 0) {
return undefined;
}
source_a = SRS_WGS84_SEMIMAJOR;
source_es = SRS_WGS84_ESQUARED;
}
var dest_a = dest.a;
var dest_b = dest.b;
var dest_es = dest.es;
if (dest.datum_type === PJD_GRIDSHIFT) {
dest_a = SRS_WGS84_SEMIMAJOR;
dest_b = SRS_WGS84_SEMIMINOR;
dest_es = SRS_WGS84_ESQUARED;
}
// Do we need to go through geocentric coordinates?
if (source_es === dest_es && source_a === dest_a && !checkParams(source.datum_type) && !checkParams(dest.datum_type)) {
return point;
}
// Convert to geocentric coordinates.
point = geodeticToGeocentric(point, source_es, source_a);
// Convert between datums
if (checkParams(source.datum_type)) {
point = geocentricToWgs84(point, source.datum_type, source.datum_params);
}
if (checkParams(dest.datum_type)) {
point = geocentricFromWgs84(point, dest.datum_type, dest.datum_params);
}
point = geocentricToGeodetic(point, dest_es, dest_a, dest_b);
if (dest.datum_type === PJD_GRIDSHIFT) {
var destGridShiftResult = applyGridShift(dest, true, point);
if (destGridShiftResult !== 0) {
return undefined;
}
}
return point;
};
function applyGridShift(source, inverse, point) {
if (source.grids === null || source.grids.length === 0) {
console.log('Grid shift grids not found');
return -1;
}
var input = {x: -point.x, y: point.y};
var output = {x: Number.NaN, y: Number.NaN};
var attemptedGrids = [];
for (var i = 0; i < source.grids.length; i++) {
var grid = source.grids[i];
attemptedGrids.push(grid.name);
if (grid.isNull) {
output = input;
break;
}
if (grid.grid === null) {
if (grid.mandatory) {
console.log("Unable to find mandatory grid '" + grid.name + "'");
return -1;
}
continue;
}
var subgrid = grid.grid.subgrids[0];
// skip tables that don't match our point at all
var epsilon = (Math.abs(subgrid.del[1]) + Math.abs(subgrid.del[0])) / 10000.0;
var minX = subgrid.ll[0] - epsilon;
var minY = subgrid.ll[1] - epsilon;
var maxX = subgrid.ll[0] + (subgrid.lim[0] - 1) * subgrid.del[0] + epsilon;
var maxY = subgrid.ll[1] + (subgrid.lim[1] - 1) * subgrid.del[1] + epsilon;
if (minY > input.y || minX > input.x || maxY < input.y || maxX < input.x ) {
continue;
}
output = applySubgridShift(input, inverse, subgrid);
if (!isNaN(output.x)) {
break;
}
}
if (isNaN(output.x)) {
console.log("Failed to find a grid shift table for location '"+
-input.x * R2D + " " + input.y * R2D + " tried: '" + attemptedGrids + "'");
return -1;
}
point.x = -output.x;
point.y = output.y;
return 0;
}
function applySubgridShift(pin, inverse, ct) {
var val = {x: Number.NaN, y: Number.NaN};
if (isNaN(pin.x)) { return val; }
var tb = {x: pin.x, y: pin.y};
tb.x -= ct.ll[0];
tb.y -= ct.ll[1];
tb.x = adjust_lon(tb.x - Math.PI) + Math.PI;
var t = nadInterpolate(tb, ct);
if (inverse) {
if (isNaN(t.x)) {
return val;
}
t.x = tb.x - t.x;
t.y = tb.y - t.y;
var i = 9, tol = 1e-12;
var dif, del;
do {
del = nadInterpolate(t, ct);
if (isNaN(del.x)) {
console.log("Inverse grid shift iteration failed, presumably at grid edge. Using first approximation.");
break;
}
dif = {x: tb.x - (del.x + t.x), y: tb.y - (del.y + t.y)};
t.x += dif.x;
t.y += dif.y;
} while (i-- && Math.abs(dif.x) > tol && Math.abs(dif.y) > tol);
if (i < 0) {
console.log("Inverse grid shift iterator failed to converge.");
return val;
}
val.x = adjust_lon(t.x + ct.ll[0]);
val.y = t.y + ct.ll[1];
} else {
if (!isNaN(t.x)) {
val.x = pin.x + t.x;
val.y = pin.y + t.y;
}
}
return val;
}
function nadInterpolate(pin, ct) {
var t = {x: pin.x / ct.del[0], y: pin.y / ct.del[1]};
var indx = {x: Math.floor(t.x), y: Math.floor(t.y)};
var frct = {x: t.x - 1.0 * indx.x, y: t.y - 1.0 * indx.y};
var val= {x: Number.NaN, y: Number.NaN};
var inx;
if (indx.x < 0 || indx.x >= ct.lim[0]) {
return val;
}
if (indx.y < 0 || indx.y >= ct.lim[1]) {
return val;
}
inx = (indx.y * ct.lim[0]) + indx.x;
var f00 = {x: ct.cvs[inx][0], y: ct.cvs[inx][1]};
inx++;
var f10= {x: ct.cvs[inx][0], y: ct.cvs[inx][1]};
inx += ct.lim[0];
var f11 = {x: ct.cvs[inx][0], y: ct.cvs[inx][1]};
inx--;
var f01 = {x: ct.cvs[inx][0], y: ct.cvs[inx][1]};
var m11 = frct.x * frct.y, m10 = frct.x * (1.0 - frct.y),
m00 = (1.0 - frct.x) * (1.0 - frct.y), m01 = (1.0 - frct.x) * frct.y;
val.x = (m00 * f00.x + m10 * f10.x + m01 * f01.x + m11 * f11.x);
val.y = (m00 * f00.y + m10 * f10.y + m01 * f01.y + m11 * f11.y);
return val;
}
var adjust_axis = function(crs, denorm, point) {
var xin = point.x,
yin = point.y,
zin = point.z || 0.0;
var v, t, i;
var out = {};
for (i = 0; i < 3; i++) {
if (denorm && i === 2 && point.z === undefined) {
continue;
}
if (i === 0) {
v = xin;
if ("ew".indexOf(crs.axis[i]) !== -1) {
t = 'x';
} else {
t = 'y';
}
}
else if (i === 1) {
v = yin;
if ("ns".indexOf(crs.axis[i]) !== -1) {
t = 'y';
} else {
t = 'x';
}
}
else {
v = zin;
t = 'z';
}
switch (crs.axis[i]) {
case 'e':
out[t] = v;
break;
case 'w':
out[t] = -v;
break;
case 'n':
out[t] = v;
break;
case 's':
out[t] = -v;
break;
case 'u':
if (point[t] !== undefined) {
out.z = v;
}
break;
case 'd':
if (point[t] !== undefined) {
out.z = -v;
}
break;
default:
//console.log("ERROR: unknow axis ("+crs.axis[i]+") - check definition of "+crs.projName);
return null;
}
}
return out;
};
var toPoint = function (array){
var out = {
x: array[0],
y: array[1]
};
if (array.length>2) {
out.z = array[2];
}
if (array.length>3) {
out.m = array[3];
}
return out;
};
var checkSanity = function (point) {
checkCoord(point.x);
checkCoord(point.y);
};
function checkCoord(num) {
if (typeof Number.isFinite === 'function') {
if (Number.isFinite(num)) {
return;
}
throw new TypeError('coordinates must be finite numbers');
}
if (typeof num !== 'number' || num !== num || !isFinite(num)) {
throw new TypeError('coordinates must be finite numbers');
}
}
function checkNotWGS(source, dest) {
return (
(source.datum.datum_type === PJD_3PARAM || source.datum.datum_type === PJD_7PARAM || source.datum.datum_type === PJD_GRIDSHIFT) && dest.datumCode !== 'WGS84') ||
((dest.datum.datum_type === PJD_3PARAM || dest.datum.datum_type === PJD_7PARAM || dest.datum.datum_type === PJD_GRIDSHIFT) && source.datumCode !== 'WGS84');
}
function transform(source, dest, point, enforceAxis) {
var wgs84;
if (Array.isArray(point)) {
point = toPoint(point);
} else {
// Clone the point object so inputs don't get modified
point = {
x: point.x,
y: point.y,
z: point.z,
m: point.m
};
}
var hasZ = point.z !== undefined;
checkSanity(point);
// Workaround for datum shifts towgs84, if either source or destination projection is not wgs84
if (source.datum && dest.datum && checkNotWGS(source, dest)) {
wgs84 = new Projection('WGS84');
point = transform(source, wgs84, point, enforceAxis);
source = wgs84;
}
// DGR, 2010/11/12
if (enforceAxis && source.axis !== 'enu') {
point = adjust_axis(source, false, point);
}
// Transform source points to long/lat, if they aren't already.
if (source.projName === 'longlat') {
point = {
x: point.x * D2R,
y: point.y * D2R,
z: point.z || 0
};
} else {
if (source.to_meter) {
point = {
x: point.x * source.to_meter,
y: point.y * source.to_meter,
z: point.z || 0
};
}
point = source.inverse(point); // Convert Cartesian to longlat
if (!point) {
return;
}
}
// Adjust for the prime meridian if necessary
if (source.from_greenwich) {
point.x += source.from_greenwich;
}
// Convert datums if needed, and if possible.
point = datum_transform(source.datum, dest.datum, point);
if (!point) {
return;
}
// Adjust for the prime meridian if necessary
if (dest.from_greenwich) {
point = {
x: point.x - dest.from_greenwich,
y: point.y,
z: point.z || 0
};
}
if (dest.projName === 'longlat') {
// convert radians to decimal degrees
point = {
x: point.x * R2D,
y: point.y * R2D,
z: point.z || 0
};
} else { // else project
point = dest.forward(point);
if (dest.to_meter) {
point = {
x: point.x / dest.to_meter,
y: point.y / dest.to_meter,
z: point.z || 0
};
}
}
// DGR, 2010/11/12
if (enforceAxis && dest.axis !== 'enu') {
return adjust_axis(dest, true, point);
}
if (!hasZ) {
delete point.z;
}
return point;
}
var wgs84 = Projection('WGS84');
function transformer(from, to, coords, enforceAxis) {
var transformedArray, out, keys;
if (Array.isArray(coords)) {
transformedArray = transform(from, to, coords, enforceAxis) || {x: NaN, y: NaN};
if (coords.length > 2) {
if ((typeof from.name !== 'undefined' && from.name === 'geocent') || (typeof to.name !== 'undefined' && to.name === 'geocent')) {
if (typeof transformedArray.z === 'number') {
return [transformedArray.x, transformedArray.y, transformedArray.z].concat(coords.splice(3));
} else {
return [transformedArray.x, transformedArray.y, coords[2]].concat(coords.splice(3));
}
} else {
return [transformedArray.x, transformedArray.y].concat(coords.splice(2));
}
} else {
return [transformedArray.x, transformedArray.y];
}
} else {
out = transform(from, to, coords, enforceAxis);
keys = Object.keys(coords);
if (keys.length === 2) {
return out;
}
keys.forEach(function (key) {
if ((typeof from.name !== 'undefined' && from.name === 'geocent') || (typeof to.name !== 'undefined' && to.name === 'geocent')) {
if (key === 'x' || key === 'y' || key === 'z') {
return;
}
} else {
if (key === 'x' || key === 'y') {
return;
}
}
out[key] = coords[key];
});
return out;
}
}
function checkProj(item) {
if (item instanceof Projection) {
return item;
}
if (item.oProj) {
return item.oProj;
}
return Projection(item);
}
function proj4$1(fromProj, toProj, coord) {
fromProj = checkProj(fromProj);
var single = false;
var obj;
if (typeof toProj === 'undefined') {
toProj = fromProj;
fromProj = wgs84;
single = true;
} else if (typeof toProj.x !== 'undefined' || Array.isArray(toProj)) {
coord = toProj;
toProj = fromProj;
fromProj = wgs84;
single = true;
}
toProj = checkProj(toProj);
if (coord) {
return transformer(fromProj, toProj, coord);
} else {
obj = {
forward: function (coords, enforceAxis) {
return transformer(fromProj, toProj, coords, enforceAxis);
},
inverse: function (coords, enforceAxis) {
return transformer(toProj, fromProj, coords, enforceAxis);
}
};
if (single) {
obj.oProj = toProj;
}
return obj;
}
}
/**
* UTM zones are grouped, and assigned to one of a group of 6
* sets.
*
* {int} @private
*/
var NUM_100K_SETS = 6;
/**
* The column letters (for easting) of the lower left value, per
* set.
*
* {string} @private
*/
var SET_ORIGIN_COLUMN_LETTERS = 'AJSAJS';
/**
* The row letters (for northing) of the lower left value, per
* set.
*
* {string} @private
*/
var SET_ORIGIN_ROW_LETTERS = 'AFAFAF';
var A = 65; // A
var I = 73; // I
var O = 79; // O
var V = 86; // V
var Z = 90; // Z
var mgrs = {
forward: forward$1,
inverse: inverse$1,
toPoint: toPoint$1
};
/**
* Conversion of lat/lon to MGRS.
*
* @param {object} ll Object literal with lat and lon properties on a
* WGS84 ellipsoid.
* @param {int} accuracy Accuracy in digits (5 for 1 m, 4 for 10 m, 3 for
* 100 m, 2 for 1000 m or 1 for 10000 m). Optional, default is 5.
* @return {string} the MGRS string for the given location and accuracy.
*/
function forward$1(ll, accuracy) {
accuracy = accuracy || 5; // default accuracy 1m
return encode(LLtoUTM({
lat: ll[1],
lon: ll[0]
}), accuracy);
}
/**
* Conversion of MGRS to lat/lon.
*
* @param {string} mgrs MGRS string.
* @return {array} An array with left (longitude), bottom (latitude), right
* (longitude) and top (latitude) values in WGS84, representing the
* bounding box for the provided MGRS reference.
*/
function inverse$1(mgrs) {
var bbox = UTMtoLL(decode(mgrs.toUpperCase()));
if (bbox.lat && bbox.lon) {
return [bbox.lon, bbox.lat, bbox.lon, bbox.lat];
}
return [bbox.left, bbox.bottom, bbox.right, bbox.top];
}
function toPoint$1(mgrs) {
var bbox = UTMtoLL(decode(mgrs.toUpperCase()));
if (bbox.lat && bbox.lon) {
return [bbox.lon, bbox.lat];
}
return [(bbox.left + bbox.right) / 2, (bbox.top + bbox.bottom) / 2];
}
/**
* Conversion from degrees to radians.
*
* @private
* @param {number} deg the angle in degrees.
* @return {number} the angle in radians.
*/
function degToRad(deg) {
return (deg * (Math.PI / 180.0));
}
/**
* Conversion from radians to degrees.
*
* @private
* @param {number} rad the angle in radians.
* @return {number} the angle in degrees.
*/
function radToDeg(rad) {
return (180.0 * (rad / Math.PI));
}
/**
* Converts a set of Longitude and Latitude co-ordinates to UTM
* using the WGS84 ellipsoid.
*
* @private
* @param {object} ll Object literal with lat and lon properties
* representing the WGS84 coordinate to be converted.
* @return {object} Object literal containing the UTM value with easting,
* northing, zoneNumber and zoneLetter properties, and an optional
* accuracy property in digits. Returns null if the conversion failed.
*/
function LLtoUTM(ll) {
var Lat = ll.lat;
var Long = ll.lon;
var a = 6378137.0; //ellip.radius;
var eccSquared = 0.00669438; //ellip.eccsq;
var k0 = 0.9996;
var LongOrigin;
var eccPrimeSquared;
var N, T, C, A, M;
var LatRad = degToRad(Lat);
var LongRad = degToRad(Long);
var LongOriginRad;
var ZoneNumber;
// (int)
ZoneNumber = Math.floor((Long + 180) / 6) + 1;
//Make sure the longitude 180.00 is in Zone 60
if (Long === 180) {
ZoneNumber = 60;
}
// Special zone for Norway
if (Lat >= 56.0 && Lat < 64.0 && Long >= 3.0 && Long < 12.0) {
ZoneNumber = 32;
}
// Special zones for Svalbard
if (Lat >= 72.0 && Lat < 84.0) {
if (Long >= 0.0 && Long < 9.0) {
ZoneNumber = 31;
}
else if (Long >= 9.0 && Long < 21.0) {
ZoneNumber = 33;
}
else if (Long >= 21.0 && Long < 33.0) {
ZoneNumber = 35;
}
else if (Long >= 33.0 && Long < 42.0) {
ZoneNumber = 37;
}
}
LongOrigin = (ZoneNumber - 1) * 6 - 180 + 3; //+3 puts origin
// in middle of
// zone
LongOriginRad = degToRad(LongOrigin);
eccPrimeSquared = (eccSquared) / (1 - eccSquared);
N = a / Math.sqrt(1 - eccSquared * Math.sin(LatRad) * Math.sin(LatRad));
T = Math.tan(LatRad) * Math.tan(LatRad);
C = eccPrimeSquared * Math.cos(LatRad) * Math.cos(LatRad);
A = Math.cos(LatRad) * (LongRad - LongOriginRad);
M = a * ((1 - eccSquared / 4 - 3 * eccSquared * eccSquared / 64 - 5 * eccSquared * eccSquared * eccSquared / 256) * LatRad - (3 * eccSquared / 8 + 3 * eccSquared * eccSquared / 32 + 45 * eccSquared * eccSquared * eccSquared / 1024) * Math.sin(2 * LatRad) + (15 * eccSquared * eccSquared / 256 + 45 * eccSquared * eccSquared * eccSquared / 1024) * Math.sin(4 * LatRad) - (35 * eccSquared * eccSquared * eccSquared / 3072) * Math.sin(6 * LatRad));
var UTMEasting = (k0 * N * (A + (1 - T + C) * A * A * A / 6.0 + (5 - 18 * T + T * T + 72 * C - 58 * eccPrimeSquared) * A * A * A * A * A / 120.0) + 500000.0);
var UTMNorthing = (k0 * (M + N * Math.tan(LatRad) * (A * A / 2 + (5 - T + 9 * C + 4 * C * C) * A * A * A * A / 24.0 + (61 - 58 * T + T * T + 600 * C - 330 * eccPrimeSquared) * A * A * A * A * A * A / 720.0)));
if (Lat < 0.0) {
UTMNorthing += 10000000.0; //10000000 meter offset for
// southern hemisphere
}
return {
northing: Math.round(UTMNorthing),
easting: Math.round(UTMEasting),
zoneNumber: ZoneNumber,
zoneLetter: getLetterDesignator(Lat)
};
}
/**
* Converts UTM coords to lat/long, using the WGS84 ellipsoid. This is a convenience
* class where the Zone can be specified as a single string eg."60N" which
* is then broken down into the ZoneNumber and ZoneLetter.
*
* @private
* @param {object} utm An object literal with northing, easting, zoneNumber
* and zoneLetter properties. If an optional accuracy property is
* provided (in meters), a bounding box will be returned instead of
* latitude and longitude.
* @return {object} An object literal containing either lat and lon values
* (if no accuracy was provided), or top, right, bottom and left values
* for the bounding box calculated according to the provided accuracy.
* Returns null if the conversion failed.
*/
function UTMtoLL(utm) {
var UTMNorthing = utm.northing;
var UTMEasting = utm.easting;
var zoneLetter = utm.zoneLetter;
var zoneNumber = utm.zoneNumber;
// check the ZoneNummber is valid
if (zoneNumber < 0 || zoneNumber > 60) {
return null;
}
var k0 = 0.9996;
var a = 6378137.0; //ellip.radius;
var eccSquared = 0.00669438; //ellip.eccsq;
var eccPrimeSquared;
var e1 = (1 - Math.sqrt(1 - eccSquared)) / (1 + Math.sqrt(1 - eccSquared));
var N1, T1, C1, R1, D, M;
var LongOrigin;
var mu, phi1Rad;
// remove 500,000 meter offset for longitude
var x = UTMEasting - 500000.0;
var y = UTMNorthing;
// We must know somehow if we are in the Northern or Southern
// hemisphere, this is the only time we use the letter So even
// if the Zone letter isn't exactly correct it should indicate
// the hemisphere correctly
if (zoneLetter < 'N') {
y -= 10000000.0; // remove 10,000,000 meter offset used
// for southern hemisphere
}
// There are 60 zones with zone 1 being at West -180 to -174
LongOrigin = (zoneNumber - 1) * 6 - 180 + 3; // +3 puts origin
// in middle of
// zone
eccPrimeSquared = (eccSquared) / (1 - eccSquared);
M = y / k0;
mu = M / (a * (1 - eccSquared / 4 - 3 * eccSquared * eccSquared / 64 - 5 * eccSquared * eccSquared * eccSquared / 256));
phi1Rad = mu + (3 * e1 / 2 - 27 * e1 * e1 * e1 / 32) * Math.sin(2 * mu) + (21 * e1 * e1 / 16 - 55 * e1 * e1 * e1 * e1 / 32) * Math.sin(4 * mu) + (151 * e1 * e1 * e1 / 96) * Math.sin(6 * mu);
// double phi1 = ProjMath.radToDeg(phi1Rad);
N1 = a / Math.sqrt(1 - eccSquared * Math.sin(phi1Rad) * Math.sin(phi1Rad));
T1 = Math.tan(phi1Rad) * Math.tan(phi1Rad);
C1 = eccPrimeSquared * Math.cos(phi1Rad) * Math.cos(phi1Rad);
R1 = a * (1 - eccSquared) / Math.pow(1 - eccSquared * Math.sin(phi1Rad) * Math.sin(phi1Rad), 1.5);
D = x / (N1 * k0);
var lat = phi1Rad - (N1 * Math.tan(phi1Rad) / R1) * (D * D / 2 - (5 + 3 * T1 + 10 * C1 - 4 * C1 * C1 - 9 * eccPrimeSquared) * D * D * D * D / 24 + (61 + 90 * T1 + 298 * C1 + 45 * T1 * T1 - 252 * eccPrimeSquared - 3 * C1 * C1) * D * D * D * D * D * D / 720);
lat = radToDeg(lat);
var lon = (D - (1 + 2 * T1 + C1) * D * D * D / 6 + (5 - 2 * C1 + 28 * T1 - 3 * C1 * C1 + 8 * eccPrimeSquared + 24 * T1 * T1) * D * D * D * D * D / 120) / Math.cos(phi1Rad);
lon = LongOrigin + radToDeg(lon);
var result;
if (utm.accuracy) {
var topRight = UTMtoLL({
northing: utm.northing + utm.accuracy,
easting: utm.easting + utm.accuracy,
zoneLetter: utm.zoneLetter,
zoneNumber: utm.zoneNumber
});
result = {
top: topRight.lat,
right: topRight.lon,
bottom: lat,
left: lon
};
}
else {
result = {
lat: lat,
lon: lon
};
}
return result;
}
/**
* Calculates the MGRS letter designator for the given latitude.
*
* @private
* @param {number} lat The latitude in WGS84 to get the letter designator
* for.
* @return {char} The letter designator.
*/
function getLetterDesignator(lat) {
//This is here as an error flag to show that the Latitude is
//outside MGRS limits
var LetterDesignator = 'Z';
if ((84 >= lat) && (lat >= 72)) {
LetterDesignator = 'X';
}
else if ((72 > lat) && (lat >= 64)) {
LetterDesignator = 'W';
}
else if ((64 > lat) && (lat >= 56)) {
LetterDesignator = 'V';
}
else if ((56 > lat) && (lat >= 48)) {
LetterDesignator = 'U';
}
else if ((48 > lat) && (lat >= 40)) {
LetterDesignator = 'T';
}
else if ((40 > lat) && (lat >= 32)) {
LetterDesignator = 'S';
}
else if ((32 > lat) && (lat >= 24)) {
LetterDesignator = 'R';
}
else if ((24 > lat) && (lat >= 16)) {
LetterDesignator = 'Q';
}
else if ((16 > lat) && (lat >= 8)) {
LetterDesignator = 'P';
}
else if ((8 > lat) && (lat >= 0)) {
LetterDesignator = 'N';
}
else if ((0 > lat) && (lat >= -8)) {
LetterDesignator = 'M';
}
else if ((-8 > lat) && (lat >= -16)) {
LetterDesignator = 'L';
}
else if ((-16 > lat) && (lat >= -24)) {
LetterDesignator = 'K';
}
else if ((-24 > lat) && (lat >= -32)) {
LetterDesignator = 'J';
}
else if ((-32 > lat) && (lat >= -40)) {
LetterDesignator = 'H';
}
else if ((-40 > lat) && (lat >= -48)) {
LetterDesignator = 'G';
}
else if ((-48 > lat) && (lat >= -56)) {
LetterDesignator = 'F';
}
else if ((-56 > lat) && (lat >= -64)) {
LetterDesignator = 'E';
}
else if ((-64 > lat) && (lat >= -72)) {
LetterDesignator = 'D';
}
else if ((-72 > lat) && (lat >= -80)) {
LetterDesignator = 'C';
}
return LetterDesignator;
}
/**
* Encodes a UTM location as MGRS string.
*
* @private
* @param {object} utm An object literal with easting, northing,
* zoneLetter, zoneNumber
* @param {number} accuracy Accuracy in digits (1-5).
* @return {string} MGRS string for the given UTM location.
*/
function encode(utm, accuracy) {
// prepend with leading zeroes
var seasting = "00000" + utm.easting,
snorthing = "00000" + utm.northing;
return utm.zoneNumber + utm.zoneLetter + get100kID(utm.easting, utm.northing, utm.zoneNumber) + seasting.substr(seasting.length - 5, accuracy) + snorthing.substr(snorthing.length - 5, accuracy);
}
/**
* Get the two letter 100k designator for a given UTM easting,
* northing and zone number value.
*
* @private
* @param {number} easting
* @param {number} northing
* @param {number} zoneNumber
* @return the two letter 100k designator for the given UTM location.
*/
function get100kID(easting, northing, zoneNumber) {
var setParm = get100kSetForZone(zoneNumber);
var setColumn = Math.floor(easting / 100000);
var setRow = Math.floor(northing / 100000) % 20;
return getLetter100kID(setColumn, setRow, setParm);
}
/**
* Given a UTM zone number, figure out the MGRS 100K set it is in.
*
* @private
* @param {number} i An UTM zone number.
* @return {number} the 100k set the UTM zone is in.
*/
function get100kSetForZone(i) {
var setParm = i % NUM_100K_SETS;
if (setParm === 0) {
setParm = NUM_100K_SETS;
}
return setParm;
}
/**
* Get the two-letter MGRS 100k designator given information
* translated from the UTM northing, easting and zone number.
*
* @private
* @param {number} column the column index as it relates to the MGRS
* 100k set spreadsheet, created from the UTM easting.
* Values are 1-8.
* @param {number} row the row index as it relates to the MGRS 100k set
* spreadsheet, created from the UTM northing value. Values
* are from 0-19.
* @param {number} parm the set block, as it relates to the MGRS 100k set
* spreadsheet, created from the UTM zone. Values are from
* 1-60.
* @return two letter MGRS 100k code.
*/
function getLetter100kID(column, row, parm) {
// colOrigin and rowOrigin are the letters at the origin of the set
var index = parm - 1;
var colOrigin = SET_ORIGIN_COLUMN_LETTERS.charCodeAt(index);
var rowOrigin = SET_ORIGIN_ROW_LETTERS.charCodeAt(index);
// colInt and rowInt are the letters to build to return
var colInt = colOrigin + column - 1;
var rowInt = rowOrigin + row;
var rollover = false;
if (colInt > Z) {
colInt = colInt - Z + A - 1;
rollover = true;
}
if (colInt === I || (colOrigin < I && colInt > I) || ((colInt > I || colOrigin < I) && rollover)) {
colInt++;
}
if (colInt === O || (colOrigin < O && colInt > O) || ((colInt > O || colOrigin < O) && rollover)) {
colInt++;
if (colInt === I) {
colInt++;
}
}
if (colInt > Z) {
colInt = colInt - Z + A - 1;
}
if (rowInt > V) {
rowInt = rowInt - V + A - 1;
rollover = true;
}
else {
rollover = false;
}
if (((rowInt === I) || ((rowOrigin < I) && (rowInt > I))) || (((rowInt > I) || (rowOrigin < I)) && rollover)) {
rowInt++;
}
if (((rowInt === O) || ((rowOrigin < O) && (rowInt > O))) || (((rowInt > O) || (rowOrigin < O)) && rollover)) {
rowInt++;
if (rowInt === I) {
rowInt++;
}
}
if (rowInt > V) {
rowInt = rowInt - V + A - 1;
}
var twoLetter = String.fromCharCode(colInt) + String.fromCharCode(rowInt);
return twoLetter;
}
/**
* Decode the UTM parameters from a MGRS string.
*
* @private
* @param {string} mgrsString an UPPERCASE coordinate string is expected.
* @return {object} An object literal with easting, northing, zoneLetter,
* zoneNumber and accuracy (in meters) properties.
*/
function decode(mgrsString) {
if (mgrsString && mgrsString.length === 0) {
throw ("MGRSPoint coverting from nothing");
}
var length = mgrsString.length;
var hunK = null;
var sb = "";
var testChar;
var i = 0;
// get Zone number
while (!(/[A-Z]/).test(testChar = mgrsString.charAt(i))) {
if (i >= 2) {
throw ("MGRSPoint bad conversion from: " + mgrsString);
}
sb += testChar;
i++;
}
var zoneNumber = parseInt(sb, 10);
if (i === 0 || i + 3 > length) {
// A good MGRS string has to be 4-5 digits long,
// ##AAA/#AAA at least.
throw ("MGRSPoint bad conversion from: " + mgrsString);
}
var zoneLetter = mgrsString.charAt(i++);
// Should we check the zone letter here? Why not.
if (zoneLetter <= 'A' || zoneLetter === 'B' || zoneLetter === 'Y' || zoneLetter >= 'Z' || zoneLetter === 'I' || zoneLetter === 'O') {
throw ("MGRSPoint zone letter " + zoneLetter + " not handled: " + mgrsString);
}
hunK = mgrsString.substring(i, i += 2);
var set = get100kSetForZone(zoneNumber);
var east100k = getEastingFromChar(hunK.charAt(0), set);
var north100k = getNorthingFromChar(hunK.charAt(1), set);
// We have a bug where the northing may be 2000000 too low.
// How
// do we know when to roll over?
while (north100k < getMinNorthing(zoneLetter)) {
north100k += 2000000;
}
// calculate the char index for easting/northing separator
var remainder = length - i;
if (remainder % 2 !== 0) {
throw ("MGRSPoint has to have an even number \nof digits after the zone letter and two 100km letters - front \nhalf for easting meters, second half for \nnorthing meters" + mgrsString);
}
var sep = remainder / 2;
var sepEasting = 0.0;
var sepNorthing = 0.0;
var accuracyBonus, sepEastingString, sepNorthingString, easting, northing;
if (sep > 0) {
accuracyBonus = 100000.0 / Math.pow(10, sep);
sepEastingString = mgrsString.substring(i, i + sep);
sepEasting = parseFloat(sepEastingString) * accuracyBonus;
sepNorthingString = mgrsString.substring(i + sep);
sepNorthing = parseFloat(sepNorthingString) * accuracyBonus;
}
easting = sepEasting + east100k;
northing = sepNorthing + north100k;
return {
easting: easting,
northing: northing,
zoneLetter: zoneLetter,
zoneNumber: zoneNumber,
accuracy: accuracyBonus
};
}
/**
* Given the first letter from a two-letter MGRS 100k zone, and given the
* MGRS table set for the zone number, figure out the easting value that
* should be added to the other, secondary easting value.
*
* @private
* @param {char} e The first letter from a two-letter MGRS 100´k zone.
* @param {number} set The MGRS table set for the zone number.
* @return {number} The easting value for the given letter and set.
*/
function getEastingFromChar(e, set) {
// colOrigin is the letter at the origin of the set for the
// column
var curCol = SET_ORIGIN_COLUMN_LETTERS.charCodeAt(set - 1);
var eastingValue = 100000.0;
var rewindMarker = false;
while (curCol !== e.charCodeAt(0)) {
curCol++;
if (curCol === I) {
curCol++;
}
if (curCol === O) {
curCol++;
}
if (curCol > Z) {
if (rewindMarker) {
throw ("Bad character: " + e);
}
curCol = A;
rewindMarker = true;
}
eastingValue += 100000.0;
}
return eastingValue;
}
/**
* Given the second letter from a two-letter MGRS 100k zone, and given the
* MGRS table set for the zone number, figure out the northing value that
* should be added to the other, secondary northing value. You have to
* remember that Northings are determined from the equator, and the vertical
* cycle of letters mean a 2000000 additional northing meters. This happens
* approx. every 18 degrees of latitude. This method does *NOT* count any
* additional northings. You have to figure out how many 2000000 meters need
* to be added for the zone letter of the MGRS coordinate.
*
* @private
* @param {char} n Second letter of the MGRS 100k zone
* @param {number} set The MGRS table set number, which is dependent on the
* UTM zone number.
* @return {number} The northing value for the given letter and set.
*/
function getNorthingFromChar(n, set) {
if (n > 'V') {
throw ("MGRSPoint given invalid Northing " + n);
}
// rowOrigin is the letter at the origin of the set for the
// column
var curRow = SET_ORIGIN_ROW_LETTERS.charCodeAt(set - 1);
var northingValue = 0.0;
var rewindMarker = false;
while (curRow !== n.charCodeAt(0)) {
curRow++;
if (curRow === I) {
curRow++;
}
if (curRow === O) {
curRow++;
}
// fixing a bug making whole application hang in this loop
// when 'n' is a wrong character
if (curRow > V) {
if (rewindMarker) { // making sure that this loop ends
throw ("Bad character: " + n);
}
curRow = A;
rewindMarker = true;
}
northingValue += 100000.0;
}
return northingValue;
}
/**
* The function getMinNorthing returns the minimum northing value of a MGRS
* zone.
*
* Ported from Geotrans' c Lattitude_Band_Value structure table.
*
* @private
* @param {char} zoneLetter The MGRS zone to get the min northing for.
* @return {number}
*/
function getMinNorthing(zoneLetter) {
var northing;
switch (zoneLetter) {
case 'C':
northing = 1100000.0;
break;
case 'D':
northing = 2000000.0;
break;
case 'E':
northing = 2800000.0;
break;
case 'F':
northing = 3700000.0;
break;
case 'G':
northing = 4600000.0;
break;
case 'H':
northing = 5500000.0;
break;
case 'J':
northing = 6400000.0;
break;
case 'K':
northing = 7300000.0;
break;
case 'L':
northing = 8200000.0;
break;
case 'M':
northing = 9100000.0;
break;
case 'N':
northing = 0.0;
break;
case 'P':
northing = 800000.0;
break;
case 'Q':
northing = 1700000.0;
break;
case 'R':
northing = 2600000.0;
break;
case 'S':
northing = 3500000.0;
break;
case 'T':
northing = 4400000.0;
break;
case 'U':
northing = 5300000.0;
break;
case 'V':
northing = 6200000.0;
break;
case 'W':
northing = 7000000.0;
break;
case 'X':
northing = 7900000.0;
break;
default:
northing = -1.0;
}
if (northing >= 0.0) {
return northing;
}
else {
throw ("Invalid zone letter: " + zoneLetter);
}
}
function Point(x, y, z) {
if (!(this instanceof Point)) {
return new Point(x, y, z);
}
if (Array.isArray(x)) {
this.x = x[0];
this.y = x[1];
this.z = x[2] || 0.0;
} else if(typeof x === 'object') {
this.x = x.x;
this.y = x.y;
this.z = x.z || 0.0;
} else if (typeof x === 'string' && typeof y === 'undefined') {
var coords = x.split(',');
this.x = parseFloat(coords[0], 10);
this.y = parseFloat(coords[1], 10);
this.z = parseFloat(coords[2], 10) || 0.0;
} else {
this.x = x;
this.y = y;
this.z = z || 0.0;
}
console.warn('proj4.Point will be removed in version 3, use proj4.toPoint');
}
Point.fromMGRS = function(mgrsStr) {
return new Point(toPoint$1(mgrsStr));
};
Point.prototype.toMGRS = function(accuracy) {
return forward$1([this.x, this.y], accuracy);
};
var C00 = 1;
var C02 = 0.25;
var C04 = 0.046875;
var C06 = 0.01953125;
var C08 = 0.01068115234375;
var C22 = 0.75;
var C44 = 0.46875;
var C46 = 0.01302083333333333333;
var C48 = 0.00712076822916666666;
var C66 = 0.36458333333333333333;
var C68 = 0.00569661458333333333;
var C88 = 0.3076171875;
var pj_enfn = function(es) {
var en = [];
en[0] = C00 - es * (C02 + es * (C04 + es * (C06 + es * C08)));
en[1] = es * (C22 - es * (C04 + es * (C06 + es * C08)));
var t = es * es;
en[2] = t * (C44 - es * (C46 + es * C48));
t *= es;
en[3] = t * (C66 - es * C68);
en[4] = t * es * C88;
return en;
};
var pj_mlfn = function(phi, sphi, cphi, en) {
cphi *= sphi;
sphi *= sphi;
return (en[0] * phi - cphi * (en[1] + sphi * (en[2] + sphi * (en[3] + sphi * en[4]))));
};
var MAX_ITER = 20;
var pj_inv_mlfn = function(arg, es, en) {
var k = 1 / (1 - es);
var phi = arg;
for (var i = MAX_ITER; i; --i) { /* rarely goes over 2 iterations */
var s = Math.sin(phi);
var t = 1 - es * s * s;
//t = this.pj_mlfn(phi, s, Math.cos(phi), en) - arg;
//phi -= t * (t * Math.sqrt(t)) * k;
t = (pj_mlfn(phi, s, Math.cos(phi), en) - arg) * (t * Math.sqrt(t)) * k;
phi -= t;
if (Math.abs(t) < EPSLN) {
return phi;
}
}
//..reportError("cass:pj_inv_mlfn: Convergence error");
return phi;
};
// Heavily based on this tmerc projection implementation
// https://github.com/mbloch/mapshaper-proj/blob/master/src/projections/tmerc.js
function init$2() {
this.x0 = this.x0 !== undefined ? this.x0 : 0;
this.y0 = this.y0 !== undefined ? this.y0 : 0;
this.long0 = this.long0 !== undefined ? this.long0 : 0;
this.lat0 = this.lat0 !== undefined ? this.lat0 : 0;
if (this.es) {
this.en = pj_enfn(this.es);
this.ml0 = pj_mlfn(this.lat0, Math.sin(this.lat0), Math.cos(this.lat0), this.en);
}
}
/**
Transverse Mercator Forward - long/lat to x/y
long/lat in radians
*/
function forward$2(p) {
var lon = p.x;
var lat = p.y;
var delta_lon = adjust_lon(lon - this.long0);
var con;
var x, y;
var sin_phi = Math.sin(lat);
var cos_phi = Math.cos(lat);
if (!this.es) {
var b = cos_phi * Math.sin(delta_lon);
if ((Math.abs(Math.abs(b) - 1)) < EPSLN) {
return (93);
}
else {
x = 0.5 * this.a * this.k0 * Math.log((1 + b) / (1 - b)) + this.x0;
y = cos_phi * Math.cos(delta_lon) / Math.sqrt(1 - Math.pow(b, 2));
b = Math.abs(y);
if (b >= 1) {
if ((b - 1) > EPSLN) {
return (93);
}
else {
y = 0;
}
}
else {
y = Math.acos(y);
}
if (lat < 0) {
y = -y;
}
y = this.a * this.k0 * (y - this.lat0) + this.y0;
}
}
else {
var al = cos_phi * delta_lon;
var als = Math.pow(al, 2);
var c = this.ep2 * Math.pow(cos_phi, 2);
var cs = Math.pow(c, 2);
var tq = Math.abs(cos_phi) > EPSLN ? Math.tan(lat) : 0;
var t = Math.pow(tq, 2);
var ts = Math.pow(t, 2);
con = 1 - this.es * Math.pow(sin_phi, 2);
al = al / Math.sqrt(con);
var ml = pj_mlfn(lat, sin_phi, cos_phi, this.en);
x = this.a * (this.k0 * al * (1 +
als / 6 * (1 - t + c +
als / 20 * (5 - 18 * t + ts + 14 * c - 58 * t * c +
als / 42 * (61 + 179 * ts - ts * t - 479 * t))))) +
this.x0;
y = this.a * (this.k0 * (ml - this.ml0 +
sin_phi * delta_lon * al / 2 * (1 +
als / 12 * (5 - t + 9 * c + 4 * cs +
als / 30 * (61 + ts - 58 * t + 270 * c - 330 * t * c +
als / 56 * (1385 + 543 * ts - ts * t - 3111 * t)))))) +
this.y0;
}
p.x = x;
p.y = y;
return p;
}
/**
Transverse Mercator Inverse - x/y to long/lat
*/
function inverse$2(p) {
var con, phi;
var lat, lon;
var x = (p.x - this.x0) * (1 / this.a);
var y = (p.y - this.y0) * (1 / this.a);
if (!this.es) {
var f = Math.exp(x / this.k0);
var g = 0.5 * (f - 1 / f);
var temp = this.lat0 + y / this.k0;
var h = Math.cos(temp);
con = Math.sqrt((1 - Math.pow(h, 2)) / (1 + Math.pow(g, 2)));
lat = Math.asin(con);
if (y < 0) {
lat = -lat;
}
if ((g === 0) && (h === 0)) {
lon = 0;
}
else {
lon = adjust_lon(Math.atan2(g, h) + this.long0);
}
}
else { // ellipsoidal form
con = this.ml0 + y / this.k0;
phi = pj_inv_mlfn(con, this.es, this.en);
if (Math.abs(phi) < HALF_PI) {
var sin_phi = Math.sin(phi);
var cos_phi = Math.cos(phi);
var tan_phi = Math.abs(cos_phi) > EPSLN ? Math.tan(phi) : 0;
var c = this.ep2 * Math.pow(cos_phi, 2);
var cs = Math.pow(c, 2);
var t = Math.pow(tan_phi, 2);
var ts = Math.pow(t, 2);
con = 1 - this.es * Math.pow(sin_phi, 2);
var d = x * Math.sqrt(con) / this.k0;
var ds = Math.pow(d, 2);
con = con * tan_phi;
lat = phi - (con * ds / (1 - this.es)) * 0.5 * (1 -
ds / 12 * (5 + 3 * t - 9 * c * t + c - 4 * cs -
ds / 30 * (61 + 90 * t - 252 * c * t + 45 * ts + 46 * c -
ds / 56 * (1385 + 3633 * t + 4095 * ts + 1574 * ts * t))));
lon = adjust_lon(this.long0 + (d * (1 -
ds / 6 * (1 + 2 * t + c -
ds / 20 * (5 + 28 * t + 24 * ts + 8 * c * t + 6 * c -
ds / 42 * (61 + 662 * t + 1320 * ts + 720 * ts * t)))) / cos_phi));
}
else {
lat = HALF_PI * sign(y);
lon = 0;
}
}
p.x = lon;
p.y = lat;
return p;
}
var names$3 = ["Fast_Transverse_Mercator", "Fast Transverse Mercator"];
var tmerc = {
init: init$2,
forward: forward$2,
inverse: inverse$2,
names: names$3
};
var sinh = function(x) {
var r = Math.exp(x);
r = (r - 1 / r) / 2;
return r;
};
var hypot = function(x, y) {
x = Math.abs(x);
y = Math.abs(y);
var a = Math.max(x, y);
var b = Math.min(x, y) / (a ? a : 1);
return a * Math.sqrt(1 + Math.pow(b, 2));
};
var log1py = function(x) {
var y = 1 + x;
var z = y - 1;
return z === 0 ? x : x * Math.log(y) / z;
};
var asinhy = function(x) {
var y = Math.abs(x);
y = log1py(y * (1 + y / (hypot(1, y) + 1)));
return x < 0 ? -y : y;
};
var gatg = function(pp, B) {
var cos_2B = 2 * Math.cos(2 * B);
var i = pp.length - 1;
var h1 = pp[i];
var h2 = 0;
var h;
while (--i >= 0) {
h = -h2 + cos_2B * h1 + pp[i];
h2 = h1;
h1 = h;
}
return (B + h * Math.sin(2 * B));
};
var clens = function(pp, arg_r) {
var r = 2 * Math.cos(arg_r);
var i = pp.length - 1;
var hr1 = pp[i];
var hr2 = 0;
var hr;
while (--i >= 0) {
hr = -hr2 + r * hr1 + pp[i];
hr2 = hr1;
hr1 = hr;
}
return Math.sin(arg_r) * hr;
};
var cosh = function(x) {
var r = Math.exp(x);
r = (r + 1 / r) / 2;
return r;
};
var clens_cmplx = function(pp, arg_r, arg_i) {
var sin_arg_r = Math.sin(arg_r);
var cos_arg_r = Math.cos(arg_r);
var sinh_arg_i = sinh(arg_i);
var cosh_arg_i = cosh(arg_i);
var r = 2 * cos_arg_r * cosh_arg_i;
var i = -2 * sin_arg_r * sinh_arg_i;
var j = pp.length - 1;
var hr = pp[j];
var hi1 = 0;
var hr1 = 0;
var hi = 0;
var hr2;
var hi2;
while (--j >= 0) {
hr2 = hr1;
hi2 = hi1;
hr1 = hr;
hi1 = hi;
hr = -hr2 + r * hr1 - i * hi1 + pp[j];
hi = -hi2 + i * hr1 + r * hi1;
}
r = sin_arg_r * cosh_arg_i;
i = cos_arg_r * sinh_arg_i;
return [r * hr - i * hi, r * hi + i * hr];
};
// Heavily based on this etmerc projection implementation
// https://github.com/mbloch/mapshaper-proj/blob/master/src/projections/etmerc.js
function init$3() {
if (!this.approx && (isNaN(this.es) || this.es <= 0)) {
throw new Error('Incorrect elliptical usage. Try using the +approx option in the proj string, or PROJECTION["Fast_Transverse_Mercator"] in the WKT.');
}
if (this.approx) {
// When '+approx' is set, use tmerc instead
tmerc.init.apply(this);
this.forward = tmerc.forward;
this.inverse = tmerc.inverse;
}
this.x0 = this.x0 !== undefined ? this.x0 : 0;
this.y0 = this.y0 !== undefined ? this.y0 : 0;
this.long0 = this.long0 !== undefined ? this.long0 : 0;
this.lat0 = this.lat0 !== undefined ? this.lat0 : 0;
this.cgb = [];
this.cbg = [];
this.utg = [];
this.gtu = [];
var f = this.es / (1 + Math.sqrt(1 - this.es));
var n = f / (2 - f);
var np = n;
this.cgb[0] = n * (2 + n * (-2 / 3 + n * (-2 + n * (116 / 45 + n * (26 / 45 + n * (-2854 / 675 ))))));
this.cbg[0] = n * (-2 + n * ( 2 / 3 + n * ( 4 / 3 + n * (-82 / 45 + n * (32 / 45 + n * (4642 / 4725))))));
np = np * n;
this.cgb[1] = np * (7 / 3 + n * (-8 / 5 + n * (-227 / 45 + n * (2704 / 315 + n * (2323 / 945)))));
this.cbg[1] = np * (5 / 3 + n * (-16 / 15 + n * ( -13 / 9 + n * (904 / 315 + n * (-1522 / 945)))));
np = np * n;
this.cgb[2] = np * (56 / 15 + n * (-136 / 35 + n * (-1262 / 105 + n * (73814 / 2835))));
this.cbg[2] = np * (-26 / 15 + n * (34 / 21 + n * (8 / 5 + n * (-12686 / 2835))));
np = np * n;
this.cgb[3] = np * (4279 / 630 + n * (-332 / 35 + n * (-399572 / 14175)));
this.cbg[3] = np * (1237 / 630 + n * (-12 / 5 + n * ( -24832 / 14175)));
np = np * n;
this.cgb[4] = np * (4174 / 315 + n * (-144838 / 6237));
this.cbg[4] = np * (-734 / 315 + n * (109598 / 31185));
np = np * n;
this.cgb[5] = np * (601676 / 22275);
this.cbg[5] = np * (444337 / 155925);
np = Math.pow(n, 2);
this.Qn = this.k0 / (1 + n) * (1 + np * (1 / 4 + np * (1 / 64 + np / 256)));
this.utg[0] = n * (-0.5 + n * ( 2 / 3 + n * (-37 / 96 + n * ( 1 / 360 + n * (81 / 512 + n * (-96199 / 604800))))));
this.gtu[0] = n * (0.5 + n * (-2 / 3 + n * (5 / 16 + n * (41 / 180 + n * (-127 / 288 + n * (7891 / 37800))))));
this.utg[1] = np * (-1 / 48 + n * (-1 / 15 + n * (437 / 1440 + n * (-46 / 105 + n * (1118711 / 3870720)))));
this.gtu[1] = np * (13 / 48 + n * (-3 / 5 + n * (557 / 1440 + n * (281 / 630 + n * (-1983433 / 1935360)))));
np = np * n;
this.utg[2] = np * (-17 / 480 + n * (37 / 840 + n * (209 / 4480 + n * (-5569 / 90720 ))));
this.gtu[2] = np * (61 / 240 + n * (-103 / 140 + n * (15061 / 26880 + n * (167603 / 181440))));
np = np * n;
this.utg[3] = np * (-4397 / 161280 + n * (11 / 504 + n * (830251 / 7257600)));
this.gtu[3] = np * (49561 / 161280 + n * (-179 / 168 + n * (6601661 / 7257600)));
np = np * n;
this.utg[4] = np * (-4583 / 161280 + n * (108847 / 3991680));
this.gtu[4] = np * (34729 / 80640 + n * (-3418889 / 1995840));
np = np * n;
this.utg[5] = np * (-20648693 / 638668800);
this.gtu[5] = np * (212378941 / 319334400);
var Z = gatg(this.cbg, this.lat0);
this.Zb = -this.Qn * (Z + clens(this.gtu, 2 * Z));
}
function forward$3(p) {
var Ce = adjust_lon(p.x - this.long0);
var Cn = p.y;
Cn = gatg(this.cbg, Cn);
var sin_Cn = Math.sin(Cn);
var cos_Cn = Math.cos(Cn);
var sin_Ce = Math.sin(Ce);
var cos_Ce = Math.cos(Ce);
Cn = Math.atan2(sin_Cn, cos_Ce * cos_Cn);
Ce = Math.atan2(sin_Ce * cos_Cn, hypot(sin_Cn, cos_Cn * cos_Ce));
Ce = asinhy(Math.tan(Ce));
var tmp = clens_cmplx(this.gtu, 2 * Cn, 2 * Ce);
Cn = Cn + tmp[0];
Ce = Ce + tmp[1];
var x;
var y;
if (Math.abs(Ce) <= 2.623395162778) {
x = this.a * (this.Qn * Ce) + this.x0;
y = this.a * (this.Qn * Cn + this.Zb) + this.y0;
}
else {
x = Infinity;
y = Infinity;
}
p.x = x;
p.y = y;
return p;
}
function inverse$3(p) {
var Ce = (p.x - this.x0) * (1 / this.a);
var Cn = (p.y - this.y0) * (1 / this.a);
Cn = (Cn - this.Zb) / this.Qn;
Ce = Ce / this.Qn;
var lon;
var lat;
if (Math.abs(Ce) <= 2.623395162778) {
var tmp = clens_cmplx(this.utg, 2 * Cn, 2 * Ce);
Cn = Cn + tmp[0];
Ce = Ce + tmp[1];
Ce = Math.atan(sinh(Ce));
var sin_Cn = Math.sin(Cn);
var cos_Cn = Math.cos(Cn);
var sin_Ce = Math.sin(Ce);
var cos_Ce = Math.cos(Ce);
Cn = Math.atan2(sin_Cn * cos_Ce, hypot(sin_Ce, cos_Ce * cos_Cn));
Ce = Math.atan2(sin_Ce, cos_Ce * cos_Cn);
lon = adjust_lon(Ce + this.long0);
lat = gatg(this.cgb, Cn);
}
else {
lon = Infinity;
lat = Infinity;
}
p.x = lon;
p.y = lat;
return p;
}
var names$4 = ["Extended_Transverse_Mercator", "Extended Transverse Mercator", "etmerc", "Transverse_Mercator", "Transverse Mercator", "tmerc"];
var etmerc = {
init: init$3,
forward: forward$3,
inverse: inverse$3,
names: names$4
};
var adjust_zone = function(zone, lon) {
if (zone === undefined) {
zone = Math.floor((adjust_lon(lon) + Math.PI) * 30 / Math.PI) + 1;
if (zone < 0) {
return 0;
} else if (zone > 60) {
return 60;
}
}
return zone;
};
var dependsOn = 'etmerc';
function init$4() {
var zone = adjust_zone(this.zone, this.long0);
if (zone === undefined) {
throw new Error('unknown utm zone');
}
this.lat0 = 0;
this.long0 = ((6 * Math.abs(zone)) - 183) * D2R;
this.x0 = 500000;
this.y0 = this.utmSouth ? 10000000 : 0;
this.k0 = 0.9996;
etmerc.init.apply(this);
this.forward = etmerc.forward;
this.inverse = etmerc.inverse;
}
var names$5 = ["Universal Transverse Mercator System", "utm"];
var utm = {
init: init$4,
names: names$5,
dependsOn: dependsOn
};
var srat = function(esinp, exp) {
return (Math.pow((1 - esinp) / (1 + esinp), exp));
};
var MAX_ITER$1 = 20;
function init$6() {
var sphi = Math.sin(this.lat0);
var cphi = Math.cos(this.lat0);
cphi *= cphi;
this.rc = Math.sqrt(1 - this.es) / (1 - this.es * sphi * sphi);
this.C = Math.sqrt(1 + this.es * cphi * cphi / (1 - this.es));
this.phic0 = Math.asin(sphi / this.C);
this.ratexp = 0.5 * this.C * this.e;
this.K = Math.tan(0.5 * this.phic0 + FORTPI) / (Math.pow(Math.tan(0.5 * this.lat0 + FORTPI), this.C) * srat(this.e * sphi, this.ratexp));
}
function forward$5(p) {
var lon = p.x;
var lat = p.y;
p.y = 2 * Math.atan(this.K * Math.pow(Math.tan(0.5 * lat + FORTPI), this.C) * srat(this.e * Math.sin(lat), this.ratexp)) - HALF_PI;
p.x = this.C * lon;
return p;
}
function inverse$5(p) {
var DEL_TOL = 1e-14;
var lon = p.x / this.C;
var lat = p.y;
var num = Math.pow(Math.tan(0.5 * lat + FORTPI) / this.K, 1 / this.C);
for (var i = MAX_ITER$1; i > 0; --i) {
lat = 2 * Math.atan(num * srat(this.e * Math.sin(p.y), - 0.5 * this.e)) - HALF_PI;
if (Math.abs(lat - p.y) < DEL_TOL) {
break;
}
p.y = lat;
}
/* convergence failed */
if (!i) {
return null;
}
p.x = lon;
p.y = lat;
return p;
}
var names$7 = ["gauss"];
var gauss = {
init: init$6,
forward: forward$5,
inverse: inverse$5,
names: names$7
};
function init$5() {
gauss.init.apply(this);
if (!this.rc) {
return;
}
this.sinc0 = Math.sin(this.phic0);
this.cosc0 = Math.cos(this.phic0);
this.R2 = 2 * this.rc;
if (!this.title) {
this.title = "Oblique Stereographic Alternative";
}
}
function forward$4(p) {
var sinc, cosc, cosl, k;
p.x = adjust_lon(p.x - this.long0);
gauss.forward.apply(this, [p]);
sinc = Math.sin(p.y);
cosc = Math.cos(p.y);
cosl = Math.cos(p.x);
k = this.k0 * this.R2 / (1 + this.sinc0 * sinc + this.cosc0 * cosc * cosl);
p.x = k * cosc * Math.sin(p.x);
p.y = k * (this.cosc0 * sinc - this.sinc0 * cosc * cosl);
p.x = this.a * p.x + this.x0;
p.y = this.a * p.y + this.y0;
return p;
}
function inverse$4(p) {
var sinc, cosc, lon, lat, rho;
p.x = (p.x - this.x0) / this.a;
p.y = (p.y - this.y0) / this.a;
p.x /= this.k0;
p.y /= this.k0;
if ((rho = Math.sqrt(p.x * p.x + p.y * p.y))) {
var c = 2 * Math.atan2(rho, this.R2);
sinc = Math.sin(c);
cosc = Math.cos(c);
lat = Math.asin(cosc * this.sinc0 + p.y * sinc * this.cosc0 / rho);
lon = Math.atan2(p.x * sinc, rho * this.cosc0 * cosc - p.y * this.sinc0 * sinc);
}
else {
lat = this.phic0;
lon = 0;
}
p.x = lon;
p.y = lat;
gauss.inverse.apply(this, [p]);
p.x = adjust_lon(p.x + this.long0);
return p;
}
var names$6 = ["Stereographic_North_Pole", "Oblique_Stereographic", "Polar_Stereographic", "sterea","Oblique Stereographic Alternative","Double_Stereographic"];
var sterea = {
init: init$5,
forward: forward$4,
inverse: inverse$4,
names: names$6
};
function ssfn_(phit, sinphi, eccen) {
sinphi *= eccen;
return (Math.tan(0.5 * (HALF_PI + phit)) * Math.pow((1 - sinphi) / (1 + sinphi), 0.5 * eccen));
}
function init$7() {
this.coslat0 = Math.cos(this.lat0);
this.sinlat0 = Math.sin(this.lat0);
if (this.sphere) {
if (this.k0 === 1 && !isNaN(this.lat_ts) && Math.abs(this.coslat0) <= EPSLN) {
this.k0 = 0.5 * (1 + sign(this.lat0) * Math.sin(this.lat_ts));
}
}
else {
if (Math.abs(this.coslat0) <= EPSLN) {
if (this.lat0 > 0) {
//North pole
//trace('stere:north pole');
this.con = 1;
}
else {
//South pole
//trace('stere:south pole');
this.con = -1;
}
}
this.cons = Math.sqrt(Math.pow(1 + this.e, 1 + this.e) * Math.pow(1 - this.e, 1 - this.e));
if (this.k0 === 1 && !isNaN(this.lat_ts) && Math.abs(this.coslat0) <= EPSLN) {
this.k0 = 0.5 * this.cons * msfnz(this.e, Math.sin(this.lat_ts), Math.cos(this.lat_ts)) / tsfnz(this.e, this.con * this.lat_ts, this.con * Math.sin(this.lat_ts));
}
this.ms1 = msfnz(this.e, this.sinlat0, this.coslat0);
this.X0 = 2 * Math.atan(this.ssfn_(this.lat0, this.sinlat0, this.e)) - HALF_PI;
this.cosX0 = Math.cos(this.X0);
this.sinX0 = Math.sin(this.X0);
}
}
// Stereographic forward equations--mapping lat,long to x,y
function forward$6(p) {
var lon = p.x;
var lat = p.y;
var sinlat = Math.sin(lat);
var coslat = Math.cos(lat);
var A, X, sinX, cosX, ts, rh;
var dlon = adjust_lon(lon - this.long0);
if (Math.abs(Math.abs(lon - this.long0) - Math.PI) <= EPSLN && Math.abs(lat + this.lat0) <= EPSLN) {
//case of the origine point
//trace('stere:this is the origin point');
p.x = NaN;
p.y = NaN;
return p;
}
if (this.sphere) {
//trace('stere:sphere case');
A = 2 * this.k0 / (1 + this.sinlat0 * sinlat + this.coslat0 * coslat * Math.cos(dlon));
p.x = this.a * A * coslat * Math.sin(dlon) + this.x0;
p.y = this.a * A * (this.coslat0 * sinlat - this.sinlat0 * coslat * Math.cos(dlon)) + this.y0;
return p;
}
else {
X = 2 * Math.atan(this.ssfn_(lat, sinlat, this.e)) - HALF_PI;
cosX = Math.cos(X);
sinX = Math.sin(X);
if (Math.abs(this.coslat0) <= EPSLN) {
ts = tsfnz(this.e, lat * this.con, this.con * sinlat);
rh = 2 * this.a * this.k0 * ts / this.cons;
p.x = this.x0 + rh * Math.sin(lon - this.long0);
p.y = this.y0 - this.con * rh * Math.cos(lon - this.long0);
//trace(p.toString());
return p;
}
else if (Math.abs(this.sinlat0) < EPSLN) {
//Eq
//trace('stere:equateur');
A = 2 * this.a * this.k0 / (1 + cosX * Math.cos(dlon));
p.y = A * sinX;
}
else {
//other case
//trace('stere:normal case');
A = 2 * this.a * this.k0 * this.ms1 / (this.cosX0 * (1 + this.sinX0 * sinX + this.cosX0 * cosX * Math.cos(dlon)));
p.y = A * (this.cosX0 * sinX - this.sinX0 * cosX * Math.cos(dlon)) + this.y0;
}
p.x = A * cosX * Math.sin(dlon) + this.x0;
}
//trace(p.toString());
return p;
}
//* Stereographic inverse equations--mapping x,y to lat/long
function inverse$6(p) {
p.x -= this.x0;
p.y -= this.y0;
var lon, lat, ts, ce, Chi;
var rh = Math.sqrt(p.x * p.x + p.y * p.y);
if (this.sphere) {
var c = 2 * Math.atan(rh / (2 * this.a * this.k0));
lon = this.long0;
lat = this.lat0;
if (rh <= EPSLN) {
p.x = lon;
p.y = lat;
return p;
}
lat = Math.asin(Math.cos(c) * this.sinlat0 + p.y * Math.sin(c) * this.coslat0 / rh);
if (Math.abs(this.coslat0) < EPSLN) {
if (this.lat0 > 0) {
lon = adjust_lon(this.long0 + Math.atan2(p.x, - 1 * p.y));
}
else {
lon = adjust_lon(this.long0 + Math.atan2(p.x, p.y));
}
}
else {
lon = adjust_lon(this.long0 + Math.atan2(p.x * Math.sin(c), rh * this.coslat0 * Math.cos(c) - p.y * this.sinlat0 * Math.sin(c)));
}
p.x = lon;
p.y = lat;
return p;
}
else {
if (Math.abs(this.coslat0) <= EPSLN) {
if (rh <= EPSLN) {
lat = this.lat0;
lon = this.long0;
p.x = lon;
p.y = lat;
//trace(p.toString());
return p;
}
p.x *= this.con;
p.y *= this.con;
ts = rh * this.cons / (2 * this.a * this.k0);
lat = this.con * phi2z(this.e, ts);
lon = this.con * adjust_lon(this.con * this.long0 + Math.atan2(p.x, - 1 * p.y));
}
else {
ce = 2 * Math.atan(rh * this.cosX0 / (2 * this.a * this.k0 * this.ms1));
lon = this.long0;
if (rh <= EPSLN) {
Chi = this.X0;
}
else {
Chi = Math.asin(Math.cos(ce) * this.sinX0 + p.y * Math.sin(ce) * this.cosX0 / rh);
lon = adjust_lon(this.long0 + Math.atan2(p.x * Math.sin(ce), rh * this.cosX0 * Math.cos(ce) - p.y * this.sinX0 * Math.sin(ce)));
}
lat = -1 * phi2z(this.e, Math.tan(0.5 * (HALF_PI + Chi)));
}
}
p.x = lon;
p.y = lat;
//trace(p.toString());
return p;
}
var names$8 = ["stere", "Stereographic_South_Pole", "Polar Stereographic (variant B)"];
var stere = {
init: init$7,
forward: forward$6,
inverse: inverse$6,
names: names$8,
ssfn_: ssfn_
};
/*
references:
Formules et constantes pour le Calcul pour la
projection cylindrique conforme à axe oblique et pour la transformation entre
des systèmes de référence.
http://www.swisstopo.admin.ch/internet/swisstopo/fr/home/topics/survey/sys/refsys/switzerland.parsysrelated1.31216.downloadList.77004.DownloadFile.tmp/swissprojectionfr.pdf
*/
function init$8() {
var phy0 = this.lat0;
this.lambda0 = this.long0;
var sinPhy0 = Math.sin(phy0);
var semiMajorAxis = this.a;
var invF = this.rf;
var flattening = 1 / invF;
var e2 = 2 * flattening - Math.pow(flattening, 2);
var e = this.e = Math.sqrt(e2);
this.R = this.k0 * semiMajorAxis * Math.sqrt(1 - e2) / (1 - e2 * Math.pow(sinPhy0, 2));
this.alpha = Math.sqrt(1 + e2 / (1 - e2) * Math.pow(Math.cos(phy0), 4));
this.b0 = Math.asin(sinPhy0 / this.alpha);
var k1 = Math.log(Math.tan(Math.PI / 4 + this.b0 / 2));
var k2 = Math.log(Math.tan(Math.PI / 4 + phy0 / 2));
var k3 = Math.log((1 + e * sinPhy0) / (1 - e * sinPhy0));
this.K = k1 - this.alpha * k2 + this.alpha * e / 2 * k3;
}
function forward$7(p) {
var Sa1 = Math.log(Math.tan(Math.PI / 4 - p.y / 2));
var Sa2 = this.e / 2 * Math.log((1 + this.e * Math.sin(p.y)) / (1 - this.e * Math.sin(p.y)));
var S = -this.alpha * (Sa1 + Sa2) + this.K;
// spheric latitude
var b = 2 * (Math.atan(Math.exp(S)) - Math.PI / 4);
// spheric longitude
var I = this.alpha * (p.x - this.lambda0);
// psoeudo equatorial rotation
var rotI = Math.atan(Math.sin(I) / (Math.sin(this.b0) * Math.tan(b) + Math.cos(this.b0) * Math.cos(I)));
var rotB = Math.asin(Math.cos(this.b0) * Math.sin(b) - Math.sin(this.b0) * Math.cos(b) * Math.cos(I));
p.y = this.R / 2 * Math.log((1 + Math.sin(rotB)) / (1 - Math.sin(rotB))) + this.y0;
p.x = this.R * rotI + this.x0;
return p;
}
function inverse$7(p) {
var Y = p.x - this.x0;
var X = p.y - this.y0;
var rotI = Y / this.R;
var rotB = 2 * (Math.atan(Math.exp(X / this.R)) - Math.PI / 4);
var b = Math.asin(Math.cos(this.b0) * Math.sin(rotB) + Math.sin(this.b0) * Math.cos(rotB) * Math.cos(rotI));
var I = Math.atan(Math.sin(rotI) / (Math.cos(this.b0) * Math.cos(rotI) - Math.sin(this.b0) * Math.tan(rotB)));
var lambda = this.lambda0 + I / this.alpha;
var S = 0;
var phy = b;
var prevPhy = -1000;
var iteration = 0;
while (Math.abs(phy - prevPhy) > 0.0000001) {
if (++iteration > 20) {
//...reportError("omercFwdInfinity");
return;
}
//S = Math.log(Math.tan(Math.PI / 4 + phy / 2));
S = 1 / this.alpha * (Math.log(Math.tan(Math.PI / 4 + b / 2)) - this.K) + this.e * Math.log(Math.tan(Math.PI / 4 + Math.asin(this.e * Math.sin(phy)) / 2));
prevPhy = phy;
phy = 2 * Math.atan(Math.exp(S)) - Math.PI / 2;
}
p.x = lambda;
p.y = phy;
return p;
}
var names$9 = ["somerc"];
var somerc = {
init: init$8,
forward: forward$7,
inverse: inverse$7,
names: names$9
};
var TOL = 1e-7;
function isTypeA(P) {
var typeAProjections = ['Hotine_Oblique_Mercator','Hotine_Oblique_Mercator_Azimuth_Natural_Origin'];
var projectionName = typeof P.PROJECTION === "object" ? Object.keys(P.PROJECTION)[0] : P.PROJECTION;
return 'no_uoff' in P || 'no_off' in P || typeAProjections.indexOf(projectionName) !== -1;
}
/* Initialize the Oblique Mercator projection
------------------------------------------*/
function init$9() {
var con, com, cosph0, D, F, H, L, sinph0, p, J, gamma = 0,
gamma0, lamc = 0, lam1 = 0, lam2 = 0, phi1 = 0, phi2 = 0, alpha_c = 0;
// only Type A uses the no_off or no_uoff property
// https://github.com/OSGeo/proj.4/issues/104
this.no_off = isTypeA(this);
this.no_rot = 'no_rot' in this;
var alp = false;
if ("alpha" in this) {
alp = true;
}
var gam = false;
if ("rectified_grid_angle" in this) {
gam = true;
}
if (alp) {
alpha_c = this.alpha;
}
if (gam) {
gamma = (this.rectified_grid_angle * D2R);
}
if (alp || gam) {
lamc = this.longc;
} else {
lam1 = this.long1;
phi1 = this.lat1;
lam2 = this.long2;
phi2 = this.lat2;
if (Math.abs(phi1 - phi2) <= TOL || (con = Math.abs(phi1)) <= TOL ||
Math.abs(con - HALF_PI) <= TOL || Math.abs(Math.abs(this.lat0) - HALF_PI) <= TOL ||
Math.abs(Math.abs(phi2) - HALF_PI) <= TOL) {
throw new Error();
}
}
var one_es = 1.0 - this.es;
com = Math.sqrt(one_es);
if (Math.abs(this.lat0) > EPSLN) {
sinph0 = Math.sin(this.lat0);
cosph0 = Math.cos(this.lat0);
con = 1 - this.es * sinph0 * sinph0;
this.B = cosph0 * cosph0;
this.B = Math.sqrt(1 + this.es * this.B * this.B / one_es);
this.A = this.B * this.k0 * com / con;
D = this.B * com / (cosph0 * Math.sqrt(con));
F = D * D -1;
if (F <= 0) {
F = 0;
} else {
F = Math.sqrt(F);
if (this.lat0 < 0) {
F = -F;
}
}
this.E = F += D;
this.E *= Math.pow(tsfnz(this.e, this.lat0, sinph0), this.B);
} else {
this.B = 1 / com;
this.A = this.k0;
this.E = D = F = 1;
}
if (alp || gam) {
if (alp) {
gamma0 = Math.asin(Math.sin(alpha_c) / D);
if (!gam) {
gamma = alpha_c;
}
} else {
gamma0 = gamma;
alpha_c = Math.asin(D * Math.sin(gamma0));
}
this.lam0 = lamc - Math.asin(0.5 * (F - 1 / F) * Math.tan(gamma0)) / this.B;
} else {
H = Math.pow(tsfnz(this.e, phi1, Math.sin(phi1)), this.B);
L = Math.pow(tsfnz(this.e, phi2, Math.sin(phi2)), this.B);
F = this.E / H;
p = (L - H) / (L + H);
J = this.E * this.E;
J = (J - L * H) / (J + L * H);
con = lam1 - lam2;
if (con < -Math.pi) {
lam2 -=TWO_PI;
} else if (con > Math.pi) {
lam2 += TWO_PI;
}
this.lam0 = adjust_lon(0.5 * (lam1 + lam2) - Math.atan(J * Math.tan(0.5 * this.B * (lam1 - lam2)) / p) / this.B);
gamma0 = Math.atan(2 * Math.sin(this.B * adjust_lon(lam1 - this.lam0)) / (F - 1 / F));
gamma = alpha_c = Math.asin(D * Math.sin(gamma0));
}
this.singam = Math.sin(gamma0);
this.cosgam = Math.cos(gamma0);
this.sinrot = Math.sin(gamma);
this.cosrot = Math.cos(gamma);
this.rB = 1 / this.B;
this.ArB = this.A * this.rB;
this.BrA = 1 / this.ArB;
if (this.no_off) {
this.u_0 = 0;
} else {
this.u_0 = Math.abs(this.ArB * Math.atan(Math.sqrt(D * D - 1) / Math.cos(alpha_c)));
if (this.lat0 < 0) {
this.u_0 = - this.u_0;
}
}
F = 0.5 * gamma0;
this.v_pole_n = this.ArB * Math.log(Math.tan(FORTPI - F));
this.v_pole_s = this.ArB * Math.log(Math.tan(FORTPI + F));
}
/* Oblique Mercator forward equations--mapping lat,long to x,y
----------------------------------------------------------*/
function forward$8(p) {
var coords = {};
var S, T, U, V, W, temp, u, v;
p.x = p.x - this.lam0;
if (Math.abs(Math.abs(p.y) - HALF_PI) > EPSLN) {
W = this.E / Math.pow(tsfnz(this.e, p.y, Math.sin(p.y)), this.B);
temp = 1 / W;
S = 0.5 * (W - temp);
T = 0.5 * (W + temp);
V = Math.sin(this.B * p.x);
U = (S * this.singam - V * this.cosgam) / T;
if (Math.abs(Math.abs(U) - 1.0) < EPSLN) {
throw new Error();
}
v = 0.5 * this.ArB * Math.log((1 - U)/(1 + U));
temp = Math.cos(this.B * p.x);
if (Math.abs(temp) < TOL) {
u = this.A * p.x;
} else {
u = this.ArB * Math.atan2((S * this.cosgam + V * this.singam), temp);
}
} else {
v = p.y > 0 ? this.v_pole_n : this.v_pole_s;
u = this.ArB * p.y;
}
if (this.no_rot) {
coords.x = u;
coords.y = v;
} else {
u -= this.u_0;
coords.x = v * this.cosrot + u * this.sinrot;
coords.y = u * this.cosrot - v * this.sinrot;
}
coords.x = (this.a * coords.x + this.x0);
coords.y = (this.a * coords.y + this.y0);
return coords;
}
function inverse$8(p) {
var u, v, Qp, Sp, Tp, Vp, Up;
var coords = {};
p.x = (p.x - this.x0) * (1.0 / this.a);
p.y = (p.y - this.y0) * (1.0 / this.a);
if (this.no_rot) {
v = p.y;
u = p.x;
} else {
v = p.x * this.cosrot - p.y * this.sinrot;
u = p.y * this.cosrot + p.x * this.sinrot + this.u_0;
}
Qp = Math.exp(-this.BrA * v);
Sp = 0.5 * (Qp - 1 / Qp);
Tp = 0.5 * (Qp + 1 / Qp);
Vp = Math.sin(this.BrA * u);
Up = (Vp * this.cosgam + Sp * this.singam) / Tp;
if (Math.abs(Math.abs(Up) - 1) < EPSLN) {
coords.x = 0;
coords.y = Up < 0 ? -HALF_PI : HALF_PI;
} else {
coords.y = this.E / Math.sqrt((1 + Up) / (1 - Up));
coords.y = phi2z(this.e, Math.pow(coords.y, 1 / this.B));
if (coords.y === Infinity) {
throw new Error();
}
coords.x = -this.rB * Math.atan2((Sp * this.cosgam - Vp * this.singam), Math.cos(this.BrA * u));
}
coords.x += this.lam0;
return coords;
}
var names$10 = ["Hotine_Oblique_Mercator", "Hotine Oblique Mercator", "Hotine_Oblique_Mercator_Azimuth_Natural_Origin", "Hotine_Oblique_Mercator_Two_Point_Natural_Origin", "Hotine_Oblique_Mercator_Azimuth_Center", "Oblique_Mercator", "omerc"];
var omerc = {
init: init$9,
forward: forward$8,
inverse: inverse$8,
names: names$10
};
function init$10() {
//double lat0; /* the reference latitude */
//double long0; /* the reference longitude */
//double lat1; /* first standard parallel */
//double lat2; /* second standard parallel */
//double r_maj; /* major axis */
//double r_min; /* minor axis */
//double false_east; /* x offset in meters */
//double false_north; /* y offset in meters */
//the above value can be set with proj4.defs
//example: proj4.defs("EPSG:2154","+proj=lcc +lat_1=49 +lat_2=44 +lat_0=46.5 +lon_0=3 +x_0=700000 +y_0=6600000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs");
if (!this.lat2) {
this.lat2 = this.lat1;
} //if lat2 is not defined
if (!this.k0) {
this.k0 = 1;
}
this.x0 = this.x0 || 0;
this.y0 = this.y0 || 0;
// Standard Parallels cannot be equal and on opposite sides of the equator
if (Math.abs(this.lat1 + this.lat2) < EPSLN) {
return;
}
var temp = this.b / this.a;
this.e = Math.sqrt(1 - temp * temp);
var sin1 = Math.sin(this.lat1);
var cos1 = Math.cos(this.lat1);
var ms1 = msfnz(this.e, sin1, cos1);
var ts1 = tsfnz(this.e, this.lat1, sin1);
var sin2 = Math.sin(this.lat2);
var cos2 = Math.cos(this.lat2);
var ms2 = msfnz(this.e, sin2, cos2);
var ts2 = tsfnz(this.e, this.lat2, sin2);
var ts0 = tsfnz(this.e, this.lat0, Math.sin(this.lat0));
if (Math.abs(this.lat1 - this.lat2) > EPSLN) {
this.ns = Math.log(ms1 / ms2) / Math.log(ts1 / ts2);
}
else {
this.ns = sin1;
}
if (isNaN(this.ns)) {
this.ns = sin1;
}
this.f0 = ms1 / (this.ns * Math.pow(ts1, this.ns));
this.rh = this.a * this.f0 * Math.pow(ts0, this.ns);
if (!this.title) {
this.title = "Lambert Conformal Conic";
}
}
// Lambert Conformal conic forward equations--mapping lat,long to x,y
// -----------------------------------------------------------------
function forward$9(p) {
var lon = p.x;
var lat = p.y;
// singular cases :
if (Math.abs(2 * Math.abs(lat) - Math.PI) <= EPSLN) {
lat = sign(lat) * (HALF_PI - 2 * EPSLN);
}
var con = Math.abs(Math.abs(lat) - HALF_PI);
var ts, rh1;
if (con > EPSLN) {
ts = tsfnz(this.e, lat, Math.sin(lat));
rh1 = this.a * this.f0 * Math.pow(ts, this.ns);
}
else {
con = lat * this.ns;
if (con <= 0) {
return null;
}
rh1 = 0;
}
var theta = this.ns * adjust_lon(lon - this.long0);
p.x = this.k0 * (rh1 * Math.sin(theta)) + this.x0;
p.y = this.k0 * (this.rh - rh1 * Math.cos(theta)) + this.y0;
return p;
}
// Lambert Conformal Conic inverse equations--mapping x,y to lat/long
// -----------------------------------------------------------------
function inverse$9(p) {
var rh1, con, ts;
var lat, lon;
var x = (p.x - this.x0) / this.k0;
var y = (this.rh - (p.y - this.y0) / this.k0);
if (this.ns > 0) {
rh1 = Math.sqrt(x * x + y * y);
con = 1;
}
else {
rh1 = -Math.sqrt(x * x + y * y);
con = -1;
}
var theta = 0;
if (rh1 !== 0) {
theta = Math.atan2((con * x), (con * y));
}
if ((rh1 !== 0) || (this.ns > 0)) {
con = 1 / this.ns;
ts = Math.pow((rh1 / (this.a * this.f0)), con);
lat = phi2z(this.e, ts);
if (lat === -9999) {
return null;
}
}
else {
lat = -HALF_PI;
}
lon = adjust_lon(theta / this.ns + this.long0);
p.x = lon;
p.y = lat;
return p;
}
var names$11 = [
"Lambert Tangential Conformal Conic Projection",
"Lambert_Conformal_Conic",
"Lambert_Conformal_Conic_1SP",
"Lambert_Conformal_Conic_2SP",
"lcc",
"Lambert Conic Conformal (1SP)",
"Lambert Conic Conformal (2SP)"
];
var lcc = {
init: init$10,
forward: forward$9,
inverse: inverse$9,
names: names$11
};
function init$11() {
this.a = 6377397.155;
this.es = 0.006674372230614;
this.e = Math.sqrt(this.es);
if (!this.lat0) {
this.lat0 = 0.863937979737193;
}
if (!this.long0) {
this.long0 = 0.7417649320975901 - 0.308341501185665;
}
/* if scale not set default to 0.9999 */
if (!this.k0) {
this.k0 = 0.9999;
}
this.s45 = 0.785398163397448; /* 45 */
this.s90 = 2 * this.s45;
this.fi0 = this.lat0;
this.e2 = this.es;
this.e = Math.sqrt(this.e2);
this.alfa = Math.sqrt(1 + (this.e2 * Math.pow(Math.cos(this.fi0), 4)) / (1 - this.e2));
this.uq = 1.04216856380474;
this.u0 = Math.asin(Math.sin(this.fi0) / this.alfa);
this.g = Math.pow((1 + this.e * Math.sin(this.fi0)) / (1 - this.e * Math.sin(this.fi0)), this.alfa * this.e / 2);
this.k = Math.tan(this.u0 / 2 + this.s45) / Math.pow(Math.tan(this.fi0 / 2 + this.s45), this.alfa) * this.g;
this.k1 = this.k0;
this.n0 = this.a * Math.sqrt(1 - this.e2) / (1 - this.e2 * Math.pow(Math.sin(this.fi0), 2));
this.s0 = 1.37008346281555;
this.n = Math.sin(this.s0);
this.ro0 = this.k1 * this.n0 / Math.tan(this.s0);
this.ad = this.s90 - this.uq;
}
/* ellipsoid */
/* calculate xy from lat/lon */
/* Constants, identical to inverse transform function */
function forward$10(p) {
var gfi, u, deltav, s, d, eps, ro;
var lon = p.x;
var lat = p.y;
var delta_lon = adjust_lon(lon - this.long0);
/* Transformation */
gfi = Math.pow(((1 + this.e * Math.sin(lat)) / (1 - this.e * Math.sin(lat))), (this.alfa * this.e / 2));
u = 2 * (Math.atan(this.k * Math.pow(Math.tan(lat / 2 + this.s45), this.alfa) / gfi) - this.s45);
deltav = -delta_lon * this.alfa;
s = Math.asin(Math.cos(this.ad) * Math.sin(u) + Math.sin(this.ad) * Math.cos(u) * Math.cos(deltav));
d = Math.asin(Math.cos(u) * Math.sin(deltav) / Math.cos(s));
eps = this.n * d;
ro = this.ro0 * Math.pow(Math.tan(this.s0 / 2 + this.s45), this.n) / Math.pow(Math.tan(s / 2 + this.s45), this.n);
p.y = ro * Math.cos(eps) / 1;
p.x = ro * Math.sin(eps) / 1;
if (!this.czech) {
p.y *= -1;
p.x *= -1;
}
return (p);
}
/* calculate lat/lon from xy */
function inverse$10(p) {
var u, deltav, s, d, eps, ro, fi1;
var ok;
/* Transformation */
/* revert y, x*/
var tmp = p.x;
p.x = p.y;
p.y = tmp;
if (!this.czech) {
p.y *= -1;
p.x *= -1;
}
ro = Math.sqrt(p.x * p.x + p.y * p.y);
eps = Math.atan2(p.y, p.x);
d = eps / Math.sin(this.s0);
s = 2 * (Math.atan(Math.pow(this.ro0 / ro, 1 / this.n) * Math.tan(this.s0 / 2 + this.s45)) - this.s45);
u = Math.asin(Math.cos(this.ad) * Math.sin(s) - Math.sin(this.ad) * Math.cos(s) * Math.cos(d));
deltav = Math.asin(Math.cos(s) * Math.sin(d) / Math.cos(u));
p.x = this.long0 - deltav / this.alfa;
fi1 = u;
ok = 0;
var iter = 0;
do {
p.y = 2 * (Math.atan(Math.pow(this.k, - 1 / this.alfa) * Math.pow(Math.tan(u / 2 + this.s45), 1 / this.alfa) * Math.pow((1 + this.e * Math.sin(fi1)) / (1 - this.e * Math.sin(fi1)), this.e / 2)) - this.s45);
if (Math.abs(fi1 - p.y) < 0.0000000001) {
ok = 1;
}
fi1 = p.y;
iter += 1;
} while (ok === 0 && iter < 15);
if (iter >= 15) {
return null;
}
return (p);
}
var names$12 = ["Krovak", "krovak"];
var krovak = {
init: init$11,
forward: forward$10,
inverse: inverse$10,
names: names$12
};
var mlfn = function(e0, e1, e2, e3, phi) {
return (e0 * phi - e1 * Math.sin(2 * phi) + e2 * Math.sin(4 * phi) - e3 * Math.sin(6 * phi));
};
var e0fn = function(x) {
return (1 - 0.25 * x * (1 + x / 16 * (3 + 1.25 * x)));
};
var e1fn = function(x) {
return (0.375 * x * (1 + 0.25 * x * (1 + 0.46875 * x)));
};
var e2fn = function(x) {
return (0.05859375 * x * x * (1 + 0.75 * x));
};
var e3fn = function(x) {
return (x * x * x * (35 / 3072));
};
var gN = function(a, e, sinphi) {
var temp = e * sinphi;
return a / Math.sqrt(1 - temp * temp);
};
var adjust_lat = function(x) {
return (Math.abs(x) < HALF_PI) ? x : (x - (sign(x) * Math.PI));
};
var imlfn = function(ml, e0, e1, e2, e3) {
var phi;
var dphi;
phi = ml / e0;
for (var i = 0; i < 15; i++) {
dphi = (ml - (e0 * phi - e1 * Math.sin(2 * phi) + e2 * Math.sin(4 * phi) - e3 * Math.sin(6 * phi))) / (e0 - 2 * e1 * Math.cos(2 * phi) + 4 * e2 * Math.cos(4 * phi) - 6 * e3 * Math.cos(6 * phi));
phi += dphi;
if (Math.abs(dphi) <= 0.0000000001) {
return phi;
}
}
//..reportError("IMLFN-CONV:Latitude failed to converge after 15 iterations");
return NaN;
};
function init$12() {
if (!this.sphere) {
this.e0 = e0fn(this.es);
this.e1 = e1fn(this.es);
this.e2 = e2fn(this.es);
this.e3 = e3fn(this.es);
this.ml0 = this.a * mlfn(this.e0, this.e1, this.e2, this.e3, this.lat0);
}
}
/* Cassini forward equations--mapping lat,long to x,y
-----------------------------------------------------------------------*/
function forward$11(p) {
/* Forward equations
-----------------*/
var x, y;
var lam = p.x;
var phi = p.y;
lam = adjust_lon(lam - this.long0);
if (this.sphere) {
x = this.a * Math.asin(Math.cos(phi) * Math.sin(lam));
y = this.a * (Math.atan2(Math.tan(phi), Math.cos(lam)) - this.lat0);
}
else {
//ellipsoid
var sinphi = Math.sin(phi);
var cosphi = Math.cos(phi);
var nl = gN(this.a, this.e, sinphi);
var tl = Math.tan(phi) * Math.tan(phi);
var al = lam * Math.cos(phi);
var asq = al * al;
var cl = this.es * cosphi * cosphi / (1 - this.es);
var ml = this.a * mlfn(this.e0, this.e1, this.e2, this.e3, phi);
x = nl * al * (1 - asq * tl * (1 / 6 - (8 - tl + 8 * cl) * asq / 120));
y = ml - this.ml0 + nl * sinphi / cosphi * asq * (0.5 + (5 - tl + 6 * cl) * asq / 24);
}
p.x = x + this.x0;
p.y = y + this.y0;
return p;
}
/* Inverse equations
-----------------*/
function inverse$11(p) {
p.x -= this.x0;
p.y -= this.y0;
var x = p.x / this.a;
var y = p.y / this.a;
var phi, lam;
if (this.sphere) {
var dd = y + this.lat0;
phi = Math.asin(Math.sin(dd) * Math.cos(x));
lam = Math.atan2(Math.tan(x), Math.cos(dd));
}
else {
/* ellipsoid */
var ml1 = this.ml0 / this.a + y;
var phi1 = imlfn(ml1, this.e0, this.e1, this.e2, this.e3);
if (Math.abs(Math.abs(phi1) - HALF_PI) <= EPSLN) {
p.x = this.long0;
p.y = HALF_PI;
if (y < 0) {
p.y *= -1;
}
return p;
}
var nl1 = gN(this.a, this.e, Math.sin(phi1));
var rl1 = nl1 * nl1 * nl1 / this.a / this.a * (1 - this.es);
var tl1 = Math.pow(Math.tan(phi1), 2);
var dl = x * this.a / nl1;
var dsq = dl * dl;
phi = phi1 - nl1 * Math.tan(phi1) / rl1 * dl * dl * (0.5 - (1 + 3 * tl1) * dl * dl / 24);
lam = dl * (1 - dsq * (tl1 / 3 + (1 + 3 * tl1) * tl1 * dsq / 15)) / Math.cos(phi1);
}
p.x = adjust_lon(lam + this.long0);
p.y = adjust_lat(phi);
return p;
}
var names$13 = ["Cassini", "Cassini_Soldner", "cass"];
var cass = {
init: init$12,
forward: forward$11,
inverse: inverse$11,
names: names$13
};
var qsfnz = function(eccent, sinphi) {
var con;
if (eccent > 1.0e-7) {
con = eccent * sinphi;
return ((1 - eccent * eccent) * (sinphi / (1 - con * con) - (0.5 / eccent) * Math.log((1 - con) / (1 + con))));
}
else {
return (2 * sinphi);
}
};
/*
reference
"New Equal-Area Map Projections for Noncircular Regions", John P. Snyder,
The American Cartographer, Vol 15, No. 4, October 1988, pp. 341-355.
*/
var S_POLE = 1;
var N_POLE = 2;
var EQUIT = 3;
var OBLIQ = 4;
/* Initialize the Lambert Azimuthal Equal Area projection
------------------------------------------------------*/
function init$13() {
var t = Math.abs(this.lat0);
if (Math.abs(t - HALF_PI) < EPSLN) {
this.mode = this.lat0 < 0 ? this.S_POLE : this.N_POLE;
}
else if (Math.abs(t) < EPSLN) {
this.mode = this.EQUIT;
}
else {
this.mode = this.OBLIQ;
}
if (this.es > 0) {
var sinphi;
this.qp = qsfnz(this.e, 1);
this.mmf = 0.5 / (1 - this.es);
this.apa = authset(this.es);
switch (this.mode) {
case this.N_POLE:
this.dd = 1;
break;
case this.S_POLE:
this.dd = 1;
break;
case this.EQUIT:
this.rq = Math.sqrt(0.5 * this.qp);
this.dd = 1 / this.rq;
this.xmf = 1;
this.ymf = 0.5 * this.qp;
break;
case this.OBLIQ:
this.rq = Math.sqrt(0.5 * this.qp);
sinphi = Math.sin(this.lat0);
this.sinb1 = qsfnz(this.e, sinphi) / this.qp;
this.cosb1 = Math.sqrt(1 - this.sinb1 * this.sinb1);
this.dd = Math.cos(this.lat0) / (Math.sqrt(1 - this.es * sinphi * sinphi) * this.rq * this.cosb1);
this.ymf = (this.xmf = this.rq) / this.dd;
this.xmf *= this.dd;
break;
}
}
else {
if (this.mode === this.OBLIQ) {
this.sinph0 = Math.sin(this.lat0);
this.cosph0 = Math.cos(this.lat0);
}
}
}
/* Lambert Azimuthal Equal Area forward equations--mapping lat,long to x,y
-----------------------------------------------------------------------*/
function forward$12(p) {
/* Forward equations
-----------------*/
var x, y, coslam, sinlam, sinphi, q, sinb, cosb, b, cosphi;
var lam = p.x;
var phi = p.y;
lam = adjust_lon(lam - this.long0);
if (this.sphere) {
sinphi = Math.sin(phi);
cosphi = Math.cos(phi);
coslam = Math.cos(lam);
if (this.mode === this.OBLIQ || this.mode === this.EQUIT) {
y = (this.mode === this.EQUIT) ? 1 + cosphi * coslam : 1 + this.sinph0 * sinphi + this.cosph0 * cosphi * coslam;
if (y <= EPSLN) {
return null;
}
y = Math.sqrt(2 / y);
x = y * cosphi * Math.sin(lam);
y *= (this.mode === this.EQUIT) ? sinphi : this.cosph0 * sinphi - this.sinph0 * cosphi * coslam;
}
else if (this.mode === this.N_POLE || this.mode === this.S_POLE) {
if (this.mode === this.N_POLE) {
coslam = -coslam;
}
if (Math.abs(phi + this.lat0) < EPSLN) {
return null;
}
y = FORTPI - phi * 0.5;
y = 2 * ((this.mode === this.S_POLE) ? Math.cos(y) : Math.sin(y));
x = y * Math.sin(lam);
y *= coslam;
}
}
else {
sinb = 0;
cosb = 0;
b = 0;
coslam = Math.cos(lam);
sinlam = Math.sin(lam);
sinphi = Math.sin(phi);
q = qsfnz(this.e, sinphi);
if (this.mode === this.OBLIQ || this.mode === this.EQUIT) {
sinb = q / this.qp;
cosb = Math.sqrt(1 - sinb * sinb);
}
switch (this.mode) {
case this.OBLIQ:
b = 1 + this.sinb1 * sinb + this.cosb1 * cosb * coslam;
break;
case this.EQUIT:
b = 1 + cosb * coslam;
break;
case this.N_POLE:
b = HALF_PI + phi;
q = this.qp - q;
break;
case this.S_POLE:
b = phi - HALF_PI;
q = this.qp + q;
break;
}
if (Math.abs(b) < EPSLN) {
return null;
}
switch (this.mode) {
case this.OBLIQ:
case this.EQUIT:
b = Math.sqrt(2 / b);
if (this.mode === this.OBLIQ) {
y = this.ymf * b * (this.cosb1 * sinb - this.sinb1 * cosb * coslam);
}
else {
y = (b = Math.sqrt(2 / (1 + cosb * coslam))) * sinb * this.ymf;
}
x = this.xmf * b * cosb * sinlam;
break;
case this.N_POLE:
case this.S_POLE:
if (q >= 0) {
x = (b = Math.sqrt(q)) * sinlam;
y = coslam * ((this.mode === this.S_POLE) ? b : -b);
}
else {
x = y = 0;
}
break;
}
}
p.x = this.a * x + this.x0;
p.y = this.a * y + this.y0;
return p;
}
/* Inverse equations
-----------------*/
function inverse$12(p) {
p.x -= this.x0;
p.y -= this.y0;
var x = p.x / this.a;
var y = p.y / this.a;
var lam, phi, cCe, sCe, q, rho, ab;
if (this.sphere) {
var cosz = 0,
rh, sinz = 0;
rh = Math.sqrt(x * x + y * y);
phi = rh * 0.5;
if (phi > 1) {
return null;
}
phi = 2 * Math.asin(phi);
if (this.mode === this.OBLIQ || this.mode === this.EQUIT) {
sinz = Math.sin(phi);
cosz = Math.cos(phi);
}
switch (this.mode) {
case this.EQUIT:
phi = (Math.abs(rh) <= EPSLN) ? 0 : Math.asin(y * sinz / rh);
x *= sinz;
y = cosz * rh;
break;
case this.OBLIQ:
phi = (Math.abs(rh) <= EPSLN) ? this.lat0 : Math.asin(cosz * this.sinph0 + y * sinz * this.cosph0 / rh);
x *= sinz * this.cosph0;
y = (cosz - Math.sin(phi) * this.sinph0) * rh;
break;
case this.N_POLE:
y = -y;
phi = HALF_PI - phi;
break;
case this.S_POLE:
phi -= HALF_PI;
break;
}
lam = (y === 0 && (this.mode === this.EQUIT || this.mode === this.OBLIQ)) ? 0 : Math.atan2(x, y);
}
else {
ab = 0;
if (this.mode === this.OBLIQ || this.mode === this.EQUIT) {
x /= this.dd;
y *= this.dd;
rho = Math.sqrt(x * x + y * y);
if (rho < EPSLN) {
p.x = this.long0;
p.y = this.lat0;
return p;
}
sCe = 2 * Math.asin(0.5 * rho / this.rq);
cCe = Math.cos(sCe);
x *= (sCe = Math.sin(sCe));
if (this.mode === this.OBLIQ) {
ab = cCe * this.sinb1 + y * sCe * this.cosb1 / rho;
q = this.qp * ab;
y = rho * this.cosb1 * cCe - y * this.sinb1 * sCe;
}
else {
ab = y * sCe / rho;
q = this.qp * ab;
y = rho * cCe;
}
}
else if (this.mode === this.N_POLE || this.mode === this.S_POLE) {
if (this.mode === this.N_POLE) {
y = -y;
}
q = (x * x + y * y);
if (!q) {
p.x = this.long0;
p.y = this.lat0;
return p;
}
ab = 1 - q / this.qp;
if (this.mode === this.S_POLE) {
ab = -ab;
}
}
lam = Math.atan2(x, y);
phi = authlat(Math.asin(ab), this.apa);
}
p.x = adjust_lon(this.long0 + lam);
p.y = phi;
return p;
}
/* determine latitude from authalic latitude */
var P00 = 0.33333333333333333333;
var P01 = 0.17222222222222222222;
var P02 = 0.10257936507936507936;
var P10 = 0.06388888888888888888;
var P11 = 0.06640211640211640211;
var P20 = 0.01641501294219154443;
function authset(es) {
var t;
var APA = [];
APA[0] = es * P00;
t = es * es;
APA[0] += t * P01;
APA[1] = t * P10;
t *= es;
APA[0] += t * P02;
APA[1] += t * P11;
APA[2] = t * P20;
return APA;
}
function authlat(beta, APA) {
var t = beta + beta;
return (beta + APA[0] * Math.sin(t) + APA[1] * Math.sin(t + t) + APA[2] * Math.sin(t + t + t));
}
var names$14 = ["Lambert Azimuthal Equal Area", "Lambert_Azimuthal_Equal_Area", "laea"];
var laea = {
init: init$13,
forward: forward$12,
inverse: inverse$12,
names: names$14,
S_POLE: S_POLE,
N_POLE: N_POLE,
EQUIT: EQUIT,
OBLIQ: OBLIQ
};
var asinz = function(x) {
if (Math.abs(x) > 1) {
x = (x > 1) ? 1 : -1;
}
return Math.asin(x);
};
function init$14() {
if (Math.abs(this.lat1 + this.lat2) < EPSLN) {
return;
}
this.temp = this.b / this.a;
this.es = 1 - Math.pow(this.temp, 2);
this.e3 = Math.sqrt(this.es);
this.sin_po = Math.sin(this.lat1);
this.cos_po = Math.cos(this.lat1);
this.t1 = this.sin_po;
this.con = this.sin_po;
this.ms1 = msfnz(this.e3, this.sin_po, this.cos_po);
this.qs1 = qsfnz(this.e3, this.sin_po);
this.sin_po = Math.sin(this.lat2);
this.cos_po = Math.cos(this.lat2);
this.t2 = this.sin_po;
this.ms2 = msfnz(this.e3, this.sin_po, this.cos_po);
this.qs2 = qsfnz(this.e3, this.sin_po);
this.sin_po = Math.sin(this.lat0);
this.cos_po = Math.cos(this.lat0);
this.t3 = this.sin_po;
this.qs0 = qsfnz(this.e3, this.sin_po);
if (Math.abs(this.lat1 - this.lat2) > EPSLN) {
this.ns0 = (this.ms1 * this.ms1 - this.ms2 * this.ms2) / (this.qs2 - this.qs1);
}
else {
this.ns0 = this.con;
}
this.c = this.ms1 * this.ms1 + this.ns0 * this.qs1;
this.rh = this.a * Math.sqrt(this.c - this.ns0 * this.qs0) / this.ns0;
}
/* Albers Conical Equal Area forward equations--mapping lat,long to x,y
-------------------------------------------------------------------*/
function forward$13(p) {
var lon = p.x;
var lat = p.y;
this.sin_phi = Math.sin(lat);
this.cos_phi = Math.cos(lat);
var qs = qsfnz(this.e3, this.sin_phi);
var rh1 = this.a * Math.sqrt(this.c - this.ns0 * qs) / this.ns0;
var theta = this.ns0 * adjust_lon(lon - this.long0);
var x = rh1 * Math.sin(theta) + this.x0;
var y = this.rh - rh1 * Math.cos(theta) + this.y0;
p.x = x;
p.y = y;
return p;
}
function inverse$13(p) {
var rh1, qs, con, theta, lon, lat;
p.x -= this.x0;
p.y = this.rh - p.y + this.y0;
if (this.ns0 >= 0) {
rh1 = Math.sqrt(p.x * p.x + p.y * p.y);
con = 1;
}
else {
rh1 = -Math.sqrt(p.x * p.x + p.y * p.y);
con = -1;
}
theta = 0;
if (rh1 !== 0) {
theta = Math.atan2(con * p.x, con * p.y);
}
con = rh1 * this.ns0 / this.a;
if (this.sphere) {
lat = Math.asin((this.c - con * con) / (2 * this.ns0));
}
else {
qs = (this.c - con * con) / this.ns0;
lat = this.phi1z(this.e3, qs);
}
lon = adjust_lon(theta / this.ns0 + this.long0);
p.x = lon;
p.y = lat;
return p;
}
/* Function to compute phi1, the latitude for the inverse of the
Albers Conical Equal-Area projection.
-------------------------------------------*/
function phi1z(eccent, qs) {
var sinphi, cosphi, con, com, dphi;
var phi = asinz(0.5 * qs);
if (eccent < EPSLN) {
return phi;
}
var eccnts = eccent * eccent;
for (var i = 1; i <= 25; i++) {
sinphi = Math.sin(phi);
cosphi = Math.cos(phi);
con = eccent * sinphi;
com = 1 - con * con;
dphi = 0.5 * com * com / cosphi * (qs / (1 - eccnts) - sinphi / com + 0.5 / eccent * Math.log((1 - con) / (1 + con)));
phi = phi + dphi;
if (Math.abs(dphi) <= 1e-7) {
return phi;
}
}
return null;
}
var names$15 = ["Albers_Conic_Equal_Area", "Albers", "aea"];
var aea = {
init: init$14,
forward: forward$13,
inverse: inverse$13,
names: names$15,
phi1z: phi1z
};
/*
reference:
Wolfram Mathworld "Gnomonic Projection"
http://mathworld.wolfram.com/GnomonicProjection.html
Accessed: 12th November 2009
*/
function init$15() {
/* Place parameters in static storage for common use
-------------------------------------------------*/
this.sin_p14 = Math.sin(this.lat0);
this.cos_p14 = Math.cos(this.lat0);
// Approximation for projecting points to the horizon (infinity)
this.infinity_dist = 1000 * this.a;
this.rc = 1;
}
/* Gnomonic forward equations--mapping lat,long to x,y
---------------------------------------------------*/
function forward$14(p) {
var sinphi, cosphi; /* sin and cos value */
var dlon; /* delta longitude value */
var coslon; /* cos of longitude */
var ksp; /* scale factor */
var g;
var x, y;
var lon = p.x;
var lat = p.y;
/* Forward equations
-----------------*/
dlon = adjust_lon(lon - this.long0);
sinphi = Math.sin(lat);
cosphi = Math.cos(lat);
coslon = Math.cos(dlon);
g = this.sin_p14 * sinphi + this.cos_p14 * cosphi * coslon;
ksp = 1;
if ((g > 0) || (Math.abs(g) <= EPSLN)) {
x = this.x0 + this.a * ksp * cosphi * Math.sin(dlon) / g;
y = this.y0 + this.a * ksp * (this.cos_p14 * sinphi - this.sin_p14 * cosphi * coslon) / g;
}
else {
// Point is in the opposing hemisphere and is unprojectable
// We still need to return a reasonable point, so we project
// to infinity, on a bearing
// equivalent to the northern hemisphere equivalent
// This is a reasonable approximation for short shapes and lines that
// straddle the horizon.
x = this.x0 + this.infinity_dist * cosphi * Math.sin(dlon);
y = this.y0 + this.infinity_dist * (this.cos_p14 * sinphi - this.sin_p14 * cosphi * coslon);
}
p.x = x;
p.y = y;
return p;
}
function inverse$14(p) {
var rh; /* Rho */
var sinc, cosc;
var c;
var lon, lat;
/* Inverse equations
-----------------*/
p.x = (p.x - this.x0) / this.a;
p.y = (p.y - this.y0) / this.a;
p.x /= this.k0;
p.y /= this.k0;
if ((rh = Math.sqrt(p.x * p.x + p.y * p.y))) {
c = Math.atan2(rh, this.rc);
sinc = Math.sin(c);
cosc = Math.cos(c);
lat = asinz(cosc * this.sin_p14 + (p.y * sinc * this.cos_p14) / rh);
lon = Math.atan2(p.x * sinc, rh * this.cos_p14 * cosc - p.y * this.sin_p14 * sinc);
lon = adjust_lon(this.long0 + lon);
}
else {
lat = this.phic0;
lon = 0;
}
p.x = lon;
p.y = lat;
return p;
}
var names$16 = ["gnom"];
var gnom = {
init: init$15,
forward: forward$14,
inverse: inverse$14,
names: names$16
};
var iqsfnz = function(eccent, q) {
var temp = 1 - (1 - eccent * eccent) / (2 * eccent) * Math.log((1 - eccent) / (1 + eccent));
if (Math.abs(Math.abs(q) - temp) < 1.0E-6) {
if (q < 0) {
return (-1 * HALF_PI);
}
else {
return HALF_PI;
}
}
//var phi = 0.5* q/(1-eccent*eccent);
var phi = Math.asin(0.5 * q);
var dphi;
var sin_phi;
var cos_phi;
var con;
for (var i = 0; i < 30; i++) {
sin_phi = Math.sin(phi);
cos_phi = Math.cos(phi);
con = eccent * sin_phi;
dphi = Math.pow(1 - con * con, 2) / (2 * cos_phi) * (q / (1 - eccent * eccent) - sin_phi / (1 - con * con) + 0.5 / eccent * Math.log((1 - con) / (1 + con)));
phi += dphi;
if (Math.abs(dphi) <= 0.0000000001) {
return phi;
}
}
//console.log("IQSFN-CONV:Latitude failed to converge after 30 iterations");
return NaN;
};
/*
reference:
"Cartographic Projection Procedures for the UNIX Environment-
A User's Manual" by Gerald I. Evenden,
USGS Open File Report 90-284and Release 4 Interim Reports (2003)
*/
function init$16() {
//no-op
if (!this.sphere) {
this.k0 = msfnz(this.e, Math.sin(this.lat_ts), Math.cos(this.lat_ts));
}
}
/* Cylindrical Equal Area forward equations--mapping lat,long to x,y
------------------------------------------------------------*/
function forward$15(p) {
var lon = p.x;
var lat = p.y;
var x, y;
/* Forward equations
-----------------*/
var dlon = adjust_lon(lon - this.long0);
if (this.sphere) {
x = this.x0 + this.a * dlon * Math.cos(this.lat_ts);
y = this.y0 + this.a * Math.sin(lat) / Math.cos(this.lat_ts);
}
else {
var qs = qsfnz(this.e, Math.sin(lat));
x = this.x0 + this.a * this.k0 * dlon;
y = this.y0 + this.a * qs * 0.5 / this.k0;
}
p.x = x;
p.y = y;
return p;
}
/* Cylindrical Equal Area inverse equations--mapping x,y to lat/long
------------------------------------------------------------*/
function inverse$15(p) {
p.x -= this.x0;
p.y -= this.y0;
var lon, lat;
if (this.sphere) {
lon = adjust_lon(this.long0 + (p.x / this.a) / Math.cos(this.lat_ts));
lat = Math.asin((p.y / this.a) * Math.cos(this.lat_ts));
}
else {
lat = iqsfnz(this.e, 2 * p.y * this.k0 / this.a);
lon = adjust_lon(this.long0 + p.x / (this.a * this.k0));
}
p.x = lon;
p.y = lat;
return p;
}
var names$17 = ["cea"];
var cea = {
init: init$16,
forward: forward$15,
inverse: inverse$15,
names: names$17
};
function init$17() {
this.x0 = this.x0 || 0;
this.y0 = this.y0 || 0;
this.lat0 = this.lat0 || 0;
this.long0 = this.long0 || 0;
this.lat_ts = this.lat_ts || 0;
this.title = this.title || "Equidistant Cylindrical (Plate Carre)";
this.rc = Math.cos(this.lat_ts);
}
// forward equations--mapping lat,long to x,y
// -----------------------------------------------------------------
function forward$16(p) {
var lon = p.x;
var lat = p.y;
var dlon = adjust_lon(lon - this.long0);
var dlat = adjust_lat(lat - this.lat0);
p.x = this.x0 + (this.a * dlon * this.rc);
p.y = this.y0 + (this.a * dlat);
return p;
}
// inverse equations--mapping x,y to lat/long
// -----------------------------------------------------------------
function inverse$16(p) {
var x = p.x;
var y = p.y;
p.x = adjust_lon(this.long0 + ((x - this.x0) / (this.a * this.rc)));
p.y = adjust_lat(this.lat0 + ((y - this.y0) / (this.a)));
return p;
}
var names$18 = ["Equirectangular", "Equidistant_Cylindrical", "eqc"];
var eqc = {
init: init$17,
forward: forward$16,
inverse: inverse$16,
names: names$18
};
var MAX_ITER$2 = 20;
function init$18() {
/* Place parameters in static storage for common use
-------------------------------------------------*/
this.temp = this.b / this.a;
this.es = 1 - Math.pow(this.temp, 2); // devait etre dans tmerc.js mais n y est pas donc je commente sinon retour de valeurs nulles
this.e = Math.sqrt(this.es);
this.e0 = e0fn(this.es);
this.e1 = e1fn(this.es);
this.e2 = e2fn(this.es);
this.e3 = e3fn(this.es);
this.ml0 = this.a * mlfn(this.e0, this.e1, this.e2, this.e3, this.lat0); //si que des zeros le calcul ne se fait pas
}
/* Polyconic forward equations--mapping lat,long to x,y
---------------------------------------------------*/
function forward$17(p) {
var lon = p.x;
var lat = p.y;
var x, y, el;
var dlon = adjust_lon(lon - this.long0);
el = dlon * Math.sin(lat);
if (this.sphere) {
if (Math.abs(lat) <= EPSLN) {
x = this.a * dlon;
y = -1 * this.a * this.lat0;
}
else {
x = this.a * Math.sin(el) / Math.tan(lat);
y = this.a * (adjust_lat(lat - this.lat0) + (1 - Math.cos(el)) / Math.tan(lat));
}
}
else {
if (Math.abs(lat) <= EPSLN) {
x = this.a * dlon;
y = -1 * this.ml0;
}
else {
var nl = gN(this.a, this.e, Math.sin(lat)) / Math.tan(lat);
x = nl * Math.sin(el);
y = this.a * mlfn(this.e0, this.e1, this.e2, this.e3, lat) - this.ml0 + nl * (1 - Math.cos(el));
}
}
p.x = x + this.x0;
p.y = y + this.y0;
return p;
}
/* Inverse equations
-----------------*/
function inverse$17(p) {
var lon, lat, x, y, i;
var al, bl;
var phi, dphi;
x = p.x - this.x0;
y = p.y - this.y0;
if (this.sphere) {
if (Math.abs(y + this.a * this.lat0) <= EPSLN) {
lon = adjust_lon(x / this.a + this.long0);
lat = 0;
}
else {
al = this.lat0 + y / this.a;
bl = x * x / this.a / this.a + al * al;
phi = al;
var tanphi;
for (i = MAX_ITER$2; i; --i) {
tanphi = Math.tan(phi);
dphi = -1 * (al * (phi * tanphi + 1) - phi - 0.5 * (phi * phi + bl) * tanphi) / ((phi - al) / tanphi - 1);
phi += dphi;
if (Math.abs(dphi) <= EPSLN) {
lat = phi;
break;
}
}
lon = adjust_lon(this.long0 + (Math.asin(x * Math.tan(phi) / this.a)) / Math.sin(lat));
}
}
else {
if (Math.abs(y + this.ml0) <= EPSLN) {
lat = 0;
lon = adjust_lon(this.long0 + x / this.a);
}
else {
al = (this.ml0 + y) / this.a;
bl = x * x / this.a / this.a + al * al;
phi = al;
var cl, mln, mlnp, ma;
var con;
for (i = MAX_ITER$2; i; --i) {
con = this.e * Math.sin(phi);
cl = Math.sqrt(1 - con * con) * Math.tan(phi);
mln = this.a * mlfn(this.e0, this.e1, this.e2, this.e3, phi);
mlnp = this.e0 - 2 * this.e1 * Math.cos(2 * phi) + 4 * this.e2 * Math.cos(4 * phi) - 6 * this.e3 * Math.cos(6 * phi);
ma = mln / this.a;
dphi = (al * (cl * ma + 1) - ma - 0.5 * cl * (ma * ma + bl)) / (this.es * Math.sin(2 * phi) * (ma * ma + bl - 2 * al * ma) / (4 * cl) + (al - ma) * (cl * mlnp - 2 / Math.sin(2 * phi)) - mlnp);
phi -= dphi;
if (Math.abs(dphi) <= EPSLN) {
lat = phi;
break;
}
}
//lat=phi4z(this.e,this.e0,this.e1,this.e2,this.e3,al,bl,0,0);
cl = Math.sqrt(1 - this.es * Math.pow(Math.sin(lat), 2)) * Math.tan(lat);
lon = adjust_lon(this.long0 + Math.asin(x * cl / this.a) / Math.sin(lat));
}
}
p.x = lon;
p.y = lat;
return p;
}
var names$19 = ["Polyconic", "poly"];
var poly = {
init: init$18,
forward: forward$17,
inverse: inverse$17,
names: names$19
};
/*
reference
Department of Land and Survey Technical Circular 1973/32
http://www.linz.govt.nz/docs/miscellaneous/nz-map-definition.pdf
OSG Technical Report 4.1
http://www.linz.govt.nz/docs/miscellaneous/nzmg.pdf
*/
/**
* iterations: Number of iterations to refine inverse transform.
* 0 -> km accuracy
* 1 -> m accuracy -- suitable for most mapping applications
* 2 -> mm accuracy
*/
function init$19() {
this.A = [];
this.A[1] = 0.6399175073;
this.A[2] = -0.1358797613;
this.A[3] = 0.063294409;
this.A[4] = -0.02526853;
this.A[5] = 0.0117879;
this.A[6] = -0.0055161;
this.A[7] = 0.0026906;
this.A[8] = -0.001333;
this.A[9] = 0.00067;
this.A[10] = -0.00034;
this.B_re = [];
this.B_im = [];
this.B_re[1] = 0.7557853228;
this.B_im[1] = 0;
this.B_re[2] = 0.249204646;
this.B_im[2] = 0.003371507;
this.B_re[3] = -0.001541739;
this.B_im[3] = 0.041058560;
this.B_re[4] = -0.10162907;
this.B_im[4] = 0.01727609;
this.B_re[5] = -0.26623489;
this.B_im[5] = -0.36249218;
this.B_re[6] = -0.6870983;
this.B_im[6] = -1.1651967;
this.C_re = [];
this.C_im = [];
this.C_re[1] = 1.3231270439;
this.C_im[1] = 0;
this.C_re[2] = -0.577245789;
this.C_im[2] = -0.007809598;
this.C_re[3] = 0.508307513;
this.C_im[3] = -0.112208952;
this.C_re[4] = -0.15094762;
this.C_im[4] = 0.18200602;
this.C_re[5] = 1.01418179;
this.C_im[5] = 1.64497696;
this.C_re[6] = 1.9660549;
this.C_im[6] = 2.5127645;
this.D = [];
this.D[1] = 1.5627014243;
this.D[2] = 0.5185406398;
this.D[3] = -0.03333098;
this.D[4] = -0.1052906;
this.D[5] = -0.0368594;
this.D[6] = 0.007317;
this.D[7] = 0.01220;
this.D[8] = 0.00394;
this.D[9] = -0.0013;
}
/**
New Zealand Map Grid Forward - long/lat to x/y
long/lat in radians
*/
function forward$18(p) {
var n;
var lon = p.x;
var lat = p.y;
var delta_lat = lat - this.lat0;
var delta_lon = lon - this.long0;
// 1. Calculate d_phi and d_psi ... // and d_lambda
// For this algorithm, delta_latitude is in seconds of arc x 10-5, so we need to scale to those units. Longitude is radians.
var d_phi = delta_lat / SEC_TO_RAD * 1E-5;
var d_lambda = delta_lon;
var d_phi_n = 1; // d_phi^0
var d_psi = 0;
for (n = 1; n <= 10; n++) {
d_phi_n = d_phi_n * d_phi;
d_psi = d_psi + this.A[n] * d_phi_n;
}
// 2. Calculate theta
var th_re = d_psi;
var th_im = d_lambda;
// 3. Calculate z
var th_n_re = 1;
var th_n_im = 0; // theta^0
var th_n_re1;
var th_n_im1;
var z_re = 0;
var z_im = 0;
for (n = 1; n <= 6; n++) {
th_n_re1 = th_n_re * th_re - th_n_im * th_im;
th_n_im1 = th_n_im * th_re + th_n_re * th_im;
th_n_re = th_n_re1;
th_n_im = th_n_im1;
z_re = z_re + this.B_re[n] * th_n_re - this.B_im[n] * th_n_im;
z_im = z_im + this.B_im[n] * th_n_re + this.B_re[n] * th_n_im;
}
// 4. Calculate easting and northing
p.x = (z_im * this.a) + this.x0;
p.y = (z_re * this.a) + this.y0;
return p;
}
/**
New Zealand Map Grid Inverse - x/y to long/lat
*/
function inverse$18(p) {
var n;
var x = p.x;
var y = p.y;
var delta_x = x - this.x0;
var delta_y = y - this.y0;
// 1. Calculate z
var z_re = delta_y / this.a;
var z_im = delta_x / this.a;
// 2a. Calculate theta - first approximation gives km accuracy
var z_n_re = 1;
var z_n_im = 0; // z^0
var z_n_re1;
var z_n_im1;
var th_re = 0;
var th_im = 0;
for (n = 1; n <= 6; n++) {
z_n_re1 = z_n_re * z_re - z_n_im * z_im;
z_n_im1 = z_n_im * z_re + z_n_re * z_im;
z_n_re = z_n_re1;
z_n_im = z_n_im1;
th_re = th_re + this.C_re[n] * z_n_re - this.C_im[n] * z_n_im;
th_im = th_im + this.C_im[n] * z_n_re + this.C_re[n] * z_n_im;
}
// 2b. Iterate to refine the accuracy of the calculation
// 0 iterations gives km accuracy
// 1 iteration gives m accuracy -- good enough for most mapping applications
// 2 iterations bives mm accuracy
for (var i = 0; i < this.iterations; i++) {
var th_n_re = th_re;
var th_n_im = th_im;
var th_n_re1;
var th_n_im1;
var num_re = z_re;
var num_im = z_im;
for (n = 2; n <= 6; n++) {
th_n_re1 = th_n_re * th_re - th_n_im * th_im;
th_n_im1 = th_n_im * th_re + th_n_re * th_im;
th_n_re = th_n_re1;
th_n_im = th_n_im1;
num_re = num_re + (n - 1) * (this.B_re[n] * th_n_re - this.B_im[n] * th_n_im);
num_im = num_im + (n - 1) * (this.B_im[n] * th_n_re + this.B_re[n] * th_n_im);
}
th_n_re = 1;
th_n_im = 0;
var den_re = this.B_re[1];
var den_im = this.B_im[1];
for (n = 2; n <= 6; n++) {
th_n_re1 = th_n_re * th_re - th_n_im * th_im;
th_n_im1 = th_n_im * th_re + th_n_re * th_im;
th_n_re = th_n_re1;
th_n_im = th_n_im1;
den_re = den_re + n * (this.B_re[n] * th_n_re - this.B_im[n] * th_n_im);
den_im = den_im + n * (this.B_im[n] * th_n_re + this.B_re[n] * th_n_im);
}
// Complex division
var den2 = den_re * den_re + den_im * den_im;
th_re = (num_re * den_re + num_im * den_im) / den2;
th_im = (num_im * den_re - num_re * den_im) / den2;
}
// 3. Calculate d_phi ... // and d_lambda
var d_psi = th_re;
var d_lambda = th_im;
var d_psi_n = 1; // d_psi^0
var d_phi = 0;
for (n = 1; n <= 9; n++) {
d_psi_n = d_psi_n * d_psi;
d_phi = d_phi + this.D[n] * d_psi_n;
}
// 4. Calculate latitude and longitude
// d_phi is calcuated in second of arc * 10^-5, so we need to scale back to radians. d_lambda is in radians.
var lat = this.lat0 + (d_phi * SEC_TO_RAD * 1E5);
var lon = this.long0 + d_lambda;
p.x = lon;
p.y = lat;
return p;
}
var names$20 = ["New_Zealand_Map_Grid", "nzmg"];
var nzmg = {
init: init$19,
forward: forward$18,
inverse: inverse$18,
names: names$20
};
/*
reference
"New Equal-Area Map Projections for Noncircular Regions", John P. Snyder,
The American Cartographer, Vol 15, No. 4, October 1988, pp. 341-355.
*/
/* Initialize the Miller Cylindrical projection
-------------------------------------------*/
function init$20() {
//no-op
}
/* Miller Cylindrical forward equations--mapping lat,long to x,y
------------------------------------------------------------*/
function forward$19(p) {
var lon = p.x;
var lat = p.y;
/* Forward equations
-----------------*/
var dlon = adjust_lon(lon - this.long0);
var x = this.x0 + this.a * dlon;
var y = this.y0 + this.a * Math.log(Math.tan((Math.PI / 4) + (lat / 2.5))) * 1.25;
p.x = x;
p.y = y;
return p;
}
/* Miller Cylindrical inverse equations--mapping x,y to lat/long
------------------------------------------------------------*/
function inverse$19(p) {
p.x -= this.x0;
p.y -= this.y0;
var lon = adjust_lon(this.long0 + p.x / this.a);
var lat = 2.5 * (Math.atan(Math.exp(0.8 * p.y / this.a)) - Math.PI / 4);
p.x = lon;
p.y = lat;
return p;
}
var names$21 = ["Miller_Cylindrical", "mill"];
var mill = {
init: init$20,
forward: forward$19,
inverse: inverse$19,
names: names$21
};
var MAX_ITER$3 = 20;
function init$21() {
/* Place parameters in static storage for common use
-------------------------------------------------*/
if (!this.sphere) {
this.en = pj_enfn(this.es);
}
else {
this.n = 1;
this.m = 0;
this.es = 0;
this.C_y = Math.sqrt((this.m + 1) / this.n);
this.C_x = this.C_y / (this.m + 1);
}
}
/* Sinusoidal forward equations--mapping lat,long to x,y
-----------------------------------------------------*/
function forward$20(p) {
var x, y;
var lon = p.x;
var lat = p.y;
/* Forward equations
-----------------*/
lon = adjust_lon(lon - this.long0);
if (this.sphere) {
if (!this.m) {
lat = this.n !== 1 ? Math.asin(this.n * Math.sin(lat)) : lat;
}
else {
var k = this.n * Math.sin(lat);
for (var i = MAX_ITER$3; i; --i) {
var V = (this.m * lat + Math.sin(lat) - k) / (this.m + Math.cos(lat));
lat -= V;
if (Math.abs(V) < EPSLN) {
break;
}
}
}
x = this.a * this.C_x * lon * (this.m + Math.cos(lat));
y = this.a * this.C_y * lat;
}
else {
var s = Math.sin(lat);
var c = Math.cos(lat);
y = this.a * pj_mlfn(lat, s, c, this.en);
x = this.a * lon * c / Math.sqrt(1 - this.es * s * s);
}
p.x = x;
p.y = y;
return p;
}
function inverse$20(p) {
var lat, temp, lon, s;
p.x -= this.x0;
lon = p.x / this.a;
p.y -= this.y0;
lat = p.y / this.a;
if (this.sphere) {
lat /= this.C_y;
lon = lon / (this.C_x * (this.m + Math.cos(lat)));
if (this.m) {
lat = asinz((this.m * lat + Math.sin(lat)) / this.n);
}
else if (this.n !== 1) {
lat = asinz(Math.sin(lat) / this.n);
}
lon = adjust_lon(lon + this.long0);
lat = adjust_lat(lat);
}
else {
lat = pj_inv_mlfn(p.y / this.a, this.es, this.en);
s = Math.abs(lat);
if (s < HALF_PI) {
s = Math.sin(lat);
temp = this.long0 + p.x * Math.sqrt(1 - this.es * s * s) / (this.a * Math.cos(lat));
//temp = this.long0 + p.x / (this.a * Math.cos(lat));
lon = adjust_lon(temp);
}
else if ((s - EPSLN) < HALF_PI) {
lon = this.long0;
}
}
p.x = lon;
p.y = lat;
return p;
}
var names$22 = ["Sinusoidal", "sinu"];
var sinu = {
init: init$21,
forward: forward$20,
inverse: inverse$20,
names: names$22
};
function init$22() {}
/* Mollweide forward equations--mapping lat,long to x,y
----------------------------------------------------*/
function forward$21(p) {
/* Forward equations
-----------------*/
var lon = p.x;
var lat = p.y;
var delta_lon = adjust_lon(lon - this.long0);
var theta = lat;
var con = Math.PI * Math.sin(lat);
/* Iterate using the Newton-Raphson method to find theta
-----------------------------------------------------*/
while (true) {
var delta_theta = -(theta + Math.sin(theta) - con) / (1 + Math.cos(theta));
theta += delta_theta;
if (Math.abs(delta_theta) < EPSLN) {
break;
}
}
theta /= 2;
/* If the latitude is 90 deg, force the x coordinate to be "0 + false easting"
this is done here because of precision problems with "cos(theta)"
--------------------------------------------------------------------------*/
if (Math.PI / 2 - Math.abs(lat) < EPSLN) {
delta_lon = 0;
}
var x = 0.900316316158 * this.a * delta_lon * Math.cos(theta) + this.x0;
var y = 1.4142135623731 * this.a * Math.sin(theta) + this.y0;
p.x = x;
p.y = y;
return p;
}
function inverse$21(p) {
var theta;
var arg;
/* Inverse equations
-----------------*/
p.x -= this.x0;
p.y -= this.y0;
arg = p.y / (1.4142135623731 * this.a);
/* Because of division by zero problems, 'arg' can not be 1. Therefore
a number very close to one is used instead.
-------------------------------------------------------------------*/
if (Math.abs(arg) > 0.999999999999) {
arg = 0.999999999999;
}
theta = Math.asin(arg);
var lon = adjust_lon(this.long0 + (p.x / (0.900316316158 * this.a * Math.cos(theta))));
if (lon < (-Math.PI)) {
lon = -Math.PI;
}
if (lon > Math.PI) {
lon = Math.PI;
}
arg = (2 * theta + Math.sin(2 * theta)) / Math.PI;
if (Math.abs(arg) > 1) {
arg = 1;
}
var lat = Math.asin(arg);
p.x = lon;
p.y = lat;
return p;
}
var names$23 = ["Mollweide", "moll"];
var moll = {
init: init$22,
forward: forward$21,
inverse: inverse$21,
names: names$23
};
function init$23() {
/* Place parameters in static storage for common use
-------------------------------------------------*/
// Standard Parallels cannot be equal and on opposite sides of the equator
if (Math.abs(this.lat1 + this.lat2) < EPSLN) {
return;
}
this.lat2 = this.lat2 || this.lat1;
this.temp = this.b / this.a;
this.es = 1 - Math.pow(this.temp, 2);
this.e = Math.sqrt(this.es);
this.e0 = e0fn(this.es);
this.e1 = e1fn(this.es);
this.e2 = e2fn(this.es);
this.e3 = e3fn(this.es);
this.sinphi = Math.sin(this.lat1);
this.cosphi = Math.cos(this.lat1);
this.ms1 = msfnz(this.e, this.sinphi, this.cosphi);
this.ml1 = mlfn(this.e0, this.e1, this.e2, this.e3, this.lat1);
if (Math.abs(this.lat1 - this.lat2) < EPSLN) {
this.ns = this.sinphi;
}
else {
this.sinphi = Math.sin(this.lat2);
this.cosphi = Math.cos(this.lat2);
this.ms2 = msfnz(this.e, this.sinphi, this.cosphi);
this.ml2 = mlfn(this.e0, this.e1, this.e2, this.e3, this.lat2);
this.ns = (this.ms1 - this.ms2) / (this.ml2 - this.ml1);
}
this.g = this.ml1 + this.ms1 / this.ns;
this.ml0 = mlfn(this.e0, this.e1, this.e2, this.e3, this.lat0);
this.rh = this.a * (this.g - this.ml0);
}
/* Equidistant Conic forward equations--mapping lat,long to x,y
-----------------------------------------------------------*/
function forward$22(p) {
var lon = p.x;
var lat = p.y;
var rh1;
/* Forward equations
-----------------*/
if (this.sphere) {
rh1 = this.a * (this.g - lat);
}
else {
var ml = mlfn(this.e0, this.e1, this.e2, this.e3, lat);
rh1 = this.a * (this.g - ml);
}
var theta = this.ns * adjust_lon(lon - this.long0);
var x = this.x0 + rh1 * Math.sin(theta);
var y = this.y0 + this.rh - rh1 * Math.cos(theta);
p.x = x;
p.y = y;
return p;
}
/* Inverse equations
-----------------*/
function inverse$22(p) {
p.x -= this.x0;
p.y = this.rh - p.y + this.y0;
var con, rh1, lat, lon;
if (this.ns >= 0) {
rh1 = Math.sqrt(p.x * p.x + p.y * p.y);
con = 1;
}
else {
rh1 = -Math.sqrt(p.x * p.x + p.y * p.y);
con = -1;
}
var theta = 0;
if (rh1 !== 0) {
theta = Math.atan2(con * p.x, con * p.y);
}
if (this.sphere) {
lon = adjust_lon(this.long0 + theta / this.ns);
lat = adjust_lat(this.g - rh1 / this.a);
p.x = lon;
p.y = lat;
return p;
}
else {
var ml = this.g - rh1 / this.a;
lat = imlfn(ml, this.e0, this.e1, this.e2, this.e3);
lon = adjust_lon(this.long0 + theta / this.ns);
p.x = lon;
p.y = lat;
return p;
}
}
var names$24 = ["Equidistant_Conic", "eqdc"];
var eqdc = {
init: init$23,
forward: forward$22,
inverse: inverse$22,
names: names$24
};
/* Initialize the Van Der Grinten projection
----------------------------------------*/
function init$24() {
//this.R = 6370997; //Radius of earth
this.R = this.a;
}
function forward$23(p) {
var lon = p.x;
var lat = p.y;
/* Forward equations
-----------------*/
var dlon = adjust_lon(lon - this.long0);
var x, y;
if (Math.abs(lat) <= EPSLN) {
x = this.x0 + this.R * dlon;
y = this.y0;
}
var theta = asinz(2 * Math.abs(lat / Math.PI));
if ((Math.abs(dlon) <= EPSLN) || (Math.abs(Math.abs(lat) - HALF_PI) <= EPSLN)) {
x = this.x0;
if (lat >= 0) {
y = this.y0 + Math.PI * this.R * Math.tan(0.5 * theta);
}
else {
y = this.y0 + Math.PI * this.R * -Math.tan(0.5 * theta);
}
// return(OK);
}
var al = 0.5 * Math.abs((Math.PI / dlon) - (dlon / Math.PI));
var asq = al * al;
var sinth = Math.sin(theta);
var costh = Math.cos(theta);
var g = costh / (sinth + costh - 1);
var gsq = g * g;
var m = g * (2 / sinth - 1);
var msq = m * m;
var con = Math.PI * this.R * (al * (g - msq) + Math.sqrt(asq * (g - msq) * (g - msq) - (msq + asq) * (gsq - msq))) / (msq + asq);
if (dlon < 0) {
con = -con;
}
x = this.x0 + con;
//con = Math.abs(con / (Math.PI * this.R));
var q = asq + g;
con = Math.PI * this.R * (m * q - al * Math.sqrt((msq + asq) * (asq + 1) - q * q)) / (msq + asq);
if (lat >= 0) {
//y = this.y0 + Math.PI * this.R * Math.sqrt(1 - con * con - 2 * al * con);
y = this.y0 + con;
}
else {
//y = this.y0 - Math.PI * this.R * Math.sqrt(1 - con * con - 2 * al * con);
y = this.y0 - con;
}
p.x = x;
p.y = y;
return p;
}
/* Van Der Grinten inverse equations--mapping x,y to lat/long
---------------------------------------------------------*/
function inverse$23(p) {
var lon, lat;
var xx, yy, xys, c1, c2, c3;
var a1;
var m1;
var con;
var th1;
var d;
/* inverse equations
-----------------*/
p.x -= this.x0;
p.y -= this.y0;
con = Math.PI * this.R;
xx = p.x / con;
yy = p.y / con;
xys = xx * xx + yy * yy;
c1 = -Math.abs(yy) * (1 + xys);
c2 = c1 - 2 * yy * yy + xx * xx;
c3 = -2 * c1 + 1 + 2 * yy * yy + xys * xys;
d = yy * yy / c3 + (2 * c2 * c2 * c2 / c3 / c3 / c3 - 9 * c1 * c2 / c3 / c3) / 27;
a1 = (c1 - c2 * c2 / 3 / c3) / c3;
m1 = 2 * Math.sqrt(-a1 / 3);
con = ((3 * d) / a1) / m1;
if (Math.abs(con) > 1) {
if (con >= 0) {
con = 1;
}
else {
con = -1;
}
}
th1 = Math.acos(con) / 3;
if (p.y >= 0) {
lat = (-m1 * Math.cos(th1 + Math.PI / 3) - c2 / 3 / c3) * Math.PI;
}
else {
lat = -(-m1 * Math.cos(th1 + Math.PI / 3) - c2 / 3 / c3) * Math.PI;
}
if (Math.abs(xx) < EPSLN) {
lon = this.long0;
}
else {
lon = adjust_lon(this.long0 + Math.PI * (xys - 1 + Math.sqrt(1 + 2 * (xx * xx - yy * yy) + xys * xys)) / 2 / xx);
}
p.x = lon;
p.y = lat;
return p;
}
var names$25 = ["Van_der_Grinten_I", "VanDerGrinten", "vandg"];
var vandg = {
init: init$24,
forward: forward$23,
inverse: inverse$23,
names: names$25
};
function init$25() {
this.sin_p12 = Math.sin(this.lat0);
this.cos_p12 = Math.cos(this.lat0);
}
function forward$24(p) {
var lon = p.x;
var lat = p.y;
var sinphi = Math.sin(p.y);
var cosphi = Math.cos(p.y);
var dlon = adjust_lon(lon - this.long0);
var e0, e1, e2, e3, Mlp, Ml, tanphi, Nl1, Nl, psi, Az, G, H, GH, Hs, c, kp, cos_c, s, s2, s3, s4, s5;
if (this.sphere) {
if (Math.abs(this.sin_p12 - 1) <= EPSLN) {
//North Pole case
p.x = this.x0 + this.a * (HALF_PI - lat) * Math.sin(dlon);
p.y = this.y0 - this.a * (HALF_PI - lat) * Math.cos(dlon);
return p;
}
else if (Math.abs(this.sin_p12 + 1) <= EPSLN) {
//South Pole case
p.x = this.x0 + this.a * (HALF_PI + lat) * Math.sin(dlon);
p.y = this.y0 + this.a * (HALF_PI + lat) * Math.cos(dlon);
return p;
}
else {
//default case
cos_c = this.sin_p12 * sinphi + this.cos_p12 * cosphi * Math.cos(dlon);
c = Math.acos(cos_c);
kp = c ? c / Math.sin(c) : 1;
p.x = this.x0 + this.a * kp * cosphi * Math.sin(dlon);
p.y = this.y0 + this.a * kp * (this.cos_p12 * sinphi - this.sin_p12 * cosphi * Math.cos(dlon));
return p;
}
}
else {
e0 = e0fn(this.es);
e1 = e1fn(this.es);
e2 = e2fn(this.es);
e3 = e3fn(this.es);
if (Math.abs(this.sin_p12 - 1) <= EPSLN) {
//North Pole case
Mlp = this.a * mlfn(e0, e1, e2, e3, HALF_PI);
Ml = this.a * mlfn(e0, e1, e2, e3, lat);
p.x = this.x0 + (Mlp - Ml) * Math.sin(dlon);
p.y = this.y0 - (Mlp - Ml) * Math.cos(dlon);
return p;
}
else if (Math.abs(this.sin_p12 + 1) <= EPSLN) {
//South Pole case
Mlp = this.a * mlfn(e0, e1, e2, e3, HALF_PI);
Ml = this.a * mlfn(e0, e1, e2, e3, lat);
p.x = this.x0 + (Mlp + Ml) * Math.sin(dlon);
p.y = this.y0 + (Mlp + Ml) * Math.cos(dlon);
return p;
}
else {
//Default case
tanphi = sinphi / cosphi;
Nl1 = gN(this.a, this.e, this.sin_p12);
Nl = gN(this.a, this.e, sinphi);
psi = Math.atan((1 - this.es) * tanphi + this.es * Nl1 * this.sin_p12 / (Nl * cosphi));
Az = Math.atan2(Math.sin(dlon), this.cos_p12 * Math.tan(psi) - this.sin_p12 * Math.cos(dlon));
if (Az === 0) {
s = Math.asin(this.cos_p12 * Math.sin(psi) - this.sin_p12 * Math.cos(psi));
}
else if (Math.abs(Math.abs(Az) - Math.PI) <= EPSLN) {
s = -Math.asin(this.cos_p12 * Math.sin(psi) - this.sin_p12 * Math.cos(psi));
}
else {
s = Math.asin(Math.sin(dlon) * Math.cos(psi) / Math.sin(Az));
}
G = this.e * this.sin_p12 / Math.sqrt(1 - this.es);
H = this.e * this.cos_p12 * Math.cos(Az) / Math.sqrt(1 - this.es);
GH = G * H;
Hs = H * H;
s2 = s * s;
s3 = s2 * s;
s4 = s3 * s;
s5 = s4 * s;
c = Nl1 * s * (1 - s2 * Hs * (1 - Hs) / 6 + s3 / 8 * GH * (1 - 2 * Hs) + s4 / 120 * (Hs * (4 - 7 * Hs) - 3 * G * G * (1 - 7 * Hs)) - s5 / 48 * GH);
p.x = this.x0 + c * Math.sin(Az);
p.y = this.y0 + c * Math.cos(Az);
return p;
}
}
}
function inverse$24(p) {
p.x -= this.x0;
p.y -= this.y0;
var rh, z, sinz, cosz, lon, lat, con, e0, e1, e2, e3, Mlp, M, N1, psi, Az, cosAz, tmp, A, B, D, Ee, F, sinpsi;
if (this.sphere) {
rh = Math.sqrt(p.x * p.x + p.y * p.y);
if (rh > (2 * HALF_PI * this.a)) {
return;
}
z = rh / this.a;
sinz = Math.sin(z);
cosz = Math.cos(z);
lon = this.long0;
if (Math.abs(rh) <= EPSLN) {
lat = this.lat0;
}
else {
lat = asinz(cosz * this.sin_p12 + (p.y * sinz * this.cos_p12) / rh);
con = Math.abs(this.lat0) - HALF_PI;
if (Math.abs(con) <= EPSLN) {
if (this.lat0 >= 0) {
lon = adjust_lon(this.long0 + Math.atan2(p.x, - p.y));
}
else {
lon = adjust_lon(this.long0 - Math.atan2(-p.x, p.y));
}
}
else {
/*con = cosz - this.sin_p12 * Math.sin(lat);
if ((Math.abs(con) < EPSLN) && (Math.abs(p.x) < EPSLN)) {
//no-op, just keep the lon value as is
} else {
var temp = Math.atan2((p.x * sinz * this.cos_p12), (con * rh));
lon = adjust_lon(this.long0 + Math.atan2((p.x * sinz * this.cos_p12), (con * rh)));
}*/
lon = adjust_lon(this.long0 + Math.atan2(p.x * sinz, rh * this.cos_p12 * cosz - p.y * this.sin_p12 * sinz));
}
}
p.x = lon;
p.y = lat;
return p;
}
else {
e0 = e0fn(this.es);
e1 = e1fn(this.es);
e2 = e2fn(this.es);
e3 = e3fn(this.es);
if (Math.abs(this.sin_p12 - 1) <= EPSLN) {
//North pole case
Mlp = this.a * mlfn(e0, e1, e2, e3, HALF_PI);
rh = Math.sqrt(p.x * p.x + p.y * p.y);
M = Mlp - rh;
lat = imlfn(M / this.a, e0, e1, e2, e3);
lon = adjust_lon(this.long0 + Math.atan2(p.x, - 1 * p.y));
p.x = lon;
p.y = lat;
return p;
}
else if (Math.abs(this.sin_p12 + 1) <= EPSLN) {
//South pole case
Mlp = this.a * mlfn(e0, e1, e2, e3, HALF_PI);
rh = Math.sqrt(p.x * p.x + p.y * p.y);
M = rh - Mlp;
lat = imlfn(M / this.a, e0, e1, e2, e3);
lon = adjust_lon(this.long0 + Math.atan2(p.x, p.y));
p.x = lon;
p.y = lat;
return p;
}
else {
//default case
rh = Math.sqrt(p.x * p.x + p.y * p.y);
Az = Math.atan2(p.x, p.y);
N1 = gN(this.a, this.e, this.sin_p12);
cosAz = Math.cos(Az);
tmp = this.e * this.cos_p12 * cosAz;
A = -tmp * tmp / (1 - this.es);
B = 3 * this.es * (1 - A) * this.sin_p12 * this.cos_p12 * cosAz / (1 - this.es);
D = rh / N1;
Ee = D - A * (1 + A) * Math.pow(D, 3) / 6 - B * (1 + 3 * A) * Math.pow(D, 4) / 24;
F = 1 - A * Ee * Ee / 2 - D * Ee * Ee * Ee / 6;
psi = Math.asin(this.sin_p12 * Math.cos(Ee) + this.cos_p12 * Math.sin(Ee) * cosAz);
lon = adjust_lon(this.long0 + Math.asin(Math.sin(Az) * Math.sin(Ee) / Math.cos(psi)));
sinpsi = Math.sin(psi);
lat = Math.atan2((sinpsi - this.es * F * this.sin_p12) * Math.tan(psi), sinpsi * (1 - this.es));
p.x = lon;
p.y = lat;
return p;
}
}
}
var names$26 = ["Azimuthal_Equidistant", "aeqd"];
var aeqd = {
init: init$25,
forward: forward$24,
inverse: inverse$24,
names: names$26
};
function init$26() {
//double temp; /* temporary variable */
/* Place parameters in static storage for common use
-------------------------------------------------*/
this.sin_p14 = Math.sin(this.lat0);
this.cos_p14 = Math.cos(this.lat0);
}
/* Orthographic forward equations--mapping lat,long to x,y
---------------------------------------------------*/
function forward$25(p) {
var sinphi, cosphi; /* sin and cos value */
var dlon; /* delta longitude value */
var coslon; /* cos of longitude */
var ksp; /* scale factor */
var g, x, y;
var lon = p.x;
var lat = p.y;
/* Forward equations
-----------------*/
dlon = adjust_lon(lon - this.long0);
sinphi = Math.sin(lat);
cosphi = Math.cos(lat);
coslon = Math.cos(dlon);
g = this.sin_p14 * sinphi + this.cos_p14 * cosphi * coslon;
ksp = 1;
if ((g > 0) || (Math.abs(g) <= EPSLN)) {
x = this.a * ksp * cosphi * Math.sin(dlon);
y = this.y0 + this.a * ksp * (this.cos_p14 * sinphi - this.sin_p14 * cosphi * coslon);
}
p.x = x;
p.y = y;
return p;
}
function inverse$25(p) {
var rh; /* height above ellipsoid */
var z; /* angle */
var sinz, cosz; /* sin of z and cos of z */
var con;
var lon, lat;
/* Inverse equations
-----------------*/
p.x -= this.x0;
p.y -= this.y0;
rh = Math.sqrt(p.x * p.x + p.y * p.y);
z = asinz(rh / this.a);
sinz = Math.sin(z);
cosz = Math.cos(z);
lon = this.long0;
if (Math.abs(rh) <= EPSLN) {
lat = this.lat0;
p.x = lon;
p.y = lat;
return p;
}
lat = asinz(cosz * this.sin_p14 + (p.y * sinz * this.cos_p14) / rh);
con = Math.abs(this.lat0) - HALF_PI;
if (Math.abs(con) <= EPSLN) {
if (this.lat0 >= 0) {
lon = adjust_lon(this.long0 + Math.atan2(p.x, - p.y));
}
else {
lon = adjust_lon(this.long0 - Math.atan2(-p.x, p.y));
}
p.x = lon;
p.y = lat;
return p;
}
lon = adjust_lon(this.long0 + Math.atan2((p.x * sinz), rh * this.cos_p14 * cosz - p.y * this.sin_p14 * sinz));
p.x = lon;
p.y = lat;
return p;
}
var names$27 = ["ortho"];
var ortho = {
init: init$26,
forward: forward$25,
inverse: inverse$25,
names: names$27
};
// QSC projection rewritten from the original PROJ4
// https://github.com/OSGeo/proj.4/blob/master/src/PJ_qsc.c
/* constants */
var FACE_ENUM = {
FRONT: 1,
RIGHT: 2,
BACK: 3,
LEFT: 4,
TOP: 5,
BOTTOM: 6
};
var AREA_ENUM = {
AREA_0: 1,
AREA_1: 2,
AREA_2: 3,
AREA_3: 4
};
function init$27() {
this.x0 = this.x0 || 0;
this.y0 = this.y0 || 0;
this.lat0 = this.lat0 || 0;
this.long0 = this.long0 || 0;
this.lat_ts = this.lat_ts || 0;
this.title = this.title || "Quadrilateralized Spherical Cube";
/* Determine the cube face from the center of projection. */
if (this.lat0 >= HALF_PI - FORTPI / 2.0) {
this.face = FACE_ENUM.TOP;
} else if (this.lat0 <= -(HALF_PI - FORTPI / 2.0)) {
this.face = FACE_ENUM.BOTTOM;
} else if (Math.abs(this.long0) <= FORTPI) {
this.face = FACE_ENUM.FRONT;
} else if (Math.abs(this.long0) <= HALF_PI + FORTPI) {
this.face = this.long0 > 0.0 ? FACE_ENUM.RIGHT : FACE_ENUM.LEFT;
} else {
this.face = FACE_ENUM.BACK;
}
/* Fill in useful values for the ellipsoid <-> sphere shift
* described in [LK12]. */
if (this.es !== 0) {
this.one_minus_f = 1 - (this.a - this.b) / this.a;
this.one_minus_f_squared = this.one_minus_f * this.one_minus_f;
}
}
// QSC forward equations--mapping lat,long to x,y
// -----------------------------------------------------------------
function forward$26(p) {
var xy = {x: 0, y: 0};
var lat, lon;
var theta, phi;
var t, mu;
/* nu; */
var area = {value: 0};
// move lon according to projection's lon
p.x -= this.long0;
/* Convert the geodetic latitude to a geocentric latitude.
* This corresponds to the shift from the ellipsoid to the sphere
* described in [LK12]. */
if (this.es !== 0) {//if (P->es != 0) {
lat = Math.atan(this.one_minus_f_squared * Math.tan(p.y));
} else {
lat = p.y;
}
/* Convert the input lat, lon into theta, phi as used by QSC.
* This depends on the cube face and the area on it.
* For the top and bottom face, we can compute theta and phi
* directly from phi, lam. For the other faces, we must use
* unit sphere cartesian coordinates as an intermediate step. */
lon = p.x; //lon = lp.lam;
if (this.face === FACE_ENUM.TOP) {
phi = HALF_PI - lat;
if (lon >= FORTPI && lon <= HALF_PI + FORTPI) {
area.value = AREA_ENUM.AREA_0;
theta = lon - HALF_PI;
} else if (lon > HALF_PI + FORTPI || lon <= -(HALF_PI + FORTPI)) {
area.value = AREA_ENUM.AREA_1;
theta = (lon > 0.0 ? lon - SPI : lon + SPI);
} else if (lon > -(HALF_PI + FORTPI) && lon <= -FORTPI) {
area.value = AREA_ENUM.AREA_2;
theta = lon + HALF_PI;
} else {
area.value = AREA_ENUM.AREA_3;
theta = lon;
}
} else if (this.face === FACE_ENUM.BOTTOM) {
phi = HALF_PI + lat;
if (lon >= FORTPI && lon <= HALF_PI + FORTPI) {
area.value = AREA_ENUM.AREA_0;
theta = -lon + HALF_PI;
} else if (lon < FORTPI && lon >= -FORTPI) {
area.value = AREA_ENUM.AREA_1;
theta = -lon;
} else if (lon < -FORTPI && lon >= -(HALF_PI + FORTPI)) {
area.value = AREA_ENUM.AREA_2;
theta = -lon - HALF_PI;
} else {
area.value = AREA_ENUM.AREA_3;
theta = (lon > 0.0 ? -lon + SPI : -lon - SPI);
}
} else {
var q, r, s;
var sinlat, coslat;
var sinlon, coslon;
if (this.face === FACE_ENUM.RIGHT) {
lon = qsc_shift_lon_origin(lon, +HALF_PI);
} else if (this.face === FACE_ENUM.BACK) {
lon = qsc_shift_lon_origin(lon, +SPI);
} else if (this.face === FACE_ENUM.LEFT) {
lon = qsc_shift_lon_origin(lon, -HALF_PI);
}
sinlat = Math.sin(lat);
coslat = Math.cos(lat);
sinlon = Math.sin(lon);
coslon = Math.cos(lon);
q = coslat * coslon;
r = coslat * sinlon;
s = sinlat;
if (this.face === FACE_ENUM.FRONT) {
phi = Math.acos(q);
theta = qsc_fwd_equat_face_theta(phi, s, r, area);
} else if (this.face === FACE_ENUM.RIGHT) {
phi = Math.acos(r);
theta = qsc_fwd_equat_face_theta(phi, s, -q, area);
} else if (this.face === FACE_ENUM.BACK) {
phi = Math.acos(-q);
theta = qsc_fwd_equat_face_theta(phi, s, -r, area);
} else if (this.face === FACE_ENUM.LEFT) {
phi = Math.acos(-r);
theta = qsc_fwd_equat_face_theta(phi, s, q, area);
} else {
/* Impossible */
phi = theta = 0;
area.value = AREA_ENUM.AREA_0;
}
}
/* Compute mu and nu for the area of definition.
* For mu, see Eq. (3-21) in [OL76], but note the typos:
* compare with Eq. (3-14). For nu, see Eq. (3-38). */
mu = Math.atan((12 / SPI) * (theta + Math.acos(Math.sin(theta) * Math.cos(FORTPI)) - HALF_PI));
t = Math.sqrt((1 - Math.cos(phi)) / (Math.cos(mu) * Math.cos(mu)) / (1 - Math.cos(Math.atan(1 / Math.cos(theta)))));
/* Apply the result to the real area. */
if (area.value === AREA_ENUM.AREA_1) {
mu += HALF_PI;
} else if (area.value === AREA_ENUM.AREA_2) {
mu += SPI;
} else if (area.value === AREA_ENUM.AREA_3) {
mu += 1.5 * SPI;
}
/* Now compute x, y from mu and nu */
xy.x = t * Math.cos(mu);
xy.y = t * Math.sin(mu);
xy.x = xy.x * this.a + this.x0;
xy.y = xy.y * this.a + this.y0;
p.x = xy.x;
p.y = xy.y;
return p;
}
// QSC inverse equations--mapping x,y to lat/long
// -----------------------------------------------------------------
function inverse$26(p) {
var lp = {lam: 0, phi: 0};
var mu, nu, cosmu, tannu;
var tantheta, theta, cosphi, phi;
var t;
var area = {value: 0};
/* de-offset */
p.x = (p.x - this.x0) / this.a;
p.y = (p.y - this.y0) / this.a;
/* Convert the input x, y to the mu and nu angles as used by QSC.
* This depends on the area of the cube face. */
nu = Math.atan(Math.sqrt(p.x * p.x + p.y * p.y));
mu = Math.atan2(p.y, p.x);
if (p.x >= 0.0 && p.x >= Math.abs(p.y)) {
area.value = AREA_ENUM.AREA_0;
} else if (p.y >= 0.0 && p.y >= Math.abs(p.x)) {
area.value = AREA_ENUM.AREA_1;
mu -= HALF_PI;
} else if (p.x < 0.0 && -p.x >= Math.abs(p.y)) {
area.value = AREA_ENUM.AREA_2;
mu = (mu < 0.0 ? mu + SPI : mu - SPI);
} else {
area.value = AREA_ENUM.AREA_3;
mu += HALF_PI;
}
/* Compute phi and theta for the area of definition.
* The inverse projection is not described in the original paper, but some
* good hints can be found here (as of 2011-12-14):
* http://fits.gsfc.nasa.gov/fitsbits/saf.93/saf.9302
* (search for "Message-Id: <9302181759.AA25477 at fits.cv.nrao.edu>") */
t = (SPI / 12) * Math.tan(mu);
tantheta = Math.sin(t) / (Math.cos(t) - (1 / Math.sqrt(2)));
theta = Math.atan(tantheta);
cosmu = Math.cos(mu);
tannu = Math.tan(nu);
cosphi = 1 - cosmu * cosmu * tannu * tannu * (1 - Math.cos(Math.atan(1 / Math.cos(theta))));
if (cosphi < -1) {
cosphi = -1;
} else if (cosphi > +1) {
cosphi = +1;
}
/* Apply the result to the real area on the cube face.
* For the top and bottom face, we can compute phi and lam directly.
* For the other faces, we must use unit sphere cartesian coordinates
* as an intermediate step. */
if (this.face === FACE_ENUM.TOP) {
phi = Math.acos(cosphi);
lp.phi = HALF_PI - phi;
if (area.value === AREA_ENUM.AREA_0) {
lp.lam = theta + HALF_PI;
} else if (area.value === AREA_ENUM.AREA_1) {
lp.lam = (theta < 0.0 ? theta + SPI : theta - SPI);
} else if (area.value === AREA_ENUM.AREA_2) {
lp.lam = theta - HALF_PI;
} else /* area.value == AREA_ENUM.AREA_3 */ {
lp.lam = theta;
}
} else if (this.face === FACE_ENUM.BOTTOM) {
phi = Math.acos(cosphi);
lp.phi = phi - HALF_PI;
if (area.value === AREA_ENUM.AREA_0) {
lp.lam = -theta + HALF_PI;
} else if (area.value === AREA_ENUM.AREA_1) {
lp.lam = -theta;
} else if (area.value === AREA_ENUM.AREA_2) {
lp.lam = -theta - HALF_PI;
} else /* area.value == AREA_ENUM.AREA_3 */ {
lp.lam = (theta < 0.0 ? -theta - SPI : -theta + SPI);
}
} else {
/* Compute phi and lam via cartesian unit sphere coordinates. */
var q, r, s;
q = cosphi;
t = q * q;
if (t >= 1) {
s = 0;
} else {
s = Math.sqrt(1 - t) * Math.sin(theta);
}
t += s * s;
if (t >= 1) {
r = 0;
} else {
r = Math.sqrt(1 - t);
}
/* Rotate q,r,s into the correct area. */
if (area.value === AREA_ENUM.AREA_1) {
t = r;
r = -s;
s = t;
} else if (area.value === AREA_ENUM.AREA_2) {
r = -r;
s = -s;
} else if (area.value === AREA_ENUM.AREA_3) {
t = r;
r = s;
s = -t;
}
/* Rotate q,r,s into the correct cube face. */
if (this.face === FACE_ENUM.RIGHT) {
t = q;
q = -r;
r = t;
} else if (this.face === FACE_ENUM.BACK) {
q = -q;
r = -r;
} else if (this.face === FACE_ENUM.LEFT) {
t = q;
q = r;
r = -t;
}
/* Now compute phi and lam from the unit sphere coordinates. */
lp.phi = Math.acos(-s) - HALF_PI;
lp.lam = Math.atan2(r, q);
if (this.face === FACE_ENUM.RIGHT) {
lp.lam = qsc_shift_lon_origin(lp.lam, -HALF_PI);
} else if (this.face === FACE_ENUM.BACK) {
lp.lam = qsc_shift_lon_origin(lp.lam, -SPI);
} else if (this.face === FACE_ENUM.LEFT) {
lp.lam = qsc_shift_lon_origin(lp.lam, +HALF_PI);
}
}
/* Apply the shift from the sphere to the ellipsoid as described
* in [LK12]. */
if (this.es !== 0) {
var invert_sign;
var tanphi, xa;
invert_sign = (lp.phi < 0 ? 1 : 0);
tanphi = Math.tan(lp.phi);
xa = this.b / Math.sqrt(tanphi * tanphi + this.one_minus_f_squared);
lp.phi = Math.atan(Math.sqrt(this.a * this.a - xa * xa) / (this.one_minus_f * xa));
if (invert_sign) {
lp.phi = -lp.phi;
}
}
lp.lam += this.long0;
p.x = lp.lam;
p.y = lp.phi;
return p;
}
/* Helper function for forward projection: compute the theta angle
* and determine the area number. */
function qsc_fwd_equat_face_theta(phi, y, x, area) {
var theta;
if (phi < EPSLN) {
area.value = AREA_ENUM.AREA_0;
theta = 0.0;
} else {
theta = Math.atan2(y, x);
if (Math.abs(theta) <= FORTPI) {
area.value = AREA_ENUM.AREA_0;
} else if (theta > FORTPI && theta <= HALF_PI + FORTPI) {
area.value = AREA_ENUM.AREA_1;
theta -= HALF_PI;
} else if (theta > HALF_PI + FORTPI || theta <= -(HALF_PI + FORTPI)) {
area.value = AREA_ENUM.AREA_2;
theta = (theta >= 0.0 ? theta - SPI : theta + SPI);
} else {
area.value = AREA_ENUM.AREA_3;
theta += HALF_PI;
}
}
return theta;
}
/* Helper function: shift the longitude. */
function qsc_shift_lon_origin(lon, offset) {
var slon = lon + offset;
if (slon < -SPI) {
slon += TWO_PI;
} else if (slon > +SPI) {
slon -= TWO_PI;
}
return slon;
}
var names$28 = ["Quadrilateralized Spherical Cube", "Quadrilateralized_Spherical_Cube", "qsc"];
var qsc = {
init: init$27,
forward: forward$26,
inverse: inverse$26,
names: names$28
};
// Robinson projection
// Based on https://github.com/OSGeo/proj.4/blob/master/src/PJ_robin.c
// Polynomial coeficients from http://article.gmane.org/gmane.comp.gis.proj-4.devel/6039
var COEFS_X = [
[1.0000, 2.2199e-17, -7.15515e-05, 3.1103e-06],
[0.9986, -0.000482243, -2.4897e-05, -1.3309e-06],
[0.9954, -0.00083103, -4.48605e-05, -9.86701e-07],
[0.9900, -0.00135364, -5.9661e-05, 3.6777e-06],
[0.9822, -0.00167442, -4.49547e-06, -5.72411e-06],
[0.9730, -0.00214868, -9.03571e-05, 1.8736e-08],
[0.9600, -0.00305085, -9.00761e-05, 1.64917e-06],
[0.9427, -0.00382792, -6.53386e-05, -2.6154e-06],
[0.9216, -0.00467746, -0.00010457, 4.81243e-06],
[0.8962, -0.00536223, -3.23831e-05, -5.43432e-06],
[0.8679, -0.00609363, -0.000113898, 3.32484e-06],
[0.8350, -0.00698325, -6.40253e-05, 9.34959e-07],
[0.7986, -0.00755338, -5.00009e-05, 9.35324e-07],
[0.7597, -0.00798324, -3.5971e-05, -2.27626e-06],
[0.7186, -0.00851367, -7.01149e-05, -8.6303e-06],
[0.6732, -0.00986209, -0.000199569, 1.91974e-05],
[0.6213, -0.010418, 8.83923e-05, 6.24051e-06],
[0.5722, -0.00906601, 0.000182, 6.24051e-06],
[0.5322, -0.00677797, 0.000275608, 6.24051e-06]
];
var COEFS_Y = [
[-5.20417e-18, 0.0124, 1.21431e-18, -8.45284e-11],
[0.0620, 0.0124, -1.26793e-09, 4.22642e-10],
[0.1240, 0.0124, 5.07171e-09, -1.60604e-09],
[0.1860, 0.0123999, -1.90189e-08, 6.00152e-09],
[0.2480, 0.0124002, 7.10039e-08, -2.24e-08],
[0.3100, 0.0123992, -2.64997e-07, 8.35986e-08],
[0.3720, 0.0124029, 9.88983e-07, -3.11994e-07],
[0.4340, 0.0123893, -3.69093e-06, -4.35621e-07],
[0.4958, 0.0123198, -1.02252e-05, -3.45523e-07],
[0.5571, 0.0121916, -1.54081e-05, -5.82288e-07],
[0.6176, 0.0119938, -2.41424e-05, -5.25327e-07],
[0.6769, 0.011713, -3.20223e-05, -5.16405e-07],
[0.7346, 0.0113541, -3.97684e-05, -6.09052e-07],
[0.7903, 0.0109107, -4.89042e-05, -1.04739e-06],
[0.8435, 0.0103431, -6.4615e-05, -1.40374e-09],
[0.8936, 0.00969686, -6.4636e-05, -8.547e-06],
[0.9394, 0.00840947, -0.000192841, -4.2106e-06],
[0.9761, 0.00616527, -0.000256, -4.2106e-06],
[1.0000, 0.00328947, -0.000319159, -4.2106e-06]
];
var FXC = 0.8487;
var FYC = 1.3523;
var C1 = R2D/5; // rad to 5-degree interval
var RC1 = 1/C1;
var NODES = 18;
var poly3_val = function(coefs, x) {
return coefs[0] + x * (coefs[1] + x * (coefs[2] + x * coefs[3]));
};
var poly3_der = function(coefs, x) {
return coefs[1] + x * (2 * coefs[2] + x * 3 * coefs[3]);
};
function newton_rapshon(f_df, start, max_err, iters) {
var x = start;
for (; iters; --iters) {
var upd = f_df(x);
x -= upd;
if (Math.abs(upd) < max_err) {
break;
}
}
return x;
}
function init$28() {
this.x0 = this.x0 || 0;
this.y0 = this.y0 || 0;
this.long0 = this.long0 || 0;
this.es = 0;
this.title = this.title || "Robinson";
}
function forward$27(ll) {
var lon = adjust_lon(ll.x - this.long0);
var dphi = Math.abs(ll.y);
var i = Math.floor(dphi * C1);
if (i < 0) {
i = 0;
} else if (i >= NODES) {
i = NODES - 1;
}
dphi = R2D * (dphi - RC1 * i);
var xy = {
x: poly3_val(COEFS_X[i], dphi) * lon,
y: poly3_val(COEFS_Y[i], dphi)
};
if (ll.y < 0) {
xy.y = -xy.y;
}
xy.x = xy.x * this.a * FXC + this.x0;
xy.y = xy.y * this.a * FYC + this.y0;
return xy;
}
function inverse$27(xy) {
var ll = {
x: (xy.x - this.x0) / (this.a * FXC),
y: Math.abs(xy.y - this.y0) / (this.a * FYC)
};
if (ll.y >= 1) { // pathologic case
ll.x /= COEFS_X[NODES][0];
ll.y = xy.y < 0 ? -HALF_PI : HALF_PI;
} else {
// find table interval
var i = Math.floor(ll.y * NODES);
if (i < 0) {
i = 0;
} else if (i >= NODES) {
i = NODES - 1;
}
for (;;) {
if (COEFS_Y[i][0] > ll.y) {
--i;
} else if (COEFS_Y[i+1][0] <= ll.y) {
++i;
} else {
break;
}
}
// linear interpolation in 5 degree interval
var coefs = COEFS_Y[i];
var t = 5 * (ll.y - coefs[0]) / (COEFS_Y[i+1][0] - coefs[0]);
// find t so that poly3_val(coefs, t) = ll.y
t = newton_rapshon(function(x) {
return (poly3_val(coefs, x) - ll.y) / poly3_der(coefs, x);
}, t, EPSLN, 100);
ll.x /= poly3_val(COEFS_X[i], t);
ll.y = (5 * i + t) * D2R;
if (xy.y < 0) {
ll.y = -ll.y;
}
}
ll.x = adjust_lon(ll.x + this.long0);
return ll;
}
var names$29 = ["Robinson", "robin"];
var robin = {
init: init$28,
forward: forward$27,
inverse: inverse$27,
names: names$29
};
function init$29() {
this.name = 'geocent';
}
function forward$28(p) {
var point = geodeticToGeocentric(p, this.es, this.a);
return point;
}
function inverse$28(p) {
var point = geocentricToGeodetic(p, this.es, this.a, this.b);
return point;
}
var names$30 = ["Geocentric", 'geocentric', "geocent", "Geocent"];
var geocent = {
init: init$29,
forward: forward$28,
inverse: inverse$28,
names: names$30
};
var mode = {
N_POLE: 0,
S_POLE: 1,
EQUIT: 2,
OBLIQ: 3
};
var params = {
h: { def: 100000, num: true }, // default is Karman line, no default in PROJ.7
azi: { def: 0, num: true, degrees: true }, // default is North
tilt: { def: 0, num: true, degrees: true }, // default is Nadir
long0: { def: 0, num: true }, // default is Greenwich, conversion to rad is automatic
lat0: { def: 0, num: true } // default is Equator, conversion to rad is automatic
};
function init$30() {
Object.keys(params).forEach(function (p) {
if (typeof this[p] === "undefined") {
this[p] = params[p].def;
} else if (params[p].num && isNaN(this[p])) {
throw new Error("Invalid parameter value, must be numeric " + p + " = " + this[p]);
} else if (params[p].num) {
this[p] = parseFloat(this[p]);
}
if (params[p].degrees) {
this[p] = this[p] * D2R;
}
}.bind(this));
if (Math.abs((Math.abs(this.lat0) - HALF_PI)) < EPSLN) {
this.mode = this.lat0 < 0 ? mode.S_POLE : mode.N_POLE;
} else if (Math.abs(this.lat0) < EPSLN) {
this.mode = mode.EQUIT;
} else {
this.mode = mode.OBLIQ;
this.sinph0 = Math.sin(this.lat0);
this.cosph0 = Math.cos(this.lat0);
}
this.pn1 = this.h / this.a; // Normalize relative to the Earth's radius
if (this.pn1 <= 0 || this.pn1 > 1e10) {
throw new Error("Invalid height");
}
this.p = 1 + this.pn1;
this.rp = 1 / this.p;
this.h1 = 1 / this.pn1;
this.pfact = (this.p + 1) * this.h1;
this.es = 0;
var omega = this.tilt;
var gamma = this.azi;
this.cg = Math.cos(gamma);
this.sg = Math.sin(gamma);
this.cw = Math.cos(omega);
this.sw = Math.sin(omega);
}
function forward$29(p) {
p.x -= this.long0;
var sinphi = Math.sin(p.y);
var cosphi = Math.cos(p.y);
var coslam = Math.cos(p.x);
var x, y;
switch (this.mode) {
case mode.OBLIQ:
y = this.sinph0 * sinphi + this.cosph0 * cosphi * coslam;
break;
case mode.EQUIT:
y = cosphi * coslam;
break;
case mode.S_POLE:
y = -sinphi;
break;
case mode.N_POLE:
y = sinphi;
break;
}
y = this.pn1 / (this.p - y);
x = y * cosphi * Math.sin(p.x);
switch (this.mode) {
case mode.OBLIQ:
y *= this.cosph0 * sinphi - this.sinph0 * cosphi * coslam;
break;
case mode.EQUIT:
y *= sinphi;
break;
case mode.N_POLE:
y *= -(cosphi * coslam);
break;
case mode.S_POLE:
y *= cosphi * coslam;
break;
}
// Tilt
var yt, ba;
yt = y * this.cg + x * this.sg;
ba = 1 / (yt * this.sw * this.h1 + this.cw);
x = (x * this.cg - y * this.sg) * this.cw * ba;
y = yt * ba;
p.x = x * this.a;
p.y = y * this.a;
return p;
}
function inverse$29(p) {
p.x /= this.a;
p.y /= this.a;
var r = { x: p.x, y: p.y };
// Un-Tilt
var bm, bq, yt;
yt = 1 / (this.pn1 - p.y * this.sw);
bm = this.pn1 * p.x * yt;
bq = this.pn1 * p.y * this.cw * yt;
p.x = bm * this.cg + bq * this.sg;
p.y = bq * this.cg - bm * this.sg;
var rh = hypot(p.x, p.y);
if (Math.abs(rh) < EPSLN) {
r.x = 0;
r.y = p.y;
} else {
var cosz, sinz;
sinz = 1 - rh * rh * this.pfact;
sinz = (this.p - Math.sqrt(sinz)) / (this.pn1 / rh + rh / this.pn1);
cosz = Math.sqrt(1 - sinz * sinz);
switch (this.mode) {
case mode.OBLIQ:
r.y = Math.asin(cosz * this.sinph0 + p.y * sinz * this.cosph0 / rh);
p.y = (cosz - this.sinph0 * Math.sin(r.y)) * rh;
p.x *= sinz * this.cosph0;
break;
case mode.EQUIT:
r.y = Math.asin(p.y * sinz / rh);
p.y = cosz * rh;
p.x *= sinz;
break;
case mode.N_POLE:
r.y = Math.asin(cosz);
p.y = -p.y;
break;
case mode.S_POLE:
r.y = -Math.asin(cosz);
break;
}
r.x = Math.atan2(p.x, p.y);
}
p.x = r.x + this.long0;
p.y = r.y;
return p;
}
var names$31 = ["Tilted_Perspective", "tpers"];
var tpers = {
init: init$30,
forward: forward$29,
inverse: inverse$29,
names: names$31
};
function init$31() {
this.flip_axis = (this.sweep === 'x' ? 1 : 0);
this.h = Number(this.h);
this.radius_g_1 = this.h / this.a;
if (this.radius_g_1 <= 0 || this.radius_g_1 > 1e10) {
throw new Error();
}
this.radius_g = 1.0 + this.radius_g_1;
this.C = this.radius_g * this.radius_g - 1.0;
if (this.es !== 0.0) {
var one_es = 1.0 - this.es;
var rone_es = 1 / one_es;
this.radius_p = Math.sqrt(one_es);
this.radius_p2 = one_es;
this.radius_p_inv2 = rone_es;
this.shape = 'ellipse'; // Use as a condition in the forward and inverse functions.
} else {
this.radius_p = 1.0;
this.radius_p2 = 1.0;
this.radius_p_inv2 = 1.0;
this.shape = 'sphere'; // Use as a condition in the forward and inverse functions.
}
if (!this.title) {
this.title = "Geostationary Satellite View";
}
}
function forward$30(p) {
var lon = p.x;
var lat = p.y;
var tmp, v_x, v_y, v_z;
lon = lon - this.long0;
if (this.shape === 'ellipse') {
lat = Math.atan(this.radius_p2 * Math.tan(lat));
var r = this.radius_p / hypot(this.radius_p * Math.cos(lat), Math.sin(lat));
v_x = r * Math.cos(lon) * Math.cos(lat);
v_y = r * Math.sin(lon) * Math.cos(lat);
v_z = r * Math.sin(lat);
if (((this.radius_g - v_x) * v_x - v_y * v_y - v_z * v_z * this.radius_p_inv2) < 0.0) {
p.x = Number.NaN;
p.y = Number.NaN;
return p;
}
tmp = this.radius_g - v_x;
if (this.flip_axis) {
p.x = this.radius_g_1 * Math.atan(v_y / hypot(v_z, tmp));
p.y = this.radius_g_1 * Math.atan(v_z / tmp);
} else {
p.x = this.radius_g_1 * Math.atan(v_y / tmp);
p.y = this.radius_g_1 * Math.atan(v_z / hypot(v_y, tmp));
}
} else if (this.shape === 'sphere') {
tmp = Math.cos(lat);
v_x = Math.cos(lon) * tmp;
v_y = Math.sin(lon) * tmp;
v_z = Math.sin(lat);
tmp = this.radius_g - v_x;
if (this.flip_axis) {
p.x = this.radius_g_1 * Math.atan(v_y / hypot(v_z, tmp));
p.y = this.radius_g_1 * Math.atan(v_z / tmp);
} else {
p.x = this.radius_g_1 * Math.atan(v_y / tmp);
p.y = this.radius_g_1 * Math.atan(v_z / hypot(v_y, tmp));
}
}
p.x = p.x * this.a;
p.y = p.y * this.a;
return p;
}
function inverse$30(p) {
var v_x = -1.0;
var v_y = 0.0;
var v_z = 0.0;
var a, b, det, k;
p.x = p.x / this.a;
p.y = p.y / this.a;
if (this.shape === 'ellipse') {
if (this.flip_axis) {
v_z = Math.tan(p.y / this.radius_g_1);
v_y = Math.tan(p.x / this.radius_g_1) * hypot(1.0, v_z);
} else {
v_y = Math.tan(p.x / this.radius_g_1);
v_z = Math.tan(p.y / this.radius_g_1) * hypot(1.0, v_y);
}
var v_zp = v_z / this.radius_p;
a = v_y * v_y + v_zp * v_zp + v_x * v_x;
b = 2 * this.radius_g * v_x;
det = (b * b) - 4 * a * this.C;
if (det < 0.0) {
p.x = Number.NaN;
p.y = Number.NaN;
return p;
}
k = (-b - Math.sqrt(det)) / (2.0 * a);
v_x = this.radius_g + k * v_x;
v_y *= k;
v_z *= k;
p.x = Math.atan2(v_y, v_x);
p.y = Math.atan(v_z * Math.cos(p.x) / v_x);
p.y = Math.atan(this.radius_p_inv2 * Math.tan(p.y));
} else if (this.shape === 'sphere') {
if (this.flip_axis) {
v_z = Math.tan(p.y / this.radius_g_1);
v_y = Math.tan(p.x / this.radius_g_1) * Math.sqrt(1.0 + v_z * v_z);
} else {
v_y = Math.tan(p.x / this.radius_g_1);
v_z = Math.tan(p.y / this.radius_g_1) * Math.sqrt(1.0 + v_y * v_y);
}
a = v_y * v_y + v_z * v_z + v_x * v_x;
b = 2 * this.radius_g * v_x;
det = (b * b) - 4 * a * this.C;
if (det < 0.0) {
p.x = Number.NaN;
p.y = Number.NaN;
return p;
}
k = (-b - Math.sqrt(det)) / (2.0 * a);
v_x = this.radius_g + k * v_x;
v_y *= k;
v_z *= k;
p.x = Math.atan2(v_y, v_x);
p.y = Math.atan(v_z * Math.cos(p.x) / v_x);
}
p.x = p.x + this.long0;
return p;
}
var names$32 = ["Geostationary Satellite View", "Geostationary_Satellite", "geos"];
var geos = {
init: init$31,
forward: forward$30,
inverse: inverse$30,
names: names$32,
};
var includedProjections = function(proj4){
proj4.Proj.projections.add(tmerc);
proj4.Proj.projections.add(etmerc);
proj4.Proj.projections.add(utm);
proj4.Proj.projections.add(sterea);
proj4.Proj.projections.add(stere);
proj4.Proj.projections.add(somerc);
proj4.Proj.projections.add(omerc);
proj4.Proj.projections.add(lcc);
proj4.Proj.projections.add(krovak);
proj4.Proj.projections.add(cass);
proj4.Proj.projections.add(laea);
proj4.Proj.projections.add(aea);
proj4.Proj.projections.add(gnom);
proj4.Proj.projections.add(cea);
proj4.Proj.projections.add(eqc);
proj4.Proj.projections.add(poly);
proj4.Proj.projections.add(nzmg);
proj4.Proj.projections.add(mill);
proj4.Proj.projections.add(sinu);
proj4.Proj.projections.add(moll);
proj4.Proj.projections.add(eqdc);
proj4.Proj.projections.add(vandg);
proj4.Proj.projections.add(aeqd);
proj4.Proj.projections.add(ortho);
proj4.Proj.projections.add(qsc);
proj4.Proj.projections.add(robin);
proj4.Proj.projections.add(geocent);
proj4.Proj.projections.add(tpers);
proj4.Proj.projections.add(geos);
};
proj4$1.defaultDatum = 'WGS84'; //default datum
proj4$1.Proj = Projection;
proj4$1.WGS84 = new proj4$1.Proj('WGS84');
proj4$1.Point = Point;
proj4$1.toPoint = toPoint;
proj4$1.defs = defs;
proj4$1.nadgrid = nadgrid;
proj4$1.transform = transform;
proj4$1.mgrs = mgrs;
proj4$1.version = '2.9.0';
includedProjections(proj4$1);
return proj4$1;
})));