User Tools

Site Tools


analysis:amplipex:phaseplot

Differences

This shows you the differences between two versions of the page.

Link to this comparison view

analysis:amplipex:phaseplot [2014/01/21 17:34]
mbonsma created
analysis:amplipex:phaseplot [2018/07/07 10:19]
Line 1: Line 1:
-<code matlab> 
  
-%% channels specified by Eric 
- 
- 
-cd('​C:​\Users\mvdmlab\Documents\Data\promoted\R036\R036-2013-05-23'​);​ 
-fname = '​R036-2013-05-23-pre.dat'; ​ %this is some of ERic's data 
- 
-%%  
-channels = [33,​57,​1,​25,​36,​58,​2,​26,​38,​62,​6,​30,​37,​61,​5,​29];​ % specify channels in the APX data file to load 
-%% 
-data = AMPX_loadData(fname,​channels_to_load,​2);​ % note decimation factor of 1 
- 
-%% frequency ranges of interest 
- 
-% low gamma: 45-65 Hz 
-lgamma = [45 65]; 
-stopbandgamm = 15; 
-Fs = data.hdr.Fs;​ 
- 
-%% filter for low gamma 
-Wp2 = lgamma * 2 / Fs;  
-Ws2 = [lgamma(1)-stopbandgamm lgamma(2)+stopbandgamm] * 2 / Fs; 
-[N2,Wn2] = cheb1ord( Wp2, Ws2, 3, 20);  % use a chebyshev because we want a steep cutoff and can sacrifice the passband quality 
-[b2,a2] = cheby1(N2,​2,​Wn2);​ 
- 
-[N2b,Wn2b] = buttord( Wp2, Ws2, 3, 20); % try a butterworth filter too. 
-[b2b,a2b] = butter(N2b,​Wn2b); ​ 
- 
-fvtool(b2,​a2,​b2b,​a2b) 
- 
-%use the butterworth here for a cleaner passband 
-%% filter data 
-ygamm = cell(1,​length(channels));​ 
- 
-%% 
-for ii = 1:​length(channels);​ 
-    ygamm{ii} = filtfilt(b2b,​a2b,​data.channels{ii});​ 
-end 
- 
-%% plot to check (optional) 
-figure; 
-plot(data.tvec,​data.channels{1},'​b'​);​ 
-hold on; 
-plot(data.tvec,​ygamm{1},'​r'​);​ 
- 
-plot(data.tvec,​data.channels{2},'​m'​);​ 
-plot(data.tvec,​data.channels{3},'​g'​);​ 
-plot(data.tvec,​data.channels{4},'​k'​);​ 
-set(gca,'​XLim',​[124.95 125]); 
- 
- 
-%% get phase 
-phi = cell(1,​length(channels));​ 
-for ii = 1:​length(channels);​ 
-    phi{ii} = angle(hilbert(ygamm{ii}));​ 
-end 
- 
-%% make nice fancy plot - this one plots filtered signal and phase 
-xwin = [124.8 125]; 
- 
-figure; 
-[AX,~,~] = plotyy(data.tvec,​(ygamm{1}-mean(ygamm{1})+1*2),​data.tvec,​(phi{1}-mean(phi{1})+1*20));​ 
-freq = AX(1); 
-ph = AX(2); 
-hold all; 
-set(AX,'​XLim',​xwin);​ 
-hold(AX(1),'​on'​);​ hold(AX(2),'​on'​);​ 
- 
-for ii = 2:​length(channels);​ 
-    plot(freq,​data.tvec,​ygamm{ii}-mean(ygamm{ii})+ii*2);​ 
-    plot(ph,​data.tvec,​phi{ii}-mean(phi{ii})+ii*20);​ 
-    hold(freq,'​on'​);​ hold(ph,'​on'​);​ 
-    count = ii; 
-end 
-set(AX(1),'​YLim',​[0,​count*2+2]);​ 
-set(AX(2),'​YLim',​[0,​count*20+20]);​ 
- 
-%% make data for colormap testing 
-i1 = find(data.tvec >= xwin(1)); 
- 
-%% 
-y = cell(1,​length(channels));​ 
-cdata = cell(1,​length(channels));​ 
- 
-x = data.tvec(i1(1):​i1(1000));​ 
- 
-for ii = 1:​length(channels) 
-    y{ii} = ygamm{ii}(i1(1):​i1(1000));​ 
-    cdata{ii} = phi{ii}(i1(1):​i1(1000));​ 
-end 
- 
-%% plot filtered signal with phase as a colour 
- 
-scalefactor = 2; 
- 
-figure; 
-colormap(hsv);​ 
- 
-for ii = 1:​length(channels);​ 
-    ​ 
-    scatter(x,​y{ii}-mean(y{ii})+ii*scalefactor,​10,​cdata{ii},'​fill'​);​ 
-    ​ 
-    hold on; 
-    count = ii; 
-end 
-ytick = [count:​scalefactor:​count*scalefactor];​ 
-set(gca,'​YLim',​[0,​count*scalefactor+scalefactor],'​YTick',​ytick,'​YTickLabel',​channels);​ 
- 
- 
-</​code>​ 
analysis/amplipex/phaseplot.txt ยท Last modified: 2018/07/07 10:19 (external edit)