<jupyter_start><jupyter_text>Dynamical Power Spectra (on fake data)<jupyter_code>%matplotlib inline
# import some modules
import numpy as np
import matplotlib.pyplot as plt
import stingray
# choose style of plots, `seaborn-talk` produce nice big figures
plt.style.use('seaborn-talk')<jupyter_output><empty_output><jupyter_text>Generate a fake lightcurve<jupyter_code># Array of timestamps, 10000 bins from 1s to 100s
times = np.linspace(1,100,10000)

# base component of the lightcurve, poisson-like
# the averaged count-rate is 100 counts/bin
noise = np.random.poisson(100,10000)

# time evolution of the frequency of our fake periodic signal
# the frequency changes with a sinusoidal shape around the value 24Hz
freq = 25 + 1.2*np.sin(2*np.pi*times/130)

# Our fake periodic variability with drifting frequency
# the amplitude of this variability is 10% of the base flux
var = 10*np.sin(2*np.pi*freq*times)

# The signal of our lightcurve is equal the base flux plus the variable flux
signal = noise+var
# Create the lightcurve object
lc = stingray.Lightcurve(times, signal)<jupyter_output><empty_output><jupyter_text>Visualizing the lightcurve<jupyter_code>lc.plot(labels=['Time (s)', 'Counts / bin'])
plt.title('Lightcurve')<jupyter_output><empty_output><jupyter_text>Zomming in..<jupyter_code>lc.plot(labels=['Time (s)', 'Counts / bin'], axis=[20,23,50,160])
plt.title('Zoomed in Lightcurve')<jupyter_output><empty_output><jupyter_text>A power spectrum of this lightcurve..<jupyter_code>ps = stingray.AveragedPowerspectrum(lc, segment_size=3, norm='leahy')
plt.plot(ps.freq, ps.power, label='segment size = {}s \n number of segments = {}'.format(3, int(lc.tseg/3)))
plt.title('Averaged Powerspectrum')
plt.xlabel('Frequency (Hz)')
plt.ylabel('Power')
plt.legend()<jupyter_output><empty_output><jupyter_text>It looks like we have at least 2 frequencies.  Let's look at the Dynamic Powerspectrum..<jupyter_code>dps = stingray.DynamicalPowerspectrum(lc, segment_size=3)
extent = min(dps.time), max(dps.time), min(dps.freq), max(dps.freq)
plt.imshow(dps.dyn_ps, aspect="auto", origin="lower", vmax=0.001,
           interpolation="none", extent=extent)
plt.title('Dynamic Powerspecttrum')
plt.xlabel('Time (s)')
plt.ylabel('Frequency (Hz)')
plt.colorbar(label='Power')<jupyter_output><empty_output><jupyter_text>It is actually only one feature drifiting along time   Rebinning in Frequency<jupyter_code>print("The current frequency resolution is {}".format(dps.df))<jupyter_output>The current frequency resolution is 0.33443708609271494<jupyter_text>Let's rebin to a frequency resolution of 1 Hz and using the average of the power<jupyter_code>dps.rebin_frequency(df_new=1.0, method="average")
print("The new frequency resolution is {}".format(dps.df))<jupyter_output>The new frequency resolution is 1.0<jupyter_text>Let's see how the Dynamical Powerspectrum looks now<jupyter_code>extent = min(dps.time), max(dps.time), min(dps.freq), max(dps.freq)
plt.imshow(dps.dyn_ps, origin="lower", aspect="auto",
           interpolation="none", extent=extent)
plt.colorbar()
plt.ylim(15, 30)<jupyter_output><empty_output><jupyter_text>Rebin time Let's  rebin our matrix in the time axis<jupyter_code>print("The current time resolution is {}".format(dps.dt))<jupyter_output>The current time resolution is 3.0<jupyter_text>Let's rebin to a time resolution of 4 s<jupyter_code>dps.rebin_time(dt_new=4.0, method="average")
print("The new time resolution is {}".format(dps.dt))
extent = min(dps.time), max(dps.time), min(dps.freq), max(dps.freq)
plt.imshow(dps.dyn_ps, origin="lower", aspect="auto",
           interpolation="none", extent=extent)
plt.colorbar()
plt.ylim(15,30)<jupyter_output><empty_output><jupyter_text>Let's trace that drifiting feature.<jupyter_code># By looking into the maximum power of each segment
max_pos = dps.trace_maximum()
plt.plot(dps.time, dps.freq[max_pos], color='red', alpha=1)
plt.xlabel('Time (s)')
plt.ylabel('Frequency (Hz)')
plt.title('Detected frequency drift')<jupyter_output><empty_output><jupyter_text>Overlaying this traced function with the Dynamical Powerspectrum<jupyter_code>extent = min(dps.time), max(dps.time), min(dps.freq), max(dps.freq)
plt.imshow(dps.dyn_ps, aspect="auto", origin="lower", vmax=0.001,
           interpolation="none", extent=extent, alpha=0.6)
plt.plot(dps.time, dps.freq[max_pos], color='C3', lw=5, alpha=1, label='drifiting function')

plt.ylim(15,30) # zoom-in around 24 hertz

plt.title('Overlay of Drifting fuction and Dynamic Powerspecttrum')
plt.xlabel('Time (s)')
plt.ylabel('Frequency (Hz)')
plt.colorbar(label='Power')
plt.legend()<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.13557858376511225, 'fraction_numerical': 0.02461139896373057, 'mean_word_length': 7.057391304347826, 'pattern_counts': {'":': 0, '<': 50, '<?xml version=': 0, '>': 50, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '13', 'dir': 'jupyter-structured-clean-dedup', 'id': '364890', 'license': 'MIT', 'n_tokens_mistral': 1791, 'n_tokens_neox': 1640, 'n_words': 380, 'path': 'DynamicalPowerspectrum/DynamicalPowerspectrum_tutorial_[fake_data].ipynb', 'repo_name': 'jdswinbank/notebooks', 'size': '4632'}