gen_mose_bd_2.m

该代码提取图片某一行或者某一列信息

%生成由MOSE所得到的表面数据
clc
clear all

fid = fopen('D:\KASKADE\KeVin\PS_MG\Multigrid_PS_Classics\Point_boundary_totalindex.dat','r');
P_b_index = fscanf(fid,'%f',[1,inf]);
fclose(fid);

fid = fopen('D:\KASKADE\KeVin\PS_MG\Multigrid_PS_Classics\Point_boundary_Dense_matlab.dat','r');
phantom_point_num = fscanf(fid,'%d',[1]);
phantom_element_num = fscanf(fid,'%d',[1]);
% "X", "Y", "Z", "Density"
Dense_temp = fscanf(fid,'%f',[4,inf]);
Dense_temp=Dense_temp.';
B_x=Dense_temp(1:phantom_point_num,1);
Dense_BP=Dense_temp(1:phantom_point_num,4);
fclose(fid);

%get the surface coordinates and density 
P_Surface=find(P_b_index~=0);
Dense_BP_Surface=Dense_temp(P_Surface,:);
%over

sum_FEM=0;
pp=0;
for ii=1:phantom_point_num
    sum_FEM=sum_FEM+Dense_BP(ii);
    if(abs(B_x(ii)-13)<=1 & (abs(Dense_temp(ii,2)^2+Dense_temp(ii,3)^2-100)<0.1))
        pp=pp+1;
    end
end
pp
sum_FEM;
max_FEM_total=max(Dense_BP);

Dense_x15_FEM=zeros(pp,1);
point_x15_FEM=zeros(pp,3);
jj=1;
for ii=1:phantom_point_num
    if(abs(B_x(ii)-13)<=1 & (abs(Dense_temp(ii,2)^2+Dense_temp(ii,3)^2-100)<0.1))
        Dense_x15_FEM(jj)=Dense_BP(ii);
        point_x15_FEM(jj,:)=Dense_temp(ii,1:3);
        jj=jj+1;
    end
end
point_x15_FEM
% Dense_x15_FEM_sort=zeros(size(Dense_x15_FEM),1);
% for ii=1:size(Dense_x15_FEM)
%     if(Dense_temp(ii,2)>=0 & Dense_temp(ii,3)>=0)
%         Dense_x15_FEM(jj)=Dense_BP(ii);
%         point_x15_FEM(jj,:)=Dense_temp(ii,1:3);
%         jj=jj+1;
%     end
% end
PD_x15_FEM=zeros(pp,4);
PD_x15_FEM(:,1:3)=point_x15_FEM;
PD_x15_FEM(:,4)=Dense_x15_FEM;
PD_x15_FEM;
max_FEM_x15=max(Dense_x15_FEM)
%0--180
index_x15_FEM_positive=find(PD_x15_FEM(:,3)>=0);
PD_x15_FEM_positive=PD_x15_FEM(index_x15_FEM_positive,1:4);
[PD_x15_FEM_positive_y_sort,I]=sort(PD_x15_FEM_positive(:,2));
D_x15_FEM_positive_sort=zeros(size(I));
theta_x15_FEM_positive_sort=zeros(size(I));
for ii=1:size(I)
    D_x15_FEM_positive_sort(ii)=PD_x15_FEM_positive(I(ii),4);
    temp_theta_FEM=acos(PD_x15_FEM_positive(I(ii),2)/10.0)*180/pi;
    if(temp_theta_FEM>=45)
        theta_x15_FEM_positive_sort(ii)=temp_theta_FEM-45;
    else
        theta_x15_FEM_positive_sort(ii)=temp_theta_FEM+315;
    end
