Source code for ecpi.common.sky.extented_body
'''
Manage Earth, moon, sun in field of view
@author: Colley Jean-Marc
'''
import logging
import numpy as np
from astropy.constants import R_earth
from astropy.constants import R_sun
import transforms3d.quaternions as tq
from scipy.linalg import norm
from ecpi.common.mission.attitude import xyz_to_direlev, direlev_to_xyz
from ecpi.common.instru.model_geom import InstruECLAIRs
from ecpi.common.mission.attitude import AttitudeECLAIRs
logger = logging.getLogger(__name__)
[docs]class ExtentedBodyInFov(object):
"""Informs if an arbitrary astronomical body is in the
filed of view of the ECLAIRs instrument.
"""
def __init__(self):
"""**Constructor**
"""
self._name = "TBD"
self._radius = None # km
self._position_body = None
self._instru_ecl = InstruECLAIRs()
self._ecl_att = AttitudeECLAIRs()
self._ecl_att.set_attitude_svom_quater([1.0, 0, 0, 0])
self._position = np.array([1, 0, 0])
self.omega_half = self._instru_ecl.get_half_fov()
[docs] def set_quater_svom(self, quater):
"""Set SVOM quaternion
..warning: set position satellite before
:param quater: unit quaternion
:type quater: array(float) of length 4
"""
self._ecl_att = AttitudeECLAIRs()
self._ecl_att.set_attitude_svom_quater(quater)
self._update()
[docs] def set_eclairs_attitude(self, ecl_att):
"""Set ECLAIRs attitude object
..warning: set position satellite before
:param ecl_att: ECLAIRs instrument attitude object
:type ecl_att: array(float) of length 3
"""
self._ecl_att = ecl_att
self._update()
[docs] def set_pos_sat_time(self, c_pos, c_time=None):
"""Set satellite position at considered time.
:param c_pos: satellite position [X, Y, Z] in km
:type c_pos: array(float) of length 3
:param c_time: considered time in s. Default is None.
:type c_time: float
"""
self._time = c_time
self._position = c_pos
self._update()
def _update(self):
self._compute_position_body()
self._compute_limb_body()
self._compute_pos_fov()
def _compute_position_body(self):
"""Compute satellite altitude in km.
.. note:: must be redefined in herit class with def of self._position_body
"""
self._dist_sat_body = norm(self._position_body)
def _compute_limb_body(self):
"""Compute limb angle in deg.
"""
self._limb = np.rad2deg(np.arcsin(self._radius / self._dist_sat_body))
[docs] def get_fov_direlev(self):
return self._body_direlev
[docs] def get_fov_xyz(self):
return self._body_fov_xyz
def _compute_pos_fov(self):
"""Compute Earth position in term of direction and elevation
(both in deg) in the field of view.
.. note :: apply reflection on Y axis so as to be consistent with
alG/UGTS. TODO: PB check if this also complies with CNES conventions.
"""
quater = self._ecl_att._q_j2000_instf
sat_body_unit = self._position_body / self._dist_sat_body
self._body_fov_xyz = tq.rotate_vector(sat_body_unit, tq.qconjugate(quater))
# reflect Y axis
self._body_fov_xyz[1] = -self._body_fov_xyz[1]
# convert from sat gcrs position to (theta, phi)
self._body_direlev = xyz_to_direlev(self._body_fov_xyz[np.newaxis, :], deg=True)[0]
direlev = self._body_direlev
logger.info(f'dir Earth={direlev[0]} deg | elev earth={direlev[1]} deg')
[docs] def limb_angle(self):
"""
:return: limb angle in degrees
:rtype: float
"""
return self._limb
[docs] def is_in_fov(self):
"""Return whether Earth appears in ECALIRs field of view
:return: if Earth is in FOV
:rtype: bool
"""
logger.info(f'omega half: {self.omega_half} deg')
limb_angle = self._limb
colat = 90 - self._body_direlev[1]
cond = not (colat > limb_angle + self.omega_half)
logger.info(f'{self._name} in FOV: {cond}')
return cond
[docs] def open_fov_map(self, n_size=199):
"""Create the open field of view map, ie. a detector size array
(199x199) filled with 0 and 1. Each pixel value is 0 if Earth is
present in the LoS and 1 otherwise.
:param idx_row: row index in astropy table
:type idx_row: int
:return: open fov 199x199 map
:rtype: 2D array(bool)
"""
size_fov_map = n_size
instru_ecl = self._instru_ecl
open_fov_m = np.empty(
shape=(size_fov_map, size_fov_map),
dtype=np.int16
)
limb_rad = np.deg2rad(self._limb)
logger.debug(f'limb angle: {self._limb} deg ({limb_rad} rad)')
vec_body = self._body_fov_xyz
logger.debug(f'u_body: {vec_body}')
for idx_y in range(size_fov_map):
for idx_z in range(size_fov_map):
# compute corresponding direction of a detector pixel
elev_ij, dir_ij = instru_ecl.skypix_to_elevdir(
(size_fov_map // 2) - idx_y, (size_fov_map // 2) - idx_z
)
ptg_ij = direlev_to_xyz(
np.array([dir_ij + 180, elev_ij]),
deg=True
)
# check condition and fill open fov map
# no need for normalization as direlev_to_xyz returns unit vector.
theta_s = np.arccos(np.dot(ptg_ij, vec_body))
bool_cond = theta_s > limb_rad
open_fov_m[idx_y, idx_z] = bool_cond
return np.flipud(open_fov_m)
#TODO: compute_body_frac plutot
[docs] def compute_earth_frac(self, open_fov_map, size_fov_map=199):
"""Compute earth fraction in FOV;
:param open_fov_map: sky map with directions occulted by Earth
:type open_fov_map: array(bool)
:param size_fov_map: siwe of sky map. Default is 199.
:type open_fov_map: int
:return: earth fraction
:rtype: float
"""
earth_frac = 1 - np.round(np.sum(open_fov_map) / size_fov_map**2, 3)
logger.info(f'fraction of fov covered by Earth: {earth_frac}')
return np.round(earth_frac, 3)
[docs]class EarthInFov(ExtentedBodyInFov):
"""Informs if Earth is in the
filed of view of the ECLAIRs instrument.
Inherits from ExtendedBodyinFov class.
"""
def __init__(self):
"""**Constructor**
"""
super().__init__()
self._name = "earth"
self._radius = 1e-3 * R_earth.value # km
self._position = np.array([self._radius, 0, 0])
def _compute_position_body(self):
'''
Must be redefined in herit class with def of self._position_body
'''
self._position_body = -self._position
super()._compute_position_body()
[docs]class MoonInFov(ExtentedBodyInFov):
"""Informs if Moon is in the
filed of view of the ECLAIRs instrument.
Inherits from ExtendedBodyinFov class.
"""
def __init__(self):
"""
**Constructor**
.. note :: Moon radius note implemented in astropy.constants
"""
super().__init__()
self._name = "moon"
self._radius = 1737.0 # km
self._position = np.array([self._radius, 0, 0])
def _compute_position_body(self):
'''
Must be redefined in herit class with def of self._position_body
'''
# TODO: position de la lune a t dans GCRS
self._position_body = "TDB"
super()._compute_position_body()
[docs]class SunInFov(ExtentedBodyInFov):
"""Informs if Sun is in the
filed of view of the ECLAIRs instrument.
Inherits from ExtendedBodyinFov class.
"""
def __init__(self):
"""
**Constructor**
"""
super().__init__()
self._name = "sun"
self._radius = 1e-3 * R_sun.value # km
self._position = np.array([self._radius, 0, 0])
def _compute_position_body(self):
'''
Must be redefined in herit class with def of self._position_body
'''
# TODO: position du soleil a t dans GCRS
self._position_body = "TDB"
super()._compute_position_body()