Skip to content

Commit

Permalink
bug fix: in getJ2.m, matrices were not reshape right when containing …
Browse files Browse the repository at this point in the history
…multiple trials

refactor: plotSpec.m and plotSpecGram.m now replace the various plotPsd....m and plotCoh....m
  • Loading branch information
Proulx-S committed Jul 24, 2024
1 parent 1cd42ba commit 5f462aa
Show file tree
Hide file tree
Showing 8 changed files with 475 additions and 95 deletions.
40 changes: 26 additions & 14 deletions addWin.m
Original file line number Diff line number Diff line change
Expand Up @@ -13,22 +13,34 @@ function addWin(H,gram,x,y)
Fs = gram.param.Fs;
winSz = gram.win(1) + 0.5/Fs;
if contains(H.Title.String,'timeseries')
yLim = ylim;
xLim = xlim;
if isempty(x)
x = xLim(1); x = x.*[1 1];
xx = x + [0 1].*winSz;
if size(gram.t,2)>1
y = ylim;
xx = gram.t(:,:,:,:,:,:,1);
yy = repmat(y(1),[size(xx,1) 1]);
line(xx,yy,'Color','m','linewidth',3)

xx = gram.t(:,:,:,:,:,:,end);
yy = repmat(y(2),[size(xx,1) 1]);
line(xx,yy,'Color','m','linewidth',3)
else
x = x.*[1 1];
xx = x + [-0.5 0.5].*winSz;
end
if isempty(y) || y==-inf
y = yLim(1);
elseif y==inf
y = yLim(2);
yLim = ylim;
xLim = xlim;
if isempty(x)
x = xLim(1); x = x.*[1 1];
xx = x + [0 1].*winSz;
else
x = x.*[1 1];
xx = x + [-0.5 0.5].*winSz;
end
if isempty(y) || y==-inf
y = yLim(1);
elseif y==inf
y = yLim(2);
end
y = y.*[1 1];
line(xx,y,'Color','m','linewidth',3)
end
y = y.*[1 1];
line(xx,y,'Color','m','linewidth',3)

elseif contains(H.Title.String,'spectrogram') || contains(H.Title.String,'coherogram')
im = findobj(H.Children,'type','image');
if isempty(x)
Expand Down
12 changes: 9 additions & 3 deletions getJ2.m
Original file line number Diff line number Diff line change
Expand Up @@ -12,10 +12,16 @@


tp2 = reshape( tp .* exp(-f.*t*2*pi*1i) ,[N Et*R*K*F*1*Wt]);
d2 = reshape( d - mean(d,1) ,[N E*R*1*1*V*W ]); % removing the mean here
d2 = reshape( d - mean(d,1) ,[N E*R*1*1*V*W]); % removing the mean here

d2 = permute(d2,[2 1]);
J = d2*tp2; % [V E*R*K*F*1*W]
% [E*R*1*1*V*W Et*R*K*F*1*Wt]

J = reshape(J,[V E R K F 1 W]); % [V E R K F 1 W]
J = permute(J,[6 2 3 4 5 1 7]); % [N E R K F V W]
J = reshape(J,[E*R*1*1*V*W Et R K F 1 Wt]);
J = permute(J,[2 3 4 5 6 7 1]);
J = reshape(J,[Et R K F 1 Wt E R 1 1 V W]);
J = permute(J,[13 7 2 3 4 11 6 1 5 8 9 10 12]); %[N ,E ,R ,K,F,V,W]
if any(size(J,9:15)~=1); dbstack; error('X'); end
% J = reshape(J,[V E R K F 1 W]); % [V E R K F 1 W]
% J = permute(J,[6 2 3 4 5 1 7]); % [N E R K F V W]
8 changes: 7 additions & 1 deletion plotPsd2.m
Original file line number Diff line number Diff line change
@@ -1,7 +1,9 @@
function [ax,F] = plotPsd2(volPsd,H,tWin)
if ~exist('H','var'); H = []; end
if isempty(H); H = figure('WindowStyle','docked'); end
if ~exist('tWin','var'); tWin = []; end
if ~exist('tWin','var'); tWin = []; end
if isfield(volPsd,'dsgn') && isfield(volPsd.dsgn,'f0'); f0 = volPsd.dsgn.f0; else; f0 = []; end


