MIRI.m
2.57 KB
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
%Concrete class for computing MIRI of model parameters. The class has two
%methods:
% - estimate_parampdfs that estimates the conditional upper and lower pdf
% of a parameter. It returns in output the array all_xbin_max and
% all_xbin_min containing the points of the ksdensity in the upper and
% lower case while all_ks_max and all_ks_min contains the corresponding
% probability density estimate;
% - MIRI_computation takes in input the two estimated pdfs of all
% parameters and computes MIRI for each of them. MIRI are returned in the
% array MIRI_param while CloudCondMax and CloudCondMin contain the mode of
% the parameter vector in the upper and lower case respectively.
classdef MIRI
properties
end
methods
function [all_xbin_max,all_xbin_min,all_ks_max,all_ks_min]=estimate_parampdfs(obj,p, Perturbation, XiMax, XiMin)
all_xbin_max=[];
all_xbin_min=[];
all_ks_max=[];
all_ks_min=[];
pdf_obj=pdfEstimator();
for ip=1:length(p)
%step=p(ip)*(Ncloud^(-1./3));
step=(max(p(ip).*Perturbation(:,ip)) - min(p(ip).*Perturbation(:,ip)))/size(Perturbation(:,ip),1);
BinEdgesMax=min(p(ip).*Perturbation(:,ip)):step:max(p(ip).*Perturbation(:,ip));
BinEdgesMin=min(p(ip).*Perturbation(:,ip)):step:max(p(ip).*Perturbation(:,ip));
[ks_max, xbin_max]=pdf_obj.evaluate_pdf(p(ip).*XiMax(:,ip)',BinEdgesMax);
[ks_min, xbin_min]=pdf_obj.evaluate_pdf(p(ip).*XiMin(:,ip)',BinEdgesMin);
all_xbin_max=[all_xbin_max;xbin_max];
all_xbin_min=[all_xbin_min;xbin_min];
all_ks_max=[all_ks_max;ks_max];
all_ks_min=[all_ks_min;ks_min];
end
end
function [MIRI_param, CloudCondMin, CloudCondMax]=MIRI_computation(obj, all_ks_max,all_ks_min,all_xbin_max,all_xbin_min)
len_p=size(all_ks_max,1);
MIRI_param=zeros(1,len_p);
CloudCondMin=zeros(1,len_p);
CloudCondMax=zeros(1,len_p);
for ip=1:len_p
ks_max=all_ks_max(ip,:);
ks_min=all_ks_min(ip,:);
xbin_max=all_xbin_max(ip,:);
xbin_min=all_xbin_min(ip,:);
[ks_max_m, ks_max_I]=max(ks_max);
[ks_min_m, ks_min_I]=max(ks_min);
CloudCondMin(1,ip)=xbin_min(ks_min_I);
CloudCondMax(1,ip)=xbin_max(ks_max_I);
s_Xi_ks=sum(abs(ks_max./sum(ks_max)-ks_min./sum(ks_min)));
MIRI_param(1,ip)=s_Xi_ks;
end
end
end
end