2014-05-03 04:49:35 +02:00
|
|
|
# -*- coding: utf-8 -*-
|
|
|
|
"""
|
|
|
|
Created on Fri May 2 11:48:55 2014
|
|
|
|
|
|
|
|
@author: rlabbe
|
|
|
|
"""
|
|
|
|
from __future__ import division
|
|
|
|
from __future__ import print_function
|
|
|
|
|
|
|
|
import copy
|
|
|
|
import numpy as np
|
|
|
|
import numpy.random as random
|
|
|
|
import matplotlib.pyplot as plt
|
|
|
|
|
|
|
|
''' should this be a class? seems like both sense and update are very
|
|
|
|
problem specific
|
|
|
|
'''
|
|
|
|
|
2015-04-05 21:26:21 +02:00
|
|
|
|
|
|
|
def bar_plot(pos, ylim=(0,1), title=None):
|
|
|
|
plt.cla()
|
|
|
|
ax = plt.gca()
|
|
|
|
x = np.arange(len(pos))
|
|
|
|
ax.bar(x, pos, color='#30a2da')
|
|
|
|
if ylim:
|
|
|
|
plt.ylim(ylim)
|
|
|
|
plt.xticks(x+0.4, x)
|
|
|
|
if title is not None:
|
|
|
|
plt.title(title)
|
|
|
|
|
|
|
|
|
2014-05-03 04:49:35 +02:00
|
|
|
class DiscreteBayes1D(object):
|
|
|
|
|
|
|
|
def __init__(self, world_map, belief=None):
|
|
|
|
self.world_map = copy.deepcopy(world_map)
|
|
|
|
self.N = len(world_map)
|
|
|
|
|
|
|
|
if belief is None:
|
|
|
|
# create belief, make all equally likely
|
|
|
|
self.belief = np.empty(self.N)
|
|
|
|
self.belief.fill (1./self.N)
|
|
|
|
|
|
|
|
else:
|
|
|
|
self.belief = copy.deepcopy(belief)
|
|
|
|
|
|
|
|
# This will be used to temporarily store calculations of the new
|
|
|
|
# belief. 'k' just means 'this iteration'.
|
|
|
|
self.belief_k = np.empty(self.N)
|
|
|
|
|
|
|
|
assert self.belief.shape == self.world_map.shape
|
|
|
|
|
|
|
|
|
|
|
|
def _normalize (self):
|
|
|
|
s = sum(self.belief)
|
|
|
|
self.belief = self.belief/s
|
|
|
|
|
|
|
|
def sense(self, Z, pHit, pMiss):
|
|
|
|
for i in range (self.N):
|
|
|
|
hit = (self.world_map[i] ==Z)
|
|
|
|
self.belief_k[i] = self.belief[i] * (pHit*hit + pMiss*(1-hit))
|
|
|
|
|
|
|
|
# copy results to the belief vector using swap-copy idiom
|
|
|
|
self.belief, self.belief_k = self.belief_k, self.belief
|
|
|
|
self._normalize()
|
|
|
|
|
|
|
|
def update(self, U, kernel):
|
|
|
|
N = self.N
|
|
|
|
kN = len(kernel)
|
|
|
|
width = int((kN - 1) / 2)
|
|
|
|
|
|
|
|
self.belief_k.fill(0)
|
|
|
|
|
|
|
|
for i in range(N):
|
|
|
|
for k in range (kN):
|
|
|
|
index = (i + (width-k)-U) % N
|
|
|
|
#print(i,k,index)
|
|
|
|
self.belief_k[i] += self.belief[index] * kernel[k]
|
|
|
|
|
|
|
|
# copy results to the belief vector using swap-copy idiom
|
|
|
|
self.belief, self.belief_k = self.belief_k, self.belief
|
|
|
|
|
|
|
|
def add_noise (Z, count):
|
|
|
|
n= len(Z)
|
|
|
|
for i in range(count):
|
|
|
|
j = random.randint(0,n)
|
|
|
|
Z[j] = random.randint(0,2)
|
|
|
|
|
|
|
|
|
|
|
|
def animate_three_doors (loops=5):
|
|
|
|
world = np.array([1,0,1,0,0,0,0,1,0,0,0,1,0,0,0,0,0])
|
|
|
|
#world = np.array([1,1,1,1,1])
|
|
|
|
#world = np.array([1,0,1,0,1,0])
|
|
|
|
|
|
|
|
|
|
|
|
f = DiscreteBayes1D(world)
|
|
|
|
|
|
|
|
measurements = np.tile(world,loops)
|
|
|
|
add_noise(measurements, 4)
|
|
|
|
|
|
|
|
for m in measurements:
|
|
|
|
f.sense (m, .8, .2)
|
|
|
|
f.update(1, (.05, .9, .05))
|
|
|
|
|
2015-04-05 21:26:21 +02:00
|
|
|
bar_plot(f.belief)
|
2014-05-03 04:49:35 +02:00
|
|
|
plt.pause(0.01)
|
|
|
|
|
|
|
|
|
|
|
|
def _test_filter():
|
|
|
|
def is_near_equal(a,b):
|
|
|
|
try:
|
|
|
|
assert sum(abs(a-b)) < 0.001
|
|
|
|
except:
|
|
|
|
print(a, b)
|
|
|
|
assert False
|
|
|
|
|
|
|
|
def test_update_1():
|
|
|
|
f = DiscreteBayes1D(np.array([0,0,1,0,0]), np.array([0,0,.8,0,0]))
|
|
|
|
f.update (1, [1])
|
|
|
|
is_near_equal (f.belief, np.array([0,0,0,.8,0]))
|
|
|
|
|
|
|
|
f.update (1, [1])
|
|
|
|
is_near_equal (f.belief, np.array([0,0,0,0,.8]))
|
|
|
|
|
|
|
|
f.update (1, [1])
|
|
|
|
is_near_equal (f.belief, np.array([.8,0,0,0,0]))
|
|
|
|
|
|
|
|
f.update (-1, [1])
|
|
|
|
is_near_equal (f.belief, np.array([0,0,0,0,.8]))
|
|
|
|
|
|
|
|
f.update (2, [1])
|
|
|
|
is_near_equal (f.belief, np.array([0,.8,0,0,0]))
|
|
|
|
|
|
|
|
f.update (5, [1])
|
|
|
|
is_near_equal (f.belief, np.array([0,.8,0,0,0]))
|
|
|
|
|
|
|
|
|
|
|
|
def test_undershoot():
|
|
|
|
f = DiscreteBayes1D(np.array([0,0,1,0,0]), np.array([0,0,.8,0,0]))
|
|
|
|
f.update (2, [.2, .8,0.])
|
|
|
|
is_near_equal (f.belief, np.array([0,0,0,.16,.64]))
|
|
|
|
|
|
|
|
def test_overshoot():
|
|
|
|
f = DiscreteBayes1D(np.array([0,0,1,0,0]), np.array([0,0,.8,0,0]))
|
|
|
|
f.update (2, [0, .8, .2])
|
|
|
|
is_near_equal (f.belief, np.array([.16,0,0,0,.64]))
|
|
|
|
|
|
|
|
|
|
|
|
test_update_1()
|
|
|
|
test_undershoot()
|
|
|
|
|
|
|
|
if __name__ == "__main__":
|
|
|
|
|
|
|
|
_test_filter()
|
|
|
|
|
|
|
|
|
|
|
|
|
2015-04-05 21:26:21 +02:00
|
|
|
animate_three_doors(loops=1)
|
2014-05-03 04:49:35 +02:00
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|