Spinny blob ball crucial error fix

import Image
import random
import math
import itertools
import copy



def line_pixels(ord1,ord2):
    # Bresenham's line algorithm
    # useful function for plotting a 2d line from 2 co-ordinates - ord1 to ord2
    # ord1 and ord2 must be tuples in the form (x,y)

    x0 = int(ord1[0]); y0 = int(ord1[1])
    x1 = int(ord2[0]); y1 = int(ord2[1])
    dx = abs(x1 - x0)
    dy = abs(y1 - y0)
    x, y = x0, y0
    sx = -1 if x0 > x1 else 1
    sy = -1 if y0 > y1 else 1

    pixlist = [(x0,y0)]
    if dx > dy:
        err = dx / 2.0
        while x != x1:
            pixlist.append((x, y))
            err -= dy
            if err < 0:
                y += sy
                err += dx
            x += sx
    else:
        err = dy / 2.0
        while y != y1:
            pixlist.append((x, y))
            err -= dx
            if err < 0:
                x += sx
                err += dy
            y += sy
    pixlist.append((x, y))
    return pixlist


i = 0

width = 26          # width/height/breadth of the" box"
                    # (i.e, the number of discrete steps)

hooke = 0.14444444  # a spring-like constant that scales down the force
                    # on adjacent points in the field

xsize, ysize = (400,400) # image size in pixels

points = [[[[0.0,0.0] for x in xrange(width)] for y in xrange(width)] for z in xrange(width)]
        # nested lists representing a 3d field with values
        # at discrete points, accessible by "points[z][y][x]"


for bbb in xrange(210000):


    im = Image.new("RGB",(xsize,ysize))
    loadobject = im.load()


    for z0, l0 in enumerate(points):
        for y0, c0 in enumerate(l0):
            for x0, p0 in enumerate(c0):
                # a bunch of 3d-projection type stuff, projects the 3d
                # co-ordinate onto the 2d image plane, mess around with
                # the scale_factors and k_constant to work with larger fields/
                # different levels of perspective
                scale_fact = 5
                k_constant = 0.013


                y1 = y0 - (width/2.0)
                x1 = (x0 - (width/2.0))*(1 + (k_constant*y1))
                y1 += -(z0-width/2)*(1 + (k_constant*y1))

                x1 *= scale_fact*1.5
                y1 *= scale_fact

                x1 += (xsize/2)
                y1 += (ysize/2)

                x1 = int(x1)
                y1 = int(y1)

                try:
                    """
                       these conditions determine whether a point is plotted
                       depending on the displacement value of the field. the first
                       2 frames always contain the full cube plotted as a reference

                       the second conditional means that only points in the
                       field with values greater than positive 6 are plotted,
                       negative values are ignored.
                    """
                    if i > 2:
                        if p0[0] < 6:
                            continue

                    scalecol = (100/width)

                    col = (128+(y0*scalecol),128-(x0*scalecol),128-(z0*scalecol))
                    loadobject[x1,y1] = col
                    loadobject[x1+1,y1] = col
                    loadobject[x1-1,y1] = col
                    loadobject[x1,y1+1] = col
                    loadobject[x1,y1-1] = col
                except:
                    pass

    # save the image and increment i
    im.save("test//TEST3_%05d.png" % i)
    i += 1


    #update all of the new points based on previous positions
    points = copy.deepcopy(newpoints)

    for np_z, l0 in enumerate(points):
        for np_y, c0 in enumerate(l0):
            for np_x, p0 in enumerate(c0):
                # all this stuff means that the velocity and displacement of the
                # field are calculated from current differences in poision

                vel = p0[1]
                adj_force = 0.0

                a1 = points[np_z][np_y][(np_x-1)%len(l0)][0]
                a2 = points[np_z][np_y][(np_x+1)%len(l0)][0]
                a3 = points[np_z][(np_y+1)%len(points)][(np_x)%len(points)][0]
                a4 = points[np_z][(np_y-1)%len(points)][(np_x)%len(points)][0]
                a5 = points[(np_z+1)%width][(np_y)][(np_x)][0]
                a6 = points[(np_z-1)%width][(np_y)][(np_x)][0]

                adj_force = float(-6*p0[0]) + a1 + a2 + a3 + a4 + a5 + a6

                vel += adj_force * hooke


                disp = p0[0] + vel
                newpoints[np_z][np_y][np_x] = [disp,vel]

    print i,

    newpoints = copy.deepcopy(points)

    close_edges = False  ### when set to true, this ensures that
                         ### all boundaries are set to 0 displacement
    if close_edges:
        for x0 in xrange(width):
            for y0 in xrange(width):
                for z0 in xrange(width):
                    if x0 == 0 or x0 == width-1:
                        points[x0][y0][z0] == [0.0,0.0]
                    if y0 == 0 or y0 == width-1:
                        points[x0][y0][z0] == [0.0,0.0]
                    if z0 == 0 or z0 == width-1:
                        points[x0][y0][z0] == [0.0,0.0]
    """
    for the first 199 frames there will be 2x points orbiting with a radius
    of 1/9 of the width, with a positive and negative displacement

    """
    upper_limit = 199
    if i < upper_limit:
        radius = width/9
        theta = i*(math.pi/15)
        points[int((width/2)+(math.sin(theta)*0.2*radius))][int((width/2)+(math.sin(theta)*radius))][int((width/2)+(math.cos(theta)*radius))] = [26,26]
        points[int((width/2)-(math.sin(theta)*0.2*radius))][int((width/2)-(math.sin(theta)*radius))][int((width/2)-(math.cos(theta)*radius))] = [-26,-26]