[问题] 多尺度tophat的实现问题

[Joy_King]

各位好，最近在用tophat作图像增强，看了下原理，还蛮简单的，但是在实现多尺度tophat上面遇到些问题，还望高手指点一下。
单尺度tophat和多尺度tophat我都做了，但是可能因为我的实现上面有些问题，导致我的多尺度tophat的增强效果几乎和单尺度的没多大差别，有时候甚至还不如单尺度的，而且运行时间上明显增加，我参考的文献是【Peihe Tang, Hao Liu, Kaiqiong Sun. Enhancement of coronary angiogram by estimation of local background】，附件是我的代码，顺便贴出来一下。
这个我已经折腾蛮久了，真心希望有高手指点一下哈，谢谢了。

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Single scale enhancement for pre-observation.
se = strel('disk', 11);
res = imtophat(Iori,se); res=imadjust(res,stretchlim(res),[0 1]);
figure;
subplot(1,3,1); imshow(Iori); title('Original Image');
subplot(1,3,2); imshow(Open_res); title('Background Estimation using Single Scale');
subplot(1,3,3); imshow(res); title('Top-hat Result using Single Scale');

%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Multiscale enhancement.
[nrow,ncol]=size(Iori);

% Construct different Structure Elements(SE).
se1=strel('disk',1); se2=strel('disk',3); se3=strel('disk',5);
se4=strel('disk',7); se5=strel('disk',9); se6=strel('disk',11);

% Opening operators.
Img=double(Iori);
Open0=Img;
Open1=imopen(Img,se1); Open2=imopen(Img,se2); Open3=imopen(Img,se3);
Open4=imopen(Img,se4); Open5=imopen(Img,se5); Open6=imopen(Img,se6);

% % Reconstruct
% Open1=imreconstruct(Open1,Img); Open2=imreconstruct(Open2,Img); Open3=imreconstruct(Open3,Img);
% Open4=imreconstruct(Open4,Img); Open5=imreconstruct(Open5,Img); Open6=imreconstruct(Open6,Img);

% % Difference between two adjacent opening results.
DI0=Open0-Open1; DI1=Open1-Open2; DI2=Open2-Open3;
DI3=Open3-Open4; DI4=Open4-Open5; DI5=Open5-Open6;
% DI0=Open0.*Open1; DI1=Open1.*Open2; DI2=Open2.*Open3;
% DI3=Open3.*Open4; DI4=Open4.*Open5; DI5=Open5.*Open6;

% DI0=Open1-Open0; DI1=Open2-Open0; DI2=Open3-Open0;
% DI3=Open4-Open0; DI4=Open5-Open0; DI5=Open6-Open0;

% search max abs (DD) and order num (ind)
% % chosen 1 : only one scale according to maximum difference
% -----------------------------------------------------------
bInd=ones(nrow,ncol);
DD=DI0;
ind0 = bInd;
KK=DI1>=DD; ind1 = bInd.*double(KK) + double(KK); bInd = bInd + 1; DD=max(DD,DI1);
KK=DI2>=DD; ind2 = bInd.*double(KK) + double(KK); bInd = bInd + 1; DD=max(DD,DI2);
KK=DI3>=DD; ind3 = bInd.*double(KK) + double(KK); bInd = bInd + 1; DD=max(DD,DI3);
KK=DI4>=DD; ind4 = bInd.*double(KK) + double(KK); bInd = bInd + 1; DD=max(DD,DI4);
KK=DI5>=DD; ind5 = bInd.*double(KK) + double(KK); bInd = bInd + 1; DD=max(DD,DI5);
ind=max(ind0,ind1);ind=max(ind,ind2); ind=max(ind,ind3); ind=max(ind,ind4); ind=max(ind,ind5);

%%
%Multi-scale opening for background estimation and enhancement
bk=zeros(nrow,ncol); % background
kk=zeros(nrow,ncol);
for i=1:nrow
for j=1:ncol
if ind(i,j)==0 bk(i,j) = Open0(i,j); end
if ind(i,j)==1 bk(i,j) = Open1(i,j); end
if ind(i,j)==2 bk(i,j) = Open2(i,j); end
if ind(i,j)==3 bk(i,j) = Open3(i,j); end
if ind(i,j)==4 bk(i,j) = Open4(i,j); end
if ind(i,j)==5 bk(i,j) = Open5(i,j); end
if ind(i,j)==6 bk(i,j) = Open6(i,j); end
end
end
figure;
subplot(1,3,1);imshow(Iori); title('Original Image');
subplot(1,3,2);imshow(bk,[]);title('Multiscale background estimation');
kk=Img-bk;
I2=uint8(kk-1);
I3=imadjust(I2,stretchlim(I2),[0 1]);
subplot(1,3,3);imshow(I3);title('Multiscale top-hat result');
% imwrite(I3,'img002enResult005.bmp');