Volumetric Function Plotter - 2.5 script with GUI

Hi all,
I am working on a script to take functions in the form d = f(x,y,z) and plot them in Blender 2.5. I originally wrote it as a pure Python 3.1 script in IDLE, but then built an interface for it in Blender. It still outputs to a function that must be loaded as a raw texture. I’ve included a few wavefunctions for a particle in an infinite height spherical well. It’s in two picies, I wrote it to save to the /Python31/ directory so it wouldn’t get lost, it’s a simple modification to change that, though. Also, if you want it in one piece, copy the function plotter in after the gui and remove the “from mat4d3 import *” line. I may tinker with it some more, but here it is for now:

Function plotting code:

(mat4d3.py)


# 4D Function Plotter
# For use with Blender Volumetrics (2.5 prerelease version)
# Henry Schreiner
# Version 0.3, 9-25-09
# Python 3.1

from math import *  # So math functions work
import os           # Some optional file handleing stuff
VERBOSE=True        # Extra details during operation

class mat4d:
    '''Make a 4D matrix with a given function'''

    TRANS = {'r':lambda x,y,z: (x,y,z),   # no tf
             's':lambda x,y,z:
                         (sqrt(x**2 + y**2 + z**2),        #r
                         atan2(y,x),                       #phi
                         (acos(z/sqrt(x**2 + y**2 + z**2)) #theta
                         if sqrt(x**2 + y**2 + z**2) else 0.0)),
             'c':lambda x,y,z:
                          (sqrt(x**2 + y**2), #s
                          atan2(y,x),         #phi
                          z)}                 #z}

    def ipos(self, a): return (a/(self.res-1)-.5)*2*self.rang
    
    def __init__(self,
                 func: 'Has to have 3 variables'=lambda x,y,z: x+y+z,
                 res=9, rang=5, coord='rectangular'):
        '''This uses x,y,z; or it uses r, phi, and theta if shp is True'''

        self.func=func
        self.rang=rang
        self.coord=coord
        self.res=res
        self.timeon=False
        self.M = [0]*self.res**3
        if VERBOSE: print('Set length of vector to', self.res**3)

    def make(self):
        for n in range(self.res**3):
            (x,y,z)=map(self.ipos,self.num2ind(n))
            (x,y,z)=self.TRANS[self.coord](x,y,z)
            try: self.M[n]=self.func(x,y,z)
            except ZeroDivisionError: self.M[n]=0

    def save(self, filename):
        '''This will save a .raw.gz file for Blender RAW volumetrics'''
        if os.access(filename,os.F_OK): return False # Don't overwrite!
        else:
            fob = open(filename + '.raw.gz','bw')
            fob.write(bytes(map(round,self.M)))
            fob.close()
            return True

    def num2ind(self,a):
        return (int(a/self.res**2),
                int((a%self.res**2)/self.res),
                a%self.res)

    def ind2num(self,i,j,k):
       return i*self.res**2 + j*self.res + k

    def normalize(self, absolute=False, scale =255):
        '''Normalizes to 0-scl, uses absolute values or total'''
        minimum = min(map(abs,self.M) if absolute else self.M)
        rang = max(map(abs,self.M) if absolute else self.M) - minimum
        for n in range(self.res**3):
            self.M[n]=((abs(self.M[n]) if absolute else self.M[n])
                       -minimum)/rang*scale

#################### Functions ####################
# Setting up wavefunction arrays
psi = [lambda r: 0,                                           #n=0
       lambda r: 1/(r*sqrt(2*pi))*sin(pi*r),                  #n=1
       lambda r: 1/(r*sqrt(2*pi))*sin(2*pi*r)]                #n=2

Y = [[lambda p,t: (1/(4*pi))**(1/2)],                         #lm=00
     
     [lambda p,t: (3/(4*pi))**(1/2)*cos(t),                   #lm=10
      lambda p,t: -(3/(8*pi))**(1/2)*sin(t)*e**(p*1j)],       #lm=11
     
     [lambda p,t: (5/(16*pi))**(1/2)*(3*cos(t)**2-1),         #lm=20
      lambda p,t: -(15/(8*pi))**(1/2)*sin(t)*cos(t)*e**(p*1j),#lm=21
      lambda p,t: (15/(32*pi))**(1/2)*sin(t)**2*e**(p*2j)]]   #lm=22

