V1-V4-LFPs-and-Visual-Attention/Code/processSpectralPower.m

addpath('./support_routines');
addpath(genpath('./support_tools'));

clear all; close all; clc;

sSubjectName='monkey 2';
sVisualArea='V4';
subjectInfos=getSubjectInfos(sSubjectName);
subjectData=getSubjectData(sSubjectName,sVisualArea,'ALL','LFPs','ALL');

SAMPLING_FREQ=1017.38;
NUM_TIME_SAMPLES=512;
NUM_FREQ_SAMPLES=NUM_TIME_SAMPLES/2+1;

timeAxisPostCue=linspace(0,NUM_TIME_SAMPLES/SAMPLING_FREQ,NUM_TIME_SAMPLES);
timeAxisPreStim=linspace(-NUM_TIME_SAMPLES/SAMPLING_FREQ,0,NUM_TIME_SAMPLES);

lamLayerTag={'Supra','Granr','Infra'};
gratCondTag={'RF','OUT1','OUT2'};
spectWinTag={'Full'};

plotPmtsSinglePens=0;
plotPmtsPenAvg=0;
plotPmtsPenAvgOutMean=1;
plotPmtsPenAvgColorCoded=0;


for reportPmtMinusMean=1%0:1 
    % reportPmtMinusMean allows to choose how to apply baseline normalization 
    % 0 reports (Pmt) ./ (MeanPmtbaseline) 
    % 1 reports (Pmt - MeanPmtbaseline)./(stdPmtbaseline)
    
    %% compute Pmts
    paramsMT=[];
    paramsMT.Fs=SAMPLING_FREQ;
    paramsMT.tapers = [2 3];
    paramsMt.pad = 0;
    [~,freqAxisPreStim] = mtspectrumc(timeAxisPreStim*1e-3, paramsMT);
    [~,freqAxisPostCue] = mtspectrumc(timeAxisPostCue*1e-3, paramsMT);
    
    pr=1;
    for pp=1:length(subjectData.penIDs) % LOOP OVER PENS
        if ~isempty(subjectData.LfpStruct(pp).Sorted)
            for bb=1%:4 % LOOP OVER SPECTRAL WINS
                for ll=1:3 % LOOP OVER LAMINAR LAYERS
                    for cnd=1:3 % LOOP OVER GRAT CONDITIONS
                        if bb==1
                            currLFPsPreStim=subjectData.LfpStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PreStimDataBi(:,:,end-NUM_TIME_SAMPLES+1:end);
                            currLFPsPostStim=subjectData.LfpStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PostStimDataBi(:,:,1:NUM_TIME_SAMPLES);
                            currLFPsStationary=subjectData.LfpStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).StationaryDataBi(:,:,1:NUM_TIME_SAMPLES);
                            currLFPsPostCue=subjectData.LfpStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PostCueDataBi(:,:,1:NUM_TIME_SAMPLES);
                            currLFPsPreFirstDim=subjectData.LfpStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PreFirstDimDataBi(:,:,end-NUM_TIME_SAMPLES+1:end);
                            
                            % SELECT GRATCONDS 2, 3, 4, 5, 9, 11, 15, 17, 20, 21, 22, 23, 27, 29, 33, 35
                            if ~isempty(subjectData.LfpStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PreSecondDimDataBi)
                                currNumTrials2=length(subjectData.LfpStruct(pp).Sorted.Full.(lamLayerTag{ll}).(gratCondTag{cnd}).gratCondSecondDim);
                                currCorrectSecondDimTrials=find(arrayfun(@(jj) any([2 3 4 5 20 21 22 23 9 11 15 17 27 29 33 35]==subjectData.LfpStruct(pp).Sorted.Full.(lamLayerTag{ll}).(gratCondTag{cnd}).gratCondSecondDim(jj)), 1:currNumTrials2));
                                currLFPsPreSecondDim=subjectData.LfpStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PreSecondDimDataBi(:,currCorrectSecondDimTrials,end-NUM_TIME_SAMPLES+1:end);
                            else
                                currLFPsPreSecondDim=[];
                            end
                            
                            % PROBABLY TO FIX FOR THIRD DIM
                            if ~isempty(subjectData.LfpStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PreThirdDimDataBi)
                                currNumTrials3=length(subjectData.LfpStruct(pp).Sorted.Full.(lamLayerTag{ll}).(gratCondTag{cnd}).gratCondThirdDim);
                                currCorrectThirdDimTrials=find(arrayfun(@(jj) any([2 3 4 5 20 21 22 23 9 11 15 17 27 29 33 35]==subjectData.LfpStruct(pp).Sorted.Full.(lamLayerTag{ll}).(gratCondTag{cnd}).gratCondThirdDim(jj)), 1:currNumTrials3));
                                currLFPsPreThirdDim=subjectData.LfpStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PreThirdDimDataBi(:,currCorrectThirdDimTrials,end-NUM_TIME_SAMPLES+1:end);
                            else
                                currLFPsPreThirdDim=[];
                            end
                            
                            currNumChs=size(currLFPsPreStim,1);  currNumTrials=size(currLFPsPreFirstDim,2);
                            currNumTrials2=size(currLFPsPreSecondDim,2); currNumTrials3=size(currLFPsPreThirdDim,2);
                        else
                            currLFPsPostStimZSc=subjectData.LfpStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PostStimDataBiZSc(:,:,1:NUM_TIME_SAMPLES);
                            currLFPsStationaryZSc=subjectData.LfpStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).StationaryDataBiZSc(:,:,1:NUM_TIME_SAMPLES);
                            currLFPsPostCueZSc=subjectData.LfpStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PostCueDataBiZSc(:,:,1:NUM_TIME_SAMPLES);
                            currLFPsPreFirstDimZSc=subjectData.LfpStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PreFirstDimDataBiZSc(:,:,end-NUM_TIME_SAMPLES+1:end);
                            
                            % SELECT GRATCONDS 2, 3, 4, 5, 9, 11, 15, 17, 20, 21, 22, 23, 27, 29, 33, 35
                            if ~isempty(subjectData.LfpStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PreSecondDimDataBiZSc)
                                currNumTrials2=length(subjectData.LfpStruct(pp).Sorted.Full.(lamLayerTag{ll}).(gratCondTag{cnd}).gratCondSecondDim);
                                currCorrectSecondDimTrials=find(arrayfun(@(jj) any([2 3 4 5 20 21 22 23 9 11 15 17 27 29 33 35]==subjectData.LfpStruct(pp).Sorted.Full.(lamLayerTag{ll}).(gratCondTag{cnd}).gratCondSecondDim(jj)), 1:currNumTrials2));
                                currLFPsPreSecondDimZSc=subjectData.LfpStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PreSecondDimDataBiZSc(:,currCorrectSecondDimTrials,end-NUM_TIME_SAMPLES+1:end);
                            else
                                currLFPsPreSecondDimZSc=[];
                            end
                            
                            % PROBABLY TO FIX FOR THIRD DIM
                            if ~isempty(subjectData.LfpStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PreThirdDimDataBiZSc)
                                currNumTrials3=length(subjectData.LfpStruct(pp).Sorted.Full.(lamLayerTag{ll}).(gratCondTag{cnd}).gratCondThirdDim);
                                currCorrectThirdDimTrials=find(arrayfun(@(jj) any([2 3 4 5 20 21 22 23 9 11 15 17 27 29 33 35]==subjectData.LfpStruct(pp).Sorted.Full.(lamLayerTag{ll}).(gratCondTag{cnd}).gratCondThirdDim(jj)), 1:currNumTrials3));
                                currLFPsPreThirdDimZSc=subjectData.LfpStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PreThirdDimDataBiZSc(:,currCorrectThirdDimTrials,end-NUM_TIME_SAMPLES+1:end);
                            else
                                currLFPsPreThirdDimZSc=[];
                            end
                            
                            currNumChs=size(currLFPsPreFirstDimZSc,1);  currNumTrials=size(currLFPsPreFirstDimZSc,2);
                            currNumTrials2=size(currLFPsPreSecondDimZSc,2); currNumTrials3=size(currLFPsPreThirdDimZSc,2);
                        end
                        currPmtPrSt=nan(currNumChs,currNumTrials,NUM_FREQ_SAMPLES);
                        currPmtPsSt=nan(currNumChs,currNumTrials,NUM_FREQ_SAMPLES);
                        currPmtStat=nan(currNumChs,currNumTrials,NUM_FREQ_SAMPLES);
                        currPmtPsCu=nan(currNumChs,currNumTrials,NUM_FREQ_SAMPLES);
                        currPmtPFDm=nan(currNumChs,currNumTrials,NUM_FREQ_SAMPLES);
                        
                        currPmtPostStimZSc=nan(currNumChs,currNumTrials,NUM_FREQ_SAMPLES);
                        currPmtStationaryZSc=nan(currNumChs,currNumTrials,NUM_FREQ_SAMPLES);
                        currPmtPostCueZSc=nan(currNumChs,currNumTrials,NUM_FREQ_SAMPLES);
                        currPmtPreFirstDimZSc=nan(currNumChs,currNumTrials,NUM_FREQ_SAMPLES);
                        currPmtPreSecondDimZSc=nan(currNumChs,currNumTrials2,NUM_FREQ_SAMPLES);
                        currPmtPreThirdDimZSc=nan(currNumChs,currNumTrials3,NUM_FREQ_SAMPLES);
                        
                        for ch=1:currNumChs
                            if bb==1
                                % FULL BAND CASE: Pmt of LFPs ZScored in frequency domain
                                currPmtPreStim=mtspectrumc(squeeze(currLFPsPreStim(ch,:,:))',paramsMT)';
                                currPmtPostStim=mtspectrumc(squeeze(currLFPsPostStim(ch,:,:))',paramsMT)';
                                currPmtStationary=mtspectrumc(squeeze(currLFPsStationary(ch,:,:))',paramsMT)';
                                currPmtPostCue=mtspectrumc(squeeze(currLFPsPostCue(ch,:,:))',paramsMT)';
                                currPmtPreFirstDim=mtspectrumc(squeeze(currLFPsPreFirstDim(ch,:,:))',paramsMT)';
                                currPmtPreSecondDim=mtspectrumc(squeeze(currLFPsPreSecondDim(ch,:,:))',paramsMT)';
                                currPmtPreThirdDim=mtspectrumc(squeeze(currLFPsPreThirdDim(ch,:,:))',paramsMT)';
                                currPmtPreStimRepMean=repmat(nanmean(currPmtPreStim,1),currNumTrials,1);
                                currPmtPreStimRepStd=repmat(nanstd(currPmtPreStim,[],1),currNumTrials,1);
                                
                                %TRY ALSO
                                % data=iddata(y,[],Ts)cw
                                % Py=spa(data,winsize,freq)
                                
                                if reportPmtMinusMean
                                    currPmtPostStimZSc(ch,:,:)=(currPmtPostStim-currPmtPreStimRepMean)./(currPmtPreStimRepStd+eps);
                                    currPmtStationaryZSc(ch,:,:)=(currPmtStationary-currPmtPreStimRepMean)./(currPmtPreStimRepStd+eps);
                                    currPmtPostCueZSc(ch,:,:)=(currPmtPostCue-currPmtPreStimRepMean)./(currPmtPreStimRepStd+eps);
                                    currPmtPreFirstDimZSc(ch,:,:)=(currPmtPreFirstDim-currPmtPreStimRepMean)./(currPmtPreStimRepStd+eps);
                                    currPmtPreSecondDimZSc(ch,:,:)=(currPmtPreSecondDim-currPmtPreStimRepMean(1:currNumTrials2,:))./(currPmtPreStimRepStd(1:currNumTrials2,:)+eps);
                                    currPmtPreThirdDimZSc(ch,:,:)=(currPmtPreThirdDim-currPmtPreStimRepMean(1:currNumTrials3,:))./(currPmtPreStimRepStd(1:currNumTrials3,:)+eps);
                                else
                                    currPmtPostStimZSc(ch,:,:)=(currPmtPostStim)./(currPmtPreStimRepMean+eps);
                                    currPmtStationaryZSc(ch,:,:)=(currPmtStationary)./(currPmtPreStimRepMean+eps);
                                    currPmtPostCueZSc(ch,:,:)=(currPmtPostCue)./(currPmtPreStimRepMean+eps);
                                    currPmtPreFirstDimZSc(ch,:,:)=(currPmtPreFirstDim)./(currPmtPreStimRepMean+eps);
                                    currPmtPreSecondDimZSc(ch,:,:)=(currPmtPreSecondDim(1:currNumTrials2,:))./(currPmtPreStimRepMean(1:currNumTrials2,:)+eps);
                                    currPmtPreThirdDimZSc(ch,:,:)=(currPmtPreThirdDim(1:currNumTrials3,:))./(currPmtPreStimRepMean(1:currNumTrials3,:)+eps);
                                end
                                currPmtPrSt(ch,:,:)=currPmtPreStim;
                                currPmtPsSt(ch,:,:)=currPmtPostStim;
                                currPmtStat(ch,:,:)=currPmtStationary;
                                currPmtPsCu(ch,:,:)=currPmtPostCue;
                                currPmtPFDm(ch,:,:)=currPmtPreFirstDim;
                                
                            else
                                % ALPHA, BETA, GAMMA BAND CASES: Pmt of LFPs ZScored in time domain
                                currPmtPostStimZSc(ch,:,:)=mtspectrumc(squeeze(currLFPsPostStimZSc(ch,:,:))',paramsMT)';
                                currPmtStationaryZSc(ch,:,:)=mtspectrumc(squeeze(currLFPsStationaryZSc(ch,:,:))',paramsMT)';
                                currPmtPostCueZSc(ch,:,:)=mtspectrumc(squeeze(currLFPsPostCueZSc(ch,:,:))',paramsMT)';
                                currPmtPreFirstDimZSc(ch,:,:)=mtspectrumc(squeeze(currLFPsPreFirstDimZSc(ch,:,:))',paramsMT)';
                                currPmtPreSecondDimZSc(ch,:,:)=mtspectrumc(squeeze(currLFPsPreSecondDimZSc(ch,:,:))',paramsMT)';
                                currPmtPreThirdDimZSc(ch,:,:)=mtspectrumc(squeeze(currLFPsPreThirdDimZSc(ch,:,:))',paramsMT)';
                            end
                        end
                        PmtStruct(pr).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PostStimZSc=permute(nanmean(currPmtPostStimZSc,2),[1 3 2]);
                        PmtStruct(pr).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).StationaryZSc=permute(nanmean(currPmtStationaryZSc,2),[1 3 2]);
                        PmtStruct(pr).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PostCueZSc=permute(nanmean(currPmtPostCueZSc,2),[1 3 2]);
                        PmtStruct(pr).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PreFirstDimZSc=permute(nanmean(currPmtPreFirstDimZSc,2),[1 3 2]);
                        PmtStruct(pr).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PreSecondDimZSc=permute(nanmean(currPmtPreSecondDimZSc,2),[1 3 2]);
                        PmtStruct(pr).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PreThirdDimZSc=permute(nanmean(currPmtPreThirdDimZSc,2),[1 3 2]);
                        
                        PmtStruct(pr).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PrStim=permute(nanmean(currPmtPrSt,2),[1 3 2]);
                        PmtStruct(pr).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PsStim=permute(nanmean(currPmtPsSt,2),[1 3 2]);
                        PmtStruct(pr).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).Statry=permute(nanmean(currPmtStat,2),[1 3 2]);
                        PmtStruct(pr).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PstCue=permute(nanmean(currPmtPsCu,2),[1 3 2]);
                        PmtStruct(pr).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PrFDim=permute(nanmean(currPmtPFDm,2),[1 3 2]);
                                                
                        PmtStruct(pr).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).(subjectData.penInfos(pp).grcColorOrder).PostStimZSc=permute(nanmean(currPmtPostStimZSc,2),[1 3 2]);
                        PmtStruct(pr).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).(subjectData.penInfos(pp).grcColorOrder).StationaryZSc=permute(nanmean(currPmtStationaryZSc,2),[1 3 2]);
                        PmtStruct(pr).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).(subjectData.penInfos(pp).grcColorOrder).PostCueZSc=permute(nanmean(currPmtPostCueZSc,2),[1 3 2]);
                        PmtStruct(pr).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).(subjectData.penInfos(pp).grcColorOrder).PreFirstDimZSc=permute(nanmean(currPmtPreFirstDimZSc,2),[1 3 2]);
                        PmtStruct(pr).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).(subjectData.penInfos(pp).grcColorOrder).PreSecondDimZSc=permute(nanmean(currPmtPreSecondDimZSc,2),[1 3 2]);
                        PmtStruct(pr).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).(subjectData.penInfos(pp).grcColorOrder).PreThirdDimZSc=permute(nanmean(currPmtPreThirdDimZSc,2),[1 3 2]);
                        clear -regexp ^curr % saves memory at runtime
                    end
                end
            end
            pr=pr+1;
        end
    end
    
    
    %% PLOT SINGLE PENs
    if plotPmtsSinglePens
        rgbCmap=eye(3);
        timeWinTag={'PostStim','Stationary','PostCue','PreFirstDim','PreSecondDim','PreThirdDim'};
        close all
        
        for pp=1:length(subjectData.penIDs) % LOOP OVER PENS
            for bb=1%:4 % LOOP OVER SPECTRAL WINS
                hfig=figure('units','normalized','outerposition',[0 0 1 1]);
                yaxislim=0;
                for tt=1:6
                    for ll=1:3 % LOOP OVER LAMINAR LAYERS
                        subplot(3,6,(ll-1)*6+tt)
                        for cnd=1:3 % LOOP OVER GRAT CONDITIONS
                            plotcmapdots(rgbCmap);
                            plotmsem(freqAxisPostCue,PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).([timeWinTag{tt} 'ZSc']),rgbCmap(cnd,:)); hold on;
                            title([timeWinTag{tt} ' - ' lamLayerTag{ll} ' Layer']);
                            yaxislim=max([yaxislim; max(nanmean(PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).([timeWinTag{tt} 'ZSc']),1))]);
                            if bb==1
                                xlim([0 128]);  ylim([-1 10]);
                                ylabel('Z-Scored Power');
                            elseif bb==2
                                xlim([0 20]);  %ylim([0 0.05]);
                                ylabel('Spectral Power');
                            elseif bb==3
                                xlim([10 35]); %ylim([0 0.05*2/3]);
                                ylabel('Spectral Power');
                            elseif bb==4
                                xlim([20 100]);  %ylim([0 0.05/3]);
                                ylabel('Spectral Power');
                            end
                            xlabel('Frequency [Hz]');
                        end
                    end
                end
                if bb>1
                    for yl=1:18;   subplot(3,6,yl); ylim([0 yaxislim+0.002]);  end
                else
                    for yl=1:18;   subplot(3,6,yl); ylim([-1 yaxislim+1]);  end
                end
                legend('Attend RF','Attend OUT1','Attend OUT2');
                supertitle([sSubjectName ' ' sVisualArea ' - PEN ' num2str(subjectData.penIDs(pp)) ' - LFPs Spectral Power ' spectWinTag{bb} ' Band']);
                %saveas(hfig,['./' sSubjectName ' ' sVisualArea '/PowSpects/SinglePens' sSubjectName(1) sVisualArea '-PEN' num2str(subjectData.penIDs(pp)) '-Pmt-' spectWinTag{bb} '.fig']);
                %saveas(hfig,['./' sSubjectName ' ' sVisualArea '/PowSpects/SinglePens/' sSubjectName(1) sVisualArea '-PEN' num2str(subjectData.penIDs(pp)) '-Pmt-' spectWinTag{bb} '.png']);
            end
            close all
        end
    end
    
    %% POOL CHANNELS ACROSS PENs
    
    for bb=1%:4
        for ll=1:3
            for cnd=1:3
                currPmtPostStimZSc=[];
                currPmtStationaryZSc=[];
                currPmtPostCueZSc=[];
                currPmtPreFirstDimZSc=[];
                currPmtPreSecondDimZSc=[];
                currPmtPreThirdDimZSc=[];
                
                currPmtPreStim=[];
                currPmtPostStim=[];
                currPmtStationary=[];
                currPmtPostCue=[];
                currPmtPreFirstDim=[];
                
                for pp=1:length(PmtStruct)
                    %currSz=size(PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PreFirstDimZSc);
                    %currSz2=size(PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PreSecondDimZSc);
                    %currSz3=size(PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PreThirdDimZSc);
                    currPmtPreStim=      [currPmtPreStim;         PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PrStim];
                    currPmtPostStim=     [currPmtPostStim;        PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PsStim];
                    currPmtStationary=   [currPmtStationary;      PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).Statry];
                    currPmtPostCue=      [currPmtPostCue;         PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PstCue];
                    currPmtPreFirstDim=  [currPmtPreFirstDim;     PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PrFDim];
                    
                    currPmtPostStimZSc=     [currPmtPostStimZSc;        PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PostStimZSc];
                    currPmtStationaryZSc=   [currPmtStationaryZSc;      PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).StationaryZSc];
                    currPmtPostCueZSc=      [currPmtPostCueZSc;         PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PostCueZSc];
                    currPmtPreFirstDimZSc=  [currPmtPreFirstDimZSc;     PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PreFirstDimZSc];
                    currPmtPreSecondDimZSc= [currPmtPreSecondDimZSc;    PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PreSecondDimZSc];
                    currPmtPreThirdDimZSc=  [currPmtPreThirdDimZSc;     PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PreThirdDimZSc];
                end
                PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PostStimZSc=currPmtPostStimZSc;
                PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).StationaryZSc=currPmtStationaryZSc;
                PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PostCueZSc=currPmtPostCueZSc;
                PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PreFirstDimZSc=currPmtPreFirstDimZSc;
                PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PreSecondDimZSc=currPmtPreSecondDimZSc;
                PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PreThirdDimZSc=currPmtPreThirdDimZSc;
                
                PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PrStim=currPmtPreStim;
                PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PsStim=currPmtPostStim;
                PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).Statry=currPmtStationary;
                PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PstCue=currPmtPostCue;
                PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PrFDim=currPmtPreFirstDim;
                
                clear -regexp ^curr % saves memory at runtime
            end
        end
    end
    
