/**
 * Quadrat analysis lays a set of equal-size areas(quadrat) over the study area and counts
 * the number of features in each quadrat and creates a frequency table.
 * The table lists the number of quadrats containing no features,
 * the number containing one feature, two features, and so on,
 * all the way up to the quadrat containing the most features.
 * The method then creates the frequency table for the random distribution, usually based on a Poisson distribution.
 * The method uses the distribution to calculate the probability for 0 feature occuring,
 * 1 feature occuring, 2 features, and so on,
 * and lists these probabilities in the frequency table.
 * By comparing the two frequency tables, you can see whether the features create a pattern.
 * If the table for the observed distribution has more quadrats containing many features than the
 * table for the random distribution dose, then the features create a clustered pattern.
 *
 * It is hard to judge the frequency tables are similar or different just by looking at them.
 * So, we can use serval statistical tests to find out how much the frequency tables differ.
 * We use Kolmogorov-Smirnov test.This method calculates cumulative probabilities for both distributions,
 * and then compares the cumulative probabilities at each class level and selects the largest absolute difference D.
 * Then, the test compares D to the critical value for a confidence level you specify.
 * If D is greater than the critical value, the difference between  the observed distribution and
 * the random distribution is significant. The greater the value the bigger the difference.
 *
 * Traditionally, squares are used for the shape of the quadrats, in a regular grid(square-grid).
 * Some researchers suggest that the quadrat size equal twice the size of mean area per feature,
 * which is simply the area of the study area divided by the number of features.
 *
 *
 * @name quadratAnalysis
 * @param {FeatureCollection<Point>} pointFeatureSet point set to study
 * @param {Object} [options={}] optional parameters
 * @param {bbox} [options.studyBbox] bbox representing the study area
 * @param {number} [options.confidenceLevel=20] a confidence level.
 * The unit is percentage . 5 means 95%, value must be in {@link K_TABLE}
 * @returns {Object} result {@link QuadratAnalysisResult}
 * @example
 *
 * var bbox = [-65, 40, -63, 42];
 * var dataset = turf.randomPoint(100, { bbox: bbox });
 * var result = turf.quadratAnalysis(dataset);
 *
 */