def wavefn(n,l,m):
    return lambda r,p,t: ((Y[l][m](p,t)*
		Y[l][m](p,t).conjugate()).real*
		psi[n](r)**2 if (r!=0 and r<1) else 0)

def magE(r,p,t):
    if r < 1: return 1*cos(t)**2
    else: return 0

def blenderrun(fn,reso,scale,c,finfile,abso=False):
    import time
    t = time.time()

    q=mat4d(fn,reso,scale,c)
    q.make()
    q.normalize(abso)
    q.save('/Python31/'+finfile)
    
    t = time.time()-t
    return t

GUI code:

(MatGUI.py)


# 4D Function Plotter GUI
# For use with Blender Volumetrics (2.5 prerelease version)
# Henry Schreiner
# Version 0.3, 9-25-09
# Python 3.1

from math import *
from mat4d3 import *
	
class ObjectButtonsPanel(bpy.types.Panel):
	__space_type__ = "PROPERTIES"
	__region_type__ = "WINDOW"
	__context__ = "texture"

class OBJECT_PT_hello(ObjectButtonsPanel):
	__label__ = "Volumetric Function Plotter"
	
	def draw(self, context):
		layout = self.layout
		scene=context.scene
		
		CONVCOORD = {'r':'x,y,z','s':'r,p,t','c':'r,p,z'}
		
		StringProperty= bpy.types.Scene.StringProperty
		IntProperty = bpy.types.Scene.IntProperty
		FloatProperty = bpy.types.Scene.FloatProperty
		EnumProperty = bpy.types.Scene.EnumProperty
		BoolProperty = bpy.types.Scene.BoolProperty
		
		StringProperty( attr="file_name",
				name="New Output File",
				description="Plot to filename",
				default="NewOut")
		
		StringProperty( attr="func_plot",
				name="Any Function",
				description="Function for Plotting",
				default='x+y+z')

		StringProperty( attr="func_name",
				name="Function Name",
				description="Function Name for Plotting",
				default="MagE")
				
		IntProperty( attr="res_grid",
				name="Any resolution",
				description="Resolution of Cube",
				min = 2, max = 512, default = 32)
		
		FloatProperty( attr="scale_grid",
				name="Any scale",
				description="scale of Cube",
				min = 0.0, max = 1000.0, default = 1.0)
		
		COORD = [('r','Rectangular','1'),
				 ('c','Cylindrical','2'),
				 ('s','Spherical','3')]
				
		EnumProperty( attr="coord_sys",
				name="Choose Coords",
				description="Coord system",
				items = COORD, default = 'r')
		
		FUNC =  [('norm','Define','1'),
				 ('def','Named','2'),
				 ('wav','Wavefunction','3')]
				
		EnumProperty( attr="func_sys",
				name="Choose Function",
				description="Function System",
				items = FUNC, default = 'norm')
				
		IntProperty( attr="wf_1",
				name="Wavenumber: n",
				description="First wave number",
				min = 1, max = 2, default = 2)
		
		IntProperty( attr="wf_2",
				name="Wavenumber: l",
				description="Second wave number",
				min = 0, max = 2, default = 2)
				
		IntProperty( attr="wf_3",
				name="Wavenumber: m",
				description="Third wave number",
				min = 0, max = 2, default = 2)
		
		FloatProperty( attr="last_run",
				name="It took this long",
				description="A time",
				min = 0.0, max = 1000.0, default = 0.0)
				
		BoolProperty( attr='func_abs',
				name="ABS value",
				description="Use ABS")
		#########################################
		
		row = layout.row()
		row.itemL(text="Enter A Filename:")
		row.itemR(scene, "file_name", text='')
		
		row = layout.row()
		row.itemL(text="Resolution")
		row.itemR(scene, "res_grid", text='')
		
		row = layout.row()
		row.itemL(text="Scale of Grid:")
		row.itemR(scene, "scale_grid", text='')

		row = layout.row()
		row.itemL(text="Function To Use:")
		row.itemR(scene, "func_sys", text='')
		
		if scene.func_sys != 'wav':
			row = layout.row()
			row.itemL(text="Coord System:")
			row.itemR(scene, "coord_sys", text='')
		
		if scene.func_sys=='norm':
			row = layout.row()
			row.itemL(text="Enter A Function d({}):".format(CONVCOORD[scene.coord_sys]))
			row = layout.row()
			row.itemR(scene, "func_plot", text='')
		elif scene.func_sys=='def':
			row = layout.row()
			row.itemL(text="Enter A Function Name:")
			row = layout.row()
			row.itemR(scene, "func_name", text='')
		else:
			row = layout.row()
			row.itemL(text="Wavenumbers:")
			row = layout.row()
