function vme = meanMaxScW(v, wl, wm, stepwl, percmi,percma)
% ---------------------------------------------------------------------------------------------
% meanMaxScW.m: Compute the average value of the maxima on data windows of the input vector.
% ---------------------------------------------------------------------------------------------
% Compute the mean values of the maxima of signal values in subwindows of a sliding window
% The distribution tails can be excluded.
% vme = meaMaxScW(v,wl,wm,stepwl,perci,percf)
% v = signal vector
% wl = width of the sliding window
% wm = width of the subwindows
% stepwl = sliding window step
% perci = % of min values to be excluded
% percf = % of max values to be excluded
%
% ---------------------------------------------------------------------------------------------
% Maurizio Varanini, Clinical Physiology Institute, CNR, Pisa, Italy
% For any comment or bug report, please send e-mail to: maurizio.varanini@ifc.cnr.it
% ---------------------------------------------------------------------------------------------
if(nargin<4), stepwl=1; end
if(nargin<5), percmi=5; end
if(nargin<6), percma=percmi; end
if nargin==0, meanMaxScW_test; return; end
if(percmi<0 || percma<0 || (percmi+percma)>=100), fprintf('\n"meanMaxScW": error in input parameters\n'); pause; end
if(stepwl<1), stepwl=1; end
if(isempty(v)), vme=[]; return; end
[nt,dummy]=size(v);
if(nt<dummy), v=v'; nt=dummy; end
wlm=fix(wl/wm);
ii=1+floor(wlm*percmi/100);
fi=wlm-floor(wlm*percma/100);
wl=wlm*wm;
nwl=fix((nt-wl)/stepwl)+1;
vme=zeros(nwl,1);
iwl=1;
for is=1:nwl
vmawm= max(reshape(v(iwl:iwl+wl-1),wm,wlm))';
omawl= sort(vmawm);
vme(is)=mean(omawl(ii:fi));
iwl=iwl+stepwl;
end
return
end % = function ===============================================================
% -------------------------------------------------------------------------------------------------
% Test program (simulating ecg channel selection for QRS detection)
function meanMaxScW_test
close all, clear variables
Fs=1000; vlen=60*Fs;
RRpmean=0.4;
ww=fix(RRpmean*1.5*Fs); % wide windows (each window contains at least one QRS)
wn=fix(RRpmean/5*Fs); % short windows (many windows do not contain QRSs)
wl=20*ww; %10*ww;
stepwl=0.2*wl; %1; %0.5*wl; %wl; %fix(0.01*Fs);
fprintf('ww=%d, wn=%d, wl=%d, stepwl=%d\n', ww,wn,wl,stepwl);
% building of test signal
sigType=3; % <== Test signal type
ns=2;
t=(1:vlen)'/Fs;
switch(sigType)
case 1
xb1 = abs(0*sin(t) + 3*randn(size(t)));
xb2 = abs(0*sin(t) + 3*randn(size(t)));
pp = floor(RRpmean*Fs);
tp = 1:pp:length(xb1);
xp=zeros(size(t)); xp(tp) = 3; %+1*abs(randn(size(tp)));
xs(:,1)=xb1;
xs(:,2)=0.3*xb2+xp;
case 2
triaNoise=ones(size(t)); nptr=1001; npCyc=1750; np1=npCyc-nptr;
nptrt=floor(size(t,1)/(npCyc))*(npCyc);
triaNoise(1:nptrt)=repmat([2+triang(nptr);ones(np1,1)],floor(size(t,1)/npCyc),1);
xb1 = abs(.5*abs(sin(t))+3*triaNoise .*randn(size(t)));
aecgd=ones(size(t)); npq=fix(0.1*Fs); npRR=floor(RRpmean*Fs); npic=npRR-npq;
npecg=floor(size(t,1)/(npRR))*(npRR);
aecgd(1:npecg)=repmat([2+triang(npq);ones(npic,1)],floor(size(t,1)/npRR),1);
xb2 = abs(.5*abs(sin(t))+2*aecgd .*randn(size(t)));
xs(:,1)=xb1;
xs(:,2)=xb2;
case 3
triaNoise=ones(size(t)); nptr=350; npCyc=1350; np1=npCyc-nptr;
nptrt=floor(size(t,1)/(npCyc))*(npCyc);
triaNoise(1:nptrt)=repmat([2+triang(nptr);ones(np1,1)],floor(size(t,1)/npCyc),1);
xb1 = abs(.5*abs(sin(t))+3*triaNoise .*randn(size(t)));
aecgd=ones(size(t)); npq=fix(0.1*Fs); npRR=floor(RRpmean*Fs); npic=npRR-npq;
npecg=floor(size(t,1)/(npRR))*(npRR);
aecgd(1:npecg)=repmat([2+triang(npq);ones(npic,1)],floor(size(t,1)/npRR),1);
xb2 = abs(.5*abs(sin(t))+2*aecgd .*randn(size(t)));
xs(:,1)=xb1;
xs(:,2)=xb2;
end
percmi=0;
percma=5;
for is=1:ns,
vmew(:,is)=meanMaxScW(xs(:,is),wl,ww,stepwl,percmi,percma);
vmen(:,is)=meanMaxScW(xs(:,is),wl,wn,stepwl,percmi,percma);
end
figure, set(gcf,'Color','white');
wgmima= mimaxscG(xs(Fs:end-Fs,1),0,0,0.1);
for is=1:ns,
subplot(ns,1,is), plot(t,xs(:,is));
ylim(wgmima);
if(is~=ns), set(gca,'XTickLabel',''); end
if(is==1), title('Signals'); end
end
vmtime=1/2*wl/Fs+(1:size(vmew,1))*stepwl/Fs; vmtime=vmtime';
figure, set(gcf,'Color','white');
wgmima= mimaxscG(vmew,0,0,0.1);
for is=1:ns,
subplot(ns,1,is), plot(vmtime,vmew(:,is));
ylim(wgmima);
if(is~=ns), set(gca,'XTickLabel',''); end
if(is==1), title('vmew'); end
end
figure, set(gcf,'Color','white');
wgmima= mimaxscG(vmen,0,0,0.1);
for is=1:ns,
subplot(ns,1,is), plot(vmtime,vmen(:,is));
ylim(wgmima);
if(is~=ns), set(gca,'XTickLabel',''); end
if(is==1), title(' vmen'); end
end
vqf=vmew./vmen;
%vqf=vmew-1.0*vmen;
figure, set(gcf,'Color','white');
wgmima= mimaxscG(vqf,0,0,0.1);
for is=1:ns,
subplot(ns,1,is), plot(vmtime,vqf(:,is));
ylim(wgmima);
if(is~=ns), set(gca,'XTickLabel',''); end
if(is==1), title('vqf'); end
end
figure, set(gcf,'Color','white');
wgmima= mimaxscG(xs,0,0,0.1);
for is=1:ns,
subplot(ns+1,1,is), plot(t,xs(:,is));
ylim(wgmima);
set(gca,'XTickLabel','');
end
subplot(ns+1,1,ns+1), hold on;
plot(vmtime,vqf(:,2)-vqf(:,1));
plot(vmtime,zeros(size(vmtime)),'r');
end %== function ================================================================
%