Effacer les filtres
Effacer les filtres

How do I calculate the radius of this culvert from the image hyberbola?

3 vues (au cours des 30 derniers jours)
Hi everyone,
I took a GPR (Ground-Penetrating Radar) B-scan of a roadway section with a round-shaped steel culvert in it.
I am trying to use MATLAB to detect the curve and/or calculate the approximate radius from the image, the actual diameter is 15 inches.
I tried using " imfindcircles" to do it, tried different ways such as increasing the intensity, turning to grayscale, binarizing it (in order to use objectpolarity option), as well as changing the [rmin , rmax] values, but that doesnt seem to work. The figure file, and the raw image texi file are also attached.
I_scan11 = imread("Scan11_color.png");
I_scan11_gray = im2gray(I_scan11);
I_scan11_gray1 = imadjust(I_scan11_gray);
%I_scan11_BW = I_scan11_gray > 210;
I_scan11_BW = imbinarize(I_scan11_gray);
% SE = strel("disk",14); % 20
% I_scan11_BW1 = imdilate(I_scan11_BW,SE);
% %I_scan11_BW1 = imclose(I_scan11_BW,SE);
% %I_scan11_BW1 = bwareaopen(I_scan11_BW1,5);
% imshow(I_scan11_BW1);
[center, radius] = imfindcircles(I_scan11_BW,[400,1120],"Sensitivity",0.99,"ObjectPolarity","dark")
% hold on;
% plot(center(:,1),center(:,2),'yx','LineWidth',2);
% hold off;
The following code generates the image as seen above: (exclusive for GPR)
load MECHANICAL_LAB_EXP_05_11_21_022_P_1A.txt; % Image text file derived from the device
Bscan = MECHANICAL_LAB_EXP_05_11_21_022_P_1A; % This is scan 11-1.
clear MECHANICAL_LAB_EXP_05_11_21_022_P_1A;
% Apparent/input dielectric constant = 6.0, subsurface material is
% non-magnetic ( magnetic permeability = 1)
[a, b] = size(Bscan); % a = cross-range (rx), b = range(r) (rx,r) = (a,b)
dr = (9 * 12) / (b-1); % Range resolution, in inches % Penetration depth is 9 ft
d = 4.11 * 12; % GPR cart rolling distance, in inches % the rolling distance is 4.11 ft
drx = d / (a - 1); % Cross-range resolution, in inches
r = (0:b-1) * dr; % Range, in inches
rx = (0:a-1) * drx; % Cross-range, in inches
imagesc(rx, r, Bscan'); colormap jet; colorbar; caxis([-2.5e6,3e6]);
xlabel('Cross-range, {\it r_x} (in)'); ylabel('Range, {\it r} (in)');
xlim([d-33 d]); ylim([0,22]);
set(gcf, "Position", [744 757 560 420]);
set(gca, 'fontsize', 12);
Please help me! Thanks!
  2 commentaires
DGM le 1 Mar 2022
Modifié(e) : DGM le 1 Mar 2022
I'm not sure I understand what I'm looking at. If the projection of a circular object is a hyperbola, then trying to work off of the curvature might not be as easy as trying to find the asymptotes and minimum range. I don't know if that's appropriate for the available data though. I can't really tell if the image shows a part of a circle or part of a hyperbola.
Koosha le 1 Mar 2022
Modifié(e) : Koosha le 1 Mar 2022
Since the GPR sends EM waves and a steel culvert is a conductor, we have 100% reflection and nothing penetrates below, and metallic objects (culvert, rebar..etc) are depicted as hyperbolas. The range axis depicts the B scan slice depth and the cross-range is the range I rolled the device on the ground, perpendicular to the culvert.
To summarize, the Image is showing a hyperbola, probably thats why the imfindcircles did not work.
Now I am trying to empirically estimate the culvert radius based on the projected hyperbola curve.

Connectez-vous pour commenter.

Réponse acceptée

William Rose
William Rose le 1 Mar 2022
@DGM, this is an inteesting and non-trivial problem. Here are suggesitons, based on my experience using the image processing tools in Labview to determine vl=blood vessel radius from ultrasound images. I have not used the Matlab image processing tools, so my suggesitons are not detailed.
Use the values from the b&w image. Let the user set a depth below which to analyze. Set all array values at shallower depths to zero to avoid getting confused later by the high reflactance values near the surface. In the image you provided, a dpeth of 15" or 16" might be a reasonable depth setting, since the top of the pipe seems to be just below 16". THen you find the depth coordinate for each "column" of the array which is maximum for that column. This is the vertical coordinate of the hyperbola, for that x coordinate. One you have done all that, you have one depth for every x-coordinate. The yu use least squares to find the coefficients of hte hyperbola that give the best fit to the x,y set of points. Then you convert the hyperbola coefficients to a radius.
Good luck.
  8 commentaires
Koosha le 2 Mar 2022
Thank you @William Rose ! it was a very pristine and well-explained code. I am actually busy for a few deadlines but I will come back to you if I had any problems understanding the indexing part.
Either way, you have solved my problem and thank you so much!

Connectez-vous pour commenter.

Plus de réponses (1)

Image Analyst
Image Analyst le 1 Mar 2022
Do you just want to find the peak/ridge line of the curved thing at the bottom of the grayscale image? I could probably find that and then to get the radius of curvature at any point, you just have to fit the data of the ridgeline (x,y) coordinates into a circle fitting function. Of course the radius of curvature changes as you move along the ridgeline so that's why you have to specify both the center location and the window width to get the radius of curvature just within that window.
Like @William Rose asked, please attach your actual data (not a pseudocolored image) in a .mat file with the paperclip icon.
  17 commentaires
William Rose
William Rose le 2 Mar 2022
Script finds peak points in the B-mode scan region of interest.
Script finds the hyperbola that gives the best fit to the peak points in the region of interest.
Script estimates the pipe radius and location that would produce the best fit hyperbola.
Script makes three plots, including the plot below.
Rik le 2 Mar 2022
Modifié(e) : Rik le 2 Mar 2022
[removed by request; water under the bridge]

Connectez-vous pour commenter.


En savoir plus sur Purple dans Help Center et File Exchange

Community Treasure Hunt

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

Start Hunting!

Translated by