Commit 8c0c97b7 authored by Martin Reinecke's avatar Martin Reinecke

try to fix volume factors

parent 70b5b7cd
......@@ -9,12 +9,12 @@ from matplotlib import pyplot as plt
np.random.seed(42)
if __name__ == '__main__':
s_space = ift.RGSpace([1024])
FFT = ift.FFTOperator(s_space)
h_space = FFT.target[0]
h_space = s_space.get_default_codomain()
FFT = ift.FFTOperator(h_space)
p_spec = lambda k: (1./(1+k)**2.5)
S = ift.create_power_operator(h_space, power_spectrum=p_spec)
sh = S.draw_sample()
s = FFT.adjoint_times(sh)
s = FFT(sh)
u = ift.Field(s_space, val = -12)
u.val[200] = 1
......
......@@ -18,12 +18,12 @@ def load_data(path):
def build_problem(data, alpha):
s_space = ift.RGSpace(data.shape, distances=len(data.shape) * [1])
h_space = s_space.get_default_codomain()
data = ift.Field(s_space,val=data)
FFT = ift.FFTOperator(s_space)
h_space = FFT.target[0]
FFT = ift.FFTOperator(h_space)
binbounds = ift.PowerSpace.useful_binbounds(h_space, logarithmic = False)
p_space = ift.PowerSpace(h_space, binbounds=binbounds)
initial_spectrum = ift.power_analyze(FFT(ift.log(data)), binbounds=p_space.binbounds)
initial_spectrum = ift.power_analyze(FFT.inverse_times(ift.log(data)), binbounds=p_space.binbounds)
initial_correlation = ift.create_power_operator(h_space, initial_spectrum)
initial_x = ift.Field(s_space, val=-1.)
alpha = ift.Field(s_space, val=alpha)
......@@ -47,7 +47,7 @@ def problem_iteration(energy, iterations=3):
h_space = energy.correlation.domain[0]
FFT = energy.FFT
binbounds = ift.PowerSpace.useful_binbounds(h_space, logarithmic=False)
new_power = ift.power_analyze(FFT(energy.s), binbounds=binbounds)
new_power = ift.power_analyze(FFT.inverse_times(energy.s), binbounds=binbounds)
new_correlation = ift.create_power_operator(h_space, new_power)
new_parameters = energy.parameters
new_parameters['correlation'] = new_correlation
......
......@@ -19,7 +19,7 @@ class SeparationEnergy(Energy):
self.q = parameters['q']
pos_tanh = parameters['pos_tanh']
self.S = self.FFT.adjoint * self.correlation * self.FFT
self.S = self.FFT * self.correlation * self.FFT.adjoint
self.a = pos_tanh(self.position)
self.a_p = pos_tanh.derivative(self.position)
......@@ -52,7 +52,7 @@ class SeparationEnergy(Energy):
@property
def curvature(self):
point = self.q * exp(-self.u) * self.u_p ** 2
R = self.FFT * self.s_p
R = self.FFT.inverse * self.s_p
N = self.correlation
S = DiagonalOperator(1/(point + 1/self.var_x))
return WienerFilterCurvature(R=R, N=N, S=S, inverter=self.inverter)
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