#!BPY
# vim:fileencoding=utf-8

################################################################################
# Demolition Script, Copyright (C) 2007 Kai Kostack
################################################################################
# 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. (www.gnu.org)
################################################################################

__author__  = 'Kai Kostack, Mathias Panzenböck'
__version__ = '1.0.71001 Beta-panzi2'
__url__     = 'www.KaiKostack.com'
__bpydoc__  = """\

This script simulates realistic demolition effects on meshes.
For more instructions see the Demolition manual.

If you like this software or you're using it commercially,
you might consider to credit me and to support my work with a donation. Thank you.

Please report bugs!


Known problems / todo
---------------------

- vmblur doesn't work with animated mode
- elast works not so well with external forces like gravity and wind
- self collision needs improvements
- multicore support
- unified baking system support


Additional technical notes
--------------------------

Following arrays will stay in memory in between frames:

						#	type     description
Blender.SpdBufGr1		#	[Vector] speed vector buffer (size = vertex count *object count)
Blender.SpdBufGr2		#	[Vector]
Blender.AniSpdBufGr1	#	[Vector] additional speed vector buffer for animated mode
Blender.AniSpdBufGr2	#	[Vector]
Blender.EdgBufGr1		#	[Vector] edge length/vector buffer (size = edge count *object count)
Blender.EdgBufGr2		#	[Vector]
Blender.CollGr1			#	[int]    collision flag (size = vertex count *object count),
Blender.CollGr2			#	[int]    0 = free, 1 = collision occured, 2 = fixed to animated object
Blender.CollLimGr1		#	[int]    collision counter (size = vertex count *object count)
Blender.CollLimGr2		#	[int]
Blender.IdxBufGr1		#	[int]    vertex index buffer for animated mode (size = vertex count *object count)
Blender.IdxBufGr2		#	[int]
	
"""

print "You are running Demolition version %s now..." % __version__

#__all__ = [
#	'main', 'saveIDProperties', 'loadIDProperties',
#	'setProperties', 'getProperties', 'getProperty', 'setProperty'
#]

################################################################################

import Blender
from Blender import *
from Blender.Mathutils import *
from math import *
import sys
from Blender.Ipo import OB_ROTX, OB_ROTY, OB_ROTZ, OB_LOCX, OB_LOCY, OB_LOCZ, \
                        OB_DLOCX, OB_DLOCY, OB_DLOCZ

# I import count as icount so there are no name
# conflicts with the common name 'count':
from itertools import izip, count as icount

MAX_FLOAT  = 1.7976931348623157E308
MIN_FLOAT  = -MAX_FLOAT
EPS_FLOAT  = 4.9E-324
MAX_INT    = sys.maxint
MIN_INT    = -MAX_INT - 1

VEC_MAX_FLOAT = 10000. #sqrt(MAX_FLOAT) / 3.
MIN_VEC    = (-VEC_MAX_FLOAT, -VEC_MAX_FLOAT, -VEC_MAX_FLOAT)
MAX_VEC    = (VEC_MAX_FLOAT, VEC_MAX_FLOAT, VEC_MAX_FLOAT)
MAX_STRLEN = 399

try:
	from profiler import Profiler

	if not hasattr(Blender,'DemolitionProfiler'):
		print "new profiler"
		Blender.DemolitionProfiler = Profiler()
	
	def DemolitionProfiler_clear():
		Blender.DemolitionProfiler.clear()

	def DemolitionProfiler_printInfos():
		Blender.DemolitionProfiler.printInfos()

	def DemolitionProfiler_addScope(scope):
		Blender.DemolitionProfiler.addScope(scope)

	def DemolitionProfiler_add(funct):
		return Blender.DemolitionProfiler.add(funct)
	
	DemolitionProfiler_addScope = Blender.DemolitionProfiler.add(DemolitionProfiler_addScope)
	DemolitionProfiler_add      = Blender.DemolitionProfiler.add(DemolitionProfiler_add)

except ImportError:
	def DemolitionProfiler_clear():
		pass

	def DemolitionProfiler_printInfos():
		pass

	def DemolitionProfiler_addScope(scope):
		pass

	def DemolitionProfiler_add(funct):
		return funct

################################################################################
##########################  Begin of constants part  ###########################

## crash group names
gr1name = "CrashGr1"
gr2name = "CrashGr2"
suffix = " Sim"

vmblur = 0		# 0			=> Tell me if you're doing a final render with vector motion blur, 1 = on, 0 = off
				# Vector Blur needs some special treatment, because the script is called three times per frame then.
				# Don't forget to activate it or you will get completely different results for the same settings.
				# Warning: Do not use motion blur together with animated mode. It's not compatible yet.

sleep = 1		# 1			=> Do no calculations beside speed vectors until (and including) this frame number

breaks = 1		# 1			=> Meshes breakable, 1 = on, 0 = off
limit = .5		# .5		=> Edge break limit (*100 % edge length)
limitTrans = .1	# .1		=> Edge break limit for tranparent faces (slower), 0 = disabled
delEdges = 0	# 0			=> Delete faceless (invisible) edges (slower), 1 = on, 0 = off
				# Use this if you don't want invisible connections between pieces. Recommended for glass pieces of a shattering window for instance.

extrude = 0.

iteration = 10	# 10		=> Iteration depth of main spring algo (must be a multiple of 3 should be a multiple of 3 when using vmblur)

maxfA = 999		# 999		=> Absolute maximum speed per vertex
springDamp = .5	# .5		=> Damping of all springs, 1 = no damping, <1 = more damping
				# Use this to prevent noise patterns or chaotic flipping bumps in the surface due to unstable high stresses, but try to keep it as high as possible to save iteration rate.
				
stiff = 0		# !Obsolete, use elast instead! => Create rigidity edges (very slow on first frame), >0 = count of edges, -1 = all possible edges

relDamp = 1.	# .8		=> Vertex speed damping, relative to objects or pieces origin, 1 = no damping, <1 = more damping
				# Use this to reduce wrapping around effect on collisions, the mesh appears more rigid and less bendable on lower values.
				# But you can see this also as rotation damping, any centrifugal force is going to slow down when relDamp <1.
			
selfcoll = 0	# 0			=> Self collision (slow), 1 = on, 0 = off

ground = 0		#			=> Ground collision, 1 = on, 0 = off
grndFrict = 1.5	# 1.5		=> Ground friction (speed = speed /grndFrict)
grndBounce = 0.	# 0			=> Ground rebounce factor for single vertices, 0 = disabled (vertices stopping immediately), >0 = rebounce factor

rigidBounce = 0.	# 0			=> Rigid objects rebounce factor for single vertices, 0 = disabled (vertices stopping immediately), >0 = rebounce factor

grav = Vector(0, 0, 0)		# -0.01		=> Gravity vector (BU per frame)
wind = Vector(0, 0, 0)		# 			=> Wind vector (BU per frame)
airRes = 1.				# 1.001		=> Air resistance (speed = speed /airRes)

shock = 0						#			=> Shock wave, 0 = off, >0 factor of strength
shockOrig = Vector(0, 0, 0)		# 			=> Shock wave origin, the wave begins here at starting frame
shockSpd = Vector(0, 0, 0)		# 			=> Shock wave speed
shockChaos = 0.					# 			=> Shock chaotic forces strength

gr1rigid = 0	#			=> Make group 1 undeformable, 0 = deformable, 1 = completely rigid (indestructible)
gr2rigid = 0	#			=> Make group 2 undeformable, 0 = deformable, 1 = completely rigid (indestructible)
fix = 0			# 			=> Make collided vertices stick together (use it with gr#rigid to have better deforming results)

collLimit = 999	# 999		=> Number of collisions per vertex before falling to sleep

animated = 0	# 999		=> Animated mode until (and including) this frame number, 0 = off
				# The system is trying to hold the shape, position and rotation of the original objects until this frame number.
selection = 0	# 0			=> Selected verts in original mesh will be set to free from goal if animated or set to be fixed if free, 1 = on, 0 = off
				# This is a nice option to give a collapsing building the initial impulse by setting some verts to loose (= weak).

#---------------# Vars only active when animated > 0

tolerance = .5	# .5			=> Maximum allowed distance between original vertex and deformed vertex for animated mode (in BUs)
				# When distance gets higher than this value then vertex is detached from goal and free.
goalForce =	.2	# .5		=> Factor for springs pulling verts to original verts (goal), 0 = no force, 1 = verts will be set immediately to original position after deformation

#---------------# Vars inactive when animated > 0

elast = 0.					# 5			=> Elasticity, objects rebounces, >0 = multiplicator of rebounce energy
elastVec = Vector(1,1,1)	# (1,1,1)	=> Rebounce direction, (1,1,1) = backward, (0,0,1) = only z-axis is taken into account
							# This is useful to force the colliding object into a special rebounce direction or axis.
elastPieces = 0				# 0			=> Elast speed correction of free pieces (very slow), 0 = disabled, 1 = enabled
							# Use this only when elast >0. This prevent loose parts within the same object also to rebounce when another part is colliding.
