Welcome, guest | Sign In | My Account | Store | Cart

IFS fractal dimension calculation using box counting method (Python recipe) by FB36
ActiveState Code (http://code.activestate.com/recipes/578037/)

IFS fractal dimension calculation using box-counting method.

Python, 134 lines
Download
  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
# IFS fractal dimension calculation using box-counting method
# http://en.wikipedia.org/wiki/Iterated_function_system
# http://en.wikipedia.org/wiki/Fractal_dimension
# http://en.wikipedia.org/wiki/Box-counting_dimension
# http://en.wikipedia.org/wiki/Simple_linear_regression
# FB - 20120211
import math
import random
from PIL import Image
imgx = 512
imgy = 512 # will be auto-re-adjusted according to aspect ratio of the fractal
maxIt = imgx * imgy * 2

fractalName = "Barnsley Fern"
mat=[[0.0, 0.0, 0.0, 0.16, 0.0, 0.0, 0.01],
     [0.85, 0.04, -0.04, 0.85, 0.0, 1.6, 0.85],
     [0.2, -0.26, 0.23, 0.22, 0.0, 1.6, 0.07],
     [-0.15, 0.28, 0.26, 0.24, 0.0, 0.44, 0.07]]

##fractalName = "Centipede"
##mat = [[0.824074, 0.281482, -0.212346,  0.864198, -1.882290, -0.110607, 0.787473],
##       [0.088272, 0.520988, -0.463889, -0.377778,  0.785360,  8.095795, 0.212527]]

##fractalName = "Levy C curve"
##mat = [[0.5, -0.5, 0.5, 0.5, 0.0, 0.0, 0.5],
##       [0.5, 0.5, -0.5, 0.5, 0.5, 0.5, 0.5]]

##fractalName = "Dragon curve"
##mat = [[0.5, -0.5, 0.5, 0.5, 0.0, 0.0, 0.5],
##       [-0.5, -0.5, 0.5, -0.5, 1.0, 0.0, 0.5]]

##fractalName = "Sierpinski Triangle"
##mat = [[0.5, 0.0, 0.0, 0.5, 0.0, 0.0, 0.33],
##       [0.5, 0.0, 0.0, 0.5, 0.5, 0.0, 0.33],
##       [0.5, 0.0, 0.0, 0.5, 0.0, 0.5, 0.33]]

##fractalName = "Sierpinski Carpet"
##mat = [[0.333, 0.0, 0.0, 0.333, 0.000, 0.999, 0.125],
##       [0.333, 0.0, 0.0, 0.333, 0.333, 0.999, 0.125],
##       [0.333, 0.0, 0.0, 0.333, 0.666, 0.999, 0.125],
##       [0.333, 0.0, 0.0, 0.333, 0.000, 0.666, 0.125],
##       [0.333, 0.0, 0.0, 0.333, 0.666, 0.666, 0.125],
##       [0.333, 0.0, 0.0, 0.333, 0.000, 0.333, 0.125],
##       [0.333, 0.0, 0.0, 0.333, 0.333, 0.333, 0.125],
##       [0.333, 0.0, 0.0, 0.333, 0.666, 0.333, 0.125]]

m = len(mat) # number of IFS transformations
# find xmin, xmax, ymin, ymax of the fractal using IFS algorithm
x = mat[0][4]
y = mat[0][5] 
xa = x
xb = x
ya = y
yb = y
for k in range(maxIt):
    i = random.randint(0, m - 1)
    if random.random() <= mat[i][6]:
        x0 = x * mat[i][0] + y * mat[i][1] + mat[i][4] 
        y  = x * mat[i][2] + y * mat[i][3] + mat[i][5] 
        x = x0 
        if x < xa:
            xa = x
        if x > xb:
            xb = x
        if y < ya:
            ya = y
        if y > yb:
            yb = y

imgy = int(imgy * (yb - ya) / (xb - xa)) # auto-re-adjust the aspect ratio 
image = Image.new("RGB", (imgx, imgy))
pixels = image.load()

# drawing using IFS algorithm
theColor = (255, 255, 255)
x=0.0
y=0.0 
for k in range(maxIt):
    i = random.randint(0, m - 1)
    if random.random() <= mat[i][6]:
        x0 = x * mat[i][0] + y * mat[i][1] + mat[i][4] 
        y  = x * mat[i][2] + y * mat[i][3] + mat[i][5] 
        x = x0 
        jx = int((x - xa) / (xb - xa) * (imgx - 1)) 
        jy = (imgy - 1) - int((y - ya) / (yb - ya) * (imgy - 1))
        if jx >= 0 and jx < imgx and jy >= 0 and jy < imgy:
            pixels[jx, jy] = theColor
image.save(fractalName + " fractal.png", "PNG")

# fractal dimension calculation using box-counting method
b = 2 # initial box size in pixels
f = 2 # box scaling factor
n = 3 # number of graph points for simple linear regression
gx = [] # x coordinates of graph points
gy = [] # y coordinates of graph points
for ib in range(n):
    bs = b * f ** ib # box size in pixels
    bnx = int(imgx / bs) # of boxes in x direction of image
    bny = int(imgy / bs) # of boxes in y direction of image
    boxCount = 0
    for by in range(bny):
        for bx in range(bnx):
            # if there are any pixels in the box then increase box count
            foundPixel = False
            for ky in range(bs):
                for kx in range(bs):
                    if pixels[bs * bx + kx, bs * by + ky] == theColor:
                        foundPixel = True
                        boxCount += 1
                        break
                if foundPixel:
                    break
    gx.append(math.log(1.0 / bs))
    gy.append(math.log(boxCount))
    
# simple linear regression
x_ = 0.0
y_ = 0.0
x2_ = 0.0
xy_ = 0.0
for j in range(n):
    x_ += gx[j]
    y_ += gy[j]
    x2_ += gx[j] ** 2
    xy_ += gx[j] * gy[j]
x_ = x_ / n
y_ = y_ / n
x2_ = x2_ / n
xy_ = xy_ / n
b = (xy_ - x_ * y_) / (x2_ - x_ ** 2) # slope of the regression line
# a = y_ - b * x_

print "Fractal Name: " + fractalName
print "Estimated Fractal Dimension: " + str(b)
Tags: fractal, math, mathematics