## Saturday, June 21, 2008

### Python: Isotropic fractal surface generator

`#! /usr/bin/env pythonfrom __future__ import divisionimport Imagefrom scipy import *class FractalSurface(object): '''Generate isotropic fractal surface image using spectral synthesis method [1, p.] References: 1. Yuval Fisher, Michael McGuire, The Science of Fractal Images, 1988 ''' def __init__(self,fd=2.5, N=256): self.N=N self.H=1-(fd-2); self.X=zeros((self.N,self.N),complex) self.A=zeros((self.N,self.N),complex) self.img=Image.Image() def genSurface(self): '''Spectral synthesis method ''' N=self.N; A=self.A powerr=-(self.H+1.0)/2.0 for i in range(int(N/2)+1): for j in range(int(N/2)+1): phase=2*pi*rand() if i is not 0 or j is not 0: rad=(i*i+j*j)**powerr*random.normal() else: rad=0.0 self.A[i,j]=complex(rad*cos(phase),rad*sin(phase)) if i is 0: i0=0.0 else: i0=N-i if j is 0: j0=0.0 else: j0=N-j self.A[i0,j0]=complex(rad*cos(phase),-rad*sin(phase)) self.A.imag[N/2][0]=0.0 self.A.imag[0,N/2]=0.0 self.A.imag[N/2][N/2]=0.0 for i in range(1,int(N/2)): for j in range(1,int(N/2)): phase=2*pi*rand() rad=(i*i+j*j)**powerr*random.normal() self.A[i,N-j]=complex(rad*cos(phase),rad*sin(phase)) self.A[N-i,j]=complex(rad*cos(phase),-rad*sin(phase)) itemp=fftpack.ifft2(self.A) itemp=itemp-itemp.min() self.X=itemp def genImage(self): #Aa=abs(Aa) Aa=self.X im=Aa.real/Aa.real.max()*255.0 self.img=Image.fromarray(uint8(im)) #img2=Image.fromstring("L",(N,N),uint8(im).tostring()) def showImg(self): self.img.show() def saveImg(self,fname="fs.tiff"): self.img.save(fname) def getFSimg(self): return self.imgdef main(): fs=FractalSurface() fs.genSurface() fs.genImage() fs.saveImg() fs.showImg()if __name__ == '__main__': main() `

FD=2.1

FD=2.5

FD=2.9