plota.m

matlab数字信号处理工具箱

        end;
        if nargin<4,
                alpha = 1; 
        elseif isnumeric(arg4),
                alpha = arg4;
        elseif ischar(arg4) & isempty(str2num(arg4))
                alpha = flag_implicit_significance;
        else
                alpha = arg4;
        end;
        if nargin<5,
                Y = [];
        else
                Y = arg5; 
        end;
             
        gf = arg2;
        if ~isfield(X.M,gf)
                error('PLOTA TFAR_ALL: field %s is unknown\n',gf);
        end;
        
        %ClassList = {'left','right','foot','tongue'};
        
        M   = size(X.M.AR,1);
        tmp = size(X.M.AR);
        MOP = tmp(2)/tmp(1);
        
        if ~isfield(X,'Label'),
                for k1 = 1:M,
                        Label{k1}=['# ',int2str(k1)];
                end;
        end;
        nr = ceil(sqrt(M));
        nc = ceil(M/nr);
        if (nargin>2) & ~isempty(arg3),
                hf = arg3;
        else
                for k = 1:M,
                        hf(k)=subplot(nr,nc,k);
                end;
        end;

        if isempty(Y);
                x0 = real(getfield(X.M,gf));
        else
                x0 = (real(getfield(X.M,gf)) - real(getfield(Y.M,gf)))./(real(getfield(X.SE,gf))*X.N + real(getfield(Y.SE,gf))*Y.N);
        end;
        
        clim = [min(x0(:)),max(x0(:))]
        caxis(clim);
        cm = colormap;
        for k = 1:M,
                subplot(hf(k));
                %imagesc(X.T,X.F,squeeze(X.M.logS1(k,:,:)))
        
                x = x0(k,1:length(X.F),:);
                ci = getfield(X.SE,gf)*(X.N-1);
                ci = ci(k,1:length(X.F),:);
                if alpha < .5,
                        xc = 2 + round(62*(squeeze(x)-clim(1))/diff(clim));
                        sz = size(x);
                        %x = x(:);
                        bf = prod(size(x));
                        xc(abs(x) < (ci*norminv(1-alpha/(2*bf)))) = 1;
                        %x(abs(x) < .5) = NaN;
                        %x = reshape(x,sz);
                        cm(1,:) = [1,1,1];
                        colormap(cm);
                else
                        xc = 1+round(63*(squeeze(x)-clim(1))/diff(clim));
                        colormap('default');
                end;
                x1 = reshape(cm(xc,1),size(xc));
                x2 = reshape(cm(xc,2),size(xc));
                x3 = reshape(cm(xc,3),size(xc));
                
                %h = imagesc(X.T,X.F,cat(3,x1,x2,x3)*diff(clim)+clim(1),clim);
                %imagesc(X.T,X.F,squeeze(X.M.logS1(k,:,:)))
                h = imagesc(X.T,X.F,cat(3,x1,x2,x3),clim);
        end;
        if isfield(X,'TITLE');	
                TIT = X.TITLE;
                TIT(TIT=='_')=' ';
                suptitle(TIT);
        else
                TIT = '';
        end
        
        
elseif strcmp(X.datatype,'TF-MVAR')    % logS2 and S1 
        
        %GF = {'C','DC','AR','PDC','DTF','dDTF','ffDTF','COH','pCOH','pCOH2','S','h','phaseS','phaseh','coh','logh','logS'};
        
        if nargin<2,
                arg2 = 'logS1';
        end;
        if nargin<3,
                alpha = 1; 
        elseif isnumeric(arg3),
                alpha = arg3;
        elseif isempty(str2num(arg3))
                alpha = flag_implicit_significance;
        else
                alpha = arg3;
        end;
        if nargin<4,
                Y = [];
        else
                Y = arg4; 
        end;
        
        gf = arg2;
        if ~isfield(X.M,gf)
                error('PLOTA TFMVAR_ALL: field %s is unknown\n',gf);
        end;
        
        %ClassList = {'left','right','foot','tongue'};
        
        M   = size(X.M.AR,1);
        tmp = size(X.M.AR);
        MOP = tmp(2)/tmp(1);
        
        if ~isfield(X,'Label'),
                for k1 = 1:M,
                        Label{k1}=['# ',int2str(k1)];
                end;
        end;
        nr = ceil(sqrt(M));
        nc = ceil(M/nr);
        if nargin>2,
                hf = arg3;
        else
                for k1 = 1:M,
                        for k2 = 1:M,
                                hf(k1,k2)=subplot(M,M,(k1-1)*M+k2);
                        end;
                end;
        end;
        
        if isempty(Y);
                x0 = real(getfield(X.M,gf));
        else
                x0 = (real(getfield(X.M,gf)) - real(getfield(Y.M,gf)))./(real(getfield(X.SE,gf))*X.N + real(getfield(Y.SE,gf))*Y.N);
        end;
        
        clim = [min(x0(:)),max(x0(:))];
        caxis(clim);
        cm = colormap;
        for k1 = 1:M,
                for k2 = 1:M,
                        subplot(hf(k1*M-M+k2));
                        x  = x0(k1,k2,1:length(X.F),:);
                        ci = getfield(X.SE,gf)*(X.N-1);
                        ci = ci(k1,k2,1:length(X.F),:);
                        if alpha < .5,
                                xc = 2 + round(62*(squeeze(x)-clim(1))/diff(clim));
                                sz = size(x);
                                %x = x(:);
                                bf = prod(size(x));
                                xc(abs(x) < (ci*norminv(1-alpha/(2*bf)))) = 1;
                                %x(abs(x) < .5) = NaN;
                                %x = reshape(x,sz);
                                cm(1,:) = [1,1,1];
                                colormap(cm);
                        else
                                xc = 1+round(63*(squeeze(x)-clim(1))/diff(clim));
                                colormap('default');
                        end;
                        x1 = reshape(cm(xc,1),size(xc));
                        x2 = reshape(cm(xc,2),size(xc));
                        x3 = reshape(cm(xc,3),size(xc));
                        
                        %h = imagesc(X.T,X.F,cat(3,x1,x2,x3)*diff(clim)+clim(1),clim);
                        h = imagesc(X.T,X.F,cat(3,x1,x2,x3),clim);
                        %h  = imagesc(X.T,X.F,squeeze(x),clim);
                        if k1==1, title(Label{k2}); end;
                        if k2==1, ylabel(Label{k1});end;
                end;
        end;
        %caxis = clim;
        h   = colorbar;
        %tmp = get(h,'ytick')'/64*diff(clim)+clim(1);
        %set(h,'yticklabel',num2str(tmp));
        
        if isfield(X,'TITLE');	
                TIT = X.TITLE;
                TIT(TIT=='_')=' ';
                suptitle(TIT);
        else
                TIT = '';
        end
        
        
