Matlab code writting help needed

[prhill]

I am writting a code to do the following..
I inport x y z coordinate points into the program and run my code to set a plane to the points and tell me the distance from to plain to each point. Once I know the distance I delet all points that are a certain distance away (in my case 128mm). The new set of data points is then saved to a file and imported during the second run of the program and in the end cut down even further.

As you can tell I am running a loop here, however I do not know how to write an actual loop in code to make the program delet points 128, 64, 32, 16mm...ect away. I am running thsi up to 20 times so a loop would be very useful. The first two runs of my code follow:

Like I said, basically I need help writting a loop to run instead of repeating all the same code for the "SECOND RUN". At the end of each "RUN" there is an if statement where the abs(plane2point_12(n_12))<64(hylighted to red to find easier) is the only change. I start with 128 decrease is by 1/2 each run to delet points.

I would also like to be able to change the variation at later times so if I is possable it would be great if the decrease of 1/2 could easily be changed to
1/4 or another fraction.

%****************FIRST SCAN, FIRST RUN****************
%load the spreadsheet
points_11=load('C:\Documents and Settings\prhill\My Documents\Hill\Point_Clouds\20_1_00.asc');
x_11=points_11,1);
y_11=points_11,2);
z_11=points_11,3);
%find centroid
N_11=size(points_11,1);
X_11=sum(x_11)/N_11;
Y_11=sum(y_11)/N_11;
Z_11=sum(z_11)/N_11;
centroid_11=[X_11 Y_11 Z_11];
%form basis for matrix elements
h1_11=x_11-X_11; h2_11=y_11-Y_11; h3_11=z_11-Z_11;

%form the elements
k11_11=sum(h1_11.*h1_11); k12_11=sum(h1_11.*h2_11); k13_11=sum(h1_11.*h3_11);
k21_11=sum(h1_11.*h2_11); k22_11=sum(h2_11.*h2_11); k23_11=sum(h2_11.*h3_11);
k31_11=sum(h1_11.*h3_11); k32_11=sum(h2_11.*h3_11); k33_11=sum(h3_11.*h3_11);

%form the matrix
M_11=[k11_11 k12_11 k13_11; k21_11 k22_11 k23_11; k31_11 k32_11 k33_11];

%find the direction vector and value
[vec_11,val_11]=eig(M_11);

%the first column of vec is the unit normal vector and corresponds to the
%minimum eigen value
un_vec_11=vec_11,1);
A_11=un_vec_11(1)
B_11=un_vec_11(2)
C_11=un_vec_11(3)

%calculate D using un_vec and centroid
D_11=-A_11.*centroid_11(1)-B_11.*centroid_11(2)-C_11.*centroid_11(3);

%formula for plane is A*x+B*y+C*z+D=0
%find the z's for the 4 corners of the plane
x_hi_11=max(x_11); x_lo_11=min(x_11);
y_hi_11=max(y_11); y_lo_11=min(y_11);
%z=(-D-A*x-B*y)/C
zc1_11=(-D_11-A_11*x_hi_11-B_11*y_hi_11)/C_11;
zc2_11=(-D_11-A_11*x_hi_11-B_11*y_lo_11)/C_11;
zc3_11=(-D_11-A_11*x_lo_11-B_11*y_lo_11)/C_11;
zc4_11=(-D_11-A_11*x_lo_11-B_11*y_hi_11)/C_11;

%the corner points are
cnr1_11=[x_hi_11 y_hi_11 zc1_11];
cnr2_11=[x_hi_11 y_lo_11 zc2_11];
cnr3_11=[x_lo_11 y_lo_11 zc3_11];
cnr4_11=[x_lo_11 y_hi_11 zc4_11];
p_cnrs_11=[cnr1_11; cnr2_11; cnr3_11; cnr4_11; cnr1_11];