export default function quadratAnalysis(pointFeatureSet: FeatureCollection<Point>, options: {
  studyBbox?: [number, number, number, number]
  confidenceLevel?: 20 | 15 | 10 | 5 | 2 | 1,
}): QuadratAnalysisResult {

  options = options || {};
  const studyBbox = options.studyBbox || turfBBox(pointFeatureSet);
  const confidenceLevel = options.confidenceLevel || 20;
  const points = pointFeatureSet.features;

  // create square-grid
  const numOfPoints = points.length;
  const sizeOfArea = area(bboxPolygon(studyBbox));
  const lengthOfSide = Math.sqrt((sizeOfArea / numOfPoints) * 2);
  const grid = squareGrid(studyBbox, lengthOfSide, {
    units: "meters",
  });
  const quadrats = grid.features;

  // count the number of features in each quadrat
  const quadratIdDict: {[key: string]: {box: BBox, cnt: number}} = {};
  for (let i = 0; i < quadrats.length; i++) {
    quadratIdDict[i] = {
      box: turfBBox(quadrats[i]),
      cnt: 0,
    };
  }

  let sumOfPoint = 0;
  for (const pt of points) {
    for (const key of Object.keys(quadratIdDict)) {
      const box = quadratIdDict[key].box;
      if (inBBox(getCoord(pt), box)) {
        quadratIdDict[key].cnt += 1;
        sumOfPoint += 1;
        break;
      }
    }
  }

  // the most amount of features in quadrat
  let maxCnt = 0;
  for (const key of Object.keys(quadratIdDict)) {
    const cnt = quadratIdDict[key].cnt;
    if (cnt > maxCnt) {
      maxCnt = cnt;
    }
  }

  const expectedDistribution = [];
  const numOfQuadrat = Object.keys(quadratIdDict).length;
  const lambda = sumOfPoint / numOfQuadrat;

  // get the cumulative probability of the random distribution
  let cumulativeProbility = 0.0;
  for (let x = 0; x < maxCnt + 1; x++) {
    cumulativeProbility += Math.exp(-lambda) * Math.pow(lambda, x) / factorial(x);
    expectedDistribution.push(cumulativeProbility);
  }
//.........这里部分代码省略.........

/**
 * Is Polygon2 in Polygon1
 * Only takes into account outer rings
 *
 * @private
 * @param {Geometry|Feature<Polygon>} feature1 Polygon1
 * @param {Geometry|Feature<Polygon>} feature2 Polygon2
 * @returns {boolean} true/false
 */
function isPolyInPoly(feature1, feature2) {
    var poly1Bbox = calcBbox(feature1);
    var poly2Bbox = calcBbox(feature2);
    if (!doBBoxOverlap(poly1Bbox, poly2Bbox)) {
        return false;
    }
    for (var i = 0; i < feature2.coordinates[0].length; i++) {
        if (!booleanPointInPolygon(feature2.coordinates[0][i], feature1)) {
            return false;
        }
    }
    return true;
}

function isLineInPoly(polygon, linestring) {
    var output = false;
    var i = 0;

    var polyBbox = calcBbox(polygon);
    var lineBbox = calcBbox(linestring);
    if (!doBBoxOverlap(polyBbox, lineBbox)) {
        return false;
    }
    for (i; i < linestring.coordinates.length - 1; i++) {
        var midPoint = getMidpoint(linestring.coordinates[i], linestring.coordinates[i + 1]);
        if (booleanPointInPolygon({type: 'Point', coordinates: midPoint}, polygon, { ignoreBoundary: true })) {
            output = true;
            break;
        }
    }
    return output;
}

export function isLineInPoly(polygon: Polygon, linestring: LineString) {
    let output = false;
    let i = 0;

    const polyBbox = calcBbox(polygon);
    const lineBbox = calcBbox(linestring);
    if (!doBBoxOverlap(polyBbox, lineBbox)) {
        return false;
    }
    for (i; i < linestring.coordinates.length - 1; i++) {
        const midPoint = getMidpoint(linestring.coordinates[i], linestring.coordinates[i + 1]);
        if (booleanPointInPolygon({type: "Point", coordinates: midPoint}, polygon, { ignoreBoundary: true })) {
            output = true;
            break;
        }
    }
    return output;
}

/**
 * Takes a {@link Feature} or {@link FeatureCollection} and returns the absolute center point of all features.
 *
 * @name center
 * @param {GeoJSON} geojson GeoJSON to be centered
 * @param {Object} [options={}] Optional parameters
 * @param {Object} [options.properties={}] Translate GeoJSON Properties to Point
 * @param {Object} [options.bbox={}] Translate GeoJSON BBox to Point
 * @param {Object} [options.id={}] Translate GeoJSON Id to Point
 * @returns {Feature<Point>} a Point feature at the absolute center point of all input features
 * @example
 * var features = turf.points([
 *   [-97.522259, 35.4691],
 *   [-97.502754, 35.463455],
 *   [-97.508269, 35.463245]
 * ]);
 *
 * var center = turf.center(features);
 *
 * //addToMap
 * var addToMap = [features, center]
 * center.properties['marker-size'] = 'large';
 * center.properties['marker-color'] = '#000';
 */
function center<P = Properties>(
    geojson: AllGeoJSON,
    options: {properties?: P, bbox?: BBox, id?: Id } = {}
): Feature<Point, P> {
    const ext = bbox(geojson);
    const x = (ext[0] + ext[2]) / 2;
    const y = (ext[1] + ext[3]) / 2;
    return point([x, y], options.properties, options);
}

export function isPolyInPoly(feature1: Feature<Polygon>|Polygon, feature2: Feature<Polygon>|Polygon) {
    // Handle Nulls
    if (feature1.type === "Feature" && feature1.geometry === null) { return false; }
    if (feature2.type === "Feature" && feature2.geometry === null) { return false; }

    const poly1Bbox = calcBbox(feature1);
    const poly2Bbox = calcBbox(feature2);
    if (!doBBoxOverlap(poly1Bbox, poly2Bbox)) {
        return false;
    }

    const coords = getGeom(feature2).coordinates;
    for (const ring of coords) {
        for (const coord of ring) {
            if (!booleanPointInPolygon(coord, feature1)) {
                return false;
            }
        }
    }
    return true;
}

        coords.forEach(function (coord) {
            if (autoComplete) coord = autoCompleteCoords(coord);

            // Largest LineString to be placed in the first position of the coordinates array
            if (orderCoords) {
                var area = calculateArea(turfBBox(lineString(coord)));
                if (area > largestArea) {
                    multiCoords.unshift(coord);
                    largestArea = area;
                } else multiCoords.push(coord);
            } else {
                multiCoords.push(coord);
            }
        });

    polygons.features.forEach(function (poly) {

        if (!poly.properties) {
            poly.properties = {};
        }
        var bbox = turfbbox(poly);
        var potentialPoints = rtree.search({minX: bbox[0], minY: bbox[1], maxX: bbox[2], maxY: bbox[3]});
        var values = [];
        potentialPoints.forEach(function (pt) {
            if (booleanPointInPolygon([pt.minX, pt.minY], poly)) {
                values.push(pt.property);
            }
        });

        poly.properties[outProperty] = values;
    });

/**
 * Nearest Neighbor Analysis calculates an index based the average distances
 * between points in the dataset, thereby providing inference as to whether the
 * data is clustered, dispersed, or randomly distributed within the study area.
 *
 * It returns a {@link Feature<Polygon>} of the study area, with the results of
 * the analysis attached as part of of the `nearestNeighborAnalysis` property
 * of the study area's `properties`. The attached
 * [_z_-score](https://en.wikipedia.org/wiki/Standard_score) indicates how many
 * standard deviations above or below the expected mean distance the data's
 * observed mean distance is. The more negative, the more clustered. The more
 * positive, the more evenly dispersed. A _z_-score between -2 and 2 indicates
 * a seemingly random distribution. That is, within _p_ of less than 0.05, the
 * distribution appears statistically significantly neither clustered nor
 * dispersed.
 *
 * **Remarks**
 *
 * - Though the analysis will work on any {@link FeatureCollection} type, it
 * works best with {@link Point} collections.
 *
 * - This analysis is _very_ sensitive to the study area provided.
 * If no {@link Feature<Polygon>} is passed as the study area, the function draws a box
 * around the data, which may distort the findings. This analysis works best
 * with a bounded area of interest within with the data is either clustered,
 * dispersed, or randomly distributed. For example, a city's subway stops may
 * look extremely clustered if the study area is an entire state. On the other
 * hand, they may look rather evenly dispersed if the study area is limited to
 * the city's downtown.
 *
 * **Bibliography**
 *
 * Philip J. Clark and Francis C. Evans, “Distance to Nearest Neighbor as a
 * Measure of Spatial Relationships in Populations,” _Ecology_ 35, no. 4
 * (1954): 445–453, doi:[10.2307/1931034](http://doi.org/10.2307/1931034).
 *
 * @name nearestNeighborAnalysis
 * @param {FeatureCollection<any>} dataset FeatureCollection (pref. of points) to study
 * @param {Object} [options={}] Optional parameters
 * @param {Feature<Polygon>} [options.studyArea] polygon representing the study area
 * @param {string} [options.units='kilometers'] unit of measurement for distances and, squared, area.
 * @param {Object} [options.properties={}] properties
 * @returns {Feature<Polygon>} A polygon of the study area or an approximation of one.
 * @example
 * var bbox = [-65, 40, -63, 42];
 * var dataset = turf.randomPoint(100, { bbox: bbox });
 * var nearestNeighborStudyArea = turf.nearestNeighborAnalysis(dataset);
 *
 * //addToMap
 * var addToMap = [dataset, nearestNeighborStudyArea];
 */
function nearestNeighborAnalysis(dataset: FeatureCollection<any>, options?: {
    studyArea?: Feature<Polygon>;
    units?: Units;
    properties?: Properties;
}): NearestNeighborStudyArea {
    // Optional params
    options = options || {};
    const studyArea = options.studyArea || bboxPolygon(bbox(dataset));
    const properties = options.properties || {};
    const units = options.units || 'kilometers';

    const features: Array<Feature<Point>> = [];
    featureEach(dataset, (feature) => {
        features.push(centroid(feature));
    });
    const n = features.length;
    const observedMeanDistance = features.map((feature, index) => {
        const otherFeatures = featureCollection<Point>(features.filter((f, i) => {
            return i !== index;
        }));
        // Have to add the ! to make typescript validation pass
        // see https://stackoverflow.com/a/40350534/1979085
        return distance(feature, nearestPoint(feature, otherFeatures).geometry!.coordinates, {units});
    }).reduce((sum, value) => { return sum + value; }, 0) / n;

    const populationDensity = n / convertArea(area(studyArea), 'meters', units);
    const expectedMeanDistance = 1 / (2 * Math.sqrt(populationDensity));
    const variance = 0.26136 / (Math.sqrt(n * populationDensity));
    properties.nearestNeighborAnalysis = {
        units: units,
        arealUnits: units + '²',
        observedMeanDistance: observedMeanDistance,
        expectedMeanDistance: expectedMeanDistance,
        nearestNeighborIndex: observedMeanDistance / expectedMeanDistance,
        numberOfPoints: n,
        zScore: (observedMeanDistance - expectedMeanDistance) / variance,
    };
    studyArea.properties = properties;

    return studyArea as NearestNeighborStudyArea;
}