#!/usr/bin/env python
# coding: utf-8

# # The Time Axis: GRBs Over Time
# 
# This notebook will demonstrate some of what you can do with WWT’s clock using a dataset of [gamma-ray burst](https://en.wikipedia.org/wiki/Gamma-ray_burst) events.

# ---
# 
# ## Step 1: Setting up WWT
# 
# First, we need to connect to the WWT app. If WWT isn't running yet, activate the JupyterLab command palette from the View menu, and select "WorldWide Telescope". Then all we need to do is run:

# In[ ]:


from pywwt.jupyter import connect_to_app
wwt = await connect_to_app().becomes_ready()


# Now the `wwt` variable will let us talk to the WWT app.
# 
# We'll also set up a utility function to help locate the data files stored alongside this notebook:

# In[ ]:


def datapath(*args):
    from os.path import join
    return join('data', *args)


# ---
# 
# ## Step 2: Loading tabular data with a time-series axis
# 
# Now that the viewer is set up, we need some data. First, some imports:

# In[ ]:


from astropy.coordinates import SkyCoord
from astropy.table import Table
from astropy.time import Time, TimeDelta
from astropy import units as u
from datetime import datetime


# Now, we'll load up a table of GRB data. We'll print out the first few rows to get a sense of the data structure. In this case it's pretty simple.

# In[ ]:


bursts = Table.read(datapath('grb_table_lite.ecsv'), format='ascii.ecsv')
bursts[:5]


# The WWT engine doesn't know enough to know that the Julian Date values above are actually dates (and not just some random floating-point measurements), so we convert them to the ISO8601 format:

# In[ ]:


bursts['dates_ISOT'] = Time(bursts['dates_JD'], format='jd').isot
bursts[:5]


# ---
# 
# ## Step 3: Plotting Data in (2+1)D
# 
# To start showing the data, we create a "table layer":

# In[ ]:


layer = wwt.layers.add_table_layer(
    table = bursts, 
    frame = 'Sky', 
    lon_att = 'ra',  # treat the column named "ra" as a longitude coordinate
    lat_att = 'dec',  # likewise for the column named "dec"
    size_scale = 200,
)


# Now let's tell the engine that we want to pay attention to the time axis of the table. **This will make all of the points disappear**, since the engine's clock is synchronized to your computer — it's plotting data for the current year, and our GRBs only go up to 2014.

# In[ ]:


layer.time_series = True
layer.time_att = 'dates_ISOT'


# We can check the engine's clock using the `get_current_time()` method:

# In[ ]:


wwt.get_current_time()


# Let's demonstrate that things are actually working by centering on the first GRB in both space *and* time — setting the WWT clock to just before it goes off.

# In[ ]:


first_grb = bursts[0]
start_pos = SkyCoord(first_grb['ra'], first_grb['dec'], unit=u.deg)
start_time = Time(first_grb['dates_ISOT'], format='isot') - TimeDelta(5 * u.second)

wwt.center_on_coordinates(start_pos)
wwt.set_current_time(start_time)


# You should see a white dot appear in the screen center after five seconds!
# 
# The WWT timeseries support has the concept of a "decay time". Points that appear on the screen slowly fade over the decay timescale — according to the WWT simulation clock, that is, not wallclock time. The default is 16 days. That's a bit short for visualizing GRBs, so let's make it longer:

# In[ ]:


layer.time_decay = 1 * u.yr


# Now, let's speed time up to see the GRBs pop off over time. If you find the right patch of the sky (namely, the ecliptic), you might also see the Sun and planets moving around at enhanced speed.

# In[ ]:


wwt.play_time(rate=3e6)


# Once you've had your fill, the following function call will reset the clock to run at real time:

# In[ ]:


wwt.play_time()


# How much time passed inside WWT?

# In[ ]:


wwt.get_current_time() - start_time


# ---
# 
# ## Credits
# 
# This notebook was prepared by:
# 
# - O. Justin Otor
# - Peter K. G. Williams
# 
# With thanks to Mark SubbaRao for creating the first version of this demonstration.