%    save(['./' sSubjectName ' ' sVisualArea '/PowSpects/' sSubjectName(1) sVisualArea '-PowSpectsAllTimeWins.mat' ],'PoolPmtStruct')
%     %% SAME AS ABOVE BUT RF COLOR-SPECIFIC
%     grcColorLabel={'RGB','GRB','BGR'};
%     f
%     %penColorsRFlabels={};
%     %penColorsRFlabnum=nan(1,length(PmtStruct));
%     %for pp=1:length(PmtStruct)
%     %    penColorsRF{pp}=subjectInfos.penInfos(pp).grcColorOrder(1);
%     %    penColorsRFlabnum(pp)=strcmpi(subjectInfos.penInfos(pp).grcColorOrder(1),'R')+...
%     %                          strcmpi(subjectInfos.penInfos(pp).grcColorOrder(1),'G')*2+...
%     %                          strcmpi(subjectInfos.penInfos(pp).grcColorOrder(1),'B')*3;
%     %end
%     
%     for bb=1%:4
%         for ll=1:3
%             for cnd=1:3
%                 for grc=1:3
%                 currPmtPostStimZSc=[];
%                 currPmtStationaryZSc=[];
%                 currPmtPostCueZSc=[];
%                 currPmtPreFirstDimZSc=[];
%                 currPmtPreSecondDimZSc=[];
%                 currPmtPreThirdDimZSc=[];
%                 for pp=1:length(PmtStruct)
%                     %if strcmpi(subjectData.penInfos(pp).grcColorOrder(1),grcColorLabel{grc})
%                     %currSz=size(PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PreFirstDimZSc);
%                     %currSz2=size(PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PreSecondDimZSc);
%                     %currSz3=size(PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PreThirdDimZSc);
%                     if isfield(PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}),grcColorLabel{grc})
%                     currPmtPostStimZSc=[currPmtPostStimZSc; PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).(grcColorLabel{grc}).PostStimZSc];
%                     currPmtStationaryZSc=[currPmtStationaryZSc; PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).(grcColorLabel{grc}).StationaryZSc];
%                     currPmtPostCueZSc=[currPmtPostCueZSc; PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).(grcColorLabel{grc}).PostCueZSc];
%                     currPmtPreFirstDimZSc=[currPmtPreFirstDimZSc; PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).(grcColorLabel{grc}).PreFirstDimZSc];
%                     currPmtPreSecondDimZSc=[currPmtPreSecondDimZSc; PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).(grcColorLabel{grc}).PreSecondDimZSc];
%                     currPmtPreThirdDimZSc=[currPmtPreThirdDimZSc; PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).(grcColorLabel{grc}).PreThirdDimZSc];
%                     end
%                 end
%                 PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).(grcColorLabel{grc}).PostStimZSc=currPmtPostStimZSc;
%                 PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).(grcColorLabel{grc}).StationaryZSc=currPmtStationaryZSc;
%                 PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).(grcColorLabel{grc}).PostCueZSc=currPmtPostCueZSc;
%                 PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).(grcColorLabel{grc}).PreFirstDimZSc=currPmtPreFirstDimZSc;
%                 PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).(grcColorLabel{grc}).PreSecondDimZSc=currPmtPreSecondDimZSc;
%                 PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).(grcColorLabel{grc}).PreThirdDimZSc=currPmtPreThirdDimZSc;
%                 
%                 clear -regexp ^curr % saves memory at runtime
%                 end
%             end
%         end
%     end
%     
%     clearvars curr*






    
    %% SAME AS ABOVE BUT COLOR-SPECIFIC
    
    grcColorLabel3={'RGB','RBG','GRB','GBR','BRG','BGR'};
    
    for bb=1%:4
        for ll=1:3
            for cnd=1:3
                for grc=1:6%3
                currPmtPostStimZSc=[];
                currPmtStationaryZSc=[];
                currPmtPostCueZSc=[];
                currPmtPreFirstDimZSc=[];
                currPmtPreSecondDimZSc=[];
                currPmtPreThirdDimZSc=[];
                for pp=1:length(PmtStruct)
                    %if strcmpi(subjectData.penInfos(pp).grcColorOrder,grcColorLabel3{grc})
                       %strcmpi(subjectInfos.penInfos(pp).grcColorOrder(1),grcColorLabel{grc})
                    %currSz=size(PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PreFirstDimZSc);
                    %currSz2=size(PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PreSecondDimZSc);
                    %currSz3=size(PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).PreThirdDimZSc);
                    if isfield(PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}),grcColorLabel3{grc})
                    currPmtPostStimZSc=     [currPmtPostStimZSc;        PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).(grcColorLabel3{grc}).PostStimZSc];
                    currPmtStationaryZSc=   [currPmtStationaryZSc;      PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).(grcColorLabel3{grc}).StationaryZSc];
                    currPmtPostCueZSc=      [currPmtPostCueZSc;         PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).(grcColorLabel3{grc}).PostCueZSc];
                    currPmtPreFirstDimZSc=  [currPmtPreFirstDimZSc;     PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).(grcColorLabel3{grc}).PreFirstDimZSc];
                    currPmtPreSecondDimZSc= [currPmtPreSecondDimZSc;    PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).(grcColorLabel3{grc}).PreSecondDimZSc];
                    currPmtPreThirdDimZSc=  [currPmtPreThirdDimZSc;     PmtStruct(pp).Sorted.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).(grcColorLabel3{grc}).PreThirdDimZSc];
                    end
                end
                PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).(grcColorLabel3{grc}).PostStimZSc=currPmtPostStimZSc;
                PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).(grcColorLabel3{grc}).StationaryZSc=currPmtStationaryZSc;
                PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).(grcColorLabel3{grc}).PostCueZSc=currPmtPostCueZSc;
                PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).(grcColorLabel3{grc}).PreFirstDimZSc=currPmtPreFirstDimZSc;
                PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).(grcColorLabel3{grc}).PreSecondDimZSc=currPmtPreSecondDimZSc;
                PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).(grcColorLabel3{grc}).PreThirdDimZSc=currPmtPreThirdDimZSc;
                
                clear -regexp ^curr % saves memory at runtime
                end
            end
        end
    end
    
    
    clearvars curr*
    
    %% PLOT MEAN ACROSS PENS
    if plotPmtsPenAvg
        timeWinTag={'PostStim','Stationary','PostCue','PreFirstDim','PreSecondDim','PreThirdDim'};
        rgbCmap=eye(3);
        close all
        for bb=1%:4
            hfig=figure('units','normalized','outerposition',[0 0 1 1]);
            for tt=1:6
                for ll=1:3
                    currSpectPowRF=PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).RF.([timeWinTag{tt} 'ZSc']);
                    currSpectPowOUT=nanmean(cat(4,PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).OUT1.([timeWinTag{tt} 'ZSc']),...
                        PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).OUT2.([timeWinTag{tt} 'ZSc'])),4);
                    numPooledChs=size(currSpectPowRF,1);
                    subplot(3,6,(ll-1)*6+tt);
                    plotcmapdots([eye(3); .5 .5 .5]);
                    if reportPmtMinusMean
                        plot(freqAxisPostCue,zeros(1,length(freqAxisPostCue)),'k');
                    end
                    for cnd=1:(1+2*double(tt<6)) % Pre3rdDim is just RF
                        plotmsem(freqAxisPostCue,PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).([timeWinTag{tt} 'ZSc']),rgbCmap(cnd,:)); hold on;
                    end
                    title([timeWinTag{tt} ' - ' lamLayerTag{ll} ' Layer']);
                    if bb==1
                        xlim([0 100]);
                        if strcmp(sVisualArea,'V1')
                            if reportPmtMinusMean
                                ylim([-1 8]);
                            else
                                ylim([0 8]);
                            end
                        else
                            if reportPmtMinusMean
                                ylim([-1 2.5]);
                            else
                                ylim([0 2.5]);
                            end
                        end
                    elseif bb==2
                        xlim([0 20]);  ylim([0 0.025]);
                        ylabel('Spectral Power');
                    elseif bb==3
                        xlim([10 35]); ylim([0 0.015]);
                        ylabel('Spectral Power');
                    elseif bb==4
                        xlim([20 100]);  ylim([0 0.006]);
                        ylabel('Spectral Power');
                    end
                    if ll==3; xlabel('Frequency [Hz]'); end;
                    if tt==1; ylabel('Baseline-Normalized Power');
                        if strcmp(sVisualArea,'V1')
                            text(3,7.5,['n=' num2str(numPooledChs) ';']);
                        else
                            text(3,2.3,['n=' num2str(numPooledChs) ';']);
                        end
                    end
                    if tt<6
                        freqIndices2Plot=find(freqAxisPostCue<=100);
                        pValuesVec=nan(1,length(freqIndices2Plot));
                        for ff=1:length(freqIndices2Plot)
                            pValuesVec(ff)=signrank(currSpectPowRF(:,ff),currSpectPowOUT(:,ff));
                        end
                        [~,~,~,pValuesVecFDR]=fdr_bh(pValuesVec);
                        pValuesBar05FDR=nan(1,length(freqIndices2Plot));
                        if reportPmtMinusMean
                            pValuesBar05FDR(pValuesVecFDR<=.05)=-1;
                        else
                            pValuesBar05FDR(pValuesVecFDR<=.05)=0;
                        end
                        plot(freqAxisPostCue(freqIndices2Plot),pValuesBar05FDR,'color',[.5 .5 .5],'linewidth',2.5);
                    end
                    set(gca,'XTick',0:20:100);%set(gca,'XTickLabel',{'0','10','20','30','40','50','60','70','80','90','100'});
                end
            end
            legend('RF','OUT1','OUT2','p\leq 0.05');
            supertitle([sSubjectName ' ' sVisualArea ' - Avg across PENs - LFPs Spectral Power ' spectWinTag{bb} ' Band']);
            
            %if reportPmtMinusMean
            %    saveas(hfig,['./Summary/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-MinusMean-' spectWinTag{bb} '.png']);
            %    saveas(hfig,['./Summary/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-MinusMean-' spectWinTag{bb} '.fig']);
            %    saveas(hfig,['./Summary/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-MinusMean-' spectWinTag{bb} '.svg']);
            %    saveas(hfig,['./' sSubjectName ' ' sVisualArea '/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-MinusMean-' spectWinTag{bb} '.fig']);
            %    saveas(hfig,['./' sSubjectName ' ' sVisualArea '/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-MinusMean-' spectWinTag{bb} '.png']);
            %else
            %    saveas(hfig,['./Summary/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-' spectWinTag{bb} '.png']);
            %    saveas(hfig,['./Summary/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-' spectWinTag{bb} '.fig']);
            %    saveas(hfig,['./Summary/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-' spectWinTag{bb} '.svg']);
            %    saveas(hfig,['./' sSubjectName ' ' sVisualArea '/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-' spectWinTag{bb} '.fig']);
            %    saveas(hfig,['./' sSubjectName ' ' sVisualArea '/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-' spectWinTag{bb} '.png']);
            %end
       end
    end
    
    %%
    %% PLOT MEAN ACROSS PENS OUT=OUT1/2+OUT2/2
    if plotPmtsPenAvgOutMean
        timeWinTag={'PostStim','Stationary','PostCue','PreFirstDim','PreSecondDim','PreThirdDim'};
        rgbCmap=eye(3);
        close all
        for bb=1%:4
            hfig=figure('units','normalized','outerposition',[0 0 1 1]);
            for tt=1:6
                for ll=1:3
                    currSpectPowRF=PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).RF.([timeWinTag{tt} 'ZSc']);
                    currSpectPowOUT=nanmean(cat(4,PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).OUT1.([timeWinTag{tt} 'ZSc']),...
                        PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).OUT2.([timeWinTag{tt} 'ZSc'])),4);
                    numPooledChs=size(currSpectPowRF,1);
                    subplot(3,6,(ll-1)*6+tt);
                    plotcmapdots([1 0 0; 0 0 1; .5 .5 .5]);
                    if reportPmtMinusMean
                        plot(freqAxisPostCue,zeros(1,length(freqAxisPostCue)),'k');
                    end
                    plotmsem(freqAxisPostCue,currSpectPowRF,'r');
                    if tt<6
                        plotmsem(freqAxisPostCue,currSpectPowOUT,'b');
                    end
                    %for cnd=1:(1+2*double(tt<6)) % Pre3rdDim is just RF
                    %    plotmsem(freqAxisPostCue,PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).([timeWinTag{tt} 'ZSc']),rgbCmap(cnd,:)); hold on;
                    %end
                    title([timeWinTag{tt} ' - ' lamLayerTag{ll} ' Layer']);
                    if bb==1
                        xlim([0 100]);
                        if strcmp(sVisualArea,'V1')
                            if reportPmtMinusMean
                                ylim([-1 8]);
                            else
                                ylim([0 8]);
                            end
                        else
                            if reportPmtMinusMean
                                ylim([-1 2.5]);
                            else
                                ylim([0 2.5]);
                            end
                        end
                    elseif bb==2
                        xlim([0 20]);  ylim([0 0.025]);
                        ylabel('Spectral Power');
                    elseif bb==3
                        xlim([10 35]); ylim([0 0.015]);
                        ylabel('Spectral Power');
                    elseif bb==4
                        xlim([20 100]);  ylim([0 0.006]);
                        ylabel('Spectral Power');
                    end
                    if ll==3; xlabel('Frequency [Hz]'); end;
                    if tt==1; ylabel('Baseline-Normalized Power');
                        if strcmp(sVisualArea,'V1')
                            text(3,7.5,['n=' num2str(numPooledChs) ';']);
                        else
                            text(3,2.3,['n=' num2str(numPooledChs) ';']);
                        end
                    end
                    if tt<6
                        freqIndices2Plot=find(freqAxisPostCue<=100);
                        pValuesVec=nan(1,length(freqIndices2Plot));
                        for ff=1:length(freqIndices2Plot)
                            pValuesVec(ff)=signrank(currSpectPowRF(:,ff),currSpectPowOUT(:,ff));
                        end
                        [~,~,~,pValuesVecFDR]=fdr_bh(pValuesVec);
                        pValuesBar05FDR=nan(1,length(freqIndices2Plot));
                        if reportPmtMinusMean
                            pValuesBar05FDR(pValuesVecFDR<=.05)=-1;
                        else
                            pValuesBar05FDR(pValuesVecFDR<=.05)=0;
                        end
                        plot(freqAxisPostCue(freqIndices2Plot),pValuesBar05FDR,'color',[0 0 0],'linewidth',2.5);
                    end
                    set(gca,'XTick',0:20:100);%set(gca,'XTickLabel',{'0','10','20','30','40','50','60','70','80','90','100'});
                end
            end
            legend('RF','OUT','p\leq 0.05');
            supertitle([sSubjectName ' ' sVisualArea ' - Avg across PENs - LFPs Spectral Power ' spectWinTag{bb} ' Band']);
            
            if reportPmtMinusMean
            %    saveas(hfig,['./Summary/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg_OUT-Pmt-MinusMean-' spectWinTag{bb} '.png']);
            %    saveas(hfig,['./Summary/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg_OUT-Pmt-MinusMean-' spectWinTag{bb} '.fig']);
            %    saveas(hfig,['./Summary/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg_OUT-Pmt-MinusMean-' spectWinTag{bb} '.svg']);
            %    saveas(hfig,['./' sSubjectName ' ' sVisualArea '/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-MinusMean-' spectWinTag{bb} '.fig']);
            %    saveas(hfig,['./' sSubjectName ' ' sVisualArea '/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-MinusMean-' spectWinTag{bb} '.png']);
            else
            %    saveas(hfig,['./Summary/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-' spectWinTag{bb} '.png']);
            %    saveas(hfig,['./Summary/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-' spectWinTag{bb} '.fig']);
            %    saveas(hfig,['./Summary/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-' spectWinTag{bb} '.svg']);
            %    saveas(hfig,['./' sSubjectName ' ' sVisualArea '/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-' spectWinTag{bb} '.fig']);
            %    saveas(hfig,['./' sSubjectName ' ' sVisualArea '/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-' spectWinTag{bb} '.png']);
            end
       end
    end
    
    
        %% PLOT MEAN ACROSS PENS OVERLAY FIRST DIM OUT=OUT1/2+OUT2/2, OVERLAYED PREFIRSTDIM
    if plotPmtsPenAvgOutMean
        timeWinTag={'PostStim','Stationary','PostCue','PreFirstDim','PreSecondDim','PreThirdDim'};
        rgbCmap=eye(3);
        for bb=1%:4
            hfig=figure('units','normalized','outerposition',[0 0 1 1]);
            for tt=1:3%6
                for ll=1:3
                    currSpectPowRF=PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).RF.([timeWinTag{tt} 'ZSc']);
                    currSpectPowOUT=nanmean(cat(4,PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).OUT1.([timeWinTag{tt} 'ZSc']),...
                        PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).OUT2.([timeWinTag{tt} 'ZSc'])),4);
                    numPooledChs=size(currSpectPowRF,1);
                    subplot(3,3,(ll-1)*3+tt);
                    plotcmapdots([1 0 0; 0 0 1; .5 .5 .5]);
                    if reportPmtMinusMean
                        plot(freqAxisPostCue,zeros(1,length(freqAxisPostCue)),'k');
                    end
                    plotmsem(freqAxisPostCue,currSpectPowRF,'r');
                    if tt<6
                        plotmsem(freqAxisPostCue,currSpectPowOUT,'b');
                    end
                    currSpectPowPFDRF=PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).RF.PreFirstDimZSc;
                    currSpectPowPFDOUT=nanmean(cat(4,PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).OUT1.PreFirstDimZSc,...
                        PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).OUT2.PreFirstDimZSc),4);
                    
                    plotmsem(freqAxisPostCue,currSpectPowPFDRF,'r:');
                    plotmsem(freqAxisPostCue,currSpectPowPFDOUT,'b:');
                    
                    %for cnd=1:(1+2*double(tt<6)) % Pre3rdDim is just RF
                    %    plotmsem(freqAxisPostCue,PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).([timeWinTag{tt} 'ZSc']),rgbCmap(cnd,:)); hold on;
                    %end
                    title([timeWinTag{tt} ' - ' lamLayerTag{ll} ' Layer']);
                    if bb==1
                        xlim([0 100]);
                        if strcmp(sVisualArea,'V1')
                            if reportPmtMinusMean
                                if strcmpi(sSubjectName,'Wyman')
                                ylim([-1 8]);
                                else
                                    ylim([-.5 4]);
                                end
                            else
                                ylim([0 8]);
                            end
                        else
                            if reportPmtMinusMean
                                ylim([-1 2.5]);
                            else
                                ylim([0 2.5]);
                            end
                        end
                    elseif bb==2
                        xlim([0 20]);  ylim([0 0.025]);
                        ylabel('Spectral Power');
                    elseif bb==3
                        xlim([10 35]); ylim([0 0.015]);
                        ylabel('Spectral Power');
                    elseif bb==4
                        xlim([20 100]);  ylim([0 0.006]);
                        ylabel('Spectral Power');
                    end
                    if ll==3; xlabel('Frequency [Hz]'); end;
                    if tt==1; ylabel('Baseline-Normalized Power');
                        if strcmp(sVisualArea,'V1')
                            if strcmpi(sSubjectName,'Wyman')
                            text(3,7.5,['n=' num2str(numPooledChs) ';']);
                            else
                            text(3,3.5,['n=' num2str(numPooledChs) ';']);
                            end
                        else
                            text(3,2.3,['n=' num2str(numPooledChs) ';']);
                        end
                    end
                    if tt<6
                        freqIndices2Plot=find(freqAxisPostCue<=100);
                        pValuesVec=nan(1,length(freqIndices2Plot));
                        for ff=1:length(freqIndices2Plot)
                            pValuesVec(ff)=signrank(currSpectPowRF(:,ff),currSpectPowOUT(:,ff));
                        end
                        [~,~,~,pValuesVecFDR]=fdr_bh(pValuesVec);
                        pValuesBar05FDR=nan(1,length(freqIndices2Plot));
                        if reportPmtMinusMean
                            pValuesBar05FDR(pValuesVecFDR<=.05)=-1;
                        else
                            pValuesBar05FDR(pValuesVecFDR<=.05)=0;
                        end
                        plot(freqAxisPostCue(freqIndices2Plot),pValuesBar05FDR,'color',[0 0 0],'linewidth',2.5);
                    end
                    set(gca,'XTick',0:20:100);%set(gca,'XTickLabel',{'0','10','20','30','40','50','60','70','80','90','100'});
                end
            end
            legend('RF','OUT','p\leq 0.05');
            supertitle([sSubjectName ' ' sVisualArea ' - Avg across PENs - LFPs Spectral Power ' spectWinTag{bb} ' Band']);
            
            %if reportPmtMinusMean
            %    saveas(hfig,['./Summary/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg_OUT-Pmt-MinusMean-' spectWinTag{bb} '.png']);
            %    saveas(hfig,['./Summary/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg_OUT-Pmt-MinusMean-' spectWinTag{bb} '.fig']);
            %    saveas(hfig,['./Summary/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg_OUT-Pmt-MinusMean-' spectWinTag{bb} '.svg']);
            %    saveas(hfig,['./' sSubjectName ' ' sVisualArea '/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-MinusMean-' spectWinTag{bb} '.fig']);
            %    saveas(hfig,['./' sSubjectName ' ' sVisualArea '/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-MinusMean-' spectWinTag{bb} '.png']);
            %else
            %    saveas(hfig,['./Summary/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-' spectWinTag{bb} '.png']);
            %    saveas(hfig,['./Summary/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-' spectWinTag{bb} '.fig']);
            %    saveas(hfig,['./Summary/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-' spectWinTag{bb} '.svg']);
            %    saveas(hfig,['./' sSubjectName ' ' sVisualArea '/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-' spectWinTag{bb} '.fig']);
            %    saveas(hfig,['./' sSubjectName ' ' sVisualArea '/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-' spectWinTag{bb} '.png']);
            %end
       end
    end
    
    
    %% PLOT MEAN ACROSS PENS FIRST DIM OUT=OUT1/2+OUT2/2, LAYERS POOLED, OVERLAY PREFIRSTDIM/POSTSTIM
    if plotPmtsPenAvgOutMean
        timeWinTag={'PostStim','Stationary','PostCue','PreFirstDim','PreSecondDim','PreThirdDim'};
        rgbCmap=eye(3);
        for bb=1%:4
            hfig=figure('units','normalized','outerposition',[0 0 1 1]);
            for tt=2%1:3%6
                %for ll=1:3
                    currSpectPowRF=[PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{1}).RF.([timeWinTag{tt} 'ZSc']);...
                                    PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{2}).RF.([timeWinTag{tt} 'ZSc']);...
                                    PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{3}).RF.([timeWinTag{tt} 'ZSc'])];
                    currSpectPowOUT=[nanmean(cat(4,PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{1}).OUT1.([timeWinTag{tt} 'ZSc']),...
                                                   PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{1}).OUT2.([timeWinTag{tt} 'ZSc'])),4);
                                     nanmean(cat(4,PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{2}).OUT1.([timeWinTag{tt} 'ZSc']),...
                                                   PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{2}).OUT2.([timeWinTag{tt} 'ZSc'])),4);
                                     nanmean(cat(4,PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{3}).OUT1.([timeWinTag{tt} 'ZSc']),...
                                                   PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{3}).OUT2.([timeWinTag{tt} 'ZSc'])),4)];
                    numPooledChs=size(currSpectPowRF,1);
                    %subplot(3,6,(ll-1)*6+tt);
                    %plotcmapdots([1 0 0; 0 0 1; .5 .5 .5]);
                    if reportPmtMinusMean
                        plot(freqAxisPostCue,zeros(1,length(freqAxisPostCue)),'k');
                    end
                    plotmsem(freqAxisPostCue,currSpectPowRF,'r:');
                    %if tt<6
                    plotmsem(freqAxisPostCue,currSpectPowOUT,'b:');
                    %end
                    currSpectPowPFDRF=[PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{1}).RF.PreFirstDimZSc;...
                                       PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{2}).RF.PreFirstDimZSc;...
                                       PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{3}).RF.PreFirstDimZSc];
                    currSpectPowPFDOUT=[nanmean(cat(4,PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{1}).OUT1.PreFirstDimZSc,...
                                                      PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{1}).OUT2.PreFirstDimZSc),4);...
                                        nanmean(cat(4,PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{2}).OUT1.PreFirstDimZSc,...
                                                      PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{2}).OUT2.PreFirstDimZSc),4);...
                                        nanmean(cat(4,PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{3}).OUT1.PreFirstDimZSc,...
                                                      PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{3}).OUT2.PreFirstDimZSc),4)];
                    
                    plotmsem(freqAxisPostCue,currSpectPowPFDRF,'r');
                    plotmsem(freqAxisPostCue,currSpectPowPFDOUT,'b');
                    
                    
                    %for cnd=1:(1+2*double(tt<6)) % Pre3rdDim is just RF
                    %    plotmsem(freqAxisPostCue,PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).([timeWinTag{tt} 'ZSc']),rgbCmap(cnd,:)); hold on;
                    %end
                    title([timeWinTag{tt} ' - ' lamLayerTag{ll} ' Layer']);
                    if bb==1
                        xlim([0 100]);
                        if strcmp(sVisualArea,'V1')
                            if reportPmtMinusMean
                                if strcmpi(sSubjectName,'Wyman')
                                ylim([-1 8]);
                                else
                                    ylim([-1 6]);
                                end
                            else
                                ylim([0 8]);
                            end
                        else
                            if reportPmtMinusMean
                                ylim([-1 2.5]);
                            else
                                ylim([0 2.5]);
                            end
                        end
                    elseif bb==2
                        xlim([0 20]);  ylim([0 0.025]);
                        ylabel('Spectral Power');
                    elseif bb==3
                        xlim([10 35]); ylim([0 0.015]);
                        ylabel('Spectral Power');
                    elseif bb==4
                        xlim([20 100]);  ylim([0 0.006]);
                        ylabel('Spectral Power');
                    end
                    if ll==3; xlabel('Frequency [Hz]'); end;
                    if tt==1; ylabel('Baseline-Normalized Power');
                        if strcmp(sVisualArea,'V1')
                            if strcmpi(sSubjectName,'Wyman')
                            text(3,7.5,['n=' num2str(numPooledChs) ';']);
                            else
                            text(3,3.5,['n=' num2str(numPooledChs) ';']);
                            end
                        else
                            text(3,2.3,['n=' num2str(numPooledChs) ';']);
                        end
                    end
                    if tt<6
                        freqIndices2Plot=find(freqAxisPostCue<=100);
                        pValuesVec=nan(1,length(freqIndices2Plot));
                        for ff=1:length(freqIndices2Plot)
                            pValuesVec(ff)=signrank(currSpectPowPFDRF(:,ff),currSpectPowPFDOUT(:,ff));
                        end
                        [~,~,~,pValuesVecFDR]=fdr_bh(pValuesVec);
                        pValuesBar05FDR=nan(1,length(freqIndices2Plot));
                        if reportPmtMinusMean
                            pValuesBar05FDR(pValuesVecFDR<=.05)=-.5;
                        else
                            pValuesBar05FDR(pValuesVecFDR<=.05)=0;
                        end
                        plot(freqAxisPostCue(freqIndices2Plot),pValuesBar05FDR,'color',[0 0 0],'linewidth',2.5);
                    end
                    set(gca,'XTick',0:20:100);%set(gca,'XTickLabel',{'0','10','20','30','40','50','60','70','80','90','100'});
                %end
            end
            %legend('RF','OUT','p\leq 0.05');
            supertitle([sSubjectName ' ' sVisualArea ' - Avg across PENs - LFPs Spectral Power ' spectWinTag{bb} ' Band']);
            
            %if reportPmtMinusMean
            %    saveas(hfig,['./Summary/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg_OUT-Pmt-MinusMean-' spectWinTag{bb} '.png']);
            %    saveas(hfig,['./Summary/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg_OUT-Pmt-MinusMean-' spectWinTag{bb} '.fig']);
            %    saveas(hfig,['./Summary/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg_OUT-Pmt-MinusMean-' spectWinTag{bb} '.svg']);
            %    saveas(hfig,['./' sSubjectName ' ' sVisualArea '/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-MinusMean-' spectWinTag{bb} '.fig']);
            %    saveas(hfig,['./' sSubjectName ' ' sVisualArea '/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-MinusMean-' spectWinTag{bb} '.png']);
            %else
            %    saveas(hfig,['./Summary/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-' spectWinTag{bb} '.png']);
            %    saveas(hfig,['./Summary/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-' spectWinTag{bb} '.fig']);
            %    saveas(hfig,['./Summary/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-' spectWinTag{bb} '.svg']);
            %    saveas(hfig,['./' sSubjectName ' ' sVisualArea '/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-' spectWinTag{bb} '.fig']);
            %    saveas(hfig,['./' sSubjectName ' ' sVisualArea '/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-' spectWinTag{bb} '.png']);
            %end
       end
    end
    
    
    %% PLOT MEAN ACROSS PENS COLOR CODED
    if plotPmtsPenAvgColorCoded
        timeWinTag={'PostStim','Stationary','PostCue','PreFirstDim','PreSecondDim','PreThirdDim'};
        rgbCmap=eye(3);
        close all
        for bb=1%:4
            hfig=figure('units','normalized','outerposition',[0 0 1 1]);
            for tt=[1 3 4 5]%1:4%6
                for ll=1:3
                    currSpectPowRF=PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).RF.([timeWinTag{tt} 'ZSc']);
                    currSpectPowOUT=nanmean(cat(4,PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).OUT1.([timeWinTag{tt} 'ZSc']),...
                        PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).OUT2.([timeWinTag{tt} 'ZSc'])),4);
                    numPooledChs=size(currSpectPowRF,1);
                    subplot(3,4,(ll-1)*4+tt-(tt>2));%subplot(3,6,(ll-1)*6+tt);
                    plotcmapdots([eye(3)]); plot(-inf,-inf,'-','linewidth',2,'color','k'); plot(-inf,-inf,':','linewidth',2,'color','k');
                    if reportPmtMinusMean
                        plot(freqAxisPostCue,zeros(1,length(freqAxisPostCue)),'k');
                    end
                    %for cnd=1:(1+2*double(tt<6)) % Pre3rdDim is just RF
                    %    plotmsem(freqAxisPostCue,PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).([timeWinTag{tt} 'ZSc']),rgbCmap(cnd,:)); hold on;
                    %end
                    %for cnd=1:(1+2*double(tt<6)) % RF is RED
                    %    PmtRFgrcRed=PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).(gratCondTag{cnd}).(grcColorLabel3{1}).([timeWinTag{tt} 'ZSc']);
                    %plotmsem(freqAxisPostCue,,rgbCmap(cnd,:)); hold on;
                    
                    % RF case when RF is R/G/B and OUT case when RF is R/G/B
                    PmtRFgrcR=[PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).RF.RGB.([timeWinTag{tt} 'ZSc']);...;
                               PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).RF.RBG.([timeWinTag{tt} 'ZSc'])];
                           
                    PmtRFgrcG=[PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).RF.GRB.([timeWinTag{tt} 'ZSc']);...
                               PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).RF.GBR.([timeWinTag{tt} 'ZSc'])];
                           
                    PmtRFgrcB=[PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).RF.BRG.([timeWinTag{tt} 'ZSc']);...
                               PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).RF.BGR.([timeWinTag{tt} 'ZSc'])];
                           
                    PmtOUTgrcR=[ 1/2*PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).OUT1.RGB.([timeWinTag{tt} 'ZSc'])+...
                                 1/2*PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).OUT2.RGB.([timeWinTag{tt} 'ZSc']);...
                                 1/2*PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).OUT1.RBG.([timeWinTag{tt} 'ZSc'])+...
                                 1/2*PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).OUT2.RBG.([timeWinTag{tt} 'ZSc'])];
                             
                    PmtOUTgrcG=[ 1/2*PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).OUT1.GRB.([timeWinTag{tt} 'ZSc'])+...
                                 1/2*PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).OUT2.GRB.([timeWinTag{tt} 'ZSc']);...
                                 1/2*PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).OUT1.GBR.([timeWinTag{tt} 'ZSc'])+...
                                 1/2*PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).OUT2.GBR.([timeWinTag{tt} 'ZSc'])];
                             
                    PmtOUTgrcB=[ 1/2*PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).OUT1.BRG.([timeWinTag{tt} 'ZSc'])+...
                                 1/2*PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).OUT2.BRG.([timeWinTag{tt} 'ZSc']);...
                                 1/2*PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).OUT1.BGR.([timeWinTag{tt} 'ZSc'])+...
                                 1/2*PoolPmtStruct.(spectWinTag{bb}).(lamLayerTag{ll}).OUT2.BGR.([timeWinTag{tt} 'ZSc'])];
                    
                             if 1
                                 plotmsem(freqAxisPostCue,PmtRFgrcR,[1 0 0],.001); hold on;
                                 plotmsem(freqAxisPostCue,PmtRFgrcG,[0 1 0],.001); hold on;
                                 plotmsem(freqAxisPostCue,PmtRFgrcB,[0 0 1],.001); hold on;
                                 if tt<6
                                     plotmsem(freqAxisPostCue,PmtOUTgrcR,[1 0 0],.001,':'); hold on;
                                     plotmsem(freqAxisPostCue,PmtOUTgrcG,[0 1 0],.001,':'); hold on;
                                     plotmsem(freqAxisPostCue,PmtOUTgrcB,[0 0 1],.001,':'); hold on;
                                 end
                             else % Peak Normalization
                                 plotmsem(freqAxisPostCue,PmtRFgrcR/max(nanmean(PmtRFgrcR,1)),[1 0 0],.001); hold on;
                                 plotmsem(freqAxisPostCue,PmtRFgrcG/max(nanmean(PmtRFgrcG,1)),[0 1 0],.001); hold on;
                                 plotmsem(freqAxisPostCue,PmtRFgrcB/max(nanmean(PmtRFgrcB,1)),[0 0 1],.001); hold on;
                                 if tt<6
                                     plotmsem(freqAxisPostCue,PmtOUTgrcR/max(nanmean(PmtOUTgrcR,1)),[1 0 0],.001,':'); hold on;
                                     plotmsem(freqAxisPostCue,PmtOUTgrcG/max(nanmean(PmtOUTgrcG,1)),[0 1 0],.001,':'); hold on;
                                     plotmsem(freqAxisPostCue,PmtOUTgrcB/max(nanmean(PmtOUTgrcB,1)),[0 0 1],.001,':'); hold on;
                                 end
                             end
                    %plotmsem(freqAxisPostCue,[PmtRFgrcR; PmtRFgrcR; PmtRFgrcR],[1 0 0],.001); hold on;
                    %plotmsem(freqAxisPostCue,[PmtOUTgrcR; PmtOUTgrcR; PmtOUTgrcR],[0 0 1],.001); hold on;
                                         
                    title([timeWinTag{tt} ' - ' lamLayerTag{ll} ' Layer']);
                    
                    if bb==1
                        xlim([0 100]);
                        if strcmp(sVisualArea,'V1') && strcmpi(sSubjectName,'Wyman')
                            ylimtop=22;
                        elseif strcmp(sVisualArea,'V1') && strcmpi(sSubjectName,'Taylor')
                            ylimtop=6;
                        elseif strcmp(sVisualArea,'V4') && strcmpi(sSubjectName,'Wyman')
                            ylimtop=2;
                        elseif strcmp(sVisualArea,'V4') && strcmpi(sSubjectName,'Taylor')
                            ylimtop=2;
                        end
                        %ylimtop=2;
                        if reportPmtMinusMean
                            ylim([-1 ylimtop]);
                        else
                            ylim([0 ylimtop]);
                        end
                    elseif bb==2
                        xlim([0 20]);  ylim([0 0.025]);
                        ylabel('Spectral Power');
                    elseif bb==3
                        xlim([10 35]); ylim([0 0.015]);
                        ylabel('Spectral Power');
                    elseif bb==4
                        xlim([20 100]);  ylim([0 0.006]);
                        ylabel('Spectral Power');
                    end
                    if ll==3; xlabel('Frequency [Hz]'); end;
                    if tt==1; ylabel('Baseline-Normalized Power');
                        text(3,ylimtop*.95,['n=' num2str(size(PmtRFgrcR,1)) ';'],'color','r');
                        text(3,ylimtop*.85,['n=' num2str(size(PmtRFgrcG,1)) ';'],'color','g');
                        text(3,ylimtop*.75,['n=' num2str(size(PmtRFgrcB,1)) ';'],'color','b');
                    end
                    %if tt<6
                    %    freqIndices2Plot=find(freqAxisPostCue<=100);
                    %    pValuesVec=nan(1,length(freqIndices2Plot));
                    %    for ff=1:length(freqIndices2Plot)
                    %        pValuesVec(ff)=signrank(currSpectPowRF(:,ff),currSpectPowOUT(:,ff));
                    %    end
                    %    [~,~,~,pValuesVecFDR]=fdr_bh(pValuesVec);
                    %    pValuesBar05FDR=nan(1,length(freqIndices2Plot));
                    %    if reportPmtMinusMean
                    %        pValuesBar05FDR(pValuesVecFDR<=.05)=-1;
                    %    else
                    %        pValuesBar05FDR(pValuesVecFDR<=.05)=0;
                    %    end
                    %    plot(freqAxisPostCue(freqIndices2Plot),pValuesBar05FDR,'color',[.5 .5 .5],'linewidth',2.5);
                    %end
                    set(gca,'XTick',0:20:100);%set(gca,'XTickLabel',{'0','10','20','30','40','50','60','70','80','90','100'});
                end
            end
            legend('RF is red','RF is green','RF is blue','Attend RF','Attend OUT');
            supertitle([sSubjectName ' ' sVisualArea ' - Avg across PENs - LFPs Spectral Power ' spectWinTag{bb} ' Band']);
            
            if reportPmtMinusMean
            %    saveas(hfig,['./Summary/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-MinusMean-' spectWinTag{bb} '_COLORCODED.png']);
            %    saveas(hfig,['./Summary/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-MinusMean-' spectWinTag{bb} '.fig']);
            %    saveas(hfig,['./Summary/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-MinusMean-' spectWinTag{bb} '.svg']);
            %    saveas(hfig,['./' sSubjectName ' ' sVisualArea '/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-MinusMean-' spectWinTag{bb} '.fig']);
            %    saveas(hfig,['./' sSubjectName ' ' sVisualArea '/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-MinusMean-' spectWinTag{bb} '.png']);
            else
            %    saveas(hfig,['./Summary/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-' spectWinTag{bb} '_COLORCODED.png']);
            %    saveas(hfig,['./Summary/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-' spectWinTag{bb} '.fig']);
            %    saveas(hfig,['./Summary/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-' spectWinTag{bb} '.svg']);
            %    saveas(hfig,['./' sSubjectName ' ' sVisualArea '/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-' spectWinTag{bb} '.fig']);
            %    saveas(hfig,['./' sSubjectName ' ' sVisualArea '/PowSpects/' sSubjectName(1) sVisualArea '-PEN_Avg-Pmt-' spectWinTag{bb} '.png']);
            end
       end
    end
end