首页 > 学术百科

学习记录——纵横断面面积计算

学习记录——纵横断⾯⾯积计算

最近在逐步学习MatLab，同时因为⼀些活动，要进⾏相关代码的编写，⽬前做了⼀个半成品出来，在这⼉记录⼀下学习的过程吧。

代码中主要以⼤量的循环为主，因为对于MatLab的基础函数以及对于程序编写中参数类型的不熟悉，导致了很多的冗余运算，今后还要将努⼒学习。

clc

clear

close all

%% 读取所有点信息并展⽰

Indata=importdata('');

% fid=fopen('','r');

% data=fscanf(fid,'%s%d,%f,%f,%f\n');

Allpoint_data=Indata.data;

Azimuth_data(3:4);

Point_data(5:end);

H0_Location=data(1));

A=isstrprop(H0_Location,'digit');

H0=str2num(H0_Location(A))/100000;

n=length(Allpoint_data); %存储数组总长，即点的总个数

X=Allpoint_data(:,1);

Y=Allpoint_data(:,2);

H=Allpoint_data(:,3);

plot(X,Y,'*y')

%% 提取K0/K1/K2的位置

SearchTool=strfind(Point_name,'K');

index=1;

for i=1:length(Allpoint_data)

if isempty(SearchTool{i}) == 0

key_point_Location(index,1)=i;

有袋目index=index+1;

end

%% 读取关键点的坐标

Key_point=zeros(length(key_point_Location),3);

for i=1:length(key_point_Location)

Key_point(i,:)=Allpoint_data(key_point_Location(i),:);

end

hold on

plot(Key_point(:,1),Key_point(:,2),'-b')

%% 计算纵断⾯

% 1.计算纵断⾯的长度

D=zeros(length(Key_point)-1,1);

for i=1:length(Key_point)-1

地质学刊D(i)=sqrt((Key_point(i+1,1)-Key_point(i,1))^2+(Key_point(i+1,2)-Key_point(i,2))^2);

end

% 2.计算内插点的平⾯坐标

%（1）先计算内插点总个数

%（2）建⽴⽅位⾓矩阵Azimuth_Key_point，存储所有⽅位⾓信息

%（3）根据循环计算内插点的平⾯坐标

deltC=10;

%内插点个数

C_num=ceil(sum(D)/10);

%⽅位⾓矩阵

Azimuth_Key_point=zeros(length(Key_point)-1,1);

for i=1:length(Key_point)-1

%该处⽅位⾓计算是按横坐标为X，纵坐标为Y计算

deltX=Key_point(i+1,1)-Key_point(i,1);

deltY=Key_point(i+1,2)-Key_point(i,2);

if deltX>0 && deltY>0

Azimuth_Key_point(i)=atan(deltX/deltY);

else if (deltX>0 && deltY<0) || (deltX<0 && deltY<0)

Azimuth_Key_point(i)=atan(deltX/deltY)+pi;

else

Azimuth_Key_point(i)=atan(deltX/deltY)+2*pi;

end

Azimuth_Key_point=rad2deg(Azimuth_Key_point);

%计算内插点坐标

Index1=1;

for i=1:C_num

L=i*deltC;

%判别内插点所在位置

if L<D(Index1)

Interpolation_point(i,1)=Key_point(Index1,1)+L*sind(Azimuth_Key_point(Index1));

Interpolation_point(i,2)=Key_point(Index1,2)+L*cosd(Azimuth_Key_point(Index1));

else if L>D(Index1)

%每当进⼊新的Kj-Kj+1时，判别内插点是否在在关键点连线之外

if length(D)>=Index1+1

if L<sum(D(1:Index1+1))

D0=sqrt((Key_point(Index1+1,1)-Key_point(1,1))^2+(Key_point(Index1+1,2)-Key_point(1,2))^2); Interpolation_point(i,1)=Key_point(Index1+1,1)+(L-D0)*sind(Azimuth_Key_point(Index1+1));

Interpolation_point(i,2)=Key_point(Index1+1,2)+(L-D0)*cosd(Azimuth_Key_point(Index1+1));

else

Index1=Index1+1;

end

else if length(D)<Index1+1

break

end

能源污染

end

hold on

plot(Interpolation_point(:,1),Interpolation_point(:,2),'*k');

% 3.计算内插点⾼程

least_num=5;

for j=1:length(Interpolation_point)

for i=1:n

Dip(i,1)=sqrt((Interpolation_point(j,1)-X(i))^2+(Interpolation_point(j,2)-Y(i))^2);

end

[B,I]=sort(Dip);

I1=I(1:5);

Interpolation_point_height(j,1)=sum(H(I1(:))./B(1:5))/sum(1./B(1:5));

end

%4.计算纵断⾯⾯积

%（1）先将所有关键点放⼊内插点矩阵中的对应位置

%（2）按⾯积公式计算

Dc=D./deltC;

Dc=floor(Dc);

Vertical_section_point=zeros(length(Key_point)+length(Interpolation_point),3);

Vertical_section_point(1,1)=Key_point(1,1);

Vertical_section_point(1,2)=Key_point(1,2);

Vertical_section_point(1,3)=Key_point(1,3);

Index2=1;

Index3=1;

for i=2:length(Vertical_section_point)

if i~=sum(Dc(1:Index2))+1+Index2

Vertical_section_point(i,1)=Interpolation_point(Index3,1);

Vertical_section_point(i,2)=Interpolation_point(Index3,2);

Vertical_section_point(i,3)=Interpolation_point_height(Index3,1);

Index3=Index3+1;

else

if Index2<=length(Dc)

Index2=Index2+1;

Vertical_section_point(i,1)=Key_point(Index2,1);

Vertical_section_point(i,2)=Key_point(Index2,2);

Vertical_section_point(i,3)=Key_point(Index2,3);

else

break

end

end罗尔中值定理

% 每个梯形计算⾯积

S=zeros(length(Vertical_section_point)-1,1);

for i=1:length(Vertical_section_point)-1

deltL=sqrt((Vertical_section_point(i,1)-Vertical_section_point(i+1,1))^2+(Vertical_section_point(i,2)-Vertical_section_point(i+1,2))^2); S(i,1)=deltL*(Vertical_section_point(i,3)+Vertical_section_point(i+1,3)-2*H0)/2;

end

%计算总⾯积

OutS=sum(S);

%% 计算横断⾯

% 1.计算横断⾯中⼼点

Central_point=zeros(length(Key_point)-1,2);

for i=1:length(Key_point)-1

Central_point(i,1)=(Key_point(i,1)+Key_point(i+1,1))/2;

Central_point(i,2)=(Key_point(i,2)+Key_point(i+1,2))/2;

end

hold on

plot(Central_point(:,1),Central_point(:,2),'*r');

%求得每根过中点垂直于关键点连线的直线斜率

Mid_K=zeros(length(Key_point)-1,1);

B=zeros(length(Key_point)-1,1);

for i=1:length(Key_point)-1

Mid_K(i)=-1*(Key_point(i+1,1)-Key_point(i,1))/(Key_point(i+1,2)-Key_point(i,2));

B(i)=Central_point(i,2)-Central_point(i,1)*Mid_K(i);

end

%计算内插点坐标

deltCS=5;

Interpolation_point_CS=zeros(ceil(25*2/deltCS*length(Central_point)),3);

%(1)根据斜率计算⽅位⾓

Azimuth_Central_point=zeros(length(Mid_K),1);

for i=1:length(Central_point)

Azimuth_Central_point(i)=atan(Mid_K(i));

end

% Azimuth_Central_point=rad2deg( Azimuth_Central_point);

%（2）计算两边的内插点坐标

for j=1:length(Azimuth_Central_point)

K_par=(2*j-2)*25/deltCS;

K_par1=(2*j-1)*25/deltCS;

Index4=1;

for i=1:deltCS

%判断⽅位⾓所在象限

if Azimuth_Central_point(j)>0 && Azimuth_Central_point(j)<pi/2

Interpolation_point_CS(K_par+Index4,1)=Central_point(j,1)+i*deltCS*cos(Azimuth_Central_point(j));

Interpolation_point_CS(K_par+Index4,2)=Central_point(j,2)+i*deltCS*sin(Azimuth_Central_point(j));

Interpolation_point_CS(K_par1+Index4,1)=Central_point(j,1)-i*deltCS*cos(Azimuth_Central_point(j));

Interpolation_point_CS(K_par1+Index4,2)=Central_point(j,2)-i*deltCS*sin(Azimuth_Central_point(j));

else if Azimuth_Central_point(j)>pi/2 && Azimuth_Central_point(j)<pi

Interpolation_point_CS(K_par+Index4,1)=Central_point(j,1)+i*deltCS*cos(Azimuth_Central_point(j));

Interpolation_point_CS(K_par+Index4,2)=Central_point(j,2)-i*deltCS*sin(Azimuth_Central_point(j));