end
D_x15_FEM_positive_sort;
size(PD_x15_FEM_positive_y_sort)
PD_x15_FEM_positive_y_sort;
%180--360
index_x15_FEM_negative=find(PD_x15_FEM(:,3)<0);
PD_x15_FEM_negative=PD_x15_FEM(index_x15_FEM_negative,1:4);
[PD_x15_FEM_negative_y_sort,K]=sort(PD_x15_FEM_negative(:,2),'descend');
D_x15_FEM_negative_sort=zeros(size(K));
theta_x15_FEM_negative_sort=zeros(size(K));
for ii=1:size(K)
    D_x15_FEM_negative_sort(ii)=PD_x15_FEM_negative(K(ii),4);
    tenp_theta_FEM2=acos(PD_x15_FEM_negative(K(ii),2)/10.0)*180/pi;
    theta_x15_FEM_negative_sort(ii)=135+tenp_theta_FEM2;
end

theta_x15_FEM_negative_sort2=sort(theta_x15_FEM_negative_sort,'ascend');
D_x15_FEM_negative_sort;
size(PD_x15_FEM_negative_y_sort)
PD_x15_FEM_negative_y_sort;

theta_x15_FEM_sort=[theta_x15_FEM_positive_sort;theta_x15_FEM_negative_sort2]
D_x15_FEM_sort=[D_x15_FEM_positive_sort; D_x15_FEM_negative_sort];

y_FEM=[PD_x15_FEM_positive_y_sort;PD_x15_FEM_negative_y_sort];
thetayD_FEM=[theta_x15_FEM_sort,D_x15_FEM_sort]
[thetay_FEM_sort, order_FEM]=sort(thetayD_FEM(:,1), 'ascend');
thetayD_FEM_new=zeros(size(thetayD_FEM(:,1)),2);
for hh=1:size(thetayD_FEM(:,1))
    thetayD_FEM_new(hh,:)=thetayD_FEM(order_FEM(hh),:);
end

size(PD_x15_FEM(:,3));
hold on, plot(thetayD_FEM_new(:,1),thetayD_FEM_new(:,2), 'b -','LineWidth',2)
FEM_curve=zeros(size(PD_x15_FEM(:,3)),2);
FEM_curve=[thetayD_FEM_new(:,1),thetayD_FEM_new(:,2)]

