vb_channel_statics

 Estimate mean amplitude by using histgram (for outlier exclusion)
 and calculate amplitude max to find bad channel & trials
 amplitude max for temporal difference is also calculated
   [ymax1,ymax2,ystd1,ystd2] = vb_channel_statics(X)
 X(ch,t,n)  : signal for ch-th channel at time t & n-th trial
 ymax1(ch,n) : amplitude max for each channel & trial
             = max_t (abs(X(ch,t,n))) 
 ymax2(ch,n) :amplitude max of temporal difference
 ystd1   : mean amplitude  (: sqrt(<X(t)^2>))
 ystd2   : mean difference (: sqrt(<dX(t)^2>))
 
 2008-6-1 Masa-aki Sato

 Copyright (C) 2011, ATR All Rights Reserved.
 License : New BSD License(see VBMEG_LICENSE.txt)

0001 function    [ymax1,ymax2,ystd1,ystd2] = vb_channel_statics(data,sw)
0002 % Estimate mean amplitude by using histgram (for outlier exclusion)
0003 % and calculate amplitude max to find bad channel & trials
0004 % amplitude max for temporal difference is also calculated
0005 %   [ymax1,ymax2,ystd1,ystd2] = vb_channel_statics(X)
0006 % X(ch,t,n)  : signal for ch-th channel at time t & n-th trial
0007 % ymax1(ch,n) : amplitude max for each channel & trial
0008 %             = max_t (abs(X(ch,t,n)))
0009 % ymax2(ch,n) :amplitude max of temporal difference
0010 % ystd1   : mean amplitude  (: sqrt(<X(t)^2>))
0011 % ystd2   : mean difference (: sqrt(<dX(t)^2>))
0012 %
0013 % 2008-6-1 Masa-aki Sato
0014 %
0015 % Copyright (C) 2011, ATR All Rights Reserved.
0016 % License : New BSD License(see VBMEG_LICENSE.txt)
0017 
0018 [Nch, T, Ntry] = size(data);
0019 
0020 % first order difference on DIM=2 (time)
0021 ydif = abs(diff(data,1,2));
0022 
0023 ymax1 = zeros(Nch, Ntry);
0024 ymax2 = zeros(Nch, Ntry);
0025 
0026 ystd1 = zeros(Nch, Ntry);
0027 ystd2 = zeros(Nch, Ntry);
0028 
0029 if nargin==1, sw = 1; end;
0030 
0031 switch    sw
0032 case    1
0033     for n=1:Nch
0034         [ymax1(n,:), ystd] = ...
0035             vb_get_max_std_ratio(squeeze(data(n,:,:)));
0036         ystd1(n,:) = ystd;
0037         
0038         [ymax2(n,:), ystd] = ...
0039             vb_get_max_std_ratio(squeeze(ydif(n,:,:)));
0040         ystd2(n,:) = ystd;
0041     end
0042 case    2
0043     for n=1:Nch
0044         [ymax1(n,:), ystd1(n,:)] = ...
0045             vb_get_max_med_ratio(squeeze(data(n,:,:)));
0046     
0047         [ymax2(n,:), ystd2(n,:)] = ...
0048             vb_get_max_med_ratio(squeeze(ydif(n,:,:)));
0049     end
0050 end
0051 %ratio = max(ratio1,ratio2);
0052 
0053 return
0054 %%% END %%%
0055 
0056 
0057 xlist = 2:0.1:15;
0058 Nhist = length(xlist);
0059 
0060 
0061 rmax1 = zeros(Nch, Ntry);
0062 rmax2 = zeros(Nch, Ntry);
0063 
0064 rhist1 = zeros(1,Nhist);
0065 rhist2 = zeros(1,Nhist);
0066 
0067 yrate1 = zeros(Nch,1);
0068 yrate2 = zeros(Nch,1);
0069 
0070 subplot(NY,NX,nfig); nfig = nfig + 1;
0071 plot(ystd1)
0072 %plot(yrate1)
0073 hold on
0074 xlim([0 Nch])
0075 plot(ystd2,'-r')
0076 %plot(yrate2,'-r')
0077 
0078 title('Mean amplitude for each channel')
0079 
0080 
0081 return
0082 
0083 
0084 NX=2;NY=2;
0085 nfig = 1;
0086 for n=1:NY
0087     subplot(NY,NX,nfig); nfig=nfig+1;
0088     plot(xh(:,n),yh(:,n))
0089     hold on
0090     plot([ymax(n) ymax(n)],[0 yh(1,n)],'-r')
0091     plot([ymed(n) ymed(n)],[0 yh(1,n)],'-b')
0092     
0093     subplot(NY,NX,nfig); nfig=nfig+1;
0094     plot(xd(:,n),yd(:,n))
0095     hold on
0096     plot([dmax(n) dmax(n)],[0 yd(1,n)],'-r')
0097     plot([dmed(n) dmed(n)],[0 yd(1,n)],'-b')
0098 end