Commit 90c4a2ca by Christian Rohlfing

### - moved find_intervals function from Notebook to ient.nb

 %% Cell type:code id: tags:  python # Copyright 2019 Institut für Nachrichtentechnik, RWTH Aachen University %matplotlib widget import ipywidgets as widgets from ipywidgets import interact, interactive, fixed, Layout from IPython.display import clear_output, display, HTML from scipy import signal # convolution from scipy.signal import find_peaks # peak finder from ient_nb.ient_plots import * from ient_nb.ient_signals import * def find_intervals(s, thresh, delta): """ Find intervals of signal s by searching for delta-functions in the second derivative of s Parameters ---------- s : array_like The signal thresh : float Threshold for delta search delta : float Sampling period Returns ------- intervals_s, peaks, dd : *intervals* are the intervals, *peaks* the found peaks and *dd* the second derivative of s. """ # derivative of derivative shows delta functions at discontinuities dd = np.diff(np.diff(s, prepend=0), prepend=0)/delta**2 # find delta functions peaks, _ = find_peaks(np.abs(dd), prominence=thresh) # return interval limits in seconds intervals = np.round(tau[peaks]*10)/10 return intervals, peaks, dd  %% Cell type:markdown id: tags:
# Demonstrator Faltung Zum Starten: Im Menü: Run Run All Cells auswählen. ## Einleitung Im Folgenden wird das Faltungsintegral $$g(t) = s(t)\ast h(t) = \int\limits_{-\infty}^{\infty} s(\tau) h(t-\tau) \,\mathrm{d}\tau$$ betrachtet. %% Cell type:markdown id: tags: ## Demo Wähle $s(t)$ und $h(t)$ sowie jeweils Verschiebung $t_0$ und Dehnungsfaktor $T$ für beide Signale: $s\left(\frac{t-t_0}{T}\right)$ und $h\left(\frac{t-t_0}{T}\right)$. Zusätzlich zu Elementarsignalen kann auch eine frei definierbare Funktion $s_0(t)$ zur Faltung verwendet werden. %% Cell type:code id: tags:  python s_0 = lambda t: rect(t/2-1/2)*(-t)  %% Cell type:code id: tags:  python def convolution(s, h): # Convolve s and h numerically return signal.convolve(s(tau), h(tau), mode='same')*deltat (tau,deltat) = np.linspace(-15, 15, 10001, retstep=True) # tau axis interval_cnt = 0; # Interval counter for stepwise plot # Signals signal_types = { 'Rechteck' : rect, 'Dreieck' : tri, 'Sprungfunktion' : unitstep, 'si-Funktion' : lambda t: si(t*np.pi), 'Exponentialimpuls' : lambda t: unitstep(t)*np.exp(-t), 'Gauß-Signal' : gauss, 'Doppelrechteck' : lambda t: rect(t*2+0.5)-rect(t*2-0.5), 'Rampe' : lambda t: t*rect(t-0.5), 'Versch. Rechteck' : lambda t: -rect(t-0.5), 'Eigene Kreation s0(t)' : s_0, } # Plot fig0, axs0 = plt.subplots(1, 2, figsize=(12,2)) def update_signals(s_type, s_T, s_t0, h_type, h_T, h_t0): # Get s(t) and h(t) global s, h, gt, intervals_g # reused in second interactive plot s = lambda tau: signal_types[s_type]((tau-s_t0)/s_T); # signal on time axis (shifted by t0 and stretched by T) h = lambda tau: signal_types[h_type]((tau-h_t0)/h_T); #intervals_s = np.array(signal_intervals[w_s_type.value])*s_T + s_t0; # get signal intervals and map to time axis (shifted by t0 and stretched by T) #intervals_h = np.array(signal_intervals[w_h_type.value])*h_T + h_t0; intervals_s,_,_ = find_intervals(s(tau), 0.49*np.max(np.abs(s(tau)))/deltat, deltat) intervals_h,_,_ = find_intervals(h(tau), 0.49*np.max(np.abs(h(tau)))/deltat, deltat) # Get g(t) = s(t) \ast h(t) gt = convolution(s, h) # numerical convolution intervals_g = np.unique(np.tile(intervals_h, len(intervals_s)) + np.repeat(intervals_s, len(intervals_h))) # intervals of g(t) if intervals_g.size == 0: # intervals_g = np.concatenate([intervals_s, intervals_h]) # Plot if not axs0[0].lines: # plot s(t) and g(t) ax = axs0[0]; ax.plot(tau, s(tau), 'rwth'); ax.set_xlabel(r'$\rightarrow t$'); ax.set_ylabel(r'$\uparrow s(t)$') ax.set_xlim([-2.9, 2.9]); ax.set_ylim([-1.19, 1.19]); ient_axis(ax); ient_grid(ax); ax = axs0[1]; ax.plot(tau, h(tau), 'rwth'); ax.set_xlabel(r'$\rightarrow t$'); ax.set_ylabel(r'$\uparrow h(t)$') ax.set_xlim(axs0[0].get_xlim()); ax.set_ylim(axs0[0].get_ylim()); ient_axis(ax); ient_grid(ax); fig0.tight_layout(); else: # update lines axs0[0].lines[0].set_ydata(s(tau)); axs0[1].lines[0].set_ydata(h(tau)); try: # if convolution figure is already opened, update s(tau) if axs[0].lines: axs[0].lines[1].set_ydata(s(tau)); ient_update_ylim(axs[0], np.concatenate((h(tau), s(tau))), 0.19, 5); ient_update_ylim(axs[1], gt, 0.19, 5); update_plot(-2) # update convolution plot update_plot_intervals() # update interval lines except: pass ient_update_ylim(axs0[0], s(tau), 0.19, 5); ient_update_ylim(axs0[1], h(tau), 0.19, 5); # Widgets w_s_type = widgets.Dropdown(options=list(signal_types.keys()), description=r'Wähle $s(t)$:') w_s_T = widgets.FloatSlider(min=0.5, max=4, value=1, step=.1, description=r'Dehnung T', style=ient_wdgtl_style) w_s_t0 = s_t0=widgets.FloatSlider(min=-2, max=2, value=0, step=.1, description=r'Verschiebung $t_0$', style=ient_wdgtl_style) w_h_type = widgets.Dropdown(options=list(signal_types.keys()), description=r'Wähle $h(t)$:') w_h_T = widgets.FloatSlider(min=0.5, max=4, value=1, step=.1, description=r'Dehnung T', style=ient_wdgtl_style) w_h_t0 = s_t0=widgets.FloatSlider(min=-2, max=2, value=0, step=.1, description=r'Verschiebung $t_0$', style=ient_wdgtl_style) w = widgets.interactive(update_signals, s_type=w_s_type, s_T=w_s_T, s_t0=w_s_t0, h_type=w_h_type, h_T=w_h_T, h_t0 = w_h_t0) display(widgets.HBox((widgets.VBox((w_s_type, w_s_T, w_s_t0)), widgets.VBox((w_h_type, w_h_T, w_h_t0), layout=Layout(margin='0 0 0 100px'))))); w.update();  %% Cell type:markdown id: tags: ... und betrachte hier das Faltungsergebnis: %% Cell type:code id: tags:  python fig, axs = plt.subplots(2, 1, figsize=(12, 12/ient_fig_aspect)) # gridspec_kw = {'width_ratios':[3, 1]} t_w = widgets.FloatSlider(min=-4, max=4, value=-2, step=.1, description='Verschiebung $t$', style=ient_wdgtl_style) intervals_gt = np.array([t_w.min, t_w.max]) @widgets.interact(t=t_w, show_integrand=widgets.Checkbox(value=True, description='Zeige Integrand', style=ient_wdgtl_style)) def update_plot(t, show_integrand=True): global interval_cnt, intervals_gt t_ind = np.where(tau>=t); t_ind = t_ind[0][0]; g_plot = gt.copy(); g_plot[t_ind:] = np.nan; # hide g(t') with t'>t sh = s(tau)*h(t-tau) # integrand interval_cnt = np.argwhere((intervals_gt[1:] > t) & (intervals_gt[:-1] <= t)); # interval interval_cnt = 0 if interval_cnt.size == 0 else interval_cnt[0][0]; # check for end of interval if not axs[0].lines: # Call plot() and decorate axes. Usually, these functions take some processing time ax = axs[0]; ax.plot(tau, h(t-tau), 'rwth', label=r'$h(t-\tau)$'); # plot h(t-tau) ax.plot(tau, s(tau), 'grun', label=r'$s(\tau)$'); # plot s(tau) ax.plot(tau, sh, '--', color='orange', lw=1, label=r'$s(\tau)h(t-\tau)$'); # plot integrand ient_annotate_xtick(ax, r'$t$', t, -0.1, 'rwth', 15); # mark t on tau axis ax.fill_between(tau, 0, sh, facecolor="none", hatch="//", edgecolor='k', linewidth=0.0); # hatch common area ax.set_xlabel(r'$\rightarrow \tau$'); ax.set_xlim([-4.2,4.2]); ient_update_ylim(ax, np.concatenate((h(tau), s(tau))), 0.19, 5); ax.legend(); ient_grid(ax); ient_axis(ax); ax = axs[1]; ax.plot(tau, g_plot); # plot g(t) ax.plot([t, t], [0, gt[t_ind]], 'ko--', lw=1); ax.set_xlabel(r'$\rightarrow t$'); ax.set_ylabel(r'$\uparrow g(t)=s(t)\ast h(t)$'); ax.set_xlim(axs[0].get_xlim()); ient_update_ylim(ax, gt, 0.19, 5); ient_grid(ax); ient_axis(ax); fig.tight_layout(); #plt.subplots_adjust(wspace=.1,hspace=.1) else: # Replace only xdata and ydata since plt.plot() takes longer time ax = axs[0]; ax.lines[0].set_ydata(h(t-tau)); ax.lines[2].set_ydata(sh); # update signals ax.texts[0].set_x(t); ax.lines[3].set_xdata([t,t]) # update labels if ax.collections: ax.collections[0].remove(); # update integrand if show_integrand: ax.fill_between(tau,0,sh, facecolor="none", hatch="//", edgecolor='k', linewidth=0.0); ax = axs[1]; ax.lines[0].set_ydata(g_plot); # update signals ax.lines[1].set_xdata([t, t]); ax.lines[1].set_ydata([0, gt[t_ind]]); # update labels axs[0].lines[2].set_visible(show_integrand) ax_intervals = axs[1].twinx(); ax_intervals.set_visible(False) def update_plot_intervals(): # Update interval lines ax_intervals.clear(); ax_intervals.axis('off'); for x in intervals_g: ax_intervals.axvline(x, color=rwth_colors['rot'],linestyle='--',lw=1) try: toggle_stepwise(stepwise.value) except: pass update_plot_intervals()  %% Cell type:markdown id: tags: Verschiebung $t$ kann automatisch abgespielt werden. Eine schrittweise Betrachtung ist ebenfalls möglich. %% Cell type:code id: tags:  python def update_t_slider(ti): global interval_cnt if interval_cnt > len(intervals_gt)-2: interval_cnt = 0 tmin = intervals_gt[interval_cnt]; tmax = intervals_gt[interval_cnt+1]; t_w.value = tmin + ti/100*(tmax-tmin) # Update float slider def reset_t_slider(b): global interval_cnt play._playing = False interval_cnt = 0 intslider.value = 0 stepwise.value = False update_t_slider(0) def update_stepwise(args): visible = args['new'] toggle_stepwise(visible) def toggle_stepwise(visible): global intervals_gt ax_intervals.set_visible(visible) if visible: intervals_gt = np.hstack([t_w.min, intervals_g, t_w.max ]) else: intervals_gt = np.array([t_w.min, t_w.max]) # Widget play = widgets.Play(value=0, min=0, max=100,step=10, description="Press play", show_repeat=False) stepwise =widgets.Checkbox(value=False, description='Schrittweise') reset = widgets.Button(description="Reset"); reset.on_click(reset_t_slider); intslider = widgets.IntSlider(description='Fortschritt') # Dummy slider in integer format, to be mapped to float slider status = widgets.Label(value='Nothing') widgets.jslink((play, 'value'), (intslider, 'value')) # Link dummy slider to player #def update_t_limits(args): # status.value = str(np.round(intslider.value)) + " " + str(t_w.value) + " " + str(interval_cnt) + " " + str(args['old']) + " " + str(args['new']); #play.observe(update_t_limits, '_playing') stepwise.observe(update_stepwise, 'value') widgets.interactive(update_t_slider, ti = intslider); intslider.layout.display = 'none'; widgets.HBox([play, stepwise, reset, intslider])#, status  %% Cell type:markdown id: tags: ### Aufgaben: * Bewege den Schieberegler für $t$ und betrachte das entstehende Faltungsintegral. Wie sind die zugehörigen Integralsgrenzen und welche Intervalle (vgl. Notebook zur Faltung) sind zu beobachten? * Wähle zwei Rechtecke unterschiedlicher Breite aus. Wie sieht das entstehende Signal aus? Wie breit ist es? Was passiert, wenn eins der Rechtecke um $t_0$ verschoben wird? * Welche Höhe bei $t=0$ hat das Resultat der Faltung $g(t) = \mathrm{rect}\left(\frac{t}{2}\right)\ast \mathrm{rect}\left(\frac{t}{2}\right)$? Überprüfe die Überlegungen mit Hilfe der entsprechenden Funktionen in der Demo. * Gilt das Kommutativgesetz $s(t) \ast h(t) = h(t) \ast s(t)$? * Wie sieht das Faltungsergebnis zweier si-Funktionen aus? Wie das Ergebnis zweier Gaußfunktionen? * Reale Rechteckimpulse weisen nur eine endliche Flankensteilheit auf. Diese können beispielsweise mit $s(t)=\mathrm{rect}(t)*\mathrm{rect}(t/T)$ oder $s(t)=\mathrm{rect}(t)*\Lambda(\frac{t}{T})$ beschrieben werden. Betrachte diese Fälle für $|T|<\frac{1}{2}$. Wie hängen Gesamtdauer und Dauer der Anstiegsflanke von $T$ ab? %% Cell type:markdown id: tags: This notebook is provided as [Open Educational Resource](https://en.wikipedia.org/wiki/Open_educational_resources) (OER). Feel free to use the notebook for your own purposes. The code is licensed under the [MIT license](https://opensource.org/licenses/MIT). Please attribute the work as follows: *Christian Rohlfing, Übungsbeispiele zur Vorlesung "Grundgebiete der Elektrotechnik 3 - Signale und Systeme"*, gehalten von Jens-Rainer Ohm, 2019, Institut für Nachrichtentechnik, RWTH Aachen University. *Christian Rohlfing, Übungsbeispiele zur Vorlesung "Grundgebiete der Elektrotechnik 3 - Signale und Systeme"*, gehalten von Jens-Rainer Ohm, 2020, Institut für Nachrichtentechnik, RWTH Aachen University. ... ...