appl121105_iotorus.observer¶
Observer related classes.
It could be too much, but there are three classes.
VirtualObserverVirtualSpacecraftVirtualSensor
They are usually used for calculating the fov direction considering ram velocity.
Spacecraft is located at 100000 km in x.
It moves toward y direction with a speed of 30 km/s.
>>> spacecraft = VirtualSpacecraft([100000, 0., 0.], [0, 30, 0])
There is a sensor, looking toward -z, and the detecting
velocity is 100 km/s
>>> sensor = VirtualSensor([0, 0, -1.5], 100.)
Note that the “absolute value” of FOV is meaningless, and inside the class, it is normalized.
>>> print(sensor.fov[2])
-1.0
>>> print(sensor.get_velocity_vector())
[ -0. -0. 100.]
Then, you may instance them into a observer class.
>>> observer = VirtualObserver(spacecraft, sensor)
To obtained the ram-velocity corrected velocity, you can use
VirtualObserver.get_observe_velocity() method.
>>> v_part = observer.get_observe_velocity()
>>> print(v_part[1]) # Vy is 30 km/s, because of ram velocity
30.0
>>> print(v_part[2]) # Vz is 100 km/s, because of particle velocity
100.0
''' Observer related classes.
It could be too much, but there are three classes.
* :class:`VirtualObserver`
* :class:`VirtualSpacecraft`
* :class:`VirtualSensor`
They are usually used for calculating the fov direction
considering ram velocity.
Spacecraft is located at 100000 km in ``x``.
It moves toward ``y`` direction with a speed of 30 km/s.
>>> spacecraft = VirtualSpacecraft([100000, 0., 0.], [0, 30, 0])
There is a sensor, looking toward ``-z``, and the detecting
velocity is 100 km/s
>>> sensor = VirtualSensor([0, 0, -1.5], 100.)
Note that the "absolute value" of FOV is meaningless, and
inside the class, it is normalized.
>>> print(sensor.fov[2])
-1.0
>>> print(sensor.get_velocity_vector())
[ -0. -0. 100.]
Then, you may instance them into a observer class.
>>> observer = VirtualObserver(spacecraft, sensor)
To obtained the ram-velocity corrected velocity, you can use
:meth:`VirtualObserver.get_observe_velocity` method.
>>> v_part = observer.get_observe_velocity()
>>> print(v_part[1]) # Vy is 30 km/s, because of ram velocity
30.0
>>> print(v_part[2]) # Vz is 100 km/s, because of particle velocity
100.0
'''
import numpy as np
class VirtualObserver:
''' Virtual observer.
Observer object.
'''
def __init__(self, spacecraft, sensor):
''' Setup the observer.
:param spacecraft: :class:`VirtualSpacecraft` object.
:param sensor: :class:`VirtualSensor` object. '''
self.spacecraft = spacecraft
self.sensor = sensor
def get_observe_velocity(self):
r''' Return the velocity vector of the observing ENA in JEqS system
:returns: np.array with (3,) shape. Velocity of ENA in JEqS system.
'''
v_noram = self.sensor.get_velocity_vector()
return v_noram + self.spacecraft.vel
class VirtualSpacecraft:
''' Virtual spacecraft. Attributes of position, velocity, and attitude.
'''
def __init__(self, pos, vel, primary_axis=None, secondary_axis=None):
''' Set up the spacecraft.
The unit of pos and vel should be in km and km/s.
primary_axis and secondary_axis define the angles.
They are normalized first so that uses do not need to normalize,
but they should be non-zero and not be parallel.
third_axis is first determined by the cross product of primary and secondary axes.
secondary_axis is now re-defined from the cross product of third and first axes.
primary_axis is 0 degrees in elevation toward 180 degrees.
secondary_axis is 0 degrees in azimutha and third_axis is 90 degrees in azimuth.
Both secondary_axis and third_axis should have 90 degrees in elevation angle.
If they are not specified, the primary_axis is taken as -pos and secondary_axis is vel.
Note that either secondary_axis or velocity of spacecraft must be specified.
:param pos: Position of observer, in km. (3,) shaped.
:param vel: Velocity of observer, in km/s. (3,) shaped.
:keyword primary_axis:
'''
self.pos = np.array(pos)
self.vel = np.array(vel)
self._setup_attitude(primary_axis, secondary_axis)
def _setup_attitude(self, primary_axis, secondary_axis):
if primary_axis == None:
primary_axis = -np.array(self.pos)
if secondary_axis == None:
secondary_axis = np.array(self.vel)
primary_axis = np.array(primary_axis)
secondary_axis = np.array(secondary_axis)
lp = np.sqrt((primary_axis * primary_axis).sum())
ls = np.sqrt((secondary_axis * secondary_axis).sum())
if lp == 0:
raise ValueError('Primary axis should not be zero. (%g %g %g)' % (primary_axis[0], primary_axis[1], primary_axis[2]))
if ls == 0:
raise ValueError('Secondary axis should not be zero (if velocity of spacecraft is 0, secondary_axis should be given explicitly). (%g %g %g)' % (secondary_axis[0], secondary_axis[1], secondary_axis[2]))
pax = primary_axis / lp
sax = secondary_axis / ls
tax = np.cross(pax, sax)
tax = tax / np.sqrt((tax * tax).sum())
sax = np.cross(tax, pax)
self.primary_axis = pax
self.secondary_axis = sax
self.third_axis = tax
def get_angles(self, kvector):
''' Return the angle corresponding to kvector in JEqS system.
:returns: (theta, phi) in radians.
See :meth:`get_vector` for usage.
'''
kvector = np.array(kvector)
n_kvector = kvector / np.sqrt((kvector ** 2).sum())
theta = np.arccos(n_kvector.dot(self.primary_axis))
phi = -np.arctan2(n_kvector.dot(self.third_axis), n_kvector.dot(self.secondary_axis))
# Minus is to rotate in a natural way (left-hand system). Just consider
# if you are the spacecraft looking toward primary_axis, and secondary_axis to the right.
# Natural roation (toward up) to be positive phi, the minus is needed.
return theta, phi
def get_vector(self, theta, phi):
''' Returns the vector corresponding to elevation and azimuth angle.
:param theta: Elevation angle in radian.
:param phi: Aximuth angle in radian.
:returns: (3,) shaped np.array as a vector in JEqS frame corresponding to the given angles.
Assume the spacecraft at (100, 0, 0) with velocity (0, 3, 0).
In this case, the default primary axis is (-1, 0, 0) and
the secondary axis (0, 1, 0).
The third axis is derived to be (0, 0, -1).
Theta is the angle from (-1, 0, 0).
Phi is the angle from (0, 1, 0) toward -(0, 0, -1) [minus here!!].
>>> sc = VirtualSpacecraft([100, 0, 0], [0, 3, 0])
You may get the vector theta=0, phi=0. Expectation is [-1, 0, 0].
>>> print(sc.get_vector(0, 0))
[-1. 0. 0.]
Reverting is also available.
>>> print(sc.get_angles([-1., 0, 0.]))
(0.0, -0)
If theta = 180 deg, the vector should be [1, 0, 0].
>>> v = sc.get_vector(np.pi, 0); print('%.2f %.2f %.2f' % (v[0], v[1], v[2]))
1.00 0.00 0.00
>>> t, p = np.rad2deg(sc.get_angles([1, 0, 0])); print('%.2f %.2f' % (t, p))
180.00 -0.00
If theta=90, x is always 0. For azimuth=0, it returns secondary_axis, [0, 1, 0].
>>> v = sc.get_vector(np.pi / 2., 0); print('%.2f %.2f %.2f' % (v[0], v[1], v[2]))
-0.00 1.00 0.00
>>> t, p = np.rad2deg(sc.get_angles([0, 1, 0])); print('%.2f %.2f' % (t, p))
90.00 -0.00
For azimuth 90, it returns minus of third_axis, i.e. [0, 0, 1]
>>> v = sc.get_vector(np.pi / 2., np.pi / 2.); print('%.2f %.2f %.2f' % (v[0], v[1], v[2]))
-0.00 0.00 1.00
>>> t, p = np.rad2deg(sc.get_angles([0, 0, 1])); print('%.2f %.2f' % (t, p))
90.00 90.00
You can give specific primary and secondary axes.
>>> sc = VirtualSpacecraft([100, 0, 0], [0, 3, 0], secondary_axis=[-100, 3, 0])
>>> v = sc.get_vector(np.pi / 2., np.pi / 2.); print('%.2f %.2f %.2f' % (v[0], v[1], v[2]))
-0.00 0.00 1.00
>>> t, p = np.rad2deg(sc.get_angles([0, 0, 1])); print('%.2f %.2f' % (t, p))
90.00 90.00
'''
p = np.cos(theta)
s = np.sin(theta) * np.cos(-phi)
t = np.sin(theta) * np.sin(-phi)
return self.primary_axis * p + self.secondary_axis * s + self.third_axis * t
class VirtualSensor:
''' Virtual sensor. Attribute of FOV vector and velocity. Spacecraft coordinate system.
'''
def __init__(self, fov, velocity):
''' Set up the sensor.
An optimal sensor that has fov vector and velocity.
:param fov: Field of view direction, in (3,) shape.
This will be normalized.
:param velocity: Velocity to observe in km/s.
'''
fov = np.array(fov)
lenfov = np.sqrt((fov * fov).sum())
self.fov = fov / lenfov
self.velocity = velocity
def get_velocity_vector(self):
''' Return a velocity vector that sensor detects in S/C frame.
:param velocity: Scalar. Positive (usually). Unit in km/s.
:returns: Velocity vector, (3,) shape.
``-velocity X fov``.
'''
return -self.velocity * self.fov
import unittest
import doctest
def doctests():
return unittest.TestSuite((
doctest.DocTestSuite(),
))
if __name__ == '__main__':
unittest.main(defaultTest='doctests')