📄 unimod.m
字号:
function [conc]=unimod(conc,rmod,cmod)
[ns,nc]=size(conc);
% 1) look for the maximum
for j=1:nc,
[y,imax(j)]=max(conc(:,j));
end
% 2) force unimodality shape
for j=1:nc,
rmax=conc(imax(j),j);
k=imax(j);
% disp('maximum at point');disp(k)
% 2a) discard left maxima (tolerance rmod)
while k>1,
k=k-1;
if conc(k,j)<=rmax,
rmax=conc(k,j);
else,
rmax2=rmax*rmod;
if conc(k,j)>rmax2,
% disp('no left unimodality in point: ');disp(k);
% pause
if cmod==0,conc(k,j)=1.0E-30;end
if cmod==1,conc(k,j)=conc(k+1,j);end
if cmod==2,
if rmax>0,
conc(k,j)=(conc(k,j)+conc(k+1,j))/2;
conc(k+1,j)=conc(k,j);
k=k+2;
else
conc(k,j)=0;
end
end
rmax=conc(k,j);
end
end
end
% 2b) discard right maxima (tolerance rmod)
rmax=conc(imax(j),j);
k=imax(j);
while k<ns,
k=k+1;
if conc(k,j)<=rmax,
rmax=conc(k,j);
else,
rmax2=rmax*rmod;
if conc(k,j)>rmax2,
if cmod==0,conc(k,j)=1.0E-30;end
if cmod==1,conc(k,j)=conc(k-1,j);end
if cmod==2,
% disp('no right unimodality in point: ');
% disp([k,conc(k,j),rmax])
% disp('rmax= ');disp(rmax);
% pause
if rmax>0,
conc(k,j)=(conc(k,j)+conc(k-1,j))/2;
conc(k-1,j)=conc(k,j);
k=k-2;
else
conc(k,j)=0;
end
end
rmax=conc(k,j);
end
end
end
end
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -