JS 离散点生成等高线图的几种方式
最近公司有需求，需要根据一组离散点来生成等高线图。我上网搜索了一些资料，然后自己也试了一下。我总结了三种实现方式：
-
turf + d3
-
turf + vjmap + d3
-
observable Plot (推荐)
接下来我详细介绍一下这三种不同的实现方式。
首先给大家看下我的数据结构：
| X | Y | Z |
|---|---|---|
| 43 | 190 | 56.86 |
| 51 | 423 | 55.06 |
| 58 | 171 | 55.84 |
| 58 | 478 | 74.29 |
| 63 | 50 | 56.56 |
| … | … | … |
turf + d3
基本思路：
① 将数据处理成地理格式
{
"type": "Feature",
"properties": {
"x": 565,
"y": 293,
"value": 96.9
},
"geometry": {
"type": "Point",
"coordinates": [
200.88888888888889,
104.17777777777778
]
}
}
② 通过 turf.interpolate (基于 IDW（反距离权重）)进行数据插值
③ 通过 turf.isobands，将插值数据计算出等边线，数据格式：
{
"type": "Feature",
"properties": {
"fill": "rgb(7, 117, 243)",
"value": "0-25"
},
"geometry": {
"type": "MultiPolygon",
"coordinates": [
[
[
[
257.14193519824784,
255.1806484413806
],
[
255.62813160239872,
255.8773511961486
],
[
257.14193519824784,
256.5056402838994
],
[
259.39336680305524,
255.8773511961486
],
[
257.14193519824784,
255.1806484413806
]
]
],
[
[
[
139.29531358573226,
219.80851022549035
],
[
138.98122511088943,
220.06321454991826
],
[
139.29531358573226,
220.31375106192775
],
[
139.70261352090898,
220.06321454991826
],
[
139.29531358573226,
219.80851022549035
]
]
],
[
[
[
205.77289500817696,
116.74637282143532
],
[
205.63455078271272,
116.8342324519602
],
[
205.77289500817696,
116.92095334929672
],
[
205.99880882366955,
116.8342324519602
],
[
205.77289500817696,
116.74637282143532
]
]
]
]
}
}
④ 通过 d3 生成svg图：
代码如下(Vue3)：
<template>
<div id="hello" style="height: 500px;width: 700px;">
<button class="but" @click="setPoint">显示统计点</button>
</div>
</template>
<script>
import * as turf from "@turf/turf"
import * as d3 from "d3"
import $ from "jquery"
export default {
name: 'turfContours',
data() {
return {
h: 440,
w: 400,
margin: 40,
hasPoint: true,
data: []
}
},
mounted() {
new Array(100).fill(0).forEach(e => {
let X=Math.floor(Math.random()*1000);
let Y=Math.floor(Math.random()*1000);
let Z=(Math.random()*100).toFixed(2);
this.data.push({
x: X,
y: Y,
z: Z
})
});
this.draw();
},
methods: {
draw() {
// 坐标系长高
let axX = this.w - this.margin * 2;
let axY = this.h - this.margin * 3;
const svg = d3.select('#hello')
.append("svg")
.attr("viewBox", [0, 0, this.w, this.h]);
// 图表
let g = svg.append('g')
.attr("transform", 'translate(' + this.margin + ',' + this.margin + ')');
// 值域区间
let val = d3.extent(this.data, d => d.z)
let valS = d3.nice(val[0], val[1], 10)
let breaks = d3.range(5).map(i => {
return d3.quantile(valS, i * 0.25)
})
// 添加颜色图例
let color = this.colorLine(svg, breaks);
// 坐标轴
let x = d3.scaleLinear()
.domain(d3.extent(this.data, (d) => d.x))
.nice() // 对齐
.rangeRound([0, axX]) // 起止位置
let xAxis = d3.axisTop(x)
let y = d3.scaleLinear()
.domain(d3.extent(this.data, (d) => d.y))
.nice()
.rangeRound([0, axY])
let yAxis = d3.axisLeft(y)
g.append("g")
.attr("transform", "translate(0, 0)")
.attr("class", "axisX")
.call(xAxis)
.selectAll("text")
.style("font-size", "4px")
g.append("g").attr("transform", "translate(0,0)")
.attr("class", "axisY")
.call(yAxis)
.selectAll("text")
.style("font-size", "4px")
// 根据坐标系,添加两端点
// 目的：端点没有值，插值的时候可能会在端点处显示有问题
let xS = d3.selectAll(".axisX text")._groups[0]
let yS = d3.selectAll(".axisY text")._groups[0]
let maxXs = d3.max(xS, d => Number(d.innerHTML))
let minXs = d3.min(xS, d => Number(d.innerHTML))
let maxYs = d3.max(yS, d => Number(d.innerHTML))
let minYs = d3.min(yS, d => Number(d.innerHTML))
let mPoint = [] // 对角顶点
let minDis = -1;
let maxDis = -1;
this.data.forEach(i => {
let distanceMin = turf.distance(turf.point([minXs, minYs]), turf.point([i.x, i.y]), {units: 'degrees'});
if (minDis === -1 || distanceMin < minDis) {
minDis = distanceMin;
mPoint[0] = {x: minXs, y: minYs, z: i.z}
}
let distanceMax = turf.distance(turf.point([i.x, i.y]), turf.point([maxXs, maxYs]), {units: 'degrees'});
if (maxDis === -1 || distanceMax < maxDis) {
maxDis = distanceMax;
mPoint[1] = {x: maxXs, y: maxYs, z: i.z}
}
})
// 原有点阵
// coordinates 通过比例计算出点对应在坐标系的位置
let featuresPoint = this.data.map(i => {
return {
type: "Feature",
properties: {
x: i.x,
y: i.y,
value: i.z
},
geometry: {
type: "Point",
coordinates: [(i.x - minXs) / (maxXs - minXs) * axX, (i.y - minYs) / (maxYs - minYs) * axY]
}
}
}
)
this.data = this.data.concat(mPoint)
// 计算点坐标
let features = this.data.map(i => {
return {
type: "Feature",
properties: {
x: i.x,
y: i.y,
value: i.z
},
geometry: {
type: "Point",
coordinates: [(i.x - minXs) / (maxXs - minXs) * axX, (i.y - minYs) / (maxYs - minYs) * axY]
}
}
}
)
// 离散点生成等值线图
let isobands = this.grid(features, color, breaks)
// 路径
g.append("g")
.selectAll("path")
.data(isobands.features)
.enter()
.append("path")
.attr("stroke", "white")
.attr("stroke-width", 0.5)
.attr("fill", d => d.properties.fill)
.attr("d", d3.geoPath());
this.showPoint(g, featuresPoint)
},
grid(features, color, breaks) {
let points = turf.featureCollection(features);
// 基于 IDW（反距离权重）算法的将数据插值为格点
let grid = turf.interpolate(points, 0.05, {
gridType: "points",
property: "value",
units: "degrees",
weight: 1
});
let isobands_options = {
zProperty: "value",
breaksProperties: []
};
breaks.map((v, i) => {
isobands_options.breaksProperties.push({'fill': color[i]});
return v
})
let isobands = turf.isobands(grid, breaks, isobands_options);
return isobands
},
colorLine(svg, breaks) {
// 显示渐变矩形条
let colorB = '#f50b36'
let colorA = '#0775f3'
let compute = d3.interpolate(colorA, colorB) // 返回一个函数
let color = [];
svg.selectAll("rect")
.data(d3.range(breaks.length - 1))
.enter()
.append("rect")
.attr("x", (d, i) => {
return this.w - (this.margin + 16 * (breaks.length - i));
})
.attr("y", 5)
.attr("width", 16)
.attr("height", 15)
.style("fill", (d, i) => {
let c = compute(i / (breaks.length - 1));
color.push(c)
return c;
});
// 数据初值
svg.selectAll("text")
.data(breaks)
.enter()
.append("text")
.attr("x", (d, i) => {
return this.w - (this.margin + 16 * (breaks.length - i) + 3);
})
.attr("y", 25)
.attr("font-size", 5)
.text(d => d)
return color;
},
showPoint(g, featuresPoint) {
// 展示点
g.append("g")
.selectAll("circle")
.data(featuresPoint)
.enter()
.append("circle")
.attr("cx", d => d.geometry.coordinates[0])
.attr("cy", d => d.geometry.coordinates[1])
.attr("r", 3)
.attr("class", 'point')
.attr("stroke", 'black')
.attr("stroke-width", 0.3)
.style("fill", "rgb(242,253,2)");
// 展示text
g.selectAll(".point")
.on("mouseover", function (d, i) {
g.append("text")
.attr("x", i.geometry.coordinates[0] + 5)
.attr("y", i.geometry.coordinates[1] + 2)
.attr("class", "tpo")
.text('(' + i.properties.x + ',' + i.properties.y + ',' + i.properties.value + ')')
.attr("font-size", 6)
d3.select(this)
.style("fill", "rgb(228,2,253)");
})
.on("mouseout", function (d, i) {
g.selectAll(".tpo")
.remove()
d3.select(this)
.style("fill", "rgb(242,253,2)");
})
},
setPoint() {
if (!this.hasPoint) {
this.hasPoint = true;
$('.but').text('隐藏统计点');
d3.select('#hello').selectAll(".point").attr("display", true);
} else {
this.hasPoint = false;
$('.but').text('显示统计点');
d3.select('#hello').selectAll(".point").attr("display", "none");
}
}
}
}
</script>
<!-- Add "scoped" attribute to limit CSS to this component only -->
<style scoped>
.but {
position: absolute;
top: 20px;
left: 70px;
}
</style>
注：
这种方式的优缺点 如下：
优点：自主可控，可以对svg图进行操作。
缺点：算法效果不理想，没有将值很好的区分；计算量大了之后速度会很慢。
turf + vjmap + d3
基本思路：
这个方案大体跟上一个方案类似，主要是针对上一种方案的缺点进行了优化。通过Vjmap的等高线图生成算法进行插值计算，生成等边线。
① 将数据处理成地理格式
{
"type": "Feature",
"properties": {
"x": 565,
"y": 293,
"value": 96.9
},
"geometry": {
"type": "Point",
"coordinates": [
200.88888888888889,
104.17777777777778
]
}
}
② 通过 vjmap.WorkerProxy(vjmap.vectorContour) (Vjmap等高线图)进行数据插值，得出等边线。
③ 通过 d3 生成svg图：
代码如下(Vue2)：
<template>
<div id="hello" style="height: 500px;width: 700px;">
<button class="but" @click="setPoint">显示统计点</button>
</div>
</template>
<script>
import * as d3 from "d3"
import $ from "jquery"
import vjmap from "vjmap";
export default {
name: 'turfContours',
data() {
return {
h: 440,
w: 400,
margin: 40,
hasPoint: true,
data: []
}
},
mounted() {
new Array(100).fill(0).forEach(e => {
let X=Math.floor(Math.random()*1000);
let Y=Math.floor(Math.random()*1000);
let Z=(Math.random()*100).toFixed(2);
this.data.push({
x: X,
y: Y,
z: Z
})
});
this.draw();
},
methods: {
async draw() {
// 坐标系长高
let axX = this.w - this.margin * 2;
let axY = this.h - this.margin * 3;
const svg = d3.select('#hello')
.append("svg")
.attr("viewBox", [0, 0, this.w, this.h]);
// 图表
let g = svg.append('g')
.attr("transform", 'translate(' + this.margin + ',' + this.margin + ')');
// 统计区间
let val = d3.extent(this.data, d => d.z)
let valS = d3.nice(val[0], val[1], 10)
let breaks = d3.range(5).map(i => {
return d3.quantile(valS, i * 0.25)
})
// 渐变颜色
let color = this.colorLine(svg, breaks);
// 坐标轴
let x = d3.scaleLinear()
.domain(d3.extent(this.data, (d) => d.x))
.nice() // 对齐
.rangeRound([0, axX]) // 起止位置
let xAxis = d3.axisTop(x)
let y = d3.scaleLinear()
.domain(d3.extent(this.data, (d) => d.y))
.nice()
.rangeRound([0, axY])
let yAxis = d3.axisLeft(y)
g.append("g")
.attr("transform", "translate(0, 0)")
.attr("class", "axisX")
.call(xAxis)
.selectAll("text")
.style("font-size", "4px")
g.append("g").attr("transform", "translate(0,0)")
.attr("class", "axisY")
.call(yAxis)
.selectAll("text")
.style("font-size", "4px")
// 根据坐标系,添加两端点
let xS = d3.selectAll(".axisX text")._groups[0]
let yS = d3.selectAll(".axisY text")._groups[0]
let maxXs = d3.max(xS, d => Number(d.innerHTML))
let minXs = d3.min(xS, d => Number(d.innerHTML))
let maxYs = d3.max(yS, d => Number(d.innerHTML))
let minYs = d3.min(yS, d => Number(d.innerHTML))
let mPoint = [] // 对角顶点
let minDis = -1;
let maxDis = -1;
this.data.forEach(i => {
let distanceMin = turf.distance(turf.point([minXs, minYs]), turf.point([i.x, i.y]), {units: 'degrees'});
if (minDis === -1 || distanceMin < minDis) {
minDis = distanceMin;
mPoint[0] = {x: minXs, y: minYs, z: i.z}
}
let distanceMax = turf.distance(turf.point([i.x, i.y]), turf.point([maxXs, maxYs]), {units: 'degrees'});
if (maxDis === -1 || distanceMax < maxDis) {
maxDis = distanceMax;
mPoint[1] = {x: maxXs, y: maxYs, z: i.z}
}
})
// 原有点阵
let featuresPoint = this.data.map(i => {
return {
type: "Feature",
properties: {
x: i.x,
y: i.y,
value: i.z
},
geometry: {
type: "Point",
coordinates: [(i.x - minXs) / (maxXs - minXs) * axX, (i.y - minYs) / (maxYs - minYs) * axY]
}
}
}
)
this.data = this.data.concat(mPoint)
// 计算点坐标
let features = this.data.map(i => {
return {
type: "Feature",
properties: {
x: i.x,
y: i.y,
value: i.z
},
geometry: {
type: "Point",
coordinates: [(i.x - minXs) / (maxXs - minXs) * axX, (i.y - minYs) / (maxYs - minYs) * axY]
}
}
}
)
// 离散点生成等值线图
let isobands = await this.grid2(features, "value", breaks, color)
// 路径
g.append("g")
.selectAll("path")
.data(isobands.features)
.enter()
.append("path")
.attr("stroke", "white")
.attr("stroke-width", 0.5)
.attr("fill", d => d.properties.fill)
.attr("d", d3.geoPath());
this.showPoint(g, featuresPoint)
},
async grid2(dataset, propField, contours, color) {
// 唯杰地图
let createContourWorker = vjmap.WorkerProxy(vjmap.vectorContour);
let points = turf.featureCollection(dataset);
let {grid, contour, variogram} = await createContourWorker(points, propField, contours, {
model: 'exponential', // 'exponential','gaussian','spherical'，三选一，默认exponential
sigma2: 0, // sigma2是σ²，对应高斯过程的方差参数，也就是这组数据z的距离，方差参数σ²的似然性反映了高斯过程中的误差，并应手动设置。一般设置为 0 ，其他数值设了可能会出空白图
alpha: 50, // [如果绘制不出来，修改此值，可以把此值改小] Alpha α对应方差函数的先验值，此参数可能控制钻孔扩散范围,越小范围越大,少量点效果明显，但点多了且分布均匀以后改变该数字即基本无效果了，默认设置为100
// extent: extent // 如果要根据数据范围自动生成此范围，则无需传此参数
}, []);
contour.features.forEach(a => {
let val = a.properties.value
let ind = color.length - 1
contours.forEach((k, i) => {
if (k <= val) {
ind = i
}
})
a.properties.fill = color[ind]
})
return contour;
},
colorLine(svg, breaks) {
// 显示渐变矩形条
let colorB = '#f50b36'
let colorA = '#0775f3'
let compute = d3.interpolate(colorA, colorB) // 返回一个函数
let color = [];
svg.selectAll("rect")
.data(d3.range(breaks.length - 1))
.enter()
.append("rect")
.attr("x", (d, i) => {
return this.w - (this.margin + 16 * (breaks.length - i));
})
.attr("y", 5)
.attr("width", 16)
.attr("height", 15)
.style("fill", (d, i) => {
let c = compute(i / (breaks.length - 1));
color.push(c)
return c;
});
// 数据初值
svg.selectAll("text")
.data(breaks)
.enter()
.append("text")
.attr("x", (d, i) => {
return this.w - (this.margin + 16 * (breaks.length - i) + 3);
})
.attr("y", 25)
.attr("font-size", 5)
.text(d => d)
return color;
},
showPoint(g, featuresPoint) {
// 展示点
g.append("g")
.selectAll("circle")
.data(featuresPoint)
.enter()
.append("circle")
.attr("cx", d => d.geometry.coordinates[0])
.attr("cy", d => d.geometry.coordinates[1])
.attr("r", 3)
.attr("class", 'point')
.attr("stroke", 'black')
.attr("stroke-width", 0.3)
.style("fill", "rgb(242,253,2)");
// 展示text
g.selectAll(".point")
.on("mouseover", function (d, i) {
g.append("text")
.attr("x", i.geometry.coordinates[0] + 5)
.attr("y", i.geometry.coordinates[1] + 2)
.attr("class", "tpo")
.text('(' + i.properties.x + ',' + i.properties.y + ',' + i.properties.value + ')')
.attr("font-size", 6)
d3.select(this)
.style("fill", "rgb(228,2,253)");
})
.on("mouseout", function (d, i) {
g.selectAll(".tpo")
.remove()
d3.select(this)
.style("fill", "rgb(242,253,2)");
})
},
setPoint() {
if (!this.hasPoint) {
this.hasPoint = true;
$('.but').text('隐藏统计点');
d3.select('#hello').selectAll(".point").attr("display", true);
} else {
this.hasPoint = false;
$('.but').text('显示统计点');
d3.select('#hello').selectAll(".point").attr("display", "none");
}
}
}
}
</script>
<!-- Add "scoped" attribute to limit CSS to this component only -->
<style scoped>
.but {
position: absolute;
top: 20px;
left: 70px;
}
</style>
注：
这种方式的优缺点 如下：
优点：自主可控，可以对svg图进行操作；相对于第一种方式，等高线的计算更加准确；
缺点：计算量大了之后速度会很慢（2000个点需要计算30~40s）
observable Plot (推荐)
这个方案是最近才有的，因为官方给的案例5月份才有。之前 d3的等高线demo一直不是我想要的处理方式，直到最近看到的这个神器 observable Plot，发现简直是鬼斧神工。
基本思路：
①直接照抄 observable Plot 的API吧 (Plot)
效果展示：
代码如下(Vue3)：
<template>
<div id="hello" style="height: 500px;width: 700px;">
</div>
</template>
<script>
import * as d3 from "d3"
import * as Plot from "@observablehq/plot";
import $ from "jquery"
export default {
name: "observable",
data() {
return {
data: []
}
},
mounted() {
new Array(100).fill(0).forEach(e => {
let X=Math.floor(Math.random()*1000);
let Y=Math.floor(Math.random()*1000);
let Z=Number((Math.random()*100).toFixed(2));
this.data.push({
x: X,
y: Y,
z: Z
})
});
this.draw();
},
methods: {
draw() {
let ps = Plot.plot({
color: {legend: true,scheme: "sinebow"},
marks: [
Plot.contour(this.data, {x: "x", y: "y", fill: "z", blur: 5}),
Plot.dot(this.data, {x: "x", y: "y", channels: {z: "z"}, tip: true}),
Plot.frame(),
Plot.text(["图1"], {frameAnchor: "top-right",dy: -15,fontSize:15,fontWeight:900}),
]
})
$("#hello").append(ps);
},
}
}
</script>
<style scoped>
</style>
注：
这种方式的优缺点 如下：
缺点：需要按照封装好的API进行绘制，有的时候达不到想要的效果(eg：颜色图例，svg图自定义)
优点：算法一流，等高线的计算更加准确，快速（真NB，5K个点1s出结果）；
支持原创！！
欢迎加入讨论！！！
573
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
