2018-10-03 18:59:01 +00:00
|
|
|
# setup: import numpy ; a = numpy.array([ [i/10., i/10., i/20.] for i in range(44440)], dtype=numpy.double) # noqa
|
|
|
|
# run: hyantes(0, 0, 90, 90, 1, 100, 80, 80, a)
|
2018-04-05 22:07:33 +00:00
|
|
|
|
2018-10-03 18:59:01 +00:00
|
|
|
# pythran export hyantes(float, float, float, float, float, float, int,
|
|
|
|
# int, float[][])
|
2018-04-05 22:07:33 +00:00
|
|
|
import numpy as np
|
2018-10-03 18:59:01 +00:00
|
|
|
|
|
|
|
|
2018-04-05 22:07:33 +00:00
|
|
|
def hyantes(xmin, ymin, xmax, ymax, step, range_, range_x, range_y, t):
|
2018-10-03 18:59:01 +00:00
|
|
|
X, Y = t.shape
|
|
|
|
pt = np.zeros((X, Y))
|
2018-04-05 22:07:33 +00:00
|
|
|
for i in range(X):
|
|
|
|
for j in range(Y):
|
|
|
|
for k in t:
|
2018-10-03 18:59:01 +00:00
|
|
|
tmp = (6368. * np.arccos(np.cos(xmin + step * i)
|
|
|
|
* np.cos(k[0]) * np.cos((ymin + step * j) - k[1])
|
|
|
|
+ np.sin(xmin + step * i) * np.sin(k[0])))
|
2018-04-05 22:07:33 +00:00
|
|
|
if tmp < range_:
|
2018-10-03 18:59:01 +00:00
|
|
|
pt[i, j] += k[2] / (1 + tmp)
|
2018-04-05 22:07:33 +00:00
|
|
|
return pt
|