Skip to content

Instantly share code, notes, and snippets.

@migurski
Created May 14, 2012 03:29
Show Gist options
  • Save migurski/2691614 to your computer and use it in GitHub Desktop.
Save migurski/2691614 to your computer and use it in GitHub Desktop.
Lambert projection for Modest Maps
var LambertProjection = function(lat0, lon0, lat1, lat2, xmin, ymin, xmax, ymax)
{
var pi = Math.PI, ln = Math.log, pow = Math.pow,
sin = Math.sin, cos = Math.cos, tan = Math.tan,
atan = Math.atan, sqrt = Math.sqrt;
function sec(t) { return 1 / cos(t); }
function cot(t) { return 1 / tan(t); }
function deg2rad(deg) { return pi * deg / 180; }
var phi0 = deg2rad(lat0),
lam0 = deg2rad(lon0),
phi1 = deg2rad(lat1),
phi2 = deg2rad(lat2);
this.lam0 = lam0;
// Equation 4
this.n = ln(cos(phi1) * sec(phi2));
this.n /= ln(tan(pi/4 + phi2/2) * cot(pi/4 + phi1/2));
// Equation 3
this.F = cos(phi1) * pow(tan(pi/4 + phi1/2), this.n) / this.n;
// Equation 6
this.P0 = this.F * pow(cot(pi/4 + phi0/2), this.n);
//
var tx = 1 / (xmax - xmin),
dx = .5 - (xmin/2 + xmax/2) * tx;
var ty = -tx,
dy = .5 - (ymin/2 + ymax/2) * ty;
var t = new MM.Transformation(tx, 0, dx, 0, ty, dy);
MM.Projection.call(this, 0, t);
};
LambertProjection.prototype = {
rawProject: function(point)
{
var pi = Math.PI, ln = Math.log, pow = Math.pow,
sin = Math.sin, cos = Math.cos, tan = Math.tan,
atan = Math.atan, sqrt = Math.sqrt;
function cot(t) { return 1 / tan(t); }
var lam = point.x, // longitude already in radians
phi = point.y; // latitude already in radians
// Equation 5
var P = this.F * pow(cot(pi/4 + phi/2), this.n);
// Equations 1 & 2
var x = P * sin(this.n * (lam - this.lam0)),
y = this.P0 - P * cos(this.n * (lam - this.lam0));
return {x: x * 6378137., y: y * 6378137.};
},
rawUnproject: function(point)
{
var pi = Math.PI, ln = Math.log, pow = Math.pow,
sin = Math.sin, cos = Math.cos, tan = Math.tan,
atan = Math.atan, sqrt = Math.sqrt;
function sgn(x)
{
if(x >= 0) {
return 1;
} else{
return -1;
}
}
var x = point.x / 6378137.,
y = point.y / 6378137.;
// Equation 9
var P = sgn(this.n) * sqrt(pow(x, 2) + pow((this.P0 - y), 2));
// Equation 10
var theta = atan(x / (this.P0 - y));
// Equation 7
var phi = 2 * atan(pow(this.F/P, 1/this.n)) - pi/2;
// Equation 8
var lam = this.lam0 + theta/this.n;
return {x: lam, y: phi};
}
};
MM.extend(LambertProjection, MM.Projection);
// adjust for latitude.
var scale = Math.cos(Math.PI * 37.5/180);
// echo the usual +/- 20m extent of spherical mercator.
var xmin = scale * -6378137 * Math.PI,
ymin = scale * -6378137 * Math.PI,
xmax = scale * 6378137 * Math.PI,
ymax = scale * 6378137 * Math.PI;
var lcc = new LambertProjection(37.5, -96, 29.5, 45.5, xmin, ymin, xmax, ymax);
var loc = {y: Math.PI * 48.38544/180, x: Math.PI * -124.72916/180};
var proj = lcc.rawProject(loc);
var loc2 = lcc.rawUnproject(proj);
console.log([loc, proj, loc2]);
var loc3 = {lat: 48.38544, lon: -124.72916};
var coord = lcc.locationCoordinate(loc3);
var loc4 = lcc.coordinateLocation(coord);
console.log([loc3, coord, loc4]);
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment