hello i need clearer circle after performing the below code help me plz

1 visualización (últimos 30 días)
clear all
close all
[FileName,PathName] = uigetfile('*.*','Select the image file');
nomfich=[PathName,FileName];
ima=imread(nomfich);
HSV=rgb2hsv(ima);
H=HSV(:,:,1);
S=HSV(:,:,2);
V=HSV(:,:,3);
%dilate image
SE = strel('line',1,1);
H = imdilate(H,SE);
S = imdilate(S,SE);
V = imdilate(V,SE);
%performing dwt2 on h,s,v channels
[hca1,hch1,hcv1,hcd1] = dwt2(H,'bior1.3');
[sca1,sch1,scv1,scd1] = dwt2(S,'bior1.3');
[vca1,vch1,vcv1,vcd1] = dwt2(V,'bior1.3');
%replacing ll(hac1) with adding all other three values(hch1,hcv1,hcd1). sililarly for s and v channls ll band is replaed with new value
LH=hch1+hcv1+hcd1;
LS=sch1+scv1+scd1;
LV=vch1+vcv1+vcd1;
reconstructing image using replaced ll value
a=idwt2(LH,hch1,hcv1,hcd1,'bior1.3');
%figure,imshow(a);
b=idwt2(LS,sch1,scv1,scd1,'bior1.3');
%figure,imshow(b);
c=idwt2(LV,vch1,vcv1,vcd1,'bior1.3');
%figure,imshow(c);
%hsvimg=a+b+c;
%combing all the three result into hsv1 image for display
HSV1(:,:,1)=a;
HSV1(:,:,2)=b;
HSV1(:,:,3)=c;
figure,imshow(HSV1);
title('reconstructed HSV image');
rgbimg=hsv2rgb(HSV1);
figure,imshow(rgbimg);
%imhist(rgbimg);
title('reconstructed rgb image');
grayimg=rgb2gray(rgbimg);
figure,imshow(grayimg);
%imhist(grayimg);
%perfotming otsu
level = graythresh(grayimg);
BW = im2bw(grayimg,level);
imshow(BW);

Respuesta aceptada

John BG
John BG el 9 de Abr. de 2018
Editada: John BG el 10 de Abr. de 2018
Hi Shivu Shetty
to acquire center and surface points on the image, use the attached ginputRed.m that has red crosshairs, as the standard ginput.m does not allow to change color, and there is no switch in figure, axes or plot handles to change crosshair color.
More images uploaded as soon as 10 uploads per day cap reloaded.
1.-
acquiring one of your start images:
clear all;close all;clc
A=imread('001.jpg');
A1=A(:,:,1); % A2=A(:,:,2);A3=A(:,:,3);
hf1=figure(1);imshow(A1);ax1=hf1.CurrentAxes;
% 002
2.-
adjusting contrast
A02=imadjust(A1,[.81;.89],[0;1]);
hf2=figure(2);h_im=imshow(A02);ax2=hf2.CurrentAxes;
% 003
3.-
trying surf
dx1=1;
[X,Y]=meshgrid([1:dx1:size(A1,2)],[1:1:size(A1,1)]);
surf(X,Y,255-A02,'EdgeColor','none')
% 004-1
.
.
surf has all rings clipping the roof of 255, which is good to count tree rings.
surf(255-A02) better than surf(A02)
.
.
yet one wonders if the rings could be a bit thicker to count them and if it is possible to have all rings with the highest amplitude, to avoid false counts, or missing rings
4.-
trying interpolation in case there's need for more samples
some ring peaks look too sharp with only one sample and many rings do really have a low level:
dx2=.1;
[Xq,Yq]=meshgrid([1:dx2:size(A1,2)],[1:dx2:size(A1,1)]);
A02q=interp2(X,Y,double(255-A02),Xq,Yq);
hf3=figure(3);surf(Xq,Yq,255-A02q,'EdgeColor','none');
comparing,
% 005
% 006
.
.
no discernible difference, still some peaks to sharp, some other peaks too low
5.-
let's carry on with fence-like rings as peaks, better than narrow trenches in A02
A03=A02q;
% 004-2
6.-
thresholding
there are some really low yet valid peaks too close to false peaks
A03(A03>=125)=255;A03(A03<125)=0;
7.-
trying single straight cross section, centre to outer crust
hf4=figure(4);hsf2=surf(Xq,Yq,A03,'EdgeColor','none');
ax4=hf4.CurrentAxes
ax4.DataAspectRatio=[1 6 1]
ax4.DataAspectRatio=[1 1 1]
campos(ax4,[250 250 3664])
how to fix ginput figure axes and plot handles not having switch for crosshair colour:
.
p=ginputRed(2)
p=floor(p)
nx=floor(linspace(p(1,1),p(2,1),max(abs(p(1,1)-p(2,1)),abs(p(1,2)-p(2,2)))))
ny=floor(linspace(p(1,2),p(2,2),max(abs(p(1,1)-p(2,1)),abs(p(1,2)-p(2,2)))))
hold(ax4,'on');plot(ax4,nx,ny,'r-');plot(ax4,[p(1,1) p(2,1)],[p(1,2) p(2,2)],'ro')
% 007
.
.
L=[];
for k=1:1:numel(nx)
L=[L A1(ny(k),nx(k))]
end
figure(5);plot(L);grid on
figure(6);plot(diff(L));grid on
.
cross sections shown in 005 and 006 above
it still doesn't get the count to the comfortable place where rings are clear tall peaks and everything else far enough below to avoid confusion.
Yet, surf has coloured the rings clearly yellow, away to everything else in magenta, let's exploit this
8.-
Isolating yellow pixels
saveas(hf4,'tree_rings_yellow.jpg')
%%%%
close all;clear all;clc
A20=imread('tree_rings_yellow.jpg');
hf1=figure(1);imshow(A20);ax1=hf1.CurrentAxes
% 008
.
.
A13=A20(:,:,1)-A20(:,:,3); % A23=A20(:,:,2)-A20(:,:,3);
[x13,y13]=find(A13>30);
% figure(2);imshow(A)
hold(ax1,'all');figure(2);h1=plot(ax1,y13,x13,'r.')
% 008
% ax1.DataAspectRatio=[9 9 1]
% axis(ax1,[791 1352 732 1201])
axis(ax1,[325 684 349 686])
the previous was marking ring pixels red, just to make sure marking the right pixels
9.-
solving the too short peaks, we need tall peaks against everything else for such purpose let's binarize
% BW better contrast
B=zeros(size(A20(:,:,1)));
for k=1:1:numel(x13)
B(x13(k),y13(k))=1;
end
hf2=figure(2);imshow(B);
ax2=hf2.CurrentAxes
% axis(ax2,[791 1352 732 1201])
axis(ax2,[325 684 349 686])
% 009
10.-
Now let's use flexible cross section, let's zig-zag, no need to count rings while pulling the cross section line,
but let's make sure that each straight section of cross section clearly gets one or more tree ring and if not just for clarity, that new sections start in between cross sections
p2=ginputRed
p2=floor(p2)
n=[p2(1,1) p2(1,2)];
for k=1:1:size(p2,1)-1
nx_section=floor(linspace(p2(k,1),p2(k+1,1),max(abs(p2(k,1)-p2(k+1,1)),abs(p2(k,2)-p2(k+1,2)))));
ny_section=floor(linspace(p2(k,2),p2(k+1,2),max(abs(p2(k,1)-p2(k+1,1)),abs(p2(k,2)-p2(k+1,2)))));
n=[n;nx_section' ny_section'];
end
n(1,:)=[];
hold(ax2,'all');
plot(ax2,p2(:,1),p2:,2),'r-','LineWidth',4);
plot(ax2,[p2(1,1) p2(end,1)],[p2(1,2) p2(end,2)],'ro','LineWidth',4);
% 010
variable peaks of interest amplitude solved, now all rings are 1s, and all space between rings are 0s
L=[];
for s=1:1:size(n,1)
L=[L B(n(s,2),n(s,1))];
end
mind the gap attempting the reading loop with B(n(s,1),n(s,2)) a result that looks ok is reached but neither the amount of rings nor the initial values of the cross section are correct.
figure(3);plot(L','r','LineWidth',3);grid on
% 011
.
counting peaks
.
[pks,locs]=findpeaks(L)
pks =
Columns 1 through 9
1 1 1 1 1 1 1 1 1
Columns 10 through 17
1 1 1 1 1 1 1 1
locs =
Columns 1 through 9
60 118 165 210 260 297 319 338 364
Columns 10 through 17
375 382 397 411 421 434 442 453
.
And the result is:
.
amount_rings=numel(pks)
=
17
.
Shivu
if you find this answer useful would you please be so kind to consider marking my answer as Accepted Answer?
To any other reader, if you find this answer useful please consider clicking on the thumbs-up vote link
thanks in advance for time and attention
John BG
.
  3 comentarios
shivu shetty
shivu shetty el 11 de Abr. de 2018
i dont have high resolution images. is there any other way to count those rings?
John BG
John BG el 11 de Abr. de 2018
Editada: John BG el 11 de Abr. de 2018
the narrowest case where we can tell that there's ring is when
'.. 101010101 ..'
how can one tell there's a ring if the very first data we get hands on is:
' .. 111111 ..'
.
?
The answer is, we cannot tell, at least not with such low amount of data.
No matter how smart a worker is, if the boss asks a worker to guess data, there are 6 possible ways for a worker to address such difficult query:
1.- Gamble one's job on the line by supplying an answer where there's only way for guessing.
Not recommended.
2.- Have a go and try your boss gambling his / her job by asking for a pair of dice, signed and approved, so when the hoax hits the fan, you can bring forward evidence that the decision to gamble was your boss's not yours.
3.- Convince some one else to take the hit for you, But again not recommended, because the worker is trying to do what the boss did to the worker.
The chances are the whoever the worker (you) find(s), the boss will have arranged in advance to pretend to help you.
But right before the critical moment of bringing forward evidence to the boss of the boss, everyone will suddenly point at you, and the mob will recklessly push the worker under the train.
It's a Microsoft world, we only live in it:The boss can always read your emails, you can't read the boss emails unless he/she wants to or pretends to.
4.- Convince your boss to give you a better camera, with more pixels in the CCD so you can get at least get
' .. 1010101 ..'
across the narrowest cross section. This is a common practice in huge budget projects: when not sure delay, spend more time and resources hopping that somehow a solution will fall down from a tree branch, like Newton's apple, needless to recall that such story was popular fiction that ended up in the books: Newton based his initial studies dropping apples 'from shoulders of giants' meaning from the roofs of high buildings with known heights.
5.- Bring in your own far better camera, that you have purchased with your own money, with the best of intentions.
I once did this one, in Spain for the contractor of a mobile operator that told me to assess base stations and told me to use my own camera when I asked for a company camera:
After solving the technical problem I almost end up having to leave my camera in the premises when the bad boss had a go on 'health & safety' and 'company regulations'.
Then is when the Union comes in really really handy.
Never ever do this in security jobs, especially if you have signed down compliance to any kind of government act.
6.- Tell boss to do it him/herself. Or tell the boss of the boss that the boss is an incompetent.
In either situation it's worth having the hand ready on the eject handle, and be certain you did all the pre-flight checks that included parachute, map, compass, knife, rations, .. in case one has to walk to next job.
Regards
John BG

Iniciar sesión para comentar.

Más respuestas (1)

Image Analyst
Image Analyst el 8 de Abr. de 2018
Contact the Laboratory of Tree Ring Research at the University of Arizona. http://ltrr.arizona.edu/ They started the field and are probably still the world's leader. They can probably point you to their best publications. You might even be able to obtain code from them.

Community Treasure Hunt

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

Start Hunting!

Translated by