%------------read SG data-------------
% fid_SG = fopen('D:\KASKADE\KeVin\PS_MG\Multigrid_PS\AMGM\2S\Point_boundary_totalindex_2S_SG.dat','r');
% P_b_index_SG = fscanf(fid_SG,'%f',[1,inf]);
% fclose(fid_SG);
% 
% fid = fopen('D:\KASKADE\KeVin\PS_MG\Multigrid_PS\AMGM\2S\Point_boundary_Dense_2S_SG_matlab.dat','r');
% phantom_point_num_SG = fscanf(fid,'%d',[1]);
% phantom_element_num_SG = fscanf(fid,'%d',[1]);
% % "X", "Y", "Z", "Density"
% Dense_temp_SG = fscanf(fid,'%f',[4,inf]);
% Dense_temp_SG=Dense_temp_SG.';
% B_x_SG=Dense_temp_SG(1:phantom_point_num_SG,1);
% Dense_BP_SG=Dense_temp(1:phantom_point_num_SG,4);
% fclose(fid);
% 
% sum_FEM_SG=0;
% pp_SG=0;
% for ii_SG=1:phantom_point_num_SG
%     sum_FEM_SG=sum_FEM_SG+Dense_BP_SG(ii_SG);
%     if(abs(B_x_SG(ii_SG)-13)<=1 & (abs(Dense_temp_SG(ii_SG,2)^2+Dense_temp_SG(ii_SG,3)^2-100)<0.1))
%         pp_SG=pp_SG+1;
%     end
% end
% pp_SG
% sum_FEM_SG;
% max_FEM_total_SG=max(Dense_BP_SG);
% 
% Dense_x15_FEM_SG=zeros(pp_SG,1);
% point_x15_FEM_SG=zeros(pp_SG,3);
% jj_SG=1;
% for ii_SG=1:phantom_point_num_SG
%     if(abs(B_x_SG(ii_SG)-13)<=1 & (abs(Dense_temp_SG(ii_SG,2)^2+Dense_temp_SG(ii_SG,3)^2-100)<0.1))
%         Dense_x15_FEM_SG(jj_SG)=Dense_BP_SG(ii_SG);
%         point_x15_FEM_SG(jj_SG,:)=Dense_temp_SG(ii_SG,1:3);
%         jj_SG=jj_SG+1;
%     end
% end
% point_x15_FEM_SG
% % Dense_x15_FEM_sort=zeros(size(Dense_x15_FEM),1);
% % for ii=1:size(Dense_x15_FEM)
% %     if(Dense_temp(ii,2)>=0 & Dense_temp(ii,3)>=0)
% %         Dense_x15_FEM(jj)=Dense_BP(ii);
% %         point_x15_FEM(jj,:)=Dense_temp(ii,1:3);
% %         jj=jj+1;
% %     end
% % end
% PD_x15_FEM_SG=zeros(pp_SG,4);
% PD_x15_FEM_SG(:,1:3)=point_x15_FEM_SG;
% PD_x15_FEM_SG(:,4)=Dense_x15_FEM_SG;
% PD_x15_FEM_SG;
% max_FEM_x15_SG=max(Dense_x15_FEM_SG)
% %0--180
% index_x15_FEM_positive_SG=find(PD_x15_FEM_SG(:,3)>=0);
% PD_x15_FEM_positive_SG=PD_x15_FEM_SG(index_x15_FEM_positive_SG,1:4);
% [PD_x15_FEM_positive_y_sort_SG,I_SG]=sort(PD_x15_FEM_positive_SG(:,2));
% D_x15_FEM_positive_sort_SG=zeros(size(I_SG));
% theta_x15_FEM_positive_sort_SG=zeros(size(I_SG));
% for ii=1:size(I_SG)
%     D_x15_FEM_positive_sort_SG(ii)=PD_x15_FEM_positive_SG(I_SG(ii),4);
%     temp_theta_FEM_SG=acos(PD_x15_FEM_positive_SG(I_SG(ii),2)/10.0)*180/pi;
%     if(temp_theta_FEM_SG>=45)
%         theta_x15_FEM_positive_sort_SG(ii)=temp_theta_FEM_SG-45;
%     else
%         theta_x15_FEM_positive_sort_SG(ii)=temp_theta_FEM_SG+315;
%     end
% end
% D_x15_FEM_positive_sort_SG;
% size(PD_x15_FEM_positive_y_sort_SG)
% PD_x15_FEM_positive_y_sort_SG;
% %180--360
% index_x15_FEM_negative_SG=find(PD_x15_FEM_SG(:,3)<0);
% PD_x15_FEM_negative_SG=PD_x15_FEM_SG(index_x15_FEM_negative_SG,1:4);
% [PD_x15_FEM_negative_y_sort_SG,K_SG]=sort(PD_x15_FEM_negative_SG(:,2),'descend');
% D_x15_FEM_negative_sort_SG=zeros(size(K_SG));
% theta_x15_FEM_negative_sort_SG=zeros(size(K_SG));
% for ii=1:size(K_SG)
%     D_x15_FEM_negative_sort_SG(ii)=PD_x15_FEM_negative_SG(K_SG(ii),4);
%     tenp_theta_FEM2_SG=acos(PD_x15_FEM_negative_SG(K_SG(ii),2)/10.0)*180/pi;
%     theta_x15_FEM_negative_sort_SG(ii)=135+tenp_theta_FEM2_SG;
% end
% 
% theta_x15_FEM_negative_sort2_SG=sort(theta_x15_FEM_negative_sort_SG,'ascend');
% D_x15_FEM_negative_sort_SG;
% size(PD_x15_FEM_negative_y_sort_SG)
% PD_x15_FEM_negative_y_sort_SG;
% 
% theta_x15_FEM_sort_SG=[theta_x15_FEM_positive_sort_SG;theta_x15_FEM_negative_sort2_SG]
% D_x15_FEM_sort_SG=[D_x15_FEM_positive_sort_SG; D_x15_FEM_negative_sort_SG];
% 
% y_FEM_SG=[PD_x15_FEM_positive_y_sort_SG;PD_x15_FEM_negative_y_sort_SG];
% thetayD_FEM_SG=[theta_x15_FEM_sort_SG,D_x15_FEM_sort_SG]
% [thetay_FEM_sort_SG, order_FEM_SG]=sort(thetayD_FEM_SG(:,1), 'ascend');
% thetayD_FEM_new_SG=zeros(size(thetayD_FEM_SG(:,1)),2);
% for hh=1:size(thetayD_FEM_SG(:,1))
%     thetayD_FEM_new_SG(hh,:)=thetayD_FEM_SG(order_FEM_SG(hh),:);
% end
% 
% size(PD_x15_FEM_SG(:,3));
% hold on, plot(thetayD_FEM_new_SG(:,1),thetayD_FEM_new_SG(:,2),  'g -','LineWidth',2)
%SG data over

