Changed Data to Survey as per Issue #61

SimPEG/Inversion.py

@@ -32,24 +32,6 @@ class BaseInversion(object):
             self.opt.printers.insert(2,IterationPrinters.phi_d)
             self.opt.printers.insert(3,IterationPrinters.phi_m)
 
-    #TODO: Move this to the data class?
-    @property
-    def phi_d_target(self):
-        """
-        target for phi_d
-
-        By default this is the number of data.
-
-        Note that we do not set the target if it is None, but we return the default value.
-        """
-        if getattr(self, '_phi_d_target', None) is None:
-            return self.data.dobs.size #
-        return self._phi_d_target
-
-    @phi_d_target.setter
-    def phi_d_target(self, value):
-        self._phi_d_target = value
-
     @Utils.timeIt
     def run(self, m0):
         """run(m0)
@@ -70,50 +52,3 @@ class BaseInversion(object):
             **finish** is called at the end of the optimization.
         """
         pass
-
-
-    def save(self, group):
-        group.attrs['phi_d'] = self.phi_d
-        group.attrs['phi_m'] = self.phi_m
-        group.setArray('m', self.m)
-        group.setArray('dpred', self.dpred)
-
-
-
-# class Inversion(Cooling, Remember, BaseInversion):
-
-#     maxIter = 10
-#     name = "SimPEG Inversion"
-
-#     def __init__(self, prob, reg, opt, data, **kwargs):
-#         BaseInversion.__init__(self, prob, reg, opt, data, **kwargs)
-
-#         self.stoppers.append(StoppingCriteria.phi_d_target_Inversion)
-
-#         if StoppingCriteria.phi_d_target_Minimize not in self.opt.stoppers:
-#             self.opt.stoppers.append(StoppingCriteria.phi_d_target_Minimize)
-
-# class TimeSteppingInversion(Remember, BaseInversion):
-#     """
-#         A slightly different view on regularization parameters,
-#         let Beta be viewed as 1/dt, and timestep by updating the
-#         reference model every optimization iteration.
-#     """
-#     maxIter = 1
-#     name = "Time-Stepping SimPEG Inversion"
-
-#     def __init__(self, prob, reg, opt, data, **kwargs):
-#         BaseInversion.__init__(self, prob, reg, opt, data, **kwargs)
-
-#         self.stoppers.append(StoppingCriteria.phi_d_target_Inversion)
-
-#         if StoppingCriteria.phi_d_target_Minimize not in self.opt.stoppers:
-#             self.opt.stoppers.append(StoppingCriteria.phi_d_target_Minimize)
-
-#     def _startup_TimeSteppingInversion(self, m0):
-
-#         def _doEndIteration_updateMref(self, xt):
-#             if self.debug: 'Updating the reference model.'
-#             self.parent.reg.mref = self.xc
-
-#         self.opt.hook(_doEndIteration_updateMref, overwrite=True)

SimPEG/ObjFunction.py

@@ -1,7 +1,7 @@
-import Utils, Parameters, numpy as np, scipy.sparse as sp
+import Utils, Parameters, Survey, Problem, numpy as np, scipy.sparse as sp
 
 class BaseObjFunction(object):
-    """BaseObjFunction(data, reg, **kwargs)"""
+    """BaseObjFunction(forward, reg, **kwargs)"""
 
     __metaclass__ = Utils.SimPEGMetaClass
 
@@ -10,6 +10,9 @@ class BaseObjFunction(object):
     debug   = False  #: Print debugging information
     counter = None   #: Set this to a SimPEG.Utils.Counter() if you want to count things
 
+    surveyPair  = Survey.BaseSurvey
+    problemPair = Problem.BaseProblem
+
     name = 'Base Objective Function'   #: Name of the objective function
 
     u_current = None   #: The most current evaluated field
@@ -31,18 +34,19 @@ class BaseObjFunction(object):
     def objFunc(self): return self
     @property
     def opt(self): return getattr(self.parent,'opt',None)
-    @property
-    def prob(self): return self.data.prob
-    @property
-    def mesh(self): return self.data.prob.mesh
-    @property
-    def model(self): return self.data.prob.model
 
 
-    def __init__(self, data, reg, **kwargs):
+    def __init__(self, forward, reg, **kwargs):
         Utils.setKwargs(self, **kwargs)
 
-        self.data = data
+        assert forward.ispaired, 'The forward problem and survey must be paired.'
+        if isinstance(forward, self.surveyPair):
+            self.survey = forward
+            self.prob = forward.prob
+        elif isinstance(forward, self.problemPair):
+            self.prob = forward
+            self.survey = forward.survey
+
 
         self.reg = reg
         self.reg.parent = self
@@ -73,13 +77,13 @@ class BaseObjFunction(object):
         self.u_current = None
         self.m_current = m
 