Interpolation_point_CS(K_par1+Index4,1)=Central_point(j,1)-i*deltCS*cos(Azimuth_Central_point(j));

Interpolation_point_CS(K_par1+Index4,2)=Central_point(j,2)+i*deltCS*sin(Azimuth_Central_point(j));索尼xperia mt27i

else if Azimuth_Central_point(j)>pi && Azimuth_Central_point(j)<3*pi/2

Interpolation_point_CS(K_par+Index4,1)=Central_point(j,1)-i*deltCS*cos(Azimuth_Central_point(j));

Interpolation_point_CS(K_par+Index4,2)=Central_point(j,2)-i*deltCS*sin(Azimuth_Central_point(j));

Interpolation_point_CS(K_par1+Index4,1)=Central_point(j,1)+i*deltCS*cos(Azimuth_Central_point(j)); Interpolation_point_CS(K_par1+Index4,2)=Central_point(j,2)+i*deltCS*sin(Azimuth_Central_point(j));

else if Azimuth_Central_point(j)>pi/2 && Azimuth_Central_point(j)<pi

Interpolation_point_CS(K_par+Index4,1)=Central_point(j,1)-i*deltCS*cos(Azimuth_Central_point(j)); Interpolation_point_CS(K_par+Index4,2)=Central_point(j,2)+i*deltCS*sin(Azimuth_Central_point(j)); Interpolation_point_CS(K_par1+Index4,1)=Central_point(j,1)+i*deltCS*cos(Azimuth_Central_point(j)); Interpolation_point_CS(K_par1+Index4,2)=Central_point(j,2)-i*deltCS*sin(Azimuth_Central_point(j));

end

Index4=Index4+1;

end

hold on

plot(Interpolation_point_CS(1:10,1),Interpolation_point_CS(1:10,2),'-b')

plot(Interpolation_point_CS(11:20,1),Interpolation_point_CS(11:20,2),'-b')

% plot(Interpolation_point_CS([1,6],1),Interpolation_point_CS([1,6],2),'*k')

% 将所有横断⾯上的点依次放⼊矩阵Cross_section_point

Cross_section_point=zeros(length(Interpolation_point_CS)+length(Central_point),3);

Index5=1;

Index6=1;

for i=1:length(Cross_section_point)

K_par2=25/deltCS+1+(Index6-1)*(25/deltCS*2+1);

if i~= K_par2

Cross_section_point(i,1)=Interpolation_point_CS(Index5,1);

Cross_section_point(i,2)=Interpolation_point_CS(Index5,2);

Cross_section_point(i,3)=Interpolation_point_CS(Index5,3);

Index5=Index5+1;

else

Cross_section_point(i,1)=Central_point(Index6,1);

Cross_section_point(i,2)=Central_point(Index6,2);

Cross_section_point(i,3)=0;

Index6=Index6+1;

end

hold on

plot(Cross_section_point(:,1),Cross_section_point(:,2),'*k')

%计算内插点⾼程

for j=1:length(Cross_section_point)

for i=1:n

Dip_CS(i,1)=sqrt((Cross_section_point(j,1)-X(i))^2+(Cross_section_point(j,2)-Y(i))^2);

end

[B_CS,I_CS]=sort(Dip_CS);

I1_CS=I_CS(1:5);

Cross_section_point(j,3)=sum(H(I1_CS(:))./B_CS(1:5))/sum(1./B_CS(1:5));

end

%计算横断⾯⾯积

S_CS=zeros((length(Cross_section_point)-length(Central_point))/2,length(Central_point));

Index8=1;

for i=1:length(Central_point)

Index7=1;

for j=1:(length(Cross_section_point)-length(Central_point))/2

K_par3=(i-1)*(25/deltCS*2+1);

deltD=sqrt((Cross_section_point(j+ K_par3,1)-Cross_section_point(j+ K_par3+1,1))^2+(Cross_section_point(j+ K_par3,2)-Cross_section_point(j+ K_par3+1, S_CS(Index7,Index8)=(Cross_section_point(j+ K_par3,3)+Cross_section_point(j+ K_par3+1,3)-2*H0)/2*deltD;寇世远

Index7=Index7+1;

end

Index8=Index8+1;

end

本文发布于:2024-09-23 21:26:58，感谢您对本站的认可！

本文链接：https://www.17tex.com/xueshu/207816.html

上一篇：故可以在材料的横断面。上进行研究...

下一篇：隧道横断面测量