氣象netCDF數據可視化分析
版權聲明：本文為博主原創文章，遵循 CC 4.0 BY-SA 版權協議，轉載請附上原文出處鏈接和本聲明。
前言
- NetCDF(network Common Data Form)網絡通用數據格式是由美國大學大氣研究協會（University Corporation for Atmospheric Research，UCAR)的Unidata項目科學家針對科學數據的特點開發的，是一種面向數組型並適於網絡共享的數據的描述和編碼標准。
對程序員來說，它和zip、jpeg、bmp文件格式類似，都是一種文件格式的標准。netcdf文件開始的目的是用於存儲氣象科學中的數據，現在已經成為許多數據采集軟件的生成文件的格式。
特點：NetCDF文件是自描述的二進制數據格式，即自帶描述屬性信息。通常包含了變量、維度和屬性，變量包含了維度、屬性（如數據單位）信息及變量的值。維度部分記錄的是每個變量的維度名及長度。屬性部分包含了一些額外信息，比如文件創建者等。
- 很多工具都可以處理NetCDF文件，比如MATLAB，Python，Java，NCL，GrADS，CDO，NCO，Panoply，ArcMap等等。NetCDF文件數據下載地址
- 這里主要講一下如何利用D3在前端處理NetCDF文件進行可視化分析。
- 核心代碼如下
<script>
//--------------------------------------
// 縮放控制
function zoomed() {
var transform = d3.event.transform;
projection.scale(scale * transform.k);
updatePaths(svg);
}
//--------------------------------------
function dragstarted() {
v0 = versor.cartesian(projection.invert(d3.mouse(this)));
r0 = projection.rotate();
q0 = versor(r0);
}
//--------------------------------------
// 拖拽控制
function dragged(d) {
var v1 = versor.cartesian(projection.rotate(r0).invert(d3.mouse(this))),
q1 = versor.multiply(q0, versor.delta(v0, v1)),
r1 = versor.rotation(q1);
projection.rotate(r1);
updatePaths(svg);
}
//--------------------------------------
function updatePaths(svg) {
svg.forEach(function(e) {
e.selectAll('path.contours').attr("d", geoPath);
e.selectAll('path.graticule').attr('d', geoPath);
e.selectAll('path.land').attr('d', geoPath);
});
}
//--------------------------------------
function createMap(id, values, range) {
var svg = d3.select('body').select(id).append('svg')
.attr('width', width)
.attr('height', height);
var group = svg.append("g").datum([]);
var extent = d3.extent(values);
// 顏色插值
var color = d3.scaleSequential(d3.interpolatePlasma)
//.domain(d3.extent(values));
.domain(range);
// console.log(d3.extent(values));
// 生成等值線
var contours = d3.contours()
.thresholds(d3.range(Math.floor(extent[0]/delta)*delta, Math.ceil(extent[1]/delta)*delta, delta))
.smooth(true)
.size([isize, jsize]);
// 對生成的等值線進行填色
group
//.attr("class", "contour-stroke")
.selectAll("path")
.data(contours(values).map(invert))
.enter().append("path")
.attr('class', 'contours')
.attr("fill", function(d) { return color(d.value); })
.attr("d", geoPath);
group.append('path')
.datum(graticule)
.attr('class', 'graticule')
.attr('d', geoPath);
group.append("path")
.datum(world)
.attr("class", "land")
.attr("d", geoPath);
// zoom on svg; drag on group
group.call(d3.drag().on('start', dragstarted)
.on('drag', dragged));
svg.call(d3.zoom().on('zoom', zoomed));
return svg;
}
//==========================================
function invert(d) {
var shared = {};
var p = {
type: "Polygon",
coordinates: d3.merge(d.coordinates.map(function(polygon) {
return polygon.map(function(ring) {
return ring.map(function(point) {
return [point[0] / isize * 360 - 180, 90 - point[1] / jsize * 180];
}).reverse();
});
}))
};
// Record the y-intersections with the antimeridian.
p.coordinates.forEach(function(ring) {
ring.forEach(function(p) {
if (p[0] === -180 || p[0] === 180) {
shared[p[1]] |= p[0] === -180 ? 1 : 2;
}
});
});
// Offset any unshared antimeridian points to prevent their stitching.
p.coordinates.forEach(function(ring) {
ring.forEach(function(p) {
if ((p[0] === -180 || p[0] === 180) && shared[p[1]] !== 3) {
p[0] = p[0] === -180 ? -179.9995 : 179.9995;
}
});
});
p = d3.geoStitch(p);
// If the MultiPolygon is empty, treat it as the Sphere.
return p.coordinates.length
? {type: "Polygon", coordinates: p.coordinates, value: d.value}
: {type: "Sphere", value: d.value};
}
//==========================================
function reverseVar(values) {
values = nj.array(values).reshape(jsize,isize);
values = values.slice([null, null, -1],null);
values = values.flatten().tolist();
return values;
}
//==========================================
var svg = [];
var world;
var graticule;
var width = 400,
height = 400,
scale = 200,
origin = {x: 55, y: -40};
var v0, // Mouse position in Cartesian coordinates at start of drag gesture.
r0, // Projection rotation as Euler angles at start.
q0; // Projection rotation as versor at start.
// 正交投影
var projection = d3.geoOrthographic()
.scale(scale)
.translate([width/2, height/2])
.rotate([origin.x, origin.y])
.center([0, 0]);
// 確定投影坐標系
var geoPath = d3.geoPath()
.projection(projection);
var min = -12;
var max = 12;
var delta = 2;
var nbLevels = Math.abs(max-min)/delta + 1;
var color = d3.scaleSequential(d3.interpolatePlasma)
.domain([min,max]);
//==========================================
var urlpath = "navy_winds_2.nc"
var reader;
var isize, jsize;
// 讀取netCDF文件數據
var oReq = new XMLHttpRequest();
oReq.open("GET", urlpath, true);
oReq.responseType = "blob";
oReq.onload = function(oEvent) {
var blob = oReq.response;
reader_url = new FileReader();
reader_url.onload = function(e) {
//====================================================================================
reader = new netcdfjs(this.result);
isize = reader.dimensions[0].size;
jsize = reader.dimensions[1].size;
var dim0Name = reader.dimensions[0].name;
var dim1Name = reader.dimensions[1].name;
axis0 = reader.getDataVariable(dim0Name);
axis1 = reader.getDataVariable(dim1Name);
var valuesVar1 = reader.getDataVariable('UWND');
valuesVar1 = reverseVar(valuesVar1);
var valuesVar2 = reader.getDataVariable('VWND');
valuesVar2 = reverseVar(valuesVar2);
range = [-12, 12];
d3.json("world-110m.json", function(error, worldJSON) {
if (error) throw error;
world = topojson.feature(worldJSON, worldJSON.objects.land);
graticule = d3.geoGraticule();
svg1 = createMap("#map1", valuesVar1, range);
svg.push(svg1);
svg2 = createMap("#map2", valuesVar2, range);
svg.push(svg2);
svgLegend = d3.select("#legend").append('svg')
.attr('width', 60)
.attr('height', 800);
svgLegend.append("g").attr("class", "legendLinear");
var legendLinear = d3.legendColor()
.shapeWidth(15)
.shapeHeight(15)
.shapePadding(1)
.cells(nbLevels)
.orient('vertical')
.ascending(true)
.labelAlign('start')
.scale(color);
svgLegend.select(".legendLinear")
.call(legendLinear);
});
//====================================================================================
}
reader_url.readAsArrayBuffer(blob);
}
oReq.send(); //start process
</script>
風場數據可視化結果圖