{ "cells": [ { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "%matplotlib inline\n", "from matplotlib import pyplot as plt\n", "from skimage.color import rgb2hed, hed2rgb\n", "import matplotlib as mp\n", "import numpy as np\n", "import math as mt\n", "\n", "from pylab import *\n", "from scipy import *\n", "\n", "from ipywidgets import *\n", "\n", "%matplotlib inline \n", "\n", "# ROZWAZANY PRZEDZIAL ([T0, T1])\n", "signalScale = [0.0, 4.0]\n", "# WYSWIETLANA SKALA NA WYKRESIE (DZIEDZINA CZASU)\n", "timePlotXScale = [0.0, 4.0]\n", "# WYSWIETLANA SKALA Y (DZIEDZINA CZASU)\n", "timePloyYScale = [-3.0, 3.0] \n", "\n", "def changeX(left = 10, right = 20, speed = 10.0):\n", " def dynamicX(x):\n", " result = x\n", " if x >= left and x <= right:\n", " angle = (x - left)/(right-left)*(mt.pi)\n", " shift = (mt.cos(angle + mt.pi) + 1.0) / 2.0\n", " result = result + shift * speed\n", " if x > right:\n", " result = result + 1.0 * speed\n", " return result\n", " return dynamicX\n", "\n", "def changeY(left = 5, right = 15, amplitude = 1.0):\n", " def dynamicY(y, x):\n", " result = y\n", " if x >= left and x <= right:\n", " angle = (x - left)/(right-left)*(2.0*mt.pi)\n", " shift = (mt.cos(angle + mt.pi) + 1.0) / 2.0\n", " result = y + shift * amplitude\n", " return result\n", " return dynamicY\n", "\n", "\n", "\n", "# KLASA SCENARIO\n", "# PO PROSTU LISTA SIGNALOW \n", "class Scenario(Widget):\n", " def __init__(self, signals):\n", " self.signals = signals\n", " \n", "# KLASA SYGNAL - PROSTY(SIN)\n", "class Signal:\n", " def __init__(self, amplitude = 1.0, frequency = 0.5, phase = 0.0, yTranslation = 0.0, color='blue', \n", " changeX=changeX(-1,-1,0), changeY=changeY(-1,-1,0)):\n", " self.amplitude = amplitude\n", " self.frequency = frequency\n", " self.phase = phase\n", " self.yTranslation = yTranslation\n", " self.color = color\n", " self.changeX = changeX\n", " self.changeY = changeY\n", " \n", " #LISTA X, Y DLA ZADANEGO PRZEDZIALU I CZESTOTLIWOSCI PROBKOWANIA\n", " def getDataForRange(self, scale, samplingFrequency):\n", " x = np.arange(scale[0], scale[1], 1.0 / samplingFrequency)\n", " fs = lambda a: (self.amplitude * mt.sin(self.phase + self.changeX(a) * self.frequency * 2.0 * mt.pi) + self.yTranslation)\n", " y = [self.changeY(fs(a), a) for a in x ]\n", " return x, y\n", "\n", "# WYLICZA WYNIKOWY SYGNAL ZE SKLADOWYCH\n", "def getFinalSignal(signals, samplingFrequency):\n", " x = []\n", " y = []\n", " for s in signals:\n", " data = s.getDataForRange(signalScale, samplingFrequency)\n", " if len(x) == 0: x = data[0]\n", " if not y: y = [0] * len(x)\n", " for i in range(len(x)):\n", " y[i] = y[i] + data[1][i]\n", " return x, y\n", " \n", "# --- RYSOWANIE WSZYSTKIEGO\n", "# --- SINGLE - NA WEJSCIU PARAMETRY DLA POJEDYNCZEGO SYGNALU\n", "# --- ALL - NA WEJSCIU SCENARIO\n", "\n", "\n", "\n", "def plotSingle(\n", " amplitude=1.0, \n", " frequency=0.5, \n", " phase=0.0, \n", " yTranslation=0.0,\n", " samplingFrequency=10.0, \n", " showComponents=True, \n", " showFinal=False, \n", " showFrequencyDomain=False, \n", " showSymmetrical=False,\n", " showInverseFFT=False):\n", " \n", " scenario = Scenario([Signal(amplitude, frequency, phase, yTranslation)])\n", " plotAll(scenario, samplingFrequency, showComponents, showFinal, showFrequencyDomain, showSymmetrical, showInverseFFT)\n", " \n", "def plotAllWithNoiseCancel(\n", " scenario = Scenario([Signal(amplitude=1.0, frequency=0.5, phase=0.0, yTranslation=0.0, color='blue')]),\n", " samplingFrequency=10.0, \n", " noiseThreshold=10000.0,\n", " showComponents=True, \n", " showFinal=False, \n", " showFrequencyDomain=False, \n", " showSymmetrical=False, \n", " showInverseFFT=False,\n", " showSpectrogram=False,\n", " windowSize=2.0,\n", " windowJump=0.5):\n", " \n", " signals = scenario.signals\n", " \n", " # --- CONSTRUCT FIGURE AND SUBPLOTS ------------------------\n", " if showFrequencyDomain or showInverseFFT or showSpectrogram:\n", " column_width_pt = 1000.0\n", " else:\n", " column_width_pt = 500.0\n", " \n", " pt_per_inch = 72\n", " size = column_width_pt / pt_per_inch;\n", " height = 0.55\n", " \n", " total = 1 + sum([showFrequencyDomain, showInverseFFT, showSpectrogram])\n", " if total in [2, 3]: height = 0.37\n", "\n", " fig = plt.figure(1, figsize=(size, height * size))\n", " sub = fig.add_subplot(1, total, 1)\n", " if showFrequencyDomain:\n", " subFreq = fig.add_subplot(1, total, 2)\n", " if showInverseFFT: \n", " subIFFT = fig.add_subplot(1, total, 2 if not showFrequencyDomain else 3)\n", "\n", " if showSpectrogram:\n", " subSpectr = fig.add_subplot(1, total, total)\n", " # ------------------------------------------------------\n", " # --- TIME DOMAIN --------------------------------------\n", " \n", " sub.set_xlabel('Time [s]')\n", " sub.set_ylabel('y = signal(x)')\n", " sub.grid('on', axis='both', color='gray', linewidth=1.25)\n", " \n", " sub.set_xlim(timePlotXScale)\n", " sub.set_ylim(timePloyYScale)\n", " \n", " if showComponents:\n", " for s,_ in enumerate(signals):\n", " data = signals[s].getDataForRange(signalScale, samplingFrequency)\n", " #ŻEBY BYLO PRZEJRZYSCIEJ\n", " if samplingFrequency < 8.0:\n", " sub.plot(data[0], data[1], marker='o', linestyle='---', color = signals[s].color)\n", " else:\n", " sub.plot(data[0], data[1], linestyle='-', color = signals[s].color)\n", " \n", " dataFinal = getFinalSignal(signals, samplingFrequency)\n", " if showFinal:\n", " dataFinal = getFinalSignal(signals, samplingFrequency)\n", " sub.plot(dataFinal[0], dataFinal[1], linestyle='-', linewidth = 2.0, color = 'red')\n", " \n", " # ------------------------------------------------------\n", " # --- FREQUENCY DOMAIN --------------------------------------\n", " \n", " ind = np.arange(len(dataFinal[0])) # the x locations for the groups \n", " signal1 = np.fft.fft(dataFinal[1])\n", " signal2 = abs(signal1) \n", "\n", " ffty = signal2\n", " freq = [v * samplingFrequency / len(dataFinal[0]) for v in ind]\n", " \n", " # USUWAMY PASMO\n", " signal3 = []\n", " for i in range(len(freq)):\n", " if noiseThreshold <= samplingFrequency / 2:\n", " if freq[i] <= samplingFrequency / 2 and freq[i] >= noiseThreshold:\n", " ffty[i] = 0.0\n", " signal1[i] = 0.0\n", " if freq[i] >= samplingFrequency / 2 and freq[i] <= samplingFrequency - noiseThreshold:\n", " ffty[i] = 0.0\n", " signal1[i] = 0.0\n", " \n", " if showSymmetrical == False:\n", " freq = freq[:len(signal2)//2]\n", " ffty = ffty[:len(signal2)//2]\n", " \n", " if showFrequencyDomain == True:\n", " \n", " subFreq.set_xlabel('Frequency [Hz]')\n", " subFreq.set_ylabel('Amplitude')\n", " subFreq.grid('on', axis='both', color='gray', linewidth=1.25)\n", " \n", " #TRANSFORM\n", " subFreq.set_ylim([0.0, 80.0])\n", " white_space = freq[-1]/len(freq)*2\n", " subFreq.set_xlim([-white_space, freq[-1]+white_space])\n", " barFreq = subFreq.stem(freq, ffty, '-*', color='blue')\n", " \n", " # ------------------------------------------------------\n", " # --- INVERSE FFT --------------------------------------\n", " if showInverseFFT == True:\n", " subIFFT.set_xlabel('Time [s]')\n", " subIFFT.set_ylabel('y = signal(x)')\n", " subIFFT.grid('on', axis='both', color='gray', linewidth=1.25)\n", " \n", " subIFFT.set_xlim(timePlotXScale)\n", " subIFFT.set_ylim(timePloyYScale)\n", " \n", " signalIFFT = np.fft.ifft(signal1)\n", " \n", " \n", " subIFFT.plot(dataFinal[0], np.real(signalIFFT), linestyle='-', linewidth = 2.0, color = 'red')\n", " \n", " # ------------------------------------------------------\n", " # --- SHOW SPECTROGRAM --------------------------------------\n", " if showSpectrogram == True:\n", " subSpectr.set_xlabel('Time [s]')\n", " subSpectr.set_ylabel('Frequency [Hz]')\n", " subSpectr.grid('on', axis='both', color='gray', linewidth=1.25)\n", " \n", " xlabel = []\n", " data = []\n", " \n", " start = 0;\n", " startInDomain = 0\n", " windowSizeInDomain = windowSize * samplingFrequency\n", " jumpInDomain = windowJump * samplingFrequency\n", " \n", " while startInDomain + windowSizeInDomain <= len(dataFinal[0]):\n", " \n", " xlabel.append(start)\n", " portion = []\n", " for i in range(int(startInDomain), int(startInDomain + windowSizeInDomain)):\n", " portion.append(dataFinal[1][i])\n", " \n", " toAdd = abs(np.fft.fft(portion))\n", " if showSymmetrical == False:\n", " toAdd = toAdd[0:int(len(toAdd)/ 2)]\n", " \n", " data.append(toAdd)\n", " \n", " start += windowJump\n", " startInDomain += jumpInDomain\n", " \n", " \n", " npdata = np.array(data)\n", " npdata = np.transpose(npdata)\n", " \n", " subSpectr.set_xlim(0, len(npdata[0]))\n", " subSpectr.set_ylim(0, len(npdata))\n", " \n", " heatmap = subSpectr.pcolor(npdata)\n", " fig.colorbar(heatmap)\n", " \n", " plt.show()" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": true }, "outputs": [], "source": [ "def plotSpectrogram(scenario = Scenario([Signal(amplitude = 1.0, frequency = 0.5, phase = 0.0, yTranslation = 0.0, color='blue')]),\n", " samplingFrequency = 10.0, \n", " showComponents = True, \n", " showFinal=False, \n", " showFrequencyDomain=False, \n", " showSymmetrical=False, \n", " showSpectrogram=False,\n", " windowSize = 2.0,\n", " windowJump = 0.5):\n", " \n", " plotAllWithNoiseCancel(scenario, samplingFrequency, 1000000.0, showComponents, showFinal, \n", " showFrequencyDomain, showSymmetrical, False, showSpectrogram,\n", " windowSize,\n", " windowJump)\n", " \n", " \n", "def plotAllWithNoiseCancelWithoutSpectrogram(scenario = Scenario([Signal(amplitude = 1.0, frequency = 0.5, phase = 0.0, yTranslation = 0.0, color='blue')]),\n", " samplingFrequency = 10.0, \n", " noiseThreshold = 10000.0,\n", " showComponents = True, \n", " showFinal=False, \n", " showFrequencyDomain=False, \n", " showSymmetrical=False, \n", " showInverseFFT=False):\n", " \n", " plotAllWithNoiseCancel(scenario, samplingFrequency, noiseThreshold, showComponents, showFinal, \n", " showFrequencyDomain, showSymmetrical, showInverseFFT, showSpectrogram=False,\n", " windowSize = 2.0,\n", " windowJump = 0.5)\n", " \n", "# --- WA - BEZ CZESCI PARAMETROW TRUE/FALSE - INACZEJ NIE WIEM JAK WYWALIC WIDGET Z INTERACT\n", "# --- WIFFT - BEZ PARAMETROW IFFT - INACZEJ NIE WIEM JAK WYWALIC WIDGET Z INTERACT\n", "def plotAll(scenario = Scenario([Signal(amplitude = 1.0, frequency = 0.5, phase = 0.0, yTranslation = 0.0, color='blue')]),\n", " samplingFrequency = 10.0, \n", " showComponents = True, \n", " showFinal=False, \n", " showFrequencyDomain=False, \n", " showSymmetrical=False, \n", " showInverseFFT=False): \n", " \n", " plotAllWithNoiseCancel(scenario, samplingFrequency, 10000.0, showComponents, showFinal, \n", " showFrequencyDomain, showSymmetrical, showInverseFFT, showSpectrogram=False,\n", " windowSize = 10,\n", " windowJump = 3)\n", "\n", "def plotSingleWithoutAll( \n", " amplitude = 1.0, \n", " frequency = 0.5, \n", " phase = 0.0, \n", " yTranslation = 0.0,\n", " samplingFrequency = 10.0):\n", " \n", " scenario = Scenario([Signal(amplitude, frequency, phase, yTranslation)])\n", " plotAll(scenario, samplingFrequency)\n", "\n", "def plotSingleWIFFT(\n", " amplitude = 1.0, \n", " frequency = 0.5, \n", " phase = 0.0, \n", " yTranslation = 0.0,\n", " samplingFrequency = 10.0, \n", " showComponents = True, \n", " showFinal=False, \n", " showFrequencyDomain=False, \n", " showSymmetrical=True):\n", " \n", " scenario = Scenario([Signal(amplitude, frequency, phase, yTranslation)])\n", " plotAll(scenario, samplingFrequency, showComponents, showFinal, showFrequencyDomain, showSymmetrical)\n", " \n", "def plotAllWithoutAll( \n", " scenario = Scenario([Signal(amplitude = 1.0, frequency = 0.5, phase = 0.0, yTranslation = 0.0, color='blue')]),\n", " samplingFrequency = 10.0, \n", " showComponents = True, \n", " showFinal=False):\n", " \n", " plotAll(scenario, samplingFrequency,showComponents, showFinal) " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# SYGNAL W DZIEDZINIE CZASU" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Podstawowe cechy, czestotliwosc probkowania, nyquista" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "interact(plotSingleWithoutAll, amplitude=(0.5,2.5,0.1),\n", " frequency=(0.1, 2.0,0.1), phase=(-2*mt.pi, 2*mt.pi, 0.25), \n", " yTranslation=(-2.5, 2.5, 0.5), samplingFrequency=(0.25,10,0.25));\n", "\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Sumowanie sygnałow" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "sc1 = Scenario([\n", " Signal(amplitude = 1.0, frequency = 0.5, phase = 0.0, yTranslation = 0.0, color = 'blue'),\n", " Signal(amplitude = 2.0, frequency = 0.25, phase = 0.0, yTranslation = 0.0, color = 'green') ])\n", "sc2 = Scenario([\n", " Signal(amplitude = 1.0, frequency = 0.5, phase = 0.0, yTranslation = 0.0, color = 'blue'),\n", " Signal(amplitude = 1.25, frequency = 0.5, phase = 0.0, yTranslation = 0.0, color = 'green') ])\n", "sc3 = Scenario([\n", " Signal(amplitude = 1.0, frequency = 0.5, phase = 0.0, yTranslation = 0.0, color = 'blue'),\n", " Signal(amplitude = 1.0, frequency = 0.5, phase = mt.pi, yTranslation = 0.0, color = 'green') ])\n", "\n", "\n", "\n", "scenarios = {\"1. Przyklad\": sc1 , \"2. Wzmocnienie\" : sc2, \"3. Tlumienie\" : sc3}\n", "interact(plotAllWithoutAll, scenario=scenarios, samplingFrequency=(0.25,10,0.25));" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# SYGNAŁ W DZIEDZINIE CZĘSTOTLIWOŚCI" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## FFT" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false, "scrolled": true }, "outputs": [], "source": [ "#TODO: Po co tu showFinal?\n", "interact(plotSingleWIFFT, amplitude=(0.5,2.5,0.1),\n", " frequency=(0.05, 2.0, 0.05), phase = (-2*mt.pi, 2*mt.pi, 0.25), \n", " yTranslation=(-2.5, 2.5, 0.5), samplingFrequency=(0.25, 10, 0.25));" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Suma sygnalow w dziedzinie czasu a efekt w dziedzinie czestotliwosci oraz odwrotne FFT" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "sc1 = Scenario([\n", " Signal(amplitude = 1.0, frequency = 0.5, phase = 0.0, yTranslation = 0.0, color = 'blue'),\n", " Signal(amplitude = 2.0, frequency = 0.25, phase = 0.0, yTranslation = 0.0, color = 'green') ])\n", "sc2 = Scenario([\n", " Signal(amplitude = 1.0, frequency = 0.5, phase = 0.0, yTranslation = 0.0, color = 'blue'),\n", " Signal(amplitude = 1.25, frequency = 0.5, phase = 0.0, yTranslation = 0.0, color = 'green') ])\n", "sc3 = Scenario([\n", " Signal(amplitude = 1.0, frequency = 0.5, phase = 0.0, yTranslation = 0.0, color = 'blue'),\n", " Signal(amplitude = 1.0, frequency = 0.5, phase = mt.pi, yTranslation = 0.0, color = 'green') ])\n", "\n", "scenarios = {\"1. Przyklad\": sc1 , \"2. Wzmocnienie\" : sc2, \"3. Tlumienie\" : sc3}\n", "interact(plotAll, scenario=scenarios, samplingFrequency=(0.25,10,0.25));\n", "\n", "#TODO: Po co tutaj inverseFFT?" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Proste filtrowanie (FFT i iFFT)" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "l1 = [\n", " Signal(amplitude=1.0, frequency=0.5, phase=0.0, yTranslation = 0.0, color = 'blue'),\n", " Signal(amplitude=1.5, frequency=0.25, phase=0.5 * mt.pi, yTranslation = 0.0, color = 'green'),\n", " Signal(amplitude=0.5, frequency=0.4, phase=1.25 * mt.pi, yTranslation = 0.0, color = 'orange')]\n", "\n", "l2 = [\n", " Signal(amplitude=1.0, frequency=0.5, phase=0.0, yTranslation = 0.0, color = 'blue'),\n", " Signal(amplitude=1.5, frequency=0.25, phase=0.5 * mt.pi, yTranslation = 0.0, color = 'green'),\n", " Signal(amplitude=0.5, frequency=0.4, phase=1.25 * mt.pi, yTranslation = 0.0, color = 'orange')]\n", "\n", "for i in range(50):\n", " amplitude = 0.05 + random.random() * 0.15\n", " #Round - Uniknięcie błędów numerycznych, co sie stanie po usunieciu round?\n", " # frequency = round(1.0 + random.random() * 3.0, 1)\n", " frequency = round(1.0 + random.random() * 3.0, 1) \n", " phase = random.random() * mt.pi\n", " color = 'black'\n", " l2.append(Signal(amplitude, frequency, phase, 0.0, color))\n", "\n", "sc1 = Scenario(l1)\n", "sc2 = Scenario(l2)\n", "\n", "scenarios = {\"1. Przyklad\": sc1, \"2. Noise\": sc2 }\n", "interact(plotAllWithNoiseCancelWithoutSpectrogram, scenario=scenarios, noiseThreshold=6.0, samplingFrequency=(0.25,10,0.25));" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Pojedynczy sygnal - zmienny w czasie" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "sc1 = Scenario([\n", " Signal(amplitude = 1.0, frequency = 0.5, phase = 0.0, yTranslation = 0.0, color = 'blue', \n", " changeX=changeX(5,15,10))])\n", "sc2 = Scenario([\n", " Signal(amplitude = 1.0, frequency = 0.5, phase = 0.0, yTranslation = 0.0, color = 'blue', \n", " changeY=changeY(5,15,1))])\n", "sc3 = Scenario([\n", " Signal(amplitude = 1.0, frequency = 0.5, phase = 0.0, yTranslation = 0.0, color = 'blue', \n", " changeX=changeX(10,20,10), changeY=changeY(0,10,1))])\n", "\n", "\n", "scenarios = {\"1. Zmienna czestotliwosc\": sc1 , \"2. Zmienna amplituda\" : sc2, \"3. Mix\" : sc3}\n", "interact(plotAll, scenario=scenarios, samplingFrequency=(0.25,10,0.25));" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Spectrogram" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "\n", "sc1 = Scenario([\n", " Signal(amplitude = 1.0, frequency = 0.5, phase = 0.0, yTranslation = 0.0, color = 'blue', \n", " changeX=changeX(5,15,10))])\n", "sc2 = Scenario([\n", " Signal(amplitude = 1.0, frequency = 0.5, phase = 0.0, yTranslation = 0.0, color = 'blue', \n", " changeY=changeY(5,15,1))])\n", "sc3 = Scenario([\n", " Signal(amplitude = 1.0, frequency = 0.5, phase = 0.0, yTranslation = 0.0, color = 'blue', \n", " changeX=changeX(10,20,10), changeY=changeY(0,10,1))])\n", "\n", "scenarios = {\"1. Czestotliwosc\": sc1, \"2. Aplituda\": sc2, \"3. Mix\": sc3 }\n", "interact(plotSpectrogram, scenario=scenarios, samplingFrequency=(0.25,10,0.25), windowSize=(4.0,10.0,0.25), \n", " windowJump=(0.25,2.0,0.25));\n" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.5.0" } }, "nbformat": 4, "nbformat_minor": 0 }