From 3dfcefa3efe579a44275680ee191042cd07cd4ef Mon Sep 17 00:00:00 2001 From: John ZuHone Date: Wed, 29 Dec 2021 13:12:19 -0500 Subject: [PATCH 1/9] Add a feature that accounts for Earth illumination --- xija/component/heat.py | 12 ++++++++++++ 1 file changed, 12 insertions(+) diff --git a/xija/component/heat.py b/xija/component/heat.py index c89ac4a..be2178f 100644 --- a/xija/component/heat.py +++ b/xija/component/heat.py @@ -544,6 +544,7 @@ class EarthHeat(PrecomputedHeatPower): def __init__(self, model, node, orbitephem0_x, orbitephem0_y, orbitephem0_z, + solarephem0_x, solarephem0_y, solarephem0_z, aoattqt1, aoattqt2, aoattqt3, aoattqt4, k=1.0): ModelComponent.__init__(self, model) @@ -551,6 +552,9 @@ def __init__(self, model, node, self.orbitephem0_x = self.model.get_comp(orbitephem0_x) self.orbitephem0_y = self.model.get_comp(orbitephem0_y) self.orbitephem0_z = self.model.get_comp(orbitephem0_z) + self.solarephem0_x = self.model.get_comp(solarephem0_x) + self.solarephem0_y = self.model.get_comp(solarephem0_y) + self.solarephem0_z = self.model.get_comp(solarephem0_z) self.aoattqt1 = self.model.get_comp(aoattqt1) self.aoattqt2 = self.model.get_comp(aoattqt2) self.aoattqt3 = self.model.get_comp(aoattqt3) @@ -613,12 +617,15 @@ def dvals(self): # Collect individual MSIDs for use in calc_earth_vis() ephem_xyzs = [getattr(self, 'orbitephem0_{}'.format(x)) for x in ('x', 'y', 'z')] + solar_xyzs = [getattr(self, 'solarephem0_{}'.format(x)) + for x in ('x', 'y', 'z')] aoattqt_1234s = [getattr(self, 'aoattqt{}'.format(x)) for x in range(1, 5)] # Note: the copy() here is so that the array becomes contiguous in # memory and allows numba to run faster (and avoids NumbaPerformanceWarning: # np.dot() is faster on contiguous arrays). ephems = np.array([x.dvals for x in ephem_xyzs]).transpose().copy() + solars = np.array([x.dvals for x in solar_xyzs]).transpose().copy() q_atts = np.array([x.dvals for x in aoattqt_1234s]).transpose() # Q_atts can have occasional bad values, maybe because the @@ -640,6 +647,11 @@ def dvals(self): else: self.calc_earth_vis_from_taco(ephems, q_atts) + cos = np.sum(ephems*solars, axis=1) + es = np.sum(ephems*ephems, axis=1)*np.sum(solars*solars, axis=1) + cos /= np.sqrt(es) + self._dvals *= 0.5*(1.0+cos) + self.put_cache() return self._dvals From 3f78f9ac58f498e86a84634734898bc2179d88d5 Mon Sep 17 00:00:00 2001 From: John ZuHone Date: Tue, 25 Jan 2022 14:13:10 -0500 Subject: [PATCH 2/9] Make this calculation optional --- xija/component/heat.py | 41 +++++++++++++++++++++++++++++------------ 1 file changed, 29 insertions(+), 12 deletions(-) diff --git a/xija/component/heat.py b/xija/component/heat.py index be2178f..3b8a2e0 100644 --- a/xija/component/heat.py +++ b/xija/component/heat.py @@ -544,21 +544,31 @@ class EarthHeat(PrecomputedHeatPower): def __init__(self, model, node, orbitephem0_x, orbitephem0_y, orbitephem0_z, - solarephem0_x, solarephem0_y, solarephem0_z, aoattqt1, aoattqt2, aoattqt3, aoattqt4, - k=1.0): + solarephem0_x=None, solarephem0_y=None, + solarephem0_z=None, k=1.0): ModelComponent.__init__(self, model) self.node = self.model.get_comp(node) self.orbitephem0_x = self.model.get_comp(orbitephem0_x) self.orbitephem0_y = self.model.get_comp(orbitephem0_y) self.orbitephem0_z = self.model.get_comp(orbitephem0_z) - self.solarephem0_x = self.model.get_comp(solarephem0_x) - self.solarephem0_y = self.model.get_comp(solarephem0_y) - self.solarephem0_z = self.model.get_comp(solarephem0_z) self.aoattqt1 = self.model.get_comp(aoattqt1) self.aoattqt2 = self.model.get_comp(aoattqt2) self.aoattqt3 = self.model.get_comp(aoattqt3) self.aoattqt4 = self.model.get_comp(aoattqt4) + if solarephem0_x is None: + self.solarephem0_x = None + else: + self.solarephem0_x = self.model.get_comp(solarephem0_x) + if solarephem0_y is None: + self.solarephem0_y = None + else: + self.solarephem0_y = self.model.get_comp(solarephem0_y) + if solarephem0_z is None: + self.solarephem0_z = None + else: + self.solarephem0_z = self.model.get_comp(solarephem0_z) + self.n_mvals = 1 self.add_par('k', k, min=0.0, max=2.0) @@ -617,15 +627,12 @@ def dvals(self): # Collect individual MSIDs for use in calc_earth_vis() ephem_xyzs = [getattr(self, 'orbitephem0_{}'.format(x)) for x in ('x', 'y', 'z')] - solar_xyzs = [getattr(self, 'solarephem0_{}'.format(x)) - for x in ('x', 'y', 'z')] aoattqt_1234s = [getattr(self, 'aoattqt{}'.format(x)) for x in range(1, 5)] # Note: the copy() here is so that the array becomes contiguous in # memory and allows numba to run faster (and avoids NumbaPerformanceWarning: # np.dot() is faster on contiguous arrays). ephems = np.array([x.dvals for x in ephem_xyzs]).transpose().copy() - solars = np.array([x.dvals for x in solar_xyzs]).transpose().copy() q_atts = np.array([x.dvals for x in aoattqt_1234s]).transpose() # Q_atts can have occasional bad values, maybe because the @@ -647,10 +654,20 @@ def dvals(self): else: self.calc_earth_vis_from_taco(ephems, q_atts) - cos = np.sum(ephems*solars, axis=1) - es = np.sum(ephems*ephems, axis=1)*np.sum(solars*solars, axis=1) - cos /= np.sqrt(es) - self._dvals *= 0.5*(1.0+cos) + # This next bit optionally checks to see if the solar ephemeris + # was passed in, and if it was it computes the fraction of the + # Earth's surface that is illuminated by the Sun. + solar_xyzs = [getattr(self, 'solarephem0_{}'.format(x)) + for x in ('x', 'y', 'z')] + + if set(solar_xyzs) != {None}: + + solars = np.array([x.dvals for x in solar_xyzs]).transpose().copy() + + cos = np.sum(ephems*solars, axis=1) + es = np.sum(ephems*ephems, axis=1)*np.sum(solars*solars, axis=1) + cos /= np.sqrt(es) + self._dvals *= 0.5*(1.0+cos) self.put_cache() From b6b6dff63ec921ecef4e2e0b5406c5ecefd2bcb3 Mon Sep 17 00:00:00 2001 From: John ZuHone Date: Tue, 15 Mar 2022 11:55:27 -0400 Subject: [PATCH 3/9] Better check for not supplying solarephem components correctly --- xija/component/heat.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/xija/component/heat.py b/xija/component/heat.py index 3b8a2e0..7a4cc3a 100644 --- a/xija/component/heat.py +++ b/xija/component/heat.py @@ -659,9 +659,9 @@ def dvals(self): # Earth's surface that is illuminated by the Sun. solar_xyzs = [getattr(self, 'solarephem0_{}'.format(x)) for x in ('x', 'y', 'z')] - - if set(solar_xyzs) != {None}: - + + if not None in solar_xyzs: + solars = np.array([x.dvals for x in solar_xyzs]).transpose().copy() cos = np.sum(ephems*solars, axis=1) From 508777674cdc8ea9b01d84a55aa38e217f286f38 Mon Sep 17 00:00:00 2001 From: John ZuHone Date: Mon, 21 Mar 2022 14:56:27 -0400 Subject: [PATCH 4/9] Small code simplifications --- xija/component/heat.py | 5 ++--- 1 file changed, 2 insertions(+), 3 deletions(-) diff --git a/xija/component/heat.py b/xija/component/heat.py index 7a4cc3a..fcc250a 100644 --- a/xija/component/heat.py +++ b/xija/component/heat.py @@ -656,9 +656,8 @@ def dvals(self): # This next bit optionally checks to see if the solar ephemeris # was passed in, and if it was it computes the fraction of the - # Earth's surface that is illuminated by the Sun. - solar_xyzs = [getattr(self, 'solarephem0_{}'.format(x)) - for x in ('x', 'y', 'z')] + # Earth's surface that is illuminated by the Sun. + solar_xyzs = [getattr(self, f'solarephem0_{x}') for x in 'xyz'] if not None in solar_xyzs: From aac2e151c41e95438f64358c8c9e717596f34aac Mon Sep 17 00:00:00 2001 From: John ZuHone Date: Mon, 21 Mar 2022 14:57:16 -0400 Subject: [PATCH 5/9] Improved syntax --- xija/component/heat.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/xija/component/heat.py b/xija/component/heat.py index fcc250a..c2baa60 100644 --- a/xija/component/heat.py +++ b/xija/component/heat.py @@ -659,7 +659,7 @@ def dvals(self): # Earth's surface that is illuminated by the Sun. solar_xyzs = [getattr(self, f'solarephem0_{x}') for x in 'xyz'] - if not None in solar_xyzs: + if None not in solar_xyzs: solars = np.array([x.dvals for x in solar_xyzs]).transpose().copy() From b2a2b79ee2a74fa255c64e80dcb698e735f3a13d Mon Sep 17 00:00:00 2001 From: John ZuHone Date: Mon, 21 Mar 2022 14:59:36 -0400 Subject: [PATCH 6/9] Keep the earth phase around and add a plotting option --- xija/component/heat.py | 14 +++++++++++++- 1 file changed, 13 insertions(+), 1 deletion(-) diff --git a/xija/component/heat.py b/xija/component/heat.py index c2baa60..dafac80 100644 --- a/xija/component/heat.py +++ b/xija/component/heat.py @@ -666,7 +666,8 @@ def dvals(self): cos = np.sum(ephems*solars, axis=1) es = np.sum(ephems*ephems, axis=1)*np.sum(solars*solars, axis=1) cos /= np.sqrt(es) - self._dvals *= 0.5*(1.0+cos) + self.earth_phase = 0.5*(1.0+cos) + self._dvals *= self.earth_phase self.put_cache() @@ -718,6 +719,17 @@ def update(self): ) self.tmal_floats = () + def plot_earth_phase__time(self, fig, ax): + lines = ax.get_lines() + if not lines: + plot_cxctime(self.model.times, self.earth_phase, ls='-', + color='#386cb0', fig=fig, ax=ax) + ax.grid() + ax.set_title('Earth Phase') + ax.set_ylabel('Earth Phase') + else: + lines[0].set_data(self.model_plotdate, self.dvals) + def plot_data__time(self, fig, ax): lines = ax.get_lines() if not lines: From 9e871000ebaee13a0215a2a8a529b8aa287e5ac8 Mon Sep 17 00:00:00 2001 From: John ZuHone Date: Mon, 21 Mar 2022 16:37:52 -0400 Subject: [PATCH 7/9] Make this two parameters--one proportional to the Earth phase and one not --- xija/component/heat.py | 13 ++++++++----- 1 file changed, 8 insertions(+), 5 deletions(-) diff --git a/xija/component/heat.py b/xija/component/heat.py index dafac80..a227605 100644 --- a/xija/component/heat.py +++ b/xija/component/heat.py @@ -546,7 +546,7 @@ def __init__(self, model, node, orbitephem0_x, orbitephem0_y, orbitephem0_z, aoattqt1, aoattqt2, aoattqt3, aoattqt4, solarephem0_x=None, solarephem0_y=None, - solarephem0_z=None, k=1.0): + solarephem0_z=None, k=1.0, k2=1.0): ModelComponent.__init__(self, model) self.node = self.model.get_comp(node) self.orbitephem0_x = self.model.get_comp(orbitephem0_x) @@ -568,9 +568,12 @@ def __init__(self, model, node, self.solarephem0_z = None else: self.solarephem0_z = self.model.get_comp(solarephem0_z) - + self.use_earth_phase = np.all([getattr(self, f"solarephem0_{ax}") is not None + for ax in "xyz"]) self.n_mvals = 1 self.add_par('k', k, min=0.0, max=2.0) + self.add_par('k2', k2, min=0.0, max=2.0) + self.earth_phase = 1.0 def calc_earth_vis_from_grid(self, ephems, q_atts): import astropy_healpix @@ -659,15 +662,13 @@ def dvals(self): # Earth's surface that is illuminated by the Sun. solar_xyzs = [getattr(self, f'solarephem0_{x}') for x in 'xyz'] - if None not in solar_xyzs: - + if self.use_earth_phase: solars = np.array([x.dvals for x in solar_xyzs]).transpose().copy() cos = np.sum(ephems*solars, axis=1) es = np.sum(ephems*ephems, axis=1)*np.sum(solars*solars, axis=1) cos /= np.sqrt(es) self.earth_phase = 0.5*(1.0+cos) - self._dvals *= self.earth_phase self.put_cache() @@ -713,6 +714,8 @@ def get_cached(self): def update(self): self.mvals = self.k * self.dvals + if self.use_earth_phase: + self.mvals += self.k2 * self.earth_phase * self.dvals self.tmal_ints = (tmal.OPCODES['precomputed_heat'], self.node.mvals_i, # dy1/dt index self.mvals_i, # mvals with precomputed heat input From 1ce7f6aeb814fe79103f1b4ce67a0e5b63139c5e Mon Sep 17 00:00:00 2001 From: John ZuHone Date: Mon, 21 Mar 2022 16:40:34 -0400 Subject: [PATCH 8/9] Add note about TWG --- xija/component/heat.py | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/xija/component/heat.py b/xija/component/heat.py index a227605..6763f48 100644 --- a/xija/component/heat.py +++ b/xija/component/heat.py @@ -659,7 +659,9 @@ def dvals(self): # This next bit optionally checks to see if the solar ephemeris # was passed in, and if it was it computes the fraction of the - # Earth's surface that is illuminated by the Sun. + # Earth's surface that is illuminated by the Sun. Originally + # discussed at: + # https://occweb.cfa.harvard.edu/twiki/bin/view/ThermalWG/MeetingNotes2022x03x03 solar_xyzs = [getattr(self, f'solarephem0_{x}') for x in 'xyz'] if self.use_earth_phase: From 29e34c76d9eba1537da7205c908fdb92e6b71cd5 Mon Sep 17 00:00:00 2001 From: John ZuHone Date: Mon, 28 Mar 2022 10:58:02 -0400 Subject: [PATCH 9/9] Make this an optional parameter --- xija/component/heat.py | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/xija/component/heat.py b/xija/component/heat.py index 6763f48..c8dd356 100644 --- a/xija/component/heat.py +++ b/xija/component/heat.py @@ -546,7 +546,7 @@ def __init__(self, model, node, orbitephem0_x, orbitephem0_y, orbitephem0_z, aoattqt1, aoattqt2, aoattqt3, aoattqt4, solarephem0_x=None, solarephem0_y=None, - solarephem0_z=None, k=1.0, k2=1.0): + solarephem0_z=None, k=1.0, k2=None): ModelComponent.__init__(self, model) self.node = self.model.get_comp(node) self.orbitephem0_x = self.model.get_comp(orbitephem0_x) @@ -572,7 +572,8 @@ def __init__(self, model, node, for ax in "xyz"]) self.n_mvals = 1 self.add_par('k', k, min=0.0, max=2.0) - self.add_par('k2', k2, min=0.0, max=2.0) + if k2 is not None: + self.add_par('k2', k2, min=0.0, max=2.0) self.earth_phase = 1.0 def calc_earth_vis_from_grid(self, ephems, q_atts):