gr1spd = Vector(0, 0, 0) 	#			=> Initial speed vector for crash group 1
gr2spd = Vector(0, 0, 0)	#			=> Initial speed vector for crash group 2

gr1rot = 0.					#			=> Initial rotation impulse for crash group 1 (degree per frame)
gr1rotAxis = 'x'			#			=> Initial rotation axis for crash group 1, 'X' 'Y' or 'Z'
gr2rot = 0.					#			=> Initial rotation impulse for crash group 2 (degree per frame)
gr2rotAxis = 'x'			#			=> Initial rotation axis for crash group 2, 'X' 'Y' or 'Z'

###########################  End of constants part  ############################
################################################################################

class PropertyType(object):
	__slots__ = ('__weakref__','group','deprecated')

	def __init__(self,group='General',deprecated=False):
		self.group      = str(group)
		self.deprecated = bool(deprecated)

	def cast(self,value):
		raise NotImplementedError("The cast method has to be overloaded by a subclass.")

class StrProperty(PropertyType):
	__slots__ = ('length')

	def __init__(self,length=MAX_STRLEN,group='General',deprecated=False):
		PropertyType.__init__(self,group,deprecated)
		self.length = int(length)

	def cast(self,value):
		if type(value) is not str:
			value = str(value)

		if len(value) > self.length:
			raise ValueError("string is to long (max length = %d): %r" % (self.length, value))

		return value

class GroupProperty(StrProperty):
	__slots__ = ()

class BoolProperty(PropertyType):
	__slots__ = ()

	def cast(self,value):
		return bool(value)

class VectorProperty(PropertyType):
	__slots__ = ('min','max')

	def __init__(self,min=MIN_VEC,max=MAX_VEC,group='General',deprecated=False):
		PropertyType.__init__(self,group,deprecated)
		self.min = min
		self.max = max
	
	def cast(self,value):
		tp = type(value)

		if tp is not Vector:
			value = Vector(value)

		if value.x < self.min[0] or value.x > self.max[0] or \
		   value.y < self.min[1] or value.y > self.max[1] or \
		   value.z < self.min[2] or value.z > self.max[2]:
			raise ValueError(
				"vector %r not in range [%r, %r]." %
				(value, self.min, self.max))

		return value

class IntProperty(PropertyType):
	__slots__ = ('min','max')

	def __init__(self,min=MIN_INT,max=MAX_INT,group='General',deprecated=False):
		PropertyType.__init__(self,group,deprecated)
		self.min = int(min)
		self.max = int(max)
		
	def cast(self,value):
		tp = type(value)

		if tp is float:
			value = int(value)
		elif tp not in (int, long):
			raise TypeError("illegal type %s." % tp.__name__)

		if value < self.min or value > self.max:
			raise ValueError("value %r not in range [%r, %r]." % (value, self.min, self.max))

		return value

class FloatProperty(PropertyType):
	__slots__ = ('min','max')

	def __init__(self,min=-1.0,max=1.0,group='General',deprecated=False):
		PropertyType.__init__(self,group,deprecated)
		self.min = float(min)
		self.max = float(max)
		
	def cast(self,value):
		tp = type(value)

		if tp in (int, long):
			value = float(value)
		elif tp is not float:
			raise TypeError("illegal type %s." % tp.__name__)

		if value < self.min or value > self.max:
			raise ValueError("value %r not in range [%r, %r]." % (value, self.min, self.max))

		return value

class SetProperty(PropertyType):
	__slots__ = ('values')

	def __init__(self,values,group='General',deprecated=False):
		PropertyType.__init__(self,group,deprecated)
		self.values = values
	
	def cast(self,value):
		if value not in self.values:
			raise ValueError("illegal value %r." % value)

		return value

propertiesDef = {
	'gr1name':     GroupProperty(group='Names'),
	'gr2name':     GroupProperty(group='Names'),
	'suffix':      StrProperty(group='Names'),
	'vmblur':      BoolProperty(),
	'sleep':       IntProperty(1,999),
	'breaks':      BoolProperty(),
	'limit':       FloatProperty(0.0,10.0),
	'limitTrans':  FloatProperty(0.0,10.0),
	'delEdges':    BoolProperty(),
	'iteration':   IntProperty(1,150),
#	'extrude':     FloatProperty(-10.,10.),
	'maxfA':       FloatProperty(0.,999.),
	'springDamp':  FloatProperty(0.0,1.0),
	'stiff':       IntProperty(-1,10000,deprecated=True),
	'relDamp':     FloatProperty(0.0,1.0),
	'selfcoll':    BoolProperty(),
	'ground':      BoolProperty(group='Ground'),
	'grndFrict':   FloatProperty(1.,5.,group='Ground'),
	'grndBounce':  FloatProperty(0.,2.,group='Ground'),
	'rigidBounce': FloatProperty(0.,2.,group='Rigid'),
	'grav':        VectorProperty((-10.,-10.,-10.),(10.,10.,10.),group='Environment'),
	'wind':        VectorProperty((-10.,-10.,-10.),(10.,10.,10.),group='Environment'),
	'airRes':      FloatProperty(1.,5.,group='Environment'),
	'shock':       FloatProperty(0.,2.,group='Shock'),
	'shockOrig':   VectorProperty(group='Shock'),
	'shockSpd':    VectorProperty((-10.,-10.,-10.),(10.,10.,10.),group='Shock'),
	'shockChaos':  FloatProperty(0.,10.,group='Shock'),
	'gr1rigid':    BoolProperty(group='Rigid'),
	'gr2rigid':    BoolProperty(group='Rigid'),
	'fix':         BoolProperty(),
	'collLimit':   IntProperty(0,1000),
	'gr1spd':      VectorProperty((-10.,-10.,-10.),(10.,10.,10.),group='Forces'),
	'gr2spd':      VectorProperty((-10.,-10.,-10.),(10.,10.,10.),group='Forces'),
	'gr1rot':      FloatProperty(-180.,180.,group='Forces'),
	'gr1rotAxis':  SetProperty({'x':'X Axis','y':'Y Axis','z':'Z Axis'},group='Forces'),
	'gr2rot':      FloatProperty(-180.,180.,group='Forces'),
	'gr2rotAxis':  SetProperty({'x':'X Axis','y':'Y Axis','z':'Z Axis'},group='Forces'),
	'animated':    IntProperty(0,999,group='Animated'),
	'selection':   BoolProperty(),
	'tolerance':   FloatProperty(0.,10.,group='Animated'),
	'goalForce':   FloatProperty(0.,1.,group='Animated'),
	'elast':       FloatProperty(0.,50.,group='Elast'),
	'elastVec':    VectorProperty((0.,0.,0.),(1.,1.,1.),group='Elast'),
	'elastPieces': BoolProperty(group='Elast')
}

##################
##################
def corrections():
	
	global vmblur, limit, limitTrans, iteration, maxfA, springDamp, elast, \
	       relDamp, gr1spd, gr2spd, gr1rot, rigidBounce, elastVec, sleep, \
	       animated, elastPieces, gr2rot, \
	       corLimit, corLimitTrans, corIteration, corMaxfA, corSpringDamp, \
		   corRelDamp, corGr1spd, corGr2spd, corGr1rot, corRigidBounce, \
	       corElastVec, corElastPieces, corGr2rot

	# set properties, if changed:
	if vmblur:
		corLimit      = limit / 3.
		corLimitTrans = limitTrans / 3.
		corIteration  = int(round(iteration /3.))
		corMaxfA      = maxfA / 3
		corSpringDamp = (springDamp + 1) *.5
		corElast      = elast / 3.
		corRelDamp    = (relDamp + 1) *.5
		corGr1spd     = gr1spd / 3
		corGr2spd     = gr2spd / 3
		corGr1rot     = gr1rot / 3
		corGr2rot     = gr2rot / 3
	else:
		corLimit      = limit
		corLimitTrans = limitTrans
		corIteration  = iteration
		corMaxfA      = maxfA
		corSpringDamp = springDamp
		corElast      = elast
		corRelDamp    = relDamp
		corGr1spd     = gr1spd
		corGr2spd     = gr2spd
		corGr1rot     = gr1rot
		corGr2rot     = gr2rot
	
	corRigidBounce = rigidBounce * elastVec
	corElastVec    = elastVec    * corElast
	
	if sleep < 1:
		sleep = 1
	
	if animated:
		corElast = 0
		corElastPieces = False
	else:
		corElastPieces = elastPieces

##############################
# correct the preset values: #
corrections()

###############################
###############################
def saveIDProperties(obj=None):
	
	if obj is None:
		obj = Scene.GetCurrent()
	
	properties = {}
	
	for name, value in getProperties().iteritems():
		if type(value) is Vector:
			properties[name] = list(value)
		else:
			properties[name] = value
	
	obj.properties['Demolition'] = properties

