I would like to represent in 3D the result of the function f(X,Y,Z) = (X²(Y+Z) + Y²(X-Z) +Z²(-X-Y) -21) with scatter points. Any advice how to do this in GNUPLOT is welcome (in -100 to +100 in each direction).
I am trying to see the "0" place at X Y and Z integer. The size of the points would be the value returned by the function.
Integer values [-100:100] in 3 dimensions is 8 million points. If you were to draw anything at all for each of those, the resulting display would be be filled at every pixel and not useful as a visualization tool. Coding by size or color at each point will not help. Now the number of points you are interested in finding may constitute a manageably small subset that can be visualized as 3D scatter plot in gnuplot, but you will have to select those points rather than plotting everything.
The simple approach would be to loop over x, y, z to calculate your function, then if the value is what you want (f(x,y,z)==0?) write out that [x,y,z] triplet to a file. If the number of points found is reasonable then you can visualize in gnuplot with the commands below. I use random points between [-1:1] as an example.
unset border
set xzeroaxis; set yzeroaxis; set zzeroaxis
set tics nomirror axis
set xyplane at 0
splot "zeros" using 1:2:3 with points pt 7 ps 0.4
(pointtype 7 is a solid circle, pointsize 0.4 makes them smaller)
Answer #2: Unless the density function is zero in most places, it fills a solid volume. In order to visualize this volume you will need to select specific points or take slices or something else. Reasonable options depend on the structure of your function. Gnuplot does now offer several approaches. See this example of using 2D slices taken from the current demo set voxel.dem. This approach makes sense for a smooth density function, but probably not for a function where the interesting bits are discrete points rather than regions of space.
thanks. I will follow later.
i just wrote a script which fullfill a discrete representation
so far enough for me with matplotlib (I was no so familiar with gnuplot).
Also: representing a function "F(X,Y,Z) - A = Result" at (X,Y,Z) all with integers in 3D. Result = 0? point colored black. Else colored.
import matplotlib as mpl
from mpl_toolkits.mplot3d.proj3d import proj_transform
import matplotlib.pyplot as plt
from matplotlib.widgets import Button
import numpy as np
mpl.use('tkagg')
def distance(point, event):
plt.sca(ax) # <------------------ introduce this one !!!!!!!!!!!!!!!!!!!!!!!!!!!
x2, y2, _ = proj_transform(point[0], point[1], point[2], plt.gca().get_proj())
x3, y3 = ax.transData.transform((x2, y2))
return np.sqrt ((x3 - event.x)**2 + (y3 - event.y)**2)
def calcClosestDatapoint(X, event):
distances = [distance(X[i, 0:3], event) for i in range(Sol)]
return np.argmin(distances)
#
def annotatePlot(X, index):
global last_mark, generated_labels
if activated_labelling:
x2, y2, _ = proj_transform(X[index, 0], X[index, 1], X[index, 2], ax.get_proj())
last_mark = plt.annotate(generated_labels[index],
xy = (x2, y2), xytext = (-20, 20), textcoords = 'offset points', ha = 'right', va = 'bottom',
bbox = dict(boxstyle = 'round,pad=0.5', fc = 'yellow', alpha = 0.5),
arrowprops = dict(arrowstyle = '->', connectionstyle = 'arc3,rad=0'))
fig.canvas.draw()
#
def onMouseMotion(event):
global Coord
if activated_labelling:
closestIndex = calcClosestDatapoint(Coord, event)
last_mark.remove()
annotatePlot(Coord, closestIndex)
def show_on(event):
global activated_labelling, last_mark,pid,mid
if activated_labelling == False:
activated_labelling = True
x2, y2, _ = proj_transform(Coord[0,0], Coord[0,1], Coord[0,2], ax.get_proj())
last_mark = plt.annotate("3D measurement on " + generated_labels[0],
xy = (x2, y2), xytext = (-20, 20), textcoords = 'offset points', ha = 'right', va = 'bottom',
bbox = dict(boxstyle = 'round,pad=0.5', fc = 'yellow', alpha = 0.5),
arrowprops = dict(arrowstyle = '->', connectionstyle = 'arc3,rad=0'))
mid = fig.canvas.mpl_connect('motion_notify_event', onMouseMotion)
#
def show_off(event):
global activated_labelling
'''
deactivate the persistent XYZ position labels at the grafic
'''
if activated_labelling:
activated_labelling = False
last_mark.remove()
fig.canvas.draw()
fig.canvas.mpl_disconnect(mid)
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
#ax = fig.gca(projection='3d')
activated_labelling = False
Wide = 100
Minimum = -50
ScanLimit = 3 # searching between o and 3; 4 and 5 are no solutions
Search = 45
Coord=[]
values=[]
generated_labels = []
#
XMin = 0
XMax = 0
YMin = 0
YMax = 0
ZMin = 0
ZMax = 0
# count the solutions found in the scan area defined above
Sol=0
for i in range(Wide+1):
for j in range(Wide+1):
for k in range(Wide+1):
########################################################################
########################################################################
####
#### THIS IS THE POLYNOM TO BE REPRESENTED
####
param_dens = ((i+Minimum)**3)+((j+Minimum)**3)+((k+Minimum)**3) -Search
if abs(param_dens) <= abs(ScanLimit):
Coord.append([i+Minimum,j+Minimum,k+Minimum])
if ScanLimit !=0:
values.append([abs(param_dens)])
labelling = "value {}\nin X:{} Y:{} Z:{}".format(Search+param_dens,i+Minimum,j+Minimum,k+Minimum)
generated_labels.append(labelling)
print(labelling+"\n")
# increase the number indicating the solutions found
Sol +=1
# for centering the window
if XMin > i+Minimum:
XMin = i+Minimum
if YMin > j+Minimum:
YMin = j+Minimum
if ZMin > k+Minimum:
ZMin = k+Minimum
if XMax < i+Minimum:
XMax = i+Minimum
if YMax < j+Minimum:
YMax = j+Minimum
if ZMax < k+Minimum:
ZMax = k+Minimum
print('######################################################')
print('## statistics / move this to a parallel search engine?')
print('## search ')
print("## total solution %d for searching center %d" % (Sol,Search))
print("## from %d to %d" % (Search-ScanLimit,Search+ScanLimit))
print("## from %d to %d" % (Minimum,Wide+Minimum))
print('##')
print('#######################################################')
#
values = np.array(values, dtype='int64')
Coord = np.array(Coord, dtype='int64')
#
if ScanLimit !=0:
cmap = plt.cm.jet # define the colormap
# extract all colors from the .jet map
cmaplist = [cmap(i) for i in range(cmap.N)]
# force the first color entry to be black
cmaplist[0] = (0, 0, 0, 1.0)
# create the new map
cmap = mpl.colors.LinearSegmentedColormap.from_list('Custom cmap', cmaplist, cmap.N)
# define the bins and normalize
bounds = np.linspace(0, ScanLimit, ScanLimit+1)
norm = mpl.colors.BoundaryNorm(bounds, cmap.N)
# create a second axes for the colorbar
ax2 = fig.add_axes([0.95, 0.1, 0.03, 0.8])
cb = mpl.colorbar.ColorbarBase(ax2, cmap=cmap, norm=norm,
spacing='proportional', ticks=bounds, boundaries=bounds, format='%1i')
#
ax.set_xlim3d(XMin-5, XMax+5)
ax.set_ylim3d(YMin-5, YMax+5)
ax.set_zlim3d(ZMin-5, ZMax+5)
#
ax.set_xlabel('X X')
ax.set_ylabel('Y Y')
ax.set_zlabel('Z Z')
ax.set_aspect(aspect=1)
# extract the scatterplot drawing in a separate function so we ca re-use the code
def draw_scatterplot():
if ScanLimit !=0:
ax.scatter3D(Coord[:,0], Coord[:,1], Coord[:,2], s=20, c=values[:,0], cmap=cmap, norm=norm)
else:
ax.scatter3D(Coord[:,0], Coord[:,1], Coord[:,2], s=20, c='green')
# draw the initial scatterplot
draw_scatterplot()
# create the "on" button, and place it somewhere on the screen
ax_on = plt.axes([0.0, 0.0, 0.1, 0.05])
button_on = Button(ax_on, 'on')
#
ax_off = plt.axes([0.12, 0.0, 0.1, 0.05])
button_off = Button(ax_off, 'off')
#
#ax_off = plt.axes([0.24, 0.0, 0.1, 0.05])
#button_off = Button(ax_off, 'off')
# link the event handler function to the click event on the button
button_on.on_clicked(show_on)
button_off.on_clicked(show_off)
#fig.colorbar(img)
plt.show()