index.js 6.0 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160
  1. 'use strict';
  2. var meta = require('@turf/meta');
  3. var invariant = require('@turf/invariant');
  4. var helpers = require('@turf/helpers');
  5. var centerMean = require('@turf/center-mean');
  6. var pointsWithinPolygon = require('@turf/points-within-polygon');
  7. var ellipse = require('@turf/ellipse');
  8. function _interopDefaultLegacy (e) { return e && typeof e === 'object' && 'default' in e ? e : { 'default': e }; }
  9. var centerMean__default = /*#__PURE__*/_interopDefaultLegacy(centerMean);
  10. var pointsWithinPolygon__default = /*#__PURE__*/_interopDefaultLegacy(pointsWithinPolygon);
  11. var ellipse__default = /*#__PURE__*/_interopDefaultLegacy(ellipse);
  12. /**
  13. * Takes a {@link FeatureCollection} and returns a standard deviational ellipse,
  14. * also known as a “directional distribution.” The standard deviational ellipse
  15. * aims to show the direction and the distribution of a dataset by drawing
  16. * an ellipse that contains about one standard deviation’s worth (~ 70%) of the
  17. * data.
  18. *
  19. * This module mirrors the functionality of [Directional Distribution](http://desktop.arcgis.com/en/arcmap/10.3/tools/spatial-statistics-toolbox/directional-distribution.htm)
  20. * in ArcGIS and the [QGIS Standard Deviational Ellipse Plugin](http://arken.nmbu.no/~havatv/gis/qgisplugins/SDEllipse/)
  21. *
  22. * **Bibliography**
  23. *
  24. * • Robert S. Yuill, “The Standard Deviational Ellipse; An Updated Tool for
  25. * Spatial Description,” _Geografiska Annaler_ 53, no. 1 (1971): 28–39,
  26. * doi:{@link https://doi.org/10.2307/490885|10.2307/490885}.
  27. *
  28. * • Paul Hanly Furfey, “A Note on Lefever’s “Standard Deviational Ellipse,”
  29. * _American Journal of Sociology_ 33, no. 1 (1927): 94—98,
  30. * doi:{@link https://doi.org/10.1086/214336|10.1086/214336}.
  31. *
  32. *
  33. * @name standardDeviationalEllipse
  34. * @param {FeatureCollection<Point>} points GeoJSON points
  35. * @param {Object} [options={}] Optional parameters
  36. * @param {string} [options.weight] the property name used to weight the center
  37. * @param {number} [options.steps=64] number of steps for the polygon
  38. * @param {Object} [options.properties={}] properties to pass to the resulting ellipse
  39. * @returns {Feature<Polygon>} an elliptical Polygon that includes approximately 1 SD of the dataset within it.
  40. * @example
  41. *
  42. * var bbox = [-74, 40.72, -73.98, 40.74];
  43. * var points = turf.randomPoint(400, {bbox: bbox});
  44. * var sdEllipse = turf.standardDeviationalEllipse(points);
  45. *
  46. * //addToMap
  47. * var addToMap = [points, sdEllipse];
  48. *
  49. */
  50. function standardDeviationalEllipse(points, options) {
  51. // Optional params
  52. options = options || {};
  53. if (!helpers.isObject(options)) throw new Error("options is invalid");
  54. var steps = options.steps || 64;
  55. var weightTerm = options.weight;
  56. var properties = options.properties || {};
  57. // Validation:
  58. if (!helpers.isNumber(steps)) throw new Error("steps must be a number");
  59. if (!helpers.isObject(properties)) throw new Error("properties must be a number");
  60. // Calculate mean center & number of features:
  61. var numberOfFeatures = meta.coordAll(points).length;
  62. var meanCenter = centerMean__default['default'](points, { weight: weightTerm });
  63. // Calculate angle of rotation:
  64. // [X, Y] = mean center of all [x, y].
  65. // theta = arctan( (A + B) / C )
  66. // A = sum((x - X)^2) - sum((y - Y)^2)
  67. // B = sqrt(A^2 + 4(sum((x - X)(y - Y))^2))
  68. // C = 2(sum((x - X)(y - Y)))
  69. var xDeviationSquaredSum = 0;
  70. var yDeviationSquaredSum = 0;
  71. var xyDeviationSum = 0;
  72. meta.featureEach(points, function (point) {
  73. var weight = point.properties[weightTerm] || 1;
  74. var deviation = getDeviations(invariant.getCoords(point), invariant.getCoords(meanCenter));
  75. xDeviationSquaredSum += Math.pow(deviation.x, 2) * weight;
  76. yDeviationSquaredSum += Math.pow(deviation.y, 2) * weight;
  77. xyDeviationSum += deviation.x * deviation.y * weight;
  78. });
  79. var bigA = xDeviationSquaredSum - yDeviationSquaredSum;
  80. var bigB = Math.sqrt(Math.pow(bigA, 2) + 4 * Math.pow(xyDeviationSum, 2));
  81. var bigC = 2 * xyDeviationSum;
  82. var theta = Math.atan((bigA + bigB) / bigC);
  83. var thetaDeg = (theta * 180) / Math.PI;
  84. // Calculate axes:
  85. // sigmaX = sqrt((1 / n - 2) * sum((((x - X) * cos(theta)) - ((y - Y) * sin(theta)))^2))
  86. // sigmaY = sqrt((1 / n - 2) * sum((((x - X) * sin(theta)) - ((y - Y) * cos(theta)))^2))
  87. var sigmaXsum = 0;
  88. var sigmaYsum = 0;
  89. var weightsum = 0;
  90. meta.featureEach(points, function (point) {
  91. var weight = point.properties[weightTerm] || 1;
  92. var deviation = getDeviations(invariant.getCoords(point), invariant.getCoords(meanCenter));
  93. sigmaXsum +=
  94. Math.pow(
  95. deviation.x * Math.cos(theta) - deviation.y * Math.sin(theta),
  96. 2
  97. ) * weight;
  98. sigmaYsum +=
  99. Math.pow(
  100. deviation.x * Math.sin(theta) + deviation.y * Math.cos(theta),
  101. 2
  102. ) * weight;
  103. weightsum += weight;
  104. });
  105. var sigmaX = Math.sqrt((2 * sigmaXsum) / weightsum);
  106. var sigmaY = Math.sqrt((2 * sigmaYsum) / weightsum);
  107. var theEllipse = ellipse__default['default'](meanCenter, sigmaX, sigmaY, {
  108. units: "degrees",
  109. angle: thetaDeg,
  110. steps: steps,
  111. properties: properties,
  112. });
  113. var pointsWithinEllipse = pointsWithinPolygon__default['default'](
  114. points,
  115. helpers.featureCollection([theEllipse])
  116. );
  117. var standardDeviationalEllipseProperties = {
  118. meanCenterCoordinates: invariant.getCoords(meanCenter),
  119. semiMajorAxis: sigmaX,
  120. semiMinorAxis: sigmaY,
  121. numberOfFeatures: numberOfFeatures,
  122. angle: thetaDeg,
  123. percentageWithinEllipse:
  124. (100 * meta.coordAll(pointsWithinEllipse).length) / numberOfFeatures,
  125. };
  126. theEllipse.properties.standardDeviationalEllipse = standardDeviationalEllipseProperties;
  127. return theEllipse;
  128. }
  129. /**
  130. * Get x_i - X and y_i - Y
  131. *
  132. * @private
  133. * @param {Array} coordinates Array of [x_i, y_i]
  134. * @param {Array} center Array of [X, Y]
  135. * @returns {Object} { x: n, y: m }
  136. */
  137. function getDeviations(coordinates, center) {
  138. return {
  139. x: coordinates[0] - center[0],
  140. y: coordinates[1] - center[1],
  141. };
  142. }
  143. module.exports = standardDeviationalEllipse;
  144. module.exports.default = standardDeviationalEllipse;