####################################################
####################################################
def loadIDProperties(obj=None,errorOnMissing=False):
	
	if obj is None:
		obj = Scene.GetCurrent()
	
	try:
		properties = obj.properties['Demolition']

	except KeyError:
		if errorOnMissing:
			raise
		else:
			return

	setProperties(save=False,**dict((key,properties[key]) for key in properties if key in propertiesDef))

##########################################
##########################################
def setProperties(save=True,**properties):
	
	for name, value in properties.iteritems():
		__setProperty(name, value)
	
	corrections()

	if save:
		saveIDProperties()

####################
####################
def getProperties():
		
	return dict((name, globals()[name]) for name in propertiesDef)

######################
######################
def getProperty(name):
	
	if name not in propertiesDef:
		raise NameError('Unknown property: ' + name)
	
	return globals()[name]

########################################
########################################
def setProperty(name, value, save=True):
	__setProperty(name, value)
	corrections()

	if save:
		saveIDProperties()

###############################
###############################
def __setProperty(name, value):
	
	propertyType = propertiesDef.get(name,None)
	
	if propertyType is None:
		raise NameError('Unknown property: ' + name)
	
	try:
		value = propertyType.cast(value)
	
	except (ValueError, TypeError, KeyError), e:
		raise TypeError(
			'illegal value %r for property %r.\noriginal error was: %s' %
			(value, name, str(e)))

	globals()[name] = value

########################################
# load privioulsy saved ID properties: #
loadIDProperties()

####################################
####################################
def createStiffnessEdges(objO, obj):
	
	if stiff < 0:
		edgesO      = objO.getData(0,1)
		edges       = obj.getData(0,1)
		edgesOlen   = len(edgesO)
		extend      = edges.extend
		edgesOrange = xrange(edgesOlen)

		print "Creating %d stiffness edges (complete)..." %len(edges)
		for i in edgesOrange:
			if not i%10: print '\r%d' %i
			iedgeO = edgesO[i]
			iedge  = edges[i]
			ik0, ik1 = iedgeO.key[0:2]
			for j in edgesOrange:
				jk0, jk1 = edgesO[j].key[0:2]
				jedge    = edges[j]
				if ik0 == jk0:
					extend(iedge.v2, jedge.v2)    #   add a single edge
				if ik0 == jk1:
					extend(iedge.v2, jedge.v1)    #   add a single edge
			for j in edgesOrange:
				jk0, jk1 = edgesO[j].key[0:2]
				jedge    = edges[j]
				if ik1 == jk0:
					extend(iedge.v1, jedge.v2)    #   add a single edge
				if ik1 == jk1:
					extend(iedge.v1, jedge.v1)    #   add a single edge
		print '\r%d' %(i +1)
	
	if stiff > 0:
		edgesO      = objO.getData(0,1).edges
		edges       = obj.getData(0,1).edges
		edgesOlen   = len(edgesO)
		extend      = edges.extend
		edgesOrange = xrange(edgesOlen)

		print "Creating %d stiffness edges..." %stiff
		for c in xrange(stiff):
			if not c%10: print '\r%d' %c
			i = int(Rand(0, edgesOlen))
			iedgeO = edgesO[i]
			iedge  = edges[i]
			ik0, ik1 = iedgeO.key[0:2]
			for j in edgesOrange:
				jk0, jk1 = edgesO[j].key[0:2]
				jedge    = edges[j]
				if ik0 == jk0:
					extend(iedge.v2, jedge.v2)    #   add a single edge
				if ik0 == jk1:
					extend(iedge.v2, jedge.v1)    #   add a single edge
			for j in edgesOrange:
				jk0, jk1 = edgesO[j].key[0:2]
				jedge    = edges[j]
				if ik1 == jk0:
					extend(iedge.v1, jedge.v2)    #   add a single edge
				if ik1 == jk1:
					extend(iedge.v1, jedge.v1)    #   add a single edge
		print '\r%d' %(c +1)

#		for i in xrange(stiff):
#			j = len(me.verts)
#			me.edges.extend(me.verts[int(Rand(0, j))], me.verts[int(Rand(0, j))])    #   add a single

# debug code
#	Mesh.Mode(Mesh.SelectModes.EDGE)
#	for i in xrange(len(meO.edges), len(me.edges)):
#		me.edges[i].sel = 0

####################################
####################################
def copyGroupsAndObjects(grpO, grp):

	for obj in grpO.objects:	#	unlink all none-meshes from group
		if obj.type != 'Mesh':
			grp.objects.unlink(obj)

	scn = Scene.GetCurrent()

	objlist = []
	for objO in grpO.objects:
		scn.objects.selected = []	#	unselect all
		objO.sel = 1
		name = obj.name
		Object.Duplicate(1, 1, 1, 1, 1, 1, 1, 1, 1, 1) 	#	mesh, surface, curve, text, metaball, armature, lamp, material, texture, ipo
		objO.drawType = Object.DrawTypes.WIRE
		objL = Object.GetSelected()
		obj  = objL[0]
		#	if len(objL) > 0: 
		#	else: print "Error: Object not duplicated."
		grpO.objects.unlink(obj)
		objlist.append(obj)
		me = obj.getData(0,1)
		obj.name = name +suffix
		me.name  = name +suffix
		me.update()
		ipo = obj.ipo
		if ipo:
			if ipo[OB_ROTX]: ipo[OB_ROTX] = None
			if ipo[OB_ROTY]: ipo[OB_ROTY] = None
			if ipo[OB_ROTZ]: ipo[OB_ROTZ] = None
			if ipo[OB_LOCX]: ipo[OB_LOCX] = None
			if ipo[OB_LOCY]: ipo[OB_LOCY] = None
			if ipo[OB_LOCZ]: ipo[OB_LOCZ] = None
			if ipo[OB_DLOCX]: ipo[OB_DLOCX] = None
			if ipo[OB_DLOCY]: ipo[OB_DLOCY] = None
			if ipo[OB_DLOCZ]: ipo[OB_DLOCZ] = None
		remove    = obj.modifiers.remove
		modifiers = obj.modifiers
		while len(modifiers) > 0:
			remove(modifiers[0])
		Mesh.Mode(Mesh.SelectModes.FACE)
		for f in me.faces:
			f.sel = 1
#		me.quadToTriangle(0)
		if not animated:
			co = Vector(objO.LocX +objO.dLocX, objO.LocY +objO.dLocY, objO.LocZ +objO.dLocZ)
			for vert in me.verts:
				vert.co += co
		obj.LocX = 0
		obj.LocY = 0
		obj.LocZ = 0
		if stiff:
			createStiffnessEdges(objO, obj)
	grp.objects = objlist
	print "%d original meshes copied from crash group \"%s\"." %(len(grpO.objects), grp.name)

			
######################################
######################################
def reCopyGroupsAndObjects(grpO, grp):

	for objO, obj in izip(grpO.objects, grp.objects):
		me = obj.getData(0,1)
		me.getFromObject(objO.name)
		Mesh.Mode(Mesh.SelectModes.FACE)
		for f in me.faces:
			f.sel = 1
		me.quadToTriangle(0)
		if not animated:
			co = Vector(objO.LocX +objO.dLocX, objO.LocY +objO.dLocY, objO.LocZ +objO.dLocZ)

			for vert in me.verts:
				vert.co += co

		if stiff:
			createStiffnessEdges(objO, obj)
		print "%d meshes from crash group \"%s\" refreshed." %(len(grpO.objects), grp.name)

####################################
####################################
def createBufferArrays(grp, grpnum):
	
	assert grpnum in (1, 2)
	
	SpdBufGr, AniSpdBufGr, EdgBufGr, CollGr, CollLimGr, IdxBufGr = \
	[], [], [], [], [], []

	# set Blender attributes to the values of the locals ending with 'Gr':
	for name, buf in locals().iteritems():
		if name.endswith('Gr'):
			setattr(Blender, '%s%d' % (name, grpnum), buf)
	
	grSpd = globals()['gr%dspd' % grpnum]
	for obj in grp.objects:
		me    = obj.getData(0,1)
		verts = me.verts
		for i in xrange(len(verts)):
			SpdBufGr.append(Vector(grSpd))
			AniSpdBufGr.append(Vector(0, 0, 0))
			CollLimGr.append(0)
			if animated:
				CollGr.append(2)
				IdxBufGr.append(i)
			else:
				CollGr.append(0)
		for edge in me.edges:
			EdgBufGr.append(verts[edge.key[0]].co -verts[edge.key[1]].co)
	print 'Buffers for crash group "%s" created.' %grp.name

############################################################################################
############################################################################################
def modifySelectedVerts(grpO, grp, BlenderSpdBufGr, BlenderCollGr, BlenderCollLimGr, grSpd):

	Mesh.Mode(Mesh.SelectModes.VERTEX)
	k = 0
	for objO in grpO.objects:
		try:    meO = Mesh.Get("$SwapMesh$")
		except:	meO = Mesh.New("$SwapMesh$")
		meO.getFromObject(objO.name)
		if animated:
			for vert in meO.verts:
				if vert.sel:
					BlenderCollGr[k] = 0
				else:
					BlenderSpdBufGr[k] = Vector(0,0,0)
				k += 1
		else:
			for vert in meO.verts:
				if vert.sel:
					BlenderCollGr[k] = 1
					BlenderCollLimGr[k] = collLimit +1
					BlenderSpdBufGr[k] = Vector(grSpd)
				k += 1
	
