behavior_description = np.random.randint(0, 3, (20, 4))
quality = np.random.random(20)
 
map_size = behavior_description.max(0) - behavior_description.min(0) + 1
 
height_size = 1
width_size = 1
height_idx = []
width_idx = []
for idx, s in enumerate(map_size):
	if height_size > width_size:
    	width_idx.append(idx)
    	width_size *= s
    else:
    	height_idx.append(idx)
		height_size *= s
 
height_ticks = map_size[height_idx]
height = np.product(height_ticks)
width_ticks = map_size[width_idx]
width = np.product(width_ticks)
bd_map = np.zeros((height, width))
 
 bd_idxs = behavior_description - behavior_description.min(0)
 
for bd_idx, qv in zip(bd_idxs, quality):
	y = 0
	for i, hi in enumerate(bd_idx[height_idx]):
    	y += hi * np.product(height_ticks[i+1:])
	x = 0
    for i, wi in enumerate(bd_idx[width_idx]):
    	x += wi * np.product(width_ticks[i+1:])
	bd_map[y][x] = qv
 
fig = plt.figure()
ax = fig.add_subplot(111)
ax.set_xticklabels([])
ax.set_yticklabels([])
n_solutions = len(quality)
search_space_size = np.product(bd_map.shape)
ax.set_xlabel(f'{n_solutions:,} / {search_space_size:,} = {n_solutions/search_space_size:.3f}')
im = ax.imshow(bd_map, cmap='jet', vmin=0, vmax=1)
fig.colorbar(im)
plt.show()
plt.close()