Gnuplot link 边缘上的所有点并填充区域
Gnuplot link all points on edges and fill the area
我有这样的 x-y 数据:
133.322 1073
399.966 1073
799.932 1073
1333.22 1073
133.322 1078
399.966 1078
799.932 1078
1333.22 1078
133.322 1083
399.966 1083
799.932 1083
1333.22 1083
133.322 1088
399.966 1088
799.932 1088
1333.22 1088
133.322 1093
399.966 1093
799.932 1093
1333.22 1093
133.322 1098
399.966 1098
799.932 1098
1333.22 1098
133.322 1103
399.966 1103
799.932 1103
1333.22 1103
1333.22 1108
我得到的是这样的:
或
我怎样才能用 Gnuplot 做到这一点?我需要以这种方式绘制许多其他数据。所以一个通用的解决方案会很棒。非常感谢。
以下是尝试围绕随机点绘制轮廓。
数学家可能有更好的程序的想法。
程序:
- 确定数据点的边界框,即左、上、右、下的极值点,通常定义四边形,但也可以是三角形或 2 边形
- 对于每个 "quadrant",即 left/top、top/right、right/bottom、bottom/left,通过确定检查其他点是否在该多边形之外两个极值点和给定点之间的曲率(或方向)。
- 根据象限用升序或降序对外部点进行排序,基本上对轮廓路径进行排序"clockwise"。问题:gnuplot 没有排序功能。您可以使用外部程序,但我通常更喜欢 gnuplot-only 解决方案,这样它会更长一些。
- 删除导致凹轮廓的所有点,这可能需要多次迭代。
我很乐意看到更简单的解决方案。欢迎提出改进建议。
代码:
### find and fill outline of random distributed points
reset session
# create some random test data
set samples 50
set table $Data
plot '+' u (invnorm(rand(0))*100):(invnorm(rand(0))*100) w table
unset table
# if you have a file, skip the above test data section
# and uncomment the following 3 lines and put your filename
#set table $Data
# plot "myData.dat" u 1:2 w table
#unset table
# get the coordinates of the leftmost, topmost, rightmost, and bottommost points
# in case of multiple points with leftmost x-value take the one with the topmost y-value
# similar of the other extreme points
#
# function to initialize extreme point cordinates
Init(x,y) = (Lx=column(x),Ly=column(y), \
Tx=column(x),Ty=column(y), \
Rx=column(x),Ry=column(y), \
Bx=column(x),By=column(y) )
# find extreme point coordinates
GetExtremes(x,y) = (Lx>column(x) ? (Lx=column(x),Ly=column(y)) : \
(Lx==column(x) && Ly<column(y) ? Ly=column(y) : NaN), \
Ty<column(y) ? (Tx=column(x),Ty=column(y)) : \
(Ty==column(y) && Tx<column(x) ? Tx=column(x) : NaN), \
Rx<column(x) ? (Rx=column(x),Ry=column(y)) : \
(Rx==column(x) && Ry>column(y) ? Ry=column(y) : NaN), \
By>column(y) ? (Bx=column(x),By=column(y)) : \
(By==column(y) && Bx>column(x) ? Bx=column(x) : NaN))
set table $Dummy
plot n=0 $Data u (n==0?(Init(1,2),n=n+1):GetExtremes(1,2)) w table
unset table
# put extremal points into array 1=left,2=top,3=right,4=bottom
# and 5=left again for closed line of the extremal polygon
array EPx[5] = [Lx,Tx,Rx,Bx,Lx]
array EPy[5] = [Ly,Ty,Ry,By,Ly]
# (re-)define sgn() function such that sgn(NaN) gives NaN (gnuplot would give 0)
mysgn(x) = x==x ? sgn(x) : NaN
# function to determine orientation of the curves through 3 points A,B,C
# output: -1=clockwise; +1=counterclockwise (mathematical positive); NaN if one point is NaN
Orientation(xa,ya,xb,yb,xc,yc) = mysgn((xb-xa)*(yc-ya)-(xc-xa)*(yb-ya))
# determine outside points for all 4 "quadrants"
# a point is outside
set datafile missing NaN # important for removing points from outline
array SinglePoint[1] # dummy array for plotting a single point
array EmptyLines[2] # dummy array for plotting two emtpy lines
array OutsidePoints[4] # number of outside points for each "quadrant"
set table $OutsidePoints
do for [n=1:4] {
i=0
plot SinglePoint u (EPx[n]):(EPy[n]) w table
plot $Data u 1:(Orientation(EPx[n],EPy[n],,,EPx[n%4+1],EPy[n%4+1]) < 0 ? (i=i+1,):NaN) w table
plot SinglePoint u (EPx[n%4+1]):(EPy[n%4+1]) w table
plot EmptyLines u ("") w table
OutsidePoints[n] = i+2
}
unset table
# sort the outside points of the 4 "quadrants" by in- or decreasing x
# since gnuplot has no built-in sort function it's a bit lengthy
AllOutsidePoints = (sum [i=1:4] OutsidePoints[i])-3 # reduce by 3 double extremal points
array Xall[AllOutsidePoints]
array Yall[AllOutsidePoints]
idx = 0
do for [n=1:4] {
array X[OutsidePoints[n]] # initialize array
array Y[OutsidePoints[n]] # initialize array
set table $Dummy
plot $OutsidePoints u (X[[=10=]+1]=, Y[[=10=]+1]=) index n-1 w table
unset table
# Bubblesort, inefficient but simple
SortDirection = n<3 ? +1 : -1 # n=1,2: +1=ascending, n=3,4: -1=descending
do for [j=OutsidePoints[n]:2:-1] {
do for [i=1:j-1] {
if ( X[i]*SortDirection > X[i+1]*SortDirection) {
tmp=X[i]; X[i]=X[i+1]; X[i+1]=tmp;
tmp=Y[i]; Y[i]=Y[i+1]; Y[i+1]=tmp;
}
}
}
# append array to datablock
set print $Outline append
do for [i=1:|X|-(n==4?0:1)] { print sprintf("%g %g",X[i],Y[i]) }
set print
# create an array for all sorted outside datapoints
last = |X|-(n==4?0:1)
do for [i=1:last] {
Xall[idx+i]=X[i]
Yall[idx+i]=Y[i]
}
idx=idx+last
}
# function checks convexity: result: >0: concave, 1=convex
# Orientation: -1=clockwise, 0=straight, +1=counterclockwise, NaN=undefined point
# function actually doesn't depend on n, just to shorten the function call later
CheckConvexity(n) = (Convex=1,sum [i=2:|Xall|-1] ( idx0=i-1, idx1=i, idx2=i+1, Convex=Convex && \
(Orientation(Xall[idx0],Yall[idx0],Xall[idx1],Yall[idx1],Xall[idx2],Yall[idx2])<0)),Convex)
# put array to datablock
set table $Convex
plot Xall u (Xall[[=10=]+1]):(Yall[[=10=]+1]) w table
unset table
# remove concave points (could take several iterations)
Count=0
while (!CheckConvexity(0)) {
Count = Count+1
print sprintf("Iteration %d ...",Count)
# set concave points to NaN
do for [i=2:|Xall|-1] {
idx0=i-1; idx1=i; idx2=i+1
if (Orientation(Xall[idx0],Yall[idx0],Xall[idx1],Yall[idx1],Xall[idx2],Yall[idx2])>=0) {
Yall[idx1] = NaN
}
}
# remove NaN points by plotting it to datablock
set table $Convex
plot Xall u (Xall[[=10=]+1]):(Yall[[=10=]+1]) w table
unset table
# create new X,Y array with reduced size
array Xall[|$Convex|]
array Yall[|$Convex|]
# put datablock into array again for next iteration if necessary
set table $Dummy
plot $Convex u (Xall[[=10=]+1]=):(Yall[[=10=]+1]=) w table
unset table
}
print sprintf("Convex after %d iterations.",Count)
set object 1 rect from Lx,By to Rx,Ty dt 3
set offsets graph 0.01, graph 0.01, graph 0.01, graph 0.01
set key out top center horizontal
plot \
keyentry w l ls 0 ti "Bounding box", \
$Outline u 1:2 w filledcurves closed lc rgb "green" fs solid 0.1 not noautoscale, \
EPx u (EPx[[=10=]+1]):(EPy[[=10=]+1]) w l lc rgb "black" ti "Extremal polygon", \
$Data u 1:2 w p pt 7 lc rgb "blue" ti "Points inside polygon", \
$OutsidePoints u 1:2 w p pt 7 lc rgb "red" ti "Points outside polygon", \
Xall u (Xall[[=10=]+1]):(Yall[[=10=]+1]) w l lw 1 lc rgb "green" ti "Concave outline", \
$Convex u 1:2 w lp pt 7 lw 2 lc rgb "orange" ti "Convex outline", \
### end of code
结果:
加法:(第二道工序)
我同意@Christoph 的观点,gnuplot 可能不是执行此操作的最佳工具,但您仍然可以这样做:-)。这是一个较短的过程。
程序:
搜索最底部的点,如果有多个点的 y 值最底部,则取 x 坐标最右边的那个点。
从这一点计算与所有其他点的角度 a
。找到与目标角度差值最小的点,开始时为aTar=0
。如果有几个相同的最小差角的点取距离最大的点dMax
。将此点添加到大纲中。
将新的目标角度设置为前一个点和新找到的点之间的线的角度。
再次转到 2,直到您得到第一个点 BxStart,ByStart
。
我希望我的描述在某种程度上是可以理解的。尽管代码比第一个版本短,但它可能效率较低,因为它必须多次(即轮廓点的数量)遍历所有点。最坏的情况是已经在圆上对齐的点。
代码:
### fill outline of random points
reset session
# create some test data
set samples 50
set table $Data
plot '+' u (invnorm(rand(0))*100):(invnorm(rand(0))*100) w table
unset table
### if you have a file, uncomment the following lines and put your filename
#set table $Data
# plot "myData.dat" u 1:2 w table
#unset table
# Angle by dx,dy (range: -90°<= angle < 270°), NaN if dx==dy==0
set angle degrees
Angle(dx,dy) = dx==dx && dy==dy ? dx==0 ? dy==0 ? NaN : sgn(dy)*90 : dx<0 ? 180+atan(dy/dx) : atan(dy/dx) : NaN
# Angle by dx,dy (range: 0°<= angle < 360°), NaN if dx==dy==0
Angle360(dx,dy) = (a_tmp=Angle(dx,dy), a_tmp==a_tmp) ? a_tmp-360*floor(a_tmp/360.) : NaN
# get the coordinates of the bottommost point
# in case of multiple points with bottommost y-coordinate take the one with the rightmost x-coordinate
Init(x,y) = (Bx=column(x),By=column(y))
GetExtremes(x,y) = (By>column(y) ? (Bx=column(x),By=column(y)) : \
(By==column(y) && Bx>column(x) ? Bx=column(x) : NaN))
set table $Dummy
plot $Data u ([=11=]==0?Init(1,2):GetExtremes(1,2)) w table
unset table
print sprintf("Bottommost point: %g,%g:", Bx,By)
# Distance
Dist(x0,y0,x1,y1) = sqrt((x1-x0)**2 + (y1-y0)**2)
# function for getting the next outline point
GetNext(x,y) = (a=Angle360(column(x)-Bx,column(y)-By), aDiff=(a<aTar?a+360-aTar:a-aTar),\
d=Dist(Bx,column(x),By,column(y)), \
aMin>aDiff ? (aNext=a, aMin=aDiff,dMax=d,xNext=column(x),yNext=column(y)) : \
aMin==aDiff ? dMax<d ? (dMax=d, xNext=column(x),yNext=column(y)) : NaN : NaN)
BxStart=Bx; ByStart=By
set print $Outline append
print sprintf("% 9g % 9g",BxStart,ByStart) # add starting point to outline
aTar=0
while 1 { # endless loop
aMin=360; dMax=0
set table $Dummy
plot $Data u (GetNext(1,2)) w table
unset table
print sprintf("% 9g % 9g",xNext,yNext)
Bx=xNext; By=yNext; aTar=aNext
if (xNext==BxStart && yNext==ByStart) { break } # exit loop
}
set print
plot $Data u 1:2 w p pt 7 ti "Datapoints",\
$Outline u 1:2 w l lc rgb "red" ti "Convex Outline"
### end of code
结果:
我有这样的 x-y 数据:
133.322 1073
399.966 1073
799.932 1073
1333.22 1073
133.322 1078
399.966 1078
799.932 1078
1333.22 1078
133.322 1083
399.966 1083
799.932 1083
1333.22 1083
133.322 1088
399.966 1088
799.932 1088
1333.22 1088
133.322 1093
399.966 1093
799.932 1093
1333.22 1093
133.322 1098
399.966 1098
799.932 1098
1333.22 1098
133.322 1103
399.966 1103
799.932 1103
1333.22 1103
1333.22 1108
我得到的是这样的:
或
我怎样才能用 Gnuplot 做到这一点?我需要以这种方式绘制许多其他数据。所以一个通用的解决方案会很棒。非常感谢。
以下是尝试围绕随机点绘制轮廓。 数学家可能有更好的程序的想法。
程序:
- 确定数据点的边界框,即左、上、右、下的极值点,通常定义四边形,但也可以是三角形或 2 边形
- 对于每个 "quadrant",即 left/top、top/right、right/bottom、bottom/left,通过确定检查其他点是否在该多边形之外两个极值点和给定点之间的曲率(或方向)。
- 根据象限用升序或降序对外部点进行排序,基本上对轮廓路径进行排序"clockwise"。问题:gnuplot 没有排序功能。您可以使用外部程序,但我通常更喜欢 gnuplot-only 解决方案,这样它会更长一些。
- 删除导致凹轮廓的所有点,这可能需要多次迭代。
我很乐意看到更简单的解决方案。欢迎提出改进建议。
代码:
### find and fill outline of random distributed points
reset session
# create some random test data
set samples 50
set table $Data
plot '+' u (invnorm(rand(0))*100):(invnorm(rand(0))*100) w table
unset table
# if you have a file, skip the above test data section
# and uncomment the following 3 lines and put your filename
#set table $Data
# plot "myData.dat" u 1:2 w table
#unset table
# get the coordinates of the leftmost, topmost, rightmost, and bottommost points
# in case of multiple points with leftmost x-value take the one with the topmost y-value
# similar of the other extreme points
#
# function to initialize extreme point cordinates
Init(x,y) = (Lx=column(x),Ly=column(y), \
Tx=column(x),Ty=column(y), \
Rx=column(x),Ry=column(y), \
Bx=column(x),By=column(y) )
# find extreme point coordinates
GetExtremes(x,y) = (Lx>column(x) ? (Lx=column(x),Ly=column(y)) : \
(Lx==column(x) && Ly<column(y) ? Ly=column(y) : NaN), \
Ty<column(y) ? (Tx=column(x),Ty=column(y)) : \
(Ty==column(y) && Tx<column(x) ? Tx=column(x) : NaN), \
Rx<column(x) ? (Rx=column(x),Ry=column(y)) : \
(Rx==column(x) && Ry>column(y) ? Ry=column(y) : NaN), \
By>column(y) ? (Bx=column(x),By=column(y)) : \
(By==column(y) && Bx>column(x) ? Bx=column(x) : NaN))
set table $Dummy
plot n=0 $Data u (n==0?(Init(1,2),n=n+1):GetExtremes(1,2)) w table
unset table
# put extremal points into array 1=left,2=top,3=right,4=bottom
# and 5=left again for closed line of the extremal polygon
array EPx[5] = [Lx,Tx,Rx,Bx,Lx]
array EPy[5] = [Ly,Ty,Ry,By,Ly]
# (re-)define sgn() function such that sgn(NaN) gives NaN (gnuplot would give 0)
mysgn(x) = x==x ? sgn(x) : NaN
# function to determine orientation of the curves through 3 points A,B,C
# output: -1=clockwise; +1=counterclockwise (mathematical positive); NaN if one point is NaN
Orientation(xa,ya,xb,yb,xc,yc) = mysgn((xb-xa)*(yc-ya)-(xc-xa)*(yb-ya))
# determine outside points for all 4 "quadrants"
# a point is outside
set datafile missing NaN # important for removing points from outline
array SinglePoint[1] # dummy array for plotting a single point
array EmptyLines[2] # dummy array for plotting two emtpy lines
array OutsidePoints[4] # number of outside points for each "quadrant"
set table $OutsidePoints
do for [n=1:4] {
i=0
plot SinglePoint u (EPx[n]):(EPy[n]) w table
plot $Data u 1:(Orientation(EPx[n],EPy[n],,,EPx[n%4+1],EPy[n%4+1]) < 0 ? (i=i+1,):NaN) w table
plot SinglePoint u (EPx[n%4+1]):(EPy[n%4+1]) w table
plot EmptyLines u ("") w table
OutsidePoints[n] = i+2
}
unset table
# sort the outside points of the 4 "quadrants" by in- or decreasing x
# since gnuplot has no built-in sort function it's a bit lengthy
AllOutsidePoints = (sum [i=1:4] OutsidePoints[i])-3 # reduce by 3 double extremal points
array Xall[AllOutsidePoints]
array Yall[AllOutsidePoints]
idx = 0
do for [n=1:4] {
array X[OutsidePoints[n]] # initialize array
array Y[OutsidePoints[n]] # initialize array
set table $Dummy
plot $OutsidePoints u (X[[=10=]+1]=, Y[[=10=]+1]=) index n-1 w table
unset table
# Bubblesort, inefficient but simple
SortDirection = n<3 ? +1 : -1 # n=1,2: +1=ascending, n=3,4: -1=descending
do for [j=OutsidePoints[n]:2:-1] {
do for [i=1:j-1] {
if ( X[i]*SortDirection > X[i+1]*SortDirection) {
tmp=X[i]; X[i]=X[i+1]; X[i+1]=tmp;
tmp=Y[i]; Y[i]=Y[i+1]; Y[i+1]=tmp;
}
}
}
# append array to datablock
set print $Outline append
do for [i=1:|X|-(n==4?0:1)] { print sprintf("%g %g",X[i],Y[i]) }
set print
# create an array for all sorted outside datapoints
last = |X|-(n==4?0:1)
do for [i=1:last] {
Xall[idx+i]=X[i]
Yall[idx+i]=Y[i]
}
idx=idx+last
}
# function checks convexity: result: >0: concave, 1=convex
# Orientation: -1=clockwise, 0=straight, +1=counterclockwise, NaN=undefined point
# function actually doesn't depend on n, just to shorten the function call later
CheckConvexity(n) = (Convex=1,sum [i=2:|Xall|-1] ( idx0=i-1, idx1=i, idx2=i+1, Convex=Convex && \
(Orientation(Xall[idx0],Yall[idx0],Xall[idx1],Yall[idx1],Xall[idx2],Yall[idx2])<0)),Convex)
# put array to datablock
set table $Convex
plot Xall u (Xall[[=10=]+1]):(Yall[[=10=]+1]) w table
unset table
# remove concave points (could take several iterations)
Count=0
while (!CheckConvexity(0)) {
Count = Count+1
print sprintf("Iteration %d ...",Count)
# set concave points to NaN
do for [i=2:|Xall|-1] {
idx0=i-1; idx1=i; idx2=i+1
if (Orientation(Xall[idx0],Yall[idx0],Xall[idx1],Yall[idx1],Xall[idx2],Yall[idx2])>=0) {
Yall[idx1] = NaN
}
}
# remove NaN points by plotting it to datablock
set table $Convex
plot Xall u (Xall[[=10=]+1]):(Yall[[=10=]+1]) w table
unset table
# create new X,Y array with reduced size
array Xall[|$Convex|]
array Yall[|$Convex|]
# put datablock into array again for next iteration if necessary
set table $Dummy
plot $Convex u (Xall[[=10=]+1]=):(Yall[[=10=]+1]=) w table
unset table
}
print sprintf("Convex after %d iterations.",Count)
set object 1 rect from Lx,By to Rx,Ty dt 3
set offsets graph 0.01, graph 0.01, graph 0.01, graph 0.01
set key out top center horizontal
plot \
keyentry w l ls 0 ti "Bounding box", \
$Outline u 1:2 w filledcurves closed lc rgb "green" fs solid 0.1 not noautoscale, \
EPx u (EPx[[=10=]+1]):(EPy[[=10=]+1]) w l lc rgb "black" ti "Extremal polygon", \
$Data u 1:2 w p pt 7 lc rgb "blue" ti "Points inside polygon", \
$OutsidePoints u 1:2 w p pt 7 lc rgb "red" ti "Points outside polygon", \
Xall u (Xall[[=10=]+1]):(Yall[[=10=]+1]) w l lw 1 lc rgb "green" ti "Concave outline", \
$Convex u 1:2 w lp pt 7 lw 2 lc rgb "orange" ti "Convex outline", \
### end of code
结果:
加法:(第二道工序)
我同意@Christoph 的观点,gnuplot 可能不是执行此操作的最佳工具,但您仍然可以这样做:-)。这是一个较短的过程。
程序:
搜索最底部的点,如果有多个点的 y 值最底部,则取 x 坐标最右边的那个点。
从这一点计算与所有其他点的角度
a
。找到与目标角度差值最小的点,开始时为aTar=0
。如果有几个相同的最小差角的点取距离最大的点dMax
。将此点添加到大纲中。将新的目标角度设置为前一个点和新找到的点之间的线的角度。
再次转到 2,直到您得到第一个点
BxStart,ByStart
。
我希望我的描述在某种程度上是可以理解的。尽管代码比第一个版本短,但它可能效率较低,因为它必须多次(即轮廓点的数量)遍历所有点。最坏的情况是已经在圆上对齐的点。
代码:
### fill outline of random points
reset session
# create some test data
set samples 50
set table $Data
plot '+' u (invnorm(rand(0))*100):(invnorm(rand(0))*100) w table
unset table
### if you have a file, uncomment the following lines and put your filename
#set table $Data
# plot "myData.dat" u 1:2 w table
#unset table
# Angle by dx,dy (range: -90°<= angle < 270°), NaN if dx==dy==0
set angle degrees
Angle(dx,dy) = dx==dx && dy==dy ? dx==0 ? dy==0 ? NaN : sgn(dy)*90 : dx<0 ? 180+atan(dy/dx) : atan(dy/dx) : NaN
# Angle by dx,dy (range: 0°<= angle < 360°), NaN if dx==dy==0
Angle360(dx,dy) = (a_tmp=Angle(dx,dy), a_tmp==a_tmp) ? a_tmp-360*floor(a_tmp/360.) : NaN
# get the coordinates of the bottommost point
# in case of multiple points with bottommost y-coordinate take the one with the rightmost x-coordinate
Init(x,y) = (Bx=column(x),By=column(y))
GetExtremes(x,y) = (By>column(y) ? (Bx=column(x),By=column(y)) : \
(By==column(y) && Bx>column(x) ? Bx=column(x) : NaN))
set table $Dummy
plot $Data u ([=11=]==0?Init(1,2):GetExtremes(1,2)) w table
unset table
print sprintf("Bottommost point: %g,%g:", Bx,By)
# Distance
Dist(x0,y0,x1,y1) = sqrt((x1-x0)**2 + (y1-y0)**2)
# function for getting the next outline point
GetNext(x,y) = (a=Angle360(column(x)-Bx,column(y)-By), aDiff=(a<aTar?a+360-aTar:a-aTar),\
d=Dist(Bx,column(x),By,column(y)), \
aMin>aDiff ? (aNext=a, aMin=aDiff,dMax=d,xNext=column(x),yNext=column(y)) : \
aMin==aDiff ? dMax<d ? (dMax=d, xNext=column(x),yNext=column(y)) : NaN : NaN)
BxStart=Bx; ByStart=By
set print $Outline append
print sprintf("% 9g % 9g",BxStart,ByStart) # add starting point to outline
aTar=0
while 1 { # endless loop
aMin=360; dMax=0
set table $Dummy
plot $Data u (GetNext(1,2)) w table
unset table
print sprintf("% 9g % 9g",xNext,yNext)
Bx=xNext; By=yNext; aTar=aNext
if (xNext==BxStart && yNext==ByStart) { break } # exit loop
}
set print
plot $Data u 1:2 w p pt 7 ti "Datapoints",\
$Outline u 1:2 w l lc rgb "red" ti "Convex Outline"
### end of code
结果: