How to find the inclination of the back bone in the attached image?

1 view (last 30 days)
Hello Everyone... Greetings...
I am working with some medical imaging problems.
I want to rotate the given image after finding the inclination angle,
So that backbone looks straight and performs further steps.
I have attached the image for the reference.
The image is high resolution with low quality.
Can anyone suggest, how to find the inclination of the backbone?
Thank you.

Accepted Answer

Ryan Comeau
Ryan Comeau on 19 May 2020
Edited: Ryan Comeau on 20 May 2020
Hello,
I am thinking the best way to do this will be analog. What i mean is just choose two sets of coordinates, compute their slope and apply the rotation yourself. It will be the least painful if you don't need to repeat this frequently.
Doing this autonomously will be quite a bit more difficult. There are ways to get it done, but guarenteeing good performance over your analog work will be difficult. Here is the way to do this analog:
image=imread('path/to/image');
%you choose 2 coordinate sets on the image
coord_1=[x1,y1];
coord_2=[x2,y2];
slope=(coord_1(2)-coord_2(2))/(coord_1(1)-coord_2(1)); %simple slope calculation
%now we convert this to degrees.
%assuming 90 is vertical and 0 is horizontal
slope_deg=atand(slope);
rotation_needed=90-slope_deg;
%perform the rotation
image=imrotate(image,rotation_needed);
Now to do this autonomouly will require some work. You will need to apply an image processing technique to segment the objects present, and identify the spine. You'll need to play around with the methods to see which one works the best. Here are some ideas that I have:
  1. Use a hough transform,this algorithm identifies lines in an image(see: https://www.mathworks.com/help/images/ref/hough.html). Locate the longest line(which should be the spine, calculate its slope and perform the rotation). Try processing without the image pre processing ans see if you can detect the spine.
%code for the Hough transform
%this block gives the endpoints of the lines it locates.
%try the Hough transform without the pre processing as well.
image=imread('path/to/image/')
image=rescale(image,0,1);%helps in binarization
BW_im=imbinarize(image);
BW_finel=bwskel(BW_im); %create image skeleton
[H,theta,rho] = hough(BW_final);
P = houghpeaks(H,45,'threshold',ceil(0.1*max(H(:))));
lines = houghlines(BW_im,theta,rho,P,'FillGap',5,'MinLength',7);
%lines is a struct, stick this in a for loop to extract
for i=1:length(lines)
xy(i,:)=[lines(i).point1,lines(i).point2];
end
%parse for the longest one to identify the spine.
2. Use Haar like features
All of these methods except 4, will require you to parse the output and figure out the slope of the spine, then apply the rotation. A Haar like feature(just google these), can be used to identify objects. Here, we will construct one for the spine, and rotate it in intervals of 5 degrees(you can change this). The Haar like features will convolute the image, and whichever angle yields the highest value(from the convolution) will give a +-2.5 degrees on your spine angle. The following block of code explains how to use a Haar like feature and the logic to determine the rotation.
image=imread('path/to/image');
image=rescale(image,0,1); %helps with speed and other things
spine_average_width=30; %this is in pixels, change for your need(tune it)
spine_average_height=100;% in pixels change for your need(tune it)
spine_haar_feature(1:spine_average_height,1:2*spine_average_width)=zeros();
spine_haar_feature(1:spine_average_height,...
spine_average_width/2:(3/2)*spine_average_width)=ones();
%we have created a binari image with a vertical bar. this is our Haar like feature
%you'll need to add a little routine to make sure the rotation is upwards
for i=1:18 %starting from vertical going horizontal
spine_angle(i)=max(conv2(image,imrotate(spine_haar_feature,-i.*5)));
end
[~,index]=max(spine_angle(:));
angle=index*5; %because 5 degrees per iteration
image=imrotate(image,angle);
I hope this answer is not overbearing, i just like medical imaging problems and think they're fun to solve.
Hope this helps,
RC
  2 Comments
Mrutyunjaya Hiremath
Mrutyunjaya Hiremath on 19 May 2020
Edited: Mrutyunjaya Hiremath on 19 May 2020
Hello Ryan Comeau,
Thanks for your time and answer.
I love the way of your coding and content presentation. Really it's nice.
I've worked out the same manual method (like, your first method).
close all;
[fileName folderName] = uigetfile({'*.jpg'}, 'Select an Image for Classification', '..\Dataset\Test');
imgTest = imread([folderName fileName]);
[~,name,~]= fileparts(fileName);
figure,imshow(imgTest);title('Query Image');
N = 5;
pts = readPoints(imgInput, N);
% fit straight line
x = pts(1,:); y = pts(2,:);
coeffs = polyfit(x, y, 1);
slope = coeffs(1);
intercept = coeffs(2);
% Get fitted values
fittedX = linspace(min(x), max(x), 200);
fittedY = polyval(coeffs, fittedX);
hold on;
plot(fittedX, fittedY, 'b-', 'LineWidth', 2);
% slope and angle of line
x1 = fittedX(1);
y1 = fittedY(1);
x2 = fittedX(end);
y2 = fittedY(end);
slope = (y2 - y1) ./ (x2 - x1); % just for better accuracy
angle = atand(slope);
% rotate image
if (angle < 0)
angle = -angle;
imgRotate = imrotate(imgInput, 90-angle, 'bicubic');
else
imgRotate = imrotate(imgInput, -(90-angle), 'bicubic');
end
figure, imshow(imgRotate);
% get center part of image
X1 = size(imgRotate,1)*cosd(angle);
Y1 = size(imgRotate,2)*cosd(angle);
W = size(imgRotate,1) - 2*Y1;
H = size(imgRotate,2) - 2*X1;
imgCrop = imcrop(imgRotate, [X1, Y1, H, W]);
figure, imshow(imgCrop);
Got good reuslts. like
But, the automotive method is good (the reason is, there are 100's of images :).
let me try your 3rd and 4th method.
Once again thank you, for the answer.

Sign in to comment.

More Answers (1)

Image Analyst
Image Analyst on 20 May 2020
I disagree with the accepted answer. Perhaps he's never heard of the radon transform - the Nobel prize winning medical imaging algorithm that is at the heart of CT and MRI. You can do this automatically just by taking the radon transform and looking for the peak.
See my attached demo. Adapt as needed by applying the FAQ: FAQ How_can_I_process_a_sequence_of_files? for your hundreds of images.
  4 Comments
Mrutyunjaya Hiremath
Mrutyunjaya Hiremath on 20 May 2020
Hello Ryan and Image Analyst,
Really it's a good discussion and explanation.
From our side, elderly people say, ‘we lose nothing, even when we are quarreling with good and knowledgeable people’. That’s true. J.
Learned a lot from this.
Thank you both of you.
Image Analyst
Image Analyst on 22 Aug 2020
Girish, start a new question - let's not hijack this one. Post your image there and a definition of what direction you consider to be the "height".

Sign in to comment.

Products


Release

R2013a

Community Treasure Hunt

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

Start Hunting!