[Gate-users] Post processing of projection

[Khalid Hussain]

Dear Gate users;

I want to calculate the fwhm of the SPECT system projection of point
source. I tried the below Matlab code but it is not giving the
satisfactory result of full-width-at half maximum.

The projection *image.jpg *and *image.mat* file is attached.

May I request you to please share your experience with the calculation of
FWHM?

Thank you so much.


V=max(max(img));
[M,N]=find(img==V);
[fwhmx, fwhmy, fwhm] =FWHM_calculation(img,N, M,pixelSize);
function [fwhm_x,fwhm_y,fwhmT]=FWHM_calculation(Image,cx,cy,pixelSize)
% this function is used to calculate the full width half maximum of image
% image: the input image, including a point source or square source in the
% center
% fwhm : the output value for fwhm
CC=ceil(cx);
CR=ceil(cy);
N=25;
Column1=CC-N:CC+N ;
Row1=CR-N:CR+N;
ROIimage=Image(Row1,Column1);
profile_x=sum(ROIimage,1);
X=1:1:2*N+1;
V=max(profile_x);
plot(X,profile_x,'r-x','LineWidth',0.7,'MarkerSize',5)
hold on
[fitresult1, gof1] =createFitSPECTfwhm2(X,profile_x,V,N);
rmse=gof1.rmse;
fwhm_x=2*sqrt(log(2))*pixelSize*fitresult1.c1
profile_y=sum(ROIimage,2);
plot(X,profile_y,'b-o','LineWidth',0.5,'MarkerSize',3)
V1=max(profile_y);
[fitresult2, gof2] =createFitSPECTfwhm2(X,profile_y,V1,N);
fwhm_y=2*sqrt(log(2))*pixelSize*fitresult2.c1
profile_yt=profile_y';
profile=(profile_x+ profile_yt)/2;
plot(X,profile,'k-*','LineWidth',0.7,'MarkerSize',5)
V=max(profile);
[fitresult2, gof2] =createFitSPECTfwhm2(X,profile,V,N);
fwhmT=2*sqrt(log(2))*pixelSize*fitresult2.c1
xlabel('Pixel Index')
ylabel('Counts')
title('Profile along Yaxis','fontSize',9)
legend('Profile along X-axis','Profile along Y-axis','Average Profile')
function [fitresult, gof] = createFitSPECTfwhm2(X, profile,maxValue,Xcenter)
%CREATEFIT(X,PROFILE_X)
% Create a fit.
%
% Data for 'untitled fit 1' fit:
% X Input : X
% Y Output: profile_x
% Output:
% fitresult : a fit object representing the fit.
% gof : structure with goodness-of fit info.
%
% See also FIT, CFIT, SFIT.
% Auto-generated by MATLAB on 26-Mar-2014 11:10:32
%% Fit: 'untitled fit 1'.
[xData, yData] = prepareCurveData( X, profile );
% Set up fittype and options.
ft = fittype( 'gauss1' );
% ft = fittype( 'gauss2' );
opts = fitoptions( ft );
opts.Display = 'Off';
% opts.Lower = [-Inf -Inf 0];
opts.Lower = [-Inf -Inf 0];
opts.StartPoint = [maxValue Xcenter 1];
opts.Upper = [Inf Inf Inf];
% opts.Upper = [Inf Inf Inf];
% Fit model to data.
[fitresult, gof] = fit( xData, yData, ft, opts );
% Plot fit with data.
figure( 'Name', 'Gaussian fit 1' );
h = plot( fitresult, xData, yData );
legend( h, 'profile_x vs. X', 'Gaussian fit 1', 'Location', 'NorthEast' );
% Label axes
xlabel( 'X' );
ylabel( 'profile_x' );


[image: pointSource.jpg]


[image: fitFunc.jpg]



I also asked this question
<https://www.mathworks.com/matlabcentral/answers/855595-how-to-calculate-image-resolution-full-width-at-half-maximum-of-a-point-source?s_tid=mlc_ans_email_view#answer_724450>
on mathwork.com but yet no acceptable solution.

Thank you for your coopeartion.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.opengatecollaboration.org/pipermail/gate-users/attachments/20210617/a5cba671/attachment-0001.html>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: pointSource.jpg
Type: image/jpeg
Size: 10247 bytes
Desc: not available
URL: <http://lists.opengatecollaboration.org/pipermail/gate-users/attachments/20210617/a5cba671/attachment-0002.jpg>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: fitFunc.jpg
Type: image/jpeg
Size: 34033 bytes
Desc: not available
URL: <http://lists.opengatecollaboration.org/pipermail/gate-users/attachments/20210617/a5cba671/attachment-0003.jpg>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: img.mat
Type: application/octet-stream
Size: 694002 bytes
Desc: not available
URL: <http://lists.opengatecollaboration.org/pipermail/gate-users/attachments/20210617/a5cba671/attachment-0001.obj>