Source code for ecpi.scripts.ecpi_simu

#!/bin/env python3
"""
"""
import gc
import os
import os.path as osp
import numpy as np
from astropy.io import fits
from astropy.time import Time, TimeDelta

from ecpi.simu.lib.eclairs_channel import SimuEclairsEnergyChannel


[docs]class EcpiSimu: """ """ def __init__(self, duration, pos, n_file, t_0, t_ref=Time('2017-01-01T00:00:00.000')): """**constructor** :param duration: exposure time in minutes. :type duration: float :param pos: ECLAIRs' position :type pos: [float, float, float] :param n_file: number of files to create. :type n_file: int :param t_0: starting time of the observation. Format ISOT :type t_0: astropy.time.Time object """ self.t_0 = t_0 self.t_ref = t_ref self.pos = pos self.duration = float(duration) self.n_file = n_file self.ptg_ori = [0,0,0] self.quater = None self.secl = None
[docs] def set_quater(self, p_quater): self.quater = p_quater
[docs] def create_obj_simu(self, use_irf=False, cxb_mode='no_cxb'): """Initiate simulation :param use_irf: whether or not to use IRF. Default is False :type use_irf: boolean :param cxb_mode: switch for the CXB model if any. Available choices are: 'flat_no_spec', 'flat_spec', 'shapebased_spec' 'shapedbased_spec_earth' and 'no_cxb'. Default is 'no_cxb'. :type cxb_mode: str """ assert hasattr(self, 'quater') # self.att = AttitudeECLAIRs() # self.att.set_attitude_instru_quater(self.quater) self.secl = SimuEclairsEnergyChannel() self.secl.context.set_quaternion_svom(self.quater) self.secl.add_src_point_with_catalog_swift(self.att, use_irf) # self.secl.context.set_quaternion_svom(self.quater) # self.secl.context.set_pos_sat(self.pos) # self.secl.set_src_point_with_custom_catalog(path_cat, use_irf=use_irf) if cxb_mode is 'flat': self.secl.add_src_cxb_flat() elif cxb_mode is 'shape_moretti': self.secl.add_src_cxb() elif cxb_mode is 'flat_moretti': self.secl.add_src_cxb_flat_moretti_spectrum() elif cxb_mode is 'shape_moretti_earth': assert self.pos is not None self.secl.add_src_cxb_earth(self.pos, self.quater)
[docs] def create_event_files(self, dest_dir, sim_arf_mean=True): """Run simulation. :param dest_dir: destination directory. :type dest_dir: str :param sim_arf_mean: whether or not to simulate ARF. Default is True. :type sim_arf: boolean """ assert hasattr(self, 'secl') dtf = (self.duration * 60) / self.n_file # s dt = TimeDelta(dtf, format='sec') times = [(self.t_0 - self.t_ref + s * dt).value for s in range(self.n_file)] for tag, t_ in enumerate(times): print(f"processing {tag+1}th time segment...") t_ = Time(t_, format='mjd', scale='tt').mjd * 24 * 3600 self.secl.context.set_duration(dtf) self.secl.context.set_t_start(t_) self.secl.reset_chan_shadows() self.secl.simu_src() # instrument effect if sim_arf_mean: self.secl.simu_arf_mean() self.secl.simu_noise_poisson() # create events self.secl.create_evts() # writing events file name_file = f"ECL-EVT-SEC-{int(t_)}-DC1-APC.fits" self.secl._evts.write_L1(osp.join(dest_dir, name_file)) del self.secl._evts.tb_evts del self.secl._evts gc.collect()
[docs] def create_attitude_file(self, dest_dir): """Create .fits file with a single line that corresponds to attitude during DC1. :param dest_dir: destination directory. :type dest_dir: str """ ra_ptg, dec_ptg = self.att.ptg_instru(deg=True)[0] ori_ptg = self.att.ori_instru(deg=True) dt = self.t_0 - self.t_ref t_new = Time(dt.value, format='mjd', scale='tt') tstart = t_new.mjd * 24 * 3600 self.ptgori = [ra_ptg, dec_ptg, ori_ptg] times = np.array([tstart]) obsids = np.array([0x01]) velocities = np.array([[0, 0, 0]]) pointings = np.array([self.ptgori]) quaters = np.array([self.quater]) status = np.array([0xAA], dtype=np.uint8) flagaavs = np.array([0x01]) origins = np.array([0x00]) create_fake_attitude( times, quaters, velocities, obsids, status, flagaavs, origins, pointings, dest_dir )
[docs] def create_orbit_file(self, dest_dir): """Create .fits file with a single line that corresponds to orbit during DC1. :param dest_dir: destination directory. :type dest_dir: str """ dt = self.t_0 - self.t_ref t_new = Time(dt.value, format='mjd', scale='tt') tstart = t_new.mjd * 24 * 3600 times = np.array([tstart]) obsids = np.array([0x01]) positions= np.array([self.pos]) velocities = np.array([[0, 0, 0]]) flagaavs = np.array([0x01]) origins = np.array([0x00]) Z = np.zeros(1) create_fake_orbit( times, obsids, positions, velocities, flagaavs, origins, Z, dest_dir )
[docs]def main(): duration = 1 # s pos = [1, 1, 1] # km n_file = 1 t_0 = Time('2017-01-02T00:00:00.000') quater = [1, 0, 0, 0] path_cat = osp.join(os.getenv("HOME"), "Bureau", "catalogs", "Crab_cat_new.fits") dest_dir = osp.join(os.getenv("HOME"), "Bureau") sdc1 = EcpiSimu(duration, pos, n_file, t_0) sdc1.set_quater(quater) sdc1.create_obj_simu() sdc1.create_attitude_file(dest_dir) sdc1.create_orbit_file(dest_dir)
if __name__ == "__main__": main()