File wellmet/simplex.py changed (mode: 100644) (index 92688f5..52521fb) |
... |
... |
from . import sball |
8 |
8 |
from . import f_models |
from . import f_models |
9 |
9 |
from .candynodes import CandyNodes |
from .candynodes import CandyNodes |
10 |
10 |
from scipy import spatial |
from scipy import spatial |
11 |
|
from scipy.spatial.distance import cdist |
|
12 |
11 |
from scipy import stats |
from scipy import stats |
13 |
12 |
|
|
14 |
13 |
import quadpy |
import quadpy |
|
... |
... |
class _Triangulation: |
254 |
253 |
|
|
255 |
254 |
sx.sample_box = sample_box |
sx.sample_box = sample_box |
256 |
255 |
sx.tri_space = tri_space |
sx.tri_space = tri_space |
|
256 |
|
sx.PDF = sample_box.pdf(tri_space) |
|
257 |
|
sx.failsi = sample_box.failsi |
257 |
258 |
sx.issi = issi #č ISSI potřebujem pro tri_estimation |
sx.issi = issi #č ISSI potřebujem pro tri_estimation |
258 |
259 |
|
|
259 |
260 |
if weighting_space is None: |
if weighting_space is None: |
|
... |
... |
class _Triangulation: |
377 |
378 |
|
|
378 |
379 |
def _tri_update(sx, tri_plan): |
def _tri_update(sx, tri_plan): |
379 |
380 |
#č separujeme, abychom vědeli, kolik času žere samotný QHull |
#č separujeme, abychom vědeli, kolik času žere samotný QHull |
380 |
|
sx.tri.add_points(tri_plan[sx.tri.npoints:]) |
|
|
381 |
|
sx.tri.add_points(tri_plan[sx.tri.npoints:]) |
|
382 |
|
#č ale jo, nejen QHull |
|
383 |
|
sx.PDF = sx.sample_box.pdf(sx.tri_space) |
|
384 |
|
sx.failsi = sx.sample_box.failsi |
381 |
385 |
|
|
382 |
386 |
|
|
383 |
387 |
|
|
|
... |
... |
class BetterCubatureIntegration(_Triangulation): |
1072 |
1076 |
sx.simplex_stats.pop(simplex, None) |
sx.simplex_stats.pop(simplex, None) |
1073 |
1077 |
|
|
1074 |
1078 |
if sx.on_delete_simplex is not None: |
if sx.on_delete_simplex is not None: |
1075 |
|
#č zpatky do ndarray... |
|
1076 |
|
sx.on_delete_simplex(indices=np.array(simplex)) |
|
|
1079 |
|
#č Bacha, teď posíláme tuple. |
|
1080 |
|
#č Kód na té straně úplně stejně bude vyhazovat |
|
1081 |
|
#č záznamy ze slovníků |
|
1082 |
|
sx.on_delete_simplex(simplex) |
1077 |
1083 |
|
|
1078 |
1084 |
|
|
1079 |
1085 |
#č vyhodil jsem simplex_id'y |
#č vyhodil jsem simplex_id'y |
|
... |
... |
class BetterCubatureIntegration(_Triangulation): |
1103 |
1109 |
|
|
1104 |
1110 |
def _integrate_mixed_simplex(sx, indices): |
def _integrate_mixed_simplex(sx, indices): |
1105 |
1111 |
|
|
1106 |
|
f = sx.f |
|
1107 |
|
vertices = f[indices] |
|
1108 |
|
vertices_model = getattr(vertices, sx.tri_space) |
|
1109 |
|
PDF = vertices.pdf(sx.tri_space) |
|
1110 |
|
failsi = sx.sample_box.failsi[indices] |
|
|
1112 |
|
vertices_model = sx.tri.points[indices] |
|
1113 |
|
PDF = sx.PDF[indices] |
|
1114 |
|
failsi = sx.failsi[indices] |
1111 |
1115 |
|
|
1112 |
1116 |
# quadpy |
# quadpy |
1113 |
1117 |
def _get_pdf(x): |
def _get_pdf(x): |
1114 |
|
dmv = 1 / cdist(x.T, vertices_model) |
|
|
1118 |
|
dmv = 1 / spatial.distance.cdist(x.T, vertices_model) |
1115 |
1119 |
dmw = dmv * PDF |
dmw = dmv * PDF |
1116 |
1120 |
|
|
1117 |
1121 |
# same as np.sum(dmv[:, failsi], axis=1) / np.sum(dmv, axis=1) |
# same as np.sum(dmv[:, failsi], axis=1) / np.sum(dmv, axis=1) |
|
... |
... |
class BetterCubatureIntegration(_Triangulation): |
1120 |
1124 |
|
|
1121 |
1125 |
# side effect |
# side effect |
1122 |
1126 |
if sx.on_add_simplex is not None: |
if sx.on_add_simplex is not None: |
1123 |
|
nodes = f.new_sample(x.T, sx.tri_space) |
|
|
1127 |
|
nodes = sx.f.new_sample(x.T, sx.tri_space) |
1124 |
1128 |
fx = nodes.pdf(sx.tri_space) |
fx = nodes.pdf(sx.tri_space) |
1125 |
1129 |
# -1=outside, 0=success, 1=failure, 2=mix |
# -1=outside, 0=success, 1=failure, 2=mix |
1126 |
1130 |
args = {'event_id':2, 'indices':indices} |
args = {'event_id':2, 'indices':indices} |
1127 |
1131 |
sx.on_add_simplex(CandyNodes(nodes, args, pfv=pfv, pfw=pfw)) |
sx.on_add_simplex(CandyNodes(nodes, args, pfv=pfv, pfw=pfw)) |
1128 |
1132 |
else: |
else: |
1129 |
|
fx = f.sample_pdf(x.T, sx.tri_space) |
|
|
1133 |
|
fx = sx.f.sample_pdf(x.T, sx.tri_space) |
1130 |
1134 |
|
|
1131 |
1135 |
return fx, pfv, pfw |
return fx, pfv, pfw |
1132 |
1136 |
|
|
|
... |
... |
class BetterCubatureIntegration(_Triangulation): |
1153 |
1157 |
|
|
1154 |
1158 |
def _integrate_red_simplex(sx, indices): |
def _integrate_red_simplex(sx, indices): |
1155 |
1159 |
|
|
1156 |
|
f = sx.f |
|
1157 |
|
vertices_model = getattr(f, sx.tri_space)[indices] |
|
|
1160 |
|
vertices_model = sx.tri.points[indices] |
1158 |
1161 |
|
|
1159 |
1162 |
# quadpy |
# quadpy |
1160 |
1163 |
def _get_pdf(x): |
def _get_pdf(x): |
1161 |
1164 |
|
|
1162 |
1165 |
# side effect |
# side effect |
1163 |
1166 |
if sx.on_add_simplex is not None: |
if sx.on_add_simplex is not None: |
1164 |
|
nodes = f.new_sample(x.T, sx.tri_space) |
|
|
1167 |
|
nodes = sx.f.new_sample(x.T, sx.tri_space) |
1165 |
1168 |
fx = nodes.pdf(sx.tri_space) |
fx = nodes.pdf(sx.tri_space) |
1166 |
1169 |
# -1=outside, 0=success, 1=failure, 2=mix |
# -1=outside, 0=success, 1=failure, 2=mix |
1167 |
1170 |
args = {'event_id':1, 'indices':indices} |
args = {'event_id':1, 'indices':indices} |
1168 |
1171 |
sx.on_add_simplex(CandyNodes(nodes, args)) |
sx.on_add_simplex(CandyNodes(nodes, args)) |
1169 |
1172 |
else: |
else: |
1170 |
|
fx = f.sample_pdf(x.T, sx.tri_space) |
|
|
1173 |
|
fx = sx.f.sample_pdf(x.T, sx.tri_space) |
1171 |
1174 |
|
|
1172 |
1175 |
return fx |
return fx |
1173 |
1176 |
|
|