############################################################
############################################################
def rotationImpulse(grp, BlenderSpdBufGr, grRot, grRotAxis):

	# make local for speed:
	_Vector            = Vector
	_RotationMatrix    = RotationMatrix
	_TranslationMatrix = TranslationMatrix
	
	sleep_1  = sleep + 1
	curframe = Get('curframe')
	if curframe <= sleep_1:
		m = 0
		for obj in grp.objects:
		 	me = obj.getData(0,1)

# center of boundary box
			verts     = me.verts
			verts_len = len(verts)
			
			vec1 = verts[0].co
			vec2 = verts[0].co
			for i in xrange(1,verts_len):
				vert = verts[i]
				if vert.co < vec1: vec1 = vert.co
				if vert.co > vec2: vec2 = vert.co
			origin = (vec1 +vec2) /2

# other method, average value of all vertices
#			vec = Vector(0, 0, 0)
#			for vert in vert:
#				vec += vert.co
#			origin = vec /verts_len

			for i in xrange(verts_len):
				vert = verts[i]
				vec  = origin -vert.co
				mat  = _RotationMatrix(grRot, 4, grRotAxis)
				matT = _TranslationMatrix(vec)
				matTotal = mat * matT
				matTotal.invert()
				if curframe == sleep_1:
					BlenderSpdBufGr[m +i] -= vert.co - (origin +_Vector(matTotal[3][:3]))
				else:
					vert.co = (origin +_Vector(matTotal[3][:3]))
			m += 1

###########################################################################################
###########################################################################################
def animateFromSource(grpO, grp, BlenderAniSpdBufGr, BlenderCollGr, BlenderIdxBufGr, time):
	
	# make local for speed:
	_Vector            = Vector
	_TranslationMatrix = TranslationMatrix
	_RotationMatrix    = RotationMatrix
	_goalForce         = goalForce
	inv_goalForce      = 1 -_goalForce
	prevtime           = time - 1
	
	m = 0
	for objO, obj in izip(grpO.objects,grp.objects):
		me = obj.getData(0,1)
		try:    meO = Mesh.Get("$SwapMesh$")
		except:	meO = Mesh.New("$SwapMesh$")
		meO.getFromObject(objO.name)
		origin = _Vector(objO.loc) +_Vector(objO.dloc)
		objO_ipo = objO.ipo
		if objO_ipo:
			icuRotX = objO_ipo[OB_ROTX]
			icuRotY = objO_ipo[OB_ROTY]
			icuRotZ = objO_ipo[OB_ROTZ]
			if icuRotX:
				rotX = -icuRotX[prevtime]
				rotX *= 10
			if icuRotY:
				rotY = -icuRotY[prevtime]
				rotY *= 10
			if icuRotZ:
				rotZ = -icuRotZ[prevtime]
				rotZ *= 10
			icuLocX = objO_ipo[OB_LOCX]
			icuLocY = objO_ipo[OB_LOCY]
			icuLocZ = objO_ipo[OB_LOCZ]
			if icuLocX:
				locXl = icuLocX[prevtime]
				origin.x = float(locXl)
			if icuLocY:
				locYl = icuLocY[prevtime]
				origin.y = float(locYl)
			if icuLocZ:
				locZl = icuLocZ[prevtime]
				origin.z = float(locZl)
			icuDLocX = objO_ipo[OB_DLOCX]
			icuDLocY = objO_ipo[OB_DLOCY]
			icuDLocZ = objO_ipo[OB_DLOCZ]
			if icuDLocX:
				dLocXl = icuDLocX[prevtime]
				origin.x += dLocXl
			if icuDLocY:
				dLocYl = icuDLocY[prevtime]
				origin.y += dLocYl
			if icuDLocZ:
				dLocZl = icuDLocZ[prevtime]
				origin.z += dLocZl
		buf = []
		buf_append = buf.append
		verts = me.verts
		for i in xrange(len(verts)):
			if BlenderCollGr[m +i] >= 2:
				co = verts[i].co
				buf_append(_Vector(co))
				co += origin
				if objO_ipo:
					if icuRotX:
						matT = _TranslationMatrix(origin -co)
						mat = _RotationMatrix(rotX, 4, 'X')
						matTotal = mat * matT
						matTotal.invert()
						co[:] = origin +_Vector(matTotal[3][:3])
					if icuRotY:
						matT = _TranslationMatrix(origin -co)
						mat = _RotationMatrix(rotY, 4, 'Y')
						matTotal = mat * matT
						matTotal.invert()
						co[:] = origin +_Vector(matTotal[3][:3])
					if icuRotZ:
						matT = _TranslationMatrix(origin -co)
						mat = _RotationMatrix(rotZ, 4, 'Z')
						matTotal = mat * matT
						matTotal.invert()
						co[:] = origin +_Vector(matTotal[3][:3])
		bufL = []
		bufL_append = bufL.append
		j = 0
		for i in xrange(len(verts)):
			if BlenderCollGr[m +i] >= 2:
				co = verts[i].co
				bufL_append(_Vector(co))
				co[:] = buf[j]
				j += 1

		origin = _Vector(objO.loc) +_Vector(objO.dloc)
		if objO_ipo:
			icuRotX = objO_ipo[OB_ROTX]
			icuRotY = objO_ipo[OB_ROTY]
			icuRotZ = objO_ipo[OB_ROTZ]
			if icuRotX:
				rotX = -icuRotX[time]
				rotX *= 10
			if icuRotY:
				rotY = -icuRotY[time]
				rotY *= 10
			if icuRotZ:
				rotZ = -icuRotZ[time]
				rotZ *= 10
			icuLocX = objO_ipo[OB_LOCX]
			icuLocY = objO_ipo[OB_LOCY]
			icuLocZ = objO_ipo[OB_LOCZ]
			if icuLocX:
				locXl = icuLocX[time]
				origin.x = float(locXl)
			if icuLocY:
				locYl = icuLocY[time]
				origin.y = float(locYl)
			if icuLocZ:
				locZl = icuLocZ[time]
				origin.z = float(locZl)
			icuDLocX = objO_ipo[OB_DLOCX]
			icuDLocY = objO_ipo[OB_DLOCY]
			icuDLocZ = objO_ipo[OB_DLOCZ]
			if icuDLocX:
				dLocXl = icuDLocX[time]
				origin.x += dLocXl
			if icuDLocY:
				dLocYl = icuDLocY[time]
				origin.y += dLocYl
			if icuDLocZ:
				dLocZl = icuDLocZ[time]
				origin.z += dLocZl
		j = 0
		vertsO = meO.verts
		origin_x = origin.x
		origin_y = origin.y
		origin_z = origin.z
		for i in xrange(len(verts)):
			m_i = m +i
			if BlenderCollGr[m_i] >= 2:
				co = verts[i].co

				coO = vertsO[BlenderIdxBufGr[m_i]].co
				co.x = (co.x *inv_goalForce) +(coO.x *_goalForce) +origin_x
				co.y = (co.y *inv_goalForce) +(coO.y *_goalForce) +origin_y
				co.z = (co.z *inv_goalForce) +(coO.z *_goalForce) +origin_z
				if objO_ipo:
					if icuRotX:
						matT = _TranslationMatrix(origin -co)
						mat = _RotationMatrix(rotX, 4, 'X')
						matTotal = mat * matT
						matTotal.invert()
						co[:] = origin +_Vector(matTotal[3][:3])
					if icuRotY:
						matT = _TranslationMatrix(origin -co)
						mat = _RotationMatrix(rotY, 4, 'Y')
						matTotal = mat * matT
						matTotal.invert()
						co[:] = origin +_Vector(matTotal[3][:3])
					if icuRotZ:
						matT = _TranslationMatrix(origin -co)
						mat = _RotationMatrix(rotZ, 4, 'Z')
						matTotal = mat * matT
						matTotal.invert()
						co[:] = origin +_Vector(matTotal[3][:3])
				BlenderAniSpdBufGr[m_i] = co -bufL[j]
				j += 1
		j = 0
		for i in xrange(len(verts)):
			if BlenderCollGr[m +i] >= 2:
				verts[i].co = _Vector(bufL[j])
				j += 1
		m += len(verts)