elseif strcmp(X.datatype,'EDF'),
        data = arg2;
        [nr,nc]=size(data);
        for k = 1:nc,
                subplot(nc,1,k);
                %t = (X.Block.number(3)+1:X.Block.number(4))/X.AS.MAXSPR*X.Dur;
                
                t = X.FILE.POS*X.Dur+(1-nr:0)/X.AS.MAXSPR*X.Dur/size(X.SIE.T,2)*size(X.SIE.T,1);
                plot(t, data(:,k));
                ylabel(deblank(X.Label(X.SIE.ChanSelect(k),:)));    
        end;
        xlabel('t [s]')
        
        
elseif strcmpi(X.datatype,'pfurtscheller_spectral_difference'),
        nc = ceil(sqrt(X.NS));
        nr = ceil(X.NS/nc);
        nch = size(X.AR,1)/X.NS;
        f = (0:.1:X.SampleRate/2)';
        H = zeros(length(f),X.NC+1);
        for k1=1:nc,
                for k2=1:nr,
                        c = k1+(k2-1)*nc;
                        if nargin>1,
                                H = X.S(:,c+X.NS*(0:X.NC));        
                                F = 0:size(X.S,1)-1;
                        else
                                for k3 = 1:X.NC+1;
                                        ix = c + X.NS*(k3-1);
                                        [H(:,k3), F] = freqz(sqrt(X.PE(ix,end)/X.SampleRate),ar2poly(X.AR(ix,:)),f,X.SampleRate);
                                end
                        end;
                        subplot(nc,nr,c);
                        semilogy(F,abs(H),'-');
                        legend({'ref','1','2'});
                        ylabel(sprintf('%s/[%s]^{1/2}',X.PhysDim,X.samplerate_units));
                        v=axis;v(2:4)=[max(F),1e-2,10];axis(v);
                        %hold on;
                        grid on;
                        if isfield(X,'Label');
                                if iscell(X.Label)
                                        title(X.Label{c});
                                else
                                        title(X.Label(c,:));
                                end;
                        else
                                title(['channel # ',int2str(c)]);
                        end;
                end
        end;
        
elseif strcmpi(X.datatype,'spectrum'),
        if nargin>1,
                Mode=arg2;
        else
                Mode='log';
        end;
        
        if ~isfield(X,'samplerate_units')
                X.samplerate_units = 'Hz';    
        end;
        if ~isfield(X,'PhysDim')
                X.PhysDim = '[1]';    
        end;
        if ~isfield(X,'QUANT')
                X.QUANT = 0;    
        end;
        
        [n,p] = size(X.AR);
        H=[]; F=[];
        for k=1:size(X.AR,1);
                [h,f] = freqz(sqrt(X.PE(k,size(X.AR,2)+1)/(X.SampleRate*2*pi)),ar2poly(X.AR(k,:)),(0:64*p)/(128*p)*X.SampleRate',X.SampleRate);
                H(:,k)=h(:);F(:,k)=f(:);
        end;
        if strcmp(lower(Mode),'log')
                semilogy(F,abs(H),'-',[0,X.SampleRate/2]',[1;1]*X.QUANT/sqrt(12*X.SampleRate),'k:');
                ylabel(sprintf('%s/[%s]^{1/2}',X.PhysDim,X.samplerate_units));
                
        elseif strcmp(lower(Mode),'log2')
                semilogy(F,real(H).^2+imag(H).^2,'-',[0,X.SampleRate/2]',[1;1]*X.QUANT.^2/(12*X.SampleRate),'k:');
                ylabel(sprintf('[%s]^2/%s',X.PhysDim,X.samplerate_units));
                
        elseif strcmp(lower(Mode),'lin')
                plot(F,abs(H),'-',[0,X.SampleRate/2]',[1;1]*X.QUANT/sqrt(12*X.SampleRate),'k:');
                ylabel(sprintf('%s/[%s]^{1/2}',X.PhysDim,X.samplerate_units));
                
        elseif strcmp(lower(Mode),'lin2')
                plot(F,real(H).^2+imag(H).^2,'-',[0,X.SampleRate/2]',[1;1]*X.QUANT.^2/(12*X.SampleRate),'k:');
                ylabel(sprintf('[%s]^2/%s',X.PhysDim,X.samplerate_units));
        end;
        xlabel(sprintf('f [%s]',X.samplerate_units));
        
        
elseif strcmp(X.datatype,'confusion'),
        if nargin>1,
                [kap,sd,H,z,OA,SA]=kappa(X.data);
                fprintf(1,'%s\n',repmat('-',1,8*(size(X.data,1)+1)));
                fprintf(1,'Kappa = %5.3f