描述
Marching squares是计算机图形学中一个用于矩阵网格点数据生成等值面的一个算法。
其主要思想是将每个2*2网格，其中4个顶点网格点根据等值线阈值处理为0，1标量。
1 代表数值大于阈值 0 代表数值小于阈值
这样处理之后有16中情况，2*2网格点根据左上角顺时针编号。
对于对角线这种情形，依据对角线大于两边处理。
算法主要流程：
1. 遍历每个小网格，从16中固定情形中选择划线类别
2.利用线性插值，结合网格点数值找寻交点
3. 划线
代码
#coding:utf-8
from pylab import *
from matplotlib.ticker import MultipleLocator, FormatStrFormatter
class PlotDemo:
def __init__(self,m,n):
self._net = None
self.m = m
self.n = n
def _set_default_figure(self,m,n):
ax = subplot(111) #注意:一般都在ax中设置,不再plot中设置
#设置主刻度标签的位置,标签文本的格式
xmajorLocator = MultipleLocator(10) #将x主刻度标签设置为20的倍数
xmajorFormatter = FormatStrFormatter('%1.1f') #设置x轴标签文本的格式
xminorLocator = MultipleLocator(1) #将x轴次刻度标签设置为5的倍数
ymajorLocator = MultipleLocator(10) #将y轴主刻度标签设置为0.5的倍数
ymajorFormatter = FormatStrFormatter('%1.1f') #设置y轴标签文本的格式
yminorLocator = MultipleLocator(1) #将x轴次刻度标签设置为5的倍数
ax.set_xlim(0,m)
ax.set_ylim(0,n)
ax.xaxis.set_major_locator(xmajorLocator)
ax.xaxis.set_major_formatter(xmajorFormatter)
ax.yaxis.set_major_locator(ymajorLocator)
ax.yaxis.set_major_formatter(ymajorFormatter)
#显示次刻度标签的位置,没有标签文本
ax.xaxis.set_minor_locator(xminorLocator)
ax.yaxis.set_minor_locator(yminorLocator)
ax.xaxis.grid(True, which='minor') #x坐标轴的网格使用主刻度
ax.yaxis.grid(True, which='minor') #y坐标轴的网格使用次刻度
def show_source(self):
self._set_default_figure(self.m, self.n)
if self._net is None:
raise Exception
net = self._net.net_info
shape = net.shape
xlen = shape[0]
ylen = shape[1]
for i in range(xlen):
for j in range(ylen):
if net[i][j] > 0:
plot(i,j,'g+')
show()
def show_contour(self):
self._set_default_figure(self.m, self.n)
net = self._net.net_info
shape = net.shape
xlen = shape[0]
ylen = shape[1]
for i in range(xlen):
for j in range(ylen):
if net[i][j] > 0:
plot(i,j,'g+')
net = self._net
utils = MarchSquareUtlis(net)
lines = utils.trancing_contours()
width,height = net.net_info.shape
arr = net.net_info
idx = 0
for i in range(width-1):
for j in range(height-1):
x,y = i,j
count,v1,v2,v3,v4,v5,v6,v7,v8 = lines[idx]
idx = idx + 1
if count == 0:
continue
if count == 1:
x1 = x + v1
y1 = y + v2
x2 = x + v3
y2 = y + v4
plot(x1,y1,x2,y2,'ro')
if count == 2:
x1 = x + v1
y1 = y + v2
x2 = x + v3
y2 = y + v4
plot(x1,y1,x2,y2,'ro')
x1 = x + v5
y1 = y + v6
x2 = x + v7
y2 = y + v8
plot(x1,y1,x2,y2,'ro')
show()
def set_net_info(self, net_info):
self._net = net_info
class RandomGenNet(object):
def __init__(self):
self.arr = self._gen_random()
def __init__(self,m,n):
self.arr = self._gen_random(m,n)
@property
def net_info(self):
return self.arr
def _gen_random(self,m=100,n=100):
return np.zeros((m,n),dtype='double')
def add_circle(self,center_x, center_y, radius,val):
r = int(radius + 0.5)
for x in range(center_x - r, center_x + r):
rx = np.abs(center_x - x)
ry = int( np.sqrt(radius * radius - rx * rx))
for y in range(ry):
self.arr[x][center_y - y] = val
self.arr[x][center_y + y] = val
return
def add_retangle(self,lf_up_x, lf_up_y, width, height,val):
for i in range(width):
for j in range(height):
self.arr[lf_up_x + i][lf_up_y + j] = val
#clock-wise, top-left|top-right|bottom-right|bottom_left
def get_retangle_bit(v1, v2, v3, v4):
return v1 << 3 | v2 << 2 | v3 << 1 | v4
#shift relative to top-left
def get_retangle_shift(bitval):
if bitval == 0 or bitval == 15:
return (0,None,None,None,None,None,None,None,None)
if bitval == 1 or 14:
return (1,0,0.5,0.5,1,None,None,None,None)
if bitval == 2 or bitval == 13:
return (1,0.5,1,1,0.5,None,None,None,None)
if bitval == 3 or bitval == 12:
return (1,0,0.5,1,0.5,None,None,None,None)
if bitval == 4 or bitval == 11:
return (1,0,0.5,1,0.5,None,None,None,None)
if bitval == 5:
return (2,0,0.5,0.5,0,0.5,1,1,0.5)
if bitval == 6 or bitval == 9:
return (1,0.5,0,0.5,1,None,None,None,None)
if bitval == 7 or bit == 8:
return (1,0,0.5,0.5,0,None,None,None,None)
if bitval == 10:
return (2,0,0.5,0.5,1,0.5,0,1,0.5)
class MarchSquareUtlis(object):
def __init__(self, net):
self.net = net
def trancing_contours(self):
ret = []
width,height = self.net.net_info.shape
arr = self.net.net_info
for i in range(width-1):
for j in range(height-1):
v1 = int(arr[i][j])
v2 = int(arr[i + 1][j])
v3 = int(arr[i + 1][j + 1])
v4 = int(arr[i][j + 1])
bitv = get_retangle_bit(v1,v2,v3,v4)
net_shift = get_retangle_shift(bitv)
ret.append(net_shift)
return ret
demo = PlotDemo(100,100)
netinfo = RandomGenNet(100,100)
netinfo.add_circle(20,20,10,1)
netinfo.add_circle(70,50,20,1)
netinfo.add_retangle(20,45,40,20,1)
netinfo.add_retangle(70,50,10,10,1)
demo.set_net_info(netinfo)
#demo.show_source()
demo.show_contour()
print get_retangle_bit(1,1,1,0)
print range(3)
参考
Marching_squares_wiki