#		   row.itemL(text="n:")
			row.itemR(scene, "wf_1", text='')
#		   row.itemL(text="l:")
			row.itemR(scene, "wf_2", text='')
#		   row.itemL(text="m:")
			row.itemR(scene, "wf_3", text='')
			
	
		row = layout.row()
		
		if scene.func_sys != 'wav':
			row.itemL(text="Use Absolute Values:")
			row.itemR(scene, "func_abs", text='')
			
		row.itemO('plot4d', text='Run')

		if scene.last_run>0:
			row = layout.row()
			row.itemL(text="Last run time: "+
					  format(scene.last_run, '.2f')+
					  ' seconds.')


class SCENE_OT_runplotter4d(bpy.types.Operator):
	"""Runs Mat4D!"""
	__idname__ = "plot4d"

	def invoke(self, context, event):
		coordlets = {'r':'x,y,z','s':'r,p,t','c':'r,p,z'}
		scene = context.scene
		
		c = scene.coord_sys
		abso = scene.func_abs
		
		if scene.func_sys == 'norm':
			fn = scene.func_plot
			fn = fn.replace('^','**')
			fn = 'lambda {}: '.format(coordlets[c])+fn
			print(fn)
			fn = eval(fn)
		elif scene.func_sys == 'def':
			fn = eval(scene.func_name)
		else:
			c='s'
			abso=False
			fn = wavefn(scene.wf_1,scene.wf_2,scene.wf_3)
		reso = scene.res_grid
		scale = scene.scale_grid

		
		finfile = scene.file_name

		context.scene.last_run = blenderrun(fn,reso,scale,c,finfile,abso)
		
		return('FINISHED',)
	
bpy.ops.add(SCENE_OT_runplotter4d)
bpy.types.register(OBJECT_PT_hello)

http://art.schreinerbooks.com/art/BlendVolGui.png

I can’t edit it because it’s too long… That’s odd. I submitted it fine. Here’s the last part, with a corrected picture:

http://schreinerbooks.com/art/BlendVolGui.png

Video examples:

Youtube vid 1
Youtube vid 2

Oh, and it will show up in the texture tab.

Here’s my blend file, it doesn’t have the shadowing that was in my vids, though.
EDIT: Improved starting, just push Alt+P on the left hand window, then you are ready to go! Don’t think you even need to save the .py file(s). With the current color ramp, you are looking at an isosurface. Turn off or adjust to taste. Push ‘Run’ to make the file (is stored in /Python31/) and then render.
here

Woooh… nice !

Man, I have lots of catchup work to do : 2.5, python 3.1.
Is coding a script for 2.5 with py3.1 very different than a script for 2.4x with py2.6 ?

Interesting idea !

very cool :slight_smile:

Very nice visualisation of a wave functions, and seeing this parameter
func: ‘Has to have 3 variables’=lambda x,y,z: x+y+z,

finally got me reading some of the python 3.0 documentation. Annotations!

Thank you. Well, 3.1 is similar to 2.6, just has some nice improvements under the hood (Unicode strings, format member function on strings, etc), one or two ‘major’ changes (print is a function! Help! Actually, it may be nicer as a function…). I think the byte code output was easier in 3.1, but to be honest, I’m not sure I ever had a binary output on a 2.6 script.

Blender 2.5 is really nice to work with for GUI programming, took me one night when I couldn’t sleep to do that part, and I’d never done it before in 2.5. I may have made a simple GUI for 2.4x before, for an atom coordinate reader for a solid state research internship, but that’s all I’ve done in 2.4x GUI’s, so I can’t vouch for difficulty there. It looks way better in 2.5, though, all nicely integrated.

Thanks! :smiley:

Yep! :wink:

