ephemeris plugin: custom phase-wrapping (#44)

* implement "wrap_at" functionality
* set initial phase-viewer limits based on wrap_at
* test coverage
* add test-coverage for non-zero dpdt
* update example notebook to use wrap_at instead of t0

Author: kecnry

lcviz/plugins/ephemeris/ephemeris.py

@@ -18,9 +18,10 @@
 
 __all__ = ['Ephemeris']
 
-_default_t0 = 0
-_default_period = 1
-_default_dpdt = 0
+_default_t0 = 0.0
+_default_period = 1.0
+_default_dpdt = 0.0
+_default_wrap_at = 1.0
 
 
 @tray_registry('ephemeris', label="Ephemeris")
@@ -39,6 +40,8 @@ class Ephemeris(PluginTemplateMixin, DatasetSelectMixin):
       Period of the ephemeris, defined at ``t0``.
     * :attr:`dpdt`:
       First derivative of the period of the ephemeris.
+    * :attr:`wrap_at`:
+      Phase at which to wrap (phased data will encompass the range 1-wrap_at to wrap_at).
     * :meth:`ephemeris`
     * :meth:`ephemerides`
     * :meth:`update_ephemeris`
@@ -68,6 +71,7 @@ class Ephemeris(PluginTemplateMixin, DatasetSelectMixin):
     period_step = Float(0.1).tag(sync=True)
     dpdt = FloatHandleEmpty(_default_dpdt).tag(sync=True)
     dpdt_step = Float(0.1).tag(sync=True)
+    wrap_at = FloatHandleEmpty(_default_wrap_at).tag(sync=True)
 
     # PERIOD FINDING
     method_items = List().tag(sync=True)
@@ -83,6 +87,7 @@ def __init__(self, *args, **kwargs):
 
         self._ignore_ephem_change = False
         self._ephemerides = {}
+        self._prev_wrap_at = _default_wrap_at
 
         self.component = EditableSelectPluginComponent(self,
                                                        mode='component_mode',
@@ -108,7 +113,7 @@ def __init__(self, *args, **kwargs):
 
     @property
     def user_api(self):
-        expose = ['component', 'period', 'dpdt', 't0',
+        expose = ['component', 'period', 'dpdt', 't0', 'wrap_at',
                   'ephemeris', 'ephemerides',
                   'update_ephemeris', 'create_phase_viewer',
                   'add_component', 'remove_component', 'rename_component',
@@ -136,6 +141,12 @@ def phase_viewer_ids(self):
     def phase_viewer_id(self):
         return self._phase_viewer_id(self.component_selected)
 
+    @property
+    def phase_viewer(self):
+        if not self.phase_viewer_exists:
+            return None
+        return self.app.get_viewer(self.phase_viewer_id)
+
     @property
     def ephemerides(self):
         return self._ephemerides
@@ -150,17 +161,19 @@ def _times_to_phases_callable(self, component):
             t0 = self.t0
             period = self.period
             dpdt = self.dpdt
+            wrap_at = self.wrap_at
         else:
             ephem = self.ephemerides.get(component, {})
             t0 = ephem.get('t0', _default_t0)
             period = ephem.get('period', _default_period)
             dpdt = ephem.get('dpdt', _default_dpdt)
+            wrap_at = ephem.get('wrap_at', _default_wrap_at)
 
         def _callable(times):
             if dpdt != 0:
-                return np.mod(1./dpdt * np.log(1 + dpdt/period*(times-t0)), 1.0)  # noqa
+                return np.mod(1./dpdt * np.log(1 + dpdt/period*(times-t0)) + (1-wrap_at), 1.0) - (1-wrap_at)  # noqa
             else:
-                return np.mod((times-t0)/period, 1.0)
+                return np.mod((times-t0)/period + (1-wrap_at), 1.0) - (1-wrap_at)
 
         return _callable
 
@@ -242,7 +255,8 @@ def create_phase_viewer(self):
         if self.phase_comp_lbl not in [comp.label for comp in dc[0].components]:
             self.update_ephemeris()  # calls _update_all_phase_arrays
 
-        if not self.phase_viewer_exists:
+        create_phase_viewer = not self.phase_viewer_exists
+        if create_phase_viewer:
             # TODO: stack horizontally by default?
             self.app._on_new_viewer(NewViewerMessage(PhaseScatterView, data=None, sender=self.app),
                                     vid=phase_viewer_id, name=phase_viewer_id)
@@ -254,6 +268,8 @@ def create_phase_viewer(self):
                 self.app.set_data_visibility(phase_viewer_id, data.label, visible == 'visible')
 
         pv = self.app.get_viewer(phase_viewer_id)
+        if create_phase_viewer:
+            pv.state.x_min, pv.state.x_max = (self.wrap_at-1, self.wrap_at)
         pv.state.x_att = self.phase_cids[self.component_selected]
         return pv
 
@@ -314,17 +330,18 @@ def _change_component(self, *args):
         self.t0 = ephem.get('t0', self.t0)
         self.period = ephem.get('period', self.period)
         self.dpdt = ephem.get('dpdt', self.dpdt)
+        self.wrap_at = ephem.get('wrap_at', self.wrap_at)
 
         # if this is a new component, update those default values back to the dictionary
-        self.update_ephemeris(t0=self.t0, period=self.period, dpdt=self.dpdt)
+        self.update_ephemeris(t0=self.t0, period=self.period, dpdt=self.dpdt, wrap_at=self.wrap_at)
         self._ignore_ephem_change = False
         if ephem:
             # if there were any changes applied by accessing the dictionary,
             # then we need to update phasing, etc (since we set _ignore_ephem_change)
             # otherwise, this is a new component and there is no need.
             self._ephem_traitlet_changed()
 
-    def update_ephemeris(self, component=None, t0=None, period=None, dpdt=None):
+    def update_ephemeris(self, component=None, t0=None, period=None, dpdt=None, wrap_at=None):
         """
         Update the ephemeris for a given component.
 
@@ -337,7 +354,9 @@ def update_ephemeris(self, component=None, t0=None, period=None, dpdt=None):
         period : float, optional
             value of period to replace
         dpdt : float, optional
-            value of period to replace
+            value of dpdt to replace
+        wrap_at : float, optional
+            value of wrap_at to replace
 
         Returns
         -------
@@ -350,7 +369,7 @@ def update_ephemeris(self, component=None, t0=None, period=None, dpdt=None):
             raise ValueError(f"component must be one of {self.component.choices}")
 
         existing_ephem = self._ephemerides.get(component, {})
-        for name, value in {'t0': t0, 'period': period, 'dpdt': dpdt}.items():
+        for name, value in {'t0': t0, 'period': period, 'dpdt': dpdt, 'wrap_at': wrap_at}.items():
             if value is not None:
                 existing_ephem[name] = value
                 if component == self.component_selected:
@@ -360,11 +379,11 @@ def update_ephemeris(self, component=None, t0=None, period=None, dpdt=None):
         self._update_all_phase_arrays(component=component)
         return existing_ephem
 
-    @observe('period', 'dpdt', 't0')
+    @observe('period', 'dpdt', 't0', 'wrap_at')
     def _ephem_traitlet_changed(self, event={}):
         if self._ignore_ephem_change:
             return
-        for value in (self.period, self.dpdt, self.t0):
+        for value in (self.period, self.dpdt, self.t0, self.wrap_at):
             if not isinstance(value, (int, float)):
                 return
         if self.period <= 0:
@@ -383,6 +402,17 @@ def round_to_1(x):
         else:
             self._update_all_phase_arrays(component=self.component_selected)
 
+        # update zoom-limits if wrap_at was changed
+        if event.get('name') == 'wrap_at':
+            old = event.get('old') if event.get('old') != '' else self._prev_wrap_at
+            if event.get('new') != '':
+                pvs = self.phase_viewer.state
+                delta_phase = event.get('new') - old
+                pvs.x_min, pvs.x_max = pvs.x_min + delta_phase, pvs.x_max + delta_phase
+                # we need to cache the old value since it could become a string
+                # if the widget is cleared
+                self._prev_wrap_at = event.get('new')
+
         # update step-sizes
         self.period_step = round_to_1(self.period/5000)
         self.dpdt_step = max(round_to_1(abs(self.period * self.dpdt)/1000) if self.dpdt != 0 else 0,

lcviz/plugins/ephemeris/ephemeris.vue

@@ -74,6 +74,20 @@
         ></v-text-field>
       </v-row>
 
+      <v-row>
+        <v-text-field
+          ref="wrap_at"
+          type="number"
+          label="Wrapping phase"
+          v-model.number="wrap_at"
+          :step="0.1"
+          type="number"
+          :hint="'Phased data will encompass the range '+wrap_at_range+'.'"
+          persistent-hint
+          :rules="[() => wrap_at!=='' || 'This field is required']"
+        ></v-text-field>
+      </v-row>
+
       <j-plugin-section-header>Period Finding/Refining</j-plugin-section-header>
 
       <plugin-dataset-select
@@ -135,3 +149,14 @@
 
   </j-tray-plugin>
 </template>
+
+<script>
+module.exports = {
+  computed: {
+    wrap_at_range() {
+      const lower = this.wrap_at - 1
+      return '('+lower.toFixed(2)+', '+this.wrap_at.toFixed(2)+')'
+    },
+  }
+};
+</script>

lcviz/tests/test_plugin_ephemeris.py

@@ -1,3 +1,4 @@
+
 def test_docs_snippets(helper, light_curve_like_kepler_quarter):
     lcviz, lc = helper, light_curve_like_kepler_quarter
 
@@ -29,6 +30,12 @@ def test_plugin_ephemeris(helper, light_curve_like_kepler_quarter):
     ephem._obj.vue_period_halve()
     assert ephem.period == 3.14
 
+    pv = ephem._obj.phase_viewer
+    # original limits are set to 0->1 (technically 1-phase_wrap -> phase_wrap)
+    assert (pv.state.x_min, pv.state.x_max) == (0.0, 1.0)
+    ephem.wrap_at = 0.5
+    assert (pv.state.x_min, pv.state.x_max) == (-0.5, 0.5)
+
     ephem.add_component('custom component')
     assert not ephem._obj.phase_viewer_exists
     ephem.create_phase_viewer()
@@ -60,3 +67,6 @@ def test_plugin_ephemeris(helper, light_curve_like_kepler_quarter):
     assert ephem._obj.method_err == ''
     ephem._obj.vue_adopt_period_at_max_power()
     assert ephem.period != 2
+
+    # test coverage for non-zero dpdt
+    ephem.dpdt = 0.00001

lcviz/viewers.py

@@ -113,6 +113,7 @@ def _set_plot_x_axes(self, dc, component_labels, light_curve):
             xlabel = f'{str(x_unit.physical_type).title()} ({x_unit})'
 
         self.figure.axes[0].label = xlabel
+        self.figure.axes[0].num_ticks = 5
 
     def _set_plot_y_axes(self, dc, component_labels, light_curve):
         self.state.y_att = dc[0].components[component_labels.index('flux')]
@@ -142,6 +143,8 @@ def _set_plot_y_axes(self, dc, component_labels, light_curve):
         self.figure.axes[0].tick_format = 'g'
         self.figure.axes[1].tick_format = 'g'
 
+        self.figure.axes[1].num_ticks = 5
+
     def _expected_subset_layer_default(self, layer_state):
         super()._expected_subset_layer_default(layer_state)
 
@@ -212,6 +215,7 @@ def _set_plot_x_axes(self, dc, component_labels, light_curve):
         # setting of y_att will be handled by ephemeris plugin
         self.state.x_att = dc[0].components[component_labels.index(f'phase:{self.ephemeris_component}')]  # noqa
         self.figure.axes[0].label = 'phase'
+        self.figure.axes[0].num_ticks = 5
 
     def times_to_phases(self, times):
         ephem = self.jdaviz_helper.plugins.get('Ephemeris', None)

notebooks/LCvizExample.ipynb

@@ -119,8 +119,8 @@
     "    (morris2017_epoch - reference_time).to_value(time_coordinates.unit) % eph.period\n",
     ")\n",
     "\n",
-    "# offset the phase axis so mid-transit occurs at phase=0.5:\n",
-    "eph.t0 += 0.5 * eph.period"
+    "# offset the wrapping phase so the transit (at phase 0) displays at center\n",
+    "eph.wrap_at = 0.5"
    ]
   },
   {