I am trying to reconstruct an algorithm belong to this paper:
Decomposition of biospeckle images in temporary spectral bands
Here is an explanation of the algorithm:
We recorded a sequence of N successive speckle images with a sampling
frequency fs. In this way it was possible to observe how a pixel
evolves through the N images. That evolution can be treated as a time
series and can be processed in the following way: Each signal
corresponding to the evolution of every pixel was used as input to a
bank of filters. The intensity values were previously divided by their
temporal mean value to minimize local differences in reflectivity or
illumination of the object. The maximum frequency that can be
adequately analyzed is determined by the sampling theorem and s half
of sampling frequency fs. The latter is set by the CCD camera, the
size of the image, and the frame grabber. The bank of filters is
outlined in Fig. 1. In our case, ten 5° order Butterworth11 filters
were used, but this number can be varied according to the required
discrimination. The bank was implemented in a computer using MATLAB
software. We chose the Butter-worth filter because, in addition to its
simplicity, it is maximally flat. Other filters, an infinite impulse
response, or a finite impulse response could be used. By means of this
bank of filters, ten corresponding signals of each filter of each
temporary pixel evolution were obtained as output. Average energy Eb
in each signal was then calculated:
where pb(n) is the intensity of the filtered pixel in the nth image
for filter b divided by its mean value and N is the total number of
images. In this way, en values of energy for each pixel were obtained,
each of hem belonging to one of the frequency bands in Fig. 1. With
these values it is possible to build ten images of the active object,
each one of which shows how much energy of time-varying speckle there
is in a certain frequency band. False color assignment to the gray
levels in the results would help in discrimination.
and here is my MATLAB code base on that :
clear all
for i=0:39
str = num2str(i);
str1 = strcat(str,'.mat');
load(str1);
D{i+1}=A;
end
new_max = max(max(A));
new_min = min(min(A));
for i=20:180
for j=20:140
ts = [];
for k=1:40
ts = [ts D{k}(i,j)]; %%% kth image pixel i,j --- ts is time series
end
ts = double(ts);
temp = mean(ts);
ts = ts-temp;
ts = ts/temp;
N = 5; % filter order
W = [0.00001 0.05;0.05 0.1;0.1 0.15;0.15 0.20;0.20 0.25;0.25 0.30;0.30 0.35;0.35 0.40;0.40 0.45;0.45 0.50];
N1 = 5;
for ind = 1:10
Wn = W(ind,:);
[B,A] = butter(N1,Wn);
ts_f(ind,:) = filter(B,A,ts);
end
for ind=1:10
imag_test1{ind}(i,j) =sum((ts_f(ind,:)./mean(ts_f(ind,:))).^2);
end
end
end
for i=1:10
temp_imag = imag_test1{i}(:,:);
x=isnan(temp_imag);
temp_imag(x)=0;
temp_imag=medfilt2(temp_imag);
t_max = max(max(temp_imag));
t_min = min(min(temp_imag));
temp_imag = (temp_imag-t_min).*(double(new_max-new_min)/double(t_max-t_min))+double(new_min);
imag_test2{i}(:,:) = temp_imag;
end
for i=1:10
A=imag_test2{i}(:,:);
B=A/max(max(A));
B=histeq(B);
figure,imshow(B)
colorbar
end
but I am not getting the same result as paper. has anybody has aby idea why? or where I have gone wrong?
Refrence Link to the paper