-        u = self.data.prob.fields(m)
+        u = self.prob.fields(m)
         self.u_current = u
 
         phi_d = self.dataObj(m, u=u)
         phi_m = self.reg.modelObj(m)
 
-        self.dpred = self.data.dpred(m, u=u)  # This is a cheap matrix vector calculation.
+        self.dpred = self.survey.dpred(m, u=u)  # This is a cheap matrix vector calculation.
 
         self.phi_d, self.phi_d_last  = phi_d, self.phi_d
         self.phi_m, self.phi_m_last  = phi_m, self.phi_m
@@ -124,7 +128,7 @@ class BaseObjFunction(object):
             u is the field of interest; d_obs is the observed data; and W is the weighting matrix.
         """
         # TODO: ensure that this is a data is vector and Wd is a matrix.
-        R = self.data.residualWeighted(m, u=u)
+        R = self.survey.residualWeighted(m, u=u)
         return 0.5*np.vdot(R, R)
 
     @Utils.timeIt
@@ -160,11 +164,11 @@ class BaseObjFunction(object):
                 \\frac{\partial \mu_\\text{data}}{\partial \mathbf{m}} = \mathbf{J}^\\top \mathbf{W \circ R}
 
         """
-        if u is None: u = self.data.prob.fields(m)
+        if u is None: u = self.prob.fields(m)
 
-        R = self.data.residualWeighted(m, u=u)
+        R = self.survey.residualWeighted(m, u=u)
 
-        dmisfit = self.data.prob.Jtvec(m, self.data.Wd * R, u=u)
+        dmisfit = self.prob.Jtvec(m, self.survey.Wd * R, u=u)
 
         return dmisfit
 
@@ -204,12 +208,12 @@ class BaseObjFunction(object):
                 \\frac{\partial^2 \mu_\\text{data}}{\partial^2 \mathbf{m}} = \mathbf{J}^\\top \mathbf{W \circ W J}
 
         """
-        if u is None: u = self.data.prob.fields(m)
+        if u is None: u = self.prob.fields(m)
 
-        R = self.data.residualWeighted(m, u=u)
+        R = self.survey.residualWeighted(m, u=u)
 
         # TODO: abstract to different norms a little cleaner.
         #                                                     \/ it goes here. in l2 it is the identity.
-        dmisfit = self.data.prob.Jtvec_approx(m, self.data.Wd * self.data.Wd * self.data.prob.Jvec_approx(m, v, u=u), u=u)
+        dmisfit = self.prob.Jtvec_approx(m, self.survey.Wd * self.survey.Wd * self.prob.Jvec_approx(m, v, u=u), u=u)
 
         return dmisfit

SimPEG/Parameters.py

@@ -41,7 +41,7 @@ class Parameter(object):
     @property
     def reg(self): return self.parent.reg
     @property
-    def data(self): return self.parent.data
+    def survey(self): return self.parent.survey
     @property
     def prob(self): return self.parent.prob
     @property
@@ -131,13 +131,13 @@ class BetaEstimate(Parameter):
             :return: beta0
         """
         objFunc  = self.parent
-        data     = objFunc.data
+        survey     = objFunc.survey
 
         m = objFunc.m_current
         u = objFunc.u_current
 
         if u is None:
-            u = data.prob.fields(m)
+            u = survey.prob.fields(m)
 
         x0 = np.random.rand(*m.shape)
         t = x0.dot(objFunc.dataObj2Deriv(m,x0,u=u))

SimPEG/Problem.py

@@ -1,4 +1,4 @@
-import Utils, Data, numpy as np, scipy.sparse as sp
+import Utils, Survey, numpy as np, scipy.sparse as sp
 import Model
 
 class BaseProblem(object):
@@ -38,40 +38,42 @@ class BaseProblem(object):
 
     counter = None   #: A SimPEG.Utils.Counter object
 
-    dataPair = Data.BaseData
+    surveyPair = Survey.BaseSurvey
     modelPair = Model.BaseModel
 
-    def __init__(self, mesh, model, **kwargs):
+    def __init__(self, model, **kwargs):
         Utils.setKwargs(self, **kwargs)
-        self.mesh = mesh
         assert (isinstance(model, self.modelPair) or
             isinstance(model, Model.ComboModel) and isinstance(model.models[0], self.modelPair)
             ), "Model object must be an instance of a %s class."%(self.modelPair.__name__)
         self.model = model
 
     @property
-    def data(self):
+    def mesh(self): return self.model.mesh
+
+    @property
+    def survey(self):
         """
-        The data object for this problem.
+        The survey object for this problem.
         """
-        return getattr(self, '_data', None)
+        return getattr(self, '_survey', None)
 
     def pair(self, d):
-        """Bind a data to this problem instance using pointers."""
-        assert isinstance(d, self.dataPair), "Data object must be an instance of a %s class."%(self.dataPair.__name__)
+        """Bind a survey to this problem instance using pointers."""
+        assert isinstance(d, self.surveyPair), "Data object must be an instance of a %s class."%(self.surveyPair.__name__)
         if d.ispaired:
-            raise Exception("The data object is already paired to a problem. Use data.unpair()")
-        self._data = d
+            raise Exception("The survey object is already paired to a problem. Use survey.unpair()")
+        self._survey = d
         d._prob = self
 
     def unpair(self):
-        """Unbind a data from this problem instance."""
+        """Unbind a survey from this problem instance."""
         if not self.ispaired: return
-        self.data._prob = None
-        self._data = None
+        self.survey._prob = None
+        self._survey = None
 
     @property
-    def ispaired(self): return self.data is not None
+    def ispaired(self): return self.survey is not None
 
     @Utils.timeIt
     def Jvec(self, m, v, u=None):
@@ -156,25 +158,26 @@ class BaseProblem(object):
         """
         pass
 
-    def createSyntheticData(self, m, std=0.05, u=None, **geometry_kwargs):
+    #TODO: Rename and refactor to createSyntheticData
+    def createSyntheticSurvey(self, m, std=0.05, u=None, **geometry_kwargs):
         """
-            Create synthetic data given a model, and a standard deviation.
+            Create synthetic survey given a model, and a standard deviation.
 
             :param numpy.array m: geophysical model
             :param numpy.array std: standard deviation
-            :rtype: numpy.array, numpy.array
-            :return: dobs, Wd
+            :rtype: SurveyObject
+            :return: survey
 
             Returns the observed data with random Gaussian noise
             and Wd which is the same size as data, and can be used to weight the inversion.
         """
-        data = self.dataPair(mtrue=m, **geometry_kwargs)
-        data.pair(self)
-        data.dtrue = data.dpred(m, u=u)
-        noise = std*abs(data.dtrue)*np.random.randn(*data.dtrue.shape)
-        data.dobs = data.dtrue+noise
-        data.std = data.dobs*0 + std
-        return data
+        survey = self.surveyPair(mtrue=m, **geometry_kwargs)
+        survey.pair(self)
+        survey.dtrue = survey.dpred(m, u=u)
+        noise = std*abs(survey.dtrue)*np.random.randn(*survey.dtrue.shape)
+        survey.dobs = survey.dtrue+noise
+        survey.std = survey.dobs*0 + std
+        return survey
 
 
 

SimPEG/Regularization.py

@@ -46,7 +46,7 @@ class BaseRegularization(object):
     @property
     def prob(self): return self.parent.prob
     @property
-    def data(self): return self.parent.data
+    def survey(self): return self.parent.survey
     @property
     def mesh(self): return self.model.mesh
 

SimPEG/Survey.py

@@ -1,8 +1,8 @@
 import Utils, numpy as np
 
 
-class BaseData(object):
-    """Data holds the observed data, and the standard deviations."""
+class BaseSurvey(object):
+    """Survey holds the observed data, and the standard deviations."""
 
     __metaclass__ = Utils.SimPEGMetaClass
 
@@ -19,25 +19,25 @@ class BaseData(object):
     @property
     def prob(self):
         """
-        The geophysical problem that explains this data, use::
+        The geophysical problem that explains this survey, use::
 
-            data.pair(prob)
+            survey.pair(prob)
         """
         return getattr(self, '_prob', None)
 
     def pair(self, p):
-        """Bind a problem to this data instance using pointers"""
-        assert hasattr(p, 'dataPair'), "Problem must have an attribute 'dataPair'."
-        assert isinstance(self, p.dataPair), "Problem requires data object must be an instance of a %s class."%(p.dataPair.__name__)
+        """Bind a problem to this survey instance using pointers"""
+        assert hasattr(p, 'surveyPair'), "Problem must have an attribute 'surveyPair'."
+        assert isinstance(self, p.surveyPair), "Problem requires survey object must be an instance of a %s class."%(p.surveyPair.__name__)
         if p.ispaired:
-            raise Exception("The problem object is already paired to a data. Use prob.unpair()")
+            raise Exception("The problem object is already paired to a survey. Use prob.unpair()")
         self._prob = p
-        p._data = self
+        p._survey = self
 
     def unpair(self):
-        """Unbind a problem from this data instance"""
+        """Unbind a problem from this survey instance"""
         if not self.ispaired: return
-        self.prob._data = None
+        self.prob._survey = None
         self._prob = None
 
     @property
@@ -108,7 +108,7 @@ class BaseData(object):
         """
             Data weighting matrix. This is a covariance matrix used in::
 
-                def data.residualWeighted(m,u=None):
+                def residualWeighted(m,u=None):
                     return self.Wd*self.residual(m, u=u)
 
             By default, this is based on the norm of the data plus a noise floor.