switch class(H)
case 'matlab.graphics.layout.TiledChartLayout'
Expand Down Expand Up @@ -52,6 +54,10 @@
addFreq([],volPsd.psdTrialGram.onsetList,volPsd.psdTrialGram.durList)
end

if ~isempty(f0)
xline(f0,'--b')
end

if isfield(volPsd,'psdTrialGramMD') && ~isempty(volPsd.psdTrialGramMD) && ~isempty(tWin)
t = volPsd.psdTrialGramMD.t(1,1,1,1,1,1,:,1);
f = volPsd.psdTrialGramMD.f(1,1,1,1,:,1,1,1);
Expand Down
104 changes: 104 additions & 0 deletions plotSpec.m
Original file line number Diff line number Diff line change
@@ -0,0 +1,104 @@
function [ax,F] = plotSpec(volPsd,metricLabel,H,tWin)
if ~exist('H','var'); H = []; end
if isempty(H); H = figure('WindowStyle','docked'); end

if isfield(volPsd,'dsgn') && isfield(volPsd.dsgn,'f0'); f0 = volPsd.dsgn.f0; else; f0 = []; end

if ~exist('tWin','var'); tWin = []; end
if ~exist('metricLabel','var'); metricLabel = []; end
if isempty(metricLabel); metricLabel = 'psd'; end % 'psd' 'coh'


switch class(H)
case 'matlab.graphics.layout.TiledChartLayout'
F = H.Parent;
case 'matlab.ui.Figure'
F = H;
otherwise
end

figure(F);
ax = {};
ax{end+1} = nexttile;


%% Select approrpiate data
switch metricLabel
case 'psd'
mt = volPsd.psd;
vec = mt.PSD;
label = 'spatially averaged spectrum';
case 'coh'
mt = volPsd.svd;
vec = mt.COH;
label = 'coherence spectrum';
otherwise
dbstack; error('code trhat')
end




%%% average across space
vec = squeeze(mean(vec,6));
f = squeeze(mt.f);
% f = squeeze(volPsd.f);
% psd = squeeze(mean(volPsd.vec,2));
plot(f,vec,'k')
grid on
axis tight
xlabel('f (Hz)')
switch metricLabel
case 'psd'
ylabel('psd')
ax{end}.YScale = 'log';
case 'coh'
ylabel('coherence')
otherwise
dbstack; error('code trhat')
end


K = mt.K;
T = mt.T;
[TW,W] = K2W(T,K,0);

paramStr = ['(K=' num2str(K) '; 2W=' num2str(W*2,'%0.4f') 'Hz; T=' num2str(T,'%0.2f') 'sec); TW=' num2str(TW)];
title(['spatially averaged spectrum ' paramStr])

yLim = vec(f>0.01);
yLim = [min(yLim) max(yLim)];
ylim(yLim)
xlim([0 f(end)])

addW([],volPsd.psd)



if isfield(volPsd,'psdTrialGram') && ~isempty(volPsd.psdTrialGram)
addFreq([],volPsd.psdTrialGram.onsetList,volPsd.psdTrialGram.durList)
end

if ~isempty(f0)
xline(f0,'--b')
end

if isfield(volPsd,'psdTrialGramMD') && ~isempty(volPsd.psdTrialGramMD) && ~isempty(tWin)
dbstack; error('double-check that');
t = volPsd.psdTrialGramMD.t(1,1,1,1,1,1,:,1);
f = volPsd.psdTrialGramMD.f(1,1,1,1,:,1,1,1);
if tWin==inf
psd = mean(mean(volPsd.psdTrialGramMD.vec.psdPC(:,:,:,:,:,:,:,:),6),7);
else
[~,wInd] = min(abs(t-tWin));
psd = mean(volPsd.psdTrialGramMD.vec.psdPC(:,:,:,:,:,:,wInd,:),6);
end
plot(squeeze(f),squeeze(psd),'--k')
end


ax = [ax{:}];




Loading

0 comments on commit 5f462aa

Please sign in to comment.