- Restructured project

crust1/__init__.py

@@ -0,0 +1,3 @@
+#from .models import CrustModel
+
+#__all__ = ['CrustModel']

model/crust1.py → crust1/models.py

@@ -1,8 +1,11 @@
 import numpy as np
 from math import floor
 
+import inspect
+import os
 
-class crustModel:
+
+class CrustModel:
     """
     Top level model object to retreive information from the LLNL Crust 1.0
     model.
@@ -25,12 +28,20 @@ class crustModel:
     Names of the nine possible layers in the model
     """
 
-    def __init__(self):
+    def __init__(self, data_path=os.path.dirname(inspect.stack()[0][1])  + '/data/'):
+        """
+        CrustModel constructor
+
+        Parameters
+        ----------
+        * data_path: Path to where the data is located. default is ./data/. Use this if you want to load external data files.
+        """
+
         # Read in data files
-        self.vp = np.loadtxt('crust1.vp')
-        self.vs = np.loadtxt('crust1.vs')
-        self.rho = np.loadtxt('crust1.rho')
-        self.bnds = np.loadtxt('crust1.bnds')
+        self.vp = np.loadtxt(data_path + 'crust1.vp')
+        self.vs = np.loadtxt(data_path + 'crust1.vs')
+        self.rho = np.loadtxt(data_path + 'crust1.rho')
+        self.bnds = np.loadtxt(data_path + 'crust1.bnds')
 
         # Reshape to a lon,lat,layer grid. The 0,0 index value
         # is at 90 south and 180 latitude.
@@ -39,10 +50,14 @@ def __init__(self):
         self.rho = self.rho.reshape((180, 360, 9))
         self.bnds = self.bnds.reshape((180, 360, 9))
 
+        # Constains a list with all the layers names.
         self.layer_names = ["water", "ice", "upper_sediments",
                             "middle_sediments", "lower_sediments",
                             "upper_crust", "middle_crust", "lower_crust",
                             "mantle"]
+
+        # Contains the list of the data that each layer has..
+        self.layer_data = ['vp','vs','rho','layer_thickness','bnd']
 
     def _get_index(self, lat, lon):
         """
@@ -77,7 +92,7 @@ def _get_index(self, lat, lon):
 
         return int(ilat), int(ilon)
 
-    def get_point(self, lat, lon):
+    def get_point(self, lat, lon, include_no_thickness=False):
         """
         Returns a model for a given latitude and longitude. Note that the model
         is only defined on a 1 degree grid starting at 89.5 and -179.5.
@@ -89,12 +104,15 @@ def get_point(self, lat, lon):
 
         lat : flaot
         Longitude of interest
+
+        include_no_thickness: bool
+        Whether to include layers without thickness or not. Default: False
 
         Returns
         -------
         model_layers : dict
-        Dictionary of layers with the keys as layer names and the values as
-        a list of vp, vs, density, layer thickness, and the top of the layer
+        Dictionary of layers with the keys as layer names and a dictionary with 
+        the values of vp, vs, density, layer thickness, and the top of the layer
         with respect to sea level.
         """
 
@@ -108,14 +126,28 @@ def get_point(self, lat, lon):
         model_layers = dict()
 
         for i, layer in enumerate(self.layer_names):
-            vp = self.vp[ilat, ilon][i]
-            vs = self.vs[ilat, ilon][i]
-            rho = self.rho[ilat, ilon][i]
-            bnd = self.bnds[ilat, ilon][i]
+
+            # Read the layer thickness
             layer_thickness = thickness[i]
-
-            # If the layer has thickness or is the mantle, write it
-            if layer_thickness >= 0.01 or layer == "mantle":
-                model_layers[layer] = [vp, vs, rho, layer_thickness, bnd]
+
+            # If no thickness (include_no_thickness=False), check if the layer has thickness (>= 0.01)
+            # or if it's the mantle, write it
+            if ((include_no_thickness == False) & (layer_thickness >= 0.01)) or layer == "mantle":
+
+                # Read each parameter
+                vp = self.vp[ilat, ilon][i]
+                vs = self.vs[ilat, ilon][i]
+                rho = self.rho[ilat, ilon][i]
+                bnd = self.bnds[ilat, ilon][i]
+
+                #####################################################################
+                # Make sure that KEYS defines here, are the same as self.layer_data #
+                #####################################################################
+                model_layers[layer] = {
+                        'vp': vp, # 
+                        'vs': vs,
+                        'rho': rho,
+                        'layer_thickness':layer_thickness,
+                        'bnd': bnd}
 
         return model_layers

setup.py

@@ -0,0 +1,12 @@
+from distutils.core import setup
+
+setup(name='crust1',
+      version='1.0.1',
+      description='A python package to retreive information from the LLNL Crust 1.0 model.',
+      author='John Leeman, Juan Manuel Haedo',
+      author_email='john@leemangeophysical.com, juanu@juanu.com.ar',
+      url='https://github.com/jrleeman/Crust1.0/',
+      packages=['crust1'], # Package name
+      package_data={'crust1': ['data/*.*']}, # Include data files
+      include_package_data=True
+     )