update the new value

1 view (last 30 days)
Jing Hua
Jing Hua on 31 Jul 2012
Hi
From the beginning of my code, I have original x_cell and an y_cell, by performing a series of calculation, some of sub-cells in x_cell and y_cell are deleted.
Now I want to use the new x_cell and y_cell to do the same calculation to get a better approximation.
when I get the new x_cell and new y_cell, matlab shows the dimension error (because some arrays are deleted in x_cell and y_cell).
How do I update the new x_cell and y_cell in some stage of the code?
Code: clc, clear all %Load all the pictures in a file files = dir('*.PNG'); num_files = numel(files); images = cell(1, num_files); %Convert them into binary image and define their boundaries black_white_pic=cell(1, num_files); boundary=cell(1, num_files); %All variables needed for the whole program num_ele_boundary=zeros(1,num_files); %Big loop containing information needed for ABAQUS
for k = 1:num_files
images{k} = imread(files(k).name);
black_white_pic{k}=im2bw(images{k},graythresh(images{k}));
boundary{k}=bwboundaries(black_white_pic{k});
num_ele_boundary(k)=length(boundary{k});
%%Generate x & y coordinate
for Q=1:length(boundary{k})
x_cell{k}{Q}=smooth(boundary{k}{Q}(:,2),'lowess');
y_cell{k}{Q}=-smooth(boundary{k}{Q}(:,1),'lowess');
%Eliminate the noise
actual_area{k}{Q}=0.5*abs(sum(x_cell{k}{Q}(1:(end-1)).*y_cell{k}{Q}(2:end)-y_cell{k}{Q}(1:(end-1)).*x_cell{k}{Q}(2:end)));
deleted_positions_1{k}=cellfun(@(x) x<15,actual_area{k},'un',0);
deleted_positions{k}=find(cell2mat(deleted_positions_1{k})==1);
for position_6=1:length(deleted_positions{k})
x_cell{k}{deleted_positions{k}(position_6)}={};
y_cell{k}{deleted_positions{k}(position_6)}={};
end
end
for q=1:length(x_cell{k})
x_cell{k}(cellfun(@isempty,x_cell{k}))= [];
y_cell{k}(cellfun(@isempty,y_cell{k}))= [];
end
%%Find ellipses (centre, semi major and minor axis & orientation)
M_1_0{k}=cellfun(@(v)v*2, x_cell{k}(:),'un',0);
M_1{k}=cellfun(@times, y_cell{k}(:),M_1_0{k}(:), 'un',0);
M_2{k}=cellfun(@(v) v.^2, y_cell{k}(:), 'un',0);
M_3{k}=cellfun(@(v) v*2, x_cell{k}(:),'un',0);
M_4{k}=cellfun(@(v) v*2, y_cell{k}(:),'un',0);
M_5{k}=cellfun(@(v) v.^0, x_cell{k}(:), 'un', 0);
M_6_1{k}=cellfun(@(v) v.^2, x_cell{k}(:), 'un',0);
M_6{k}=cellfun(@(v) v./-1, M_6_1{k}(:), 'un', 0);
M{k}=strcat(M_1{k}, M_2{k}, M_3{k}, M_4{k}, M_5{k});
M_a{k}=cellfun(@pinv, M{k}(:),'un',0);
e{k}=cellfun(@mtimes, M_a{k}(:),M_6{k}(:),'un',0);
for j=1:length(e{k})
e_a{k}(j)=cat(2,(e{k}(j)));
end
e_a_1{k}=cell2mat(e_a{k});
end
%%Regroup the 'e ' matrix
for i=1:5
for p=1:num_files
z{i}{p}=e_a_1{p}(i,:);
end
end
b=z{1};c=z{2};d=z{3};f=z{4};g=z{5};
for N=1:num_files
%%Centre of x and y coordinates
delta_1{N}=b{N}.^2;
delta{N}=-c{N}+delta_1{N};
x_centre_1{N}=c{N}.*d{N};
x_centre_2{N}=b{N}.*f{N};
x_centre_3{N}=x_centre_1{N}-x_centre_2{N};
x_centre{N}=x_centre_3{N}./delta{N};
y_centre_1{N}=f{N};
y_centre_2{N}=b{N}.*d{N};
y_centre_3{N}=y_centre_1{N}-y_centre_2{N};
y_centre{N}=y_centre_3{N}./delta{N};
%%Length of semi-major, semi-minor axis
nom_1{N}=f{N}.^2;
nom_2_1{N}=d{N}.^2;
nom_2{N}=c{N}.*nom_2_1{N};
nom_3_1{N}=b{N}.^2;
nom_3{N}=g{N}.*nom_3_1{N};
nom_4_1{N}=2*b{N}.*d{N};
nom_4{N}=nom_4_1{N}.*f{N};
nom_5{N}=c{N}.*g{N};
nom{N}=2*(+nom_1{N}+nom_2{N}+nom_3{N}-nom_4{N}-nom_5{N});
s_1{N}=4*nom_3_1{N}(:); s_2{N}=1-c{N}(:); s_3{N}=s_2{N}(:).^2; s_4{N}=s_1{N}(:)./s_3{N}(:); s_5{N}=1+s_4{N}(:); s{N}=sqrt(s_5{N}(:));
a_prime_1{N}=c{N}(:)-1; a_prime_2{N}=a_prime_1{N}(:).*s{N}(:); a_prime_3{N}=c{N}(:)+1; a_prime_4{N}=a_prime_2{N}(:)-a_prime_3{N}(:); a_prime_5{N}=delta{N}(:).*a_prime_4{N}(:); a_prime_6{N}=nom{N}(:)./a_prime_5{N}(:); a_prime{N}=sqrt(a_prime_6{N}(:)); b_prime_1{N}=-c{N}(:)+1; b_prime_2{N}=b_prime_1{N}(:).*s{N}(:); b_prime_3{N}=c{N}(:)+1; b_prime_4{N}=b_prime_2{N}(:)-b_prime_3{N}(:); b_prime_5{N}=delta{N}(:).*b_prime_4{N}(:); b_prime_6{N}=nom{N}(:)./b_prime_5{N}(:); b_prime{N}=sqrt(b_prime_6{N}(:)); semimajor_axis{N}=max(a_prime{N},b_prime{N}); semiminor_axis{N}=min(a_prime{N},b_prime{N}); %find the orientation of the ellipse (positive anti-clockwise) for T=1: length(x_cell{N}) if (b{N}(T)==0) & (1<c) angle{N}(T)=0 end if (b{N}(T)==0) & (1>c) angle{N}(T)=pi/2 end if 1<abs(c{N}(T)) angle{N}(T)=0.5*(acot((1-c{N}(T))/(2*b{N}(T)))); end if 1>c{N}(T) angle{N}(T)=pi/2+0.5*(acot((1-c{N}(T))/(2*b{N}(T)))); end end
%%Error Based on Area for P=1:length(x_cell{N}) actual_area{N}{P}=0.5*abs(sum(x_cell{N}{P}(1:(end-1)).*y_cell{N}{P}(2:end)-y_cell{N}{P}(1:(end-1)).*x_cell{N}{P}(2:end))); model_area{N}{P}=pi*semimajor_axis{N}(P)*semiminor_axis{N}(P); error_based_on_area{N}{P}=((model_area{N}{P}-actual_area{N}{P}))/actual_area{N}{P}*100; grouped_error{N}=cell2mat(error_based_on_area{N});
% for graph=1:length(x_cell{1}) % plot(x_cell{3}{graph},y_cell{3}{graph}); % axis equal % hold on % end
%Address first complication-two or more ellipse interval_1=6; %This is the first number to be adjusted (interval should be a even number) %position_1 is the position of the first point position_1{N}{P}=1:interval_1:((round(((length(x_cell{N}{P}))/interval_1)))*interval_1-interval_1); %position_2 is the position of the second point position_2{N}{P}=(1+interval_1):interval_1:((round(((length(x_cell{N}{P}))/interval_1)))*interval_1); %position_3 is the position of the point is the middle of position_1 and 2 position_3{N}{P}=(position_1{N}{P}+position_2{N}{P})./2;
%perimeter vector perimeter_vector_x_dis{N}{P}=x_cell{N}{P}(position_2{N}{P})-x_cell{N}{P}(position_1{N}{P}); perimeter_vector_y_dis{N}{P}=y_cell{N}{P}(position_2{N}{P})-y_cell{N}{P}(position_1{N}{P}); perimeter_vector_x_1{N}{P}=x_cell{N}{P}(position_1{N}{P}); perimeter_vector_x_2{N}{P}=x_cell{N}{P}(position_2{N}{P}); perimeter_vector_y_1{N}{P}=y_cell{N}{P}(position_1{N}{P}); perimeter_vector_y_2{N}{P}=y_cell{N}{P}(position_2{N}{P}); perimeter_vector{N}{P}=[perimeter_vector_x_dis{N}{P} perimeter_vector_y_dis{N}{P}];
%centre vector centre_vector_x_1{N}{P}=x_centre{N}(P); centre_vector_x_2{N}{P}=x_cell{N}{P}(position_3{N}{P}); centre_vector_y_1{N}{P}=y_centre{N}(P); centre_vector_y_2{N}{P}=y_cell{N}{P}(position_3{N}{P}); centre_vector_x_dis{N}{P}=centre_vector_x_2{N}{P}-centre_vector_x_1{N}{P}; centre_vector_y_dis{N}{P}=centre_vector_y_2{N}{P}-centre_vector_y_1{N}{P}; centre_vector{N}{P}=[centre_vector_x_dis{N}{P} centre_vector_y_dis{N}{P}];
%x-axis vector x_axis_vector=[1 0];
%angle between two vectors dot_product_1{N}{P}=perimeter_vector{N}{P}(:,1).*centre_vector{N}{P}(:,1)+perimeter_vector{N}{P}(:,2).*centre_vector{N}{P}(:,2); arc_first_angle{N}{P} = dot_product_1{N}{P}./sqrt((centre_vector{N}{P}(:,1)).^2+(centre_vector{N}{P}(:,2)).^2)./sqrt((perimeter_vector{N}{P}(:,1)).^2+(perimeter_vector{N}{P}(:,2)).^2); first_angle{N}{P}=acos(arc_first_angle{N}{P})*180/pi;
%angle between second vector and x-axis dot_product_2{N}{P}=x_axis_vector(1)*centre_vector{N}{P}(:,1)+x_axis_vector(2)*centre_vector{N}{P}(:,2); arc_second_angle{N}{P} = dot_product_2{N}{P}./sqrt((centre_vector{N}{P}(:,1)).^2+(centre_vector{N}{P}(:,2)).^2)/sqrt((x_axis_vector(1))^2+(x_axis_vector(2))^2); second_angle{N}{P}=acos(arc_second_angle{N}{P})*180/pi;
%the difference between two angles in between two vectors for position_4=1:(length(first_angle{N}{P})-1); difference_in_y_angle{N}{P}(position_4)=-first_angle{N}{P}(position_4+1)+first_angle{N}{P}(position_4); end
threshold=50; %this is the second parameter needs to be adjusted (check if results is consistence with image) critical_points_position{N}{P}=find(difference_in_y_angle{N}{P}>threshold); number_of_critical_points{N}{P}=length(critical_points_position{N}{P}); critical_coordinates_position{N}{P}=(critical_points_position{N}{P}-1).*interval_1; critical_point_coordinates_x{N}{P}=x_cell{N}{P}(critical_coordinates_position{N}{P}); critical_point_coordinates_y{N}{P}=y_cell{N}{P}(critical_coordinates_position{N}{P}); critical_point{N}{P}=cat(2,critical_point_coordinates_x{N}{P},critical_point_coordinates_y{N}{P});
%linear interpolation if number_of_critical_points{N}{P}==2 linear_slop{N}{P}=(critical_point{N}{P}(1,2)-critical_point{N}{P}(2,2))/(critical_point{N}{P}(1,1)-critical_point{N}{P}(2,1)); lineat_constant{N}{P}=linear_slop{N}{P}*critical_point{N}{P}(1,1)+critical_point{N}{P}(1,2);
%Get the additional coordinates interval_2=1; %this is the third parameter needs to be adjusted (this is the interval between two consecutive additional x_coordinates additional_x_coordinates{N}{P}=critical_point_coordinates_x{N}{P}(1):interval_2:critical_point_coordinates_x{N}{P}(2); additional_y_coordinates{N}{P}=additional_x_coordinates{N}{P}.*linear_slop{N}{P}+lineat_constant{N}{P};
end
%delete the corresponding coordinates in x_cell & y_cell deleted_position_3{N}=critical_points_position{N}; emptyIndex = cellfun(@isempty,deleted_position_3{N}); nonemptyIndex=~emptyIndex; deleted_position_3{N}(emptyIndex) = {0}; deleted_position_3{N}(nonemptyIndex) = {1}; deleted_positions_2nd{N}=find(cell2mat(deleted_position_3{N})==1); for position_7=1:length(deleted_positions_2nd{N}) x_cell{N}{deleted_positions_2nd{N}(position_7)}={}; y_cell{N}{deleted_positions_2nd{N}(position_7)}={}; end x_cell{N}(cellfun(@isempty,x_cell{N}))= []; y_cell{N}(cellfun(@isempty,y_cell{N}))= []; end end
  3 Comments
Jan
Jan on 31 Jul 2012
@Jing Hua: This is a big piece of code. It is unlikely that somebody will analyze it to suggest the wanted modifications - even if you will format it correctly (please read the instructions in the "About Matlab Answers" link again). I stop reading the code whenever I find a "clear all".
I suggest to post the relevant part only and explain the question again with more details: "How do I update the new x_cell and y_cell in some stage of the code?"

Sign in to comment.

Answers (0)

Tags

Community Treasure Hunt

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

Start Hunting!

Translated by