Это довольно сложная задача с моим кодом. Прежде всего, код, который я помещаю здесь, не работает, потому что я использую лист Excel (но я рад отправить его по электронной почте, если люди хотят попробовать использовать мой код).
У меня есть лист Excel с данными о поперечных волокнах в микроскопическом изображении, которое я взял. Информация в основном: location
раздела, area
, angle
вращения.
Отсюда я вычисляю угол ориентации Phi и Gamma. После этого я использую функцию рассеяния для построения точки разных цветов для каждого значения угла Phi. Я использую постоянный цвет для диапазона 10 градусов. Что дает мне такую картину:
Теперь моя цель - рассчитать площадь каждой однородной области. Поэтому я ищу способ построения, чтобы сказать все точки в области -10 +10 (сейчас я делаю 20 градусов, но буду делать 10 после). Я использовал внешний вид, и я получаю следующее:
Белый соответствует тому, где точки находятся в пределах диапазона, который я выбрал. После этого я использую панель инструментов в MATLAB для преобразования каждой точки в пиксель. Поэтому я получаю черный фон с множеством белых пикселей, затем я использую imdilate
для создания кругов, заполнения отверстий и выделения каждого региона определенным цветом. Наконец, я использую границы функций и патча, чтобы создать каждую границу и заполнить их цветом. И я получаю изображение следующим образом:
Это то, что я хочу, и я могу получить область каждого региона и общую площадь (я использовал порог, чтобы отбросить небольшие области). Затем я запускаю код несколько раз для каждого региона, и я использую imfuse
, чтобы вернуть их вместе и посмотреть, как он выглядит.
ПРОБЛЕМА, они перекрываются довольно много, и это потому, что в моих данных есть некоторые ошибки, и поэтому некоторые синие точки будут красными и так далее. Поэтому я хочу запустить код один раз, а затем, когда я перезапущу его с другим диапазоном, он делает то же самое, но не принимает во внимание значение, когда уже что-то закладывается раньше.
Я попытался сделать это, выполнив один раз, сохранив матрицу bw4 и добавив условие при построении черно-белой картинки, сказав, что если Phi находится в моем диапазоне. И здесь нет белого цвета, тогда вы можете поместить белый, иначе он черный. Но, похоже, это не работает.
Я понимаю, что это довольно сложная вещь для объяснения, но я был бы признателен за любые идеи и открывал бы чат по электронной почте или по-другому. Теперь я помещаю полный код, и я могу отправить вам свой лист Excel, если вы хотите запустить его на своем компьютере и убедиться сами.
clearvars -except data colheaders bw4
close all
clc
%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%% CHANGE DATA FOR EACH SAMPLE %%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
cd 'C:\Users\dkarta\Desktop\Sample 12\12.6'
data=xlsread('Sample12_6res.xlsx');
cd 'C:\Users\dkarta\Documents\MATLAB'
%data=Sample121res; % Data name
imax=length(data); % Numbers of rows in data sheet
y=11900; % Number of pixels in the y on image j
%%
data(:,15)=data(:,9)*pi/180; % Convers Column 9 (angle of rotation) in rads
data(:,16)=y-data(:,6); % Reset the Y coordinate axis to bottom left
delta = 0 : 0.01 : 2*pi; % Angle in paramteric equations
theta=45*pi/180; % Sample cutting angle in rads
%AA=[data(:,5)' data(:,16)' phi']
% Define colors
beta=acos(data(1:imax,8)./data(1:imax,7));%./acos(0);
phi=atan(sin(beta).*cos(data(1:imax,15))./(sin(theta)*sin(beta).*sin(data(1:imax,15))+cos(theta)*cos(beta)))/(pi/2);
phi2=phi/2+1/2; % Scales in plane angle phi between 0 and 1
gamma=atan((cos(theta)*sin(beta).*sin(data(1:imax,15))-sin(theta)*cos(beta))./...
(sin(theta)*sin(beta).*sin(data(1:imax,15))+cos(theta)*cos(beta)))/(pi/2);
gamma2=gamma+1/2; % Scales out of plane angle gamma between 0 and 1
%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%% MESHGRID AND COLOURMAP %%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%
x1=data(1:imax,5);
y1=data(1:imax,16);
z1=phi*90;
z2=gamma*90;
n=300;
%Create regular grid across data space
[X,Y] = meshgrid(linspace(min(x1),max(x1),n), linspace(min(y1),max(y1),n));
% Creating a colormap with 10 degree constant colors
map4=[0 0 1;0 1/3 1;0 2/3 1; 0 1 1;0 1 2/3;0 1 1/3;0 1 0;1/3 1 0;2/3 1 0;1 1 0;1 0.75 0;1 0.5 0;1 0.25 0;1 0 0;0.75 0 0.25;0.5 0 0.5;0.25 0 0.75; 0 0 1];
Colormap4=colormap(map4);
h=colorbar;
caxis([-90 90])
set(h, 'YTick', [-90:10:90])
%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%% PLOT USING SCATTER - ISOLATE SOME REGIONS %%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
a=-10; % Lower boundary for angle interval
b=10; % Upper boundary for angle interval
c=z1>a & z1 < b;
c=c.*1;
%j=1;
y1=(y1-min(y1)+1);
y2=max(y1)-y1+1;
[X1,Y1]=meshgrid(1:500,1:500);
griddata(x1,y2,c,X1,Y1);
clear c1
for i=1:imax
if z1(i)< b && z1(i)> a %&& bw4(round(y1(i)),round(x1(i))) == 0
c(i) = 1;
c1(round(y2(i)),round(x1(i)))=1;
else
c(i)= 0;
c1(round(y2(i)),round(x1(i)))=0;
end
end
C=[c c c];
%c(find(c==0)) = NaN;
%contourf(X,Y,griddata(x1,y1,c,X,Y),100,'EdgeColor', 'None')
figure(1), scatter(x1,y1,3,z1,'filled');
axis equal
axis ([0 8000 0 12000])
axis off
figure(2), scatter(x1,y1,3,C,'filled');
axis equal
axis ([0 8000 0 12000])
axis off
se=strel('disk',50,8);
bw2=imdilate(c1,se);
bw4=bwlabel(bw2);
bw3=imfill(bw4,'holes');
max(bw4(:));
figure(3),imshow(c1,'InitialMagnification', 10);
figure(4), imshow(bw2,'InitialMagnification', 10);
figure(5), imshow(bw3,'InitialMagnification', 10);
figure(6),imshow(label2rgb(bw4),'InitialMagnification', 10);
k=ones(max(bw4(:)),1);
clear bw5
for i=1:length(x1)
if bw3(round(y2(i)),round(x1(i))) ~= 0
m=bw3(round(y2(i)),round(x1(i)));
bw5{m}(k(m),1)=x1(i); bw5{m}(k(m),2)=y2(i);
k(m)=k(m)+1;
end
end
figure(7), imshow(~c1,'InitialMagnification', 10);
hold on
for i=1:max(bw4(:))
%scatter(bw5{i}(:,1),bw5{i}(:,2))
j = boundary(bw5{i}(:,1),bw5{i}(:,2),0.5);
%poly=convhull(bw5{i}(:,1),bw5{i}(:,2));
%plot(bw5{i}(poly,1),bw5{i}(poly,2)), title('convhull')
if polyarea(bw5{i}(j,1),bw5{i}(j,2))> 10^5;
patch(bw5{i}(j,1),bw5{i}(j,2),'r'), title('boundary')
indexminy(i)=find(min(bw5{i}(:,2)) == bw5{i}(:,2));
indexminx(i)=find(min(bw5{i}(:,1)) == bw5{i}(:,1));
indexmaxy(i)=find(max(bw5{i}(:,2)) == bw5{i}(:,2));
indexmaxx(i)=find(max(bw5{i}(:,1)) == bw5{i}(:,1));
%xmin = bw5{i}(indexminx); xmax = bw5{i}(indexmaxx);
%ymin = bw5{i}(indexminy); ymax = bw5{i}(indexmaxy);
str=[(indexminx(i)+indexmaxx(i))/2,(indexminy(i)+indexmaxy(i))/2,'Region no.',num2str(i)];
text((min(x1(i))+max(x1(i)))/2,(min(y1(i))+max(y1(i)))/2,str)
polya(i)=polyarea(bw5{i}(j,1),bw5{i}(j,2));
end
end
spolya=sum(polya(:))
print -dpng -r500 B
Просто чтобы показать вам больше фотографий, когда я соберу несколько из них:
И когда я плаваю:
Как вы видите, они перекрываются, чего я не хочу, поэтому я хочу, чтобы каждое изображение, которое я создаю, "знало", что я делаю в предыдущих запусках, чтобы оно не перекрывалось. Я хочу получить процентную область каждого региона, и если они перекрываются, я не могу использовать фактическую общую площадь моего образца, и результаты неверны.