########################################################
########################################################
def reAnimateFromSource(grpO, grp, BlenderCollGr, time):
	
	# make local for speed:
	_Vector            = Vector
	_TranslationMatrix = TranslationMatrix
	_RotationMatrix    = RotationMatrix
	prevtime           = time - 1
	
	m = 0
	for objO, obj in izip(grpO.objects, grp.objects):
		me = obj.getData(0,1)
		origin = _Vector(objO.loc) +_Vector(objO.dloc)
		objO_ipo = objO.ipo
		if objO_ipo:
			icuRotX = objO_ipo[OB_ROTX]
			icuRotY = objO_ipo[OB_ROTY]
			icuRotZ = objO_ipo[OB_ROTZ]
			if icuRotX:
				rotX = -icuRotX[prevtime]
				rotX *= 10
			if icuRotY:
				rotY = -icuRotY[prevtime]
				rotY *= 10
			if icuRotZ:
				rotZ = -icuRotZ[prevtime]
				rotZ *= 10
			icuLocX = objO_ipo[OB_LOCX]
			icuLocY = objO_ipo[OB_LOCY]
			icuLocZ = objO_ipo[OB_LOCZ]
			if icuLocX:
				locXl = icuLocX[prevtime]
				origin.x = float(locXl)
			if icuLocY:
				locYl = icuLocY[prevtime]
				origin.y = float(locYl)
			if icuLocZ:
				locZl = icuLocZ[prevtime]
				origin.z = float(locZl)
			icuDLocX = objO_ipo[OB_DLOCX]
			icuDLocY = objO_ipo[OB_DLOCY]
			icuDLocZ = objO_ipo[OB_DLOCZ]
			if icuDLocX:
				dLocXl = icuDLocX[prevtime]
				origin.x += dLocXl
			if icuDLocY:
				dLocYl = icuDLocY[prevtime]
				origin.y += dLocYl
			if icuDLocZ:
				dLocZl = icuDLocZ[prevtime]
				origin.z += dLocZl
		verts = me.verts
		for i in xrange(len(verts)):
			if BlenderCollGr[m +i] >= 2:
				co = verts[i].co
				if objO_ipo:
					if icuRotZ:
						matT = _TranslationMatrix(origin -co)
						mat = _RotationMatrix(-rotZ, 4, 'Z')
						matTotal = mat * matT
						matTotal.invert()
						co[:] = origin +_Vector(matTotal[3][:3])
					if icuRotY:
						matT = _TranslationMatrix(origin -co)
						mat = _RotationMatrix(-rotY, 4, 'Y')
						matTotal = mat * matT
						matTotal.invert()
						co[:] = origin +_Vector(matTotal[3][:3])
					if icuRotX:
						matT = _TranslationMatrix(origin -co)
						mat = _RotationMatrix(-rotX, 4, 'X')
						matTotal = mat * matT
						matTotal.invert()
						co[:] = origin +_Vector(matTotal[3][:3])
				co -= origin
		m += len(verts)

########################################################################
########################################################################
def checkToleranceFromSource(grpO, grp, BlenderCollGr, BlenderIdxBufGr):
	
	m = 0
	for objO, obj in izip(grpO.objects, grp.objects):
		me = obj.getData(0,1)
		try:    meO = Mesh.Get("$SwapMesh$")
		except:	meO = Mesh.New("$SwapMesh$")
		meO.getFromObject(objO.name)
		verts  = me.verts
		vertsO = meO.verts
		for i in xrange(len(verts)):
			dif = verts[i].co -vertsO[BlenderIdxBufGr[m +i]].co
			if BlenderCollGr[m +i] >= 2:
				if dif.length > tolerance:
					BlenderCollGr[m +i] = 3
		m += len(verts)

############################################
############################################
def gravity(BlenderSpdBufGr, BlenderCollGr):
	
	# make local for speed:
	_grav = grav

	for i in xrange(len(BlenderSpdBufGr)):
		if BlenderCollGr[i] == 0:
			BlenderSpdBufGr[i] += _grav

##########################################
##########################################
def winds(BlenderSpdBufGr, BlenderCollGr):

	# make local for speed:
	_wind = wind

	for i in xrange(len(BlenderSpdBufGr)):
		if BlenderCollGr[i] == 0:
			BlenderSpdBufGr[i] += _wind

######################################################
######################################################
def shockImpulse(grp, BlenderSpdBufGr, BlenderCollGr):

	# make local for speed:
	_Vector     = Vector
	_Rand       = Rand
	_shockChaos = shockChaos
	_shock      = shock
	_shockSpd   = shockSpd

	fac = abs(Get('curframe') -Get('staframe'))
	min = Vector(_shockSpd)
	min.x = (min.x *fac) +shockOrig.x -(_shockSpd.x /2)
	min.y = (min.y *fac) +shockOrig.y -(_shockSpd.y /2)
	min.z = (min.z *fac) +shockOrig.z -(_shockSpd.z /2)
	max = min +_shockSpd
	shockSpd_shock = _shockSpd *_shock
	j = 0
	for obj in grp.objects:
	 	verts = obj.getData(0,1).verts
		for i in xrange(len(verts)):
			if BlenderCollGr[j] != 1:
				vec = verts[i].co
				if ((vec[0] > min[0] and vec[0] < max[0]) or (vec[0] < min[0] and vec[0] > max[0])) \
				or ((vec[1] > min[1] and vec[1] < max[1]) or (vec[1] < min[1] and vec[1] > max[1])) \
				or ((vec[2] > min[2] and vec[2] < max[2]) or (vec[2] < min[2] and vec[2] > max[2])):
					BlenderSpdBufGr[j] += shockSpd_shock + \
						_Vector(
							_Rand(-_shockChaos, _shockChaos),
							_Rand(-_shockChaos, _shockChaos),
							_Rand(-_shockChaos, _shockChaos))
			j += 1

##################################################
##################################################
def airResistanceAndMaximumSpeed(BlenderSpdBufGr):
	
	# make local for speed:
	_maxfA = corMaxfA

	for spd in BlenderSpdBufGr:
		spd /= airRes
		if spd.length > _maxfA:
			spd *= _maxfA /spd.length

###########################################################################
###########################################################################
def countingPiecesForElasticitySpeedCorrection(grp, bufPiece, bufPieceOfs):
	
	pieces = 1
	pieceOfs = []
	for obj in grp.objects:
	 	me = obj.getData(0,1)
		bufKey0 = [edge.key[0] for edge in me.edges]
		bufKey1 = [edge.key[1] for edge in me.edges]
		subIdx = 1
		while subIdx < len(me.edges) -1:
			q = True
			edgeIdx1 = 0
			while edgeIdx1 < subIdx: # max = length -2
				key0idx1b = bufKey0[edgeIdx1]
				key1idx1b = bufKey1[edgeIdx1]
				edgeIdx2 = subIdx +1
				while edgeIdx2 < len(me.edges):
					# check if edges connected and swap them (sorting)
					if key0idx1b == bufKey1[edgeIdx2] or key0idx1b == bufKey0[edgeIdx2] or \
					   key1idx1b == bufKey0[edgeIdx2] or key1idx1b == bufKey1[edgeIdx2]:
						buf = bufKey0[subIdx]
						bufKey0[subIdx] = bufKey0[edgeIdx2]
						bufKey0[edgeIdx2] = buf
						buf = bufKey1[subIdx]
						bufKey1[subIdx] = bufKey1[edgeIdx2]
						bufKey1[edgeIdx2] = buf
						edgeIdx1 = subIdx
						q = False
						break
					edgeIdx2 += 1
				edgeIdx1 += 1
			if q:	# one loose piece complete
				pieces += 1
				pieceOfs.append(subIdx)
			subIdx += 1
		pieceOfs.append(subIdx)
		
		j = 0
		for i in xrange(0, len(me.edges)):
			if bufKey0[i] not in bufPiece: bufPiece.append(bufKey0[i])
			if bufKey1[i] not in bufPiece: bufPiece.append(bufKey1[i])
			if i == pieceOfs[j] -1:
				bufPieceOfs.append(len(bufPiece))
				if j < len(pieceOfs) -1: j += 1
# debug code
#		Mesh.Mode(Mesh.SelectModes.VERTEX)
#		for i in xrange(bufPieceOfs[1], bufPieceOfs[2]):
#			me.verts[bufPiece[i]].sel = 0
	print "Crash group \"%s\" contains %d loose piece(s)." %(grp.name, pieces)

######################################################################
######################################################################
def measureSpeed(grp, BlenderSpdBufGr, bufPiece, bufPieceOfs, vecOld):

	# make local for speed improvement:
	_elastPieces = corElastPieces
	_Vector      = Vector

	m             = 0
	j_range       = xrange(len(bufPieceOfs) -1)
	vecOld_append = vecOld.append

	for obj in grp.objects:
	 	me = obj.getData(0,1)
		verts_len = len(me.verts)
		
		if _elastPieces:
			for j in j_range:
				vO = vecOld[j]
				vecold_append(_Vector(0, 0, 0))
				for i in xrange(bufPieceOfs[j], bufPieceOfs[j +1]):
					vO += BlenderSpdBufGr[m +bufPiece[i]]
				vO /= (bufPieceOfs[j +1] -bufPieceOfs[j])
		else:
			vecOld_append(_Vector(0, 0, 0))
			vecOld0 = vecOld[0]
			for i in xrange(verts_len):
				vecOld0 += BlenderSpdBufGr[m +i]
			vecOld0 /= verts_len
		
		m += verts_len
	
##############################################################################
##############################################################################
def springBodyAlgorithm(grp, BlenderSpdBufGr, BlenderEdgBufGr, BlenderCollGr):
	
	# make often used globals local for speed:
	_springDamp   = corSpringDamp
	_springDamp_2 = _springDamp * 2

	m = 0
	k = 0
	for obj in grp.objects:
	 	me  = obj.getData(0,1)
		buf = [Vector(vert.co) for vert in me.verts]
		
		edges       = me.edges
		edges_len   = len(edges)
		edges_range = xrange(edges_len)
		verts       = me.verts
		verts_len   = len(verts)

		for it in xrange(corIteration):
			for j in edges_range:
				edge    = edges[j]
				edge_v1 = edge.v1
				edge_v2 = edge.v2
				m_edge_v1_idx = m +edge_v1.index
				m_edge_v2_idx = m +edge_v2.index
				dif = (edge_v1.co +BlenderSpdBufGr[m_edge_v1_idx]) \
				     -(edge_v2.co +BlenderSpdBufGr[m_edge_v2_idx])
				if dif.length > 0:
					fac = (BlenderEdgBufGr[k +j].length /dif.length) -1
				else:
					fac = 0
				cv1 = BlenderCollGr[m_edge_v1_idx]
				cv2 = BlenderCollGr[m_edge_v2_idx]
				if cv1 != 1 and cv2 != 1:
					edge_v1.co += dif *fac *_springDamp
					edge_v2.co -= dif *fac *_springDamp
				elif cv1 != 1 and cv2 == 1:
					edge_v1.co += dif *fac *_springDamp_2
				elif cv1 == 1 and cv2 != 1:
					edge_v2.co -= dif *fac *_springDamp_2
		
		for i in xrange(verts_len):
			if BlenderCollGr[m +i] != 1:
				BlenderSpdBufGr[m +i] -= buf[i] -verts[i].co
			verts[i].co = buf[i]

		m += verts_len
		k += edges_len

#####################################################################################
#####################################################################################
def correctSpeed(grp, BlenderSpdBufGr, bufPiece, bufPieceOfs, vecOld, BlenderCollGr):

	# make local for speed:
	_elastPieces = corElastPieces
	_Vector      = Vector
	_relDamp     = corRelDamp
	_elastVec    = corElastVec

	m = 0
	for obj in grp.objects:
	 	me        = obj.getData(0,1)
		verts     = me.verts
		verts_len = len(verts)
		
		if _elastPieces:
			for j in xrange(len(bufPieceOfs) -1):
				vecNew = _Vector(0, 0, 0)
				for i in xrange(bufPieceOfs[j], bufPieceOfs[j +1]):
					vecNew += BlenderSpdBufGr[m +bufPiece[i]]
				vecNew /= bufPieceOfs[j +1] -bufPieceOfs[j]
				dif = vecOld[j] -vecNew
				dif[0] *= _elastVec[0]
				dif[1] *= _elastVec[1]
				dif[2] *= _elastVec[2]
				for i in xrange(bufPieceOfs[j], bufPieceOfs[j +1]):
					m_bufP = m +bufPiece[i]
					BlenderSpdBufGr[m_bufP] -= dif
					difV = vecNew -BlenderSpdBufGr[m_bufP]
					if BlenderCollGr[m +i] != 1: BlenderSpdBufGr[m_bufP] = vecNew -(difV *_relDamp)
		else:
			m_verts_range = xrange(m, m + verts_len)

			vecNew = _Vector(0, 0, 0)
			for i in m_verts_range:
				vecNew += BlenderSpdBufGr[i]

			vecNew /= verts_len
			dif = vecOld[0] -vecNew
			dif[0] *= _elastVec[0]
			dif[1] *= _elastVec[1]
			dif[2] *= _elastVec[2]

			for i in m_verts_range:
				BlenderSpdBufGr[i] -= dif
				difV = vecNew -BlenderSpdBufGr[i]
				if BlenderCollGr[i] != 1: BlenderSpdBufGr[i] = vecNew -(difV *_relDamp)
						
		m += verts_len

##############################################################
##############################################################
def groundCollisionCheck(grp, BlenderSpdBufGr, BlenderCollGr):
	
	# make local for speed:
	_grndBounce = grndBounce
	_grndFrict  = grndFrict

	j = 0
	for obj in grp.objects:
	 	me = obj.getData(0,1)

		for vert in me.verts:
			if vert.co.z < 0:
		 		vert.co.z += abs(BlenderSpdBufGr[j].z)
				BlenderSpdBufGr[j].x /= _grndFrict
				BlenderSpdBufGr[j].y /= _grndFrict
				if _grndBounce:
					BlenderSpdBufGr[j].z = abs(BlenderSpdBufGr[j].z) *_grndBounce
				else:
					BlenderSpdBufGr[j].z = 0
				if BlenderCollGr[j] == 0: BlenderCollGr[j] = 1
# debug code
#				Mesh.Mode(Mesh.SelectModes.VERTEX)
#				vert.sel = 0
			j += 1
					
#############################################################################################################################################################################################################
#############################################################################################################################################################################################################
def geometryCollisionCheck(grp1, grp2, BlenderSpdBufGr1, BlenderSpdBufGr2, SpdBufTmp1, SpdBufTmp2, BlenderCollGr1, BlenderCollGr2, BlenderCollLimGr1, BlenderCollLimGr2, gr1Rigid, gr1Spd, gr2Rigid, gr2Spd):
	
	# make often used globals local for speed:
	_Intersect     = Intersect
	_rigidBounce_x = corRigidBounce.x
	_rigidBounce_y = corRigidBounce.y
	_rigidBounce_z = corRigidBounce.z

	k = 0
	for obj2 in grp2.objects:
	 	me2 = obj2.getData(0,1)
		l = 0
		for obj1 in grp1.objects:
		 	me1 = obj1.getData(0,1)
			me1_verts       = me1.verts
			me1_verts_range = xrange(len(me1_verts))

			for face in me2.faces:
				fverts = face.verts
				v0idx = fverts[0].index
				v1idx = fverts[1].index
				v2idx = fverts[2].index
				v0co = fverts[0].co
				v1co = fverts[1].co
				v2co = fverts[2].co
				k_v0idx = k + v0idx
				k_v1idx = k + v1idx
				k_v2idx = k + v2idx
				for i in me1_verts_range:
					l_i   = l +i
					me1vi = me1_verts[i]
					spd1  = BlenderSpdBufGr1[l_i]
					
					if BlenderCollGr1[l_i] != 1:
						vec = _Intersect(v0co, v1co, v2co, spd1, me1vi.co, 1)
						if vec:
							min = me1vi.co
							max = me1vi.co +spd1
							if ((vec[0] > min[0] and vec[0] < max[0]) or (vec[0] < min[0] and vec[0] > max[0])) \
							or ((vec[1] > min[1] and vec[1] < max[1]) or (vec[1] < min[1] and vec[1] > max[1])) \
							or ((vec[2] > min[2] and vec[2] < max[2]) or (vec[2] < min[2] and vec[2] > max[2])):
								if gr1Rigid:
									me1vi.co = vec -spd1 /100
									vec.x = spd1.x *_rigidBounce_x
									vec.y = spd1.y *_rigidBounce_y
									vec.z = spd1.z *_rigidBounce_z
									SpdBufTmp1[l_i] = gr1Spd -vec
								else:
									if gr2Rigid:
										vec.x = spd1.x *_rigidBounce_x
										vec.y = spd1.y *_rigidBounce_y
										vec.z = spd1.z *_rigidBounce_z
										SpdBufTmp2[k_v0idx] = gr2Spd +vec
										SpdBufTmp2[k_v1idx] = gr2Spd +vec
										SpdBufTmp2[k_v2idx] = gr2Spd +vec
									else:
										vec = spd1 /2
										SpdBufTmp1[l_i] = vec + \
											((BlenderSpdBufGr2[k_v0idx] +
											  BlenderSpdBufGr2[k_v1idx] +
											  BlenderSpdBufGr2[k_v2idx]) /6)
										SpdBufTmp2[k_v0idx] = BlenderSpdBufGr2[k_v0idx] /2 +vec
										SpdBufTmp2[k_v1idx] = BlenderSpdBufGr2[k_v1idx] /2 +vec
										SpdBufTmp2[k_v2idx] = BlenderSpdBufGr2[k_v2idx] /2 +vec
									if BlenderCollGr2[k_v0idx] == 0: BlenderCollGr2[k_v0idx] = 1
									if BlenderCollGr2[k_v1idx] == 0: BlenderCollGr2[k_v1idx] = 1
									if BlenderCollGr2[k_v2idx] == 0: BlenderCollGr2[k_v2idx] = 1
									BlenderCollLimGr2[k_v0idx] += 1
									BlenderCollLimGr2[k_v1idx] += 1
									BlenderCollLimGr2[k_v2idx] += 1
# debug code
#									Mesh.Mode(Mesh.SelectModes.VERTEX)
#									me2.verts[v0idx].sel = 0
#									me2.verts[v1idx].sel = 0
#									me2.verts[v2idx].sel = 0
#								Mesh.Mode(Mesh.SelectModes.VERTEX)
#								me1_verts[i].sel = 0
								if BlenderCollGr1[l_i] == 0: BlenderCollGr1[l_i] = 1
								BlenderCollLimGr1[l_i] += 1
			l += len(me1_verts)
		k += len(me2.verts)
	
############################################
############################################
def applySpeedVectors(grp, BlenderSpdBufGr):
	
	j = 0
	for obj in grp.objects:
	 	me = obj.getData(0,1)
		for vert in me.verts:
			vert.co += BlenderSpdBufGr[j]
			j += 1

##############################################################################################################
##############################################################################################################
def edgeStressTest(grp, BlenderSpdBufGr, BlenderAniSpdBufGr, BlenderEdgBufGr, BlenderCollGr, BlenderIdxBufGr):
	
	# make often used globals local for speed:
	_delEdges   = delEdges
	_limit      = corLimit
	_limitTrans = corLimitTrans

	k = 0
	m = 0
	for obj in grp.objects:
		me    = obj.getData(0,1)
		faces = me.faces
		edges = me.edges
		edges_delete = edges.delete

		if _limitTrans:
			transpVerts = []
			transpVerts_append = transpVerts.append
			mat = me.materials
			if mat:
				for fac in faces:
					if mat[fac.mat].alpha < 1:
						fac_verts = fac.verts
						transpVerts_append(fac_verts[0].index)
						transpVerts_append(fac_verts[1].index)
						transpVerts_append(fac_verts[2].index)
# replaced, code caused memory leak on older blender versions
#			if me.materials:
#				for l in xrange(len(faces)):
#					mat = me.materials[me.faces[l].mat]
#					if mat.alpha < 1:
#						transpVerts.append(me.faces[l].verts[0].index)
#						transpVerts.append(me.faces[l].verts[1].index)
#						transpVerts.append(me.faces[l].verts[2].index)
		j     = 0
		count = 0
		while j < len(edges):
			lim       = _limit
			edge      = edges[j]
			edge_key  = edge.key
			edge_key0 = edge_key[0]
			edge_key1 = edge_key[1]
			
			if _delEdges:
				q = True
				for fac in faces:
					fac_verts = fac.verts
					if fac_verts[0].index == edge_key0 or \
					   fac_verts[1].index == edge_key0 or \
					   fac_verts[2].index == edge_key0:
						q = False
						break
				if not q:
					q = True
					for fac in me.faces:
						if fac.verts[0].index == edge_key1 or \
						   fac.verts[1].index == edge_key1 or \
						   fac.verts[2].index == edge_key1:
							q = False
							break
			else:
				q = False
			if not q and _limitTrans:
				if edge_key0 in transpVerts and edge_key1 in transpVerts:
					lim = _limitTrans
			if q or abs(edge.length - BlenderEdgBufGr[k].length) > lim:
				q1 = 0
				q2 = 0
				for n in xrange(len(edges)):
					nedge_key = edges[n].key
					if nedge_key[0] == edge_key0 or nedge_key[1] == edge_key0: q1 += 1
					if nedge_key[0] == edge_key1 or nedge_key[1] == edge_key1: q2 += 1
				if q1 == 1:
					if m +edge_key0 < len(BlenderSpdBufGr):
						del BlenderSpdBufGr[m +edge_key0]
						del BlenderCollGr[m +edge_key0]
						if animated:
							del BlenderAniSpdBufGr[m +edge_key0]
							del BlenderIdxBufGr[m +edge_key0]
				if q2 == 1:
					if m +edge_key1 < len(BlenderSpdBufGr):
						del BlenderSpdBufGr[m +edge_key1]
						del BlenderCollGr[m +edge_key1]
						if animated:
							del BlenderAniSpdBufGr[m +edge_key1]
							del BlenderIdxBufGr[m +edge_key1]
				del BlenderEdgBufGr[k]
				edges_delete(j)
				j -= 1
				edge      = edges[j]
				edge_key  = edge.key
				edge_key0 = edge_key[0]
				edge_key1 = edge_key[1]
				k -= 1
				count += 1
			j += 1
			k += 1
		me.update() # "this method may be only temporary and may be removed in future (blender) releases."
		m += len(me.verts)
		if count > 0:
			print "%d edges from object \"%s\" overstressed and removed." %(count, obj.name)
	

###########
###########
###########
def main():

	grp1o = Group.Get(gr1name)
	grp2o = Group.Get(gr2name)
	q = False
	try:
		grp1 = Group.Get(gr1name +suffix)
		grp2 = Group.Get(gr2name +suffix)
	except:
		q = True

#--------------------------------------------

	if Get('curframe') == Get('staframe') and q:	# new if no crash group exists
		print "======== NEW SIMULATION ========"
		DemolitionProfiler_clear()

###########################
##	copy groups and objects

		grp1 = Group.New(gr1name +suffix)
		copyGroupsAndObjects(grp1o, grp1)
		grp2 = Group.New(gr2name +suffix)
		copyGroupsAndObjects(grp2o, grp2)

########################
##	create buffer arrays

		createBufferArrays(grp1, 1)
		createBufferArrays(grp2, 2)

		if animated:
			animateFromSource(grp1o, grp1, Blender.AniSpdBufGr1, Blender.CollGr1, Blender.IdxBufGr1, Get('curtime'))
			animateFromSource(grp2o, grp2, Blender.AniSpdBufGr2, Blender.CollGr2, Blender.IdxBufGr2, Get('curtime'))
		if selection:
			modifySelectedVerts(grp1o, grp1, Blender.SpdBufGr1, Blender.CollGr1, Blender.CollLimGr1, corGr1spd)
			modifySelectedVerts(grp2o, grp2, Blender.SpdBufGr2, Blender.CollGr2, Blender.CollLimGr2, corGr2spd)

#--------------------------------------------

	if Get('curframe') == Get('staframe') and not q:	# new and a crash group exists
		DemolitionProfiler_printInfos()
		DemolitionProfiler_clear()
		print "======== RESTART SIMULATION ========"

#############################
##	recopy groups and objects

		reCopyGroupsAndObjects(grp1o, grp1)
		reCopyGroupsAndObjects(grp2o, grp2)
		
########################
##	create buffer arrays

		createBufferArrays(grp1, 1)
		createBufferArrays(grp2, 2)

		if animated:
			animateFromSource(grp1o, grp1, Blender.AniSpdBufGr1, Blender.CollGr1, Blender.IdxBufGr1, Get('curtime'))
			animateFromSource(grp2o, grp2, Blender.AniSpdBufGr2, Blender.CollGr2, Blender.IdxBufGr2, Get('curtime'))
		if selection:
			modifySelectedVerts(grp1o, grp1, Blender.SpdBufGr1, Blender.CollGr1, Blender.CollLimGr1, corGr1spd)
			modifySelectedVerts(grp2o, grp2, Blender.SpdBufGr2, Blender.CollGr2, Blender.CollLimGr2, corGr2spd)

#--------------------------------------------

	if Get('curframe') > Get('staframe'):	# frame processing
		print "Processing frame %d..." %Get('curframe')

##	recalc normals
		for obj in grp1.objects:
		 	me = obj.getData(0,1)
			me.calcNormals()
		for obj in grp2.objects:
		 	me = obj.getData(0,1)
			me.calcNormals()
		
##############################################
##	modifying speed vectors - rotation impulse
##	modifying speed vectors - gravity
##	modifying speed vectors - winds
##	modifying speed vectors - air resistance and maximum speed
##	backup speed vectors for later comparison

		### make local for speed: ###
		BlenderCollGr1    = Blender.CollGr1
		BlenderCollGr2    = Blender.CollGr2
		BlenderCollLimGr1 = Blender.CollLimGr1
		BlenderCollLimGr2 = Blender.CollLimGr2
		_collLimit        = collLimit

		if Get('curframe') > animated:
			if corGr1rot != 0: rotationImpulse(grp1, Blender.SpdBufGr1, corGr1rot, gr1rotAxis)
			if corGr2rot != 0: rotationImpulse(grp2, Blender.SpdBufGr2, corGr2rot, gr2rotAxis)
		else:
			reAnimateFromSource(grp1o, grp1, BlenderCollGr1, Get('curtime'))
			reAnimateFromSource(grp2o, grp2, BlenderCollGr2, Get('curtime'))
			checkToleranceFromSource(grp1o, grp1, BlenderCollGr1, Blender.IdxBufGr1)
			checkToleranceFromSource(grp2o, grp2, BlenderCollGr2, Blender.IdxBufGr2)
			animateFromSource(grp1o, grp1, Blender.AniSpdBufGr1, BlenderCollGr1, Blender.IdxBufGr1, Get('curtime'))
			animateFromSource(grp2o, grp2, Blender.AniSpdBufGr2, BlenderCollGr2, Blender.IdxBufGr2, Get('curtime'))
			for i, coll1, spd1, aniSpd1 in izip(icount(), BlenderCollGr1, Blender.SpdBufGr1, Blender.AniSpdBufGr1):
				if coll1 >= 2: spd1 += aniSpd1
				# 'coll1 = 0' would only change the local variable coll1 but not the corresponding list item:
				if coll1 == 3: BlenderCollGr1[i] = 0
			for i, coll2, spd2, aniSpd2 in izip(icount(), BlenderCollGr2, Blender.SpdBufGr2, Blender.AniSpdBufGr2):
				if coll2 >= 2: spd2 += aniSpd2
				if coll2 == 3: BlenderCollGr2[i] = 0
		if Get('curframe') == animated +1:
			for i, coll1, spd1, aniSpd1 in izip(icount(), BlenderCollGr1, Blender.SpdBufGr1, Blender.AniSpdBufGr1):
				if coll1 == 2:
					spd1 += aniSpd1
					BlenderCollGr1[i] = 0
			for i, coll2, spd2, aniSpd2 in izip(icount(), BlenderCollGr2, Blender.SpdBufGr2, Blender.AniSpdBufGr2):
				if coll2 == 2:
					spd2 += aniSpd2
					BlenderCollGr2[i] = 0
		if not gr1rigid:
			gravity(Blender.SpdBufGr1, BlenderCollGr1)
			winds(Blender.SpdBufGr1, BlenderCollGr1)
			if shock: shockImpulse(grp1, Blender.SpdBufGr1, Blender.CollGr1)
		if not gr2rigid:
			gravity(Blender.SpdBufGr2, BlenderCollGr2)
			winds(Blender.SpdBufGr2, BlenderCollGr2)
			if shock: shockImpulse(grp2, Blender.SpdBufGr2, Blender.CollGr2)