%distance from points to plane
%use formula (Ax+By+Cz+D)/sqrt(A^2+B^2+C^2)
%the denominator = 1, just evaluate the numerator
plane2point_11=[x_11 y_11 z_11 ones(size(x_11,1),1)]*[A_11 B_11 C_11 D_11]'
m_11=1;
for n_11=1 : size(plane2point_11,1);
if abs(plane2point_11(n_11))<128
cull1_11(m_11,=[x_11(n_11) y_11(n_11) z_11(n_11)];
m_11=m_11+1;
end
end

fid_11=fopen('C:\MATLAB701\work\Point_Clouds11','w')
fprintf(fid_11,'%6.3f %6.3f %6.3f\n',cull1_11');
fclose(fid_11)

%****************FIRST SCAN,SECOND RUN****************
%load the spreadsheet
points_12=load('C:\MATLAB701\work\Point_Clouds11');
x_12=points_12,1);
y_12=points_12,2);
z_12=points_12,3);
%find centroid
N_12=size(points_12,1);
X_12=sum(x_12)/N_12;
Y_12=sum(y_12)/N_12;
Z_12=sum(z_12)/N_12;
centroid_12=[X_12 Y_12 Z_12];
%form basis for matrix elements
h1_12=x_12-X_12; h2_12=y_12-Y_12; h3_12=z_12-Z_12;

%form the elements
k11_12=sum(h1_12.*h1_12); k12_12=sum(h1_12.*h2_12); k13_12=sum(h1_12.*h3_12);
k21_12=sum(h1_12.*h2_12); k22_12=sum(h2_12.*h2_12); k23_12=sum(h2_12.*h3_12);
k31_12=sum(h1_12.*h3_12); k32_12=sum(h2_12.*h3_12); k33_12=sum(h3_12.*h3_12);

%form the matrix
M_12=[k11_12 k12_12 k13_12; k21_12 k22_12 k23_12; k31_12 k32_12 k33_12];

%find the direction vector and value
[vec_12,val_12]=eig(M_12);

%the first column of vec is the unit normal vector and corresponds to the
%minimum eigen value
un_vec_12=vec_12,1);
A_12=un_vec_12(1)
B_12=un_vec_12(2)
C_12=un_vec_12(3)

%calculate D using un_vec and centroid
D_12=-A_12.*centroid_12(1)-B_12.*centroid_12(2)-C_12.*centroid_12(3);

%formula for plane is A*x+B*y+C*z+D=0
%find the z's for the 4 corners of the plane
x_hi_12=max(x_12); x_lo_12=min(x_12);
y_hi_12=max(y_12); y_lo_12=min(y_12);
%z=(-D-A*x-B*y)/C
zc1_12=(-D_12-A_12*x_hi_12-B_12*y_hi_12)/C_12;
zc2_12=(-D_12-A_12*x_hi_12-B_12*y_lo_12)/C_12;
zc3_12=(-D_12-A_12*x_lo_12-B_12*y_lo_12)/C_12;
zc4_12=(-D_12-A_12*x_lo_12-B_12*y_hi_12)/C_12;

%the corner points are
cnr1_12=[x_hi_12 y_hi_12 zc1_12];
cnr2_12=[x_hi_12 y_lo_12 zc2_12];
cnr3_12=[x_lo_12 y_lo_12 zc3_12];
cnr4_12=[x_lo_12 y_hi_12 zc4_12];
p_cnrs_12=[cnr1_12; cnr2_12; cnr3_12; cnr4_12; cnr1_12];

%distance from points to plane
%use formula (Ax+By+Cz+D)/sqrt(A^2+B^2+C^2)
%the denominator = 1, just evaluate the numerator
plane2point_12=[x_12 y_12 z_12 ones(size(x_12,1),1)]*[A_12 B_12 C_12 D_12]'
m_12=1;
for n_12=1 : size(plane2point_12,1);
if abs(plane2point_12(n_12))<64
cull1_12(m_12,=[x_12(n_12) y_12(n_12) z_12(n_12)];
m_12=m_12+1;
end
end

fid_12=fopen('C:\MATLAB701\work\Point_Clouds12','w')
fprintf(fid_12,'%6.3f %6.3f %6.3f\n',cull1_12');
fclose(fid_12)

 

[Greg Truby]

I'm not even sure this is Excel. I am, however pretty sure it's English. So moving to an English forum.

 

[Joe4]

This was also posted in the Excel forum earlier and I moved it to the same place.
http://www.mrexcel.com/forum/showthread.php?t=332345

Per forum rules, please do NOT start new posts on the same questions. Such threads will be locked.