forked from rougier/ten-rules
-
Notifications
You must be signed in to change notification settings - Fork 0
/
projections.py
122 lines (98 loc) · 4.02 KB
/
projections.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
# -*- coding: utf-8 -*-
# -----------------------------------------------------------------------------
# Copyright INRIA
# Contributors: Wahiba Taouali ([email protected])
# Nicolas P. Rougier ([email protected])
#
# This software is governed by the CeCILL license under French law and abiding
# by the rules of distribution of free software. You can use, modify and/ or
# redistribute the software under the terms of the CeCILL license as circulated
# by CEA, CNRS and INRIA at the following URL
# http://www.cecill.info/index.en.html.
#
# As a counterpart to the access to the source code and rights to copy, modify
# and redistribute granted by the license, users are provided only with a
# limited warranty and the software's author, the holder of the economic
# rights, and the successive licensors have only limited liability.
#
# In this respect, the user's attention is drawn to the risks associated with
# loading, using, modifying and/or developing or reproducing the software by
# the user in light of its specific status of free software, that may mean that
# it is complicated to manipulate, and that also therefore means that it is
# reserved for developers and experienced professionals having in-depth
# computer knowledge. Users are therefore encouraged to load and test the
# software's suitability as regards their requirements in conditions enabling
# the security of their systems and/or data to be ensured and, more generally,
# to use and operate it in the same conditions as regards security.
#
# The fact that you are presently reading this means that you have had
# knowledge of the CeCILL license and that you accept its terms.
# -----------------------------------------------------------------------------
import os
import numpy as np
from parameters import *
def cartesian_to_polar(x, y):
''' Cartesian to polar coordinates. '''
rho = np.sqrt(x**2+y**2)
theta = np.arctan2(y,x)
return rho,theta
def polar_to_cartesian(rho, theta):
''' Polar to cartesian coordinates. '''
x = rho*np.cos(theta)
y = rho*np.sin(theta)
return x,y
def polar_to_logpolar(rho, theta):
''' Polar to logpolar coordinates. '''
# Shift in the SC mapping function in deg
A = 3.0
# Collicular magnification along u axe in mm/rad
Bx = 1.4
# Collicular magnification along v axe in mm/rad
By = 1.8
xmin, xmax = 0.0, 4.80743279742
ymin, ymax = -2.76745559565, 2.76745559565
rho = rho*90.0
x = Bx*np.log(np.sqrt(rho*rho+2*A*rho*np.cos(theta)+A*A)/A)
y = By*np.arctan(rho*np.sin(theta)/(rho*np.cos(theta)+A))
x = (x-xmin)/(xmax-xmin)
y = (y-ymin)/(ymax-ymin)
return x, y
def retina_projection(Rs=retina_shape, Ps=projection_shape):
'''
Compute the projection indices from retina to colliculus
Parameters
----------
Rs : (int,int)
Half-retina shape
Ps : (int,int)
Retina projection shape (might be different from colliculus)
'''
filename = "retina (%d,%d) - colliculus (%d,%d).npy" % (Rs[0],Rs[1],Ps[0],Ps[1])
if os.path.exists(filename):
return np.load(filename)
s = 4
rho = ((np.logspace(start=0, stop=1, num=s*Rs[1],base=10)-1)/9.)
theta = np.linspace(start=-np.pi/2,stop=np.pi/2, num=s*Rs[0])
rho = rho.reshape((s*Rs[1],1))
rho = np.repeat(rho,s*Rs[0], axis=1)
theta = theta.reshape((1,s*Rs[0]))
theta = np.repeat(theta,s*Rs[1], axis=0)
y,x = polar_to_cartesian(rho,theta)
xmin,xmax = x.min(), x.max()
x = (x-xmin)/(xmax-xmin)
ymin,ymax = y.min(), y.max()
y = (y-ymin)/(ymax-ymin)
P = np.zeros((Ps[0],Ps[1],2), dtype=int)
xi = np.rint(x*(Rs[0]-1)).astype(int)
yi = np.rint((0.0+1.0*y)*(Rs[1]-1)).astype(int)
yc,xc = polar_to_logpolar(rho,theta)
xmin,xmax = xc.min(), xc.max()
xc = (xc-xmin)/(xmax-xmin)
ymin,ymax = yc.min(), yc.max()
yc = (yc-ymin)/(ymax-ymin)
xc = np.rint(xc*(Ps[0]-1)).astype(int)
yc = np.rint((.0+yc*1.0)*(Ps[1]-1)).astype(int)
P[xc,yc,0] = xi
P[xc,yc,1] = yi
np.save(filename, P)
return P