Commit 939c6c7d authored by Philipp Arras's avatar Philipp Arras

Some pep8

parent d5563c04
Pipeline #24109 failed with stage
in 4 minutes and 10 seconds
......@@ -22,9 +22,9 @@ from ..minimization.energy import Energy
class NoiseEnergy(Energy):
def __init__(self, position, d, m, D, t, HarmonicTransform, Instrument, nonlinearity,
alpha, q, Projection, munit=1., sunit=1., dunit=1., samples=3, sample_list=None,
inverter=None):
def __init__(self, position, d, m, D, t, HarmonicTransform, Instrument,
nonlinearity, alpha, q, Projection, munit=1., sunit=1.,
dunit=1., samples=3, sample_list=None, inverter=None):
super(NoiseEnergy, self).__init__(position=position)
self.m = m
self.D = D
......@@ -53,14 +53,15 @@ class NoiseEnergy(Energy):
self.sample_list = sample_list
self.inverter = inverter
A = Projection.adjoint_times(munit * exp(.5*self.t)) # unit: munit
A = Projection.adjoint_times(munit * exp(.5 * self.t)) # unit: munit
map_s = self.ht(A * m)
self._gradient = None
for sample in self.sample_list:
map_s = self.ht(A * sample)
residual = self.d - self.Instrument(sunit * self.nonlinearity(map_s))
residual = self.d - \
self.Instrument(sunit * self.nonlinearity(map_s))
lh = .5 * residual.vdot(self.N.inverse_times(residual))
grad = -.5 * self.N.inverse_times(residual.conjugate() * residual)
......@@ -71,11 +72,13 @@ class NoiseEnergy(Energy):
self._value += lh
self._gradient += grad
self._value *= 1./len(self.sample_list)
self._value += .5 * self.position.integrate() + (self.alpha - 1.).vdot(self.position) + self.q.vdot(exp(-self.position))
self._value *= 1. / len(self.sample_list)
self._value += .5 * self.position.integrate()
self._value += (self.alpha - 1.).vdot(self.position) + \
self.q.vdot(exp(-self.position))
self._gradient *= 1./len(self.sample_list)
self._gradient += (self.alpha-0.5) - self.q * (exp(-self.position))
self._gradient *= 1. / len(self.sample_list)
self._gradient += (self.alpha - 0.5) - self.q * (exp(-self.position))
def at(self, position):
return self.__class__(
......
......@@ -20,12 +20,15 @@ from ..operators.inversion_enabler import InversionEnabler
from .response_operators import LinearizedPowerResponse
def NonlinearPowerCurvature(position, HarmonicTransform, Instrument, nonlinearity,
Projection, N, T, sample_list, inverter, munit=1., sunit=1.):
def NonlinearPowerCurvature(position, HarmonicTransform, Instrument,
nonlinearity, Projection, N, T, sample_list,
inverter, munit=1., sunit=1.):
result = None
for sample in sample_list:
LinR = LinearizedPowerResponse(Instrument, nonlinearity, HarmonicTransform, Projection, position, sample, munit, sunit)
op = LinR.adjoint*N.inverse*LinR
LinR = LinearizedPowerResponse(Instrument, nonlinearity,
HarmonicTransform, Projection, position,
sample, munit, sunit)
op = LinR.adjoint * N.inverse * LinR
result = op if result is None else result + op
result = result*(1./len(sample_list)) + T
result = result * (1. / len(sample_list)) + T
return InversionEnabler(result, inverter)
......@@ -51,7 +51,9 @@ class NonlinearPowerEnergy(Energy):
default : 3
"""
def __init__(self, position, d, N, m, D, HarmonicTransform, Instrument, nonlinearity, Projection, sigma=0., samples=3, sample_list=None, munit=1., sunit=1., inverter=None):
def __init__(self, position, d, N, m, D, HarmonicTransform, Instrument,
nonlinearity, Projection, sigma=0., samples=3,
sample_list=None, munit=1., sunit=1., inverter=None):
super(NonlinearPowerEnergy, self).__init__(position)
self.d, self.N, self.m, self.D, self.ht = d, N, m, D, HarmonicTransform
self.Instrument = Instrument
......@@ -72,16 +74,25 @@ class NonlinearPowerEnergy(Energy):
self.T = SmoothnessOperator(domain=self.position.domain[0],
strength=sigma, logarithmic=True)
A = Projection.adjoint_times(munit * exp(.5*position)) # unit: munit
A = Projection.adjoint_times(munit * exp(.5 * position)) # unit: munit
map_s = self.ht(A * m)
Tpos = self.T(position)
self._gradient = None
for sample in self.sample_list:
map_s = self.ht(A * sample)
LinR = LinearizedPowerResponse(Instrument, nonlinearity, self.ht, Projection, position, sample, munit, sunit)
residual = self.d - self.Instrument(sunit * self.nonlinearity(map_s))
LinR = LinearizedPowerResponse(
Instrument,
nonlinearity,
self.ht,
Projection,
position,
sample,
munit,
sunit)
residual = self.d - \
self.Instrument(sunit * self.nonlinearity(map_s))
lh = 0.5 * residual.vdot(self.N.inverse_times(residual))
grad = LinR.adjoint_times(self.N.inverse_times(residual))
......@@ -92,9 +103,9 @@ class NonlinearPowerEnergy(Energy):
self._value += lh
self._gradient += grad
self._value *= 1./len(self.sample_list)
self._value += 0.5*self.position.vdot(Tpos)
self._gradient *= -1./len(self.sample_list)
self._value *= 1. / len(self.sample_list)
self._value += 0.5 * self.position.vdot(Tpos)
self._gradient *= -1. / len(self.sample_list)
self._gradient += Tpos
def at(self, position):
......
......@@ -23,8 +23,8 @@ from .response_operators import LinearizedSignalResponse
class NonlinearWienerFilterEnergy(Energy):
def __init__(self, position, d, Instrument, nonlinearity, HarmonicTransform, power, N, S, sunit=1.,
inverter=None):
def __init__(self, position, d, Instrument, nonlinearity, HarmonicTransform,
power, N, S, sunit=1., inverter=None):
super(NonlinearWienerFilterEnergy, self).__init__(position=position)
self.d = d
self.sunit = sunit
......
......@@ -19,12 +19,15 @@
from ..field import exp
def LinearizedSignalResponse(Instrument, nonlinearity, HarmonicTransform, power, s, sunit):
return sunit * (Instrument * nonlinearity.derivative(s) * HarmonicTransform * power)
def LinearizedSignalResponse(Instrument, nonlinearity, HarmonicTransform, power,
s, sunit):
return sunit * (Instrument * nonlinearity.derivative(s) *
HarmonicTransform * power)
def LinearizedPowerResponse(Instrument, nonlinearity, HarmonicTransform, Projection, t, m, munit, sunit):
power = exp(0.5*t) * munit
def LinearizedPowerResponse(Instrument, nonlinearity, HarmonicTransform,
Projection, t, m, munit, sunit):
power = exp(0.5 * t) * munit
position = HarmonicTransform(Projection.adjoint_times(power) * m)
linearization = nonlinearity.derivative(position)
return sunit * (0.5 * Instrument * linearization * HarmonicTransform * m *
......
......@@ -45,7 +45,7 @@ class WienerFilterCurvature(EndomorphicOperator):
self.R = R
self.N = N
self.S = S
op = R.adjoint*N.inverse*R + S.inverse
op = R.adjoint * N.inverse * R + S.inverse
self._op = InversionEnabler(op, inverter, S.times)
@property
......@@ -91,8 +91,10 @@ class WienerFilterCurvature(EndomorphicOperator):
def generate_posterior_sample2(self):
power = self.S.diagonal()
mock_signal = Field.from_random(random_type="normal",
domain=self.S.domain, dtype=power.dtype)
mock_signal = Field.from_random(
random_type="normal",
domain=self.S.domain,
dtype=power.dtype)
mock_signal *= sqrt(power)
noise = self.N.diagonal()
......
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