# defining the center of an image and recognise the contrast disks

6 views (last 30 days)
sesilia maidelin on 10 Feb 2022
Commented: sesilia maidelin on 17 Feb 2022
so i have these images of a test object. i need to make a code that recognises the contrast disks ( the dark dots around the disk), where the disks are actually present around the whole disk as in the image of the test object ( not the x ray ). i used imcrop and ginput to crop a rectangle in the required area, but that is not good enough because on the parts where the contrast disk is not visible it goes awry also not at all standardised. i need matlab to be able to locate the center of the contrast disks around the bigger disk with minimum user input so the results are repeatable. please help! thank you
##### 2 CommentsShow NoneHide None
Simon Chan on 10 Feb 2022
Edited: Simon Chan on 10 Feb 2022
Why not use a template to locate the disks? Leeds phantom has a fixed geometry, you can first determine the entire phantom center and its angle of rotation. Finally put the template on it, then you know where are the disks.
sesilia maidelin on 16 Feb 2022
yes, i think that's what the code below is doing. Thank you!

Sign in to comment.

### Accepted Answer

DGM on 10 Feb 2022
Edited: DGM on 10 Feb 2022
Consider the example:
A = imread('https://www.mathworks.com/matlabcentral/answers/uploaded_files/890405/image.png');
B = bwareafilt(A<128,3); % pick the 3 objects with largest area
imshow(B)
% get properties
S = regionprops(B,'area','solidity','centroid','minferetproperties','boundingbox');
C = vertcat(S.Centroid);
ara = vertcat(S.Area);
sol = vertcat(S.Solidity);
mind = vertcat(S.MinFeretDiameter);
% try to discern which objects are which
[~,isring] = max((mind.*ara)./sol);
[~,isbox] = max(sol./ara);
% get relevant geometry information
Bring = S(isring).BoundingBox;
Cring = Bring(1:2) + Bring(3:4)/2; % center of ring
Cbox = C(isbox,:); % center of box
Rring = mind(isring)/2; % approx radius of ring
boxloc = (Cbox-Cring).*[1 1]; % relative position vector of box
baserot = atan2d(boxloc(2),boxloc(1)); % base rotation of entire object
% these are constants extracted from a modified copy
% radii are relative to Rring
Rbox = 0.666*Rring;
Rdot = 0.786*Rring;
rdot = 0.0541*Rring;
dotrot = 30;
% find dot locations
dotangles = linspace(dotrot,180-dotrot,9) + [0; 180] + baserot;
dotangles = flipud(reshape(dotangles.',[],1));
[dotx doty] = pol2cart(dotangles*pi/180,Rdot);
dotxy = [dotx doty] + Cring;
imshow(A); hold on
plot(dotxy(:,1),dotxy(:,2),'o')
%% try to get samples by doing a local average around each point
dotsamples = zeros(numel(dotx),1);
sc = size(A);
[xx yy] = meshgrid(1:sc(2),1:sc(1));
for k = 1:numel(dotx)
mask = ((xx-dotxy(k,1)).^2+(yy-dotxy(k,2)).^2) <= (rdot*0.5)^2;
samplesize = nnz(mask);
dotsamples(k) = sum(A(mask))/samplesize;
end
dotsamples
dotsamples = 18×1
124.4250 127.0500 129.5641 131.3846 136.0250 132.4872 134.2973 134.8649 137.3250 139.5750
This should work even if the object is rotated or moved/scaled. The thresholding might need to be adjusted if the contrast and levels change.
##### 7 CommentsShow 5 older commentsHide 5 older comments
DGM on 16 Feb 2022
Edited: DGM on 16 Feb 2022
Regarding the question you asked about how the locations are determined:
Once a new image is thresholded, there are only two objects of interest -- the dark ring near the perimeter, and the dark off-center square. The location of the phantom is represented by the center of the ring. The rotation of the phantom is determined by the location of the square with respect to the ring center.
Once the center and base rotation angle are known, the angular displacement of the dots is known. The angle between the square and the first dot is 30d. The angle between each dot is 15d. The radial position of the dots is known relative to the ring radius, so the ring radius will resolve that.
In short, the only information that's needed from the incoming image are the ring center and radius, and the square's location.
Of course, one might ask where those precalculated constants came from. I started by editing a copy of the original image to make it easier to segment.
% extract geometry constants from simplified copy
A = imread('diskxrt.png');
B = bwareafilt(A<128,3);
imshow(B)
S = regionprops(B,'area','solidity','centroid','minferetproperties','boundingbox');
C = vertcat(S.Centroid);
ara = vertcat(S.Area);
sol = vertcat(S.Solidity);
mind = vertcat(S.MinFeretDiameter);
[~,isring] = max((mind.*ara)./sol);
[~,isdot] = max(sol./ara);
isbox = 3;
Bring = S(isring).BoundingBox;
Cring = Bring(1:2) + Bring(3:4)/2; % ring center calculated as center of bbox
%Cring = C(isring,:); % ring center as calculated from image centroid
Cbox = C(isbox,:);
Rring = mind(isring)/2;
boxloc = (Cbox-Cring).*[1 -1];
baserot = atan2d(boxloc(2),boxloc(1));
% these are relative to Rring (as calc via Cring/Bring)
Rbox = norm(Cbox-Cring)/Rring % 0.6474 0.6661
Rbox = 0.6661
Rdot = norm(C(isdot,:)-Cring)/Rring % 0.7670 0.7853
Rdot = 0.7853
rdot = S(isdot).MinFeretDiameter/2/Rring % 0.0541 0.0541
rdot = 0.0541
dotloc = (C(isdot,:)-Cring).*[1 -1];
dotrot = atan2d(dotloc(2),dotloc(1)) - baserot % 30.285 29.499
dotrot = 29.4990
dotangles = linspace(dotrot,180-dotrot,9) + [0; 180]
dotangles = 2×9
29.4990 44.6242 59.7495 74.8747 90.0000 105.1253 120.2505 135.3758 150.5010 209.4990 224.6242 239.7495 254.8747 270.0000 285.1253 300.2505 315.3758 330.5010
hold on;
plot(C(isring,1),C(isring,2),'yo')
plot(C(isdot,1),C(isdot,2),'*')
plot(C(isbox,1),C(isbox,2),'s')
sesilia maidelin on 17 Feb 2022
ah i see, thank you so much for your explanation! huge appreciation!

Sign in to comment.

### Categories

Find more on Agriculture in Help Center and File Exchange

R2021a

### Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!