I’d like to see a way to see parameters for functions in blender, and if they all had annotations, there is a way to see them using .annotations. If IDLE-style tooltips were there, you wouldn’t even need that… Maybe there is a way to access that easily as it is, I couldn’t find it though.

I would really like to see the ability to animate colorbands, that would be great! Moving isosurfaces… :slight_smile:

I just tried this in a recent build of 2.5 but received an error message when I tried to run the script:

“TypeError: ‘bpy_ops_submodule’ is not callable”

which is a reference to this line from the script I think:

bpy.ops.add(SCENE_OT_runplotter4d)

I’m sure it’s something simple, but I’m not familiar with python. Can anyone suggest a solution?

Cheers.

I’ll try to take a look at it tonight after work. I’ll bet the blender-python API has changed a little. That line adds an operator.

Thanks mate, I appreciate it.

Umm, a little was an understatement. The auto-width ui, standardization to bl_ names, and adjustments to bpy are taking a little time to fix… I’ll post here when I’ve got it working again.

Okay, here you go. http://www.pasteall.org/blend/2973 This is my new, current-version-compatable script in a blend file. You have to reload the texture file (/python31/NewOut.tar.gz) every time you rerun the script. Not sure why… Anyway, let me know if you have problems, and I’ll possibly work on a version as an add-on. The new ui is nicer! It took me 1/2 the space to rewrite the important part of the draw function. :slight_smile:

Thanks alot for doing this, I really appreciate it. I’m playing with it now and it all seems to be working (although I’m still trying the different options etc).

Hopefully this is the last time you’ll need to update the script (the UI at least)!

Loving this script. I’m still playing around with the settings but I’m already getting some promising results:

http://img249.imageshack.us/img249/9274/4dwavefunctionsv0006ble.png

I hope I don’t have to again, but I love the changes they made, make my code nicer, shorter, and more uniform, so if that continues, I’ll be happy to update again. :slight_smile:

That looks excellent, love the use of color bands. Thanks for letting me see what you are doing with it, makes the work it took worthwhile! What are you graphing, if I may ask?

Edit: Oh, the center part is a center part, not something behind it… Is that nlm = 210 then?

Yes that’s just nlm=210. I’m looking to use the script to plot some ‘artists impression’ images of light coupling to nano-particles. So I need to generate a light field, a suggestive (but not necessarily accurate) wavefunction around the nano-particles, and perhaps show some coupling between the two. I’ll post any noteworthy updates.

Thanks again for the script and the update, it’s really nice. :cool:

do you have one updated py. file where we can work with it ?

and can this run from the tet edtior or may be as an addons?

keep up the good work

i like the 3D shape here very strange
any better short intro on these nano things ?

happy 2.5


# 4D Function Plotter
# For use with Blender Volumetrics (2.5 prerelease version)
# Henry Schreiner
# Version 0.3, 9-25-09
# Python 3.1

from math import *  # So math functions work
import os           # Some optional file handleing stuff
VERBOSE=True        # Extra details during operation

class mat4d:
    '''Make a 4D matrix with a given function'''

    TRANS = {'r':lambda x,y,z: (x,y,z),   # no tf
             's':lambda x,y,z:
                         (sqrt(x**2 + y**2 + z**2),        #r
                         atan2(y,x),                       #phi
                         (acos(z/sqrt(x**2 + y**2 + z**2)) #theta
                         if sqrt(x**2 + y**2 + z**2) else 0.0)),
             'c':lambda x,y,z:
                          (sqrt(x**2 + y**2), #s
                          atan2(y,x),         #phi
                          z)}                 #z}

    def ipos(self, a): return (a/(self.res-1)-.5)*2*self.rang
    
    def __init__(self,
                 func: 'Has to have 3 variables'=lambda x,y,z: x+y+z,
                 res=9, rang=5, coord='rectangular'):
        '''This uses x,y,z; or it uses r, phi, and theta if shp is True'''

        self.func=func
        self.rang=rang
        self.coord=coord
        self.res=res
        self.timeon=False
        self.M = [0]*self.res**3
        if VERBOSE: print('Set length of vector to', self.res**3)

    def make(self):
        for n in range(self.res**3):
            (x,y,z)=map(self.ipos,self.num2ind(n))
            (x,y,z)=self.TRANS[self.coord](x,y,z)
            try: self.M[n]=self.func(x,y,z)
            except ZeroDivisionError: self.M[n]=0

    def save(self, filename):
        '''This will save a .raw.gz file for Blender RAW volumetrics'''
        if os.access(filename,os.F_OK): return False # Don't overwrite!
        else:
            fob = open(filename + '.raw.gz','bw')
            fob.write(bytes(map(round,self.M)))
            fob.close()
            return True

    def num2ind(self,a):
        return (int(a/self.res**2),
                int((a%self.res**2)/self.res),
                a%self.res)

    def ind2num(self,i,j,k):
       return i*self.res**2 + j*self.res + k

    def normalize(self, absolute=False, scale =255):
        '''Normalizes to 0-scl, uses absolute values or total'''
        minimum = min(map(abs,self.M) if absolute else self.M)
        rang = max(map(abs,self.M) if absolute else self.M) - minimum
        for n in range(self.res**3):
            self.M[n]=((abs(self.M[n]) if absolute else self.M[n])
                       -minimum)/rang*scale

