Source code :: main

[Return]
[Download]#!/usr/bin/env python # -*- coding: utf-8 -*- import matplotlib matplotlib.rcParams['backend'] = 'Agg' from scipy.signal import * from numpy import * from matplotlib.pyplot import * import numpy as np from cwt import cwt from matplotlib.ticker import MultipleLocator, FormatStrFormatter,LogLocator import os matplotlib.rcParams['xtick.direction'] = 'out' matplotlib.rcParams['ytick.direction'] = 'out' matplotlib.rcParams['xtick.major.size'] = 10 matplotlib.rcParams['xtick.minor.size'] = 7 matplotlib.rcParams['ytick.major.size'] = 10 matplotlib.rcParams['ytick.minor.size'] = 7 def Spectrogram(omega0,frequenciCutOff,signal,name,start, end): print name tvec = signal[:,0] startAdv = max(start*0.8, tvec[0]) endAdv = min((end-start)*1.2+start, tvec[-1]) dt = tvec[1]-tvec[0] N = len(tvec) for i in range(3): filtered = medfilt(signal[:,1], floor(N/100)+1) res = abs(signal[:,1] - filtered) ind = argsort(res) n_out = floor(N/400) ind_2 = in1d(range(N), ind[-n_out:]) & (res > 0.2*amax(filtered)) signal[ ind_2, 1] = filtered[ind_2] #remove "n_out" most distant points from moving median spec, scale = cwt(signal[:,1], dt, dj=0.05, wf='morlet', p=omega0, extmethod='none', extlength='powerof2') # approximate scales through frequencies freq = (omega0 + np.sqrt(2.0 + omega0**2))/(4*np.pi * scale[1:]) contrast = 10 #spec = spec[:,(tvec > startAdv )*( tvec <endAdv)] pole = log(1+contrast*np.abs(spec)) pole = pole[freq > frequenciCutOff, : ] freq = freq[freq > frequenciCutOff] fig = figure() ax = fig.add_axes([0.1, 0.1, 0.8, 0.85]) img = ax.imshow(pole, extent=[startAdv*1000,endAdv*1000 ,freq[-1], freq[0]], aspect='auto') ax.set_yscale('log', nonposy='clip') minorLocator = MultipleLocator(1) ax.xaxis.set_minor_locator(minorLocator) ax.plot([start*1000, (start)*1000], [amin(freq), amax(freq)], 'w--') ax.plot([end*1000, (end)*1000], [amin(freq), amax(freq)], 'w--') ax.axis([startAdv*1000,endAdv*1000,amin(freq), amax(freq)]) ax.set_xlabel('time [ms]') ax.set_ylabel('Frequency [Hz]') savefig('spectrogram_'+name+'.png') clf() signal= signal[(tvec > startAdv )*( tvec<endAdv),1 ] plot(1000*tvec[(tvec > startAdv )*( tvec<endAdv) ],signal,'k',linewidth=0.1) #axvline(linewidth=4) minimum = -max(abs(amin(signal)), abs(amax(signal))) maximum = -minimum axis([startAdv*1000,endAdv*1000,minimum, maximum]) plot([start*1000, start*1000], [minimum, maximum], 'r--') plot([end*1000, end*1000], [minimum, maximum], 'r--') xlabel('time [ms]') ylabel('Intensity [a.u.]') #show() savefig('signal_'+name+'.png') clf() start = loadtxt('PlasmaStart') end = loadtxt('PlasmaEnd') omega0 = 15 #čím větší tím lepší frekvenční rozlišení a opačně frequenciCutOff = 900 signal = loadtxt('NIdata.lvm') for i in range(1,size(signal,1)): try: name = str(i) Spectrogram(omega0,frequenciCutOff,signal[:,[0,i]],name,start,end) except: print 'Failed mirnov' raise os.system('convert -resize 150 spectrogram_1.png icon.png') #try: #start = 0 #end = 50 #omega0 = 30 #frequenciCutOff = 900 #signal = loadtxt('comedi_spectra.txt') #signal[:,1] -=mean(signal[:,1]) #name = 'comedi' #Spectrogram(omega0,frequenciCutOff,signal,name,start,end) #except: #print 'Failed comedi' #raise [Return]

Navigation