%%%读入MOSE_OFF各个文件
fid = fopen('D:\KASKADE\KeVin\PS_MG\Multigrid_PS\AMGM\1S_OFF\surface.off','r');
file_ext = fscanf(fid,'%s',[1]);
num_data = fscanf(fid,'%d',[3,1]);
num_P=num_data(1)-1;
triangle_num = num_data(2);
dummy = fscanf(fid,'%f',[3,1]);
off_point = fscanf(fid,'%f',[3,num_P]);
off_point=off_point.';
off_triangle0 = fscanf(fid,'%f',[4,triangle_num]);
off_triangle0=off_triangle0.';
off_triangle=off_triangle0(:,2:4);
fclose(fid);


fid = fopen('D:\KASKADE\KeVin\PS_MG\Multigrid_PS\AMGM\1S\1S_Results\3DVertexFluxOutput.dat','r');
MOSE_data = fscanf(fid,'%f',[4,inf]);
MOSE_data = MOSE_data.';
%index_MOSE=MOSE_data(:,1);
points_MOSE=MOSE_data(:,1:3);
dense_MOSE=MOSE_data(:,4);
fclose(fid);

sum_MOSE=0;
qq=0;
for ii=1:length(dense_MOSE)
    sum_MOSE=sum_MOSE+dense_MOSE(ii);
    if(abs(points_MOSE(ii,1)-13)<=0.18)
        qq=qq+1;
    end
end
sum_MOSE;
max_MOSE=max(dense_MOSE);


Dense_x15_MOSE=zeros(qq,1);
point_x15_MOSE=zeros(qq,3);
jj=1;
for ii=1:length(dense_MOSE)
    if(abs(points_MOSE(ii,1)-13)<=0.18)
        Dense_x15_MOSE(jj)=MOSE_data(ii,4);
        point_x15_MOSE(jj,:)=MOSE_data(ii,1:3);
        jj=jj+1;
    end
end

PD_x15_MOSE=zeros(qq,4);
PD_x15_MOSE(:,1:3)=point_x15_MOSE;
PD_x15_MOSE(:,4)=Dense_x15_MOSE;
PD_x15_MOSE;

%0--180_MOSE
index_x15_MOSE_positive=find(PD_x15_MOSE(:,3)>=0);
PD_x15_MOSE_positive=PD_x15_MOSE(index_x15_MOSE_positive,1:4);
[PD_x15_MOSE_positive_z_sort,I_MOSE]=sort(PD_x15_MOSE_positive(:,2));
theta_x15_MOSE_positive_sort=zeros(size(I_MOSE));
D_x15_MOSE_positive_sort=zeros(size(I_MOSE));
for ii=1:size(I_MOSE)
    D_x15_MOSE_positive_sort(ii)=PD_x15_MOSE_positive(I_MOSE(ii),4);
    temp_theta=acos(PD_x15_MOSE_positive(I_MOSE(ii),2)/10)*180/pi-45;
    if(temp_theta>=0)
        theta_x15_MOSE_positive_sort(ii)=temp_theta;
    else
        theta_x15_MOSE_positive_sort(ii)=temp_theta+360;
    end
end
D_x15_MOSE_positive_sort;