#########################
##	reset collision flags
##	modifying speed vectors - ground collision check
##	modifying speed vectors - geometry collision check
##  counting pieces for elasticity speed correction
##	modifying speed vectors - spring body algorithm
##	modifying vertices - applying speed vectors
##	edge stress test
	
		if Get('curframe') > sleep:
			if Get('curframe') > animated:
				bufPiece1 = []; bufPieceOfs1 = [0]
				if corElastPieces: countingPiecesForElasticitySpeedCorrection(grp1, bufPiece1, bufPieceOfs1)
				vecOld = []
				measureSpeed(grp1, Blender.SpdBufGr1, bufPiece1, bufPieceOfs1, vecOld)
			if ground and not gr1rigid: groundCollisionCheck(grp1, Blender.SpdBufGr1, BlenderCollGr1)
			SpdBufTmp1 = [Vector(obj) for obj in Blender.SpdBufGr1]
			SpdBufTmp2 = [Vector(obj) for obj in Blender.SpdBufGr2]
			geometryCollisionCheck(grp1, grp2, Blender.SpdBufGr1, Blender.SpdBufGr2, SpdBufTmp1, SpdBufTmp2, BlenderCollGr1, BlenderCollGr2, BlenderCollLimGr1, BlenderCollLimGr2, gr2rigid, corGr2spd, gr1rigid, corGr1spd)
			if selfcoll: geometryCollisionCheck(grp1, grp1, Blender.SpdBufGr1, Blender.SpdBufGr1, SpdBufTmp1, SpdBufTmp1, BlenderCollGr1, BlenderCollGr1, BlenderCollLimGr1, BlenderCollLimGr1, gr1rigid, corGr1spd, gr1rigid, corGr1spd)
			for coll1, spd1, spdTmp1 in izip(BlenderCollGr1, Blender.SpdBufGr1, SpdBufTmp1):
				if coll1 == 1: spd1[:] = spdTmp1
			for coll2, spd2, spdTmp2 in izip(BlenderCollGr2, Blender.SpdBufGr2, SpdBufTmp2):
				if coll2 == 1: spd2[:] = spdTmp2
			if not gr1rigid: springBodyAlgorithm(grp1, Blender.SpdBufGr1, Blender.EdgBufGr1, BlenderCollGr1)
			if not fix:
				for i in xrange(len(BlenderCollGr1)):
					if BlenderCollGr1[i] == 1 and BlenderCollLimGr1[i] <= _collLimit: BlenderCollGr1[i] = 0
				for i in xrange(len(BlenderCollGr2)):
					if BlenderCollGr2[i] == 1 and BlenderCollLimGr2[i] <= _collLimit: BlenderCollGr2[i] = 0
			if Get('curframe') > animated: correctSpeed(grp1, Blender.SpdBufGr1, bufPiece1, bufPieceOfs1, vecOld, BlenderCollGr1)
		airResistanceAndMaximumSpeed(Blender.SpdBufGr1)
		applySpeedVectors(grp1, Blender.SpdBufGr1)
		if Get('curframe') > sleep:
			if Get('curframe') > animated:
				bufPiece2 = []; bufPieceOfs2 = [0]
				if corElastPieces: countingPiecesForElasticitySpeedCorrection(grp2, bufPiece2, bufPieceOfs2)
				vecOld = []
				measureSpeed(grp2, Blender.SpdBufGr2, bufPiece2, bufPieceOfs2, vecOld)
			if ground and not gr2rigid: groundCollisionCheck(grp2, Blender.SpdBufGr2, BlenderCollGr2)
			SpdBufTmp1 = [Vector(obj) for obj in Blender.SpdBufGr1]
			SpdBufTmp2 = [Vector(obj) for obj in Blender.SpdBufGr2]
			geometryCollisionCheck(grp2, grp1, Blender.SpdBufGr2, Blender.SpdBufGr1, SpdBufTmp2, SpdBufTmp1, BlenderCollGr2, BlenderCollGr1, BlenderCollLimGr2, BlenderCollLimGr1, gr1rigid, corGr1spd, gr2rigid, corGr2spd)
			if selfcoll: geometryCollisionCheck(grp2, grp2, Blender.SpdBufGr2, Blender.SpdBufGr2, SpdBufTmp2, SpdBufTmp2, BlenderCollGr2, BlenderCollGr2, BlenderCollLimGr2, BlenderCollLimGr2, gr2rigid, corGr2spd, gr2rigid, corGr2spd)
			for coll1, spd1, spdTmp1, in izip(BlenderCollGr1, Blender.SpdBufGr1, SpdBufTmp1):
				if coll1 == 1: spd1[:] = spdTmp1
			for coll2, spd2, spdTmp2, in izip(BlenderCollGr2, Blender.SpdBufGr2, SpdBufTmp2):
				if coll2 == 1: spd2[:] = spdTmp2
			if not gr2rigid: springBodyAlgorithm(grp2, Blender.SpdBufGr2, Blender.EdgBufGr2, BlenderCollGr2)
			if not fix:
				for i in xrange(len(BlenderCollGr1)):
					if BlenderCollGr1[i] == 1 and BlenderCollLimGr1[i] <= _collLimit: BlenderCollGr1[i] = 0
				for i in xrange(len(BlenderCollGr2)):
					if BlenderCollGr2[i] == 1 and BlenderCollLimGr2[i] <= _collLimit: BlenderCollGr2[i] = 0
			if Get('curframe') > animated: correctSpeed(grp2, Blender.SpdBufGr2, bufPiece2, bufPieceOfs2, vecOld, BlenderCollGr2)
		airResistanceAndMaximumSpeed(Blender.SpdBufGr2)
		applySpeedVectors(grp2, Blender.SpdBufGr2)
	
		if Get('curframe') > sleep:
			if breaks:
				if not gr1rigid: edgeStressTest(grp1, Blender.SpdBufGr1, Blender.AniSpdBufGr1, Blender.EdgBufGr1, BlenderCollGr1, Blender.IdxBufGr1)
				if not gr2rigid: edgeStressTest(grp2, Blender.SpdBufGr2, Blender.AniSpdBufGr2, Blender.EdgBufGr2, BlenderCollGr2, Blender.IdxBufGr2)
		if Get('curframe') <= animated:
			for coll1, spd1, aniSpd1 in izip(BlenderCollGr1, Blender.SpdBufGr1, Blender.AniSpdBufGr1):
				if coll1 >= 2: spd1 -= aniSpd1
			for coll2, spd2, aniSpd2 in izip(BlenderCollGr2, Blender.SpdBufGr2, Blender.AniSpdBufGr2):
				if coll2 >= 2: spd2 -= aniSpd2

		for obj in grp1.objects:
			me = obj.getData(0,1)
			me.update() # "this method may be only temporary and may be removed in future (blender) releases." 
		for obj in grp2.objects:
			me = obj.getData(0,1)
			me.update() # "this method may be only temporary and may be removed in future (blender) releases." 


######################################
######################################
DemolitionProfiler_addScope(globals())

try:
	import psyco
except ImportError:
	print 'The demolition script will run faster if you install psyco.'
	print 'see: http://psyco.sourceforge.net/'
else:
	DemolitionProfiler_add(psyco.full)()
	print 'Found and used psyco. Therefore the demolition script will run faster. :)'

##########################
##########################
##########################
if __name__ == '__main__':
	main()

# vim:noexpandtab:sw=4:ts=4: