function [ESN,PSN,MUN,SIN]=emlht(w,ES,PS,MU,SI,zm)
%function [ESN,PSN,MUN,SIN]=emlht(w,ES,PS,MU,SI,zm)
% updates LH subband of HMT model once
% by tying within scale for each subband.
%
% Author: H. Choi
% Last update : 12/14/1998
%
% the data structure of 2-D DWT follows the format used
% in the Rice wavelelet matlab toolbox (see README)
%
% internal variables:
% M : no. of mixture densities
% P : size of image (PxP pixels)
% level : no. of levels of HMT
%
% Input :
% w : data PxP matrix (wavelet transform of an image)
% ES : state transition matrix (MxMxPxP)
% PS : state probability matrix (MxPxP)
% MU : mean matrix (MxPxP)
% SI : variance matrix (MxPxP)
% zm : type of density functions
% zm = 1 : zero mean (do not update MU)
% zm = 0 : nonzero mean (update MU)
%
% Output:
% Updated ES PS MU SI in ESN PSN MUN SIN
M=size(ES,1);
P=size(w,1);
level=log2(P);
BE=zeros(M,P,P);
BEP=zeros(M,P,P);
BER=zeros(M,P,P);
AL=zeros(M,P,P);
P1=zeros(M,P,P);
P2=zeros(M,M,P,P);
%UP step
wtmp = repmat(w,[1 1 M]);
wtmp = shiftdim(wtmp,2);
si=1; ei=2^(level-1); sj=2^(level-1)+1; ej=P;
gtmp = gauss(wtmp,MU,SI);
scale = repmat(mean(gtmp,1),[M 1 1]);
BE(:,si:ei,sj:ej) = gtmp(:,si:ei,sj:ej)./scale(:,si:ei,sj:ej);
%clear MUtmp SItmp;
for k=level:-1:2
J=2^(k-1);J2=J*J; si = 1; ei = J; sj = J+1; ej = 2*J;
EStmp = reshape(ES(:,:,si:ei,sj:ej),M,M*J2);
if M == 2
%%%%%% For M=2 the following is faster
BEtmp = zeros(M,M*J2);
BEtmp(:,1:M:(M*J2))=reshape(BE(:,si:ei,sj:ej),M,J2);
BEtmp(:,2:M:(M*J2))=BEtmp(:,1:M:(M*J2));
else
% For general M (not equal to 2) use the following
%
BEtmp = zeros(M,M*4^(k-1));
for m=1:M
BEtmp(:,m:M:(M*4^(k-1)))=reshape(BE(:,si:ei,sj:ej,:),M,4^(k-1));
end;
end;
BEtmp = reshape(EStmp.*BEtmp,[M M J J]);
BEP(:,si:ei,sj:ej) = squeeze(sum(BEtmp,1));
sni = 1; eni = J/2; snj = J/2+1; enj = J;
%construct betachild matrix here
BCtmp = BEP(:,si:2:ei,sj:2:ej);
BCtmp = BCtmp.*BEP(:,si+1:2:ei,sj:2:ej);
BCtmp = BCtmp.*BEP(:,si:2:ei,sj+1:2:ej);
BCtmp = BCtmp.*BEP(:,si+1:2:ei,sj+1:2:ej);
scaletmp = repmat(mean(BCtmp,1),[M 1 1]);
scale(:,sni:eni,snj:enj) = scale(:,sni:eni,snj:enj).*scaletmp;
BE(:,sni:eni,snj:enj)=gtmp(:,sni:eni,snj:enj)./scale(:,sni:eni,snj:enj).*BCtmp;
%construct BE(:,pai(i),paj(j),dindex) matrix
Btmp=zeros(M,J,J);
Btmp(:,1:2:J,1:2:J)=BE(:,sni:eni,snj:enj);
Btmp(:,2:2:J,1:2:J)=BE(:,sni:eni,snj:enj);
Btmp(:,1:2:J,2:2:J)=BE(:,sni:eni,snj:enj);
Btmp(:,2:2:J,2:2:J)=BE(:,sni:eni,snj:enj);
BER(:,si:ei,sj:ej)=Btmp./BEP(:,si:ei,sj:ej);
end;
clear EStmp BEtmp BCtmp Btmp;
%DOWN step
%initialize
AL(:,1,2) = PS(:,1,2);
for k=2:level
J = 2^(k-1); J2=J*J;
si=1; ei=J; sj=J+1; ej=2*J;
sni = 1; eni = J/2; snj = J/2+1; enj = J;
Atmp=zeros(M,J,J);
Atmp(:,1:2:J,1:2:J)=AL(:,sni:eni,snj:enj);
Atmp(:,2:2:J,1:2:J)=AL(:,sni:eni,snj:enj);
Atmp(:,1:2:J,2:2:J)=AL(:,sni:eni,snj:enj);
Atmp(:,2:2:J,2:2:J)=AL(:,sni:eni,snj:enj);
Atmp = repmat(reshape(Atmp.*BER(:,si:ei,sj:ej),1,M*J2),[M 1]);
EStmp = reshape(ES(:,:,si:ei,sj:ej),M,M*J2);
ALtmp = reshape(EStmp.*Atmp,[M M J J]);
AL(:,si:ei,sj:ej) = squeeze(sum(ALtmp,2));
end;
clear Atmp EStmp ALtmp;
%compute probabilities
for k=2:level
J=2^(k-1); J2=J*J;
si=1; ei=J; sj=J+1; ej=2*J;
sni = 1; eni = J/2; snj = J/2+1; enj = J;
temp = repmat(sum(AL(:,si:ei,sj:ej).*BE(:,si:ei,sj:ej), 1),[M 1]);
P1(:,si:ei,sj:ej) = AL(:,si:ei,sj:ej).*BE(:,si:ei,sj:ej)./temp;
%compute P2
if M == 2
% For M=2 the following may be faster
BEtmp = zeros(M,M*J2);
BEtmp(:,1:M:(M*J2))=reshape(BE(:,si:ei,sj:ej),M,J2);
BEtmp(:,2:M:(M*J2))=BEtmp(:,1:M:(M*J2));
else
% For general M (not equal to 2) use the following
BEtmp = zeros(M,M*J2);
for m=1:M
BEtmp(:,m:M:(M*J2))=reshape(BE(:,si:ei,sj:ej,:),M,J2);
end;
end;
BEtmp = reshape(BEtmp,[M M J J]);
EStmp = ES(:,:,si:ei,sj:ej);
Atmp=zeros(M,J,J);
Atmp(:,1:2:J,1:2:J)=AL(:,sni:eni,snj:enj);
Atmp(:,2:2:J,1:2:J)=AL(:,sni:eni,snj:enj);
Atmp(:,1:2:J,2:2:J)=AL(:,sni:eni,snj:enj);
Atmp(:,2:2:J,2:2:J)=AL(:,sni:eni,snj:enj);
Atmp = repmat(reshape(Atmp,1,M*J2),[M 1]);
Atmp = reshape(Atmp,[M M J J]);
BERtmp = repmat(reshape(BER(:,si:ei,sj:ej),1,M*J2),[M 1]);
BERtmp = reshape(BERtmp,[M M J J]);
temp = repmat(reshape(temp,1,M*J2),[M 1]);
temp = reshape(temp, [M M J J]);
P2(:,:,si:ei,sj:ej)=BEtmp.*EStmp.*Atmp.*BERtmp./temp;
end;
P1(:,1,2)=AL(:,1,2).*BE(:,1,2)./repmat(sum(AL(:,1,2).*BE(:,1,2),1),[M 1 1]);
clear temp BEtmp EStmp Atmp BERtmp;
%PSP=PS; ESP=ES; MUP=MU; SIP=SI;
%M step
PS(:,1,2)=P1(:,1,2);
for k=2:level
J=2^(k-1); J2=J*J;
si=1; ei=J; sj=J+1; ej=2*J;
sni = 1; eni = J/2; snj = J/2+1; enj = J;
pstmp = sum(sum(P1(:,si:ei,sj:ej),3),2)/J2;
pstmp = pstmp.*(pstmp>1e-4)+1e-4*(pstmp<=1e-4);
PS(:,si:ei,sj:ej) = repmat(pstmp,[1 J J]);
if zm == 0
% do not update MU if zero mean densities
mutmp = sum(sum(wtmp(:,si:ei,sj:ej).*P1(:,si:ei,sj:ej),3),2)/J2;
MU(:,si:ei,sj:ej) = repmat(mutmp,[1 J J])./PS(:,si:ei,sj:ej);
end;
sitmp = sum(sum((wtmp(:,si:ei,sj:ej)-MU(:,si:ei,sj:ej)).^2.*P1(:,si:ei,sj:ej),3),2)/J2;
SI(:,si:ei,sj:ej) = repmat(sitmp,[1 J J])./PS(:,si:ei,sj:ej);
estmp =sum(sum(P2(:,:,si:ei,sj:ej),4),3)/J2;
ptmp = [PS(:,sni,snj)'; PS(:,sni,snj)'];
ES(:,:,si:ei,sj:ej)= repmat(estmp,[1 1 J J])./repmat(ptmp,[1 1 J J]);
end; %k
ESN=ES; PSN=PS; MUN=MU; SIN=SI;
没有合适的资源?快使用搜索试试~ 我知道了~
数学建模常用算法matlab代码.rar
共37个文件
m:29个
png:2个
tar:1个
1.该资源内容由用户上传,如若侵权请联系客服进行举报
2.虚拟产品一经售出概不退款(资源遇到问题,请及时私信上传者)
2.虚拟产品一经售出概不退款(资源遇到问题,请及时私信上传者)
版权申诉
5星 · 超过95%的资源 5 下载量 158 浏览量
2021-08-23
00:52:19
上传
评论 1
收藏 18.42MB RAR 举报
温馨提示
数学建模常用算法matlab代码.rar
资源推荐
资源详情
资源评论
收起资源包目录
数学建模常用算法matlab代码.rar (37个子文件)
隐马尔可夫预测代码(含有大量案例)
whmt1
ex512.m 511B
gauss.m 81B
hmtdeno.m 1KB
makenoise.m 1KB
up.m 976B
posthh.m 3KB
Idwt.m 802B
ex256.m 511B
barbara.png 174KB
testcon.m 187B
Iidwt.m 936B
daubcqf.m 2KB
postlh.m 3KB
Idwt2.m 3KB
hmttrain.m 2KB
lena512 256KB
psnr.m 88B
example3.m~ 954B
hmtmodel.m 2KB
hdenoise.m 4KB
lenamodel.mat 16.5MB
lena256 64KB
emhlt.m 5KB
makeh1.m 84B
vec2mat.m 1KB
Iidwt2.m 2KB
peppers.png 157KB
example1.m 508B
posthl.m 3KB
emhht.m 5KB
emlht.m 5KB
README 3KB
dfilters.m 3KB
hmtdeno0.m 1KB
example3.m 1KB
example2.m 541B
whmt1.tar 17.23MB
共 37 条
- 1
资源评论
- 南城转暖2022-02-24用户下载后在一定时间内未进行评价,系统默认好评。
- Wymundluo2023-05-11这个资源对我启发很大,受益匪浅,学到了很多,谢谢分享~
- toremimum2022-09-20终于找到了超赞的宝藏资源,果断冲冲冲,支持!
- m0_712482302022-09-24资源中能够借鉴的内容很多,值得学习的地方也很多,大家一起进步!
- lzy200506302024-04-08感谢大佬分享的资源给了我灵感,果断支持!感谢分享~
花花老爷子
- 粉丝: 77
- 资源: 1
上传资源 快速赚钱
- 我的内容管理 展开
- 我的资源 快来上传第一个资源
- 我的收益 登录查看自己的收益
- 我的积分 登录查看自己的积分
- 我的C币 登录后查看C币余额
- 我的收藏
- 我的下载
- 下载帮助
安全验证
文档复制为VIP权益,开通VIP直接复制
信息提交成功