Skip to content
Snippets Groups Projects
Commit da72aa9a authored by Fengshan Zheng's avatar Fengshan Zheng
Browse files

Reconstruction Updated

parent c89aaf69
No related branches found
No related tags found
1 merge request!12Reconstruction
...@@ -3,7 +3,7 @@ ...@@ -3,7 +3,7 @@
# Author: J. Caron # Author: J. Caron
# #
"""This module provides the :class:`~.Costfunction` class which represents a strategy to calculate """This module provides the :class:`~.Costfunction` class which represents a strategy to calculate
the so called `cost` of a threedimensional magnetization distribution.""" the so called `cost` of a three dimensional magnetization distribution."""
import logging import logging
...@@ -20,7 +20,7 @@ class Costfunction(object): ...@@ -20,7 +20,7 @@ class Costfunction(object):
Represents a strategy for the calculation of the `cost` of a 3D magnetic distribution in Represents a strategy for the calculation of the `cost` of a 3D magnetic distribution in
relation to two-dimensional phase maps. The `cost` is a measure for the difference of the relation to two-dimensional phase maps. The `cost` is a measure for the difference of the
simulated phase maps from the magnetic distributions to the given set of phase maps and from simulated phase maps from the magnetic distributions to the given set of phase maps and from
a priori knowledge represented by a :class:`~.Regularisator` object. Furthermore this class a prior knowledge represented by a :class:`~.Regularisator` object. Furthermore this class
provides convenient methods for the calculation of the derivative :func:`~.jac` or the product provides convenient methods for the calculation of the derivative :func:`~.jac` or the product
with the Hessian matrix :func:`~.hess_dot` of the costfunction, which can be used by with the Hessian matrix :func:`~.hess_dot` of the costfunction, which can be used by
optimizers. All required data should be given in a :class:`~DataSet` object. optimizers. All required data should be given in a :class:`~DataSet` object.
...@@ -41,7 +41,7 @@ class Costfunction(object): ...@@ -41,7 +41,7 @@ class Costfunction(object):
Size of the input space. Size of the input space.
Se_inv : :class:`~numpy.ndarray` (N=2), optional Se_inv : :class:`~numpy.ndarray` (N=2), optional
Inverted covariance matrix of the measurement errors. The matrix has size `m x m` with m Inverted covariance matrix of the measurement errors. The matrix has size `m x m` with m
being the length of the targetvector y. being the length of the target vector y.
""" """
......
...@@ -16,15 +16,15 @@ import logging ...@@ -16,15 +16,15 @@ import logging
import numpy as np import numpy as np
from pyramid.fielddata import VectorData from pyramid.fielddata import VectorData, ScalarData
__all__ = ['optimize_linear', 'optimize_nonlin', 'optimize_splitbregman'] __all__ = ['optimize_linear', 'optimize_linear_charge', 'optimize_nonlin', 'optimize_splitbregman']
_log = logging.getLogger(__name__) _log = logging.getLogger(__name__)
def optimize_linear(costfunction, mag_0=None, ramp_0=None, max_iter=None, verbose=False): def optimize_linear(costfunction, mag_0=None, ramp_0=None, max_iter=None, verbose=False):
"""Reconstruct a three-dimensional magnetic distribution from given phase maps via the """Reconstruct a three-dimensional magnetic distribution from given phase maps via the
conjugate gradient optimizaion method :func:`~.scipy.sparse.linalg.cg`. conjugate gradient optimization method :func:`~.scipy.sparse.linalg.cg`.
Blazingly fast for l2-based cost functions. Blazingly fast for l2-based cost functions.
Parameters Parameters
...@@ -35,11 +35,11 @@ def optimize_linear(costfunction, mag_0=None, ramp_0=None, max_iter=None, verbos ...@@ -35,11 +35,11 @@ def optimize_linear(costfunction, mag_0=None, ramp_0=None, max_iter=None, verbos
mag_0: :class:`~.VectorData` mag_0: :class:`~.VectorData`
The starting magnetisation distribution used for the reconstruction. A zero vector will be The starting magnetisation distribution used for the reconstruction. A zero vector will be
used if no VectorData object is specified. used if no VectorData object is specified.
mag_0: :class:`~.Ramp` ramp_0: :class:`~.Ramp`
The starting ramp for the reconstruction. A zero vector will be The starting ramp for the reconstruction. A zero vector will be
used if no Ramp object is specified. used if no Ramp object is specified.
max_iter : int, optional max_iter : int, optional
The maximum number of iterations for the opimization. The maximum number of iterations for the optimization.
verbose: bool, optional verbose: bool, optional
If set to True, information like a progressbar is displayed during reconstruction. If set to True, information like a progressbar is displayed during reconstruction.
The default is False. The default is False.
...@@ -78,8 +78,8 @@ def optimize_linear(costfunction, mag_0=None, ramp_0=None, max_iter=None, verbos ...@@ -78,8 +78,8 @@ def optimize_linear(costfunction, mag_0=None, ramp_0=None, max_iter=None, verbos
def optimize_linear_charge(costfunction, charge_0=None, ramp_0=None, max_iter=None, verbose=False): def optimize_linear_charge(costfunction, charge_0=None, ramp_0=None, max_iter=None, verbose=False):
"""Reconstruct a three-dimensional magnetic distribution from given phase maps via the """Reconstruct a three-dimensional charge distribution from given phase maps via the
conjugate gradient optimizaion method :func:`~.scipy.sparse.linalg.cg`. conjugate gradient optimization method :func:`~.scipy.sparse.linalg.cg`.
Blazingly fast for l2-based cost functions. Blazingly fast for l2-based cost functions.
Parameters Parameters
...@@ -87,22 +87,22 @@ def optimize_linear_charge(costfunction, charge_0=None, ramp_0=None, max_iter=No ...@@ -87,22 +87,22 @@ def optimize_linear_charge(costfunction, charge_0=None, ramp_0=None, max_iter=No
costfunction : :class:`~.Costfunction` costfunction : :class:`~.Costfunction`
A :class:`~.Costfunction` object which implements a specified forward model and A :class:`~.Costfunction` object which implements a specified forward model and
regularisator which is minimized in the optimization process. regularisator which is minimized in the optimization process.
mag_0: :class:`~.VectorData` charge_0: :class:`~.VectorData`
The starting magnetisation distribution used for the reconstruction. A zero vector will be The starting magnetisation distribution used for the reconstruction. A zero vector will be
used if no VectorData object is specified. used if no VectorData object is specified.
mag_0: :class:`~.Ramp` ramp_0: :class:`~.Ramp`
The starting ramp for the reconstruction. A zero vector will be The starting ramp for the reconstruction. A zero vector will be
used if no Ramp object is specified. used if no Ramp object is specified.
max_iter : int, optional max_iter : int, optional
The maximum number of iterations for the opimization. The maximum number of iterations for the optimization.
verbose: bool, optional verbose: bool, optional
If set to True, information like a progressbar is displayed during reconstruction. If set to True, information like a progressbar is displayed during reconstruction.
The default is False. The default is False.
Returns Returns
------- -------
magdata : :class:`~pyramid.fielddata.VectorData` elecdata : :class:`~pyramid.fielddata.ScalarData`
The reconstructed magnetic distribution as a :class:`~.VectorData` object. The reconstructed charge distribution as a :class:`~.ScalarData` object.
""" """
import jutil.cg as jcg import jutil.cg as jcg
...@@ -113,8 +113,8 @@ def optimize_linear_charge(costfunction, charge_0=None, ramp_0=None, max_iter=No ...@@ -113,8 +113,8 @@ def optimize_linear_charge(costfunction, charge_0=None, ramp_0=None, max_iter=No
# Get starting distribution vector x_0: # Get starting distribution vector x_0:
x_0 = np.empty(costfunction.n) x_0 = np.empty(costfunction.n)
if charge_0 is not None: if charge_0 is not None:
costfunction.fwd_model.magdata = charge_0 costfunction.fwd_model.elecdata = charge_0
x_0[:data_set.n] = costfunction.fwd_model.magdata.get_vector(mask=data_set.mask) x_0[:data_set.n] = costfunction.fwd_model.elecdata.get_vector(mask=data_set.mask)
if ramp_0 is not None: if ramp_0 is not None:
ramp_vec = ramp_0.param_cache.ravel() ramp_vec = ramp_0.param_cache.ravel()
else: else:
...@@ -126,10 +126,10 @@ def optimize_linear_charge(costfunction, charge_0=None, ramp_0=None, max_iter=No ...@@ -126,10 +126,10 @@ def optimize_linear_charge(costfunction, charge_0=None, ramp_0=None, max_iter=No
_log.info('Cost after optimization: {:.3e}'.format(costfunction(x_opt))) _log.info('Cost after optimization: {:.3e}'.format(costfunction(x_opt)))
# Cut ramp parameters if necessary (this also saves the final parameters in the ramp class!): # Cut ramp parameters if necessary (this also saves the final parameters in the ramp class!):
x_opt = costfunction.fwd_model.ramp.extract_ramp_params(x_opt) x_opt = costfunction.fwd_model.ramp.extract_ramp_params(x_opt)
# Create and return fitting VectorData object: # Create and return fitting ScalarData object:
mag_opt = VectorData(data_set.a, np.zeros((3,) + data_set.dim)) charge_opt = ScalarData(data_set.a, np.zeros(data_set.dim))
mag_opt.set_vector(x_opt, data_set.mask) charge_opt.set_vector(x_opt, data_set.mask)
return mag_opt return charge_opt
def optimize_nonlin(costfunction, first_guess=None): def optimize_nonlin(costfunction, first_guess=None):
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment