def check_visibility(camera, pts_w, zmin=0, zmax=100):
"""
Check if points are visible given fov of camera.
This method checks for both horizontal and vertical
fov.
camera: type Camera
"""
# Transform points in to camera's reference
# Camera: p_cw
pts_c = camera.c2w(pts_w.reshape(-1, 3))
# Determine look-at vector, and check angle
# subtended with camera's z-vector (3rd column)
z = pts_c[:,2]
v = pts_c / np.linalg.norm(pts_c, axis=1).reshape(-1, 1)
hangle, vangle = np.arctan2(v[:,0], v[:,2]), np.arctan2(-v[:,1], v[:,2])
# Provides inds mask for all points that are within fov
return np.fabs(hangle) < camera.fov[0] * 0.5 and \
np.fabs(vangle) < camera.fov[1] * 0.5 and \
z >= zmin and z <= zmax
python类arctan2()的实例源码
def setFromQTransform(self, tr):
p1 = Point(tr.map(0., 0.))
p2 = Point(tr.map(1., 0.))
p3 = Point(tr.map(0., 1.))
dp2 = Point(p2-p1)
dp3 = Point(p3-p1)
## detect flipped axes
if dp2.angle(dp3) > 0:
#da = 180
da = 0
sy = -1.0
else:
da = 0
sy = 1.0
self._state = {
'pos': Point(p1),
'scale': Point(dp2.length(), dp3.length() * sy),
'angle': (np.arctan2(dp2[1], dp2[0]) * 180. / np.pi) + da
}
self.update()
def setFromQTransform(self, tr):
p1 = Point(tr.map(0., 0.))
p2 = Point(tr.map(1., 0.))
p3 = Point(tr.map(0., 1.))
dp2 = Point(p2-p1)
dp3 = Point(p3-p1)
## detect flipped axes
if dp2.angle(dp3) > 0:
#da = 180
da = 0
sy = -1.0
else:
da = 0
sy = 1.0
self._state = {
'pos': Point(p1),
'scale': Point(dp2.length(), dp3.length() * sy),
'angle': (np.arctan2(dp2[1], dp2[0]) * 180. / np.pi) + da
}
self.update()
def great_circle_dist(p1, p2):
"""Return the distance (in km) between two points in
geographical coordinates.
"""
lon0, lat0 = p1
lon1, lat1 = p2
EARTH_R = 6372.8
lat0 = np.radians(float(lat0))
lon0 = np.radians(float(lon0))
lat1 = np.radians(float(lat1))
lon1 = np.radians(float(lon1))
dlon = lon0 - lon1
y = np.sqrt(
(np.cos(lat1) * np.sin(dlon)) ** 2
+ (np.cos(lat0) * np.sin(lat1)
- np.sin(lat0) * np.cos(lat1) * np.cos(dlon)) ** 2)
x = np.sin(lat0) * np.sin(lat1) + \
np.cos(lat0) * np.cos(lat1) * np.cos(dlon)
c = np.arctan2(y, x)
return EARTH_R * c
def get_phases(self):
sizeimg = np.real(self.imgfft).shape
mag = np.zeros(sizeimg)
for x in range(sizeimg[0]):
for y in range(sizeimg[1]):
mag[x][y] = np.arctan2(np.real(self.imgfft[x][y]), np.imag(self.imgfft[x][y]))
rpic = MyImage(mag)
rpic.limit(1)
return rpic
# int my = y-output.height/2;
# int mx = x-output.width/2;
# float angle = atan2(my, mx) - HALF_PI ;
# float radius = sqrt(mx*mx+my*my) / factor;
# float ix = map(angle,-PI,PI,input.width,0);
# float iy = map(radius,0,height,0,input.height);
# int inputIndex = int(ix) + int(iy) * input.width;
# int outputIndex = x + y * output.width;
# if (inputIndex <= input.pixels.length-1) {
# output.pixels[outputIndex] = input.pixels[inputIndex];
def circles(self, x1, y1, x2, y2, r, nmin=2):
arc = self.ctx.arc
fill = self.ctx.fill
dx = x1-x2
dy = y1-y2
dd = sqrt(dx*dx+dy*dy)
n = int(dd/self.pix)
n = n if n>nmin else nmin
a = arctan2(dy, dx)
scale = linspace(0, dd, n)
xp = x1-scale*cos(a)
yp = y1-scale*sin(a)
for x, y in zip(xp, yp):
arc(x, y, r, 0, pi*2.)
fill()
def sandstroke_non_linear(self,xys,grains=10,left=True):
pix = self.pix
rectangle = self.ctx.rectangle
fill = self.ctx.fill
dx = xys[:,2] - xys[:,0]
dy = xys[:,3] - xys[:,1]
aa = arctan2(dy,dx)
directions = column_stack([cos(aa),sin(aa)])
dd = sqrt(square(dx)+square(dy))
for i,d in enumerate(dd):
rnd = sqrt(random((grains,1)))
if left:
rnd = 1.0-rnd
for x,y in xys[i,:2] + directions[i,:]*rnd*d:
rectangle(x,y,pix,pix)
fill()
def sandstroke(self,xys,grains=10):
pix = self.pix
rectangle = self.ctx.rectangle
fill = self.ctx.fill
dx = xys[:,2] - xys[:,0]
dy = xys[:,3] - xys[:,1]
aa = arctan2(dy,dx)
directions = column_stack([cos(aa),sin(aa)])
dd = sqrt(square(dx)+square(dy))
for i,d in enumerate(dd):
for x,y in xys[i,:2] + directions[i,:]*random((grains,1))*d:
rectangle(x,y,pix,pix)
fill()
def dir_threshold(img, sobel_kernel=3, thresh=(0, np.pi/2)):
# Apply the following steps to img
# 1) Convert to grayscale
gray = cv2.cvtColor(img, cv2.COLOR_RGB2GRAY)
# 2) Take the gradient in x and y separately
sobelx = cv2.Sobel(gray, cv2.CV_64F, 1, 0, ksize=sobel_kernel)
sobely = cv2.Sobel(gray, cv2.CV_64F, 0, 1, ksize=sobel_kernel)
# 3) Take the absolute value of the x and y gradients
abs_sobelx = np.absolute(sobelx)
abs_sobely = np.absolute(sobely)
# 4) Use np.arctan2(abs_sobely, abs_sobelx) to calculate the direction of the gradient
absgraddir = np.arctan2(abs_sobely, abs_sobelx)
# 5) Create a binary mask where direction thresholds are met
binary_output = np.zeros_like(absgraddir)
binary_output[(absgraddir >= thresh[0]) & (absgraddir <= thresh[1])] = 1
# 6) Return this mask as your binary_output image
return binary_output
# Define a function that applies Sobel x and y,
# then computes the magnitude of the gradient
# and applies a threshold
def make_spiral(self, r=0.25, G=0.0001):
for k in range(10):
x = self.X[:, 0] - 0.5
y = self.X[:, 1] - 0.5
theta = np.arctan2(x, y)
ds = [r * (i + theta / (2 * np.pi)) for i in range(int(1 / r))]
alphas = [np.sqrt(x ** 2 + y ** 2) / d for d in ds]
for alpha in alphas:
d = np.concatenate([(x * (1 - alpha))[:, None], (y * (1 - alpha))[:, None]], axis=1)
f = -G * d / (d ** 2).sum(axis=1, keepdims=True)
self.X += f
self.X = np.clip(self.X, 0, 1)
rs = np.arange(0, 0.7, 0.001)
theta = 2 * np.pi * rs / r
y = rs * np.sin(theta) + 0.5
x = -rs * np.cos(theta) + 0.5
spiral = zip(x, y)
self.collection = matplotlib.collections.LineCollection([spiral], colors='k')
def ecc(self, val):
'''
We need to update the time of pericenter passage whenever the
eccentricty, longitude of pericenter, period, or time of transit
changes. See the appendix in Shields et al. (2016).
'''
self._ecc = val
fi = (3 * np.pi / 2.) - self._w
self.tperi0 = self.t0 + (self.per * np.sqrt(1. - self.ecc * self.ecc) /
(2. * np.pi) * (self.ecc * np.sin(fi) /
(1. + self.ecc * np.cos(fi)) - 2.
/ np.sqrt(1. - self.ecc * self.ecc)
* np.arctan2(np.sqrt(1. - self.ecc * self.ecc)
* np.tan(fi/2.), 1. + self.ecc)))
def preprocess_edges(self):
# Calculate the sobel edge features
denom = 3 * 255.
grey = np.sum(self.image/denom, axis=2, keepdims=False, dtype=np.float32)
sobel_x = scipy.ndimage.sobel(grey, axis=0)
sobel_y = scipy.ndimage.sobel(grey, axis=1)
self.edge_angles = np.arctan2(sobel_y, sobel_x) # Counterclockwise
self.edge_magnitues = (sobel_x ** 2 + sobel_y ** 2) ** .5
assert(self.edge_angles.dtype == np.float32)
assert(self.edge_magnitues.dtype == np.float32)
if False:
plt.figure("EDGES")
plt.subplot(1,2,1)
plt.imshow(self.edge_magnitues, interpolation='nearest')
plt.title("MAG")
plt.subplot(1,2,2)
plt.imshow(self.edge_angles, interpolation='nearest')
plt.title("ANG")
plt.show()
def azimuth(_lAz_data):
_inc = _lAz_data[0]
_lat = _lAz_data[1]
velocity_eq = _lAz_data[2]
@jit(nopython=True)
def _az_calc():
inert_az = np.arcsin(max(min(np.cos(np.deg2rad(_inc)) / np.cos(np.deg2rad(_lat)), 1), -1))
_VXRot = _lAz_data[3] * np.sin(inert_az) - velocity_eq * np.cos(np.deg2rad(_lat))
_VYRot = _lAz_data[3] * np.cos(inert_az)
return np.rad2deg(np.fmod(np.arctan2(_VXRot, _VYRot) + (2 * pi), (2 * pi)))
_az = _az_calc()
if _lAz_data[4] == "Ascending": return _az
if _lAz_data[4] == "Descending":
if _az <= 90: return 180 - _az
elif _az >= 270: return 540 - _az
def azimuth(_lAz_data):
_inc = _lAz_data[0]
_lat = _lAz_data[1]
velocity_eq = _lAz_data[2]
@jit(nopython=True)
def _az_calc():
inert_az = np.arcsin(max(min(np.cos(np.deg2rad(_inc)) / np.cos(np.deg2rad(_lat)), 1), -1))
_VXRot = _lAz_data[3] * np.sin(inert_az) - velocity_eq * np.cos(np.deg2rad(_lat))
_VYRot = _lAz_data[3] * np.cos(inert_az)
return np.rad2deg(np.fmod(np.arctan2(_VXRot, _VYRot) + (2 * pi), (2 * pi)))
_az = _az_calc()
if _lAz_data[4] == "Ascending": return _az
if _lAz_data[4] == "Descending":
if _az <= 90: return 180 - _az
elif _az >= 270: return 540 - _az
def azimuth(_lAz_data):
_inc = _lAz_data[0]
_lat = _lAz_data[1]
velocity_eq = _lAz_data[2]
@jit(nopython=True)
def _az_calc():
inert_az = np.arcsin(max(min(np.cos(np.deg2rad(_inc)) / np.cos(np.deg2rad(_lat)), 1), -1))
_VXRot = _lAz_data[3] * np.sin(inert_az) - velocity_eq * np.cos(np.deg2rad(_lat))
_VYRot = _lAz_data[3] * np.cos(inert_az)
return np.rad2deg(np.fmod(np.arctan2(_VXRot, _VYRot) + (2 * pi), (2 * pi)))
_az = _az_calc()
if _lAz_data[4] == "Ascending": return _az
if _lAz_data[4] == "Descending":
if _az <= 90: return 180 - _az
elif _az >= 270: return 540 - _az
def azimuth(_lAz_data):
_inc = _lAz_data[0]
_lat = _lAz_data[1]
velocity_eq = _lAz_data[2]
@jit(nopython=True)
def _az_calc():
inert_az = np.arcsin(max(min(np.cos(np.deg2rad(_inc)) / np.cos(np.deg2rad(_lat)), 1), -1))
_VXRot = _lAz_data[3] * np.sin(inert_az) - velocity_eq * np.cos(np.deg2rad(_lat))
_VYRot = _lAz_data[3] * np.cos(inert_az)
return np.rad2deg(np.fmod(np.arctan2(_VXRot, _VYRot) + 360, 360))
_az = _az_calc()
if _lAz_data[4] == "Ascending": return _az
if _lAz_data[4] == "Descending":
if _az <= 90: return 180 - _az
elif _az >= 270: return 540 - _az
def asSpherical(self):
'''
Computes and returns a representation of this point in spherical coordinates: (r,phi,theta).
r = radius or distance of the point from the origin.
phi = is the angle of the projection on the xy plain and the x axis
theta = is the angle with the z axis.
x = r*cos(phi)*sin(theta)
y = r*sin(phi)*sin(theta)
z = r*cos(theta)
'''
x,y,z,_ = self.asArray()
r = np.sqrt(x**2+y**2+z**2)
phi = np.arctan2(y,x)
theta = np.arctan2(np.sqrt(x**2+y**2),z)
return r,phi,theta
def __init__(self,jet,kernels,k,x,y,pt,subpixel):
self.jet = jet
self.kernels = kernels
self.k = k
self.x = x
self.y = y
re = np.real(jet)
im = np.imag(jet)
self.mag = np.sqrt(re*re + im*im)
self.phase = np.arctan2(re,im)
if subpixel:
d = np.array([[pt.X()-x],[pt.Y()-y]])
comp = np.dot(self.k,d)
self.phase -= comp.flatten()
self.jet = self.mag*np.exp(1.0j*self.phase)
def vector_angle_in_degrees(v):
'''
Given a vector, returns its elevation angle in degrees (-180..180).
>>> vector_angle_in_degrees([1, 0])
0.0
>>> vector_angle_in_degrees([1, 1])
45.0
>>> vector_angle_in_degrees([0, 1])
90.0
>>> vector_angle_in_degrees([-1, 1])
135.0
>>> vector_angle_in_degrees([-1, 0])
180.0
>>> vector_angle_in_degrees([-1, -1])
-135.0
>>> vector_angle_in_degrees([0, -1])
-90.0
>>> vector_angle_in_degrees([1, -1])
-45.0
'''
return np.arctan2(v[1], v[0]) * 180 / np.pi
def calcAngSepDeg(ra0, dec0, ra1, dec1):
'''Return the angular separation between two objects. Use the
special case of the Vincenty formula that is accurate for all
distances'''
C = numpy.pi / 180
d0 = C * dec0
d1 = C * dec1
r12 = C * (ra0 - ra1)
cd0 = numpy.cos(d0)
sd0 = numpy.sin(d0)
cd1 = numpy.cos(d1)
sd1 = numpy.sin(d1)
cr12 = numpy.cos(r12)
sr12 = numpy.sin(r12)
num = numpy.sqrt((cd0 * sr12) ** 2 + (cd1 * sd0 - sd1 * cd0 * cr12) ** 2)
den = sd0 * sd1 + cd0 * cd1 * cr12
return numpy.arctan2(num, den) / C
def nodeEpochsCalc(paramsDI,omegaDIoffset):
"""
Calculate the epochs for the Ascending and Descending nodes, might be in different orbital
periods and AN/DN might be the wrong order, but should work for plotting... I hope...
"""
taAtNodes = [-1.0*paramsDI[9]+omegaDIoffset,180.0-1.0*paramsDI[9]+omegaDIoffset]
nodeEpochs = []
for ta in taAtNodes:
if ta<0.0:
ta =ta+360.0
elif ta>360:
ta =ta-360.0
TA_s_rad = ta*(np.pi/180.0)
top = np.sqrt(1.0-paramsDI[4])*np.sin(TA_s_rad/2.0)
btm = np.sqrt(1.0+paramsDI[4])*np.cos(TA_s_rad/2.0)
ATAN_rad = np.arctan2(top, btm)
#NOTE: both math.atan2 and np.arctan2 tried with same results, both produce negatives rather than continuous 0-360
#thus, must correct for negative outputs
if ATAN_rad<0:
ATAN_rad = ATAN_rad+(2.0*np.pi)
M_s_rad = ATAN_rad*2.0-paramsDI[4]*np.sin(ATAN_rad*2.0)
delta_t = (M_s_rad*paramsDI[7]*days_per_year)/(2.0*np.pi)
nodeEpochs.append(paramsDI[5]+delta_t)
return nodeEpochs
def polar3d(x):
"""
Polar coordinate representation of a three-dimensional vector.
@param x: vector (i.e. rank one array)
@return: polar coordinates (radius and polar angles)
"""
if x.shape != (3,):
raise ValueError(x)
r = norm(x)
theta = numpy.arccos(x[2] / r)
phi = numpy.arctan2(x[1], x[0])
return numpy.array([r, theta, phi])
dtopotools_horiz_okada_and_1d.py 文件源码
项目：finite_volume_seismic_model
作者: cjvogl
项目源码
文件源码
阅读 28
收藏 0
点赞 0
评论 0
def strike_direction(x1, y1, x2, y2):
"""
Calculate strike direction between two points.
Actually calculates "initial bearing" from (x1,y1) in direction
towards (x2,y2), following
http://www.movable-type.co.uk/scripts/latlong.html
"""
x1 = x1*numpy.pi/180.
y1 = y1*numpy.pi/180.
x2 = x2*numpy.pi/180.
y2 = y2*numpy.pi/180.
dx = x2-x1
theta = numpy.arctan2(numpy.sin(dx)*numpy.cos(y2), \
numpy.cos(y1)*numpy.sin(y2) \
- numpy.sin(y1)*numpy.cos(y2)*numpy.cos(dx))
s = theta*180./numpy.pi
if s<0:
s = 360+s
return s
dtopotools_horiz_okada_and_1d.py 文件源码
项目：finite_volume_seismic_model
作者: cjvogl
项目源码
文件源码
阅读 27
收藏 0
点赞 0
评论 0
def strike_direction(x1, y1, x2, y2):
"""
Calculate strike direction between two points.
Actually calculates "initial bearing" from (x1,y1) in direction
towards (x2,y2), following
http://www.movable-type.co.uk/scripts/latlong.html
"""
x1 = x1*numpy.pi/180.
y1 = y1*numpy.pi/180.
x2 = x2*numpy.pi/180.
y2 = y2*numpy.pi/180.
dx = x2-x1
theta = numpy.arctan2(numpy.sin(dx)*numpy.cos(y2), \
numpy.cos(y1)*numpy.sin(y2) \
- numpy.sin(y1)*numpy.cos(y2)*numpy.cos(dx))
s = theta*180./numpy.pi
if s<0:
s = 360+s
return s
def value(self, xyz):
xyz = xyz.reshape(-1,3)
a = self.a
b = self.b
c = self.c
d = self.d
vec1 = xyz[b] - xyz[a]
vec2 = xyz[c] - xyz[b]
vec3 = xyz[d] - xyz[c]
cross1 = np.cross(vec2, vec3)
cross2 = np.cross(vec1, vec2)
arg1 = np.sum(np.multiply(vec1, cross1)) * \
np.sqrt(np.sum(vec2**2))
arg2 = np.sum(np.multiply(cross1, cross2))
answer = np.arctan2(arg1, arg2)
return answer
def value(self, xyz):
xyz = xyz.reshape(-1,3)
a = np.array(self.a)
b = self.b
c = self.c
d = np.array(self.d)
xyza = np.mean(xyz[a], axis=0)
xyzd = np.mean(xyz[d], axis=0)
vec1 = xyz[b] - xyza
vec2 = xyz[c] - xyz[b]
vec3 = xyzd - xyz[c]
cross1 = np.cross(vec2, vec3)
cross2 = np.cross(vec1, vec2)
arg1 = np.sum(np.multiply(vec1, cross1)) * \
np.sqrt(np.sum(vec2**2))
arg2 = np.sum(np.multiply(cross1, cross2))
answer = np.arctan2(arg1, arg2)
return answer
def value(self, xyz):
xyz = xyz.reshape(-1,3)
a = self.a
b = self.b
c = self.c
d = self.d
vec1 = xyz[b] - xyz[a]
vec2 = xyz[c] - xyz[b]
vec3 = xyz[d] - xyz[c]
cross1 = np.cross(vec2, vec3)
cross2 = np.cross(vec1, vec2)
arg1 = np.sum(np.multiply(vec1, cross1)) * \
np.sqrt(np.sum(vec2**2))
arg2 = np.sum(np.multiply(cross1, cross2))
answer = np.arctan2(arg1, arg2)
return answer
def q_to_euler(q):
"""Converts Quaternions to Euler angles.
Parameters
----------
q : array_like
Array holding Quaternions.
Returns
-------
phi : float
`phi` angle in radians.
theta :float
`theta` angle in radians.
psi : float
`psi` angle in radians.
"""
phi = np.arctan2(2*(q[0]*q[1]+q[2]*q[3]),(q[0]**2+q[3]**2-q[1]**2-q[2]**2))
theta = np.arcsin(2*(q[0]*q[2]-q[1]*q[3]))
psi = np.arctan2(2*(q[0]*q[3]+q[1]*q[2]),(q[0]**2+q[1]**2-q[2]**2-q[3]**2))
return phi, theta, psi
def _circle_reference(state,
time,
finished,
radius=None,
speed=None,
init_angle=None,
z_vel=None):
ref = StateVector()
angle = init_angle + speed / radius * time
ref.pos = array([radius * cos(angle),
radius * sin(angle),
z_vel * time])
ref.vel[:] = [-speed * sin(angle), speed * cos(angle), z_vel]
ref.euler[2] = pi + np.arctan2(state.pos[1], state.pos[0])
# reference.omega_b[2] = speed / radius
return ref
# private stationary reference function
def erf(xgoal, x):
"""
Returns error e given two states xgoal and x.
Angle differences are taken properly on SO3.
"""
e = xgoal - x
c = np.cos(x[2])
s = np.sin(x[2])
cg = np.cos(xgoal[2])
sg = np.sin(xgoal[2])
e[2] = np.arctan2(sg*c - cg*s, cg*c + sg*s)
return e
################################################# OBJECTIVES
# Initial condition and goal
