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

IFS fractal dimension calculation using box-counting method.

Python, 134 lines
 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)
 Created by FB36 on Sun, 12 Feb 2012 (MIT)