%pylab inline
import pandas as pd
import csv
from scipy.interpolate import SmoothBivariateSpline, bisplrep
from threshold_functions import *
from equivalent_ellipse import *
Populating the interactive namespace from numpy and matplotlib
x_list = list()
y_list = list()
with open('../data/cutout_x.csv', 'r') as x_csvfile:
with open('../data/cutout_y.csv', 'r') as y_csvfile:
x_reader = csv.reader(x_csvfile, delimiter=',', lineterminator='\n')
y_reader = csv.reader(y_csvfile, delimiter=',', lineterminator='\n')
for row in x_reader:
x_list += [row]
for row in y_reader:
y_list += [row]
num_cutouts = len(x_list)
x_array = [0,]*num_cutouts
y_array = [0,]*num_cutouts
for i in range(num_cutouts):
x_array[i] = array(x_list[i], dtype='float')
y_array[i] = array(y_list[i], dtype='float')
cutout = [0,]*num_cutouts
for i in range(num_cutouts):
cutout[i] = shapely_cutout(x_array[i],y_array[i])
cutout_dimensions = pd.DataFrame.from_csv('../data/cutout_dimensions.csv')
ellipse_def = transpose(array([cutout_dimensions['a'].values,
cutout_dimensions['b'].values,
cutout_dimensions['width'].values,
cutout_dimensions['length'].values,
cutout_dimensions['theta'].values]))
width_init = cutout_dimensions['width'].values
cutoutRef_init = arange(len(width_init))
sorted_ref = cutoutRef_init[argsort(width_init)]
similar_list = diff(width_init[sorted_ref]) < 0.0001
duplicates_removed_ref = sorted_ref[~similar_list]
boundBox = cutout_dimensions['box_bounds'].values
too_big = boundBox[duplicates_removed_ref] > 11
too_big_removed_ref = duplicates_removed_ref[~too_big]
final_ref = too_big_removed_ref
width = cutout_dimensions['width'].values[final_ref]
length = cutout_dimensions['length'].values[final_ref]
cutoutRef = final_ref
ratio = width/length
# plot(append(width,width),append(log2(ratio),log2(1/ratio)),'.')
def give_calc(xGrid,yGrid,width,ratio,w):
w = w/max(w)
ref = (w < 0.01)
width = delete(width,find(ref))
ratio = delete(ratio,find(ref))
w = delete(w,find(ref))
return fit_give(xGrid,log2(yGrid),
append(width,width),
append(log2(ratio),log2(1/ratio)),
zeros(len(width)*2),append(w,w))
def gap_calc(xGrid,yGrid,width,ratio):
return angle_gap(xGrid,log2(yGrid),append(width,width),
append(log2(ratio),log2(1/ratio)),1,2)
# w = ones(len(width))
# valid = (give_calc(width,ratio,width,ratio,w) < 0.15) & (gap_calc(width,ratio,width,ratio) < 130)
# curr_give = give_calc(width[valid],ratio[valid],width,ratio,w)
# curr_gap = gap_calc(width[valid],ratio[valid],width,ratio)
# curr_gap
# curr_give = give_calc(width[valid],ratio[valid],width,ratio,w)
# curr_give[2] = 0.2
# sum((curr_give/0.15)**15)
# curr_gap[2] = 100
# sum((curr_gap/130)**30)
# sum(w)
# sort(rand(100))[0:-5]
global keeping_track
def toMinimise(w):
if (sum(w) > 13):
curr_give = give_calc(width,ratio,width,ratio,w)
else:
curr_give = ones(len(w))
weights_bias = sum((2*exp(-((w-0.5)*2)**2) + 3*w - 1)/3)
too_few_bias = exp(-(sum(w)-15))
give_bias_init = 4/(1 + exp(-(curr_give-0.15)*200))
give_bias = sum(sort(give_bias_init)[0:-5])
output = weights_bias + give_bias + too_few_bias
return output
# t = linspace(0,20,1000)
# y = exp(-(t-13))
# plot(t,y)
# ylim([0,1])
# t = linspace(0,1,1000)
# y = (2*exp(-((t-0.5)*2)**2) + 3*t - 1)/3
# plot(t,y)
# t = linspace(0,1,1000)
# y = 4/(1 + exp(-(t-0.15)*200))
# t5 = 0.15
# y5 = 4/(1 + exp(-(t5-0.15)*200))
# print(y5)
# plot(t,y)
# scatter(t5,y5)
def print_fun(w, f, accepted):
global successCount_basinhopping
global minVals_basinhopping
global nSuccess_basinhopping
print("Output %.4f" % (f))
w_show = reshape(w,[3,23])
imshow(w_show, interpolation='nearest', cmap=cm.gray, vmin=0, vmax=1)
show()
bounds = ([0,1],)*len(width)
# x0 = ones(len(width))
x0 = loadtxt('x0')
keeping_track = 0
global nSuccess_basinhopping
nSuccess_basinhopping = 5
global minVals_basinhopping
minVals_basinhopping = np.empty(nSuccess_basinhopping)
minVals_basinhopping[:] = np.nan
global successCount_basinhopping
successCount_basinhopping = 0
minimizer_kwargs = {"tol":0.01,"method": 'SLSQP', "bounds": bounds}
# curr_give = give_calc(width[valid],ratio[valid],width,ratio,ones(len(width)))
# curr_give
# give_calc(width[valid],ratio[valid],width,ratio,x0)
# (give_calc(width[valid],ratio[valid],width,ratio,ones(len(width))) -
# give_calc(width[valid],ratio[valid],width,ratio,x0))
resltn = 25 # Must be odd
xRange = [min(width),max(width)]
yRange = [min(ratio),max(ratio)]
xVec = linspace(xRange[0],xRange[1],resltn)
yVec = linspace(yRange[0],yRange[1],resltn)
xMesh, yMesh = meshgrid(xVec, yVec)
giveMesh = give_calc(xMesh,yMesh,width,ratio,ones(len(width)))
while 1:
output = basinhopping(toMinimise,x0,
niter=20,minimizer_kwargs=minimizer_kwargs,
callback=print_fun)
print("Final %.4f" % (output.fun))
w_show = reshape(output.x,[3,23])
imshow(w_show, interpolation='nearest', cmap=cm.gray, vmin=0, vmax=1)
show()
w = output.x/max(output.x)
# give_pic = give_calc(width,ratio,width,ratio,w)
# print(around(give_pic,4))
saveStuffs = w > 0.5
widthSave = width[saveStuffs]
ratioSave = ratio[saveStuffs]
plot(width[~saveStuffs],ratio[~saveStuffs],'rx')
scatter(widthSave,ratioSave)
print("Currently %d cutouts require measurement:" % (sum(saveStuffs)))
newGiveMesh = give_calc(xMesh,yMesh,width,ratio,w)
contour(xMesh,yMesh,giveMesh, levels=[0.05,0.1,0.15], colors='g',alpha=0.7)
contour(xMesh,yMesh,newGiveMesh, levels=[0.05,0.1,0.15], colors='r',alpha=0.7)
title('Optimisation of where to measure')
xlabel('Width (cm)')
ylabel('Ratio')
x0 = w
# print(output)
savetxt('cutouts_to_be_measured',cutoutRef[saveStuffs])
savetxt('x0',x0)
show()
Output 71.5533
Output 74.4339
Output 72.9055
Output 71.1501
Output 70.4769
Output 74.4066
Output 76.1726
Output 71.1335
Output 69.6944
Output 67.4352
Output 77.5026
Output 74.9814
Output 78.1372
Output 73.4267
Output 73.0739
Output 86.0535
Output 69.8285
Output 86.1100
Output 70.5087
Output 76.8883
Final 62.6919
Currently 37 cutouts require measurement:
Output 69.6023
Output 74.6371
Output 78.0186
Output 75.0176
Output 73.5341
Output 67.1599
Output 71.1460
Output 75.8123
# x0 = output.x
# oldF = 200
# savetxt('oldF',oldF)
# newGiveMesh = give_calc(xMesh,yMesh,width,ratio,x0)
# contour(xMesh,yMesh,giveMesh, levels=[0.05,0.1,0.15], colors='g',alpha=0.7)
# contour(xMesh,yMesh,newGiveMesh, levels=[0.05,0.1,0.15], colors='r',alpha=0.7)
# len(width) - sum(w)
# saveStuffs = x0 > 0.5
# widthSave = width[saveStuffs]
# ratioSave = ratio[saveStuffs]
# plot(width[~saveStuffs],ratio[~saveStuffs],'rx')
# scatter(widthSave,ratioSave)
# sum(saveStuffs)
# oldF
# output.fun