Skip to content
Snippets Groups Projects
Commit 8c9884c9 authored by Philipp Arras's avatar Philipp Arras
Browse files

nifty5 -> nifty7 (4/n)

parent b3d6cc08
No related branches found
No related tags found
1 merge request!2Draft: Nifty5 to nifty7
Pipeline #107683 failed
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program. If not, see <http://www.gnu.org/licenses/>.
#
# Copyright(C) 2013-2021 Max-Planck-Society
#
# NIFTy is being developed at the Max-Planck-Institut fuer Astrophysik.
import numpy as np import numpy as np
import nifty7 as ift import nifty7 as ift
......
...@@ -32,7 +32,7 @@ data_space = R.target ...@@ -32,7 +32,7 @@ data_space = R.target
data = ift.makeField(data_space, data) data = ift.makeField(data_space, data)
ground_truth = ift.makeField(position_space, ground_truth) ground_truth = ift.makeField(position_space, ground_truth)
plot_WF('data', ground_truth, data) plot_WF('1_data', ground_truth, data)
N = ift.ScalingOperator(data_space, 0.1) N = ift.ScalingOperator(data_space, 0.1)
...@@ -45,12 +45,12 @@ S = HT @ S_h @ HT.adjoint ...@@ -45,12 +45,12 @@ S = HT @ S_h @ HT.adjoint
D_inv = S.inverse + R.adjoint @ N.inverse @ R D_inv = S.inverse + R.adjoint @ N.inverse @ R
j = (R.adjoint @ N.inverse)(data) j = (R.adjoint @ N.inverse)(data)
IC = ift.GradientNormController(iteration_limit=1, tol_abs_gradnorm=1e-7) IC = ift.GradientNormController(iteration_limit=100, tol_abs_gradnorm=1e-7)
D = ift.InversionEnabler(D_inv.inverse, IC, approximation=S) D = ift.InversionEnabler(D_inv.inverse, IC, approximation=S)
m = D(j) m = D(j)
plot_WF('result', ground_truth, data, m) plot_WF('1_result', ground_truth, data, m)
N = ift.SamplingDtypeSetter(N, np.float64) N = ift.SamplingDtypeSetter(N, np.float64)
S_h = ift.SamplingDtypeSetter(S_h, np.float64) S_h = ift.SamplingDtypeSetter(S_h, np.float64)
...@@ -59,4 +59,4 @@ Dinv = ift.WienerFilterCurvature(R, N, S, IC, IC) ...@@ -59,4 +59,4 @@ Dinv = ift.WienerFilterCurvature(R, N, S, IC, IC)
N_samples = 10 N_samples = 10
samples = [Dinv.draw_sample(from_inverse=True) + m for i in range(N_samples)] samples = [Dinv.draw_sample(from_inverse=True) + m for i in range(N_samples)]
plot_WF('result', ground_truth, data, m=m, samples=samples) plot_WF('1_result_with_uncertainty', ground_truth, data, m=m, samples=samples)
...@@ -15,25 +15,27 @@ ...@@ -15,25 +15,27 @@
# #
# NIFTy is being developed at the Max-Planck-Institut fuer Astrophysik. # NIFTy is being developed at the Max-Planck-Institut fuer Astrophysik.
import numpy as np
import nifty7 as ift import nifty7 as ift
from helpers import (checkerboard_response, generate_gaussian_data,
plot_prior_samples_2d, plot_reconstruction_2d)
position_space = ift.RGSpace(2*(256,)) position_space = ift.RGSpace(2*(256,))
harmonic_space = position_space.get_default_codomain()
HT = ift.HarmonicTransformOperator(harmonic_space, target=position_space) args = {
power_space = ift.PowerSpace(harmonic_space) 'offset_mean': 0,
'offset_std': (1e-3, 1e-6),
A = ift.SLAmplitude(
**{ # Amplitude of field fluctuations
'target': power_space, 'fluctuations': (1., 0.8), # 1.0, 1e-2
'n_pix': 64, # 64 spectral bins
# Smoothness of spectrum # Exponent of power law power spectrum component
'a': 10, # relatively high variance of spectral curvature 'loglogavgslope': (-3., 1), # -6.0, 1
'k0': .2, # quefrency mode below which cepstrum flattens
# Power-law part of spectrum # Amplitude of integrated Wiener process power spectrum component
'sm': -4, # preferred power-law slope 'flexibility': (2, 1.), # 1.0, 0.5
'sv': .6, # low variance of power-law slope
'im': -2, # y-intercept mean, in-/decrease for more/less contrast # How ragged the integrated Wiener process component is
'iv': 2. # y-intercept variance 'asperity': (0.5, 0.4) # 0.1, 0.5
}) }
signal = ift.SimpleCorrelatedField(position_space, **args)
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment