contours.js 6.8 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225
  1. import {extent, nice, thresholdSturges, ticks} from "d3-array";
  2. import {slice} from "./array.js";
  3. import ascending from "./ascending.js";
  4. import area from "./area.js";
  5. import constant from "./constant.js";
  6. import contains from "./contains.js";
  7. import noop from "./noop.js";
  8. var cases = [
  9. [],
  10. [[[1.0, 1.5], [0.5, 1.0]]],
  11. [[[1.5, 1.0], [1.0, 1.5]]],
  12. [[[1.5, 1.0], [0.5, 1.0]]],
  13. [[[1.0, 0.5], [1.5, 1.0]]],
  14. [[[1.0, 1.5], [0.5, 1.0]], [[1.0, 0.5], [1.5, 1.0]]],
  15. [[[1.0, 0.5], [1.0, 1.5]]],
  16. [[[1.0, 0.5], [0.5, 1.0]]],
  17. [[[0.5, 1.0], [1.0, 0.5]]],
  18. [[[1.0, 1.5], [1.0, 0.5]]],
  19. [[[0.5, 1.0], [1.0, 0.5]], [[1.5, 1.0], [1.0, 1.5]]],
  20. [[[1.5, 1.0], [1.0, 0.5]]],
  21. [[[0.5, 1.0], [1.5, 1.0]]],
  22. [[[1.0, 1.5], [1.5, 1.0]]],
  23. [[[0.5, 1.0], [1.0, 1.5]]],
  24. []
  25. ];
  26. export default function() {
  27. var dx = 1,
  28. dy = 1,
  29. threshold = thresholdSturges,
  30. smooth = smoothLinear;
  31. function contours(values) {
  32. var tz = threshold(values);
  33. // Convert number of thresholds into uniform thresholds.
  34. if (!Array.isArray(tz)) {
  35. const e = extent(values, finite);
  36. tz = ticks(...nice(e[0], e[1], tz), tz);
  37. while (tz[tz.length - 1] >= e[1]) tz.pop();
  38. while (tz[1] < e[0]) tz.shift();
  39. } else {
  40. tz = tz.slice().sort(ascending);
  41. }
  42. return tz.map(value => contour(values, value));
  43. }
  44. // Accumulate, smooth contour rings, assign holes to exterior rings.
  45. // Based on https://github.com/mbostock/shapefile/blob/v0.6.2/shp/polygon.js
  46. function contour(values, value) {
  47. const v = value == null ? NaN : +value;
  48. if (isNaN(v)) throw new Error(`invalid value: ${value}`);
  49. var polygons = [],
  50. holes = [];
  51. isorings(values, v, function(ring) {
  52. smooth(ring, values, v);
  53. if (area(ring) > 0) polygons.push([ring]);
  54. else holes.push(ring);
  55. });
  56. holes.forEach(function(hole) {
  57. for (var i = 0, n = polygons.length, polygon; i < n; ++i) {
  58. if (contains((polygon = polygons[i])[0], hole) !== -1) {
  59. polygon.push(hole);
  60. return;
  61. }
  62. }
  63. });
  64. return {
  65. type: "MultiPolygon",
  66. value: value,
  67. coordinates: polygons
  68. };
  69. }
  70. // Marching squares with isolines stitched into rings.
  71. // Based on https://github.com/topojson/topojson-client/blob/v3.0.0/src/stitch.js
  72. function isorings(values, value, callback) {
  73. var fragmentByStart = new Array,
  74. fragmentByEnd = new Array,
  75. x, y, t0, t1, t2, t3;
  76. // Special case for the first row (y = -1, t2 = t3 = 0).
  77. x = y = -1;
  78. t1 = above(values[0], value);
  79. cases[t1 << 1].forEach(stitch);
  80. while (++x < dx - 1) {
  81. t0 = t1, t1 = above(values[x + 1], value);
  82. cases[t0 | t1 << 1].forEach(stitch);
  83. }
  84. cases[t1 << 0].forEach(stitch);
  85. // General case for the intermediate rows.
  86. while (++y < dy - 1) {
  87. x = -1;
  88. t1 = above(values[y * dx + dx], value);
  89. t2 = above(values[y * dx], value);
  90. cases[t1 << 1 | t2 << 2].forEach(stitch);
  91. while (++x < dx - 1) {
  92. t0 = t1, t1 = above(values[y * dx + dx + x + 1], value);
  93. t3 = t2, t2 = above(values[y * dx + x + 1], value);
  94. cases[t0 | t1 << 1 | t2 << 2 | t3 << 3].forEach(stitch);
  95. }
  96. cases[t1 | t2 << 3].forEach(stitch);
  97. }
  98. // Special case for the last row (y = dy - 1, t0 = t1 = 0).
  99. x = -1;
  100. t2 = values[y * dx] >= value;
  101. cases[t2 << 2].forEach(stitch);
  102. while (++x < dx - 1) {
  103. t3 = t2, t2 = above(values[y * dx + x + 1], value);
  104. cases[t2 << 2 | t3 << 3].forEach(stitch);
  105. }
  106. cases[t2 << 3].forEach(stitch);
  107. function stitch(line) {
  108. var start = [line[0][0] + x, line[0][1] + y],
  109. end = [line[1][0] + x, line[1][1] + y],
  110. startIndex = index(start),
  111. endIndex = index(end),
  112. f, g;
  113. if (f = fragmentByEnd[startIndex]) {
  114. if (g = fragmentByStart[endIndex]) {
  115. delete fragmentByEnd[f.end];
  116. delete fragmentByStart[g.start];
  117. if (f === g) {
  118. f.ring.push(end);
  119. callback(f.ring);
  120. } else {
  121. fragmentByStart[f.start] = fragmentByEnd[g.end] = {start: f.start, end: g.end, ring: f.ring.concat(g.ring)};
  122. }
  123. } else {
  124. delete fragmentByEnd[f.end];
  125. f.ring.push(end);
  126. fragmentByEnd[f.end = endIndex] = f;
  127. }
  128. } else if (f = fragmentByStart[endIndex]) {
  129. if (g = fragmentByEnd[startIndex]) {
  130. delete fragmentByStart[f.start];
  131. delete fragmentByEnd[g.end];
  132. if (f === g) {
  133. f.ring.push(end);
  134. callback(f.ring);
  135. } else {
  136. fragmentByStart[g.start] = fragmentByEnd[f.end] = {start: g.start, end: f.end, ring: g.ring.concat(f.ring)};
  137. }
  138. } else {
  139. delete fragmentByStart[f.start];
  140. f.ring.unshift(start);
  141. fragmentByStart[f.start = startIndex] = f;
  142. }
  143. } else {
  144. fragmentByStart[startIndex] = fragmentByEnd[endIndex] = {start: startIndex, end: endIndex, ring: [start, end]};
  145. }
  146. }
  147. }
  148. function index(point) {
  149. return point[0] * 2 + point[1] * (dx + 1) * 4;
  150. }
  151. function smoothLinear(ring, values, value) {
  152. ring.forEach(function(point) {
  153. var x = point[0],
  154. y = point[1],
  155. xt = x | 0,
  156. yt = y | 0,
  157. v1 = valid(values[yt * dx + xt]);
  158. if (x > 0 && x < dx && xt === x) {
  159. point[0] = smooth1(x, valid(values[yt * dx + xt - 1]), v1, value);
  160. }
  161. if (y > 0 && y < dy && yt === y) {
  162. point[1] = smooth1(y, valid(values[(yt - 1) * dx + xt]), v1, value);
  163. }
  164. });
  165. }
  166. contours.contour = contour;
  167. contours.size = function(_) {
  168. if (!arguments.length) return [dx, dy];
  169. var _0 = Math.floor(_[0]), _1 = Math.floor(_[1]);
  170. if (!(_0 >= 0 && _1 >= 0)) throw new Error("invalid size");
  171. return dx = _0, dy = _1, contours;
  172. };
  173. contours.thresholds = function(_) {
  174. return arguments.length ? (threshold = typeof _ === "function" ? _ : Array.isArray(_) ? constant(slice.call(_)) : constant(_), contours) : threshold;
  175. };
  176. contours.smooth = function(_) {
  177. return arguments.length ? (smooth = _ ? smoothLinear : noop, contours) : smooth === smoothLinear;
  178. };
  179. return contours;
  180. }
  181. // When computing the extent, ignore infinite values (as well as invalid ones).
  182. function finite(x) {
  183. return isFinite(x) ? x : NaN;
  184. }
  185. // Is the (possibly invalid) x greater than or equal to the (known valid) value?
  186. // Treat any invalid value as below negative infinity.
  187. function above(x, value) {
  188. return x == null ? false : +x >= value;
  189. }
  190. // During smoothing, treat any invalid value as negative infinity.
  191. function valid(v) {
  192. return v == null || isNaN(v = +v) ? -Infinity : v;
  193. }
  194. function smooth1(x, v0, v1, value) {
  195. const a = value - v0;
  196. const b = v1 - v0;
  197. const d = isFinite(a) || isFinite(b) ? a / b : Math.sign(a) / Math.sign(b);
  198. return isNaN(d) ? x : x + d - 0.5;
  199. }