#################### Functions ####################
# Setting up wavefunction arrays
psi = [lambda r: 0,                                           #n=0
       lambda r: 1/(r*sqrt(2*pi))*sin(pi*r),                  #n=1
       lambda r: 1/(r*sqrt(2*pi))*sin(2*pi*r)]                #n=2

Y = [[lambda p,t: (1/(4*pi))**(1/2)],                         #lm=00
     
     [lambda p,t: (3/(4*pi))**(1/2)*cos(t),                   #lm=10
      lambda p,t: -(3/(8*pi))**(1/2)*sin(t)*e**(p*1j)],       #lm=11
     
     [lambda p,t: (5/(16*pi))**(1/2)*(3*cos(t)**2-1),         #lm=20
      lambda p,t: -(15/(8*pi))**(1/2)*sin(t)*cos(t)*e**(p*1j),#lm=21
      lambda p,t: (15/(32*pi))**(1/2)*sin(t)**2*e**(p*2j)]]   #lm=22

def wavefn(n,l,m):
    return lambda r,p,t: ((Y[l][m](p,t)*
        Y[l][m](p,t).conjugate()).real*
        psi[n](r)**2 if (r!=0 and r<1) else 0)

def magE(r,p,t):
    if r < 1: return 1*cos(t)**2
    else: return 0

def blenderrun(fn,reso,scale,c,finfile,abso=False):
    import time
    t = time.time()

    q=mat4d(fn,reso,scale,c)
    q.make()
    q.normalize(abso)
    q.save('/Python31/'+finfile)
    
    t = time.time()-t
    return t

# 4D Function Plotter GUI
# For use with Blender Volumetrics (2.5 prerelease version - Beta 2 approx)
# Henry Schreiner
# Version 0.5, 6-20-10
# Python 3.1

class ObjectButtonsPanel(bpy.types.Panel):
    bl_space_type = "PROPERTIES"
    bl_region_type = "WINDOW"
    bl_context = "texture"

