Page MenuHome

object_laplace_lightning.py

File Metadata

Author
Thomas Eldredge (teldredge)
Created
Nov 13 2013, 3:01 PM

object_laplace_lightning.py

# ***** BEGIN GPL LICENSE BLOCK *****
#
# This program is free software; you can redistribute it and/or
# modify it under the terms of the GNU General Public License
# as published by the Free Software Foundation; either version 2
# of the License, or (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program; if not, write to the Free Software Foundation,
# Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
#
# ***** END GPL LICENCE BLOCK *****
bl_info = {
"name": "Laplacian Lightning",
"author": "teldredge",
"version": (0, 2, 4),
"blender": (2, 5, 7),
"api": 36147,
"location": "View3D > ToolShelf > Laplacian Lightning",
"description": "Lightning mesh generator using laplacian growth algorithm",
"warning": "Beta",
"wiki_url": "http://www.funkboxing.com",
"tracker_url": "http://www.funkboxing.com",
"category": "Object"}
######################################################################
######################################################################
##################### BLENDER LAPLACIAN LIGHTNING ####################
############################ teldredge ###############################
######################## www.funkboxing.com ##########################
######################################################################
######################## using algorithm from ########################
######################################################################
################## FAST SIMULATION OF LAPLACIAN GROWTH ###############
#################### http://gamma.cs.unc.edu/FRAC/ ###################
######################################################################
###################### and a few ideas ideas from ####################
######################################################################
########## FAST ANIMATION OF LIGHTNING USING AN ADAPTIVE MESH ########
################ http://gamma.cs.unc.edu/FAST_LIGHTNING/ #############
######################################################################
######################################################################
"""
-----RELEASE LOG/NOTES/PONTIFICATIONS-----
v0.1.0 - 04.11.11
basic generate functions and UI
object creation report (Custom Properties: FSLG_REPORT)
v0.2.0 - 04.15.11
started spelling laplacian right.
add curve function (not in UI) ...twisting problem
classify stroke by MAIN path, h-ORDER paths, TIP paths
jitter cells for mesh creation
add materials if present
v0.2.1 - 04.16.11
mesh classification speedup
v0.2.2 - 04.21.11
fxns to write/read array to file
restrict growth to insulator cells (object bounding box)
origin/ground defineable by object
gridunit more like 'resolution'
v0.2.3 - 04.24.11
cloud attractor object (termintates loop if hit)
secondary path orders disabled in UI (set to 1)
v0.2.4 - 04.26.11
fixed object selection in UI
will not run if required object not selected
moved to view 3d > toolbox
v0.x -
voxelize by faces for arbitrary insulator shapes
fix add curve function
user defined attractor path
environment map boundary conditions - requires Eqn. 15 from FSLG...
assign wattage at each segment?
animated arcs... ionization path
default settings for -lightning, -teslacoil, -spark/arc
blender setup -> settings text file ->
c++ for growth loop (maybe classify) ->
mesh text file -> blender visualize
multiple 'MAIN' brances for non-lightning discharges
n-symmetry option, create mirror images, snowflakes, etc...
slo-mo animation option -
build vert by vert, shapekey each step...
keyframe shapekeys for each...
quick-n'dirty animation option
animate emit settings as a strike
"""
######################################################################
######################################################################
######################################################################
import bpy
import time
import random
from math import sqrt
from mathutils import Vector
import struct
import bisect
import os.path
notZero = 0.0000000001
scn = bpy.context.scene
######################################################################
########################### UTILITY FXNS #############################
######################################################################
def dist(ax, ay, az ,bx, by, bz):
xsq = pow(ax-bx, 2)
ysq = pow(ay-by, 2)
zsq = pow(az-bz, 2)
d = sqrt(xsq+ysq+zsq)
return d
def distINT(ax, ay, az ,bx, by, bz):
xsq = pow(ax-bx, 2)
ysq = pow(ay-by, 2)
zsq = pow(az-bz, 2)
d = int(sqrt(xsq+ysq+zsq))
return d
def splitList(aList, idx):
ll = []
for x in aList:
ll.append(x[idx])
return ll
def deDupe(aList):
###---THANKS TO THIS GUY - http://www.peterbe.com/plog/uniqifiers-benchmark
noDupes = []
[noDupes.append(i) for i in aList if not noDupes.count(i)]
return noDupes
def getLowHigh(aList):
tLow = aList[0]; tHigh = aList[0]
for a in aList:
if a < tLow: tLow = a
if a > tHigh: tHigh = a
return tLow, tHigh
def weightedRandomChoice(aList):
tL = []
tweight = 0
for a in range(len(aList)):
idex = a; weight = aList[a]
if weight > 0.0:
tweight += weight
tL.append((tweight, idex))
i = bisect.bisect(tL, (random.uniform(0, tweight), None))
r = tL[i][1]
return r
def getStencil3D_26(x,y,z):
nL = []
for xT in range(x-1, x+2):
for yT in range(y-1, y+2):
for zT in range(z-1, z+2):
nL.append((xT, yT, zT))
nL.remove((x,y,z))
return nL
def jitterCells(aList, jit):
j = jit/2
bList = []
for a in aList:
ax = a[0] + random.uniform(-j, j)
ay = a[1] + random.uniform(-j, j)
az = a[2] + random.uniform(-j, j)
bList.append((ax, ay, az))
return bList
def addReportProp(ob, str):
bpy.types.Object.FSLG_REPORT = bpy.props.StringProperty(
name = 'fslg_report', default = '')
ob.FSLG_REPORT = str
######################################################################
######################## VISUALIZATION FXNS ##########################
######################################################################
def writeArrayToVoxel(arr, filename):
gridS = 64
half = int(gridS/2)
bitOn = 255
aGrid = [[[0 for z in range(gridS)] for y in range(gridS)] for x in range(gridS)]
for a in arr:
try:
aGrid[a[0]+half][a[1]+half][a[2]+half] = bitOn
except:
print('Particle beyond voxel domain')
file = open(filename, "wb")
for z in range(gridS):
for y in range(gridS):
for x in range(gridS):
file.write(struct.pack('B', aGrid[x][y][z]))
file.flush()
file.close()
def writeArrayToFile(arr, filename):
file = open(filename, "w")
for a in arr:
tstr = str(a[0]) + ',' + str(a[1]) + ',' + str(a[2]) + '\n'
file.write(tstr)
file.close
def readArrayFromFile(filename):
file = open(filename, "r")
arr = []
for f in file:
pt = f[0:-1].split(',')
arr.append((int(pt[0]), int(pt[1]), int(pt[2])))
return arr
def makeMeshCube(msize):
msize = msize/2
mmesh = bpy.data.meshes.new('q')
mmesh.vertices.add(8)
mmesh.vertices[0].co = [-msize, -msize, -msize]
mmesh.vertices[1].co = [-msize, msize, -msize]
mmesh.vertices[2].co = [ msize, msize, -msize]
mmesh.vertices[3].co = [ msize, -msize, -msize]
mmesh.vertices[4].co = [-msize, -msize, msize]
mmesh.vertices[5].co = [-msize, msize, msize]
mmesh.vertices[6].co = [ msize, msize, msize]
mmesh.vertices[7].co = [ msize, -msize, msize]
mmesh.faces.add(6)
mmesh.faces[0].vertices_raw = [0,1,2,3]
mmesh.faces[1].vertices_raw = [0,4,5,1]
mmesh.faces[2].vertices_raw = [2,1,5,6]
mmesh.faces[3].vertices_raw = [3,2,6,7]
mmesh.faces[4].vertices_raw = [0,3,7,4]
mmesh.faces[5].vertices_raw = [5,4,7,6]
mmesh.update(calc_edges=True)
return(mmesh)
def writeArrayToCubes(arr, gridBU, cBOOL = False, jBOOL = False):
for a in arr:
x = a[0]; y = a[1]; z = a[2]
me = makeMeshCube(gridBU)
ob = bpy.data.objects.new('xCUBE', me)
ob.location.x = x*gridBU
ob.location.y = y*gridBU
ob.location.z = z*gridBU
if cBOOL:
###---POS+BLUE, NEG-RED, ZERO:BLACK
col = (1.0, 1.0, 1.0, 1.0)
if a[3] == 0: col = (0.0, 0.0, 0.0, 1.0)
if a[3] < 0: col = (-a[3], 0.0, 0.0, 1.0)
if a[3] > 0: col = (0.0, 0.0, a[3], 1.0)
ob.color = col
bpy.context.scene.objects.link(ob)
bpy.context.scene.update()
if jBOOL:
### CAN'T JOIN ALL CUBES TO A SINGLE MESH RIGHT... ARGH...
### bpy.ops.object.join()
for q in bpy.context.scene.objects:
q.select = False
if q.name[0:5] == 'xCUBE':
q.select = True
bpy.context.scene.objects.active = q
def addVert(ob, pt, conni = -1):
mmesh = ob.data
mmesh.vertices.add(1)
vcounti = len(mmesh.vertices)-1
mmesh.vertices[vcounti].co = [pt[0], pt[1], pt[2]]
if conni > -1:
mmesh.edges.add(1)
ecounti = len(mmesh.edges)-1
mmesh.edges[ecounti].vertices = [conni, vcounti]
mmesh.update()
def addEdge(ob, va, vb):
mmesh = ob.data
mmesh.edges.add(1)
ecounti = len(mmesh.edges)-1
mmesh.edges[ecounti].vertices = [va, vb]
mmesh.update()
def newMesh(mname):
mmesh = bpy.data.meshes.new(mname)
omesh = bpy.data.objects.new(mname, mmesh)
bpy.context.scene.objects.link(omesh)
return omesh
def writeArrayToMesh(mname, arr, gridBU, rpt = None):
mob = newMesh(mname)
mob.scale = (gridBU, gridBU, gridBU)
if rpt: addReportProp(mob, rpt)
addVert(mob, arr[0], -1)
for ai in range(1, len(arr)):
a = arr[ai]
addVert(mob, a, ai-1)
return mob
### SOME PROBLEM WITH IT ADDING (0,0,0)
def writeArrayToCurves(cname, arr, gridBU, bd = .05, rpt = None):
cur = bpy.data.curves.new('fslg_curve', 'CURVE')
cur.use_fill_front = False
cur.use_fill_back = False
cur.bevel_depth = bd
cur.bevel_resolution = 2
cob = bpy.data.objects.new(cname, cur)
cob.scale = (gridBU, gridBU, gridBU)
if rpt: addReportProp(cob, rpt)
bpy.context.scene.objects.link(cob)
cur.splines.new('BEZIER')
cspline = cur.splines[0]
div = 1 ### SPACING FOR HANDLES (2 - 1/2 WAY, 1 - NEXT BEZIER
for a in range(len(arr)):
cspline.bezier_points.add(1)
bp = cspline.bezier_points[len(cspline.bezier_points)-1]
if a-1 < 0: hL = arr[a]
else:
hx = arr[a][0] - ((arr[a][0]-arr[a-1][0]) / div)
hy = arr[a][1] - ((arr[a][1]-arr[a-1][1]) / div)
hz = arr[a][2] - ((arr[a][2]-arr[a-1][2]) / div)
hL = (hx,hy,hz)
if a+1 > len(arr)-1: hR = arr[a]
else:
hx = arr[a][0] + ((arr[a+1][0]-arr[a][0]) / div)
hy = arr[a][1] + ((arr[a+1][1]-arr[a][1]) / div)
hz = arr[a][2] + ((arr[a+1][2]-arr[a][2]) / div)
hR = (hx,hy,hz)
bp.co = arr[a]
bp.handle_left = hL
bp.handle_right = hR
def addArrayToMesh(mob, arr):
addVert(mob, arr[0], -1)
mmesh = mob.data
vcounti = len(mmesh.vertices)-1
for ai in range(1, len(arr)):
a = arr[ai]
addVert(mob, a, len(mmesh.vertices)-1)
def addMaterial(ob, matname):
mat = bpy.data.materials[matname]
ob.active_material = mat
def writeStokeToMesh(arr, MAINi, HORDERi, TIPSi, orig, gs, rpt=None):
###---JITTER CELLS
jarr = jitterCells(arr, 1)
###---MAIN BRANCH
print(' WRITING MAIN BRANCH')
llmain = []
for x in MAINi:
llmain.append(jarr[x])
mob = writeArrayToMesh('la0MAIN', llmain, gs)
mob.location = orig
#writeArrayToCurves('laMAIN', llmain, .10, .25)
###---hORDER BRANCHES
for hOi in range(len(HORDERi)):
print(' WRITING ORDER', hOi)
hO = HORDERi[hOi]
hob = newMesh('la1H'+str(hOi))
for y in hO:
llHO = []
for x in y:
llHO.append(jarr[x])
addArrayToMesh(hob, llHO)
hob.scale = (gs, gs, gs)
hob.location = orig
###---TIPS
print(' WRITING TIP PATHS')
tob = newMesh('la2TIPS')
for y in TIPSi:
llt = []
for x in y:
llt.append(jarr[x])
addArrayToMesh(tob, llt)
tob.scale = (gs, gs, gs)
tob.location = orig
###---ADD MATERIALS TO OBJECTS (IF THEY EXIST)
try:
addMaterial(mob, 'edgeMAT-h0')
addMaterial(hob, 'edgeMAT-h1')
addMaterial(tob, 'edgeMAT-h2')
print(' ADDED MATERIALS')
except: print(' MATERIALS NOT FOUND')
###---ADD GENERATION REPORT
if rpt:
addReportProp(mob, rpt)
addReportProp(hob, rpt)
addReportProp(tob, rpt)
def visualizeArray(cg, origin, gs, vm, vc, vv, rst):
###---READ/WRITE ARRAY TO FILE
#tfile = 'c:\\testarr.txt'
#writeArrayToFile(cg, tfile)
#cg = readArrayFromFile(tfile)
###---WRITE ARRAY TO MESH
if vm:
aMi, aHi, aTi = classifyStroke(cg, scn.HORDER)
print(':::WRITING TO MESH')
writeStokeToMesh(cg, aMi, aHi, aTi, origin, gs, rst)
print(':::MESH WRITTEN')
###---WRITE ARRAY TO CUBE OBJECTS
if vc:
print(' WRITING TO CUBES')
writeArrayToCubes(cg, gs, False, True)
print(' CUBITS WRITTEN')
###---WRITE ARRAY TO VOXEL DATA FILE
if vv:
print(' WRITING TO VOXELS')
fname = "FSLGvoxels.raw"
path = os.path.dirname(bpy.data.filepath)
writeArrayToVoxel(cg, path + "\\" + fname)
print(' VOXEL DATA WRITTEN')
######################################################################
########################### ALGORITHM FXNS ###########################
########################## FROM OTHER PAPER ##########################
######################### AND STUFF I MADE UP ########################
######################################################################
def buildCPGraph(arr):
###---IN -XYZ ARRAY AS BUILT BY GENERATOR
###---OUT -[(CHILDindex, PARENTindex)]
sgarr = []
sgarr.append((1, 0)) #
for ai in range(2, len(arr)):
cs = arr[ai]
cpts = arr[0:ai]
cslap = getStencil3D_26(cs[0], cs[1], cs[2])
for nc in cslap:
ct = cpts.count(nc)
if ct>0:
cti = cpts.index(nc)
sgarr.append((ai, cti))
return sgarr
def findChargePath(oc, fc, ngraph, restrict = [], partial = True):
###---oc -ORIGIN CHARGE INDEX, fc -FINAL CHARGE INDEX
###---ngraph -NODE GRAPH, restrict- INDEX OF SITES CANNOT TRAVERSE
###---partial -RETURN PARTIAL PATH IF RESTRICTION ENCOUNTERD
cList = splitList(ngraph, 0)
pList = splitList(ngraph, 1)
aRi = []
cNODE = fc
for x in range(len(ngraph)):
pNODE = pList[cList.index(cNODE)]
aRi.append(cNODE)
cNODE = pNODE
npNODECOUNT = cList.count(pNODE)
if cNODE == oc: ### STOP IF ORIGIN FOUND
aRi.append(cNODE) ### RETURN PATH
return aRi
if npNODECOUNT == 0: ### STOP IF NO PARENTS
return [] ### RETURN []
if pNODE in restrict: ### STOP IF PARENT IS IN RESTRICTION
if partial: ### RETURN PARTIAL OR []
aRi.append(cNODE)
return aRi
else: return []
def findTips(arr):
lt = []
for ai in arr[0:len(arr)-1]:
a = ai[0]
cCOUNT = 0
for bi in arr:
b = bi[1]
if a == b:
cCOUNT += 1
if cCOUNT == 0:
lt.append(a)
return lt
def findChannelRoots(path, ngraph, restrict = []):
roots = []
for ai in range(len(ngraph)):
chi = ngraph[ai][0]
par = ngraph[ai][1]
if par in path and not chi in path and \
not chi in restrict:
roots.append(par)
droots = deDupe(roots)
return droots
def findChannels(roots, tips, ngraph, restrict):
cPATHS = []
for ri in range(len(roots)):
r = roots[ri]
sL = 1
sPATHi = []
for ti in range(len(tips)):
t = tips[ti]
if t < r: continue
tPATHi = findChargePath(r, t, ngraph, restrict, False)
tL = len(tPATHi)
if tL > sL:
if countChildrenOnPath(tPATHi, ngraph) > 1:
sL = tL
sPATHi = tPATHi
tTEMP = t; tiTEMP = ti
if len(sPATHi) > 0:
print(' found path/idex from', r, ri, 'of',
len(roots), 'possible | (tip/index)', tTEMP, tiTEMP)
cPATHS.append(sPATHi)
tips.remove(tTEMP)
return cPATHS
def countChildrenOnPath(aPath, ngraph, quick = True):
###---RETURN HOW MANY BRANCHES
### COUNT WHEN NODE IS A PARENT >1 TIMES
### quick -STOP AND RETURN AFTER FIRST
cCOUNT = 0
pList = splitList(ngraph,1)
for ai in range(len(aPath)-1):
ap = aPath[ai]
pc = pList.count(ap)
if quick and pc > 1:
return pc
return cCOUNT
###---CLASSIFY CHANNELS INTO 'MAIN', 'hORDER/SECONDARY' and 'SIDE'
def classifyStroke(sarr, hORDER = 1):
print(':::CLASSIFYING STROKE')
###---BUILD CHILD/PARENT GRAPH (indexes of sarr)
sgarr = buildCPGraph(sarr)
###---FIND MAIN CHANNEL
print(' finding MAIN')
oCharge = sgarr[0][1]
fCharge = sgarr[len(sgarr)-1][0]
aMAINi = findChargePath(oCharge, fCharge, sgarr)
###---FIND TIPS
print(' finding TIPS')
aTIPSi = findTips(sgarr)
###---FIND hORDER CHANNEL ROOTS
### hCOUNT = ORDERS BEWTEEN MAIN and SIDE/TIPS
### !!!STILL BUGGY!!!
hRESTRICT = list(aMAINi) ### add to this after each time...
allHPATHSi = [] ### all hO paths, [[h0], [h1]...]
curPATHSi = [aMAINi] ### list of paths find roots on
for h in range(hORDER):
allHPATHSi.append([])
for pi in range(len(curPATHSi)): ### loop through all paths in this order
p = curPATHSi[pi]
### get roots for this path
aHROOTSi = findChannelRoots(p, sgarr, hRESTRICT)
print(' found', len(aHROOTSi), 'roots in ORDER', h, ':#paths:', len(curPATHSi))
### get channels for these roots
if len(aHROOTSi) == 0:
print('NO ROOTS FOR FOUND FOR CHANNEL')
aHPATHSi = []
continue
else:
aHPATHSiD = findChannels(aHROOTSi, aTIPSi, sgarr, hRESTRICT)
aHPATHSi = aHPATHSiD
allHPATHSi[h] += aHPATHSi
### set these channels as restrictions for next iterations
for hri in aHPATHSi:
hRESTRICT += hri
curPATHSi = aHPATHSi
###---SIDE BRANCHES, FINAL ORDER OF HEIRARCHY
### FROM TIPS THAT ARE NOT IN AN EXISTING PATH
### BACK TO ANY OTHER POINT THAT IS ALREADY ON A PATH
aDRAWNi = []
aDRAWNi += aMAINi
for oH in allHPATHSi:
for o in oH:
aDRAWNi += o
aTPATHSi = []
for a in aTIPSi:
if not a in aDRAWNi:
aPATHi = findChargePath(oCharge, a, sgarr, aDRAWNi)
aDRAWNi += aPATHi
aTPATHSi.append(aPATHi)
return aMAINi, allHPATHSi, aTPATHSi
######################################################################
########################### ALGORITHM FXNS ###########################
############################## FROM FSLG #############################
######################################################################
def getGrowthProbability(uN, aList):
###---IN: uN -USER TERM, cList -CANDIDATE SITES, oList -CANDIDATE SITE CHARGES
### OUT: LIST OF [(XYZ), POT, PROB]
cList = splitList(aList, 0)
oList = splitList(aList, 1)
Omin, Omax = getLowHigh(oList)
if Omin == Omax: Omax += notZero; Omin -= notZero
PdL = []
E = 0
E = notZero ###===DIVISOR FOR (FSLG - Eqn. 12)
for o in oList:
Uj = (o - Omin) / (Omax - Omin) ###===(FSLG - Eqn. 13)
E += pow(Uj, uN)
for oi in range(len(oList)):
o = oList[oi]
Ui = (o - Omin) / (Omax - Omin)
Pd = (pow(Ui, uN)) / E ###===(FSLG - Eqn. 12)
PdINT = Pd * 100
PdL.append(Pd)
return PdL
def updatePointCharges(p, cList, eList = []):
###---IN: pNew -NEW GROWTH CELL
### cList -OLD CANDIDATE SITES, eList -SAME
### OUT: LIST OF NEW CHARGE AT CANDIDATE SITES
r1 = 1/2 ###===(FSLG - Eqn. 10)
nOiL = []
for oi in range(len(cList)):
o = cList[oi][1]
c = cList[oi][0]
iOe = 0
#rit = distINT(c[0], c[1], c[2], p[0], p[1], p[2])
rit = dist(c[0], c[1], c[2], p[0], p[1], p[2])
iOe += (1 - (r1/rit))
Oit = o + iOe
nOiL.append((c, Oit))
return nOiL
def initialPointCharges(pList, cList, eList = []):
###---IN: p -CHARGED CELL (XYZ), cList -CANDIDATE SITES (XYZ, POT, PROB)
### OUT: cList -WITH POTENTIAL CALCULATED
r1 = 1/2 ###===(Eqn. 10)
npList = []
for p in pList:
npList.append(((p[0], p[1], p[2]), 1.0))
for e in eList:
npList.append(((e[0], e[1], e[2]), e[3]))
OiL = []
for i in cList:
Oi = 0
for j in npList:
if i != j[0]:
#rij = distINT(i[0], i[1], i[2], j[0][0], j[0][1], j[0][2])
rij = dist(i[0], i[1], i[2], j[0][0], j[0][1], j[0][2])
Oi += (1 - (r1 / rij)) * j[1] ### CHARGE INFLUENCE
OiL.append(((i[0], i[1], i[2]), Oi))
return OiL
def fakeGroundChargePlane(z, charge):
eCL = []
xy = abs(z)/2
eCL += [(0, 0, z, charge)]
eCL += [(xy, 0, z, charge)]
eCL += [(0, xy, z, charge)]
eCL += [(-xy, 0, z, charge)]
eCL += [(0, -xy, z, charge)]
return eCL
def getCellPlane(pt, l):
###---POINT, LENGTH OF A SIDE / 2
ll = []
len = int(l/2)
x = pt[0]; y = pt[1]; z = pt[2]
for xT in range(x-(len-2), x+(len-1)):
for yT in range(y-(len-2), y+(len-1)):
ll.append((xT, yT, z))
return ll
def getCandidateSites(aList, iList = []):
###---IN: aList -(X,Y,Z) OF CHARGED CELL SITES, iList -insulator sites
### OUT: CANDIDATE LIST OF GROWTH SITES [(X,Y,Z)]
cList = []
for c in aList:
tempList = getStencil3D_26(c[0], c[1], c[2])
for t in tempList:
if not t in aList and not t in iList:
cList.append(t)
ncList = deDupe(cList)
return ncList
def insulatorCells(orig, ob, gs):
###---RETURN GRID CELLS FROM BOUNDING BOX
### OF AN INSULATOR OBJECT
### (RETURNS SOLID, QUICKER IF IT WAS HULL)
cloc = ob.location
dims = ob.dimensions
cx1 = int( ((cloc[0] - (dims[0] / 2 ) ) - orig[0] ) / gs )
cx2 = int( ((cloc[0] + (dims[0] / 2 ) ) - orig[0] ) / gs )
cy1 = int( ((cloc[1] - (dims[1] / 2 ) ) - orig[1] ) / gs )
cy2 = int( ((cloc[1] + (dims[1] / 2 ) ) - orig[1] ) / gs )
cz1 = int( ((cloc[2] - (dims[2] / 2 ) ) - orig[2] ) / gs )
cz2 = int( ((cloc[2] + (dims[2] / 2 ) ) - orig[2] ) / gs )
ll = []
for x in range(cx1, cx2):
for y in range(cy1, cy2):
for z in range(cz1, cz2):
ll.append((x,y,z))
return ll
def cloudCells(orig, ob, gs, charge):
###---RETURN GRID CELLS FROM BOUNDING BOX (JUST PERMIMETER)
cloc = ob.location
dims = ob.dimensions
cx1 = int( ((cloc[0] - (dims[0] / 2 ) ) - orig[0] ) / gs )
cx2 = int( ((cloc[0] + (dims[0] / 2 ) ) - orig[0] ) / gs )
cy1 = int( ((cloc[1] - (dims[1] / 2 ) ) - orig[1] ) / gs )
cy2 = int( ((cloc[1] + (dims[1] / 2 ) ) - orig[1] ) / gs )
cz1 = int( ((cloc[2] - (dims[2] / 2 ) ) - orig[2] ) / gs )
cz2 = int( ((cloc[2] + (dims[2] / 2 ) ) - orig[2] ) / gs )
llc = []; ll = []
for x in range(cx1, cx2):
for y in range(cy1, cy2):
for z in range(cz1, cz2):
llc.append((x,y,z,charge))
ll.append((x,y,z))
return llc, ll
def obList():
ll = []
for obi in range(len(bpy.data.objects)):
ob = bpy.data.objects[obi]
ll.append((ob.name, ob.name, ob.name))
return ll
def setupObjects():
oOB = bpy.data.objects.new('ELorigin', None)
oOB.location = ((0,0,10))
bpy.context.scene.objects.link(oOB)
gOB = bpy.data.objects.new('ELground', None)
gOB.empty_draw_type = 'ARROWS'
bpy.context.scene.objects.link(gOB)
cME = makeMeshCube(1)
cOB = bpy.data.objects.new('ELcloud', cME)
cOB.location = ((-2,8,12))
bpy.context.scene.objects.link(cOB)
iME = makeMeshCube(1)
iOB = bpy.data.objects.new('ELinsulator', iME)
iOB.location = ((0,0,5))
iOB.scale = ((5,5,.1))
bpy.context.scene.objects.link(iOB)
try:
scn.OOB = 'ELorigin'
scn.GOB = 'ELground'
scn.COB = 'ELcloud'
scn.IOB = 'ELinsulator'
except: pass
def checkSettings():
check = True
if scn.OOB == "":
print('ERROR: NO ORIGIN OBJECT SELECTED')
check = False
if scn.GROUNDBOOL and scn.GOB == "":
print('ERROR: NO GROUND OBJECT SELECTED')
check = False
if scn.CLOUDBOOL and scn.COB == "":
print('ERROR: NO CLOUD OBJECT SELECTED')
check = False
if scn.IBOOL and scn.IOB == "":
print('ERROR: NO INSULATOR OBJECT SELECTED')
check = False
#make a popup here
return check
#def meshToCells(orig, ob, gs, charge)
######################################################################
############################### MAIN #################################
######################################################################
def FSLG():
print('\n<<<<<<------GO GO GADGET: FAST SIMULATION OF LAPLACIAN GROWTH!')
tc1 = time.clock()
TSTEPS = scn.TSTEPS
obORIGIN = scn.objects[scn.OOB]
obGROUND = scn.objects[scn.GOB]
scn.ORIGIN = obORIGIN.location
scn.GROUNDZ = int((obGROUND.location[2] - scn.ORIGIN[2]) / scn.GSCALE)
###======FAST SIMULATION OF LAPLACIAN GROWTH======###
###====== 1) INSERT INTIAL CHARGE POINT
cgrid = [(0, 0, 0)]
###---GROUND CHARGE CELL / INSULATOR LISTS (eChargeList/icList)
eChargeList = []; icList = []; cloudList = []
if scn.GROUNDBOOL:
eChargeList = fakeGroundChargePlane(scn.GROUNDZ, scn.GROUNDC)
if scn.CLOUDBOOL:
obCLOUD = scn.objects[scn.COB]
eChargeList, cloudList = cloudCells(scn.ORIGIN, obCLOUD, scn.GSCALE, scn.CLOUDC)
if scn.IBOOL:
obINSULATOR = scn.objects[scn.IOB]
icList = insulatorCells(scn.ORIGIN, obINSULATOR, scn.GSCALE)
###====== 2) LOCATE CANDIDATE SITES AROUND CHARGE
cSites = getCandidateSites(cgrid, icList)
###====== 3) CALC POTENTIAL AT EACH SITE (Eqn. 10)
cSites = initialPointCharges(cgrid, cSites, eChargeList)
ts = 1
while ts <= TSTEPS:
###====== 1) SELECT NEW GROWTH SITE (Eqn. 12)
###===GET PROBABILITIES AT CANDIDATE SITES
gProbs = getGrowthProbability(scn.BIGVAR, cSites)
###===CHOOSE NEW GROWTH SITE BASED ON PROBABILITIES
gSitei = weightedRandomChoice(gProbs)
gsite = cSites[gSitei][0]
###====== 2) ADD NEW POINT CHARGE AT GROWTH SITE
###===ADD NEW GROWTH CELL TO GRID
cgrid.append(gsite)
###===REMOVE NEW GROWTH CELL FROM CANDIDATE SITES
cSites.remove(cSites[gSitei])
###====== 3) UPDATE POTENTIAL AT CANDIDATE SITES (Eqn. 11)
cSites = updatePointCharges(gsite, cSites, eChargeList)
###====== 4) ADD NEW CANDIDATES SURROUNDING GROWTH SITE
###===GET CANDIDATE 'STENCIL'
ncSitesT = getCandidateSites([gsite], icList)
###===REMOVE CANDIDATES ALREADY IN CANDIDATE LIST OR CHARGE GRID
ncSites = []
cSplit = splitList(cSites, 0)
for cn in ncSitesT:
if not cn in cSplit and \
not cn in cgrid:
ncSites.append((cn, 0))
###====== 5) CALC POTENTIAL AT NEW CANDIDATE SITES (Eqn. 10)
ncSplit = splitList(ncSites, 0)
ncSites = initialPointCharges(cgrid, ncSplit, eChargeList)
###===ADD NEW CANDIDATE SITES TO CANDIDATE LIST
for ncs in ncSites:
cSites.append(ncs)
###===ITERATION COMPLETE
istr1 = ':::T-STEP: ' + str(ts) + '/' + str(TSTEPS)
istr12 = ' | GROUNDZ: ' + str(scn.GROUNDZ) + ' | '
istr2 = 'CANDS: ' + str(len(cSites)) + ' | '
istr3 = 'GSITE: ' + str(gsite)
print(istr1 + istr12 + istr2 + istr3)
ts += 1
###---EARLY TERMINATION FOR GROUND/CLOUD STRIKE
if scn.GROUNDBOOL:
if gsite[2] == scn.GROUNDZ:
ts = TSTEPS+1
print('<<<<<<------EARLY TERMINATION DUE TO GROUNDSTRIKE')
continue
if scn.CLOUDBOOL:
if gsite in cloudList:
ts = TSTEPS+1
print('<<<<<<------EARLY TERMINATION DUE TO CLOUDSTRIKE')
continue
tc2 = time.clock()
tcRUN = tc2 - tc1
print('<<<<<<------LAPLACIAN GROWTH LOOP COMPLETED: ' + str(len(cgrid)) + ' / ' + str(tcRUN)[0:5] + ' SECONDS')
print('<<<<<<------VISUALIZING DATA')
reportSTRING = getReportString(tcRUN)
###---VISUALIZE ARRAY
visualizeArray(cgrid, scn.ORIGIN, scn.GSCALE, scn.VMESH, scn.VCUBE, scn.VVOX, reportSTRING)
print('<<<<<<------COMPLETE')
######################################################################
################################ GUI #################################
######################################################################
###---NOT IN UI
bpy.types.Scene.ORIGIN = bpy.props.FloatVectorProperty(name = "origin charge")
bpy.types.Scene.GROUNDZ = bpy.props.IntProperty(name = "ground Z coordinate")
###---IN UI
bpy.types.Scene.TSTEPS = bpy.props.IntProperty(name = "iterations")
bpy.types.Scene.GSCALE = bpy.props.FloatProperty(name = "grid unit size")
bpy.types.Scene.BIGVAR = bpy.props.FloatProperty(name = "straightness")
bpy.types.Scene.GROUNDBOOL = bpy.props.BoolProperty(name = "use ground object")
bpy.types.Scene.GROUNDC = bpy.props.IntProperty(name = "ground charge")
bpy.types.Scene.CLOUDBOOL = bpy.props.BoolProperty(name = "use cloud object")
bpy.types.Scene.CLOUDC = bpy.props.IntProperty(name = "cloud charge")
bpy.types.Scene.HORDER = bpy.props.IntProperty(name = "secondary paths orders")
bpy.types.Scene.VMESH = bpy.props.BoolProperty(name = "mesh")
bpy.types.Scene.VCUBE = bpy.props.BoolProperty(name = "cubes")
bpy.types.Scene.VVOX = bpy.props.BoolProperty(name = "voxel (experimental)")
bpy.types.Scene.IBOOL = bpy.props.BoolProperty(name = "use insulator object")
bpy.types.Scene.OOB = bpy.props.StringProperty()
bpy.types.Scene.GOB = bpy.props.StringProperty()
bpy.types.Scene.COB = bpy.props.StringProperty()
bpy.types.Scene.IOB = bpy.props.StringProperty()
###---DEFAULT USER SETTINGS
scn.TSTEPS = 350
scn.HORDER = 1
scn.GSCALE = 0.15
scn.BIGVAR = 6.3
scn.GROUNDBOOL = True
scn.GROUNDC = -250
scn.CLOUDBOOL = False
scn.CLOUDC = -1
scn.VMESH = True
scn.VCUBE = False
scn.VVOX = False
scn.IBOOL = False
try:
scn.OOB = "ELorigin"
scn.GOB = "ELground"
scn.COB = "ELcloud"
scn.IOB = "ELinsulator"
except: pass
### TESTING
if False:
#if True:
scn.TSTEPS = 500
scn.GSCALE = 0.20
scn.BIGVAR = 6.3
scn.GROUNDC = -250
class runFSLGLoopOperator(bpy.types.Operator):
'''By The Mighty Hammer Of Thor!!!'''
bl_idname = "object.runfslg_operator"
bl_label = "run FSLG Loop Operator"
def execute(self, context):
if checkSettings():
FSLG()
else: pass
return {'FINISHED'}
class setupObjectsOperator(bpy.types.Operator):
'''Create Origin/Ground/Cloud/Insulator Objects'''
bl_idname = "object.setup_objects_operator"
bl_label = "Setup Objects Operator"
def execute(self, context):
setupObjects()
return {'FINISHED'}
class OBJECT_PT_fslg(bpy.types.Panel):
bl_label = "Laplacian Lightning"
bl_space_type = "VIEW_3D"
bl_region_type = "TOOLS"
bl_context = "objectmode"
def draw(self, context):
scn = context.scene
layout = self.layout
colR = layout.column()
#row1 = layout.row()
#colL = row1.column()
#colR = row1.column()
colR.prop(scn, 'TSTEPS')
colR.prop(scn, 'GSCALE')
colR.prop(scn, 'BIGVAR')
colR.operator('object.setup_objects_operator', text = 'create setup objects')
colR.label('origin object')
colR.prop_search(scn, "OOB", context.scene, "objects")
colR.prop(scn, 'GROUNDBOOL')
colR.prop_search(scn, "GOB", context.scene, "objects")
colR.prop(scn, 'GROUNDC')
colR.prop(scn, 'CLOUDBOOL')
colR.prop_search(scn, "COB", context.scene, "objects")
colR.prop(scn, 'CLOUDC')
colR.prop(scn, 'IBOOL')
colR.prop_search(scn, "IOB", context.scene, "objects")
colR.operator('object.runfslg_operator', text = 'generate lightning')
#col.prop(scn, 'HORDER')
colR.prop(scn, 'VMESH')
colR.prop(scn, 'VCUBE')
colR.prop(scn, 'VVOX')
def getReportString(rtime):
rSTRING1 = 't:' + str(scn.TSTEPS) + ',sc:' + str(scn.GSCALE)[0:4] + ',uv:' + str(scn.BIGVAR)[0:4] + ','
rSTRING2 = 'ori:' + str(scn. ORIGIN[0]) + '/' + str(scn. ORIGIN[1]) + '/' + str(scn. ORIGIN[2]) + ','
rSTRING3 = 'gz:' + str(scn.GROUNDZ) + ',gc:' + str(scn.GROUNDC) + ',rtime:' + str(int(rtime))
return rSTRING1 + rSTRING2 + rSTRING3
def register():
bpy.utils.register_class(runFSLGLoopOperator)
bpy.utils.register_class(setupObjectsOperator)
bpy.utils.register_class(OBJECT_PT_fslg)
def unregister():
bpy.utils.unregister_class(runFSLGLoopOperator)
bpy.utils.unregister_class(setupObjectsOperator)
bpy.utils.unregister_class(OBJECT_PT_fslg)
if __name__ == "__main__":
### RUN FOR TESTING
#FSLG()
### UI
register()
pass

Event Timeline