Skip to content

Instantly share code, notes, and snippets.

Show Gist options
  • Save Houdini/290833 to your computer and use it in GitHub Desktop.
Save Houdini/290833 to your computer and use it in GitHub Desktop.
import pylab
from math import exp
def scatter_plot(x, y, r, color='b'):
for i in range(len(x)):
pylab.scatter(x[i], y[i], s=r, c=color)
#x1 = [0, 1., 0.1, 2, 3, 2.3, 3.1, 1.7, 1.9]
#y1 = [0.1, 0.2, 0.8, 1, 1.5, 1.9, 2.1, 1.7, 1.9]
#x2 = [0.1, 0.5, 2.1, 2.3, 0.8, 2 , 1.2, 1 ,-0.1, 1.5, 0.1]
#y2 = [1.5, 3 , 2.9, 2.9, 1.5, 2.5, 3.2, 2.9, 1.8, 1.54, 2.1]
x1 = [-3, -4, -5]
y1 = [-0.1, 2, 0]
x2 = [1, 2, 3]
y2 = [0.5, 0, -0.9]
eps = 0.01
sigma = 1.3
def get_function_value(x, y, mx, my, sigma_this):
value = 1./(2*3.14159265*sigma_this**2)*exp(-0.5/sigma_this**2*((x-mx)**2 + (y-my)**2))
return value
x = [i*0.001-3 for i in range(6000)]
y = [get_function_value(xi,0,0,0, sigma) for xi in x]
def do_one_plot():
scatter_plot(x1, y1, 200)
scatter_plot(x2, y2, 200, 'g')
pylab.plot(x, y)
def get_coef(h, sigma):
def find_border(min_or_max, x_or_y):
return eval(min_or_max + '([' + min_or_max + '(' + x_or_y + '1),' + min_or_max + '(' + x_or_y + '2)])')
left_border = find_border('min', 'x')
right_border = find_border('max', 'x')
top_border = find_border('max', 'y')
bottom_border = find_border('min', 'y')
coef = 0
for i in range(int(round((right_border - left_border)/h))):
for j in range(int(round((top_border - bottom_border)/h))):
center_i = ((i*h + left_border) + h)/2
center_j = ((j*h + bottom_border) + h)/2
sum1 = 0
sum2 = 0
for i_point in range(len(x1)):
sum1 += get_function_value(center_i, center_j, x1[i_point], y1[i_point], sigma)
for i_point in range(len(x1)):
sum2 += get_function_value(center_i, center_j, x2[i_point], y2[i_point], sigma)
sum = min([sum1, sum2])
coef += sum*h**2
print "coef = ", coef
print "point = ",
return coef
def plot_full(step_length, all_steps):
steps = [i*step_length for i in range(1,all_steps)]
coefs = [get_coef(0.01, step) for step in steps]
max_coef_i = 0
max_coef = coefs[0]
for i in range(len(steps)):
if coefs[i] > max_coef:
max_coef = coefs[i]
max_coef_i = steps[i]
print "max_coef_i = ", max_coef_i
print "max_coef = ", max_coef
pylab.plot(steps, coefs)
do_one_plot()
#plot_full(0.05, 80)
#get_coef(0.01, 0.2)
pylab.show()
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment