In [47]:
%pylab inline
import csv
import pandas as pd
from scipy.interpolate import SmoothBivariateSpline
from scipy.optimize import minimize
from mpl_toolkits.mplot3d import Axes3D
from threshold_functions import *
In [48]:
pylab.rcParams['savefig.dpi'] = 130
In [49]:
# Define thresholds
giveThreshold = 0.25
gapThreshold = 130
In [50]:
with open('cutout_factors','r') as f:
loaded_factor = eval(f.read())
dictArray = array(list(loaded_factor.items()),float)
cutoutRef = dictArray[:,0].astype(int)
ref = argsort(cutoutRef)
factor = dictArray[:,1]
cutoutRef = cutoutRef[ref]
factor = factor[ref]
with open('circle_cutout_factors','r') as f:
loaded_circle_factor = eval(f.read())
circleDictArray = array(list(loaded_circle_factor.items()),float)
circleDiameter = circleDictArray[:,0]
circleFactor = circleDictArray[:,1]
factor = append(factor,circleFactor)
with open('custom_cutout_factors','r') as f:
loaded_custom_factor = eval(f.read())
customDimensionsDict = {'5x13':[5,13],
'5x10':[5,10],
'5x8':[5,8],
'4x13':[4,13],
'4x10':[4,10],
'4x8':[4,8],
'4x6.5':[4,6.5],
'3x13':[3,13],
'3x9':[3,9],
'3x6.5':[3,6.5],
'3x5':[3,5]}
customDataLabelList = list(customDimensionsDict.keys())
custom_factors = [loaded_custom_factor[x] for x in customDataLabelList]
custom_widths = [customDimensionsDict[x][0] for x in customDataLabelList]
custom_lengths = [customDimensionsDict[x][1] for x in customDataLabelList]
factor = append(factor,custom_factors)
In [51]:
# loaded_factor.items()
In [52]:
cutoutRef
Out[52]:
In [52]:
In [53]:
cutout_dimensions = pd.DataFrame.from_csv('new_cutoutdata.csv',index_col =None)
ref = argsort(cutout_dimensions['cutoutRef'].values.astype(int))
width = cutout_dimensions['new_width'].values[ref]
length = cutout_dimensions['new_length'].values[ref]
weight = concatenate((ones(shape(width)), 1*ones(shape(circleDiameter)), ones(shape(custom_widths))))
width = append(width,circleDiameter)
length = append(length,circleDiameter)
nanFill = empty(shape(circleDiameter))
nanFill.fill(nan)
cutoutRef = append(cutoutRef,nanFill)
width = append(width,custom_widths)
length = append(length,custom_lengths)
nanFill = empty(shape(custom_widths))
nanFill.fill(nan)
cutoutRef = append(cutoutRef,nanFill)
# editRef = cutoutRef == 14
# width[editRef] = 4.3
# length[editRef] = 4.6
ratio = width / length
In [54]:
circleDiameter
Out[54]:
In [55]:
scatter(circleDiameter,circleFactor)
Out[55]:
In [56]:
scatter(width,ratio)
title('Data so far')
xlabel('Width (cm)')
ylabel('Ratio')
Out[56]:
In [57]:
# Spline fitting
widthMirrored = np.append(width,width)
ratioLogMirrored = np.log2(np.append(width/length,length/width))
factorMirrored = np.append(factor,factor)
weightMirrored = append(weight,weight)
splineFit = SmoothBivariateSpline(widthMirrored,ratioLogMirrored,factorMirrored,w=weightMirrored,s=len(width))
In [58]:
# Create surface mesh
xVec2 = linspace(width.min(),width.max(),100)
yVec2 = linspace(ratio.min(),ratio.max(),100)
xMesh2, yMesh2 = meshgrid(xVec2, yVec2)
gapMesh2 = angle_gap(xMesh2,log2(yMesh2),widthMirrored,ratioLogMirrored,1,2)
giveMesh2 = fit_give(xMesh2,log2(yMesh2),widthMirrored,ratioLogMirrored,factorMirrored,s=len(width))
zMesh2 = splineFit.ev(xMesh2, log2(yMesh2))
ref2 = (gapMesh2>gapThreshold) | (giveMesh2>giveThreshold)
zMesh2[ref2] = nan
In [59]:
# Plot the fit
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
bot = np.amin(factor) - 0.04 * np.ptp(factor)
top = np.amax(factor) + 0.04 * np.ptp(factor)
ax.contour(xMesh2,yMesh2,gapMesh2, offset=bot, levels=[gapThreshold], colors='r',alpha=0.3)
ax.contour(xMesh2,yMesh2,giveMesh2, offset=bot, levels=[giveThreshold], colors='g',alpha=0.3)
proxyGap, = ax.plot([width.min(),width.min()],[ratio.min(),ratio.min()],[bot,bot],'r-',alpha=0.3)
proxyGive, = ax.plot([width.min(),width.min()],[ratio.min(),ratio.min()],[bot,bot],'g-',alpha=0.3)
ax.plot_surface(xMesh2,yMesh2,zMesh2,alpha=0.3,rstride=4,cstride=4)
sc = ax.scatter(width,ratio,factor)
ax.set_title('Fit displayed against width and aspect ratio')
ax.set_xlabel('Width')
ax.set_ylabel('Ratio')
ax.set_zlabel('Factor')
legend([proxyGap,proxyGive],['Gap threshold', 'Give threshold'],loc=7,prop={'size':8})
right = ratio.min()-0.01*ratio.ptp()
ax.set_ylim(right,1)
ax.set_zlim(bot,top)
ax.view_init(elev=10, azim=-80)
In [60]:
# Remove data and predict
predictedFactor = zeros(shape(width))
predictionAngleGap = zeros(shape(width))
predictionFitGive = zeros(shape(width))
for i in range(len(width)):
x = width[i]
y = log2(width[i]/length[i])
adjWidth = np.delete(width,i)
adjLength = np.delete(length,i)
adjFactor = np.delete(factor,i)
adjWeight = delete(weight,i)
adjWidthMirrored = np.append(adjWidth,adjWidth)
adjRatioLogMirrored = np.log2(np.append(adjWidth/adjLength,adjLength/adjWidth))
adjFactorMirrored = np.append(adjFactor,adjFactor)
adjWeightMirrored = append(adjWeight,adjWeight)
adjSplineFit = SmoothBivariateSpline(adjWidthMirrored,adjRatioLogMirrored,adjFactorMirrored,w=adjWeightMirrored,s=len(width))
predictedFactor[i] = adjSplineFit.ev(x,y)
predictionAngleGap[i] = angle_gap(x,y,adjWidthMirrored,adjRatioLogMirrored,1,2)
predictionFitGive[i] = fit_give(x,y,adjWidthMirrored,adjRatioLogMirrored,adjFactorMirrored,s=len(width))
difference = factor - predictedFactor
withinThresholds = ((predictionAngleGap<gapThreshold) &
(predictionFitGive<giveThreshold)
)
In [61]:
# Display table
df = pd.DataFrame(transpose([cutoutRef,width, length, factor, predictedFactor, difference,
predictionAngleGap, predictionFitGive, withinThresholds]))
df.columns = ['Reference','width', 'length', 'factor', 'Predicted Factor', 'Difference',
'Angle Gap', 'Fit Give', 'Within thresholds']
df.to_csv(path_or_buf = "cutoutdata.csv", index=False)
df
Out[61]:
In [62]:
print("There are", sum(withinThresholds), "points within give and gap thresholds."+
"\nThe standard deviation of the prediction differences for these is", std(difference[withinThresholds]))
In [63]:
# Histogram
fig = plt.figure()
ax = fig.add_subplot(111)
n, bins, patches = hist(difference[withinThresholds])
setp(patches, 'facecolor', 'g', 'alpha', 0.5)
ax.set_xlabel('Prediction differences')
ax.set_ylabel('Frequency')
ax.set_title('Histogram of prediction differences')
Out[63]:
In [64]:
# Visual representation of differences
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
ax.plot_surface(xMesh2,yMesh2,zMesh2,alpha=0.3,rstride=4,cstride=4)
diffColour = copy(difference)
diffColour[~withinThresholds] = nan
colourRange = nanmax(abs(diffColour))
sc = ax.scatter(width,ratio,factor,c=diffColour, vmin=-colourRange, vmax=colourRange)
colorbar(sc)
ax.set_title('Differences')
ax.set_xlabel('Width')
ax.set_ylabel('Ratio')
ax.set_zlabel('Factor')
ax.set_ylim(right,1)
ax.set_zlim(bot,top)
ax.view_init(elev=10, azim=-80)
In [65]:
diffTest = abs(difference)
diffTest[~withinThresholds] = 0
testRef = argmax(diffTest)
cutoutRef[testRef]
Out[65]:
In [66]:
factor[testRef]
Out[66]:
In [67]:
scatter(width,ratio)
scatter(width[testRef],ratio[testRef],c='red')
Out[67]:
In [68]:
def to_minimise(cut_increase,ref):
test_width = width[ref] + cut_increase
test_length = length[ref] + cut_increase
test_ratio = test_width / test_length
return (factor[ref] - splineFit.ev(test_width,log2(test_ratio)))**2
In [69]:
output = minimize(to_minimise,0,args=(testRef,))
output
Out[69]:
In [69]: