Homec4science
Diffusion Roughness_ENAC_Project (master) Browse

R10592/BOX COUNTING/BoxCounting.py
86c9f01a7d54master

Last Change
Raw File

BoxCounting.py
View Options

import numpy as num
import matplotlib.pyplot as plt
import math
import time
import csv
from mpl_toolkits.mplot3d import Axes3D
#latest version and best optimization method by Gianluca Costagliola
#The code calculate the number of boxes intersected for chosen scaling factors for an input file which contains x,y,z coo. in each column
#####################################
#input file contains x,y,z coordinates in each column
with open('inputfile.txt','r') as f:
lines = list(csv.reader(f, delimiter = ' ', skipinitialspace = True))
x = num.zeros(len(lines));
y = num.zeros(len(lines));
z = num.zeros(len(lines));
for i in range(len(lines)):
x[i] = lines[i][0];
y[i] = lines[i][1];
z[i] = lines[i][2];
#calculation of the centroid of the object
Cx=num.mean(x)
Cy=num.mean(y)
Cz=num.mean(z)
#initialization of new coordinates
xc = num.zeros(len(lines));
yc = num.zeros(len(lines));
zc = num.zeros(len(lines));
#Translation of the object to have the centroid at the origin (Cx = Cy = Cz = 0)
for i in range(len(lines)):
xc[i] = x[i] - Cx
yc[i] = y[i] - Cy
zc[i] = z[i] - Cz
XBoxmin=1.0000001*min(min(xc), min(yc), min(zc)); #Min. coordinate of the first box
XBoxmax=1.0000001*max(max(xc), max(yc), max(zc)); #Max. coordinate of the first box
print('number of points:')
print(len(lines))
print('Coordinates min. and max. of the first box:')
print(XBoxmin)
print(XBoxmax)
LBox = XBoxmax - XBoxmin; #Side length of the first box
print('Length of the side of the first box:')
print(LBox)
steps=18; #Number of steps for the computation
#sf=num.zeros(steps); #scaling factors
sf=[5,7,10,15,20,25,30,40,50,60,70,80,100,120,150,200,250,300];
NBox=num.zeros(steps,dtype=int); #nr of total boxes for each scaling factor
for i in range(steps):
NBox[i]=int(sf[i]**3);
print('scaling factors:')
print(sf)
print('Nr. of total boxes:')
print(NBox)
NN=num.zeros(steps,dtype=int);
t = time.time()
for n in range(0,steps): #loop scaling factor
lbox=LBox/sf[n]
array=num.zeros(NBox[n],dtype=int)
for i in range(0,len(lines)):
ii=int( (xc[i]-XBoxmin)/lbox)
ij=int( (yc[i]-XBoxmin)/lbox)
ik=int( (zc[i]-XBoxmin)/lbox)
ind=int(ii+sf[n]*ij+sf[n]*sf[n]*ik)
array[ind]=array[ind]+1
for s in range(0,NBox[n]):
if array[s]>0:
NN[n]=NN[n]+1
print('Nr of boxes intersected:')
print(NN)
elapsed = time.time() - t
print('time elapsed [s]:')
print(elapsed)
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
ax.scatter(xc, yc, zc, c='r', marker='o')
plt.show()
c4science · Help