oxjs/source/Ox/js/Geo.js

453 lines
16 KiB
JavaScript
Raw Normal View History

2011-11-05 16:46:53 +00:00
'use strict';
(function() {
// fixme: make all this work with different types of "points"
// i.e. {lat, lng}, [lat, lng]
2011-12-31 12:57:02 +00:00
function deg(point) {
return Ox.map(point, function(val) {
return Ox.mod(Ox.deg(val) + 180, 360) - 180;
});
}
function rad(point) {
2011-12-31 12:57:02 +00:00
return Ox.map(point, function(val) {
return Ox.rad(val);
});
}
function splitArea(area) {
return Ox.crossesDateline(area.sw, area.ne) ? [
{sw: area.sw, ne: {lat: area.ne.lat, lng: 180}},
{sw: {lat: area.sw.lat, lng: -180}, ne: area.ne}
] : [area];
}
/*@
2011-12-31 12:57:02 +00:00
Ox.crossesDateline <f> Returns true if a given line crosses the dateline
2012-01-06 11:57:54 +00:00
> Ox.crossesDateline({lat: 0, lng: -90}, {lat: 0, lng: 90})
false
> Ox.crossesDateline({lat: 0, lng: 90}, {lat: 0, lng: -90})
true
@*/
// FIXME: argument should be {w: ..., e: ...}
2011-12-31 12:57:02 +00:00
Ox.crossesDateline = function(pointA, pointB) {
return pointA.lng > pointB.lng;
};
/*@
Ox.getArea <f> Returns the area in square meters of a given rectancle
@*/
// FIXME: argument should be {sw: ..., ne: ...}
2011-12-31 12:57:02 +00:00
Ox.getArea = function(pointA, pointB) {
/*
2012-01-04 07:42:48 +00:00
area of a ring between two latitudes:
2 * PI * r^2 * abs(sin(latA) - sin(latB))
see http://mathforum.org/library/drmath/view/63767.html
=>
2 * Math.PI
* Math.pow(Ox.EARTH_RADIUS, 2)
* Math.abs(Math.sin(Ox.rad(latA)) - Math.sin(Ox.rad(latB)))
* Math.abs(Ox.rad(lngA) - Ox.rad(lngB))
/ (2 * Math.PI)
*/
2011-12-31 12:57:02 +00:00
if (Ox.crossesDateline(pointA, pointB)) {
pointB.lng += 360;
}
2011-12-31 12:57:02 +00:00
pointA = rad(pointA);
pointB = rad(pointB);
return Math.pow(Ox.EARTH_RADIUS, 2)
* Math.abs(Math.sin(pointA.lat) - Math.sin(pointB.lat))
* Math.abs(pointA.lng - pointB.lng);
};
/*@
Ox.getBearing <f> Returns the bearing from one point to another
> Ox.getBearing({lat: -45, lng: 0}, {lat: 45, lng: 0})
0
2012-01-06 11:57:54 +00:00
> Ox.getBearing({lat: 0, lng: -90}, {lat: 0, lng: 90})
90
@*/
2011-12-31 12:57:02 +00:00
Ox.getBearing = function(pointA, pointB) {
var pointA = rad(pointA),
pointB = rad(pointB),
x = Math.cos(pointA.lat) * Math.sin(pointB.lat)
- Math.sin(pointA.lat) * Math.cos(pointB.lat)
* Math.cos(pointB.lng - pointA.lng),
y = Math.sin(pointB.lng - pointA.lng)
* Math.cos(pointB.lat);
return (Ox.deg(Math.atan2(y, x)) + 360) % 360;
};
/*@
Ox.getCenter <f> Returns the center of a recangle on a spehre
> Ox.getCenter({lat: -45, lng: -90}, {lat: 45, lng: 90})
{lat: 0, lng: 0}
@*/
2011-12-31 12:57:02 +00:00
Ox.getCenter = function(pointA, pointB) {
var pointA = rad(pointA),
pointB = rad(pointB),
x = Math.cos(pointB.lat)
* Math.cos(pointB.lng - pointA.lng),
y = Math.cos(pointB.lat)
* Math.sin(pointB.lng - pointA.lng),
d = Math.sqrt(
2011-12-31 12:57:02 +00:00
Math.pow(Math.cos(pointA.lat) + x, 2) + Math.pow(y, 2)
),
2011-12-31 12:57:02 +00:00
lat = Math.atan2(Math.sin(pointA.lat) + Math.sin(pointB.lat), d),
lng = pointA.lng + Math.atan2(y, Math.cos(pointA.lat) + x);
return deg({lat: lat, lng: lng});
};
/*@
Ox.getCircle <f> Returns points on a circle around a given point
(center, radius, precision) -> <a> Points
center <o> Center point ({lat, lng})
radius <n> Radius in meters
precision <n> Precision (the circle will have 2^precision segments)
@*/
Ox.getCircle = function(center, radius, precision) {
return Ox.range(
0, 360, 360 / Math.pow(2, precision)
).map(function(bearing) {
return Ox.getPoint(center, radius, bearing);
});
};
/*@
Ox.getDegreesPerMeter <f> Returns degrees per meter at a given latitude
> 360 / Ox.getDegreesPerMeter(0)
Ox.EARTH_CIRCUMFERENCE
@*/
Ox.getDegreesPerMeter = function(lat) {
return 360 / Ox.EARTH_CIRCUMFERENCE / Math.cos(lat * Math.PI / 180);
};
/*@
Ox.getDistance <f> Returns the distance in meters between two points
> Ox.getDistance({lat: -45, lng: -90}, {lat: 45, lng: 90}) * 2
Ox.EARTH_CIRCUMFERENCE
@*/
2011-12-31 12:57:02 +00:00
Ox.getDistance = function(pointA, pointB) {
var pointA = rad(pointA),
pointB = rad(pointB);
return Math.acos(
2011-12-31 12:57:02 +00:00
Math.sin(pointA.lat) * Math.sin(pointB.lat)
+ Math.cos(pointA.lat) * Math.cos(pointB.lat)
* Math.cos(pointB.lng - pointA.lng)
) * Ox.EARTH_RADIUS;
};
/*@
Ox.getLatLngByXY <f> Returns lat/lng for a given x/y on a 1x1 mercator projection
> Ox.getLatLngByXY({x: 0.5, y: 0.5})
{lat: 0, lng: 0}
@*/
Ox.getLatLngByXY = function(xy) {
function getVal(val) {
return (val - 0.5) * 2 * Math.PI;
}
return {
lat: -Ox.deg(Math.atan(Ox.sinh(getVal(xy.y)))),
lng: Ox.deg(getVal(xy.x))
2011-12-31 12:57:02 +00:00
};
};
/*@
Ox.getLine <f> Returns points on a line between two points
(pointA, pointB, precision) -> <a> Points
pointA <o> Start point ({lat, lng})
pointB <o> End point ({lat, lng})
precision <n> Precision (the line will have 2^precision segments)
@*/
Ox.getLine = function(pointA, pointB, precision) {
var line = [pointA, pointB], points;
while (precision > 0) {
points = [line[0]];
Ox.loop(line.length - 1, function(i) {
points.push(
Ox.getCenter(line[i], line[i + 1]),
line[i + 1]
);
});
line = points;
precision--;
}
2011-12-31 12:57:02 +00:00
return line;
};
/*@
Ox.getMetersPerDegree <f> Returns meters per degree at a given latitude
> Ox.getMetersPerDegree(0) * 360
Ox.EARTH_CIRCUMFERENCE
@*/
Ox.getMetersPerDegree = function(lat) {
return Math.cos(lat * Math.PI / 180) * Ox.EARTH_CIRCUMFERENCE / 360;
};
2011-12-31 12:57:02 +00:00
/*@
Ox.getPoint <f> Returns a point at a given distance/bearing from a given point
> Ox.getPoint({lat: -45, lng: 0}, Ox.EARTH_CIRCUMFERENCE / 4, 0)
{lat: 45, lng: 0}
@*/
Ox.getPoint = function(point, distance, bearing) {
var pointB = {};
point = rad(point);
distance /= Ox.EARTH_RADIUS;
bearing = Ox.rad(bearing);
pointB.lat = Math.asin(
Math.sin(point.lat) * Math.cos(distance)
+ Math.cos(point.lat) * Math.sin(distance) * Math.cos(bearing)
);
pointB.lng = point.lng + Math.atan2(
Math.sin(bearing) * Math.sin(distance) * Math.cos(point.lat),
Math.cos(distance) - Math.sin(point.lat) * Math.sin(pointB.lat)
);
return deg(pointB);
};
/*@
Ox.getXYByLatLng <f> Returns x/y on a 1x1 mercator projection for a given lat/lng
> Ox.getXYByLatLng({lat: 0, lng: 0})
{x: 0.5, y: 0.5}
@*/
Ox.getXYByLatLng = function(latlng) {
function getVal(val) {
return (val / (2 * Math.PI) + 0.5)
}
return {
x: getVal(Ox.rad(latlng.lng)),
y: getVal(Ox.asinh(Math.tan(Ox.rad(-latlng.lat))))
};
};
2011-12-31 12:57:02 +00:00
/*@
Ox.isPolar <f> Returns true if a given point is outside the bounds of a mercator projection
2012-01-06 11:57:54 +00:00
> Ox.isPolar({lat: 90, lng: 0})
true
2011-12-31 12:57:02 +00:00
@*/
Ox.isPolar = function(point) {
return point.lat < Ox.MIN_LATITUDE || point.lat > Ox.MAX_LATITUDE;
};
/*@
Ox.containsArea <f> Returns true if an area contains another area
<script>
Ox.test.areas = [
{sw: {lat: -30, lng: -30}, ne: {lat: 30, lng: 30}},
{sw: {lat: -20, lng: -40}, ne: {lat: 20, lng: 40}},
{sw: {lat: -30, lng: 150}, ne: {lat: 30, lng: -150}},
{sw: {lat: 10, lng: -170}, ne: {lat: 20, lng: -160}}
];
</script>
> Ox.containsArea(Ox.test.areas[0], Ox.test.areas[1])
false
> Ox.containsArea(Ox.test.areas[2], Ox.test.areas[3])
true
@*/
2012-04-02 21:11:56 +00:00
// FIXME: Shouldn't this be rewritten as a test
// if the intersection is equal to the inner area?
Ox.containsArea = function(areaA, areaB) {
// If an area crosses the dateline,
// we split it into two parts,
// west and east of the dateline
var areas = [areaA, areaB].map(splitArea),
ret;
function contains(areaA, areaB) {
return areaA.sw.lat <= areaB.sw.lat
&& areaA.sw.lng <= areaB.sw.lng
&& areaA.ne.lat >= areaB.ne.lat
&& areaA.ne.lng >= areaB.ne.lng;
}
// For each part of the inner area, test if it
// is contained in any part of the outer area
Ox.forEach(areas[1], function(area1) {
Ox.forEach(areas[0], function(area0) {
ret = contains(area0, area1);
// Break if the outer part contains the inner part
2012-05-22 07:11:26 +00:00
ret && Ox.break();
});
// Break if no outer part contains the inner part
2012-05-22 07:11:26 +00:00
!ret && Ox.break();
});
return ret;
};
/*@
Ox.intersectAreas <f> Returns the intersection of two areas, or null
<script>
Ox.test.areas = [
{sw: {lat: -10, lng: -10}, ne: {lat: 0, lng: 0}},
{sw: {lat: 0, lng: 0}, ne: {lat: 10, lng: 10}},
{sw: {lat: -30, lng: 150}, ne: {lat: 30, lng: -150}},
{sw: {lat: 25, lng: -155}, ne: {lat: 35, lng: -145}}
];
</script>
> Ox.intersectAreas([Ox.test.areas[0], Ox.test.areas[1]])
{sw: {lat: 0, lng: 0}, ne: {lat: 0, lng: 0}}
> Ox.intersectAreas([Ox.test.areas[2], Ox.test.areas[3]])
{sw: {lat: 25, lng: -155}, ne: {lat: 30, lng: -150}}
@*/
2012-04-02 21:11:56 +00:00
// FIXME: handle the a corner case where
// two areas have two intersections
Ox.intersectAreas = function(areas) {
var intersections, ret;
// If an area crosses the dateline,
// we split it into two parts,
// west and east of the dateline
areas = areas.map(splitArea);
ret = areas[0];
function intersect(areaA, areaB) {
return areaA.sw.lat > areaB.ne.lat
|| areaA.sw.lng > areaB.ne.lng
|| areaA.ne.lat < areaB.sw.lat
|| areaA.ne.lng < areaB.sw.lng
? null : {
sw: {
lat: Math.max(areaA.sw.lat, areaB.sw.lat),
lng: Math.max(areaA.sw.lng, areaB.sw.lng)
},
ne: {
lat: Math.min(areaA.ne.lat, areaB.ne.lat),
lng: Math.min(areaA.ne.lng, areaB.ne.lng)
}
};
}
2012-05-24 09:02:59 +00:00
Ox.forEach(areas.slice(1), function(parts) {
if (ret.length == 1 && parts.length == 1) {
ret = intersect(ret[0], parts[0]);
} else {
// intersect each part of the intersection
// with all parts of the next area
intersections = Ox.compact(ret.map(function(part) {
2012-05-24 07:45:33 +00:00
return Ox.intersectAreas(parts.concat(part));
}));
ret = intersections.length == 0 ? null
: Ox.joinAreas(intersections);
}
if (ret === null) {
2012-05-22 07:11:26 +00:00
Ox.break();
} else {
ret = splitArea(ret);
}
});
return ret ? Ox.joinAreas(ret) : null;
};
/*@
Ox.joinAreas <f> Joins an array of areas
<script>
Ox.test.areas = [
{sw: {lat: -30, lng: 150}, ne: {lat: -20, lng: 160}},
{sw: {lat: -10, lng: 170}, ne: {lat: 10, lng: -170}},
{sw: {lat: 20, lng: -160}, ne: {lat: 30, lng: -150}}
];
</script>
> Ox.joinAreas(Ox.test.areas)
{sw: {lat: -30, lng: 150}, ne: {lat: 30, lng: -150}}
@*/
Ox.joinAreas = function(areas) {
// While the combined latitude is trivial (min to max), the combined longitude
// spans from the eastern to the western edge of the largest gap between areas
var ret = areas[0],
gaps = [{
sw: {lat: -90, lng: ret.ne.lng},
ne: {lat: 90, lng: ret.sw.lng}
}];
function containsGaps(area) {
2012-05-22 14:29:37 +00:00
return Ox.getIndices(gaps, function(gap) {
return Ox.containsArea({
sw: {lat: -90, lng: area.sw.lng},
ne: {lat: 90, lng: area.ne.lng}
2012-05-22 14:29:37 +00:00
}, gap);
});
}
function intersectsWithGaps(area) {
var ret = {};
gaps.forEach(function(gap, i) {
var intersection = Ox.intersectAreas([area, gap]);
if (intersection) {
ret[i] = intersection;
}
});
return ret;
}
function isContainedInGap(area) {
var ret = -1;
Ox.forEach(gaps, function(gap, i) {
if (Ox.containsArea(gap, area)) {
ret = i;
2012-05-22 07:11:26 +00:00
Ox.break();
}
});
return ret;
}
2012-05-24 09:02:59 +00:00
areas.slice(1).forEach(function(area) {
var index, indices, intersections;
if (area.sw.lat < ret.sw.lat) {
ret.sw.lat = area.sw.lat;
}
if (area.ne.lat > ret.ne.lat) {
ret.ne.lat = area.ne.lat;
}
// If the area is contained in a gap, split the gap in two
index = isContainedInGap(area);
if (index > -1) {
gaps.push({
sw: gaps[index].sw,
ne: {lat: 90, lng: area.sw.lng}
});
gaps.push({
sw: {lat: -90, lng: area.ne.lng},
ne: gaps[index].ne
});
gaps.splice(index, 1);
} else {
// If the area contains gaps, remove them
indices = containsGaps(area);
Ox.reverse(indices).forEach(function(index) {
gaps.splice(index, 1);
});
// If the area intersects with gaps, shrink them
intersections = intersectsWithGaps(area);
Ox.forEach(intersections, function(intersection, index) {
gaps[index] = {
sw: {
lat: -90,
lng: gaps[index].sw.lng == intersection.sw.lng
? intersection.ne.lng : gaps[index].sw.lng
},
ne: {
lat: 90,
lng: gaps[index].ne.lng == intersection.ne.lng
? intersection.sw.lng : gaps[index].ne.lng
}
};
});
}
});
if (gaps.length == 0) {
ret.sw.lng = -180;
ret.ne.lng = 180;
} else {
gaps.sort(function(a, b) {
return (
b.ne.lng
+ (Ox.crossesDateline(b.sw, b.ne) ? 360 : 0)
- b.sw.lng
) - (
a.ne.lng
+ (Ox.crossesDateline(a.sw, a.ne) ? 360 : 0)
- a.sw.lng
);
});
ret.sw.lng = gaps[0].ne.lng;
ret.ne.lng = gaps[0].sw.lng;
}
return ret;
};
}());