Press n or j to go to the next uncovered block, b, p or k for the previous block.
| 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 | 2x 2x 2x 2x 44x 44x 44x 44x 44x 44x 44x 44x 44x 44x 44x 44x 120x 120x 120x 44x 42x 42x 42x 42x 2x 2x 117x 284x 284x 120x 120x 120x 2x | // Adapted from https://github.com/mourner/geokdbush
// ISC License
// Copyright (c) 2017, Vladimir Agafonkin
// Permission to use, copy, modify, and/or distribute this software for any purpose
// with or without fee is hereby granted, provided that the above copyright notice
// and this permission notice appear in all copies.
// THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL WARRANTIES WITH
// REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY AND
// FITNESS. IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT,
// INDIRECT, OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS
// OF USE, DATA OR PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER
// TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR PERFORMANCE OF
// THIS SOFTWARE.
import { KDBush } from './kdbush';
import TinyQueue from './tinyqueue';
const earthRadius = 6371;
const rad = Math.PI / 180;
export function around(index: KDBush, lng: number, lat: number, maxResults = Infinity, maxDistance = Infinity) {
let maxHaverSinDist = 1;
const result = [];
Iif (maxResults === undefined) maxResults = Infinity;
if (maxDistance !== undefined) maxHaverSinDist = haverSin(maxDistance / earthRadius);
// a distance-sorted priority queue that will contain both points and kd-tree nodes
const q = new TinyQueue([], compareDist);
// an object that represents the top kd-tree node (the whole Earth)
let node = {
left: 0, // left index in the kd-tree array
right: index.ids.length - 1, // right index
axis: 0, // 0 for longitude axis and 1 for latitude axis
dist: 0, // will hold the lower bound of children's distances to the query point
minLng: -180, // bounding box of the node
minLat: -90,
maxLng: 180,
maxLat: 90
};
const cosLat = Math.cos(lat * rad);
while (node) {
const right = node.right;
const left = node.left;
if (right - left <= index.nodeSize) { // leaf node
// add all points of the leaf node to the queue
for (let i = left; i <= right; i++) {
const id = index.ids[i];
const dist = haverSinDist(lng, lat, index.coords[2 * i], index.coords[2 * i + 1], cosLat);
q.push({ id, dist });
}
} else E{ // not a leaf node (has child nodes)
const m = (left + right) >> 1; // middle index
const midLng = index.coords[2 * m];
const midLat = index.coords[2 * m + 1];
// add middle point to the queue
const id = index.ids[m];
const dist = haverSinDist(lng, lat, midLng, midLat, cosLat);
q.push({ id, dist });
const nextAxis = (node.axis + 1) % 2;
// first half of the node
const leftNode = {
left,
right: m - 1,
axis: nextAxis,
minLng: node.minLng,
minLat: node.minLat,
maxLng: node.axis === 0 ? midLng : node.maxLng,
maxLat: node.axis === 1 ? midLat : node.maxLat,
dist: 0
};
// second half of the node
const rightNode = {
left: m + 1,
right,
axis: nextAxis,
minLng: node.axis === 0 ? midLng : node.minLng,
minLat: node.axis === 1 ? midLat : node.minLat,
maxLng: node.maxLng,
maxLat: node.maxLat,
dist: 0
};
leftNode.dist = boxDist(lng, lat, cosLat, leftNode);
rightNode.dist = boxDist(lng, lat, cosLat, rightNode);
// add child nodes to the queue
q.push(leftNode);
q.push(rightNode);
}
// fetch closest points from the queue; they're guaranteed to be closer
// than all remaining points (both individual and those in kd-tree nodes),
// since each node's distance is a lower bound of distances to its children
while (q.length && q.peek().id != null) {
const candidate = q.pop()!;
Iif (candidate.dist > maxHaverSinDist) return result;
result.push(candidate.id);
if (result.length === maxResults) return result;
}
// the next closest kd-tree node
node = q.pop();
}
return result;
}
// lower bound for distance from a location to points inside a bounding box
function boxDist(lng: number, lat: number, cosLat: number, node: any) {
const minLng = node.minLng;
const maxLng = node.maxLng;
const minLat = node.minLat;
const maxLat = node.maxLat;
// query point is between minimum and maximum longitudes
Iif (lng >= minLng && lng <= maxLng) {
Iif (lat < minLat) return haverSin((lat - minLat) * rad);
Iif (lat > maxLat) return haverSin((lat - maxLat) * rad);
return 0;
}
// query point is west or east of the bounding box;
// calculate the extremum for great circle distance from query point to the closest longitude;
const haverSinDLng = Math.min(haverSin((lng - minLng) * rad), haverSin((lng - maxLng) * rad));
const extremumLat = vertexLat(lat, haverSinDLng);
// if extremum is inside the box, return the distance to it
Iif (extremumLat > minLat && extremumLat < maxLat) {
return haverSinDistPartial(haverSinDLng, cosLat, lat, extremumLat);
}
// otherwise return the distan e to one of the bbox corners (whichever is closest)
return Math.min(
haverSinDistPartial(haverSinDLng, cosLat, lat, minLat),
haverSinDistPartial(haverSinDLng, cosLat, lat, maxLat)
);
}
function compareDist(a: any, b: any) {
return a.dist - b.dist;
}
function haverSin(theta: number) {
const s = Math.sin(theta / 2);
return s * s;
}
function haverSinDistPartial(haverSinDLng: number, cosLat1: number, lat1: number, lat2: number) {
return cosLat1 * Math.cos(lat2 * rad) * haverSinDLng + haverSin((lat1 - lat2) * rad);
}
function haverSinDist(lng1: number, lat1: number, lng2: number, lat2: number, cosLat1: number) {
const haverSinDLng = haverSin((lng1 - lng2) * rad);
return haverSinDistPartial(haverSinDLng, cosLat1, lat1, lat2);
}
export function distance(lng1: number, lat1: number, lng2: number, lat2: number) {
const h = haverSinDist(lng1, lat1, lng2, lat2, Math.cos(lat1 * rad));
return 2 * earthRadius * Math.asin(Math.sqrt(h));
}
function vertexLat(lat: number, haverSinDLng: number) {
const cosDLng = 1 - 2 * haverSinDLng;
Iif (cosDLng <= 0) return lat > 0 ? 90 : -90;
return Math.atan(Math.tan(lat * rad) / cosDLng) / rad;
} |