空气污染扩散-高斯烟羽模型的克里金Kriging插值模拟地图实现

1、示例地址：GIS动态扩散模型-演示版 (simu.cool)

2、插件源码地址：https://github.com/oeo4b/kriging.js

3、安装插件

#安装克里金插值所需的插件
npm i @turf/turf# 
npm i @sakitam-gis/kriging

4、引入依赖

// 引入turf
import * as turf from "@turf/turf";

// 引入kriging.js
import kriging from "@sakitam-gis/kriging";

5、定义变量

      // 插值矢量图层组
      featureLayerGroup: undefined,
      imageLayerGroup: undefined,

6、预先定义DOM元素

  <div id="mapDiv" class="map">
    <!-- leaflet 地图容器 -->
    <canvas id="canvasMap" style="display: none"></canvas>
  </div>

7、渲染代码，使用Leaflet

 // 矢量图层组
         this.featureLayerGroup = new L.FeatureGroup()
          .addTo(this.map)
          .bringToFront();
        // 图像图层组
        this.imageLayerGroup = new L.FeatureGroup()
          .addTo(this.map)
          .bringToFront();
        // 完全透明
        var collection = turf.featureCollection(features);
        let bboxPolygon = turf.convex(collection); // 外包多边形范围
        let scope = L.geoJSON(bboxPolygon, {
          style: function () {
            return {
              fillColor: "6666ff",
              color: "red",
              weight: 2,
              opacity: 0,
              fillOpacity: 0,
            };
          },
        }).addTo(this.imageLayerGroup);

        this.map.fitBounds(scope.getBounds());
        var _self = this;
        //根据scope边界线，生成范围信息
        let xlim = [
          scope.getBounds()._southWest.lng,
          scope.getBounds()._northEast.lng,
        ];
        let ylim = [
          scope.getBounds()._southWest.lat,
          scope.getBounds()._northEast.lat,
        ];

        function loadkriging(points) {
          let canvas = document.getElementById("canvasMap");
          canvas.width = 2000;
          canvas.height = 2000;
          // 数量
          // let pointLength = points.features.length;
          let pointLength = points.length;
          let t = []; // 数值
          let x = []; // 经度
          let y = []; // 纬度
          // 加载点数过多的话，会出现卡顿
          for (let i = 0; i < pointLength; i++) {
            x.push(points[i].geometry.coordinates[0]);
            y.push(points[i].geometry.coordinates[1]);
            t.push(points[i].properties.dB);
          }

          // 克里金插值参数
          const params = {
            krigingModel: "exponential",
            krigingSigma2: 0,
            krigingAlpha: 100,
            canvasAlpha: 0.8, //canvas图层透明度
            colors: [
              "#2c75ab",
              "#a0cada",
              "#dfe0b6",
              "#dca25d",
              "#c81b1e",
            ],
            idx: [50, 55, 60, 65, 70],
          };
          // 对数据集进行训练
          let variogram = kriging.train(
            t,
            x,
            y,
            params.krigingModel,
            params.krigingSigma2,
            params.krigingAlpha
          );
          // 将插值范围封装成特定格式
          var collection = turf.featureCollection(points);
          let bboxPolygon = turf.convex(collection); // 外包多边形范围
          // 根据外包矩形范围生成外包矩形面Polygon
          // let bboxPolygon = turf.bboxPolygon(bbox);
          let positions = [];
          bboxPolygon.geometry.coordinates[0].forEach((v) => {
            positions.push([v[0], v[1]]);
          });
          // 将边界封装成特定的格式
          let range = [positions];
          // 最后一个参数是插值格点精度大小
          let grid = kriging.grid(range, variogram, (xlim[1] - xlim[0]) / 200);
          // 将得到的格网grid渲染至canvas上
          kriging.plot(
            canvas,
            grid,
            [xlim[0], xlim[1]],
            [ylim[0], ylim[1]],
            params.colors
          );
        }

        //将canvas对象转换成image的URL
        function returnImgae() {
          let mycanvas = document.getElementById("canvasMap");
          return mycanvas.toDataURL("image/png");
        }

        // 执行克里金插值函数
        loadkriging(features);

        let imageBounds = [
          [ylim[0], xlim[0]],
          [ylim[1], xlim[1]],
        ];
        L.imageOverlay(returnImgae(), imageBounds, { opacity: 0.8 }).addTo(
          this.imageLayerGroup
        );

8、离散的点如果过多，js会报异常，因此如果遇到这种情况，需要对点进行抽稀，或者修改替换Kriging.js源代码如下：

Array.prototype.rep = function (n) {
  var arrayn = new Array(n);
  var v = this[0];
  for (var i = 0; i < n; i++) {
    arrayn[i] = v;
  }
  return arrayn;
};

4、实现功能效果：