Commit 47b0ac0b authored by Lukas Platz's avatar Lukas Platz
Browse files

demo prints zero mode values of mock and reconstruction

parent a4217c4e
...@@ -76,8 +76,7 @@ if __name__ == '__main__': ...@@ -76,8 +76,7 @@ if __name__ == '__main__':
power_space_p = ift.PowerSpace(harmonic_space_p) power_space_p = ift.PowerSpace(harmonic_space_p)
power_space_e = ift.PowerSpace(harmonic_space_e) power_space_e = ift.PowerSpace(harmonic_space_e)
amp_p = ift.SLAmplitude( amp_p_dict = {
**{
'target': power_space_p, 'target': power_space_p,
'n_pix': 64, # 64 spectral bins 'n_pix': 64, # 64 spectral bins
'keys': ['tau_p', 'phi_p', 'zm_p'], 'keys': ['tau_p', 'phi_p', 'zm_p'],
...@@ -95,10 +94,9 @@ if __name__ == '__main__': ...@@ -95,10 +94,9 @@ if __name__ == '__main__':
# zero-mode # zero-mode
'za': 5., # IG alpha parameter 'za': 5., # IG alpha parameter
'zq': 0.001 # IG q parameter 'zq': 0.001 # IG q parameter
}) }
amp_e = ift.SLAmplitude( amp_e_dict = {
**{
'target': power_space_e, 'target': power_space_e,
'n_pix': 64, # 64 spectral bins 'n_pix': 64, # 64 spectral bins
'keys': ['tau_e', 'phi_e', 'zm_e'], 'keys': ['tau_e', 'phi_e', 'zm_e'],
...@@ -116,7 +114,10 @@ if __name__ == '__main__': ...@@ -116,7 +114,10 @@ if __name__ == '__main__':
# zero-mode # zero-mode
'za': 6., # IG alpha parameter 'za': 6., # IG alpha parameter
'zq': 0.002 # IG q parameter 'zq': 0.002 # IG q parameter
}) }
amp_p = ift.SLAmplitude(**amp_p_dict)
amp_e = ift.SLAmplitude(**amp_e_dict)
# Define sky operator # Define sky operator
# za=4, zq=5 leads to an amplitude median of one # za=4, zq=5 leads to an amplitude median of one
...@@ -125,10 +126,9 @@ if __name__ == '__main__': ...@@ -125,10 +126,9 @@ if __name__ == '__main__':
elif mode == 0: elif mode == 0:
amplitudes = [None, amp_e] amplitudes = [None, amp_e]
log_diffuse = ift.MultiCorrelatedField(sky_target, diff_zm_params = {'gs_a': 11., 'gs_q': 50.}
amplitudes, log_diffuse = ift.MultiCorrelatedField(sky_target, amplitudes,
gs_a=11., **diff_zm_params)
gs_q=50.)
sky = log_diffuse.exp() sky = log_diffuse.exp()
M = ift.DiagonalOperator(exposure) M = ift.DiagonalOperator(exposure)
...@@ -142,7 +142,6 @@ if __name__ == '__main__': ...@@ -142,7 +142,6 @@ if __name__ == '__main__':
data = lamb(mock_position) data = lamb(mock_position)
signal = sky(mock_position) signal = sky(mock_position)
# DEBUG # DEBUG
mock_log_diffuse_vals = log_diffuse(mock_position).to_global_data() mock_log_diffuse_vals = log_diffuse(mock_position).to_global_data()
print("std(mock_log_diffuse):", mock_log_diffuse_vals.std()) print("std(mock_log_diffuse):", mock_log_diffuse_vals.std())
...@@ -165,7 +164,7 @@ if __name__ == '__main__': ...@@ -165,7 +164,7 @@ if __name__ == '__main__':
iteration = 0 iteration = 0
# Minimize # Minimize
for i in range(3): for i in range(10):
KL = ift.MetricGaussianKL(pos, H, 10) KL = ift.MetricGaussianKL(pos, H, 10)
KL, convergence = minimizer(KL) KL, convergence = minimizer(KL)
pos = KL.position pos = KL.position
...@@ -182,7 +181,9 @@ if __name__ == '__main__': ...@@ -182,7 +181,9 @@ if __name__ == '__main__':
plot.add(cd(GR.adjoint(data)), title='Data') plot.add(cd(GR.adjoint(data)), title='Data')
plot.add(cd(reconst), title='Posterior Mean') plot.add(cd(reconst), title='Posterior Mean')
plot.add(cd(reconst - signal), title='Error') plot.add(cd(reconst - signal), title='Error')
plot.output(xsize=12, ysize=10, name=filename.format(f'reconstr-{iteration:02d}')) plot.output(xsize=12,
ysize=10,
name=filename.format(f'{iteration:02d}-reconstr'))
p_p = ift.PowerDistributor(harmonic_space_p, power_space_p) @ amp_p p_p = ift.PowerDistributor(harmonic_space_p, power_space_p) @ amp_p
p_e = ift.PowerDistributor(harmonic_space_e, power_space_e) @ amp_e p_e = ift.PowerDistributor(harmonic_space_e, power_space_e) @ amp_e
...@@ -197,13 +198,26 @@ if __name__ == '__main__': ...@@ -197,13 +198,26 @@ if __name__ == '__main__':
for sample in KL.samples: for sample in KL.samples:
sc.add(ift.power_analyze(p.force(pos + sample))) sc.add(ift.power_analyze(p.force(pos + sample)))
plot = ift.Plot() plot = ift.Plot()
plot.add([ plot.add(
sc.mean, [sc.mean, ift.power_analyze(p.force(mock_position))],
ift.power_analyze(p.force(mock_position))
],
linewidth=[1., 3.], linewidth=[1., 3.],
label=['Posterior Mean', 'Ground Truth'], label=['Posterior Mean', 'Ground Truth'],
title=f"Posterior {n} spectrum") title=f"Posterior {n} spectrum")
plot.output(xsize=12, ysize=10, name=filename.format(f'power_{n}-{iteration:02d}')) plot.output(xsize=12,
ysize=10,
name=filename.format(f'{iteration:02d}-power_{n}'))
print("Saved results as '{}'.".format(filename.format('*'))) print("Saved results as '{}'.".format(filename.format('*')))
# Calculate InverseGamma prior misfits
keys = ['zm_e', 'zm_p', 'xi_pspec_scaling']
hyperparams = [(amp_e_dict['za'], amp_e_dict['zq']),
(amp_p_dict['za'], amp_p_dict['zq']),
(diff_zm_params['gs_a'], diff_zm_params['gs_q'])]
print("InverseGamma missfits:")
for key, hpar in zip(keys, hyperparams):
op = ift.InverseGammaOperator(sky.domain[key], *hpar)
print(f"{key}:")
print(f"\tmock: {float(op(mock_position[key]).local_data):+2.3f}")
print(f"\t rec: {float(op(pos[key]).local_data):+2.3f}\n")
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment