-
Notifications
You must be signed in to change notification settings - Fork 3
/
main_topo
107 lines (96 loc) · 4.18 KB
/
main_topo
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
%==========================================================================
%~~~~~~~~~~~~~~~~~~ooooo~~~~ooooo~~~~0~~~0~~~~ooooo~~~~~~~~~~~~~~~~~~~~~~~~
%~~~~~~~~~~~~~~~~~~~~0~~~~~~~~0~~~~~~0o~~0~~~~0~~~~~~~~~~~~~~~~~~~~~~~~~~~~
%~~~~~~~~~~~~~~~~~~~~0~~~~~~~~0~~~~~~0~o~0~~~~0~oo0~~~~~~~~~~~~~~~~~~~~~~~~
%~~~~~~~~~~~~~~~~~~~~0~~~~~~~~0~~~~~~0~~o0~~~~0~~~0~~~~~~~~~~~~~~~~~~~~~~~~
%-----------------ooo0------oo0oo----0---0----0oooo------------------------
%==========================================================================
% ^---^ Meow
% (o . o)
% ||=||
%________________________UU-UU_____________________________________________
%__________________________Y_______________________________________________
% |
% J
%==========================================================================
%data input
%==========================================================================
%read models
clear;
%dname='E:\jing\percolation';
dname='\\INFPWFS803.ad.unsw.edu.au\Staff053$\z3379317\Documents\Documents\MATLAB\matlab_script_centos\percolation';
foldername='percolation';%'realisation1_500'; %for realisations
file='derivation_original_inflation_2.nc';%'realisation_mineralised_1.nc'; %for realisations
%file='substack_500_cleats_clear_mineralised.nc';
%file='substack_500.nc';
fname=fullfile(dname,foldername,file); %for realisations
%fname=fullfile(dname,file);
variable='segmented'; %for realisations
%variable='1';
mat=nc2mat(fname,variable);
%%
%this is for original image
mat(mat==0)=5;
mat(mat==2)=0;
mat=~logical(mat);
mat=permute(mat,[1 3 2]);
mat=mat(:,:,1:500);
%%
%this is for cleared original image
mat=logical(mat);
mat=permute(mat,[1 3 2]);
mat=mat(:,:,1:500);
%%
%==========================================================================
%local porosity and percolation probability
%==========================================================================
sub_length=10:20:170; %size of subsample
%sub_length=190;
n_sub=500; %number of subsamples
nbin=20;
%%local porosity and percolation probability
[ probability_poro, probability_perc, bin, pav ] = probability_poro_perc( mat,sub_length,n_sub,nbin );
%%
%==========================================================================
%percolation threshold
%==========================================================================
%mat_cleat=permute(mat,[1 3 2]); %original
matrix=dfn; %realisation
matrix=double(matrix);
%mat_cleat_3=imresize3(mat_cleat,3,'nearest');
for i=1:size(matrix,3)
matrix_3(:,:,i)=imresize(matrix(:,:,i),2,'nearest');
end
step=10;
%[ mat_deflation ] = deflation(mat_cleat_3,step ); %target phase is matris, so mat should be ~
[ mat_deflation ] = deflation(matrix_3,step ); %target phase is matris, so mat should be ~
%morphy_NC( mat_deflation,0.0165,'test.nc');
perc=[];
mat_deflation_binary=mat_deflation;
%mat_deflation_binary=permute(mat_deflation_binary,[1 3 2]); %in the direction of horizental
%mat_deflation_binary=permute(mat_deflation_binary,[3 2 1]); %in the direction of vertical
mat_deflation_binary(:,:,1:step)= [];
mat_deflation_binary(:,:,end-step:end)=[];
mat_deflation_edge=mat_deflation_binary;
poro=[];
for i=1:step
mat_deflation_binary(mat_deflation_edge==(i+1))=0;
mat_label=bwlabeln(logical(mat_deflation_binary));
cluster=intersect(unique(mat_label(:,:,1)),unique(mat_label(:,:,end)),'rows');
perc(i)=~isempty(find(cluster>0, 1));
poro(i)=sum(logical(mat_deflation_binary(:)))/numel(mat_deflation_binary);
end
%%
perc_filter=perc(find(perc));
i=length(perc_filter);
clear mat_deflation_binary
mat_critical=mat_deflation_edge;
mat_critical(mat_deflation_edge<=(i+1)& mat_deflation_edge>1)=0; %critical model
percolation_threshold=poro(i);%sum(logical(mat_critical(:)))/numel(mat_critical);
%%
%==========================================================================
%fractal dimension
%==========================================================================
for i=1:length(sub_length)
density(i) = mass_density1( mat, sub_length(i),n_sub );
end