class OBJECT_PT_hello(ObjectButtonsPanel):
    bl_label = "Volumetric Function Plotter"
    bl_default_closed = False
    
    def draw(self, context):
        
        scene=context.scene
        
        CONVCOORD = {'r':'x,y,z','s':'r,p,t','c':'r,p,z'}
        
        StringProperty= bpy.types.Scene.StringProperty
        IntProperty = bpy.types.Scene.IntProperty
        FloatProperty = bpy.types.Scene.FloatProperty
        EnumProperty = bpy.types.Scene.EnumProperty
        BoolProperty = bpy.types.Scene.BoolProperty
        
        StringProperty( attr="file_name",
                name="New Output File",
                description="Plot to filename",
                default="NewOut")
        
        StringProperty( attr="func_plot",
                name="Any Function",
                description="Function for Plotting",
                default='x+y+z')

        StringProperty( attr="func_name",
                name="Function Name",
                description="Function Name for Plotting",
                default="magE")
                
        IntProperty( attr="res_grid",
                name="Any resolution",
                description="Resolution of Cube",
                min = 2, max = 512, default = 32)
        
        FloatProperty( attr="scale_grid",
                name="Any scale",
                description="scale of Cube",
                min = 0.0, max = 1000.0, default = 1.0)
        
        COORD = [('r','Rectangular','1'),
                 ('c','Cylindrical','2'),
                 ('s','Spherical','3')]
                
        EnumProperty( attr="coord_sys",
                name="Choose Coords",
                description="Coord system",
                items = COORD, default = 'r')
        
        FUNC =  [('norm','Define','1'),
                 ('def','Named','2'),
                 ('wav','Wavefunction','3')]
                
        EnumProperty( attr="func_sys",
                name="Choose Function",
                description="Function System",
                items = FUNC, default = 'norm')
                
        IntProperty( attr="wf_1",
                name="Wavenumber: n",
                description="First wave number",
                min = 1, max = 2, default = 2)
        
        IntProperty( attr="wf_2",
                name="Wavenumber: l",
                description="Second wave number",
                min = 0, max = 2, default = 2)
                
        IntProperty( attr="wf_3",
                name="Wavenumber: m",
                description="Third wave number",
                min = 0, max = 2, default = 2)
        
        FloatProperty( attr="last_run",
                name="It took this long",
                description="A time",
                min = 0.0, max = 1000.0, default = 0.0)
                
        BoolProperty( attr='func_abs',
                name="ABS value",
                description="Use ABS")
                
        #########################################
        layout = self.layout

        narrowui = bpy.context.user_preferences.view.properties_width_check

        wide_ui = context.region.width > narrowui

        split = layout.split()        
        col = split.column()
        
        col.prop(scene, "file_name", text='Enter a filename')
        col.separator()
        sub=col.row()
        sub.prop(scene, "res_grid", text='Resolution')
        sub.prop(scene, "scale_grid", text='Scale for grid')
        col.separator()   
        col.prop(scene, "func_sys", text='Function to plot')

        if scene.func_sys != 'wav':
            col.prop(scene, "coord_sys", text='Coord System')
        
        if scene.func_sys == 'norm':
            col.prop(scene, "func_plot", text='Enter A Function d({}):'.format(CONVCOORD[scene.coord_sys]))

        elif scene.func_sys == 'def':
            col.prop(scene, "func_name", text='Enter A Function Name:')
    
        else:
            col.label(text="Wavenumbers:")
            sub=col.row()
            sub.prop(scene, "wf_1", text='n')
            sub.prop(scene, "wf_2", text='l')
            sub.prop(scene, "wf_3", text='m')


        if scene.func_sys != 'wav':
            col.prop(scene, "func_abs", text='Use Absolute Values:')

        col.operator('plot4d', text='Run')
        
        if scene.last_run>0:
            col.label(text="Last run time: "+
                      format(scene.last_run, '.2f')+
                      ' seconds.')

#            
#        row.itemO('plot4d', text='Run')
#

class SCENE_OT_runplotter4d(bpy.types.Operator):
    """Runs Mat4D!"""
    bl_label = "Plot 4D"
    bl_idname = "plot4d"
    

    def invoke(self, context, event):
        coordlets = {'r':'x,y,z','s':'r,p,t','c':'r,p,z'}
        scene = context.scene
        
        c = scene.coord_sys
        abso = scene.func_abs
        
        if scene.func_sys == 'norm':
            fn = scene.func_plot
            fn = fn.replace('^','**')
            fn = 'lambda {}: '.format(coordlets[c])+fn
            print(fn)
            fn = eval(fn)
        elif scene.func_sys == 'def':
            fn = eval(scene.func_name)
        else:
            c='s'
            abso=False
            fn = wavefn(scene.wf_1,scene.wf_2,scene.wf_3)
        reso = scene.res_grid
        scale = scene.scale_grid

        
        finfile = scene.file_name

        context.scene.last_run = blenderrun(fn,reso,scale,c,finfile,abso)
        
        return('FINISHED')
    
bpy.types.register(SCENE_OT_runplotter4d)
bpy.types.register(OBJECT_PT_hello)

Wow, that was 57 characters too long, so I had to move the following reply to a new post:

Sure, it’s in the .blend, you can run it in a text editor. It’ll take me a little to look into add-ons, and I really want to know if I can directly access the voxel values.

Manbitesdog, please do keep us updated if you have anything that you can show! I’ll be glad to help if you need something.

i tried to run this from text editor
but get error on line 107

return tgui line !

how can this be corrected ?

Thanks

Thanks mate, will do.