绘制等高线图的算法
假设我们在一个矩形范围内，有一些离散的高度数据，形如(x,y,height)这样的数据集，然后我们要得到一张这样的等高线图：我们描述一下生成等高线图的算法。一、图形概述实际地图上实用的单位是米，但我们显示时使用的单位是像素，这里面有一个转换的关系。而且，显示的图形有可能需要缩放。所以我们收集的数据，x和y统一使用相对值。我们把总长和总宽都定为1，实际的坐标相对于单位1来定。例如地图...
·
假设我们在一个矩形范围内，有一些离散的高度数据，形如(x,y,height)这样的数据集，然后我们要得到一张这样的等高线图：
我们描述一下生成等高线图的算法。
一、图形概述
实际地图上实用的单位是米，但我们显示时使用的单位是像素，这里面有一个转换的关系。而且，显示的图形有可能需要缩放。所以我们收集的数据，x和y统一使用相对值。我们把总长和总宽都定为1，实际的坐标相对于单位1来定。例如地图长8km，然后我们的坐标点x为1km，那么x就是1/8=0.125。
实际上，等高线图是分成一个一个网格的，我所使用的网格大小是10*10。每个网格里的线和颜色组成了整张的等高线图。
二、插值
我们的数据都是形如(x,y,height)这样的格式，如(0.85f, 0.05f, 0.49f)、(0.75f, 0.1f, 0.7f)、(0.5f, 0.65f, 0.04f)。
假如我们要显示的图形大小是800*600，切成网格后，就会有80*60=4800个网格。一般情况下，我们收集到的数据量远远小于这个数，所以我们就需要插值。
我们用比较简单的距离反比策略进行插值。我们知道，跟某个点的距离越近，受它的影响就越大。所以，跟点V(x1,y1)相距d的点P(x,y)，它受V的影响可按以下公式计算：
P-V=V/d
其中d=(x-x1)^2+(y-y1)^2
而P的最终值，由平面上所有已知值的点确定。
代码如下：
double D = 0;
double DV = 0;
foreach (IntMeasureData imd in measure_data)
{
double d = 1.0 / ((imd.X - i) * (imd.X - i) + (imd.Y - j) * (imd.Y - j));
D += d;
DV += imd.Z * d;
}
value = (float)(DV / D);三、Marching squares算法
这就是生成等高线图的一种算法，我们来描述一下。
我们上面提到过，把800*600的图切成了80*60个网格，而在第二步里面，我们把这些网格都填满数据了。我们假设数据如下：
图中的数据就是高度。我们现在要画一条高度值为5的线。我们把5跟表格中的数据比较，比5大的为1，用实心黑点表示，比5小的为0，化为一个二值的表格，如下图所示：
然后如果我们把每个格子取出来，就会发现有以下的16种情况：
只要我们把上图的线画出来，就会拼出一张等高线图。
四、算法细节
问题1，我们可能需要的是有颜色的图，而不是只有线条的图。
我们只要把上图中格子里的线和实心黑点围成的区域填充颜色即可。
问题2，不同的高度值，如何绘图。
先画低阈值的，再画高阈值的，一层层往上叠即可。
问题3，有没有办法把所有格子的线条围成一个区域，然后对这个区域进行一次性颜色的填充。
如果有，就不会用这种分成一格一格的方法去做了。只能一格一格地填充颜色。
五、代码解析
代码是C#、WPF。
我们先定义测量数据的数据结构：
/// <summary>
/// 测量数据
/// </summary>
public struct MeasureData
{
/// <summary>
/// 初始化测量数据
/// </summary>
/// <param name="x">坐标x</param>
/// <param name="y">坐标y</param>
/// <param name="z">高度</param>
public MeasureData(float x, float y, float z)
{
X = x;
Y = y;
Z = z;
}
/// <summary>
/// 坐标X
/// </summary>
public float X;
/// <summary>
/// 坐标Y
/// </summary>
public float Y;
/// <summary>
/// 高度
/// </summary>
public float Z;
}
internal struct IntMeasureData
{
public IntMeasureData(MeasureData md, int x_num, int y_num)
{
X = (int)(md.X * x_num);
if (X >= x_num)
{
X = x_num - 1;
}
Y = (int)(md.Y * y_num);
if (Y >= y_num)
{
Y = y_num - 1;
}
Z = md.Z;
}
public int X;
public int Y;
public float Z;
}IntMeasureData是存放在网格里的数据。
然后我们做第一步：数据插值：
/// <summary>
/// 数据插值
/// </summary>
private void InitData()
{
x_num = (int)(w / grid_w);
y_num = (int)(h / grid_h);
data = new float[x_num, y_num];
IntMeasureData[] measure_data = new IntMeasureData[HeightDots.Length];
for (int i = HeightDots.Length - 1; i >= 0; i--)
{
measure_data[i] = new IntMeasureData(HeightDots[i], x_num, y_num);
}
min = float.MaxValue;
max = float.MinValue;
for (int i = 0; i < x_num; i++)
{
for (int j = 0; j < y_num; j++)
{
float value = 0;
bool find = false;
foreach (IntMeasureData imd in measure_data)
{
if (i == imd.X && j == imd.Y)
{
value = imd.Z;
find = true;
break;
}
}
if (!find)
{
double D = 0;
double DV = 0;
foreach (IntMeasureData imd in measure_data)
{
double d = 1.0 / ((imd.X - i) * (imd.X - i) + (imd.Y - j) * (imd.Y - j));
D += d;
DV += imd.Z * d;
}
value = (float)(DV / D);
}
data[i, j] = value;
min = Math.Min(min, value);
max = Math.Max(max, value);
}
}
}/// <summary>
/// 获取某个阈值下的图形数据
/// </summary>
/// <param name="threshold">阈值</param>
/// <returns>图形数据</returns>
private List<List<PointF>> CreateMapData(float threshold)
{
byte[,] binary_data = new byte[x_num, y_num];
for (int i = 0; i < x_num; i++)
{
for (int j = 0; j < y_num; j++)
{
binary_data[i, j] = (byte)((data[i, j] >= threshold) ? 1 : 0);
}
}
List<List<PointF>> shapes = new List<List<PointF>>();
for (int i = 1; i < x_num; i++)
{
for (int j = 1; j < y_num; j++)
{
int num = (binary_data[i - 1, j - 1] << 3) + (binary_data[i, j - 1] << 2) + (binary_data[i, j] << 1) + binary_data[i - 1, j];
List<PointF> points = new List<PointF>();
switch (num)
{
case 0:
break;
case 1:
AddLeft(points, i, j, threshold);
AddLeftBottom(points, i, j);
AddBottom(points, i, j, threshold);
break;
case 2:
AddRight(points, i, j, threshold);
AddRightBottom(points, i, j);
AddBottom(points, i, j, threshold);
break;
case 3:
AddLeft(points, i, j, threshold);
AddLeftBottom(points, i, j);
AddRightBottom(points, i, j);
AddRight(points, i, j, threshold);
break;
case 4:
AddTop(points, i, j, threshold);
AddRightTop(points, i, j);
AddRight(points, i, j, threshold);
break;
case 5:
AddLeft(points, i, j, threshold);
AddRightTop(points, i, j);
AddRight(points, i, j, threshold);
AddBottom(points, i, j, threshold);
AddLeftBottom(points, i, j);
AddLeft(points, i, j, threshold);
break;
case 6:
AddTop(points, i, j, threshold);
AddRightTop(points, i, j);
AddRightBottom(points, i, j);
AddBottom(points, i, j, threshold);
break;
case 7:
AddTop(points, i, j, threshold);
AddRightTop(points, i, j);
AddRightBottom(points, i, j);
AddLeftBottom(points, i, j);
AddLeft(points, i, j, threshold);
break;
case 8:
AddTop(points, i, j, threshold);
AddLeftTop(points, i, j);
AddLeft(points, i, j, threshold);
break;
case 9:
AddTop(points, i, j, threshold);
AddLeftTop(points, i, j);
AddLeftBottom(points, i, j);
AddBottom(points, i, j, threshold);
break;
case 10:
AddTop(points, i, j, threshold);
AddLeftTop(points, i, j);
AddLeft(points, i, j, threshold);
AddBottom(points, i, j, threshold);
AddRightBottom(points, i, j);
AddRight(points, i, j, threshold);
break;
case 11:
AddTop(points, i, j, threshold);
AddLeftTop(points, i, j);
AddLeftBottom(points, i, j);
AddRightBottom(points, i, j);
AddRight(points, i, j, threshold);
break;
case 12:
AddLeft(points, i, j, threshold);
AddRight(points, i, j, threshold);
AddRightTop(points, i, j);
AddLeftTop(points, i, j);
break;
case 13:
AddRight(points, i, j, threshold);
AddRightTop(points, i, j);
AddLeftTop(points, i, j);
AddLeftBottom(points, i, j);
AddBottom(points, i, j, threshold);
break;
case 14:
AddLeft(points, i, j, threshold);
AddBottom(points, i, j, threshold);
AddRightBottom(points, i, j);
AddRightTop(points, i, j);
AddLeftTop(points, i, j);
break;
case 15:
AddLeftTop(points, i, j);
AddRightTop(points, i, j);
AddRightBottom(points, i, j);
AddLeftBottom(points, i, j);
break;
}
if (num != 0)
{
shapes.Add(points);
}
}
}
return shapes;
}private void AddLeft(List<PointF> list, int x, int y, float threshold)
{
list.Add(new PointF((x - 1) * grid_w, (y - 1 + V(data[x - 1, y - 1], data[x - 1, y], threshold)) * grid_h));
}
private void AddRight(List<PointF> list, int x, int y, float threshold)
{
list.Add(new PointF(x * grid_w, (y - 1 + V(data[x, y - 1], data[x, y], threshold)) * grid_h));
}
private void AddTop(List<PointF> list, int x, int y, float threshold)
{
list.Add(new PointF((x - 1 + V(data[x - 1, y - 1], data[x, y - 1], threshold)) * grid_w, (y - 1) * grid_h));
}
private void AddBottom(List<PointF> list, int x, int y, float threshold)
{
list.Add(new PointF((x - 1 + V(data[x - 1, y], data[x, y], threshold)) * grid_w, y * grid_h));
}
private void AddLeftTop(List<PointF> list, int x, int y)
{
list.Add(new PointF((x - 1) * grid_w, (y - 1) * grid_h));
}
private void AddRightTop(List<PointF> list, int x, int y)
{
list.Add(new PointF(x * grid_w, (y - 1) * grid_h));
}
private void AddLeftBottom(List<PointF> list, int x, int y)
{
list.Add(new PointF((x - 1) * grid_w, y * grid_h));
}
private void AddRightBottom(List<PointF> list, int x, int y)
{
list.Add(new PointF(x * grid_w, y * grid_h));
}
private float V(float min, float max, float threshold)
{
return (threshold - min) / (max - min);
}结果是我们得到了每个网格要填充颜色的端点列表。把这些端点围成一个多边形，再绘制成图形即可。
开放原子开发者工作坊旨在鼓励更多人参与开源活动，与志同道合的开发者们相互交流开发经验、分享开发心得、获取前沿技术趋势。工作坊有多种形式的开发者活动，如meetup、训练营等，主打技术交流，干货满满，真诚地邀请各位开发者共同参与！
更多推荐
10
0- 0
已为社区贡献15条内容
所有评论(0)