#!/usr/bin/env python3

#
# sampling3.py - Sampling illustration
#
# 23May20  Written by Everett Lipman
#
# Copyright (c) 2020 Everett A. Lipman.  All rights reserved.
#
# The following uses are prohibited without prior written permission:
#    * Duplication in any form
#    * Creation of derivative works
#    * Electronic posting and distribution
#
USAGE="""
usage: sampling3.py
"""
import numpy as np
import matplotlib.pyplot as plt
from matplotlib.lines import Line2D

#
# sine wave
#
cycles = 3
spc = 2.4
npts = 1200
maxt = 2*np.pi*cycles

ymax = 2.5
ymin = -ymax*3/5
dot = 6

t = np.linspace(0, maxt, npts)
y = np.cos(t)

step = int(npts/(cycles*spc))
sindex = np.arange(npts)[::step]

fig, (ax1, ax2) = plt.subplots(2, 1, sharex=True, dpi=200)

ax1.plot(t[sindex], y[sindex], 'o', markersize=dot, color='blue')
ax1.bar(t[sindex], y[sindex], width=0.07, color='black')

ax1.set_ylim(ymin, ymax)
ax1.spines['top'].set_visible(False)
ax1.spines['right'].set_visible(False)
ax1.spines['bottom'].set_visible(False)
ax1.axes.get_xaxis().set_visible(False)

# ax2.plot(t, y)
ax2.plot(t[sindex], y[sindex], 'o', markersize=dot, color='blue')
badline, = ax2.plot(t[sindex], y[sindex], color='orange')
ax2.set_ylim(ymin, ymax)
ax2.set_xlabel('$t$')

ax2.spines['top'].set_visible(False)
ax2.spines['right'].set_visible(False)

#
# centered title
#
ax1.text(0.5, 1.15, 'Convolution and signal recovery',
         fontsize=18, horizontalalignment='center', transform=ax1.transAxes)
ax1.text(0.5, 0.9,
         r'$(f*g)(t) = \int_{-\infty}^\infty f(\tau)g(t - \tau) d\tau$',
         fontsize=18, horizontalalignment='center', transform=ax1.transAxes)

ax1.text(1.83, -1.9, r'$f(\tau)$', color='blue', fontsize=16,
         horizontalalignment='left')

ax1.text(7.824, -1.9,
         r'$g(\tau) = \frac{\sin(\pi\tau/T)}{\pi\tau/T}$',
         color='indianred', fontsize=16, horizontalalignment='left')

ax2.text(1.5, 1.7, r'$(f * g)(t)$', color='green', fontsize=16,
         horizontalalignment='left')

g = np.zeros(len(t))
fig.show()
sumline = Line2D(t, g, color='green')

dt = t[sindex[1]] - t[sindex[0]]
omega = 1/dt  # for sinc function; np.sinc(x) is sin(np.pi*x)/(np.pi*x)

addline = True
for i in sindex:
   input('Press <Enter> to advance...')

   ts = t - t[i]  # shifted time
   sincvals = y[i]*np.sinc(omega*ts)  # shifted sinc function

   ax1.plot(t, sincvals)
   g = g + sincvals
   if addline:
      addline = False
      ax2.add_line(sumline)
   sumline.set_data(t, g)
   fig.show()

input('Press <Enter> to remove bad lines...')
badline.set_visible(False)
fig.canvas.draw()

input("\nPress <Enter> to exit...\n")