I am new to Python. I have two scripts for generating and plotting a 2D lattice which values represent a spatially distributed attribute in the form of a fractal image. The first script contains the function generating the fractal (orthodox name: Spectral Synthesis with Inverse Fourier Transform), while the second simply calls the function, plots it, ad writes it to a .csv file. The function relies on multiple random.gauss(0.0, sigma)
calls for obtaining normally distributed random numbers. It also calls random.random()
since it needs other random numbers. Below are the codes. Thanks to anyone who will help!
The problem: the plot always shows the same pattern, hence making me think that there is something NOT random in the code. My intention is to obtain a different fractal image every time I run the code. What makes it NOT random? Is that seed=122
in the second script?
PS: please note that since the procedure introduces a given degree of periodicity, and since I do not want that periodicity to be shown, the plotted lattice is smaller than the generated lattice. In particular, I am using size(plotted)=101x101, while size(generated)=256x256.
First script: the function itself
from __future__ import division #Avoids the floor of the mathematical result of division if the args are ints or longs
import numpy
import random
import math
import pylab
def SpectralSynthesisFM2D(max_level, sigma, H, seed=0, normalise=True, bounds=[0,1]):
"""
________________________________________________________________________
Args:
max_level : Maximum number of recursions (N = 2^max_level)
sigma : Initial standard deviation
H : Roughness constant (Hurst exponent), varies form 0.0 to 1.0
H=0.8 is a good representation of many natural phenomena (Voss, 1985)
seed : seed value for random number generator
normalise : normalizes the data using bound
bounds : used for normalization of the grid data
Result:
Output is given in the form of an array (grid) which holds surface
values for a square region.
_________________________________________________________________________
"""
N = int(2**max_level)
A = numpy.zeros((N,N), dtype = complex)
random.seed() #Seeds the random number generator
PI = 3.141592
for i in range(0,int((N-1)/2)):
for j in range(0,int((N-1)/2)):
phase = 2*PI*random.random() #/random.randrange(1,Arand)
if i != 0 or j != 0:
rad = pow((i*i + j*j),(-(H+1)/2) )*random.gauss(0.0, sigma)
else:
rad = 0.0
A[i][j] = rad*math.cos(phase) + rad*math.sin(phase)*j
if i ==0:
i0 = 0
else:
i0 = N - i
if j==0:
j0 = 0
else:
j0 = N - j
A[i0][j0] = rad * math.cos(phase) - rad*math.sin(phase)*j
for i in range(1,int((N-1)/2)):
for j in range(1,int((N-1)/2)):
phase = 2*PI*random.random() #/random.randrange(1,Arand)
rad = pow((i*i + j*j),(-(H+1)/2) )*random.gauss(0.0, sigma)
A[i][N-j] = rad * math.cos(phase) + rad* math.sin(phase)*j
A[N-i][j] = rad * math.cos(phase) - rad* math.sin(phase)*j
Grid = numpy.real(pylab.ifft2(( A ))) #Implements the Discrete Inverse Fast Fourier Transform on a 2D lattice (grid)
if(normalise):
Grid += numpy.amin(Grid)*-1 + bounds[0]
Grid = (Grid/numpy.amax(Grid)) * bounds[1]
return Grid
Second script: plotting the function and obtaining the fractal image:
from __future__ import division #Avoids the floor of the mathematical result of division if the args are ints or longs
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.colors as mc
import SSFM2D
#Parameter assignments
max_level=8 #This is the exponent controlling the grid size. In this case N=2^8=256. Use only integers.
N=2**max_level
sigma=1 #Variation for random Gauss generation (standardised normal distribution)
H=0.8 #Hurst exponent (0.8 is a recommended value for natural phenomena)
seed=122 #Setting the seed for random Gauss generation
print ('The lattice size is '+str(N)+'x'+str(N))
#Lattice initialization
Lattice=np.zeros((256,256))
#Calling Spectral fBm function
Lattice=SSFM2D.SpectralSynthesisFM2D(max_level, sigma, H, seed, normalise=True, bounds=[0,1])
#Plotting the lattice
M=np.zeros((101,101)) #Lattice to be plotted. Size(M) != size(Lattice) to minimize visual fBm periodic features.
for i in range(0,101):
for j in range(0,101):
M[i][j]=Lattice[i][j]
#Lattice statistics
print M[-101:,-101:].sum()
print M[-101:,-101:].max()
print M[-101:,-101:].min()
print M[-101:,-101:].mean()
#Writing X, Y and values to a .csv file from scratch
import numpy
import csv
with open('C:\\Users\\Francesco\\Desktop\\Python_files\\csv\\fBm_256x256.csv', 'w') as f: # Change directory if necessary
writer = csv.writer(f)
writer.writerow(['X', 'Y', 'Value'])
for (x, y), val in numpy.ndenumerate(M):
writer.writerow([x, y, val])
# Plotting - Showing interpolation of randomization
plt.imshow(M[-101:,-101:].T, origin='lower',interpolation='nearest',cmap='Blues', norm=mc.Normalize(vmin=0,vmax=M.max()))
title_string=('fBm: Inverse FFT on Spectral Synthesis')
subtitle_string=('Lattice size: 101x101')
plt.suptitle(title_string, y=0.99, fontsize=17)
plt.title(subtitle_string, fontsize=9)
plt.show()
# Makes a custom list of tick mark intervals for colour bar (assumes minimum is always zero)
numberOfTicks = 5
ticksListIncrement = M.max()/(numberOfTicks)
ticksList = []
for i in range((numberOfTicks+1)):
ticksList.append(ticksListIncrement * i)
cb=plt.colorbar(orientation='horizontal', format='%0.2f', ticks=ticksList)
cb.set_label('Values')
plt.show()
plt.xlim(0, 100)
plt.xlabel('Easting (Cells)')
plt.ylim(100, 0)
plt.ylabel('Northing (Cells)')
usually a seed
is used to initialize a random number generator. same seed
-> same sequence of random numbers (which is very userful for testing).
import random
random.seed(22)
for i in range(10):
print(random.randint(0, 100))
so no matter how many times you run that program (or even on what machine), the output will always be the same.
but when you want 'real' randomness (or at least better randomness) you better call
random.seed() # or: random.seed(None)
with no arguments. python then calls seed
using the current system time - which provides a source of some randomness.
as i do not know your library i can't tell for sure but i'd bet it works in a simliar way. you could try to set seed=None
or seed=random.randint(something)
.
in any case: be sure to read the api documentation of your library and check what seed
does.