如何计算平均径向从2D噪声功率谱1D功率谱

问题描述 投票:0回答:1

大家好我试图从图像的二维FFT计算1D功率谱。我平均水平,但通过观察它没有让我感觉到的图表做了。能否请你建议如何在2D数据将达到噪声功率谱的一维表示做径向平均。谢谢你,我会感谢您的帮助。

这里是我的代码$

 fid = fopen('C:\Users\3772khobrap\Desktop\project  related\NPS_cal_data_UB\100000006.raw','r');
img = fread(fid,[512 512],'uint16');
roi = zeros(64);
avg = zeros(64);
Ux= 0.0075;% Pixel size
Uy = 0.0075;% Pixel size
%% This block of code is subdividing imaage into smaller ROI and averaging purpose
for r = 1:8 
    r_shift = (r-1)*64;
    for c = 1:8 
        c_shift = (c-1)*64;
        for i = 1:64 
            for j = 1:64
                p = img(i+r_shift,j+c_shift);
                roi(i,j) = p;
            end
        end
        avg = avg+roi;
    end
end
avg = avg./64;
%%Actual process of NPS calculation
scale = (Ux*Uy)./(64*64);%Scaling fator as per NPS calculation formula
f_x = 1/(2*Ux);%Nyquiest frequecy along x direction
f_y = 1/(2*Ux);%Nyquiest frequecy along y direction
FFT_2d = (fftshift(fft2(avg))).^2;% Power spectrum calculation
NPS = abs(FFT_2d).*scale; %% 2D NPS 
f = linspace(-f_x,f_y,64);% x-axis  vector for 1D NPS 
X_img = linspace(-f_x,f_x,512);% X axis of NPS image
Y_img = linspace(-f_x,f_x,512);% Y axis of NPS image

figure(1)
subplot(2,2,1)
imagesc(X_img,Y_img,img)
colormap gray
xlabel('X [cm]'); ylabel('Y [cm]')
title('noise image')
subplot(2,2,2)
imagesc(f,f,log(NPS))
colormap gray
xlabel('frequency [cm^{-1}]'); ylabel('frequency [cm^{-1}]');
title('2D NPS')
subplot(2,2,3)
plot(f_p,NPS(:,32))
xlabel('frequency [cm^{-2}]'); ylabel('NPS [cm^{-2}]')
title('1D NPS from central slice')
subplot(2,2,4)
plot(f_p,mean(NPS,2))
xlabel('frequency [cm^{-2}]'); ylabel('NPS [cm^2]')
title('1D NPS along X direction')
matlab image-processing fft noise
1个回答
0
投票

您可以设定这样的功能:

function profile = radialAverage(IMG, cx, cy, w)
    % computes the radial average of the image IMG around the cx,cy point
    % w is the vector of radii starting from zero
    [a,b] = size(IMG);
    [X, Y] = meshgrid( (1:a)-cx, (1:b)-cy);
    R = sqrt(X.^2 + Y.^2);
    profile = [];
    for i = w % radius of the circle
        mask = (i-1<R & R<i+1); % smooth 1 px around the radius
        values = (1-abs(R(mask)-i)) .* double(IMG(mask)); % smooth based on distance to ring
        % values = IMG(mask); % without smooth
        profile(end+1) = mean( values(:) );
    end
end
© www.soinside.com 2019 - 2024. All rights reserved.