@@ -139,20 +139,82 @@ class BaseData(object):
         """
         return Utils.mkvc(self.Wd*self.residual(m, u=u))
 
-    @property
-    def RHS(self):
-        """
-            Source matrix.
-        """
-        return getattr(self, '_RHS', None)
-    @RHS.setter
-    def RHS(self, value):
-        self._RHS = value
-
     @property
     def isSynthetic(self):
         "Check if the data is synthetic."
-        return (self.mtrue is not None)
+        return self.mtrue is not None
+
+
+    #TODO: Move this to the data class?
+    # @property
+    # def phi_d_target(self):
+    #     """
+    #     target for phi_d
+
+    #     By default this is the number of data.
+
+    #     Note that we do not set the target if it is None, but we return the default value.
+    #     """
+    #     if getattr(self, '_phi_d_target', None) is None:
+    #         return self.data.dobs.size #
+    #     return self._phi_d_target
+
+    # @phi_d_target.setter
+    # def phi_d_target(self, value):
+    #     self._phi_d_target = value
+
+
+class BaseRxList(object):
+    """SimPEG Receiver List Object"""
+
+    locs = None   #: Locations (nRx x 3)
+
+    knownRxTypes = None  #: Set this to a list of strings to ensure that txType is known
+
+    def __init__(self, locs, rxType, **kwargs):
+        self.locs = locs
+        self.rxType = rxType
+        Utils.setKwargs(self, **kwargs)
+
+    @property
+    def rxType(self):
+        """Receiver Type"""
+        return getattr(self, '_rxType', None)
+    @rxType.setter
+    def rxType(self, value):
+        known = self.knownRxTypes
+        if known is not None:
+            assert value in known, "rxType must be in ['%s']" % ("', '".join(known))
+        self._rxType = value
+
+class BaseTx(object):
+    """SimPEG Transmitter Object"""
+
+    loc    = None #: Location [x,y,z]
+
+    rxList = None #: SimPEG Receiver List
+    rxListPair = BaseRxList
+
+    knownTxTypes = None #: Set this to a list of strings to ensure that txType is known
+
+    def __init__(self, loc, txType, rxList, **kwargs):
+        assert isinstance(rxList, self.rxListPair), 'rxList must be a %s'%self.rxListPair.__name__
+
+        self.loc    = loc
+        self.txType = txType
+        self.rxList = rxList
+        Utils.setKwargs(self, **kwargs)
+
+    @property
+    def txType(self):
+        """Transmitter Type"""
+        return getattr(self, '_txType', None)
+    @txType.setter
+    def txType(self, value):
+        known = self.knownTxTypes
+        if known is not None:
+            assert value in known, "txType must be in ['%s']" % ("', '".join(known))
+        self._txType = value
 
 if __name__ == '__main__':
     d = BaseData()

SimPEG/Utils/__init__.py

@@ -159,11 +159,11 @@ def requires(var):
 
         .. note::
 
-            To use data.%s(), SimPEG requires that a problem be bound to the data.
+            To use survey.%s(), SimPEG requires that a problem be bound to the survey.
             If a problem has not been bound, an Exception will be raised.
             To bind a problem to the Data object::
 
-                data.pair(myProblem)
+                survey.pair(myProblem)
 
             """ % f.__name__
         else:

SimPEG/__init__.py

@@ -5,7 +5,7 @@ from Solver import Solver
 import Mesh
 import Model
 import Problem
-import Data
+import Survey
 import Regularization
 import ObjFunction
 import Optimization

Tutorials/Linear.py

@@ -5,8 +5,8 @@ import matplotlib.pyplot as plt
 class LinearProblem(Problem.BaseProblem):
     """docstring for LinearProblem"""
 
-    def __init__(self, mesh, model, G, **kwargs):
-        Problem.BaseProblem.__init__(self, mesh, model, **kwargs)
+    def __init__(self, model, G, **kwargs):
+        Problem.BaseProblem.__init__(self, model, **kwargs)
         self.G = G
 
     def fields(self, m, u=None):
@@ -20,8 +20,7 @@ class LinearProblem(Problem.BaseProblem):
 
 
 def example(N):
-    h = np.ones(N)/N
-    M = Mesh.TensorMesh([h])
+    M = Mesh.TensorMesh([N])
 
     nk = 20
     jk = np.linspace(1.,20.,nk)
@@ -43,8 +42,8 @@ def example(N):
 
 
     model = Model.BaseModel(M)
-    prob = LinearProblem(M, model, G)
-    data = prob.createSyntheticData(mtrue, std=0.01)
+    prob = LinearProblem(model, G)
+    data = prob.createSyntheticSurvey(mtrue, std=0.01)
 
     return prob, data, model