Repository: hobuinc/silvimetric (language: python)

Issue:
Added FUSION metrics
Added FUSION metrics except the cover-related metrics that require a height threshold. I think I have everything in place but haven't looked at the values that come out. I suspect there are much faster ways to compute some things. For example, I compute 15 percentile height metrics using a function call for each. It is easy with numpy.percentile() to compute all with one call but things aren't set up to pull a list of metrics into TileDB. Same thing goes for the l-moments and l-moment ratios. Having each in a separate function means the l-moments are computed 7 times. Same thing could be done with a single call.


Discussion / hints:
it should be a default that the type for all metrics is `np.float32` so we don't have to copypasta that allover the place.

I also wonder if each metric shouldn't provide a ``test()`` harness that can be run by pytest to protect from dry rot. 
Is there an example of what the test() should look like? I poked around in the tests a bit but didn't see anything. 

Most of the metrics are computed in simple ways but i suppose the underlying numpy methods could change and cause problems or differences.
> Is there an example of what the test() should look like? I

I was thinking add something like this to `Metric`:

```
    def test(self, vector):
        try:
            results = self._method(vector)
        except Exception as e:
            # if we div0 or or throw or something we failed
            raise e

        # running the method over a vector of input 
        # should only return a single value that is a number (for now)
        assert isinstance(results, (int, float, complex)) and not isinstance(results, bool)
        return True
```

Then the test harness would pass in stuff that could be weird (like a vector of 0s or a string or other bad input). If any of the `test` calls raise an exception, we report to the user, etc.

        

@bmcgaughey1 I made a PR with some changes to fix some small stuff as well as adding a simple framework for testing metrics. If you want to look that over, I can merge that into this PR whenever.

One thing that I'm finding with this PR is that some of the methods end up dividing by zero with data that hasn't been sanitized. The ones I can think of off the top of my head are `m_crr` and `lmom4`.
The divide by zero makes sense. Do you think it is best to test for divide by zero or just trap the exception and return a value indicating the metric couldn't be computed?
I accepted the change to the call for np.percentile 1 to 100 instead of 1 to 99. the 100th percentile is max and it isn't used in the profilearea metric because of potential outliers. No problem leaving it at 100 but it adds time (arguable a small amount of time but time nonetheless)
> I accepted the change to the call for np.percentile 1 to 100 instead of 1 to 99. the 100th percentile is max and it isn't used in the profilearea metric because of potential outliers. No problem leaving it at 100 but it adds time (arguable a small amount of time but time nonetheless)

The range method in python returns results that are inclusive on the bottom and exclusive on the top, so `p[98]` was throwing an out of bounds error [here](https://github.com/hobuinc/silvimetric/pull/70/files#diff-4136e4fbb3c10e3226f73f1cebd4b0859755c14224f4df4f25fda191d3d69aadR250)
```
>>> print(list(range(9)))
[0, 1, 2, 3, 4, 5, 6, 7, 8]
>>> print(list(range(1,10)))
[1, 2, 3, 4, 5, 6, 7, 8, 9]
```
> The divide by zero makes sense. Do you think it is best to test for divide by zero or just trap the exception and return a value indicating the metric couldn't be computed?

I guess we should leave it up to the Metric to decide what it wants to do, but it should be a requirement that if a Metric throws for ANY input on computation, it's `default` value is written into TileDB. Then a Metric writer has the option of doing `nan` or `min/max` or `-9999` or whatever it wants. 
Still waiting for how we want to handle the height thresholds. There is a "base" threshold for all of the simple metrics, a cover threshold, and a "no threshold" case for computing metrics by height strata. The last actually needs a list of thresholds that define the height strata. Seems to me that we want these are part of init. Calls to simple metric functions could filter the data prior to the call or pass the threshold value to the function and let the filtering happen in the function. For the cover metrics, ratios of return counts above the threshold to total return count, the function will need all the points and the threshold. For strata, we might need to move ahead with the ability to have the metric function return a list of values. However, we could call the strata-metric functions with a lower and upper threshold...one call for each strata. just one more case where it could be done more efficiently (compute multiple metrics with a single function call).
So what we definitely need in order to function correctly are optional arguments for low and high thresholds for each metric.

I agree that it would be simpler and more efficient for the `Shatter` process if we allowed for metrics to output multiple results, I think the complexity is transferred to the `Extract` step, because we would have to include a lot more information in what's being passed to it in order to know how to extract that information.

The biggest problem I see with this is how to generalize this information across `Metrics`, and what is the relevance of each entry in a list passed to `Extract`? For `percentile`, each entry may correspond to it's percentage, but how can we generalize that across multiple `Metrics`? 

I don't immediately see an elegant way to solve this problem. I think we should probably do it eventually, but maybe we should just do the thresholds for now and wait to do multiple value outputs until we have a fleshed out approach to it?

What are your thoughts @bmcgaughey1 ?
The single output metric functions makes sense for now. Adding thresholds passed to each metric function also makes sense. These are set globally for all metrics (but one set (high & low) for the base metrics, another for the cover metrics, and then either a list for the strata metrics or the function is called for each strata with strata lower & upper thresholds). I would see these thresholds stored in TileDB along with the list of metrics.

There aren't many of the FUSION metrics that would have multiple return values. percentiles and l-moment metrics are the ones that stand out (just because there are easy ways to compute them). I agree that having a flexible way to handle the multiple return cases is not trivial (or at least not trivial to make it robust).
In general, I pushed the changes to add thresholds so I could get to the code while on the road this week. I didn't deal with the testing framework to add the threshold parameters (hence the failed tests). The way I added thresholds isn't correct given that the metric function gets a height threshold but applies it (if I left in the code that uses the threshold) to a different dimension of the data.
I pushed changes to add doc strings to all metric functions and remove use of height thresholds. I didn't have time to test code so not sure if tests will crash. Looks like things are failing

Fix patch:
diff --git a/.gitignore b/.gitignore
index aeecde4..85d76b8 100644
--- a/.gitignore
+++ b/.gitignore
@@ -4,9 +4,12 @@ __pycache__/
 build/
 dist/
 stats/
-*.html
 .env
 
+*.html
+**/html/**
+**/_build/**
+
 #tiledb
 **/__fragments/**
 **/__commits/**
diff --git a/README.md b/README.md
index 0bde116..d7be137 100644
--- a/README.md
+++ b/README.md
@@ -1,8 +1,7 @@
 ## SilviMetric
 
 
-SilviMetric is an open source library and set of command line utilities for extracting point cloud metrics
-into a TileDB database. See https://silvimetric.com for documentation and tutorials.
+SilviMetric is an open source library and set of command line utilities for extracting point cloud metrics into a TileDB database. See https://silvimetric.com for documentation and tutorials.
 
 [<img src="https://github.com/hobuinc/silvimetric/blob/main/docs/source/logo/Logos/PNG/SilviMeteric_Logo_2c.png?raw=true">](https://silvimetric.com/)
 
@@ -10,3 +9,26 @@ into a TileDB database. See https://silvimetric.com for documentation and tutori
 
 GitHub hosts the project at https://github.com/hobuinc/silvimetric
 
+
+### Installation
+These scripts will install `Silvimetric` dependencies as python libraries to the conda environment silvimetric.
+
+`Silvimetric` requires that we install some packages from `conda` (`TileDB`, and `PDAL`) so it's usually easier to use only `conda` to handle your environment. If this is unavailable to you, you will need to install `TileDB` and `PDAL` from source before installing the python packages that are dependent on those (`python-pdal` and `tiledb-py`).
+
+##### Pip and Conda
+
+```
+conda env create -f https://raw.githubusercontent.com/hobuinc/silvimetric/main/environment.yml
+conda activate silvimetric
+pip install silvimetric
+```
+
+##### Source
+
+```
+git clone https://github.com/hobuinc/silvimetric.git && cd silvimetric
+conda env create -f environment.yml && conda activate silvimetric
+pip install .
+# To install for development:
+# pip install -e .
+```
\ No newline at end of file
diff --git a/docs/create.sh b/docs/create.sh
new file mode 100755
index 0000000..6bf92bd
--- /dev/null
+++ b/docs/create.sh
@@ -0,0 +1,4 @@
+#!/bin/bash
+
+# This assumes you've already done pip install -r requirements.txt in docs dir
+python -m sphinx -T -b html -d _build/doctrees -D language=en source html
\ No newline at end of file
diff --git a/docs/requirements.txt b/docs/requirements.txt
index 7cf05da..a06057b 100644
--- a/docs/requirements.txt
+++ b/docs/requirements.txt
@@ -3,3 +3,7 @@ sphinxcontrib-bibtex
 sphinx-rtd-theme
 sphinxcontrib-spelling
 sphinx-notfound-page
+numpy
+tiledb
+dask
+.
\ No newline at end of file
diff --git a/docs/source/api/commands/delete.rst b/docs/source/api/commands/delete.rst
new file mode 100644
index 0000000..8e53fca
--- /dev/null
+++ b/docs/source/api/commands/delete.rst
@@ -0,0 +1,4 @@
+Delete
+======
+
+.. autofunction:: silvimetric.commands.manage.delete
\ No newline at end of file
diff --git a/docs/source/api/commands/extract.rst b/docs/source/api/commands/extract.rst
new file mode 100644
index 0000000..031f5eb
--- /dev/null
+++ b/docs/source/api/commands/extract.rst
@@ -0,0 +1,7 @@
+Extract
+=============================
+
+.. automodule:: silvimetric.commands.extract
+    :members:
+    :undoc-members:
+    :show-inheritance:
\ No newline at end of file
diff --git a/docs/source/api/commands/index.rst b/docs/source/api/commands/index.rst
new file mode 100644
index 0000000..6f491b7
--- /dev/null
+++ b/docs/source/api/commands/index.rst
@@ -0,0 +1,11 @@
+.. toctree::
+    :maxdepth: 1
+
+    initialize
+    scan
+    shatter
+    info
+    extract
+    delete
+    resume
+    restart
\ No newline at end of file
diff --git a/docs/source/api/commands/info.rst b/docs/source/api/commands/info.rst
new file mode 100644
index 0000000..4d21ef0
--- /dev/null
+++ b/docs/source/api/commands/info.rst
@@ -0,0 +1,7 @@
+Info
+===========================
+
+.. automodule:: silvimetric.commands.info
+   :members:
+   :undoc-members:
+   :show-inheritance:
\ No newline at end of file
diff --git a/docs/source/api/commands/initialize.rst b/docs/source/api/commands/initialize.rst
new file mode 100644
index 0000000..796c612
--- /dev/null
+++ b/docs/source/api/commands/initialize.rst
@@ -0,0 +1,7 @@
+Initialize
+========================
+
+.. automodule:: silvimetric.commands.initialize
+   :members:
+   :undoc-members:
+   :show-inheritance:
\ No newline at end of file
diff --git a/docs/source/api/commands/restart.rst b/docs/source/api/commands/restart.rst
new file mode 100644
index 0000000..585e6fa
--- /dev/null
+++ b/docs/source/api/commands/restart.rst
@@ -0,0 +1,4 @@
+Restart
+=======
+
+.. autofunction:: silvimetric.commands.manage.restart
\ No newline at end of file
diff --git a/docs/source/api/commands/resume.rst b/docs/source/api/commands/resume.rst
new file mode 100644
index 0000000..03a7aa3
--- /dev/null
+++ b/docs/source/api/commands/resume.rst
@@ -0,0 +1,4 @@
+Resume
+======
+
+.. autofunction:: silvimetric.commands.manage.resume
\ No newline at end of file
diff --git a/docs/source/api/commands/scan.rst b/docs/source/api/commands/scan.rst
new file mode 100644
index 0000000..7e7b4ed
--- /dev/null
+++ b/docs/source/api/commands/scan.rst
@@ -0,0 +1,7 @@
+Scan
+===========================
+
+.. automodule:: silvimetric.commands.scan
+   :members:
+   :undoc-members:
+   :show-inheritance:
\ No newline at end of file
diff --git a/docs/source/api/commands/shatter.rst b/docs/source/api/commands/shatter.rst
new file mode 100644
index 0000000..9804b84
--- /dev/null
+++ b/docs/source/api/commands/shatter.rst
@@ -0,0 +1,7 @@
+Shatter
+============================
+
+.. automodule:: silvimetric.commands.shatter
+   :members:
+   :undoc-members:
+   :show-inheritance:
\ No newline at end of file
diff --git a/docs/source/api.rst b/docs/source/api/index.rst
similarity index 55%
rename from docs/source/api.rst
rename to docs/source/api/index.rst
index 910969b..28962d2 100644
--- a/docs/source/api.rst
+++ b/docs/source/api/index.rst
@@ -4,7 +4,19 @@
 API
 ================================================================================
 
-.. autosummary::
-   :toctree: generated
 
-   silvimetric
+Resources
+---------
+
+.. toctree::
+    :maxdepth: 1
+
+    resources/index
+
+Commands
+--------
+
+.. toctree::
+    :maxdepth: 1
+
+    commands/index
\ No newline at end of file
diff --git a/docs/source/api/resources/bounds.rst b/docs/source/api/resources/bounds.rst
new file mode 100644
index 0000000..674efa2
--- /dev/null
+++ b/docs/source/api/resources/bounds.rst
@@ -0,0 +1,6 @@
+Bounds
+==========================
+
+.. automodule:: silvimetric.resources.bounds
+    :members:
+    :undoc-members:
\ No newline at end of file
diff --git a/docs/source/api/resources/config.rst b/docs/source/api/resources/config.rst
new file mode 100644
index 0000000..0210924
--- /dev/null
+++ b/docs/source/api/resources/config.rst
@@ -0,0 +1,5 @@
+Config
+==============================
+
+.. automodule:: silvimetric.resources.config
+   :members:
\ No newline at end of file
diff --git a/docs/source/api/resources/data.rst b/docs/source/api/resources/data.rst
new file mode 100644
index 0000000..7beb989
--- /dev/null
+++ b/docs/source/api/resources/data.rst
@@ -0,0 +1,7 @@
+Data
+=====================
+
+.. automodule:: silvimetric.resources.data
+   :members:
+   :undoc-members:
+   :show-inheritance:
\ No newline at end of file
diff --git a/docs/source/api/resources/entry.rst b/docs/source/api/resources/entry.rst
new file mode 100644
index 0000000..c9445ef
--- /dev/null
+++ b/docs/source/api/resources/entry.rst
@@ -0,0 +1,15 @@
+Attribute
+----------------------------------
+
+.. automodule:: silvimetric.resources.attribute
+   :members:
+   :undoc-members:
+   :show-inheritance:
+
+Metric
+----------------------------------
+
+.. automodule:: silvimetric.resources.metric
+   :members:
+   :undoc-members:
+   :show-inheritance:
diff --git a/docs/source/api/resources/extents.rst b/docs/source/api/resources/extents.rst
new file mode 100644
index 0000000..ae87329
--- /dev/null
+++ b/docs/source/api/resources/extents.rst
@@ -0,0 +1,5 @@
+Extents
+============================
+
+.. autoclass:: silvimetric.resources.extents.Extents
+    :members:
\ No newline at end of file
diff --git a/docs/source/api/resources/index.rst b/docs/source/api/resources/index.rst
new file mode 100644
index 0000000..8f01c8a
--- /dev/null
+++ b/docs/source/api/resources/index.rst
@@ -0,0 +1,10 @@
+.. toctree::
+    :maxdepth: 1
+
+    config
+    bounds
+    data
+    log
+    storage
+    entry
+    extents
\ No newline at end of file
diff --git a/docs/source/api/resources/log.rst b/docs/source/api/resources/log.rst
new file mode 100644
index 0000000..613e790
--- /dev/null
+++ b/docs/source/api/resources/log.rst
@@ -0,0 +1,7 @@
+Log
+======================
+
+.. automodule:: silvimetric.resources.log
+   :members:
+   :undoc-members:
+   :show-inheritance:
diff --git a/docs/source/api/resources/storage.rst b/docs/source/api/resources/storage.rst
new file mode 100644
index 0000000..f661c07
--- /dev/null
+++ b/docs/source/api/resources/storage.rst
@@ -0,0 +1,6 @@
+Storage
+------------------------------------
+
+.. automodule:: silvimetric.resources.storage
+    :members:
+    :undoc-members:
\ No newline at end of file
diff --git a/docs/source/cli/extract.rst b/docs/source/cli/extract.rst
index a4703a0..cf07872 100644
--- a/docs/source/cli/extract.rst
+++ b/docs/source/cli/extract.rst
@@ -16,16 +16,22 @@ Synopsis
 
 .. code-block::
 
-    Usage: silvimetric DATABASE extract [OPTIONS]
+    Usage: silvimetric [OPTIONS] extract [OPTIONS]
 
     Extract silvimetric metrics from DATABASE
 
     Options:
-    -a, --attributes ATTRS  List of attributes to include in Database
-    -m, --metrics METRICS   List of metrics to include in Database
-    --bounds BOUNDS
-    -o, --outdir TEXT       [required]
+    -a, --attributes ATTRS  List of attributes to include output
+    -m, --metrics METRICS   List of metrics to include in output
+    --bounds BOUNDS         Bounds for data to include in output
+    -o, --outdir PATH       Output directory.  [required]
     --help                  Show this message and exit.
 
+Example
+-------
+
+.. code-block::
+
+    silvimetric -d test.tdb extract -o test_tifs/
 
 .. include:: ../substitutions.txt
diff --git a/docs/source/cli/index.rst b/docs/source/cli/index.rst
index b35a57f..57f6f4e 100755
--- a/docs/source/cli/index.rst
+++ b/docs/source/cli/index.rst
@@ -4,8 +4,44 @@
 Command Line Interface
 ================================================================================
 
+
+.. only:: html
+
+    The command line interface (CLI) will facilitate your interaction with
+    Silvimetric from the terminal.
+
+.. code-block::
+
+    Usage: silvimetric [OPTIONS] COMMAND [ARGS]...
+
+    Options:
+    -d, --database PATH             Database path
+    --debug                         Changes logging level from INFO to DEBUG.
+    --log-dir TEXT                  Directory for log output
+    --progress BOOLEAN              Report progress
+    --workers INTEGER               Number of workers for Dask
+    --threads INTEGER               Number of threads per worker for Dask
+    --watch                         Open dask diagnostic page in default web
+                                    browser.
+    --dasktype [threads|processes]  What Dask uses for parallelization. For
+                                    moreinformation see here https://docs.dask.o
+                                    rg/en/stable/scheduling.html#local-threads
+    --scheduler [distributed|local|single-threaded]
+                                    Type of dask scheduler. Both are local, but
+                                    are run with different dask libraries. See
+                                    more here https://docs.dask.org/en/stable/sc
+                                    heduling.html.
+    --help                          Show this message and exit.
+
+    Commands:
+    extract     Extract silvimetric metrics from DATABASE
+    info        Retrieve information on current state of DATABASE
+    initialize  Create an empty DATABASE
+    scan        Scan point cloud and determine the optimal tile size.
+    shatter     Shatter point clouds into DATABASE cells.
+
 .. toctree::
-   :maxdepth: 2
+   :maxdepth: 1
    :caption: Commands
 
    extract
diff --git a/docs/source/cli/info.rst b/docs/source/cli/info.rst
index 4120556..d232041 100644
--- a/docs/source/cli/info.rst
+++ b/docs/source/cli/info.rst
@@ -6,7 +6,8 @@ info
 
 .. only:: html
 
-    Info provides detailed information about a SilviMetrics database
+    Info will query the database and output information about the `Shatter`
+    processes that been run as well as the database schema.
 
 .. Index:: info
 
@@ -16,13 +17,26 @@ Synopsis
 
 .. code-block::
 
-    Usage: silvimetric DATABASE info [OPTIONS]
+    Usage: silvimetric [OPTIONS] info [OPTIONS]
 
     Print info about Silvimetric database
 
     Options:
-    --history
-    --help     Show this message and exit.
+        --bounds BOUNDS                 Bounds to filter by
+        --date [%Y-%m-%d|%Y-%m-%dT%H:%M:%SZ]
+                                        Select processes with this date
+        --history                       Show the history section of the output.
+        --metadata                      Show the metadata section of the output.
+        --attributes                    Show the attributes section of the output.
+        --dates <DATETIME DATETIME>...  Select processes within this date range
+        --name TEXT                     Select processes with this name
+        --help                          Show this message and exit.
+
+Example
+-------
 
+.. code-block::
+
+    silvimetric -d test.tdb info
 
 .. include:: ../substitutions.txt
\ No newline at end of file
diff --git a/docs/source/cli/initialize.rst b/docs/source/cli/initialize.rst
index 97de933..4785829 100644
--- a/docs/source/cli/initialize.rst
+++ b/docs/source/cli/initialize.rst
@@ -20,15 +20,26 @@ Synopsis
 
 .. code-block::
 
-    Usage: silvimetric DATABASE initialize [OPTIONS] BOUNDS CRS
+    Usage: silvimetric [OPTIONS] initialize [OPTIONS]
 
     Initialize silvimetrics DATABASE
 
     Options:
+    --bounds BOUNDS         Root bounds that encapsulates all data  [required]
+    --crs CRS               Coordinate system of data  [required]
     -a, --attributes ATTRS  List of attributes to include in Database
     -m, --metrics METRICS   List of metrics to include in Database
     --resolution FLOAT      Summary pixel resolution
     --help                  Show this message and exit.
 
 
+Example
+-------
+
+.. code-block::
+
+    silvimetric --database test.tdb initialize --crs "EPSG:3857" \
+        --bounds '[300, 300, 600, 600]'
+
+
 .. include:: ../substitutions.txt
\ No newline at end of file
diff --git a/docs/source/cli/scan.rst b/docs/source/cli/scan.rst
index 3bf37ab..ae46644 100644
--- a/docs/source/cli/scan.rst
+++ b/docs/source/cli/scan.rst
@@ -6,7 +6,8 @@ scan
 
 .. only:: html
 
-    Scan provides detailed information about a SilviMetrics database
+    This will inspect a point cloud file with the database as a reference and
+    supply an estimate of what tile size you could use for `Shatter` operations.
 
 .. Index:: scan
 
@@ -16,8 +17,24 @@ Synopsis
 
 .. code-block::
 
-    Usage: silvimetric DATABASE scan [OPTIONS]
+    silvimetric [OPTIONS]  scan [OPTIONS] POINTCLOUD
 
+    Options:
+    --resolution FLOAT     Summary pixel resolution
+    --filter               Remove empty space in computation. Will take extra
+                            time.
+    --point_count INTEGER  Point count threshold.
+    --depth INTEGER        Quadtree depth threshold.
+    --bounds BOUNDS        Bounds to scan.
+    --help                 Show this message and exit.
+
+
+Usage
+-----
+
+.. code-block::
+
+    silvimetric --database test.tdb scan tests/data/test_data.copc.laz
 
 
 .. include:: ../substitutions.txt
diff --git a/docs/source/cli/shatter.rst b/docs/source/cli/shatter.rst
index c50e9c9..fbd3662 100644
--- a/docs/source/cli/shatter.rst
+++ b/docs/source/cli/shatter.rst
@@ -15,19 +15,27 @@ Synopsis
 
 .. code-block::
 
-    Usage: silvimetric DATABASE shatter [OPTIONS] POINTCLOUD
+    Usage: silvimetric [OPTIONS] shatter [OPTIONS] POINTCLOUD
 
     Insert data provided by POINTCLOUD into the silvimetric DATABASE
 
     Options:
-    --workers INTEGER
-    --bounds BOUNDS
-    --tilesize INTEGER
-    --threads INTEGER
-    --watch
-    --dasktype [cluster|threads|processes|single-threaded]
+    --bounds BOUNDS                 Bounds for data to include in processing
+    --tilesize INTEGER              Number of cells to include per tile
+    --report                        Whether or not to write a report of the
+                                    process, useful for debugging
+    --date [%Y-%m-%d|%Y-%m-%dT%H:%M:%SZ]
+                                    Date the data was produced.
+    --dates <DATETIME DATETIME>...  Date range the data was produced during
     --help                          Show this message and exit.
 
 
+Example
+-------
+
+.. code-block::
+
+    silvimetric -d test.tdb shatter --date 2023-1-1 tests/data/test_data.copc.laz
+
 
 .. include:: ../substitutions.txt
\ No newline at end of file
diff --git a/docs/source/conf.py b/docs/source/conf.py
index 8aea005..2e72f22 100644
--- a/docs/source/conf.py
+++ b/docs/source/conf.py
@@ -10,6 +10,7 @@
 import sys, os, re
 import time
 import datetime
+
 if os.environ.get('SOURCE_DATE_EPOCH'):
     year  = datetime.datetime.utcfromtimestamp(int(os.environ.get('SOURCE_DATE_EPOCH', time.gmtime()))).year
     today = datetime.datetime.utcfromtimestamp(int(os.environ.get('SOURCE_DATE_EPOCH', time.gmtime()))).strftime('%B %d, %Y')
@@ -75,3 +76,20 @@ def read_version(filename):
    'sphinx.ext.autodoc',
    'sphinx.ext.autosummary',
 ]
+
+autodoc_default_options = {
+    'member-order': 'groupwise'
+}
+
+autodoc_mock_imports = [
+    "pyproj",
+    "osgeo",
+    "distributed",
+    "dask.dataframe",
+    "scipy",
+    "dill",
+    "shapely",
+    "pdal",
+    "gdal",
+    "pandas"
+]
\ No newline at end of file
diff --git a/docs/source/index.rst b/docs/source/index.rst
index df3af13..cd59787 100644
--- a/docs/source/index.rst
+++ b/docs/source/index.rst
@@ -15,12 +15,12 @@ Find out more about SilviMetric by visiting :ref:`about`. A slide deck about
 SilviMetric is also available on `Google Slides <https://docs.google.com/presentation/d/1E561EgWwLgN5R9P0LBxuI1r7kG155u8E6-MOWpkycSM/edit?usp=sharing>`__.
 
 .. toctree::
-   :maxdepth: 2
    :caption: Contents
+   :maxdepth: 2
 
    about
    cli/index
-   api
+   api/index
    QuickStart <quickstart>
    development
    tutorial
diff --git a/docs/source/quickstart.rst b/docs/source/quickstart.rst
index 1ad7260..9b87a19 100644
--- a/docs/source/quickstart.rst
+++ b/docs/source/quickstart.rst
@@ -67,7 +67,7 @@ need a bounds and a coordinate reference system.
     .. code-block:: shell-session
 
         pdal info https://s3.amazonaws.com/hobu-lidar/autzen-classified.copc.laz  \
-            --readers.copc.resolution=10 | jq -c '.stats.bbox.native.bbox'
+            --readers.copc.resolution=1 | jq -c '.stats.bbox.native.bbox'
 
     Our boundary is emitted in expanded form.
 
@@ -124,26 +124,125 @@ need a bounds and a coordinate reference system.
     Be careful with your shell's quote escaping rules!
 
 
+Scan
+--------------------------------------------------------------------------------
+
+The `scan` command will tell us information about the pointcloud with respect
+to the database we already created, including a best guess at the correct number
+of cells per tile, or `tile size`.
+
+    .. code-block:: shell-session
+
+        silvimetric -d ${db_name} scan ${pointcloud}
+
+
+We should see output like the output below, recommending we use a `tile size`
+of 185.
+
+::
+
+    silvimetric - INFO - info:156 - Pointcloud information:
+    silvimetric - INFO - info:156 -   Storage Bounds: [635579.2, 848884.83, 639003.73, 853536.21]
+    silvimetric - INFO - info:156 -   Pointcloud Bounds: [635577.79, 848882.15, 639003.73, 853537.66]
+    silvimetric - INFO - info:156 -   Point Count: 10653336
+    silvimetric - INFO - info:156 - Tiling information:
+    silvimetric - INFO - info:156 -   Mean tile size: 91.51758793969849
+    silvimetric - INFO - info:156 -   Std deviation: 94.31396536316173
+    silvimetric - INFO - info:156 -   Recommended split size: 185
+
 Shatter
 --------------------------------------------------------------------------------
 
-We can now insert data into the SMDB
+We can now insert data into the SMDB.
+
+If we run this command without the argument `--tilesize`, `Silvimetric` will
+determine a tile size for you. The method will be the same as the `Scan` method,
+but will filter out the tiles that have no data in them.
 
    .. code-block:: shell-session
 
-     silvimetric autzen-smdb.tdb \
+     silvimetric -d autzen-smdb.tdb \
+        --threads 4 \
+        --workers 4 \
+        --watch \
         shatter \
-        https://s3.amazonaws.com/hobu-lidar/autzen-classified.copc.laz \
-        --threads 10 \
-        --watch
+        --date 2008-12-01 \
+        https://s3.amazonaws.com/hobu-lidar/autzen-classified.copc.laz
+
+If we grab the tile size from the `scan` that we ran earlier, we'll skip the
+filtering step.
+
+   .. code-block:: shell-session
+
+     silvimetric -d autzen-smdb.tdb \
+        --threads 4 \
+        --workers 4 \
+        --watch \
+        shatter \
+        --tilesize 185 \
+        --date 2008-12-01 \
+        https://s3.amazonaws.com/hobu-lidar/autzen-classified.copc.laz
 
 Extract
 --------------------------------------------------------------------------------
 
-After data is inserted, we can extract
+After data is inserted, we can extract it into different rasters. When we created
+the database we gave it a list of `Attributes` and `Metrics`. When we ran
+`Shatter`, we filled in the values for those in each cell. If we have a database
+with the `Attributes` Intensity and Z, in combination with the `Metrics` min and
+max, each cell will contain values for `min_Intensity`, `max_Intensity`,
+`min_Z`, and `max_Z`. This is also the list of available rasters we can extract.
+
+   .. code-block:: shell-session
+
+     silvimetric -d autzen-smdb.tdb extract -o output-directory
+
+Info
+--------------------------------------------------------------------------------
+
+We can query past shatter processes and the schema for the database with the
+Info call.
+
+    .. code-block:: shell-session
+
+        silvimetric -d autzen-smdb.tdb info --history
+
+This will print out a JSON object containing information about the current state
+of the database. We can find the `name` key here, which necessary for
+`Delete`, `Restart`, and `Resume`. For the following commands we will have
+copied the value of the `name` key in the variable `uuid`.
+
+Delete
+--------------------------------------------------------------------------------
+
+We can also remove a `shatter` process by using the `delete` command. This will
+remove all data associated with that shatter process from the database, but
+will leave an updated config of it in the database config should you want to
+reference it later.
+
+    .. code-block:: shell-session
+
+        silvimetric -d autzen-smdb.tdb delete --id $uuid
+
+Restart
+--------------------------------------------------------------------------------
+
+If you would like to rerun a `Shatter` process, whether or not it was previously
+finished, you can use the `restart` command. This will call the `delete` method
+and use the config from that to re-run the `shatter` process.
+
+    .. code-block:: shell-session
+
+        silvimetric -d autzen-smdb.tdb restart --id $uuid
+
+Resume
+--------------------------------------------------------------------------------
+
+If a `Shatter` process is cancelled partway through, we can pick up where we
+left off with the `Resume` command.
 
    .. code-block:: shell-session
 
-     silvimetric autzen-smdb.tdb extract -o output-directory
+        silvimetric -d autzen-smdb.tdb resume --id $uuid
 
 .. include:: ./substitutions.txt
\ No newline at end of file
diff --git a/docs/source/substitutions.txt b/docs/source/substitutions.txt
index 52763fc..09e60ea 100644
--- a/docs/source/substitutions.txt
+++ b/docs/source/substitutions.txt
@@ -13,3 +13,4 @@
 .. |CondaForge| replace:: `Conda Forge <https://conda-forge.org/docs/>`__
 .. |S3| replace:: `S3 <https://aws.amazon.com/pm/serv-s3/>`__
 .. |FUSION| replace:: `FUSION <http://forsys.sefs.uw.edu/fusion/fusionlatest.html>`__
+.. |Click| replace:: `Click <https://click.palletsprojects.com/en/8.1.x/>`__
diff --git a/docs/source/tutorial.rst b/docs/source/tutorial.rst
index 3769847..435a500 100644
--- a/docs/source/tutorial.rst
+++ b/docs/source/tutorial.rst
@@ -189,7 +189,7 @@ Example:
     $ METRIC_PATH="./path/to/python_metrics.py"
     $ silvimetric --database $DB_NAME initialize --bounds "$BOUNDS" \
           --crs "EPSG:$EPSG" \
-          -m $METRIC_PATH -m min -m max -m mean
+          -m "${METRIC_PATH},min,max,mean"
 
 .. warning::
 
diff --git a/environment.yml b/environment.yml
index 6674d0a..62257ac 100644
--- a/environment.yml
+++ b/environment.yml
@@ -3,7 +3,9 @@ channels:
   - conda-forge
 dependencies:
   - python
+  - pip
   - tiledb-py
+  - pdal
   - python-pdal
   - numpy
   - shapely
@@ -15,5 +17,5 @@ dependencies:
   - websocket-client
   - python-json-logger
   - dill
-  - pip:
-    - line-profiler
+  - pandas
+  - lmoments3
\ No newline at end of file
diff --git a/src/silvimetric/__init__.py b/src/silvimetric/__init__.py
index 1f356cc..23c8d91 100644
--- a/src/silvimetric/__init__.py
+++ b/src/silvimetric/__init__.py
@@ -1,1 +1,20 @@
-__version__ = '1.0.0'
+__version__ = '1.2.1'
+
+from .resources.bounds import Bounds
+from .resources.extents import Extents
+from .resources.storage import Storage
+from .resources.metric import Metric, run_metrics
+from .resources.metrics import grid_metrics, l_moments, percentiles, statistics, all_metrics
+from .resources.metrics import product_moments
+from .resources.log import Log
+from .resources.data import Data
+from .resources.attribute import Attribute, Pdal_Attributes, Attributes
+from .resources.config import StorageConfig, ShatterConfig, ExtractConfig
+from .resources.config import ApplicationConfig
+
+from .commands.shatter import shatter
+from .commands.extract import extract
+from .commands.info import info
+from .commands.scan import scan
+from .commands.initialize import initialize
+from .commands.manage import delete, resume, restart
diff --git a/src/silvimetric/cli/cli.py b/src/silvimetric/cli/cli.py
index 7cd246c..19d1c8b 100644
--- a/src/silvimetric/cli/cli.py
+++ b/src/silvimetric/cli/cli.py
@@ -4,35 +4,40 @@
 import logging
 
 
-from ..resources import Bounds, Log
-from ..resources import Attribute, Metric
-from ..resources import StorageConfig, ShatterConfig, ExtractConfig, ApplicationConfig
-from ..commands import shatter, extract, scan, info, initialize
+from .. import __version__
+from .. import Attribute, Metric, Bounds, Log
+from .. import StorageConfig, ShatterConfig, ExtractConfig, ApplicationConfig
+from ..commands import shatter, extract, scan, info, initialize, manage
 from .common import BoundsParamType, CRSParamType, AttrParamType, MetricParamType
-from .common import dask_handle
+from .common import dask_handle, close_dask
 
 @click.group()
 @click.option("--database", '-d', type=click.Path(exists=False), help="Database path")
-@click.option("--debug", is_flag=True, default=False, help="Print debug messages?")
-@click.option("--log-level", default="INFO", help="Log level (INFO/DEBUG)")
+@click.option("--debug", is_flag=True, default=False, help="Changes logging level from INFO to DEBUG.")
 @click.option("--log-dir", default=None, help="Directory for log output", type=str)
-@click.option("--progress", default=True, type=bool, help="Report progress")
-@click.option("--workers", type=int, default=12, help="Number of workers for Dask")
-@click.option("--threads", type=int, default=4, help="Number of threads per worker for Dask")
+@click.option("--progress", is_flag=True, default=True, type=bool, help="Report progress")
+@click.option("--workers", type=int, help="Number of workers for Dask")
+@click.option("--threads", type=int, help="Number of threads per worker for Dask")
 @click.option("--watch", is_flag=True, default=False, type=bool,
         help="Open dask diagnostic page in default web browser.")
 @click.option("--dasktype", default='processes', type=click.Choice(['threads',
         'processes']), help="What Dask uses for parallelization. For more"
         "information see here https://docs.dask.org/en/stable/scheduling.html#local-threads")
-@click.option("--scheduler", default='distributed', type=click.Choice(['distributed',
+@click.option("--scheduler", default='local', type=click.Choice(['distributed',
         'local', 'single-threaded']), help="Type of dask scheduler. Both are "
         "local, but are run with different dask libraries. See more here "
         "https://docs.dask.org/en/stable/scheduling.html.")
+@click.version_option(__version__)
 @click.pass_context
-def cli(ctx, database, debug, log_level, log_dir, progress, dasktype, scheduler,
+def cli(ctx, database, debug, log_dir, progress, dasktype, scheduler,
         workers, threads, watch):
 
     # Set up logging
+    if debug:
+        log_level = 'DEBUG'
+    else:
+        log_level = 'INFO'
+
     numeric_level = getattr(logging, log_level.upper(), None)
     if not isinstance(numeric_level, int):
         raise ValueError(f"Invalid log level: {log_level}")
@@ -42,9 +47,13 @@ def cli(ctx, database, debug, log_level, log_dir, progress, dasktype, scheduler,
             log=log,
             debug=debug,
             progress=progress,
-            scheduler=scheduler)
-    dask_handle(dasktype, scheduler, workers, threads, watch, app.log)
+            scheduler=scheduler,
+            dasktype=dasktype,
+            workers=workers,
+            threads=threads,
+            watch=watch)
     ctx.obj = app
+    ctx.call_on_close(close_dask)
 
 
 @cli.command("info")
@@ -52,6 +61,12 @@ def cli(ctx, database, debug, log_level, log_dir, progress, dasktype, scheduler,
         help="Bounds to filter by")
 @click.option("--date", type=click.DateTime(['%Y-%m-%d','%Y-%m-%dT%H:%M:%SZ']),
         help="Select processes with this date")
+@click.option("--history", type=bool, is_flag=True, default=False,
+        help="Show the history section of the output.")
+@click.option("--metadata", type=bool, is_flag=True, default=False,
+        help="Show the metadata section of the output.")
+@click.option("--attributes", type=bool, is_flag=True, default=False,
+        help="Show the attributes section of the output.")
 @click.option("--dates", type=click.Tuple([
         click.DateTime(['%Y-%m-%d','%Y-%m-%dT%H:%M:%SZ']),
         click.DateTime(['%Y-%m-%d','%Y-%m-%dT%H:%M:%SZ'])]), nargs=2,
@@ -59,17 +74,34 @@ def cli(ctx, database, debug, log_level, log_dir, progress, dasktype, scheduler,
 @click.option("--name", type=str, default=None,
         help="Select processes with this name")
 @click.pass_obj
-def info_cmd(app, bounds, date, dates, name):
+def info_cmd(app, bounds, date, dates, name, history, metadata, attributes):
     import json
     if date is not None and dates is not None:
-        app.log.warning("Both 'date' and 'dates' specified. Prioritizing 'dates'")
+        app.log.warning("Both 'date' and 'dates' specified. Prioritizing"
+            "'dates'")
 
     start_date = dates[0] if dates else date
     end_date = dates[1] if dates else date
 
     i = info.info(app.tdb_dir, bounds=bounds, start_time=start_date,
         end_time=end_date, name=name)
-    app.log.info(json.dumps(i, indent=2))
+
+    if any([history, metadata, attributes]):
+        filtered = { }
+        if history:
+            filtered['history'] = i['history']
+        if metadata:
+            filtered['metadata'] = i['metadata']
+        if attributes:
+            filtered['attributes'] = i['attributes']
+
+        app.log.info(json.dumps(filtered, indent=2))
+
+    else:
+        app.log.info(json.dumps(i, indent=2))
+        return
+
+
 
 
 @cli.command("scan")
@@ -86,9 +118,12 @@ def info_cmd(app, bounds, date, dates, name):
         help="Bounds to scan.")
 @click.pass_obj
 def scan_cmd(app, resolution, point_count, pointcloud, bounds, depth, filter):
-    """Scan point cloud and determine the optimal tile size."""
+    """Scan point cloud, output information on it, and determine the optimal
+    tile size."""
+    dask_handle(app.dasktype, app.scheduler, app.workers, app.threads,
+            app.watch, app.log)
     return scan.scan(app.tdb_dir, pointcloud, bounds, point_count, resolution,
-            depth, filter)
+            depth, filter, log=app.log)
 
 
 @cli.command('initialize')
@@ -98,22 +133,21 @@ def scan_cmd(app, resolution, point_count, pointcloud, bounds, depth, filter):
         help="Coordinate system of data")
 @click.option("--attributes", "-a", multiple=True, type=AttrParamType(),
         help="List of attributes to include in Database")
-@click.option("--metrics", "-m", multiple=True, type=MetricParamType(),
-        help="List of metrics to include in Database")
+@click.option("--metrics", "-m", type=MetricParamType(), default=[],
+        help="List of metrics to include in output, eg. '-m stats,percentiles'")
 @click.option("--resolution", type=float, default=30.0,
         help="Summary pixel resolution")
 @click.pass_obj
 def initialize_cmd(app: ApplicationConfig, bounds: Bounds, crs: pyproj.CRS,
         attributes: list[Attribute], resolution: float, metrics: list[Metric]):
-    import itertools
-    """Initialize silvimetrics DATABASE
-    """
+    """Initialize silvimetrics DATABASE"""
+
     storageconfig = StorageConfig(tdb_dir = app.tdb_dir,
             log = app.log,
             root = bounds,
             crs = crs,
             attrs = attributes,
-            metrics = list(itertools.chain(*metrics)),
+            metrics = metrics,
             resolution = resolution)
     return initialize.initialize(storageconfig)
 
@@ -134,6 +168,10 @@ def initialize_cmd(app: ApplicationConfig, bounds: Bounds, crs: pyproj.CRS,
         help="Date range the data was produced during")
 @click.pass_obj
 def shatter_cmd(app, pointcloud, bounds, report, tilesize, date, dates):
+    """Insert data provided by POINTCLOUD into the silvimetric DATABASE"""
+
+    dask_handle(app.dasktype, app.scheduler, app.workers, app.threads,
+            app.watch, app.log)
 
     if date is not None and dates is not None:
         app.log.warning("Both 'date' and 'dates' specified. Prioritizing 'dates'")
@@ -141,31 +179,33 @@ def shatter_cmd(app, pointcloud, bounds, report, tilesize, date, dates):
     if date is None and dates is None:
         raise ValueError("One of '--date' or '--dates' must be provided.")
 
-    """Insert data provided by POINTCLOUD into the silvimetric DATABASE"""
-    config = ShatterConfig(tdb_dir = app.tdb_dir,
+    config = ShatterConfig(tdb_dir=app.tdb_dir,
             date=dates if dates else tuple([date]),
-            log = app.log,
-            filename = pointcloud,
-            bounds = bounds,
+            log=app.log,
+            filename=pointcloud,
+            bounds=bounds,
             tile_size=tilesize)
 
     if report:
         if app.scheduler != 'distributed':
+            from dask.diagnostics import ProgressBar
             app.log.warning('Report option is incompatible with scheduler'
                             '{scheduler}, skipping.')
-        report_path = f'reports/{config.name}.html'
-        with performance_report(report_path) as pr:
             shatter.shatter(config)
-        print(f'Writing report to {report_path}.')
+        else:
+            report_path = f'reports/{config.name}.html'
+            with performance_report(report_path) as pr:
+                shatter.shatter(config)
+            print(f'Writing report to {report_path}.')
     else:
         shatter.shatter(config)
 
 
 @cli.command('extract')
-@click.option("--attributes", "-a", multiple=True, type=AttrParamType(),
+@click.option("--attributes", "-a", multiple=True, type=AttrParamType(), default=[],
         help="List of attributes to include output")
-@click.option("--metrics", "-m", multiple=True, type=MetricParamType(),
-        help="List of metrics to include in output")
+@click.option("--metrics", "-m", type=MetricParamType(), default=[],
+        help="List of metrics to include in output, eg. '-m stats,percentiles'")
 @click.option("--bounds", type=BoundsParamType(), default=None,
         help="Bounds for data to include in output")
 @click.option("--outdir", "-o", type=click.Path(exists=False), required=True,
@@ -176,6 +216,8 @@ def extract_cmd(app, attributes, metrics, outdir, bounds):
 
     #TODO only allow metrics and attributes to be added if they're present
     # in the storage config.
+    dask_handle(app.dasktype, app.scheduler, app.workers, app.threads,
+        app.watch, app.log)
 
     config = ExtractConfig(tdb_dir = app.tdb_dir,
             log = app.log,
@@ -185,6 +227,26 @@ def extract_cmd(app, attributes, metrics, outdir, bounds):
             bounds = bounds)
     extract.extract(config)
 
+@cli.command('delete')
+@click.option('--id', type=click.UUID, required=True,
+    help="Shatter Task UUID.")
+@click.pass_obj
+def delete_cmd(app, id):
+    manage.delete(app.tdb_dir, id)
+
+@cli.command('restart')
+@click.option('--id', type=click.UUID, required=True,
+    help="Shatter Task UUID.")
+@click.pass_obj
+def restart_cmd(app, id):
+    manage.restart(tdb_dir=app.tdb_dir, name=id)
+
+@cli.command('resume')
+@click.option('--id', type=click.UUID, required=True,
+    help="Shatter Task UUID.")
+@click.pass_obj
+def resume_cmd(app, id):
+    manage.resume(tdb_dir=app.tdb_dir, name=id)
 
 if __name__ == "__main__":
-    cli()
+    cli()
\ No newline at end of file
diff --git a/src/silvimetric/cli/common.py b/src/silvimetric/cli/common.py
index a247d0c..a21cf62 100644
--- a/src/silvimetric/cli/common.py
+++ b/src/silvimetric/cli/common.py
@@ -5,8 +5,10 @@
 import dask
 from dask.diagnostics import ProgressBar
 from dask.distributed import Client, LocalCluster
+from ..resources.metrics import l_moments, percentiles, statistics, product_moments
+from ..resources.metrics import aad, grid_metrics, all_metrics
 
-from ..resources import Bounds, Attribute, Metric, Attributes, Metrics, Log
+from .. import Bounds, Attribute, Metric, Attributes, Log
 
 
 class BoundsParamType(click.ParamType):
@@ -22,7 +24,7 @@ def convert(self, value, param, ctx):
 class CRSParamType(click.ParamType):
     name = "CRS"
 
-    def convert(self, value, param, ctx):
+    def convert(self, value, param, ctx) -> pyproj.CRS:
         try:
             crs = pyproj.CRS.from_user_input(value)
             return crs
@@ -33,46 +35,79 @@ class AttrParamType(click.ParamType):
     name="Attrs"
     #TODO add import similar to metrics
     def convert(self, value, param, ctx) -> list[Attribute]:
-        try:
-            return [Attributes[a] for a in value]
-        except Exception as e:
-            self.fail(f"{value!r} is not available in Attributes, {e}", param, ctx)
-
-class MetricParamType(click.ParamType):
-    name="Metrics"
-    def convert(self, value, param, ctx) -> list[Metric]:
-        if '.py' in value:
+        if isinstance(value, list):
             try:
-                import importlib.util
-                import os
-                from pathlib import Path
-
-                cwd = os.getcwd()
-                p = Path(cwd, value)
-                if not p.exists():
-                    self.fail("Failed to find import file for metrics at"
-                            f" {str(p)}", param, ctx)
-
-                spec = importlib.util.spec_from_file_location('user_metrics', str(p))
-                user_metrics = importlib.util.module_from_spec(spec)
-                spec.loader.exec_module(user_metrics)
-                ms = user_metrics.metrics()
+                return [Attributes[a] for a in value]
             except Exception as e:
-                self.fail(f"Failed to import metrics from {str(p)} with error {e}",
-                        param, ctx)
-
-            for m in ms:
-                if not isinstance(m, Metric):
-                    self.fail(f"Invalid Metric supplied: {m}")
-            return user_metrics.metrics()
+                self.fail(f"{value!r} is not available in Attributes, {e}", param, ctx)
+        elif isinstance(value, str):
+            return Attributes[value]
+        else:
+            self.fail(f"{value!r} is of an invalid type, {e}", param, ctx)
 
-        try:
-            return [ Metrics[value] ]
-        except Exception as e:
-            self.fail(f"{value!r} is not available in Metrics, {e}", param, ctx)
+class MetricParamType(click.ParamType):
+    name="metrics"
+    def convert(self, value, param, ctx) -> list[Metric]:
+        if value is None or not value:
+            return list(all_metrics.values())
+        parsed_values = value.split(',')
+        metrics: set[Metric] = set()
+        for val in parsed_values:
+            if '.py' in val:
+                # user imported metrics from external file
+                try:
+                    import importlib.util
+                    import os
+                    from pathlib import Path
+
+                    cwd = os.getcwd()
+                    p = Path(cwd, val)
+                    if not p.exists():
+                        self.fail("Failed to find import file for metrics at"
+                                f" {str(p)}", param, ctx)
+
+                    spec = importlib.util.spec_from_file_location('user_metrics', str(p))
+                    user_metrics = importlib.util.module_from_spec(spec)
+                    spec.loader.exec_module(user_metrics)
+                    ms = user_metrics.metrics()
+                except Exception as e:
+                    self.fail(f"Failed to import metrics from {str(p)} with error {e}",
+                            param, ctx)
+
+                for m in ms:
+                    if not isinstance(m, Metric):
+                        self.fail(f"Invalid Metric supplied: {m}")
+
+                metrics.update(list(user_metrics.metrics()))
+            else:
+                # SilviMetric defined metrics
+                try:
+                    if val == 'stats':
+                        metrics.update(list(statistics.values()))
+                    elif val == 'p_moments':
+                        metrics.update(list(product_moments.values()))
+                    elif val == 'l_moments':
+                        metrics.update(list(l_moments.values()))
+                    elif val == 'percentiles':
+                        metrics.update(list(percentiles.values()))
+                    elif val == 'aad':
+                        metrics.update(list(aad.aad.values()))
+                    elif val == 'grid_metrics':
+                        metrics.update(list(grid_metrics.values()))
+                    elif val == 'all':
+                        metrics.update(list(all_metrics.values()))
+                    else:
+                        m = all_metrics[val]
+                        if isinstance(m, Metric):
+                            metrics.add(m)
+                        else:
+                            metrics.udpate(list(m))
+                except Exception as e:
+                    self.fail(f"{val!r} is not available in Metrics", param, ctx)
+        return list(metrics)
 
 def dask_handle(dasktype: str, scheduler: str, workers: int, threads: int,
-        watch: bool, log: Log):
+        watch: bool, log: Log) -> None:
     dask_config = { }
 
     if dasktype == 'threads':
@@ -83,9 +118,6 @@ def dask_handle(dasktype: str, scheduler: str, workers: int, threads: int,
         dask_config['threads_per_worker'] = threads
 
     if scheduler == 'local':
-        if scheduler != 'distributed':
-            log.warning("Selected scheduler type does not support continuously"
-                            "updated config information.")
         # fall back to dask type to determine the scheduler type
         dask_config['scheduler'] = dasktype
         if watch:
@@ -108,9 +140,11 @@ def dask_handle(dasktype: str, scheduler: str, workers: int, threads: int,
             webbrowser.open(client.cluster.dashboard_link)
 
     elif scheduler == 'single-threaded':
-        if scheduler != 'distributed':
-            log.warning("""Selected scheduler type does not support continuously\
-                            updated config information.""")
         dask_config['scheduler'] = scheduler
 
     dask.config.set(dask_config)
+
+def close_dask() -> None:
+    client = dask.config.get('distributed.client')
+    if isinstance(client, Client):
+        client.close()
\ No newline at end of file
diff --git a/src/silvimetric/commands/extract.py b/src/silvimetric/commands/extract.py
index ceb047c..04dafeb 100644
--- a/src/silvimetric/commands/extract.py
+++ b/src/silvimetric/commands/extract.py
@@ -1,13 +1,15 @@
+from pathlib import Path
+from itertools import chain
+
+
+from typing import Union
 from osgeo import gdal, osr
+import dask
 import numpy as np
-from pathlib import Path
 import pandas as pd
-import dask
-import dask.dataframe as dd
-import dask.bag as db
 
-from ..resources import Storage, ExtractConfig, Metric, Attribute
-from ..resources import Extents
+
+from .. import Storage, Extents, ExtractConfig
 
 np_to_gdal_types = {
     np.dtype(np.byte).str: gdal.GDT_Byte,
@@ -23,17 +25,25 @@
 }
 
 def write_tif(xsize: int, ysize: int, data:np.ndarray, name: str,
-              config: ExtractConfig):
+        config: ExtractConfig) -> None:
+    """
+    Write out a raster with GDAL
+
+    :param xsize: Length of X plane.
+    :param ysize: Length of Y plane.
+    :param data: Data to write to raster.
+    :param name: Name of raster to write.
+    :param config: ExtractConfig.
+    """
     osr.UseExceptions()
     path = Path(config.out_dir) / f'{name}.tif'
     crs = config.crs
     srs = osr.SpatialReference()
     srs.ImportFromWkt(crs.to_wkt())
-    # transform = [x, res, 0, y, 0, res]
     b = config.bounds
 
     transform = [b.minx, config.resolution, 0,
-                 b.maxy, 0, -1* config.resolution]
+                 b.maxy, 0, -1*config.resolution]
 
     driver = gdal.GetDriverByName("GTiff")
     gdal_type = np_to_gdal_types[np.dtype(data.dtype).str]
@@ -45,31 +55,46 @@ def write_tif(xsize: int, ysize: int, data:np.ndarray, name: str,
     tif.FlushCache()
     tif = None
 
-def get_metrics(data_in, config: ExtractConfig):
+def get_metrics(data_in: pd.DataFrame, storage: Storage) -> Union[None, pd.DataFrame]:
+    """
+    Reruns a metric over this cell. Only called if there is overlapping data.
+
+    :param data_in: Dataframe to be rerun.
+    :param storage: Base storage object.
+    :return: Combined dict of attribute and newly derived metric data.
+    """
+
+    #TODO should just use the metric calculation methods from shatter
     if data_in is None:
         return None
 
-    # make sure it's not empty. No empty writes
-    if not np.any(data_in['count']):
-        return None
+    def expl(x):
+        return x.explode()
+
+    attrs = [a.name for a in storage.config.attrs if a.name not in ['X','Y']]
+
+    # set index so we can apply to the whole dataset without needing to skip X and Y
+    # then reset in the index because that's what metric.do expects
+    exploded = data_in.set_index(['X','Y']).apply(expl)[attrs].reset_index()
+    metric_data = dask.persist(*[ m.do(exploded) for m in storage.config.metrics ])
 
-    # doing dask compute inside the dict array because it was too fine-grained
-    # when it was outside
-    #TODO move this operation to a dask bag and convert to dataframe and then
-    # to a dict for better efficiency?
-    metric_data = {
-        f'{m.entry_name(attr.name)}': [m(cell_data) for cell_data in data_in[attr.name]]
-        for attr in config.attrs for m in config.metrics
-    }
-    # md = dask.compute(metric_data)[0]
-    data_out = data_in | metric_data
+    data_out = data_in.set_index(['X','Y']).join([m for m in metric_data])
     return data_out
 
-def handle_overlaps(config: ExtractConfig, storage: Storage, indices: np.ndarray):
-    ma_list = [ m.entry_name(a.name) for m in config.metrics for a in
-        config.attrs ]
-    att_list = [ a.name for a in config.attrs ] + [ 'count' ]
-    out_list = [ *ma_list, 'X', 'Y' ]
+def handle_overlaps(config: ExtractConfig, storage: Storage, indices: np.ndarray) -> pd.DataFrame:
+    """
+    Handle cells that have overlapping data. We have to re-run metrics over these
+    cells as there's no other accurate way to determined metric values. If there
+    are no overlaps, this will do nothing.
+
+    :param config: ExtractConfig.
+    :param storage: Database storage object.
+    :param indices: Indices with overlap.
+    :return: Dataframe of rerun data.
+    """
+
+    ma_list = storage.getDerivedNames()
+    att_list = [ a.name for a in config.attrs ]
 
     minx = indices['x'].min()
     maxx = indices['x'].max()
@@ -84,81 +109,63 @@ def handle_overlaps(config: ExtractConfig, storage: Storage, indices: np.ndarray
         att_meta[a.name] =  a.dtype
 
     with storage.open("r") as tdb:
-        data = tdb.query(attrs=ma_list, order='F', coords=True).df[minx:maxx,
-                miny:maxy]
-        data['X'] = data['X'] - data['X'].min()
-        data['Y'] = data['Y'] - data['Y'].min()
-
-        # 1. should find values that are not unique, meaning they have multiple
-        # entries
-        recarray = data.to_records()
-        xys = recarray[['X', 'Y']]
-        unq, idx, counts = np.unique(xys, return_index=True, return_counts=True)
-        redos = np.where(counts >= 2)
-
-        if not np.any(redos):
-            return data
-
-        leaves = data.loc[idx[np.where(counts == 1)]]
-
-        # 2. then should combine the values of those attribute/cell combos
-        ridx = np.unique(unq[redos])
-        rxs = list(ridx['X'])
-        rys = list(ridx['Y'])
-        if np.any(rxs):
-            redo = pd.DataFrame(tdb.query(attrs=att_list).multi_index[[rxs], [rys]])
-            redo['X'] = redo['X'] - minx
-            redo['Y'] = redo['Y'] - miny
-
-        parts = min(int(redo['X'].max() * redo['Y'].max()), 1000)
-        r = dd.from_pandas(redo, npartitions=parts, name='MetricDataFrame')
-        def squash(x):
-            def s2(vals):
-                if vals.dtype == object:
-                    return np.array([*vals], vals.dtype).reshape(-1)
-                elif vals.name == 'X' or vals.name == 'Y':
-                    return vals.values[0]
-                else:
-                    return sum(vals)
-            val = x.aggregate(s2)
-            p = pd.DataFrame(val).T
-            return p
-
-        recs = r.groupby(['X','Y']).apply(squash, meta=att_meta)
-        rs = recs.compute().to_dict(orient='list')
-
-        # 3. Should rerun the metrics over them
-        metrics = get_metrics(rs, config)
-        ms = pd.DataFrame.from_dict(metrics)[out_list]
-        return pd.concat((ms, leaves[out_list]))
-
-def extract(config: ExtractConfig):
+        # TODO this can be more efficient. Use count to find indices, then work
+        # with that smaller set from there. Working as is for now, but slow.
+        dit = tdb.query(attrs=[*att_list, *ma_list], order='F', coords=True,
+                return_incomplete=True, use_arrow=False).df[minx:maxx, miny:maxy]
+        data = pd.DataFrame()
+
+        storage.config.log.info('Collecting database information...')
+        for d in dit:
+            if data.empty:
+                data = d
+            else:
+                data = pd.concat([data, d])
+
+
+        # should find values that are not unique, meaning they have multiple entries
+        data = data.set_index(['X','Y'])
+        redo_indices = data.index[data.index.duplicated(keep='first')]
+        if redo_indices.empty:
+            return data.reset_index()
+
+        # data with overlaps
+        redo_data = data.loc[redo_indices][att_list].groupby(['X','Y']).agg(lambda x: list(chain(*x)))
+        # data that has no overlaps
+        clean_data = data.loc[data.index[~data.index.duplicated(False)]]
+
+        storage.config.log.warning('Overlapping data detected. Rerunning metrics over these cells...')
+        new_metrics = get_metrics(redo_data.reset_index(), storage)
+    return pd.concat([clean_data, new_metrics]).reset_index()
+
+def extract(config: ExtractConfig) -> None:
+    """
+    Pull data from database for each desired metric and output them to rasters
+
+    :param config: ExtractConfig.
+    """
 
     dask.config.set({"dataframe.convert-string": False})
 
-    ma_list = [ m.entry_name(a.name) for m in config.metrics for a in
-        config.attrs ]
-
-
     storage = Storage.from_db(config.tdb_dir)
+    ma_list = storage.getDerivedNames()
+    config.log.debug(f'Extracting metrics {[m for m in ma_list]}')
     root_bounds=storage.config.root
 
     e = Extents(config.bounds, config.resolution, root=root_bounds)
     i = e.get_indices()
-    minx = i['x'].min()
-    maxx = i['x'].max()
-    miny = i['y'].min()
-    maxy = i['y'].max()
-    x1 = maxx - minx + 1
-    y1 = maxy - miny + 1
+    xsize = e.x2
+    ysize = e.y2
 
+    # figure out if there are any overlaps and handle them
     final = handle_overlaps(config, storage, i).sort_values(['Y', 'X'])
 
     # output metric data to tifs
     for ma in ma_list:
-        m_data = np.full(shape=(y1,x1), fill_value=np.nan, dtype=final[ma].dtype)
+        # TODO should output in sections so we don't run into memory problems
+        m_data = np.full(shape=(ysize,xsize), fill_value=np.nan, dtype=final[ma].dtype)
         a = final[['X','Y',ma]].to_numpy()
         for x,y,md in a[:]:
             m_data[int(y)][int(x)] = md
 
-        write_tif(x1, y1, m_data, ma, config)
\ No newline at end of file
+        write_tif(xsize, ysize, m_data, ma, config)
\ No newline at end of file
diff --git a/src/silvimetric/commands/info.py b/src/silvimetric/commands/info.py
index 721b871..48dcf61 100644
--- a/src/silvimetric/commands/info.py
+++ b/src/silvimetric/commands/info.py
@@ -1,9 +1,25 @@
 from datetime import datetime
 from uuid import UUID
+from typing import Union
 
-from ..resources import Storage, Bounds
+from .. import Storage, Bounds
 
-def check_values(start_time, end_time, bounds, name):
+def check_values(start_time: datetime, end_time: datetime, bounds: Bounds,
+        name: Union[UUID, str]):
+    """
+    Validate arguments for info command.
+
+    :param start_time: Starting datetime object.
+    :param end_time: Ending datetime object.
+    :param bounds: Bounds to query by.
+    :param name: Name to query by.
+    :raises TypeError: Incorrect type of start_time argument.
+    :raises TypeError: Incorrect type of end_time argument.
+    :raises TypeError: Incorrect type of bounds argument.
+    :raises TypeError: Incorrect type of name argument.
+    :raises TypeError: Incorrect type of name argument.
+    :meta public:
+    """
     if start_time is not None and not isinstance(start_time, datetime):
         raise TypeError(f'Incorrect type of "start_time" argument.')
     if end_time is not None and not isinstance(end_time, datetime):
@@ -11,14 +27,28 @@ def check_values(start_time, end_time, bounds, name):
     if bounds is not None and not isinstance(bounds, Bounds):
         raise TypeError(f'Incorrect type of "bounds" argument.')
     if name is not None:
-        try:
-            UUID(name)
-        except:
+        if isinstance(name, UUID):
+            pass
+        elif isinstance(name, str):
+            try:
+                UUID(name)
+            except:
+                raise TypeError(f'Incorrect type of "name" argument.')
+        else:
             raise TypeError(f'Incorrect type of "name" argument.')
 
-def info(tdb_dir, start_time:datetime=None, end_time:datetime=None,
-          bounds:Bounds=None, name:str=None):
-    """Print info about Silvimetric database"""
+def info(tdb_dir:str, start_time:datetime=None, end_time:datetime=None,
+          bounds:Bounds=None, name:Union[str, UUID]=None) -> dict:
+    """
+    Collect information about database in current state
+
+    :param tdb_dir: TileDB database directory path.
+    :param start_time: Process starting time query, defaults to None
+    :param end_time: Process ending time query, defaults to None
+    :param bounds: Bounds query, defaults to None
+    :param name: Name query, defaults to None
+    :return: Returns json object containing information on database.
+    """
     check_values(start_time, end_time, bounds, name)
 
     with Storage.from_db(tdb_dir) as tdb:
@@ -29,7 +59,8 @@ def info(tdb_dir, start_time:datetime=None, end_time:datetime=None,
         if bounds is None:
             bounds = tdb.config.root
         if name is not None:
-            name = UUID(name)
+            if isinstance(name, str):
+                name = UUID(name)
 
         # We always have these
         meta = tdb.getConfig()
@@ -54,4 +85,4 @@ def info(tdb_dir, start_time:datetime=None, end_time:datetime=None,
         except KeyError:
             history = {}
 
-        return info
\ No newline at end of file
+    return info
\ No newline at end of file
diff --git a/src/silvimetric/commands/initialize.py b/src/silvimetric/commands/initialize.py
index cda309a..fccdd97 100644
--- a/src/silvimetric/commands/initialize.py
+++ b/src/silvimetric/commands/initialize.py
@@ -1,16 +1,13 @@
-from silvimetric.resources import Storage, StorageConfig
+from .. import Storage, StorageConfig
 
 def initialize(storage: StorageConfig):
     """
-    Initialize a Silvimetric TileDB instance for a given StorageConfig instance
-
-    Parameters
-    ----------
-    StorageConfig :  StorageConfig
+    Initialize a Silvimetric TileDB instance for a given StorageConfig instance.
 
+    :param storage: :class:`silvimetric.resources.config.StorageConfig`.
+    :return: :class:`silvimetric.resources.storage.Storage` database object.
     """
 
-    storage.log.debug(f"Initializing SilviMetric Database at '{storage.tdb_dir}'")
-
     s = Storage.create(storage)
+    storage.log.info(f"Initialized SilviMetric Database at '{storage.tdb_dir}'")
     return s
diff --git a/src/silvimetric/commands/manage.py b/src/silvimetric/commands/manage.py
new file mode 100644
index 0000000..cace6d4
--- /dev/null
+++ b/src/silvimetric/commands/manage.py
@@ -0,0 +1,54 @@
+from .. import Storage, ShatterConfig
+from .info import info
+from .shatter import shatter
+
+def delete(tdb_dir: str, name:str) -> ShatterConfig:
+    """
+    Delete Shatter process from database and return config for that process.
+
+    :param tdb_dir: TileDB database directory path.
+    :param name: UUID name of the Shatter process.
+    :raises KeyError: Shatter process with ID does not exist.
+    :raises ValueError: Shatter process with ID is missing a time reservation
+    :return: Config of process that was deleted.
+    """
+    res = info(tdb_dir=tdb_dir, name=name)
+
+    try:
+        config = ShatterConfig.from_dict(res['history'][0])
+    except LookupError:
+        raise KeyError(f'Shatter process with ID {name} does not exist')
+
+    try:
+        time_slot = config.time_slot
+    except KeyError:
+        raise ValueError(f'Shatter process with ID {name} is missing a time reservation.')
+
+    storage = Storage.from_db(tdb_dir)
+    return storage.delete(time_slot)
+
+def restart(tdb_dir: str, name: str) -> int:
+    """
+    Delete shatter process from database and run it again with the same config.
+
+    :param tdb_dir: TileDB database directory path.
+    :param name: UUID name of Shatter process.
+    :return: Point count of the restarted shatter process.
+    """
+
+    cfg = delete(tdb_dir, name)
+    return shatter(cfg)
+
+def resume(tdb_dir: str, name: str) -> int:
+    """
+    Resume partially completed shatter process. Process must partially completed
+    and have an already established time slot.
+
+    :param tdb_dir: TileDB database directory path.
+    :param name: UUID name of Shatter process.
+    :return: Point count of the restarted shatter process.
+    """
+    res = info(tdb_dir=tdb_dir, name=name)
+    assert len(res['history']) == 1
+    config = res['history'][0]
+    return shatter(config)
\ No newline at end of file
diff --git a/src/silvimetric/commands/scan.py b/src/silvimetric/commands/scan.py
index 66ef3de..a17987f 100644
--- a/src/silvimetric/commands/scan.py
+++ b/src/silvimetric/commands/scan.py
@@ -3,17 +3,41 @@
 import dask.bag as db
 import dask
 import math
+import json
+
+from .. import Storage, Data, Extents, Bounds, Log
+
+def scan(tdb_dir: str, pointcloud: str, bounds: Bounds, point_count:int=600000, resolution:float=100,
+        depth:int=6, filter:bool=False, log: Log=None):
+    """
+    Scan pointcloud and determine appropriate tile sizes.
+
+    :param tdb_dir: TileDB database directory.
+    :param pointcloud: Path to point cloud.
+    :param bounds: Bounding box to filter by.
+    :param point_count: Point count threshold., defaults to 600000
+    :param resolution: Resolution threshold., defaults to 100
+    :param depth: Tree depth threshold., defaults to 6
+    :param filter: Remove empty Extents. This takes longer, but is more accurage., defaults to False
+    :return: Returns list of point counts.
+    """
+
+    # TODO Scan should output other information about a file like bounds, pc, avg points per cell
+    with Storage.from_db(tdb_dir) as tdb:
 
-from ..resources import Storage, Data, Extents, Bounds
 
-def scan(tdb_dir, pointcloud, bounds, point_count=600000, resolution=100,
-        depth=6, filter=False):
+        if log is None:
+            logger = logging.getLogger('silvimetric')
+        else:
+            logger = log
+        data = Data(pointcloud, tdb.config, bounds)
 
-    logger = logging.getLogger('silvimetric')
-    with Storage.from_db(tdb_dir) as tdb:
+        thresholds = dict(thresholds=dict(resolution=resolution, point_count=point_count, depth=depth))
+        logger.debug(json.dumps(thresholds, indent=2))
 
-        data = Data(pointcloud, tdb.config, bounds)
         extents = Extents.from_sub(tdb_dir, data.bounds)
+        count = dask.delayed(data.estimate_count)(extents.bounds).persist()
+
 
         if filter:
             chunks = extents.chunk(data, resolution, point_count, depth)
@@ -21,22 +45,42 @@ def scan(tdb_dir, pointcloud, bounds, point_count=600000, resolution=100,
 
         else:
             cell_counts = extent_handle(extents, data, resolution, point_count,
-                depth)
+                depth, log)
 
+        num_cells = np.sum(cell_counts).item()
         std = np.std(cell_counts)
         mean = np.mean(cell_counts)
         rec = int(mean + std)
 
-        logger.info(f'Tiling information:')
-        logger.info(f'  Mean tile size: {mean}')
-        logger.info(f'  Std deviation: {std}')
-        logger.info(f'  Recommended split size: {rec}')
+        pc_info = dict(pc_info=dict(storage_bounds=tdb.config.root.to_json(),
+            data_bounds=data.bounds.to_json(), count=dask.compute(count)))
+        tiling_info = dict(tile_info=dict(num_cells=num_cells,
+            num_tiles=len(cell_counts), mean=mean, std_dev=std, recommended=rec))
 
-        return rec
+        final_info = pc_info | tiling_info
+        logger.info(json.dumps(final_info, indent=2))
 
+        return final_info
 
-def extent_handle(extent: Extents, data: Data, res_threshold=100,
-        pc_threshold=600000, depth_threshold=6):
+
+def extent_handle(extent: Extents, data: Data, res_threshold:int=100,
+        pc_threshold:int=600000, depth_threshold:int=6, log: Log = None) -> list[int]:
+    """
+    Recurisvely iterate through quad tree of this Extents object with given
+    threshold parameters.
+
+    :param extent: Current Extent.
+    :param data: Data object created from point cloud file.
+    :param res_threshold: Resolution threshold., defaults to 100
+    :param pc_threshold: Point count threshold., defaults to 600000
+    :param depth_threshold: Tree depth threshold., defaults to 6
+    :return: Returns list of Extents that fit thresholds.
+    """
+
+    if log is None:
+        logger = logging.getLogger('silvimetric')
+    else:
+        logger = log
 
     if extent.root is not None:
         bminx, bminy, bmaxx, bmaxy = extent.root.get()
@@ -58,12 +102,12 @@ def extent_handle(extent: Extents, data: Data, res_threshold=100,
 
     curr = db.from_delayed(tile_info(chunk, data, res_threshold, pc_threshold,
             depth_threshold))
-    logger = logging.getLogger('silvimetric')
     a = [ ]
 
     curr_depth = 0
     while curr.npartitions > 0:
-        logger.info(f'Chunking {curr.npartitions} tiles at depth {curr_depth}')
+        logger = logging.getLogger('silvimetric')
+        logger.debug(f'Chunking {curr.npartitions} tiles at depth {curr_depth}')
         n = curr.compute()
         to_add = [ x for x in n if isinstance(x, int) ]
         a = a + to_add
@@ -77,9 +121,20 @@ def extent_handle(extent: Extents, data: Data, res_threshold=100,
 
 
 @dask.delayed
-def tile_info(extent: Extents, data: Data, res_threshold=100,
-        pc_threshold=600000, depth_threshold=6, depth=0):
-
+def tile_info(extent: Extents, data: Data, res_threshold:int=100,
+        pc_threshold:int=600000, depth_threshold:int=6, depth:int=0):
+    """
+    Recursively explore current extents, use thresholds to determine when to
+    stop searching.
+
+    :param extent: Current Extent.
+    :param data: Data object created from point cloud file.
+    :param res_threshold: Resolution threshold., defaults to 100
+    :param pc_threshold: Point count threshold., defaults to 600000
+    :param depth_threshold: Tree depth threshold., defaults to 6
+    :param depth: Current Tree depth., defaults to 0
+    :return: Returns list of Extents that fit thresholds.
+    """
 
     pc = data.estimate_count(extent.bounds)
     target_pc = pc_threshold
diff --git a/src/silvimetric/commands/shatter.py b/src/silvimetric/commands/shatter.py
index c3a3389..bd9213c 100644
--- a/src/silvimetric/commands/shatter.py
+++ b/src/silvimetric/commands/shatter.py
@@ -1,160 +1,234 @@
 import numpy as np
-import gc
-import copy
 import signal
 import datetime
+import copy
+from typing import Generator
+import pandas as pd
+import tiledb
 
-import dask
-from dask.distributed import Client, as_completed, futures_of, CancelledError
+from dask.distributed import as_completed, futures_of, CancelledError
+from distributed.client import _get_global_client as get_client
+from dask.delayed import Delayed
 import dask.array as da
 import dask.bag as db
+from dask.diagnostics import ProgressBar
 
-from tiledb import SparseArray
-
-from ..resources import Extents, Storage, ShatterConfig, Data
+from .. import Extents, Storage, Data, ShatterConfig, run_metrics
 
 def get_data(extents: Extents, filename: str, storage: Storage):
-    #TODO look at making a record array here
+    """
+    Execute pipeline and retrieve point cloud data for this extent
+
+    :param extents: :class:`silvimetric.resources.extents.Extents` being operated on.
+    :param filename: Path to either PDAL pipeline or point cloud.
+    :param storage: :class:`silvimetric.resources.storage.Storage` database object.
+    :return: Point data array from PDAL.
+    """
     data = Data(filename, storage.config, bounds = extents.bounds)
+    p = data.pipeline
     data.execute()
-    return data.array
+    return p.get_dataframe(0)
+
+def arrange(points: pd.DataFrame, leaf, attrs: list[str]):
+    """
+    Arrange data to fit key-value TileDB input format.
+
+    :param data: Tuple of indices and point data array (xis, yis, data).
+    :param leaf: :class:`silvimetric.resources.extents.Extent` being operated on.
+    :param attrs: List of attribute names.
+    :raises Exception: Missing attribute error.
+    :return: None if no work is done, or a tuple of indices and rearranged data.
+    """
+    if points is None:
+        return None
+    if points.size == 0:
+        return None
 
-def cell_indices(xpoints, ypoints, x, y):
-    return da.logical_and(xpoints == x, ypoints == y)
+    points = points.loc[points.Y < leaf.bounds.maxy]
+    points = points.loc[points.Y >= leaf.bounds.miny]
+    points = points.loc[points.X >= leaf.bounds.minx]
+    points = points.loc[points.X < leaf.bounds.maxx, [*attrs, 'xi', 'yi']]
 
-def get_atts(points: np.ndarray, leaf: Extents, attrs: list[str]):
     if points.size == 0:
         return None
 
-    xis = da.floor(points[['xi']]['xi'])
-    yis = da.floor(points[['yi']]['yi'])
+    points.loc[:, 'xi'] = da.floor(points.xi)
+    # ceil for y because origin is at top left
+    points.loc[:, 'yi'] = da.ceil(points.yi)
 
-    att_view = points[:][attrs]
-    idx = leaf.get_indices()
-    l = [att_view[cell_indices(xis, yis, x, y)] for x,y in idx]
-    return l
+    return points
 
-def arrange(data: tuple[np.ndarray, np.ndarray, np.ndarray], leaf: Extents, attrs):
-    if data is None:
+def get_metrics(data_in, storage: Storage):
+    """
+    Run DataFrames through metric processes
+    """
+    if data_in is None:
         return None
 
-    di = data
-    dd = {}
-    for att in attrs:
-        try:
-            dd[att] = np.fromiter([*[np.array(col[att], col[att].dtype) for col in di], None], dtype=object)[:-1]
-        except Exception as e:
-            raise Exception(f"Missing attribute {att}: {e}")
-    counts = np.array([z.size for z in dd['Z']], np.int32)
-
-    ## remove empty indices
-    empties = np.where(counts == 0)[0]
-    dd['count'] = counts
-    idx = leaf.get_indices()
-    if bool(empties.size):
-        for att in dd:
-            dd[att] = np.delete(dd[att], empties)
-        idx = np.delete(idx, empties)
-
-    dx = idx['x']
-    dy = idx['y']
-    return (dx, dy, dd)
-
-def get_metrics(data_in, attrs: list[str], storage: Storage):
+    return run_metrics(data_in, storage.config.metrics)
+
+def agg_list(data_in):
+    """
+    Make variable-length point data attributes into lists
+    """
     if data_in is None:
         return None
 
-    ## data comes in as [dx, dy, { 'att': [data] }]
-    dx, dy, data = data_in
+    # grouped = data_in.groupby(['xi','yi'])
+    # grouped = grouped.agg(list)
+
+    # return grouped.assign(count=lambda x: [len(z) for z in grouped.Z])
+
+    # # coerce datatypes to object so we can store np arrays as cell values
+    old_dtypes = data_in.dtypes
+    xyi_dtypes = { 'xi': np.float64, 'yi': np.float64 }
+    o = np.dtype('O')
+    col_dtypes = { a: o for a in data_in.columns if a not in ['xi','yi'] }
 
-    # make sure it's not empty. No empty writes
-    if not np.any(data['count']):
+
+    coerced = data_in.astype(col_dtypes | xyi_dtypes)
+    a = coerced.groupby(['xi','yi']).agg(lambda x: np.array(x, old_dtypes[x.name]))
+    return a.assign(count=lambda x: [len(z) for z in a.Z])
+
+def join(list_data: pd.DataFrame, metric_data):
+    """
+    Join the list data and metric DataFrames together.
+    """
+    if list_data is None or metric_data is None:
         return None
 
-    # doing dask compute inside the dict array because it was too fine-grained
-    # when it was outside
-    metric_data = {
-        f'{m.entry_name(attr)}': [m(cell_data) for cell_data in data[attr]]
-        for attr in attrs for m in storage.config.metrics
-    }
-    data_out = data | metric_data
-    return (dx, dy, data_out)
+    if isinstance(metric_data, Delayed):
+        metric_data = metric_data.compute()
 
-def write(data_in, tdb):
+    return list_data.join(metric_data).reset_index()
 
-    if data_in is None:
-        return 0
+def write(data_in, storage, timestamp):
+    """
+    Write cell data to database
+
+    :param data_in: Data to be written to database.
+    :param tdb: TileDB write stream.
+    :return: Number of points written.
+    """
+
+    # data_in = data_in.reset_index()
+    data_in = data_in.rename(columns={'xi':'X','yi':'Y'})
 
-    dx, dy, dd = data_in
-    tdb[dx,dy] = dd
-    pc = int(dd['count'].sum())
+    attr_dict = {f'{a.name}': a.dtype for a in storage.config.attrs}
+    xy_dict = { 'X': data_in.X.dtype, 'Y': data_in.Y.dtype }
+    metr_dict = {f'{m.name}': m.dtype for m in storage.config.metrics}
+    dtype_dict = attr_dict | xy_dict | metr_dict
+
+    varlen_types = {a.dtype for a in storage.config.attrs}
+
+    tiledb.from_pandas(uri=storage.config.tdb_dir, sparse=True,
+        dataframe=data_in, mode='append', timestamp=timestamp,
+        column_types=dtype_dict, varlen_types=varlen_types)
+
+    pc = data_in['count'].sum().item()
     p = copy.deepcopy(pc)
-    del pc, data_in
 
+    del pc, data_in
     return p
 
-def run(leaves: db.Bag, config: ShatterConfig, storage: Storage,
-        tdb: SparseArray):
+Leaves = Generator[Extents, None, None]
+def get_processes(leaves: Leaves, config: ShatterConfig, storage: Storage) -> db.Bag:
+    """ Create dask bags and the order of operations.  """
+
+    ## Handle dask bag transitions through work states
+    attrs = [a.name for a in config.attrs]
+    timestamp = (config.time_slot, config.time_slot)
+
+    # remove any extents that have already been done, only skip if full overlap
+    leaf_bag: db.Bag = db.from_sequence(leaves)
+    if config.mbr:
+        def mbr_filter(one: Extents):
+            return all(one.disjoint_by_mbr(m) for m in config.mbr)
+        leaf_bag = leaf_bag.filter(mbr_filter)
+
+    def pc_filter(d: pd.DataFrame):
+        if d is None:
+            return False
+        return not d.empty
 
-    start_time = config.start_time
-    ## Handle a kill signal
+    points: db.Bag = leaf_bag.map(get_data, config.filename, storage)
+    arranged: db.Bag = points.map(arrange, leaf_bag, attrs)
+    filtered = arranged.filter(pc_filter)
+    metrics: db.Bag = filtered.map(run_metrics, storage.config.metrics)
+    lists: db.Bag = filtered.map(agg_list)
+    joined: db.Bag = lists.map(join, metrics)
+    writes: db.Bag = joined.map(write, storage, timestamp)
+
+    return writes
+
+def run(leaves: Leaves, config: ShatterConfig, storage: Storage) -> int:
+    """
+    Coordinate running of shatter process and handle any interruptions
+
+    :param leaves: Generator of Leaf nodes.
+    :param config: :class:`silvimetric.resources.config.ShatterConfig`
+    :param storage: :class:`silvimetric.resources.storage.Storage`
+    :return: Number of points processed.
+    """
+
+    # Process kill handler. Make sure we write out a config even if we fail.
     def kill_gracefully(signum, frame):
-        client = dask.config.get('distributed.client')
+        client = get_client()
         if client is not None:
             client.close()
         end_time = datetime.datetime.now().timestamp() * 1000
+
+        storage.consolidate_shatter(config.time_slot)
         config.end_time = end_time
-        with storage.open(timestamp=(start_time, end_time)) as a:
-            config.nonempty_domain = a.nonempty_domain()
-            config.finished=False
+        config.mbrs = storage.mbrs(config.time_slot)
+        config.finished=False
+        config.log.info('Saving config before quitting...')
 
-            config.log.info('Saving config before quitting...')
-            tdb.meta['shatter'] = str(config)
-            config.log.info('Quitting.')
+        storage.saveMetadata('shatter', str(config), config.time_slot)
+        config.log.info('Quitting.')
 
     signal.signal(signal.SIGINT, kill_gracefully)
 
-    ## Handle dask bag transitions through work states
-    attrs = [a.name for a in config.attrs]
+    processes = get_processes(leaves, config, storage)
 
-    leaf_bag: db.Bag = db.from_sequence(leaves)
-    points: db.Bag = leaf_bag.map(get_data, config.filename, storage)
-    att_data: db.Bag = points.map(get_atts, leaf_bag, attrs)
-    arranged: db.Bag = att_data.map(arrange, leaf_bag, attrs)
-    metrics: db.Bag = arranged.map(get_metrics, attrs, storage)
-    writes: db.Bag = metrics.map(write, tdb)
-
-    ## If dask is distributed, use the futures feature
-    if dask.config.get('scheduler') == 'distributed':
-        pc_futures = futures_of(writes.persist())
+        ## If dask is distributed, use the futures feature
+    dc = get_client()
+    if dc is not None:
+        pc_futures = futures_of(processes.persist())
         for batch in as_completed(pc_futures, with_results=True).batches():
             for future, pack in batch:
                 if isinstance(pack, CancelledError):
                     continue
                 for pc in pack:
                     config.point_count = config.point_count + pc
+                    del pc
 
         end_time = datetime.datetime.now().timestamp() * 1000
         config.end_time = end_time
         config.finished = True
-        a = storage.open(timestamp=(start_time, end_time))
-        return config.point_count
-
-    ## Handle non-distributed dask scenarios
-    pc = sum(writes)
-    end_time = datetime.datetime.now().timestamp() * 1000
-    with storage.open(timestamp=(start_time, end_time)) as a:
-
-        config.end_time = end_time
-        config.finished=True
-        config.point_count = pc
-        config.nonempty_domain = a.nonempty_domain()
-        tdb.meta['shatter'] = str(config)
-
-    return pc
-
-
-def shatter(config: ShatterConfig):
+    else:
+        # Handle non-distributed dask scenarios
+        with ProgressBar() as p:
+            config.point_count = sum(processes)
+
+    # modify config to reflect result of shattter process
+    config.mbr = storage.mbrs(config.time_slot)
+    config.log.debug('Saving shatter metadata')
+    config.end_time = datetime.datetime.now().timestamp() * 1000
+    config.finished = True
+
+    storage.saveMetadata('shatter', str(config), config.time_slot)
+    return config.point_count
+
+def shatter(config: ShatterConfig) -> int:
+    """
+    Handle setup and running of shatter process.
+    Will look for a config that has already been run before and needs to be
+    resumed.
+
+    :param config: :class:`silvimetric.resources.config.ShatterConfig`.
+    :return: Number of points processed.
+    """
     # get start time in milliseconds
     config.start_time = datetime.datetime.now().timestamp() * 1000
     # set up tiledb
@@ -162,27 +236,31 @@ def shatter(config: ShatterConfig):
     data = Data(config.filename, storage.config, config.bounds)
     extents = Extents.from_sub(config.tdb_dir, data.bounds)
 
+    config.log.info('Beginning shatter process...')
+    config.log.debug(f'Shatter Config: {config}')
+    config.log.debug(f'Data: {str(data)}')
+    config.log.debug(f'Extents: {str(extents)}')
 
-    with storage.open('w') as tdb:
+    if get_client() is None:
+        config.log.warning("Selected scheduler type does not support"
+            "continuously updated config information.")
 
-        ## shatter should still output a config if it's interrupted
-        # TODO add slices yet to be finished so we can pick up this process
-        # in the future
+    if not config.time_slot: # defaults to 0, which is reserved for storage cfg
+        config.time_slot = storage.reserve_time_slot()
 
-        if config.bounds is None:
-            config.bounds = extents.bounds
+    if config.bounds is None:
+        config.bounds = extents.bounds
 
-        config.log.debug('Grabbing leaf nodes...')
-        if config.tile_size is not None:
-            leaves = extents.get_leaf_children(config.tile_size)
-        else:
-            leaves = extents.chunk(data, 100)
+    config.log.debug('Grabbing leaf nodes...')
+    if config.tile_size is not None:
+        leaves = extents.get_leaf_children(config.tile_size)
+    else:
+        leaves = extents.chunk(data, 100)
 
-        # Begin main operations
-        config.log.debug('Fetching and arranging data...')
-        pc = run(leaves, config, storage, tdb)
-        config.point_count = int(pc)
+    # Begin main operations
+    config.log.debug('Fetching and arranging data...')
+    pc = run(leaves, config, storage)
 
-        config.log.debug('Saving shatter metadata')
-        tdb.meta['shatter'] = str(config)
-        return config.point_count
+    #consolidate the fragments in this time slot down to just one
+    storage.consolidate_shatter(config.time_slot)
+    return pc
\ No newline at end of file
diff --git a/src/silvimetric/resources/__init__.py b/src/silvimetric/resources/__init__.py
index a8c2dbe..e69de29 100644
--- a/src/silvimetric/resources/__init__.py
+++ b/src/silvimetric/resources/__init__.py
@@ -1,10 +0,0 @@
-__version__ = '0.0.1'
-
-from .bounds import Bounds
-from .extents import Extents
-from .storage import Storage
-from .metric import Metric, Metrics
-from .log import Log
-from .data import Data
-from .entry import Attribute, Pdal_Attributes, Attributes
-from .config import StorageConfig, ShatterConfig, ExtractConfig, ApplicationConfig
\ No newline at end of file
diff --git a/src/silvimetric/resources/array_extensions.py b/src/silvimetric/resources/array_extensions.py
new file mode 100644
index 0000000..36b07c0
--- /dev/null
+++ b/src/silvimetric/resources/array_extensions.py
@@ -0,0 +1,115 @@
+from __future__ import annotations
+
+from typing import Sequence
+
+import re
+
+import numpy as np
+import pandas as pd
+
+@pd.api.extensions.register_extension_dtype
+class AttributeDtype(pd.api.extensions.ExtensionDtype):
+    type = np.generic
+
+    def __init__(self, subtype=np.float32):
+        self._subtype = np.dtype(subtype)
+
+    def __str__(self) -> str:
+        return f'AttributeDtype[{self.subtype}]'
+
+    def __repr__(self) -> str:
+        return f'AttributeDtype[{self.subtype}]'
+
+    # TestDtypeTests
+    def __hash__(self) -> int:
+        return hash(str(self))
+
+
+    @property
+    def subtype(self):
+        return self._subtype
+
+    @property
+    def name(self):
+        return f"AttributeDtype[{self.subtype}]"
+
+    @classmethod
+    def construct_array_type(cls):
+        return AttributeArray
+
+    @classmethod
+    def construct_from_string(cls, string):
+        if string.lower() == "attribute":
+            return cls()
+        match = re.match(r"^attribute\[(\w+)\]$", string, re.IGNORECASE)
+        if match:
+            return cls(match.group(1))
+        raise TypeError(f"Cannot construct a 'AttributeDtype' from '{string}'")
+
+
+class AttributeArray(pd.api.extensions.ExtensionArray):
+    def __init__(self, arrays, dtype):
+        assert isinstance(dtype, AttributeDtype)
+        self._dtype = dtype
+        # self._data = np.array([np.array(array, dtype=dtype.subtype) for array in arrays], dtype=dtype.subtype)
+        self._data = [np.asarray(array, dtype=dtype.subtype) for array in arrays]
+
+    @classmethod
+    def _from_sequence(cls, scalars, dtype=None, copy=False):
+        if dtype is None and bool(scalars.dtype):
+            return cls(scalars, AttributeDtype(scalars.dtype))
+        elif dtype is None:
+            cls(scalars, AttributeDtype(np.generic))
+        return cls(scalars, dtype)
+
+    def __len__(self):
+        return len(self._data)
+
+    def __getitem__(self, i):
+        return self._data[i]
+
+    @property
+    def dtype(self):
+        return self._dtype
+
+    def copy(self):
+        return type(self)(self._data, self._dtype)
+
+    def take(self, indices, allow_fill=False, fill_value=None):
+        """ Take elements from an array.  """
+
+        if allow_fill and fill_value is None:
+            fill_value = self.dtype.na_value
+
+        result = pd.core.algorithms.take(self._data, indices, allow_fill=allow_fill,
+                                         fill_value=fill_value)
+        return self._from_sequence(result)
+
+    @classmethod
+    def _from_factorized(cls, uniques: np.ndarray, original: AttributeArray):
+        """ Reconstruct an AttributeArray after factorization.  """
+        return cls(uniques, AttributeDtype(uniques.dtype))
+
+    @pd.core.ops.unpack_zerodim_and_defer('__eq__')
+    def __eq__(self, other):
+        return self._apply_operator('__eq__', other, recast=False)
+
+    def nbytes(self):
+        return self._data.nbytes
+
+    def isna(self):
+        """ A 1-D array indicating if each value is missing.  """
+        return pd.isnull(self._data)
+
+    def copy(self):
+        """ Return a copy of the array.  """
+        copied = self._data.copy()
+        return type(self)(copied, self.dtype)
+
+    @classmethod
+    def _concat_same_type(cls, to_concat: Sequence[AttributeArray]) -> AttributeArray:
+        """
+        Concatenate multiple AngleArrays.
+        """
+        return cls(np.concatenate(to_concat))
+
diff --git a/src/silvimetric/resources/attribute.py b/src/silvimetric/resources/attribute.py
new file mode 100644
index 0000000..05dff52
--- /dev/null
+++ b/src/silvimetric/resources/attribute.py
@@ -0,0 +1,84 @@
+import json
+import numpy as np
+import pdal
+from tiledb import Attr
+from .array_extensions import AttributeArray, AttributeDtype
+
+class Attribute():
+    """Represents point data from a PDAL execution that has been binned, and
+    provides the information necessary to transfer that data to the database."""
+
+    def __init__(self, name: str, dtype) -> None:
+        self.name = name
+        """Name of the attribute, eg. Intensity."""
+        self.dtype: AttributeDtype
+        """SilviMetric representation of array of numpy dtype"""
+        if isinstance(dtype, AttributeDtype):
+            self.dtype = dtype
+        else:
+            # AttributeDtype takes any dtype that can be passed to np.dtype
+            try:
+                self.dtype = AttributeDtype(subtype=dtype)
+            except Exception as e:
+                raise AttributeError(f"Invalid dtype passed to Attribute: {dtype}") from e
+
+    def make_array(self, data, copy=False):
+        return AttributeArray(data=data, copy=copy)
+
+    def entry_name(self) -> str:
+        """Return TileDB attribute name."""
+        return self.name
+
+    def __eq__(self, other):
+        if self.dtype != other.dtype:
+            return False
+        elif self.name != other.name:
+            return False
+        else:
+            return True
+
+    def __hash__(self):
+        return hash(('name', self.name, 'dtype', self.dtype))
+
+
+    def schema(self) -> Attr:
+        """
+        Create the tiledb schema for this attribute.
+        :return: TileDB attribute schema
+        """
+        return Attr(name=self.name, dtype=self.dtype.subtype, var=True)
+
+    def to_json(self) -> object:
+        return {
+            'name': self.name,
+            'dtype': np.dtype(self.dtype.subtype).str
+        }
+
+    @staticmethod
+    def from_dict(data: dict) -> "Attribute":
+        """
+        Make an Attribute from a JSON like object
+        """
+        name = data['name']
+        dtype = data['dtype']
+        return Attribute(name, dtype)
+
+    @staticmethod
+    def from_string(data: str) -> "Attribute":
+        """
+        Create Attribute from string or dict version of it.
+
+        :param data: Stringified or json object of attribute
+        :raises TypeError: Incorrect type of incoming data, must be string or dict
+        :return: Return derived Attribute
+        """
+        j = json.loads(data)
+        return Attribute.from_dict(j)
+
+
+    def __repr__(self) -> str:
+        return json.dumps(self.to_json())
+
+# A list of pdal dimensions can be found here https://pdal.io/en/2.6.0/dimensions.html
+Pdal_Attributes = { d['name']: Attribute(d['name'], d['dtype']) for d in pdal.dimensions }
+Attributes = Pdal_Attributes
\ No newline at end of file
diff --git a/src/silvimetric/resources/bounds.py b/src/silvimetric/resources/bounds.py
index afd8242..d8a500c 100644
--- a/src/silvimetric/resources/bounds.py
+++ b/src/silvimetric/resources/bounds.py
@@ -1,12 +1,19 @@
 import json
 import ast
 
+#TODO should these bounds have a buffer?
 class Bounds(dict): #for JSON serializing
+    """Simple class to represent a 2 or 3-dimensional bounding box that can be
+    generated from both JSON or PDAL bounds form."""
     def __init__(self, minx: float, miny: float, maxx: float, maxy: float):
         self.minx = float(minx)
+        """minimum X Plane"""
         self.miny = float(miny)
+        """minimum Y plane"""
         self.maxx = float(maxx)
+        """maximum X plane"""
         self.maxy = float(maxy)
+        """maximum Y plane"""
 
     def __eq__(self, other):
         return  other.minx == self.minx and \
@@ -19,18 +26,24 @@ def __ne__(self, other):
     def __bool__(self):
         return True
 
+    def __repr__(self) -> str:
+        return str(self.get())
 
     @staticmethod
     def from_string(bbox_str: str):
-        """
-        Accepts bounds from strings in the form:
+        """Create Bounds object from a PDAL bounds string in the form:
 
         "([1,101],[2,102],[3,103])"
         "{\"minx\": 1,\"miny\": 2,\"maxx\": 101,\"maxy\": 102}"
         "[1,2,101,102]"
         "[1,2,3,101,102,103]"
 
+        :param bbox_str: Bounds string
+        :raises Exception: Unable to load Bounds via json or PDAL bounds type
+        :raises Exception: Bounding boxes must have either 4 or 6 elements
+        :return: Bounds object
         """
+
         try:
             bbox = json.loads(bbox_str)
         except json.decoder.JSONDecodeError as e:
@@ -67,20 +80,33 @@ def from_string(bbox_str: str):
 
 
     def get(self) -> list[float]:
-        return [self.minx, self.miny, self.maxx, self.maxy]
-
+        """Return Bounds as a list of floats
 
-    def __repr__(self) -> str:
-        return str(self.get())
+        :return: list of floats in form [minx, miny, maxx, maxy]
+        """
+        return [self.minx, self.miny, self.maxx, self.maxy]
 
 
     def to_string(self) -> str:
+        """Return string representation of Bounds
+
+        :return: string of a list of floats in form [minx, miny, maxx, maxy]
+        """
         return self.__repr__()
-    def to_json(self) -> str:
-        return json.dumps(self.get())
+
+    def to_json(self) -> list[float]:
+        """Return object as a json serializable list
+
+        :return: list of floats in form [minx, miny, maxx, maxy]
+        """
+        return self.get()
 
 
     def bisect(self):
+        """Bisects the current Bounds
+
+        :yield: 4 child bounds
+        """
         centerx = self.minx + ((self.maxx - self.minx)/ 2)
         centery = self.miny + ((self.maxy - self.miny)/ 2)
         yield from [
@@ -91,8 +117,47 @@ def bisect(self):
         ]
 
     def disjoint(self, other):
+        """Determine if two bounds are disjointed
+
+        :param other: Bounds this object is being compared to
+        :return: True if this box shares no point with the other box, otherwise False
+        """
         if other.minx > self.maxx or other.maxx < self.minx:
             return True
         if other.miny > self.maxy or other.maxy < self.miny:
             return True
         return False
+
+    def adjust_to_cell_lines(self, resolution):
+        xmindif = self.minx % resolution
+        xmaxdif = self.maxx % resolution
+        ymaxdif = self.maxy % resolution
+        ymindif = self.miny % resolution
+
+        if xmindif:
+            self.minx = self.minx - xmindif
+        if xmaxdif:
+            self.maxx = self.maxx + (resolution - xmaxdif)
+        if ymindif:
+            self.miny = self.miny - ymindif
+        if ymaxdif:
+            self.maxy = self.maxy + (resolution - ymaxdif)
+
+    @staticmethod
+    def shared_bounds(first, second):
+        """Find the Bounds that is shared between two Bounds.
+
+        :param first: First Bounds object for comparison.
+        :param second: Second Bounds object for comparison.
+        :returns: None if there is no overlap, otherwise the shared Bounds
+        """
+
+        if first.disjoint(second):
+            return None
+
+        minx = max(first.minx, second.minx)
+        maxx = min(first.maxx, second.maxx)
+        miny = max(first.miny, second.miny)
+        maxy = min(first.maxy, second.maxy)
+
+        return Bounds(minx, miny, maxx, maxy)
diff --git a/src/silvimetric/resources/config.py b/src/silvimetric/resources/config.py
index c2fba36..ce0a81e 100644
--- a/src/silvimetric/resources/config.py
+++ b/src/silvimetric/resources/config.py
@@ -13,10 +13,10 @@
 
 from .log import Log
 from .extents import Bounds
-from .metric import Metric, Metrics
-from .entry import Attribute, Attributes
-from . import constants
-from . import __version__
+from .metric import Metric
+from .metrics import grid_metrics
+from .attribute import Attribute, Attributes
+from .. import __version__
 
 
 class SilviMetricJSONEncoder(json.JSONEncoder):
@@ -25,9 +25,13 @@ def default(self, o):
 
 @dataclass(kw_only=True)
 class Config(ABC):
+    """Base config"""
     tdb_dir: str
+    """Path to TileDB directory to use."""
     log: Log = field(default_factory = lambda: Log("INFO"))
+    """Log object."""
     debug: bool = field(default=False)
+    """Debug flag."""
 
     def to_json(self):
         keys = self.__dataclass_fields__.keys()
@@ -38,6 +42,10 @@ def to_json(self):
             d['log'] = d['log'].to_json()
         return d
 
+    @classmethod
+    def from_json(self, data: str):
+        return self.from_string(json.dumps(data))
+
     @classmethod
     @abstractmethod
     def from_string(self, data: str):
@@ -51,18 +59,30 @@ def __repr__(self):
 
 @dataclass
 class StorageConfig(Config):
+    """ Config for constructing a Storage object """
     root: Bounds
+    """Root project bounding box"""
     crs: pyproj.CRS
-    resolution: float = DEFAULT_RESOLUTION
+    """Coordinate reference system, same for all data in a project"""
+    resolution: float = 30.0
+    """Resolution of cells, same for all data in a project, defaults to 30.0"""
 
     attrs: list[Attribute] = field(default_factory=lambda: [
         Attribute(a, Attributes[a].dtype)
         for a in [ 'Z', 'NumberOfReturns', 'ReturnNumber', 'Intensity' ]])
-    metrics: list[Metric] = field(default_factory=lambda: [ Metrics[m]
-                                  for m in Metrics.keys() ])
+    """List of :class:`silvimetric.resources.attribute.Attribute` attributes that
+    represent point data, defaults to Z, NumberOfReturns, ReturnNumber, Intensity"""
+    metrics: list[Metric] = field(default_factory=lambda: [ grid_metrics[m]
+                                  for m in grid_metrics.keys() ])
+    """List of :class:`silvimetric.resources.metrics.grid_metrics` grid_metrics that
+    represent derived data, defaults to values in grid_metrics object"""
     version: str = __version__
+    """Silvimetric version"""
     capacity: int = 1000000
-
+    """TileDB Capacity, defaults to 1000000"""
+    next_time_slot: int = 1
+    """Next time slot to be allocated to a shatter process. Increment after
+    use., defaults to 1"""
 
     def __post_init__(self) -> None:
 
@@ -77,12 +97,14 @@ def __post_init__(self) -> None:
             self.attrs = [Attributes[a] for a in [ 'Z', 'NumberOfReturns',
                                             'ReturnNumber', 'Intensity' ]]
         if not len(self.metrics):
-            self.metrics = [ Metrics[m] for m in Metrics.keys() ]
+            self.metrics = [ grid_metrics[m] for m in grid_metrics.keys() ]
 
         if not self.crs.is_projected:
-            raise Exception(f"Given coordinate system is not a rectilinear projected coordinate system")
+            raise Exception("Given coordinate system is not a rectilinear"
+                    " projected coordinate system")
 
-        self.metric_definitions = { m.name: str(m) for m in self.metrics}
+        self.metric_definitions = { m.name: str(m) for m in
+                self.metrics}
 
     def __eq__(self, other):
 
@@ -100,7 +122,8 @@ def to_json(self):
         d['attrs'] = [a.to_json() for a in self.attrs]
         d['metrics'] = [m.to_json() for m in self.metrics]
         d['crs'] = json.loads(self.crs.to_json())
-        d['root'] = json.loads(self.root.to_json())
+        d['root'] = self.root.to_json()
+
         return d
 
     @classmethod
@@ -108,17 +131,18 @@ def from_string(cls, data: str):
         x = json.loads(data)
         root = Bounds(*x['root'])
         if 'metrics' in x:
-            ms = [ Metric.from_string(m) for m in x['metrics']]
+            ms = [ Metric.from_dict(m) for m in x['metrics']]
         else:
-            ms = None
+            ms = [ ]
         if 'attrs' in x:
-            attrs = [ Attribute.from_string(a) for a in x['attrs']]
+            attrs = [ Attribute.from_dict(a) for a in x['attrs']]
         else:
-            attrs = None
+            attrs = [ ]
         if 'crs' in x:
             crs = pyproj.CRS.from_user_input(json.dumps(x['crs']))
         else:
             crs = None
+
         n = cls(tdb_dir = x['tdb_dir'],
                 root = root,
                 log = Log(**x['log']),
@@ -126,18 +150,36 @@ def from_string(cls, data: str):
                 attrs = attrs,
                 crs = crs,
                 metrics = ms,
-                capacity = x['capacity'])
+                capacity = x['capacity'],
+                version = x['version'])
 
         return n
 
     def __repr__(self):
-        return json.dumps(self.to_json())
+        j = self.to_json()
+        return json.dumps(j)
 
 @dataclass
 class ApplicationConfig(Config):
+    """ Base application config """
+
     debug: bool = False,
+    """Debug mode, defaults to False"""
     progress: bool = False,
+    """Should processes display progress bars, defaults to False"""
+
+    # Dask configuration
+    dasktype: str = 'processes'
+    """Dask parallelization type. For information see
+    https://docs.dask.org/en/stable/scheduling.html#local-threads """
     scheduler: str = 'distributed'
+    """Dask scheduler, defaults to 'distributed'"""
+    workers: int = 12
+    """Number of dask workers"""
+    threads: int = 4
+    """Number of threads per dask worker"""
+    watch: bool = False
+    """Open dask diagnostic page in default web browser"""
 
     def to_json(self):
         d = super().to_json()
@@ -148,31 +190,63 @@ def from_string(cls, data: str):
         x = json.loads(data)
         n = cls(tdb_dir = x['tdb_dir'],
                 debug = x['debug'],
-                progress = x['progress'])
+                progress = x['progress'],
+                dasktype = x['dasktype'],
+                scheduler = x['scheduler'],
+                workers = x['workers'],
+                threads = x['threads'],
+                watch = x['watch'])
         return n
 
     def __repr__(self):
         return json.dumps(self.to_json())
 
+Mbr = tuple[tuple[tuple[int,int], tuple[int,int]], ...]
 @dataclass
 class ShatterConfig(Config):
+    """Config for Shatter process"""
     filename: str
+    """Input filename referencing a PDAL pipeline or point cloud file."""
     date: Union[datetime, Tuple[datetime, datetime]]
+    """A date or date range representing data collection times."""
     attrs: list[Attribute] = field(default_factory=list)
+    """List of attributes to use in shatter. If this is not set it will be
+    filled by the attributes in the database instance."""
     metrics: list[Metric] = field(default_factory=list)
+    """A list of metrics to use in shatter. If this is not set it will be filled
+    by the metrics in the database instance."""
     bounds: Union[Bounds, None] = field(default=None)
+    """The bounding box of the shatter process., defaults to None"""
     name: uuid.UUID = field(default=uuid.uuid4())
+    """UUID representing this shatter process and will be generated if not
+    provided., defaults to uuid.uuid()"""
     tile_size: Union[int, None] = field(default=None)
+    """The number of cells to include in a tile., defaults to None"""
     start_time:float = 0
+    """The process starting time in seconds since Jan 1 1970., defaults to 0"""
     end_time:float = 0
+    """The process ending time in seconds since Jan 1 1970., defaults to 0"""
     point_count: int = 0
-    nonempty_domain: tuple[tuple[int, int], ...] = ()
+    """The number of points that has been processed so far., defaults to 0"""
+    mbr: Mbr = field(default_factory=tuple)
+    """The minimum bounding rectangle derived from TileDB array fragments.
+    This will be used to for resuming shatter processes and making sure it
+    doesn't repeat work., defaults to tuple()"""
     finished: bool = False
+    """Finished flag for shatter process., defaults to False"""
+    time_slot: int = 0
+    """The time slot that has been reserved for this shatter process. Will be
+    used as the timestamp in tiledb writes to better organize and manage
+    processes., defaults to 0"""
 
     def __post_init__(self) -> None:
         from .storage import Storage
         s = Storage.from_db(self.tdb_dir)
 
+        if isinstance(self.tile_size, float):
+            self.tile_size = int(self.tile_size)
+            self.log.warning(f'Truncating tile size to integer({self.tile_size})')
+            pass
         if not self.attrs:
             self.attrs = s.getAttributes()
         if not self.metrics:
@@ -184,10 +258,11 @@ def to_json(self):
         d = super().to_json()
 
         d['name'] = str(self.name)
-        d['bounds'] = json.loads(self.bounds.to_json()) if self.bounds is not None else None
+        d['time_slot'] = self.time_slot
+        d['bounds'] = self.bounds.to_json() if self.bounds is not None else None
         d['attrs'] = [a.to_json() for a in self.attrs]
         d['metrics'] = [m.to_json() for m in self.metrics]
-        d['nonempty_domain'] = [ list(a) for a in self.nonempty_domain]
+        d['mbr'] = list(self.mbr)
 
         if isinstance(self.date, tuple):
             d['date'] = [ dt.strftime('%Y-%m-%dT%H:%M:%SZ') for dt in self.date]
@@ -198,26 +273,37 @@ def to_json(self):
     @classmethod
     def from_string(cls, data: str):
         x = json.loads(data)
+        return cls.from_dict(x)
 
-        ms = list([ Metric.from_string(m) for m in x['metrics']])
-        attrs = list([ Attribute.from_string(a) for a in x['attrs']])
+    @classmethod
+    def from_dict(cls, data: dict):
+        x = data
+
+        ms = list([ Metric.from_dict(m) for m in x['metrics']])
+        attrs = list([ Attribute.from_dict(a) for a in x['attrs']])
         if isinstance(x['date'], list):
             date = tuple(( datetime.strptime(d, '%Y-%m-%dT%H:%M:%SZ') for d in x['date']))
         else:
             date = datetime.strptime(x['date'], '%Y-%m-%dT%H:%M:%SZ')
-        nonempty_domains = tuple(tuple(ned) for ned in x['nonempty_domain'])
+        mbr = tuple(tuple(tuple(mb) for mb in m) for m in x['mbr'])
         # TODO key error if these aren't there. If we're calling from_string
         # then these keys need to exist.
 
-        n = cls(tdb_dir = x['tdb_dir'],
-                filename = x['filename'],
-                attrs = attrs,
-                metrics = ms,
-                debug = x['debug'],
-                name = uuid.UUID(x['name']),
+        n = cls(tdb_dir=x['tdb_dir'],
+                filename=x['filename'],
+                attrs=attrs,
+                metrics=ms,
+                debug=x['debug'],
+                name=uuid.UUID(x['name']),
                 bounds=Bounds(*x['bounds']),
-                nonempty_domain=nonempty_domains,
-                date= date)
+                tile_size=x['tile_size'],
+                start_time=x['start_time'],
+                end_time=x['end_time'],
+                point_count=x['point_count'],
+                mbr=mbr,
+                date=date,
+                time_slot=x['time_slot'],
+                finished=x['finished'])
 
         return n
 
@@ -227,10 +313,17 @@ def __repr__(self):
 
 @dataclass
 class ExtractConfig(Config):
+    """Config for the Extract process."""
     out_dir: str
+    """The directory where derived rasters should be written."""
     attrs: list[Attribute] = field(default_factory=list)
+    """List of attributes to use in shatter. If this is not set it
+    will be filled by the attributes in the database instance."""
     metrics: list[Metric] = field(default_factory=list)
+    """A list of metrics to use in shatter. If this is not set it
+    will be filled by the metrics in the database instance."""
     bounds: Bounds = field(default=None)
+    """The bounding box of the shatter process., defaults to None"""
 
     def __post_init__(self) -> None:
         from .storage import Storage
@@ -251,20 +344,28 @@ def __post_init__(self) -> None:
     def to_json(self):
         d = super().to_json()
 
-        d['metrics'] = [a.to_json() for a in self.attrs]
+        d['attrs'] = [a.to_json() for a in self.attrs]
         d['metrics'] = [m.to_json() for m in self.metrics]
         d['crs'] = json.loads(self.crs.to_json())
-        d['bounds'] = json.loads(self.bounds.to_json())
+        d['bounds'] = self.bounds.to_json()
         return d
 
     @classmethod
     def from_string(cls, data: str):
         x = json.loads(data)
         if 'metrics' in x:
-            ms = [ Metric.from_string(m) for m in x['metrics']]
+            ms = [ Metric.from_dict(m) for m in x['metrics']]
         if 'attrs' in x:
-            attrs = [ Attribute.from_string(a) for a in x['attrs']]
-        n = cls(x['out_dir'], attrs, ms, x['debug'])
+            attrs = [ Attribute.from_dict(a) for a in x['attrs']]
+        if 'bounds' in x:
+            bounds = Bounds(*x['bounds'])
+        if 'log' in x:
+            l = x['log']
+            log = Log(l['log_level'], l['logdir'], l['logtype'], l['logfilename'])
+        else:
+            log = Log("INFO")
+        n = cls(tdb_dir=x['tdb_dir'], out_dir=x['out_dir'], attrs=attrs,
+                metrics=ms, debug=x['debug'], bounds=bounds, log=log)
 
         return n
 
diff --git a/src/silvimetric/resources/data.py b/src/silvimetric/resources/data.py
index ac2e03b..d17fd2e 100644
--- a/src/silvimetric/resources/data.py
+++ b/src/silvimetric/resources/data.py
@@ -1,6 +1,7 @@
-from . import Bounds
+from .bounds import Bounds
 from .config import StorageConfig
 import numpy as np
+import copy
 
 import pdal
 
@@ -8,25 +9,52 @@
 import json
 
 class Data:
+    """Represents a point cloud or PDAL pipeline, and performs essential operations
+    necessary to understand and execute a Shatter process."""
 
     def __init__(self,
                  filename: str,
                  storageconfig: StorageConfig,
                  bounds: Bounds = None):
         self.filename = filename
+        """Path to either PDAL pipeline or point cloud file"""
+
         self.bounds = bounds
+        """Bounds of this section of data"""
 
         self.reader_thread_count = 2
+        """Thread count for PDAL reader. Keep to 2 so we don't hog threads"""
 
-        self.storageconfig = storageconfig
         self.reader = self.get_reader()
-        self.pipeline = self.get_pipeline()
+        """PDAL reader"""
+
         if self.bounds is None:
             self.bounds = Data.get_bounds(self.reader)
 
+        # adjust bounds if necessary
+        self.bounds.adjust_to_cell_lines(storageconfig.resolution)
+        self.bounds = Bounds.shared_bounds(self.bounds, storageconfig.root)
+
+        self.storageconfig = storageconfig
+        """:class:`silvimetric.resources.StorageConfig`"""
+        self.pipeline = self.get_pipeline()
+        """PDAL pipeline"""
+
+        self.log = storageconfig.log
+
+    def to_json(self):
+        j = dict(filename=self.filename, bounds=self.bounds.get(),
+                pipeline=json.loads(self.pipeline.pipeline), is_pipeline=self.is_pipeline())
+        return j
+
+    def __repr__(self):
+        return json.dumps(self.to_json(), indent=2)
 
     def is_pipeline(self) -> bool:
-        """Does this instance represent a pdal.Pipeline or a simple filename"""
+        """Does this instance represent a pdal.Pipeline or a simple filename
+
+        :return: Return true if input is a pipeline
+        """
 
         p = pathlib.Path(self.filename)
         if p.suffix == '.json' and p.name != 'ept.json':
@@ -35,26 +63,37 @@ def is_pipeline(self) -> bool:
 
 
     def make_pipeline(self) -> pdal.Pipeline:
-        """Take a COPC or EPT endpoint and generate a PDAL pipeline for it"""
+        """Take a COPC or EPT endpoint and generate a PDAL pipeline for it
+
+        :return: Return PDAL pipeline
+        """
 
         reader = pdal.Reader(self.filename, tag='reader')
         reader._options['threads'] = self.reader_thread_count
         if self.bounds:
             reader._options['bounds'] = str(self.bounds)
-        class_zero = pdal.Filter.assign(value="Classification = 0")
-        rn = pdal.Filter.assign(value="ReturnNumber = 1 WHERE ReturnNumber < 1")
-        nor = pdal.Filter.assign(value="NumberOfReturns = 1 WHERE NumberOfReturns < 1")
-        ferry = pdal.Filter.ferry(dimensions="X=>xi, Y=>yi")
 
-        # assign_x = pdal.Filter.assign(
-        #     value=f"xi = (X - {self.root.minx}) / {resolution}")
-        # assign_y = pdal.Filter.assign(
-        #     value=f"yi = ({self.root.maxy} - Y) / {resolution}")
 
-        return reader | class_zero | rn | nor #| ferry | assign_x | assign_y #| smrf | hag
+        return reader.pipeline()
+        ## TODO
+        ## Remove these filters.assign stages for now
+        ## Waiting for a pdal/pdal-python fix to this
+        ## https://github.com/PDAL/python/issues/174
 
-    def get_pipeline(self):
-        """Fetch the pipeline for the instance"""
+        # class_zero = pdal.Filter.assign(value="Classification = 0")
+        # rn = pdal.Filter.assign(value="ReturnNumber = 1 WHERE ReturnNumber < 1")
+        # nor = pdal.Filter.assign(value="NumberOfReturns = 1 WHERE NumberOfReturns < 1")
+        # ferry = pdal.Filter.ferry(dimensions="X=>xi, Y=>yi")
+
+        # return reader | class_zero | rn | nor
+
+    def get_pipeline(self) -> pdal.Pipeline:
+        """Fetch the pipeline for the instance
+
+        :raises Exception: File type isn't COPC or EPT
+        :raises Exception: More than one reader detected
+        :return: Return PDAL pipline
+        """
 
         # If we are a pipeline, read and parse it. If we
         # aren't, go make_pipeline using some options that
@@ -84,7 +123,15 @@ def get_pipeline(self):
             # we only answer to copc or ept readers
             if stage_kind in allowed_readers:
                 if self.bounds:
-                    stage._options['bounds'] = str(self.bounds)
+                    res = self.storageconfig.resolution
+                    collar = Bounds(
+                        self.bounds.minx - res,
+                        self.bounds.miny - res,
+                        self.bounds.maxx + res,
+                        self.bounds.maxy + res
+                    )
+                    stage._options['bounds'] = str(collar)
+                    # stage._options['bounds'] = str(self.bounds)
 
             # We strip off any writers from the pipeline that were
             # given to us and drop them  on the floor
@@ -100,30 +147,46 @@ def get_pipeline(self):
         # Add xi and yi – only need this for PDAL < 2.6
         ferry = pdal.Filter.ferry(dimensions="X=>xi, Y=>yi")
         assign_x = pdal.Filter.assign(value=f"xi = (X - {self.storageconfig.root.minx}) / {resolution}")
-        assign_y = pdal.Filter.assign(value=f"yi = ({self.storageconfig.root.maxy} - Y) / {resolution}")
+        assign_y = pdal.Filter.assign(value=f"yi = (({self.storageconfig.root.maxy} - Y) / {resolution}) - 1")
+        # hag = pdal.Filter.hag_nn()
+
         stages.append(ferry)
         stages.append(assign_x)
         stages.append(assign_y)
+        # stages.append(hag)
 
         # return our pipeline
         return pdal.Pipeline(stages)
 
     def execute(self):
+        """Execute PDAL pipeline
+
+        :raises Exception: PDAL error message passed from execution
+        """
         try:
             self.pipeline.execute()
-            # self.storageconfig.log.debug(f"PDAL log: {self.pipeline.log}")
+            if self.pipeline.log and self.pipeline.log is not None:
+                self.log.debug(f"PDAL log: {self.pipeline.log}")
         except Exception as e:
+            if self.pipeline.log and self.pipeline.log is not None:
+                self.log.debug(f"PDAL log: {self.pipeline.log}")
             print(self.pipeline.pipeline, e)
             raise e
 
     def get_array(self) -> np.ndarray:
-        """Fetch the array from the execute()'d pipeline"""
+        """Fetch the array from the execute()'d pipeline
+
+        :return: get data as a numpy ndarray
+        """
         return self.pipeline.arrays[0]
     array = property(get_array)
 
     def get_reader(self) -> pdal.Reader:
-        """Grab or make the reader for this instance so we can use it to do things like
-        get the count()"""
+        """Grab or make the reader for this instance so we can use it to do things
+        like get the count()
+
+        :return: get PDAL reader for input
+        """
         if self.is_pipeline():
             p = pathlib.Path(self.filename)
             j = p.read_bytes().decode('utf-8')
@@ -140,12 +203,23 @@ def get_reader(self) -> pdal.Reader:
 
     @staticmethod
     def get_bounds(reader: pdal.Reader) -> Bounds:
+        """Get the bounding box of a point cloud from PDAL.
+
+        :param reader: PDAL Reader representing input data
+        :return: bounding box of point cloud
+        """
         p = reader.pipeline()
         qi = p.quickinfo[reader.type]
         return Bounds.from_string(json.dumps(qi['bounds']))
 
     def estimate_count(self, bounds: Bounds) -> int:
-        """For the provided bounds, estimate the maximum number of points that could be inside them for this instance."""
+        """For the provided bounds, estimate the maximum number of points that
+        could be inside them for this instance.
+
+        :param bounds: query bounding box
+        :return: estimated point count
+        """
+
         reader = self.get_reader()
         if bounds:
             reader._options['bounds'] = str(bounds)
@@ -156,8 +230,14 @@ def estimate_count(self, bounds: Bounds) -> int:
         return pc
 
     def count(self, bounds: Bounds) -> int:
-        """For the provided bounds, read and count the number of points that are inside them for this instance."""
-        reader = self.get_reader()
+        """For the provided bounds, read and count the number of points that are
+        inside them for this instance.
+
+        :param bounds: query bounding box
+        :return: point count
+        """
+
+        reader = copy.deepcopy(self.get_reader())
         if bounds:
             reader._options['bounds'] = str(bounds)
         pipeline = reader.pipeline()
diff --git a/src/silvimetric/resources/entry.py b/src/silvimetric/resources/entry.py
deleted file mode 100644
index bdbf6e6..0000000
--- a/src/silvimetric/resources/entry.py
+++ /dev/null
@@ -1,76 +0,0 @@
-import json
-import numpy as np
-import pdal
-from typing import Self, Union
-from abc import ABC, abstractmethod
-from tiledb import Attr
-
-class Entry(ABC):
-
-    def __eq__(self, other: Self):
-        return self.name == other.name and \
-            self.dtype == other.dtype and \
-            self.dependencies == other.dependencies
-
-    @abstractmethod
-    def entry_name(self) -> str:
-        raise NotImplementedError
-
-    @abstractmethod
-    def schema(self) -> Attr:
-        raise NotImplementedError
-
-    @abstractmethod
-    def to_json(self) -> object:
-        raise NotImplementedError
-
-    @staticmethod
-    @abstractmethod
-    def from_string(data: str) -> Self:
-        raise NotImplementedError
-
-    @abstractmethod
-    def __repr__(self):
-        raise NotImplementedError
-
-
-class Attribute(Entry):
-
-    def __init__(self, name: str, dtype: np.dtype, deps: list[Self]=None):
-        super().__init__()
-        self.name = name
-        self.dtype = dtype
-        self.dependencies = deps
-
-    def entry_name(self) -> str:
-        return self.name
-
-    def schema(self) -> Attr:
-        return Attr(name=self.name, dtype=self.dtype, var=True)
-
-    def to_json(self) -> object:
-        return {
-            'name': self.name,
-            'dtype': np.dtype(self.dtype).str,
-            'dependencies': self.dependencies,
-        }
-
-    @staticmethod
-    def from_string(data: Union[str, dict]) -> Self:
-        if isinstance(data, str):
-            j = json.loads(data)
-        elif isinstance(data, dict):
-            j = data
-        else:
-            raise TypeError(f'Type of data, {data} is not str or dict')
-        name = j['name']
-        dtype = j['dtype']
-        deps = j['dependencies']
-        return Attribute(name, dtype, deps)
-
-    def __repr__(self):
-        return json.dumps(self.to_json())
-
-# A list of pdal dimensions can be found here https://pdal.io/en/2.6.0/dimensions.html
-Pdal_Attributes = { d['name']: Attribute(d['name'], d['dtype']) for d in pdal.dimensions }
-Attributes = Pdal_Attributes
\ No newline at end of file
diff --git a/src/silvimetric/resources/extents.py b/src/silvimetric/resources/extents.py
index e9c164b..d3d0959 100644
--- a/src/silvimetric/resources/extents.py
+++ b/src/silvimetric/resources/extents.py
@@ -5,14 +5,18 @@
 import dask
 import dask.bag as db
 
-from typing import Self
 from math import ceil
 
+from .log import Log
 from .bounds import Bounds
 from .storage import Storage
 from .data import Data
 
+IndexDomain = tuple[float, float]
+IndexDomainList = tuple[IndexDomain, IndexDomain]
+
 class Extents(object):
+    """Handles bounds operations for point cloud data."""
 
     def __init__(self,
                  bounds: Bounds,
@@ -20,29 +24,87 @@ def __init__(self,
                  root: Bounds):
 
         self.bounds = bounds
+        """Bounding box of this section of data."""
         self.root = root
+        """Root bounding box of the database."""
 
+        # adjust bounds so they're matching up with cell lines
+        self.bounds.adjust_to_cell_lines(resolution)
         minx, miny, maxx, maxy = self.bounds.get()
 
+
         self.rangex = maxx - minx
+        """Range of X Indices"""
         self.rangey = maxy - miny
+        """Range of Y indices"""
         self.resolution = resolution
+        """Resolution of database."""
         self.cell_count = int((self.rangex * self.rangey) / self.resolution ** 2)
+        """Number of cells in this Extents"""
 
         self.x1 = math.floor((minx - self.root.minx) / resolution)
-        self.y1 = math.floor((self.root.maxy - maxy) / resolution)
+        """Minimum X index"""
         self.x2 = math.ceil((maxx - self.root.minx) / resolution)
+        """Maximum X index"""
+
+        self.y1 = math.ceil((self.root.maxy - maxy) / resolution)
+        """Minimum Y index, or maximum Y value in point cloud"""
         self.y2 = math.ceil((self.root.maxy - miny) / resolution)
+        """Maximum Y index, or minimum Y value in point cloud"""
+        self.domain: IndexDomainList = ((self.x1, self.x2), (self.y1, self.y2))
+        """Minimum bounding rectangle of this Extents"""
 
     def get_indices(self):
+        """
+        Create indices for this section of the database relative to the root
+        bounds.
+
+        :return: Indices of this bounding box
+        """
         return np.array(
             [(i,j) for i in range(self.x1, self.x2)
             for j in range(self.y1, self.y2)],
             dtype=[('x', np.int32), ('y', np.int32)]
         )
 
+    def disjoint_by_mbr(self, mbr):
+        """
+        Determine if this Extents shares any points with a minimum bounding
+        rectangle.
+
+        :param mbr: Minimum bounding rectangle as defined by TileDB.
+        :return: True if no shared points, false otherwise.
+        """
+        xs, ys = mbr
+        x1, x2 = xs
+        y1, y2 = ys
+        if x1 > self.x2 or x2 < self.x1:
+            return True
+        if y1 > self.y2 or y2 < self.y1:
+            return True
+        return False
+
+    def disjoint(self, other):
+        """
+        Determined if this Extents shares any points with another Extents object.
+
+        :param other: Extents object to compare against.
+        :return: True if no shared points, false otherwise.
+        """
+        return self.bounds.disjoint(other.bounds)
+
     def chunk(self, data: Data, res_threshold=100, pc_threshold=600000,
-            depth_threshold=6, scan=False):
+            depth_threshold=6):
+        """
+        Split up a dataset into tiles based on the given thresholds. Unlike Scan
+        this will filter out any tiles that contain no points.
+
+        :param data: Incoming Data object to oeprate on.
+        :param res_threshold: Resolution threshold., defaults to 100
+        :param pc_threshold: Point count threshold., defaults to 600000
+        :param depth_threshold: Tree depth threshold., defaults to 6
+        :return: Return list of Extents that fit the criteria
+        """
         if self.root is not None:
             bminx, bminy, bmaxx, bmaxy = self.root.get()
             r = self.root
@@ -53,6 +115,7 @@ def chunk(self, data: Data, res_threshold=100, pc_threshold=600000,
         # make bounds in scale with the desired resolution
         minx = bminx + (self.x1 * self.resolution)
         maxx = bminx + (self.x2 * self.resolution)
+
         miny = bmaxy - (self.y2 * self.resolution)
         maxy = bmaxy - (self.y1 * self.resolution)
 
@@ -62,13 +125,14 @@ def chunk(self, data: Data, res_threshold=100, pc_threshold=600000,
             self.root = chunk.bounds
 
         filtered = []
-        curr = db.from_delayed(chunk.filter(data, res_threshold, pc_threshold, depth_threshold))
-        curr_depth = 0
+        curr = db.from_delayed([ch.filter(data, res_threshold, pc_threshold,
+            depth_threshold, 1) for ch in chunk.split()])
+        curr_depth = 1
 
-        logger = logging.getLogger('silvimetric')
+        logger = data.storageconfig.log
         while curr.npartitions > 0:
 
-            logger.info(f'Filtering {curr.npartitions} tiles at depth {curr_depth}')
+            logger.debug(f'Filtering {curr.npartitions} tiles at depth {curr_depth}')
             n = curr.compute()
             to_add = [ne for ne in n if isinstance(ne, Extents)]
             if to_add:
@@ -81,27 +145,41 @@ def chunk(self, data: Data, res_threshold=100, pc_threshold=600000,
 
 
     def split(self):
+        """
+        Split this extent into 4 children along the cell lines
+
+        :return: Returns 4 child extents
+        """
         minx, miny, maxx, maxy = self.bounds.get()
 
         x_adjusted = math.floor((maxx - minx) / 2 / self.resolution)
-        y_adjusted = math.floor((maxy - miny) / 2 / self.resolution)
+        y_adjusted = math.ceil((maxy - miny) / 2 / self.resolution)
 
         midx = minx + (x_adjusted * self.resolution)
         midy = maxy - (y_adjusted * self.resolution)
 
-        return [
+        exts =  [
             Extents(Bounds(minx, miny, midx, midy), self.resolution, self.root), #lower left
             Extents(Bounds(midx, miny, maxx, midy), self.resolution, self.root), #lower right
             Extents(Bounds(minx, midy, midx, maxy), self.resolution, self.root), #top left
             Extents(Bounds(midx, midy, maxx, maxy), self.resolution, self.root)  #top right
         ]
+        return exts
 
-    # create quad tree of chunks for this bounds, run pdal quickinfo over this
-    # chunk to determine if there are any points available in this
-    # set a bottom resolution of ~1km
     @dask.delayed
-    def filter(self, data: Data, res_threshold=100, pc_threshold=600000, depth_threshold=6, depth=0) -> Self:
-
+    def filter(self, data: Data, res_threshold=100, pc_threshold=600000, depth_threshold=6, depth=0):
+        """
+        Creates quad tree of chunks for this bounds, runs pdal quickinfo over
+        this to determine if there are any points available. Uses a bottom resolution
+        of 1km.
+
+        :param data: Data object containing point cloud details.
+        :param res_threshold: Resolution threshold., defaults to 100
+        :param pc_threshold: Point count threshold., defaults to 600000
+        :param depth_threshold: Tree depth threshold., defaults to 6
+        :param depth: Current tree depth., defaults to 0
+        :return: Returns a list of Extents.
+        """
 
         pc = data.estimate_count(self.bounds)
         target_pc = pc_threshold
@@ -133,6 +211,12 @@ def filter(self, data: Data, res_threshold=100, pc_threshold=600000, depth_thres
 
 
     def _find_dims(self, tile_size):
+        """
+        Find most square-like Extents given the number of cells per tile.
+
+        :param tile_size: Number of cells per tile.
+        :return: Returns x and y coordinates in a list.
+        """
         s = math.sqrt(tile_size)
         if int(s) == s:
             return [s, s]
@@ -144,6 +228,12 @@ def _find_dims(self, tile_size):
         return [x, y]
 
     def get_leaf_children(self, tile_size):
+        """
+        Get children Extents with given number of cells per tile.
+
+        :param tile_size: Cells per tile.
+        :yield: Yield from list of child extents.
+        """
         res = self.resolution
         xnum, ynum = self._find_dims(tile_size)
 
@@ -167,13 +257,25 @@ def get_leaf_children(self, tile_size):
 
     @staticmethod
     def from_storage(tdb_dir: str):
+        """
+        Create Extents from information stored in database.
+
+        :param tdb_dir: TileDB database directory.
+        :return: Returns resulting Extents.
+        """
         storage = Storage.from_db(tdb_dir)
         meta = storage.getConfig()
         return Extents(meta.root, meta.resolution, meta.root)
 
     @staticmethod
     def from_sub(tdb_dir: str, sub: Bounds):
+        """
+        Create an Extents that is less than the overall extents of the database.
 
+        :param tdb_dir: TileDB database directory.
+        :param sub: Desired bounding box.
+        :return: Returns resulting Extents.
+        """
         storage = Storage.from_db(tdb_dir)
 
         meta = storage.getConfig()
@@ -183,6 +285,7 @@ def from_sub(tdb_dir: str, sub: Bounds):
 
         return Extents(sub, res, base)
 
+
     def __repr__(self):
         minx, miny, maxx, maxy = self.bounds.get()
         return f"([{minx:.2f},{maxx:.2f}],[{miny:.2f},{maxy:.2f}])"
diff --git a/src/silvimetric/resources/lmom4.py b/src/silvimetric/resources/lmom4.py
deleted file mode 100644
index 2f0253c..0000000
--- a/src/silvimetric/resources/lmom4.py
+++ /dev/null
@@ -1,29 +0,0 @@
-import numpy as np
-
-# adapted from: https://xiaoganghe.github.io/python-climate-visuals/chapters/data-analytics/scipy-basic.html
-# added b3 and l4 from:
-#   Vogel, R. M. and Fennessey, N. M.: L moment diagrams should replace product moment diagrams,
-#   Water Resour. Res., 29, 1745–1752, https://doi.org/10.1029/93WR00341, 1993.
-def lmom4(data):
-    # lmom4 returns the first four L-moments of data
-    # data is the 1-d array    
-    # n is the total number of points in data, j is the j_th point
-    #
-    # j range in for loops starts with 1 so we need to subtract 1 for all b# equations
-
-    n = len(data)
-    # sort in descending order
-    data = np.sort(data.reshape(n))[::-1]
-    b0 = np.mean(data)
-    b1 = np.array([(n - j - 1) * data[j] / n / (n - 1)
-                   for j in range(n)]).sum()
-    b2 = np.array([(n - j - 1) * (n - j - 2) * data[j] / n / (n - 1) / (n - 2)
-                   for j in range(n - 1)]).sum()
-    b3 = np.array([(n - j - 1) * (n - j - 2) * (n - j - 3) * data[j] / n / (n - 1) / (n - 2) / (n - 3)
-                   for j in range(n - 2)]).sum()
-    l1 = b0
-    l2 = 2 * b1 - b0
-    l3 = 6 * (b2 - b1) + b0
-    l4 = 20 * b3 - 30 * b2 + 12 * b1 - b0
-
-    return l1, l2, l3, l4
diff --git a/src/silvimetric/resources/log.py b/src/silvimetric/resources/log.py
index 0f5f00f..4b99be1 100644
--- a/src/silvimetric/resources/log.py
+++ b/src/silvimetric/resources/log.py
@@ -7,10 +7,10 @@
 """
 import logging
 import pathlib
-import os
+import sys
+
 from typing import Any
 from typing import Dict
-from typing import TYPE_CHECKING
 
 try:
     import websocket
@@ -62,15 +62,11 @@ def __init__(self,
         """
         Creates logging formatting and structure
 
-        Parameters
-        ----------
-        config:
-            Application config representing the runtime config
+        :param config: Application config representing the runtime config
         """
 
+        # need to be careful not to pull logging from previous runs
         self.logger = logging.getLogger("silvimetric")
-        self.logger.setLevel(log_level)
-        self.log_level = log_level
         self.logdir = logdir
         if logdir:
             self.logtype = 'file'
@@ -78,6 +74,13 @@ def __init__(self,
             self.logtype = logtype
         self.logfilename = logfilename
 
+        # do not recreate handlers if they're already present
+        if self.logger.handlers:
+            self.log_level = self.logger.level
+            return
+
+        self.log_level = log_level
+        self.logger.setLevel(log_level)
 
         # File Handler for Logging
         log_format = logging.Formatter(
@@ -121,7 +124,7 @@ def __init__(self,
         else:
             # if the user didn't specify a log dir, we just do the StreamHandler
             if not self.logdir:
-                log_handler = logging.StreamHandler()
+                log_handler = logging.StreamHandler(stream=sys.stdout)
                 log_handler.setFormatter(log_format)
                 self.logger.addHandler(log_handler)
 
diff --git a/src/silvimetric/resources/metric.py b/src/silvimetric/resources/metric.py
index 6fc5355..15790df 100644
--- a/src/silvimetric/resources/metric.py
+++ b/src/silvimetric/resources/metric.py
@@ -1,1090 +1,257 @@
 import json
-# import warnings
-import numpy as np
-from typing import Callable, Optional, Any, Union
-from scipy import stats
-from inspect import getsource
-from tiledb import Attr
-import dask
 import base64
-import dill
 
-from .entry import Attribute, Entry
-from . import lmom4
-from . import constants
+from typing import Callable, Optional, Any, Union, List, Self
+from uuid import uuid4
+from functools import reduce
+
+from tiledb import Attr
+import numpy as np
+import dill
+import pandas as pd
 
-# trap warnings as errors...mainly for scipy.stats.skew and scipy.stats.kurtosis
-# looks like this causes trouble for dask
-# warnings.filterwarnings("error")
+import dask
+from dask.delayed import Delayed
+from distributed import Future
+from .attribute import Attribute
 
-MetricFn = Callable[[np.ndarray, float, float, Optional[Union[Any, None]]], np.ndarray]
+MetricFn = Callable[[pd.DataFrame, Any], pd.DataFrame]
+FilterFn = Callable[[pd.DataFrame, Optional[Union[Any, None]]], pd.DataFrame]
 
 # Derived information about a cell of points
-
 ## TODO should create list of metrics as classes that derive from Metric?
-class Metric(Entry):
-    def __init__(self, name: str, method: MetricFn, dtype: np.dtype=np.float32,
-            deps: list[Attribute]=None):
-        super().__init__()
+class Metric():
+    """
+    A Metric is a TileDB entry representing derived cell data. There is a base set of
+    metrics available through Silvimetric, or you can create your own. A Metric
+    object has all the information necessary to facilitate the derivation of
+    data as well as its insertion into the database.
+    """
+    def __init__(self, name: str, dtype: np.dtype, method: MetricFn,
+            dependencies: list[Self]=[], filters: List[FilterFn]=[],
+            attributes: List[Attribute]=[]) -> None:
+
+        #TODO make deps, filters, attrs into tuples or sets, not lists so they're hashable
+
         self.name = name
+        """Metric name. eg. mean"""
+        self.dtype = np.dtype(dtype).str
+        """Numpy data type."""
+        self.dependencies = dependencies
+        """Metrics this is dependent on."""
         self._method = method
-        self.dtype = dtype
-        self.dependencies = deps
+        """The method that processes this data."""
+        self.filters = filters
+        """List of user-defined filters to perform before performing method."""
+        self.attributes = attributes
+        """List of Attributes this Metric applies to. If empty it's used for all
+        Attributes"""
 
-    def schema(self, attr: Attribute):
+    def __eq__(self, other):
+        if self.name != other.name:
+            return False
+        elif self.dtype != other.dtype:
+            return False
+        elif self.dependencies != other.dependencies:
+            return False
+        elif self._method != other._method:
+            return False
+        elif self.filters != other.filters:
+            return False
+        elif self.attributes != other.attributes:
+            return False
+        else:
+            return True
+
+    def __hash__(self):
+        return hash((
+            'name', self.name,
+            'dtype', self.dtype,
+            'dependencies', frozenset(self.dependencies),
+            'method', base64.b64encode(dill.dumps(self._method)).decode(),
+            'filters', frozenset(base64.b64encode(dill.dumps(f)).decode() for f in self.filters),
+            'attrs', frozenset(self.attributes)))
+
+    def schema(self, attr: Attribute) -> Any:
+        """
+        Create schema for TileDB creation.
+
+        :param attr: :class:`silvimetric.resources.entry.Atttribute`
+        :return: TileDB Attribute
+        """
         entry_name = self.entry_name(attr.name)
         return Attr(name=entry_name, dtype=self.dtype)
 
-    # common name, storage name
     def entry_name(self, attr: str) -> str:
+        """Name for use in TileDB and extract file generation."""
         return f'm_{attr}_{self.name}'
 
-    @dask.delayed
-    def delayed(self, data: np.ndarray) -> np.ndarray:
-        return self._method(data)
+    def sanitize_and_run(self, d, *args):
+        # args come in as a value wrapped in one index of a 2D dataframe
+        # we need to remove the wrapping and pass the values to the method
+        # to make things easier
+        a = tuple( a.values[0][0] for a in args)
+        return self._method(d, *a)
+
+    def do(self, data: pd.DataFrame, *args) -> pd.DataFrame:
+        """Run metric and filters. Use previously run metrics to avoid running
+        the same thing multiple times."""
+
+        # the index columns are determined by where this data is coming from
+        # if it has xi and yi, then it's coming from shatter
+        # if it has X and Y, then it's coming from extract as a rerun of a cell
+        if isinstance(data, Future):
+            data = data.result()
+
+        idx = ['xi','yi']
+        if any([i not in data.columns for i in idx]):
+            idx = ['X','Y']
+
+        if self.attributes:
+            attrs = [*[a.name for a in self.attributes],*idx]
+            data = data[attrs]
+
+        data = self.run_filters(data)
+        gb = data.groupby(idx)
+
+        # create map of current column name to tuple of new column name and metric method
+        cols = data.columns
+        runner = lambda d: self.sanitize_and_run(d, *args)
+
+        new_cols = {
+            c: [ (self.entry_name(c), runner) ]
+            for c in cols if c not in idx
+        }
+
+        val = gb.aggregate(new_cols, args=args)
+
+        #remove hierarchical columns
+        val.columns = val.columns.droplevel(0)
+        return val
+
+    #TODO make dict with key for each Attribute effected? {att: [fn]}
+    # for now these filters apply to all Attributes
+    def add_filter(self, fn: FilterFn, desc: str):
+        """
+        Add filter method to list of filters to run before calling main method.
+        """
+        self.filters.append(fn)
+
+    def run_filters(self, data: pd.DataFrame) -> pd.DataFrame:
+        for f in self.filters:
+            ndf = f(data)
+            #TODO should this check be here?
+            if not isinstance(ndf, pd.DataFrame):
+                raise TypeError('Filter outputs must be a DataFrame. '
+                        f'Type detected: {type(ndf)}')
+            data = ndf
+        return data
+
 
     def to_json(self) -> dict[str, any]:
         return {
             'name': self.name,
             'dtype': np.dtype(self.dtype).str,
-            'dependencies': self.dependencies,
-            'method_str': getsource(self._method),
-            'method': base64.b64encode(dill.dumps(self._method)).decode()
+            'dependencies': [d.to_json() for d in self.dependencies],
+            'method': base64.b64encode(dill.dumps(self._method)).decode(),
+            'filters': [base64.b64encode(dill.dumps(f)).decode() for f in self.filters],
+            'attributes': [a.to_json() for a in self.attributes]
         }
 
-    def from_string(data: Union[str, dict]):
-        if isinstance(data, str):
-            j = json.loads(data)
-        elif isinstance(data, dict):
-            j = data
-        name = j['name']
-        dtype = np.dtype(j['dtype'])
-        dependencies = j['dependencies']
-        method = dill.loads(base64.b64decode(j['method'].encode()))
-
-        return Metric(name=name, method=method, dtype=dtype, deps=dependencies)
-
-    def __eq__(self, other):
+    @staticmethod
+    def from_dict(data: dict) -> "Metric":
+        name = data['name']
+        dtype = np.dtype(data['dtype'])
+        method = dill.loads(base64.b64decode(data['method'].encode()))
+
+        if 'dependencies' in data.keys() and \
+                data['dependencies'] and \
+                data['dependencies'] is not None:
+            dependencies = [ Metric.from_dict(d) for d in data['dependencies'] ]
+        else:
+            dependencies = [ ]
+
+        if 'attributes' in data.keys() and \
+                data['attributes'] and \
+                data['attributes'] is not None:
+            attributes = [ Attribute.from_dict(a) for a in data['attributes']]
+        else:
+            attributes = [ ]
+
+        if 'filters' in data.keys() and \
+                data['filters'] and \
+                data['filters'] is not None:
+            filters = [ dill.loads(base64.b64decode(f)) for f in data['filters'] ]
+        else:
+            filters = [ ]
+
+        return Metric(name, dtype, method, dependencies, filters, attributes)
+
+    @staticmethod
+    def from_string(data: str) -> Self:
+        j = json.loads(data)
+        return Metric.from_dict(j)
+
+    def __eq__(self, other) -> tuple:
         return (self.name == other.name and
                 self.dtype == other.dtype and
                 self.dependencies == other.dependencies and
-                self._method == other._method)
+                self._method == other._method,
+                self.attributes == other.attributes,
+                self.filters == other.filters)
 
-    def __call__(self, data: np.ndarray) -> np.ndarray:
-        return self._method(data)
+    def __call__(self, data: pd.DataFrame) -> pd.DataFrame:
+        return self.do(data)
 
     def __repr__(self) -> str:
-        return json.dumps(self.to_json())
-
-#TODO add all metrics from https://github.com/hobuinc/silvimetric/issues/5
-
-def m_count(data):
-    return len(data)
-
-def m_mean(data):
-    """
-    Return the arithmetic mean.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    return np.mean(data)
-
-def m_mode(data):
-    """
-    Return the mode (most common value).
-
-    Mode is practically undefined for floating point values. FUSION's logic is to
-    partition data into 64 bins then find the bin with the highest count.
-    This approach has the disadvantage that the bin size varies depending on
-    the Z range.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    nbins = 64
-    maxv = np.max(data)
-    minv = np.min(data)
-    if minv == maxv:
-        return minv
-    
-    bins, binedges = np.histogram(data, bins = nbins, density = False)
-    
-    thebin = np.argmax(bins)
-    
-    # compute the height and return...nbins - 1 is to get the bottom value of the bin
-    return minv + thebin * (maxv - minv) / (nbins - 1)
-
-def m_median(data):
-    """
-    Return the median.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    return np.median(data)
-
-def m_min(data):
-    """
-    Return the minimum value.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    return np.min(data)
-
-def m_max(data):
-    """
-    Return the maximum value.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    return np.max(data)
-
-def m_stddev(data):
-    """
-    Return the standard deviation of values.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    return np.std(data)
-
-# start of new metrics to match FUSION
-def m_variance(data):
-    """
-    Return the variance of values.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    return np.var(data)
-
-def m_cv(data):
-    """
-    Return the coefficient of variation of values.
-    Standard deviation / mean.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    return np.std(data) / np.mean(data)
-
-# TODO check performance of other methods
-def m_abovemean(data):
-    """
-    Return the proportion of values above the mean.
+        return f"Metric_{self.name}"
 
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    return (data > np.mean(data)).sum() / len(data)
-
-# TODO check performance of other methods
-def m_abovemode(data):
-    """
-    Return the proportion of values above the mode.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    return (data > m_mode(data)).sum() / len(data)
-
-# both stats.skew and stats.kurtosis throw a warning:
-# Precision loss occurred in moment calculation due to catastrophic cancellation. 
-# This occurs when the data are nearly identical. Results may be unreliable.
-#
-# I think this is because values are very similar...end up close to the mean
-# GridMetrics computes these using a different formula that may also have trouble
-# with numeric stability.
-def m_skewness(data):
+def get_methods(data: Union[pd.DataFrame, Delayed], metrics: Metric | list[Metric],
+        uuid=None) -> list[Delayed]:
     """
-    Return the skewness of values.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    if len(data) < 4:
-        return NODATA
-
-    try:
-        s = stats.skew(data)
-    except:
-        s = NODATA
-
-    return s
-
-def m_kurtosis(data):
+    Create Metric dependency graph by iterating through desired metrics and
+    their dependencies, creating Delayed objects that can be run later.
     """
-    Return the kurtosis of values.
+    # identitity for this graph, can be created before or during this method
+    # call, but needs to be the same across this graph, and unique compared
+    # to other graphs
+    if uuid is None:
+        uuid = uuid4()
 
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
+    # don't duplicate a delayed object
+    if not isinstance(data, Delayed):
+        data = dask.delayed(data)
 
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
+    if isinstance(metrics, Metric):
+        metrics = [ metrics ]
 
-    if len(data) < 4:
-        return NODATA
+    # iterate through metrics and their dependencies.
+    # uuid here will help guide metrics to use the same dependency method
+    # calls from dask
+    seq = []
+    for m in metrics:
+        if not isinstance(m, Metric):
+            continue
+        ddeps = get_methods(data, m.dependencies, uuid)
+        dd = dask.delayed(m.do)(data, *ddeps,
+            dask_key_name=f'{m.name}-{str(uuid)}')
+        seq.append(dd)
 
-    try:
-        k = stats.kurtosis(data)
-    except:
-        k = NODATA
+    return seq
 
-    return k
-
-def m_aad(data):
+def run_metrics(data: pd.DataFrame, metrics: Union[Metric, list[Metric]]) -> pd.DataFrame:
     """
-    Return the average absolute deviation from the mean.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    m = m_mean(data)
-    return np.mean(np.absolute(data - m))
-
-def m_madmedian(data):
-    """
-    Return the median absolute difference from the median value.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    return stats.median_abs_deviation(data)
-
-def m_madmean(data):
-    """
-    Return the median absolute difference from the mean value.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    return stats.median_abs_deviation(data, center=np.mean)
-
-# TODO needs work
-def m_madmode(data):
-    """
-    Return the median absolute difference from the mode value.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
+    Collect Metric dependency graph and run it, then merge the results together.
     """
+    graph = get_methods(data, metrics)
 
-    m = m_mode(data)
-    return np.median(np.absolute(data - m))
-
-# TODO test various methods for interpolation=... for all percentile-related metrics
-# I think the default matches FUSION method but need to test
-def m_iq(data):
-    """
-    Return the interquartile difference.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    return stats.iqr(data)
-
-def m_90m10(data):
-    """
-    Return the 90th percentile minus the 10th percentile.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    p = np.percentile(data, [10,90])
-    return p[1] - p[0]
-
-def m_95m05(data):
-    """
-    Return the 95th percentile minnus the 5th percentile.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    p = np.percentile(data, [5,95])
-    return p[1] - p[0]
-
-def m_crr(data):
-    """
-    Return the cannopy relief ratio.
-    (mean - min) / (max - min).
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    maxv = np.max(data)
-    minv = np.min(data)
-    if minv == maxv:
-        return NODATA
-    
-    return (np.mean(data) - minv) / (maxv - minv)
-
-def m_sqmean(data):
-    """
-    Return the square root of the average squared value.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    return np.sqrt(np.mean(np.square(data)))
-
-def m_cumean(data):
-    """
-    Return the cube root of the average cubed value.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    return np.cbrt(np.mean(np.power(np.absolute(data), 3)))
-
-# TODO compute L-moments. These are done separately because we only add
-# a single element to TileDB. This is very inefficient since we have to
-# compute all L-moments at once. Ideally, we would have a single metric
-# function that returns an array with 7 values
-
-# added code to compute first 4 l-moments in lmom4.py. There is a package,
-# lmoments3 that can compute the same values but I wanted to avoid the
-# package as it has some IBM copyright requirements (need to include copyright
-# statement with derived works)
-
-# L1 is same as mean...compute using np.mean for speed
-def m_l1(data):
-    """
-    Return the first L-moment (mean).
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    if len(data) < 4:
-        return NODATA
-
-    return np.mean(data)
-
-def m_l2(data):
-    """
-    Return the second L-moment.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    if len(data) < 4:
-        return NODATA
-
-    l = lmom4(data)
-    return l[1]
-
-def m_l3(data):
-    """
-    Return the third L-moment.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    if len(data) < 4:
-        return NODATA
-
-    l = lmom4(data)
-    return l[2]
-
-def m_l4(data):
-    """
-    Return the fourth L-moment.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    if len(data) < 4:
-        return NODATA
-
-    l = lmom4(data)
-    return l[3]
-
-def m_lcv(data):
-    """
-    Return the L-moment coefficient of variation (L2 / mean).
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    if len(data) < 4:
-        return NODATA
-
-    l = lmom4(data)
-
-    if l[0] == 0.0:
-        return NODATA
-    
-    return l[1] / l[0]
-
-def m_lskewness(data):
-    """
-    Return the L-moment skewness (L3 / L2).
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    if len(data) < 4:
-        return NODATA
-
-    l = lmom4(data)
-    if l[1] == 0.0:
-        return NODATA
-    
-    return l[2] / l[1]
-
-def m_lkurtosis(data):
-    """
-    Return the L-moment kurtosis (L4 / L2).
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    if len(data) < 4:
-        return NODATA
-
-    l = lmom4(data)
-    if l[1] == 0.0:
-        return NODATA
-    
-    return l[3] / l[1]
-
-# not sure how an array of metrics can be ingested by shatter
-# so do these as 15 separate metrics. may be slower than doing all in one call
-#def m_percentiles(data):
-#    return(np.percentile(data, [1,5,10,20,25,30,40,50,60,70,75,80,90,95,99]))
-
-def m_p01(data):
-    """
-    Return the 1st percentile value.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    return(np.percentile(data, 1))
-
-def m_p05(data):
-    """
-    Return the 5th percentile value.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    return(np.percentile(data, 5))
-
-def m_p10(data):
-    """
-    Return the 10th percentile value.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    return(np.percentile(data, 10))
-
-def m_p20(data):
-    """
-    Return the 20th percentile value.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    return(np.percentile(data, 20))
-
-def m_p25(data):
-    """
-    Return the 25th percentile value.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    return(np.percentile(data, 25))
-
-def m_p30(data):
-    """
-    Return the 30th percentile value.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    return(np.percentile(data, 30))
-
-def m_p40(data):
-    """
-    Return the 40th percentile value.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    return(np.percentile(data, 40))
-
-def m_p50(data):
-    """
-    Return the 50th percentile value.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    return(np.percentile(data, 50))
-
-def m_p60(data):
-    """
-    Return the 60th percentile value.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    return(np.percentile(data, 60))
-
-def m_p70(data):
-    """
-    Return the 70th percentile value.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    return(np.percentile(data, 70))
-
-def m_p75(data):
-    """
-    Return the 75th percentile value.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    return(np.percentile(data, 75))
-
-def m_p80(data):
-    """
-    Return the 80th percentile value.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    return(np.percentile(data, 80))
-
-def m_p90(data):
-    """
-    Return the 90th percentile value.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    return(np.percentile(data, 90))
-
-def m_p95(data):
-    """
-    Return the 95th percentile value.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    return(np.percentile(data, 95))
-
-def m_p99(data):
-    """
-    Return the 99th percentile value.
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    return(np.percentile(data, 99))
-
-def m_profilearea(data):
-    """
-    Return the profile area.
-
-    Profile area was first described in 
-
-        Hu, Tianyu, Qin Ma, Yanjun Su, John J. Battles, Brandon M. Collins, 
-        Scott L. Stephens, Maggi Kelly, Qinghua Guo. 2019. A simple and 
-        integrated approach for fire severity assessment using bi-temporal 
-        airborne LiDAR data. Int J. Appl Earth Obs Geoinformation, 78 (2019): 25-38.
-
-    as the area under the height percentile profile or curve. They found the 
-    metric useful to compare pre- and post-fire canopy structure at both the 
-    individual tree and pixel scales. The implementation in CloudMetrics and 
-    GridMetrics varies from that described in Hu et al. 2019. Heights are 
-    normalized using the 99th percentile height instead of the maximum height 
-    to eliminate problems related to high outliers and the area under the 
-    percentile curve is computed directly from 1 percent slices instead of 
-    fitting a polynomial to the percentile heights and computing the area 
-    under the polynomial
-
-    Parameters
-    ----------
-    data : numpy.ndarray
-        Data for dimension
-
-    Returns
-    -------
-    float
-        Value for metric computed using values in data. Returns
-        NODATA value when metric cannot be computed
-    """
-
-    # sanity check...must have valid heights/elevations
-    if np.max(data) <= 0:
-        return NODATA
-
-
-    p = np.percentile(data, range(1, 100))
-    p0 = max(np.min(data), 0.0)
-
-    # second sanity check...99th percentile must be > 0
-    if p[98] > 0.0:
-        # compute area under normalized percentile height curve using composite trapeziod rule
-        pa = p0 / p[98]
-        for ip in p[:97]:
-            pa += 2.0 * ip / p[98]
-        pa += 1.0
+    # try returning just the graph and see if that can speed thigns up
+    # return graph
 
-        return pa * 0.5
-    else:
-        return NODATA
+    computed_list = dask.persist(*graph, optimize_graph=True)
 
-# not sure how this will be computed once height and cover thresholds are
-# sorted out
-#def m_allcover(data):
-#    return (data > coverthreshold).sum() / len(data) * 100.0
+    def merge(x, y):
+        return x.merge(y, on=['xi','yi'])
+    merged = reduce(merge, computed_list)
 
-Metrics = {
-    'count' : Metric('count', m_count),
-    'mean' : Metric('mean', m_mean),
-    'mode' : Metric('mode', m_mode),
-    'median' : Metric('median', m_median),
-    'min' : Metric('min', m_min),
-    'max' : Metric('max', m_max),
-    'stddev' : Metric('stddev', m_stddev),
-    'variance' : Metric('variance', m_variance),
-    'cv' : Metric('cv', m_cv),
-    'abovemean' : Metric('abovemean', m_abovemean),
-    'abovemode' : Metric('abovemode', m_abovemode),
-    'skewness' : Metric('skewness', m_skewness),
-    'kurtosis' : Metric('kurtosis', m_kurtosis),
-    'aad' : Metric('aad', m_aad),
-    'madmedian' : Metric('madmedian', m_madmedian),
-    'madmode' : Metric('madmode', m_madmode),
-    'iq' : Metric('iq', m_iq),
-    'crr' : Metric('crr', m_crr),
-    'sqmean' : Metric('sqmean', m_sqmean),
-    'cumean' : Metric('cumean', m_cumean),
-    'l1' : Metric('l1', m_l1),
-    'l2' : Metric('l2', m_l2),
-    'l3' : Metric('l3', m_l3),
-    'l4' : Metric('l4', m_l4),
-    'lcv' : Metric('lcv', m_lcv),
-    'lskewness' : Metric('lskewness', m_lskewness),
-    'lkurtosis' : Metric('lkurtosis', m_lkurtosis),
-    '90m10' : Metric('90m10', m_90m10),
-    '95m05' : Metric('95m05', m_95m05),
-    'p01' : Metric('p01', m_p01),
-    'p05' : Metric('p05', m_p05),
-    'p10' : Metric('p10', m_p10),
-    'p20' : Metric('p20', m_p20),
-    'p25' : Metric('p25', m_p25),
-    'p30' : Metric('p30', m_p30),
-    'p40' : Metric('p40', m_p40),
-    'p50' : Metric('p50', m_p50),
-    'p60' : Metric('p60', m_p60),
-    'p70' : Metric('p70', m_p70),
-    'p75' : Metric('p75', m_p75),
-    'p80' : Metric('p80', m_p80),
-    'p90' : Metric('p90', m_p90),
-    'p95' : Metric('p95', m_p95),
-    'p99' : Metric('p99', m_p99),
-#    'allcover' : Metric('allcover', m_allcover),
-    'profilearea' : Metric('profilearea', m_profilearea),
-}
+    return merged
diff --git a/src/silvimetric/resources/metrics/__init__.py b/src/silvimetric/resources/metrics/__init__.py
new file mode 100644
index 0000000..c53e3e0
--- /dev/null
+++ b/src/silvimetric/resources/metrics/__init__.py
@@ -0,0 +1,13 @@
+from ..metric import Metric
+from .percentiles import percentiles
+from .l_moments import l_moments
+from .stats import statistics
+from .p_moments import product_moments
+
+#TODO make each one of these have a version with the NumberOfReturns>2 filter
+
+grid_metrics: dict[str, Metric] = dict(percentiles | l_moments | statistics |
+        product_moments)
+
+all_metrics: dict[str, Metric] = dict(percentiles | l_moments | statistics |
+        product_moments)
\ No newline at end of file
diff --git a/src/silvimetric/resources/metrics/aad.py b/src/silvimetric/resources/metrics/aad.py
new file mode 100644
index 0000000..5ef2c6c
--- /dev/null
+++ b/src/silvimetric/resources/metrics/aad.py
@@ -0,0 +1,24 @@
+import numpy as np
+from scipy import stats
+from ..metric import Metric
+from .p_moments import product_moments
+
+def m_aad(data, *args):
+    mean = args[0]
+    return np.mean(np.absolute(data - mean))
+
+def m_madmedian(data, *args):
+    return stats.median_abs_deviation(data)
+
+def m_madmean(data, *args):
+    return stats.median_abs_deviation(data, center=np.mean)
+
+def m_madmode(data):
+    def mode_center(data, axis):
+        return stats.mode(data, axis=axis).mode
+    return stats.median_abs_deviation(data, center=mode_center)
+
+aad: dict[str, Metric] = {}
+aad['aad'] = Metric('aad', np.float32, m_aad, product_moments['mean'])
+aad['madmedian'] = Metric('madmedian', np.float32, m_madmedian)
+aad['madmode'] = Metric('madmode', np.float32, m_madmode)
\ No newline at end of file
diff --git a/src/silvimetric/resources/metrics/filters.py b/src/silvimetric/resources/metrics/filters.py
new file mode 100644
index 0000000..51ca10c
--- /dev/null
+++ b/src/silvimetric/resources/metrics/filters.py
@@ -0,0 +1,4 @@
+from numpy import np
+
+def f_2plus(data):
+    return data[data['HeightAboveGround'] > 2]
\ No newline at end of file
diff --git a/src/silvimetric/resources/metrics/l_moments.py b/src/silvimetric/resources/metrics/l_moments.py
new file mode 100644
index 0000000..7d927ec
--- /dev/null
+++ b/src/silvimetric/resources/metrics/l_moments.py
@@ -0,0 +1,142 @@
+from ..metric import Metric
+import numpy as np
+from .p_moments import mean
+from lmoments3 import lmom_ratios
+from lmoments3 import distr
+
+import warnings
+# suppress warnings from dividing by 0, these are handled in the metric creation
+warnings.filterwarnings(
+    action='ignore',
+    category=RuntimeWarning,
+    module='lmoments3'
+)
+
+def lmom4(data):
+    n = data.count()
+    try:
+        paras = distr.gam.lmom_fit(data)
+        if n < 1:
+            return [np.nan, np.nan, np.nan, np.nan]
+        elif n < 2:
+            return [data.mean, np.nan, np.nan, np.nan]
+        elif n < 3:
+            return [*distr.gam.lmom(nmom=2, **paras), np.nan, np.nan, np.nan]
+        elif n < 4:
+            return [*distr.gam.lmom(nmom=3, **paras), np.nan, np.nan, np.nan]
+        return [*distr.gam.lmom(nmom=4, **paras)]
+    except:
+        return [data.mean(), np.nan, np.nan, np.nan]
+
+# L1 is same as mean...compute using np.mean for speed
+def m_l1(data, *args):
+    return args[0][0]
+
+def m_l2(data, *args):
+    return args[0][1]
+
+def m_l3(data, *args):
+    return args[0][2]
+
+def m_l4(data, *args):
+    return args[0][3]
+
+def m_lcv(data, *args):
+    l: tuple[float, float, float, float] = args[0]
+
+    try:
+        return l[1] / l[0]
+    except ZeroDivisionError as e:
+        return np.nan
+
+def m_lskewness(data, *args):
+    l: tuple[float, float, float, float] = args[0]
+    try:
+        return l[2] / l[1]
+    except ZeroDivisionError as e:
+        return np.nan
+
+def m_lkurtosis(data, *args):
+    l: tuple[float, float, float, float] = args[0]
+    try:
+        return l[3] / l[1]
+    except ZeroDivisionError as e:
+        return np.nan
+
+# intermediate metric, not intended for insertion into db
+l_mom_base = Metric('lmombase', object, lmom4, [])
+
+
+l1 = Metric('l1', np.float32, m_l1, [l_mom_base])
+l2 = Metric('l2', np.float32, m_l2, [l_mom_base])
+l3 = Metric('l3', np.float32, m_l3, [l_mom_base])
+l4 = Metric('l4', np.float32, m_l4, [l_mom_base])
+lcv = Metric('lcv', np.float32, m_lcv, [l_mom_base])
+lskewness = Metric('lskewness', np.float32, m_lskewness,
+    [l_mom_base])
+lkurtosis = Metric('lkurtosis', np.float32, m_lkurtosis,
+    [l_mom_base])
+
+l_moments: dict[str, Metric] = {
+    'l1': l1,
+    'l2': l2,
+    'l3': l3,
+    'l4': l4,
+    'lcv': lcv,
+    'lskewness': lskewness,
+    'lkurtosis': lkurtosis,
+}
+
+# adapted from: https://xiaoganghe.github.io/python-climate-visuals/chapters/data-analytics/scipy-basic.html
+# added b3 and l4 from:
+#   Vogel, R. M. and Fennessey, N. M.: L moment diagrams should replace product moment diagrams,
+#   Water Resour. Res., 29, 1745–1752, https://doi.org/10.1029/93WR00341, 1993.
+# def lmom4_og(data, *args):
+#     # lmom4 returns the first four L-moments of data
+#     # data is the 1-d array
+#     # n is the total number of points in data, j is the j_th point
+#     #
+#     # j range in for loops starts with 1 so we need to subtract 1 for all b# equations
+
+#     data = data.values
+#     n = len(data)
+#     # sort in descending order
+#     data = np.sort(data.reshape(n))[::-1]
+
+#     b0 = args[0]
+#     l1: float = b0
+
+#     # cannot compute l-moments greater than the number of points
+#     if n < 2:
+#         l2 = np.nan
+#     else:
+#         b1 = np.array([(n - j - 1) * data[j] / n / (n - 1)
+#                     for j in range(n)]).sum()
+#         l2: float = 2 * b1 - b0
+
+#     if n < 3:
+#         l3 = np.nan
+#     else:
+#         b2 = np.array([(n - j - 1) * (n - j - 2) * data[j] / n / (n - 1) / (n - 2)
+#                     for j in range(n - 1)]).sum()
+#         l3: float = 6 * (b2 - b1) + b0
+
+#     if n < 4:
+#         l4 = np.nan
+#     else:
+#         b3 = np.array([(n - j - 1) * (n - j - 2) * (n - j - 3) * data[j] / n / (n - 1) / (n - 2) / (n - 3)
+#                     for j in range(n - 2)]).sum()
+#         l4: float = 20 * b3 - 30 * b2 + 12 * b1 - b0
+
+
+#     return l1, l2, l3, l4
+
+# TODO compute L-moments. These are done separately because we only add
+# a single element to TileDB. This is very inefficient since we have to
+# compute all L-moments at once. Ideally, we would have a single metric
+# function that returns an array with 7 values
+
+# added code to compute first 4 l-moments in lmom4.py. There is a package,
+# lmoments3 that can compute the same values but I wanted to avoid the
+# package as it has some IBM copyright requirements (need to include copyright
+# statement with derived works)
\ No newline at end of file
diff --git a/src/silvimetric/resources/metrics/p_moments.py b/src/silvimetric/resources/metrics/p_moments.py
new file mode 100644
index 0000000..2913c97
--- /dev/null
+++ b/src/silvimetric/resources/metrics/p_moments.py
@@ -0,0 +1,31 @@
+import numpy as np
+from scipy.stats import moment
+
+from ..metric import Metric
+
+def m_moments(data, *args):
+    mean = args[0]
+    return moment(data, center=mean, order=[2,3,4], nan_policy='omit').tolist()
+
+def m_mean(data, *args):
+    return np.mean(data)
+
+def m_variance(data, *args):
+    return args[0][0]
+
+def m_skewness(data, *args):
+    return args[0][1]
+
+def m_kurtosis(data, *args):
+    return args[0][2]
+
+mean = Metric(name='mean', dtype=np.float32, method=m_mean)
+moment_base = Metric(name='moment_base', dtype=object, method=m_moments, dependencies=[mean])
+variance = Metric(name='variance', dtype=np.float32, method=m_variance, dependencies=[moment_base])
+skewness = Metric(name='skewness', dtype=np.float32, method=m_skewness, dependencies=[moment_base])
+kurtosis = Metric(name='kurtosis', dtype=np.float32, method=m_kurtosis, dependencies=[moment_base])
+
+product_moments: dict[str, Metric] = dict(mean=mean,
+    variance=variance,
+    skewness=skewness,
+    kurtosis=kurtosis)
\ No newline at end of file
diff --git a/src/silvimetric/resources/metrics/percentiles.py b/src/silvimetric/resources/metrics/percentiles.py
new file mode 100644
index 0000000..0ae46dc
--- /dev/null
+++ b/src/silvimetric/resources/metrics/percentiles.py
@@ -0,0 +1,78 @@
+import numpy as np
+import pandas as pd
+from ..metric import Metric
+
+def percentile_base(data: pd.DataFrame):
+    return np.percentile(data, range(100)).tolist()
+
+def m_p01(data: pd.DataFrame, *args):
+    return args[0][1]
+
+def m_p05(data: pd.DataFrame, *args):
+    return args[0][5]
+
+def m_p10(data: pd.DataFrame, *args):
+    return args[0][10]
+
+def m_p20(data: pd.DataFrame, *args):
+    return args[0][20]
+
+def m_p25(data: pd.DataFrame, *args):
+    return args[0][25]
+
+def m_p30(data: pd.DataFrame, *args):
+    return args[0][30]
+
+def m_p40(data: pd.DataFrame, *args):
+    return args[0][40]
+
+def m_p50(data: pd.DataFrame, *args):
+    return args[0][50]
+
+def m_p60(data: pd.DataFrame, *args):
+    return args[0][60]
+
+def m_p70(data: pd.DataFrame, *args):
+    return args[0][70]
+
+def m_p75(data: pd.DataFrame, *args):
+    return args[0][75]
+
+def m_p80(data: pd.DataFrame, *args):
+    return args[0][80]
+
+def m_p90(data: pd.DataFrame, *args):
+    return args[0][90]
+
+def m_p95(data: pd.DataFrame, *args):
+    return args[0][95]
+
+def m_p99(data: pd.DataFrame, *args):
+    return args[0][99]
+
+def m_90m10(data, *args):
+    return args[0][90] - args[0][10]
+
+def m_95m05(data, *args):
+    return args[0][95] - args[0][5]
+
+pct_base = Metric('pct_base', object, percentile_base)
+
+percentiles: dict[str, Metric] = { }
+percentiles['p01'] = Metric('p01', np.float32, m_p01, [pct_base])
+percentiles['p05'] = Metric('p05', np.float32, m_p05, [pct_base])
+percentiles['p10'] = Metric('p10', np.float32, m_p10, [pct_base])
+percentiles['p20'] = Metric('p20', np.float32, m_p20, [pct_base])
+percentiles['p25'] = Metric('p25', np.float32, m_p25, [pct_base])
+percentiles['p30'] = Metric('p30', np.float32, m_p30, [pct_base])
+percentiles['p40'] = Metric('p40', np.float32, m_p40, [pct_base])
+percentiles['p50'] = Metric('p50', np.float32, m_p50, [pct_base])
+percentiles['p60'] = Metric('p60', np.float32, m_p60, [pct_base])
+percentiles['p70'] = Metric('p70', np.float32, m_p70, [pct_base])
+percentiles['p75'] = Metric('p75', np.float32, m_p75, [pct_base])
+percentiles['p80'] = Metric('p80', np.float32, m_p80, [pct_base])
+percentiles['p90'] = Metric('p90', np.float32, m_p90, [pct_base])
+percentiles['p95'] = Metric('p95', np.float32, m_p95, [pct_base])
+percentiles['p99'] = Metric('p99', np.float32, m_p99, [pct_base])
+percentiles['90m10'] = Metric('90m10', np.float32, m_90m10, [pct_base])
+percentiles['95m05'] = Metric('95m05', np.float32, m_95m05, [pct_base])
\ No newline at end of file
diff --git a/src/silvimetric/resources/metrics/stats.py b/src/silvimetric/resources/metrics/stats.py
new file mode 100644
index 0000000..98cc968
--- /dev/null
+++ b/src/silvimetric/resources/metrics/stats.py
@@ -0,0 +1,128 @@
+import numpy as np
+from scipy import stats
+from ..metric import Metric
+from .p_moments import mean
+from .percentiles import pct_base
+
+import warnings
+# suppress warnings from dividing by 0, these are handled in the metric creation
+warnings.filterwarnings(
+    action='ignore',
+    category=RuntimeWarning,
+    module='lmoments3'
+)
+
+def m_mode(data):
+    u, c = np.unique(data, return_counts=True)
+    i = np.where(c == c.max())
+    v = u[i[0][0]]
+    return v
+
+def m_median(data, *args):
+    return np.median(data)
+
+def m_min(data, *args):
+    return np.min(data)
+
+def m_max(data, *args):
+    return np.max(data)
+
+def m_stddev(data, *args):
+    return np.std(data)
+
+def m_cv(data, *args):
+    stddev, mean = args
+    if mean == 0:
+        return np.nan
+    return stddev / mean
+
+# TODO check performance of other methods
+def m_abovemean(data, *args):
+    mean = args[0]
+    l = len(data)
+    if l == 0:
+        return np.nan
+    return (data > mean).sum() / l
+
+# TODO check performance of other methods
+def m_abovemode(data, *args):
+    mode = args[0]
+    l = len(data)
+    if l == 0:
+        return np.nan
+    return (data > mode).sum() / l
+
+def m_iq(data):
+    return stats.iqr(data)
+
+def m_crr(data, *args):
+    mean, minimum, maximum = args
+    den = (maximum - minimum)
+    if den == 0:
+        return np.nan
+    return (mean - minimum) / den
+
+def m_sqmean(data):
+    return np.sqrt(np.mean(np.square(data)))
+
+def m_cumean(data):
+    return np.cbrt(np.mean(np.power(np.absolute(data), 3)))
+
+def m_profilearea(data, *args):
+    # sanity check...must have valid heights/elevations
+    dmax, dmin , p = args
+    if dmax <= 0:
+        return -9999.0
+
+    # p = np.percentile(data, range(1, 100))
+    p0 = max(dmin, 0.0)
+
+    # second sanity check...99th percentile must be > 0
+    p99 = p[99]
+    if p99 > 0.0:
+        # compute area under normalized percentile height curve using composite trapeziod rule
+        pcts = np.array(p[:98])
+        areas = pcts * 2 / p99
+        pa = p0/p99 + areas.sum() + 1
+
+        return pa * 0.5
+    else:
+        return -9999.0
+
+# # TODO example for cover using all returns and a height threshold
+# # the threshold must be a parameter and not hardcoded
+# def m_cover(data):
+#     threshold = 2
+#     return (data > threshold).sum() / len(data)
+
+mode = Metric('mode', np.float32, m_mode)
+median = Metric('median', np.float32, m_median)
+# TODO better names?
+sm_min = Metric('min', np.float32, m_min)
+sm_max = Metric('max', np.float32, m_max)
+stddev = Metric('stddev', np.float32, m_stddev)
+cv = Metric('cv', np.float32, m_cv, [ mean, stddev ])
+abovemean = Metric('abovemean', np.float32, m_abovemean, [ mean ])
+abovemode = Metric('abovemode', np.float32, m_abovemode, [ mode ])
+iq = Metric('iq', np.float32, m_iq)
+crr = Metric('crr', np.float32, m_crr, [ mean, sm_min, sm_max ])
+sqmean = Metric('sqmean', np.float32, m_sqmean)
+cumean = Metric('cumean', np.float32, m_cumean)
+profilearea = Metric('profilearea', np.float32, m_profilearea, [ sm_max, sm_min, pct_base ])
+
+statistics: dict[str, Metric] = dict(
+    mode=mode,
+    median=median,
+    min=sm_min,
+    max=sm_max,
+    stddev=stddev,
+    cv=cv,
+    abovemean=abovemean,
+    abovemode=abovemode,
+    iq=iq,
+    crr=crr,
+    sqmean=sqmean,
+    cumean=cumean,
+    profilearea=profilearea,
+)
+# statistics['cover'] = Metric('cover', np.float32, m_cover)
diff --git a/src/silvimetric/resources/storage.py b/src/silvimetric/resources/storage.py
index 798349b..60eb3eb 100644
--- a/src/silvimetric/resources/storage.py
+++ b/src/silvimetric/resources/storage.py
@@ -1,25 +1,23 @@
+import os
 import tiledb
 import numpy as np
 from datetime import datetime
 import pathlib
 import contextlib
 import json
-from typing import Any
+import urllib
+from typing import Generator
 
 from math import floor
 
-from .config import StorageConfig
+from .config import StorageConfig, ShatterConfig
 from .metric import Metric, Attribute
-from ..resources import Bounds
+from .bounds import Bounds
 
 class Storage:
     """ Handles storage of shattered data in a TileDB Database. """
 
     def __init__(self, config: StorageConfig, ctx:tiledb.Ctx=None):
-        # if not ctx:
-        #     self.ctx = tiledb.default_ctx()
-        # else:
-        #     self.ctx = ctx
 
         if not tiledb.object_type(config.tdb_dir) == "array":
             raise Exception(f"Given database directory '{config.tdb_dir}' does not exist")
@@ -61,6 +59,9 @@ def create( config:StorageConfig, ctx:tiledb.Ctx=None):
         if not ctx:
             ctx = tiledb.default_ctx()
 
+        # adjust cell bounds if necessary
+        config.root.adjust_to_cell_lines(config.resolution)
+
         # dims = { d['name']: d['dtype'] for d in pdal.dimensions if d['name'] in config.attrs }
         xi = floor((config.root.maxx - config.root.minx) / float(config.resolution))
         yi = floor((config.root.maxy - config.root.miny) / float(config.resolution))
@@ -72,13 +73,22 @@ def create( config:StorageConfig, ctx:tiledb.Ctx=None):
         count_att = tiledb.Attr(name="count", dtype=np.int32)
         dim_atts = [attr.schema() for attr in config.attrs]
 
-        metric_atts = [m.schema(a) for m in config.metrics for a in config.attrs]
+        metric_atts = [m.schema(a) for m in config.metrics for a in config.attrs if a in m.attributes or not m.attributes]
+
+        # Check that all attributes required for metric usage are available
+        att_list = [a.name for a in config.attrs]
+        required_atts = [d.name for m in config.metrics for d in m.dependencies
+                if isinstance(d, Attribute)]
+        for ra in required_atts:
+            if ra not in att_list:
+                raise ValueError(f'Missing required dependency, {ra}.')
 
         # allows_duplicates lets us insert multiple values into each cell,
         # with each value representing a set of values from a shatter process
         # https://docs.tiledb.com/main/how-to/performance/performance-tips/summary-of-factors#allows-duplicates
         schema = tiledb.ArraySchema(domain=domain, sparse=True,
-            attrs=[count_att, *dim_atts, *metric_atts], allows_duplicates=True)
+            attrs=[count_att, *dim_atts, *metric_atts], allows_duplicates=True,
+            capacity=1000)
         schema.check()
 
         tiledb.SparseArray.create(config.tdb_dir, schema)
@@ -87,23 +97,17 @@ def create( config:StorageConfig, ctx:tiledb.Ctx=None):
             a.meta['config'] = meta
 
         s = Storage(config, ctx)
+        s.saveConfig()
 
         return s
 
     @staticmethod
     def from_db(tdb_dir: str):
         """
-        Create Storage class from previously created storage path
-
-        Parameters
-        ----------
-        tdb_dir : str
-            Path to storage directory
+        Create Storage object from information stored in a database.
 
-        Returns
-        -------
-        Storage
-            Return the generated storage
+        :param tdb_dir: TileDB database directory.
+        :return: Returns the derived storage.
         """
         with tiledb.open(tdb_dir, 'r') as a:
             s = a.meta['config']
@@ -114,19 +118,15 @@ def from_db(tdb_dir: str):
     def saveConfig(self) -> None:
         """
         Save StorageConfig to the Database
-
         """
-        with self.open('w') as a:
-            a.meta['config'] = str(self.config)
+        with self.open('w') as w:
+            w.meta['config'] = str(self.config)
 
     def getConfig(self) -> StorageConfig:
         """
-        Get the StorageConfig currently in use by the Storage
+        Get the StorageConfig currently in use by the Storage.
 
-        Returns
-        -------
-        StorageConfig
-            StorageConfig object
+        :return: StorageConfig representing this object.
         """
         with self.open('r') as a:
             s = a.meta['config']
@@ -134,95 +134,79 @@ def getConfig(self) -> StorageConfig:
 
         return config
 
-    def getMetadata(self, key: str, default=None) -> str:
+    def getMetadata(self, key: str, timestamp: int) -> str:
         """
-        Return metadata at given key
-
-        Parameters
-        ----------
-        key : str
-            Metadata key
+        Return metadata at given key.
 
-        Returns
-        -------
-        str
-            Metadata value found in storage
+        :param key: Key to look for in metadata.
+        :param timestamp: Time stamp for querying database.
+        :return: Metadata value found in storage.
         """
-        with self.open('r') as r:
-            try:
-                val = r.meta[key]
-            except KeyError as e:
-                if default is not None:
-                    return default
-                raise(e)
+
+        with self.open('r', (timestamp, timestamp)) as r:
+            val = r.meta[f'{key}']
 
         return val
 
-    def saveMetadata(self, key: str, data: Any) -> None:
+
+    def saveMetadata(self, key: str, data: str, timestamp: int) -> None:
         """
-        Save metadata to storage
+        Save metadata to storage.
 
-        Parameters
-        ----------
-        key : str
-            Metadata key
-        data : any
-            Metadata value. Must be translatable to and from string.
+        :param key: Metadata key to save to.
+        :param data: Data to save to metadata.
         """
-        with self.open('w') as w:
-            w.meta[key] = data
+        with self.open('w', (timestamp, timestamp)) as w:
+            w.meta[f'{key}'] = data
+        return
 
     def getAttributes(self) -> list[Attribute]:
         """
-        Return list of attribute names from storage config
+        Find list of attribute names from storage config.
 
-        Returns
-        -------
-        list[str]
-            List of attribute names
+        :return: List of attribute names.
         """
         return self.getConfig().attrs
 
     def getMetrics(self) -> list[Metric]:
         """
-        Return List of metric names from storage config
+        Find List of metric names from storage config
 
-        Returns
-        -------
-        list[str]
-            List of metric names
+        :return: List of metric names.
         """
         return self.getConfig().metrics
 
+    def getDerivedNames(self) -> list[str]:
+        # if no attributes are set in the metric, use all
+        return [m.entry_name(a.name) for m in self.config.metrics
+                for a in self.config.attrs if not m.attributes
+                or a.name in [ma.name for ma in m.attributes]]
+
     @contextlib.contextmanager
-    def open(self, mode:str='r', timestamp=None) -> tiledb.SparseArray:
+    def open(self, mode:str='r', timestamp=None) -> Generator[tiledb.SparseArray, None, None]:
+        """
+        Open stream for TileDB database in given mode and at given timestamp.
+
+        :param mode: Mode to open TileDB stream in. Valid options are
+            'w', 'r', 'm', 'd'., defaults to 'r'.
+        :param timestamp: Timestamp to open database at., defaults to None.
+        :raises Exception: Incorrect Mode was given, only valid modes are 'w' and 'r'.
+        :raises Exception: Path exists and is not a TileDB array.
+        :raises Exception: Path does not exist.
+        :yield: TileDB array context manager.
         """
-        Open either a read or write stream for TileDB database
 
-        Parameters
-        ----------
-        mode : str, optional
-            Stream mode. Valid options are 'r' and 'w', by default 'r'
+        # tiledb and dask have bad interaction with opening an array if
+        # other threads present
 
-        Raises
-        ------
-        Exception
-            Incorrect Mode was given, only valid modes are 'w' and 'r'
-        Exception
-            Path exists and is not a TileDB array
-        Exception
-            Path does not exist
-        """
         if tiledb.object_type(self.config.tdb_dir) == "array":
-            if mode == 'w':
-                tdb = tiledb.open(self.config.tdb_dir, 'w', timestamp=timestamp)
-            elif mode == 'r':
-                tdb = tiledb.open(self.config.tdb_dir, 'r', timestamp=timestamp)
+            if mode in ['w', 'r', 'd', 'm']:
+                tdb = tiledb.open(self.config.tdb_dir, mode, timestamp=timestamp)
             else:
                 raise Exception(f"Given open mode '{mode}' is not valid")
         elif pathlib.Path(self.config.tdb_dir).exists():
             raise Exception(f"Path {self.config.tdb_dir} already exists and is not" +
-                            " initialized for TileDB access.")
+                    " initialized for TileDB access.")
         else:
             raise Exception(f"Path {self.config.tdb_dir} does not exist")
 
@@ -231,99 +215,143 @@ def open(self, mode:str='r', timestamp=None) -> tiledb.SparseArray:
         finally:
             tdb.close()
 
+    def reserve_time_slot(self) -> int:
+        """
+        Increment time slot in database and reserve that spot for a new
+        shatter process.
+
+        :param config: Shatter config will be written as metadata to reserve
+            time slot.
+
+        :return: Time slot.
+        """
+        with tiledb.open(self.config.tdb_dir, 'r') as r:
+            latest = json.loads(r.meta['config'])
+
+        time = latest['next_time_slot']
+        latest['next_time_slot'] = time + 1
+
+        with tiledb.open(self.config.tdb_dir, 'w') as w:
+            w.meta['config'] = json.dumps(latest)
+
+        return time
+
     def get_history(self, start_time: datetime, end_time: datetime,
                     bounds: Bounds, name:str=None):
+        """
+        Retrieve history of the database at current point in time.
 
-        from ..resources import ShatterConfig
+        :param start_time: Query parameter, starting datetime of process.
+        :param end_time: Query parameter, ending datetime of process.
+        :param bounds: Query parameter, bounds to query by.
+        :param name: Query paramter, shatter process uuid., by default None
+        :return: Returns list of array fragments that meet query parameters.
+        """
 
         af = tiledb.array_fragments(self.config.tdb_dir, True)
         m = [ ]
         for idx in range(len(af)):
 
-            # if processes are too short, begin and end times can be the same
-            # so it's necessary to adjust time bounds so we can grab information
+            time_range =  af[idx].timestamp_range
+            # all shatter processes should be input as a point in time, eg (1,1)
+            if isinstance(time_range, tuple):
+                time_range = time_range[0]
 
-            if idx == 0:
-                begin = 0
+            try:
+                s_str = self.getMetadata('shatter', time_range)
+            except KeyError:
+                continue
+            s = ShatterConfig.from_string(s_str)
+            if s.bounds.disjoint(bounds):
+                continue
+
+            # filter name
+            if name is not None and name != s.name:
+                continue
+
+            # filter dates
+            if isinstance(s.date, tuple) and len(s.date) == 2:
+                if s.date[1] < start_time or s.date[0] > end_time:
+                    continue
+            elif isinstance(s.date, tuple) and len(s.date) == 1:
+                if s.date[0] < start_time or s.date[0] > end_time:
+                    continue
             else:
-                begin = af[idx-1].timestamp_range[1]
-
-            end = af[idx].timestamp_range[1]
-
-            with self.open('r', (begin, end)) as r:
-                try:
-                    if not bool(r.meta['shatter']):
-                        continue
-
-                    s = ShatterConfig.from_string(r.meta['shatter'])
-
-                    # filter bounds
-                    if s.bounds.disjoint(bounds):
-                        continue
-
-                    # filter name
-                    if name is not None and name != s.name:
-                        continue
-
-                    # filter dates
-                    if isinstance(s.date, tuple) and len(s.date) == 2:
-                        if s.date[1] < start_time or s.date[0] > end_time:
-                            continue
-                    elif isinstance(s.date, tuple) and len(s.date) == 1:
-                        if s.date[0] < start_time or s.date[0] > end_time:
-                            continue
-                    else:
-                        if s.date < start_time or s.date > end_time:
-                            continue
-
-                    m.append(json.loads(r.meta['shatter']))
-                except KeyError:
+                if s.date < start_time or s.date > end_time:
                     continue
 
+            m.append(s.to_json())
+
         return m
 
-    #TODO what are we reading? queries are probably going to be specific
-    def read(self, xs: np.ndarray, ys: np.ndarray) -> np.ndarray:
+    def mbrs(self, proc_num: int):
         """
-        Read from the Database
-        Parameters
-        ----------
-        xs : np.ndarray
-            X index
-        ys : np.ndarray
-            Y index
+        Get minimum bounding rectangle of a given shatter process. If this process
+        has been finished and consolidated the mbr will be much less granulated
+        than if the fragments are still intact. Mbrs are represented as tuples
+        in the form of ((minx, maxx), (miny, maxy))
 
-        Returns
-        -------
-        np.ndarray
-            Items found at the indicated cell
+        :param proc_num: Process number or time slot of the shatter process.
+        :return: Returns mbrs that match the given process number.
         """
-        with self.open('r') as tdb:
-            data = tdb[xs, ys]
-        return data
+        af_all = self.get_fragments_by_time(proc_num)
+        mbrs_list = tuple(mbrs for af in af_all for mbrs in af.mbrs)
+        mbrs = tuple(tuple(tuple(a.item() for a in mb) for mb in m) for m in mbrs_list)
+        return mbrs
 
-    def write(self, xs: np.ndarray, ys: np.ndarray, data: np.ndarray) -> None:
+    def get_fragments_by_time(self, proc_num: int) -> list[tiledb.FragmentInfo]:
         """
-        Write data to TileDB database
+        Get TileDB array fragments from the time slot specified.
 
-        Parameters
-        ----------
-        xs : np.ndarray
-            X cell indices
-        ys : np.ndarray
-            Y cell indices
-        data : np.ndarray
-            Numpy object of data values for attributes in each index pairing
+        :param proc_num: Requested time slot.
+        :return: Array fragments from time slot.
+        """
+        af = tiledb.array_fragments(self.config.tdb_dir, include_mbrs=True)
+        return [a for a in af if a.timestamp_range == (proc_num, proc_num)]
+
+    def delete(self, proc_num: int) -> ShatterConfig:
+        """
+        Delete Shatter process and all associated data from database.
+
+        :param proc_num: Shatter process time slot
+        :return: Config of deleted Shatter process
         """
 
-        with self.open('w') as tdb:
-            # data = {k: v.astype(np.dtype(v.dtype)) for k,v in data.items()}
+        self.config.log.debug(f'Deleting time slot {proc_num}...')
+        with self.open('r', (proc_num, proc_num)) as r:
+            sh_cfg: ShatterConfig = ShatterConfig.from_string(r.meta['shatter'])
+            sh_cfg.mbr = ()
+            sh_cfg.finished = False
+
+        self.config.log.debug('Deleting fragments...')
+        tiledb.Array.delete_fragments(self.config.tdb_dir,
+            timestamp_start=proc_num, timestamp_end=proc_num)
+        self.config.log.debug('Rewriting config.')
+        with self.open('w', (proc_num, proc_num)) as w:
+            w.meta['shatter'] = json.dumps(sh_cfg.to_json())
+        return sh_cfg
+
+    def consolidate_shatter(self, proc_num: int, retries=0) -> None:
+        """
+        Consolidate the fragments from a shatter process into one fragment.
+        This makes the database perform better, but reduces the granularity of
+        time traveling.
 
-            # if self.config.app.debug:
-            #     tiledb.stats_reset()
-            #     tiledb.stats_enable()
-            tdb[xs, ys] = data
+        :param proc_num: Time slot associated with shatter process.
+        """
+        try:
+            afs = self.get_fragments_by_time(proc_num)
+            uris = [ os.path.split(urllib.parse.urlparse(f.uri).path)[-1] for f in afs ]
+            tiledb.consolidate(self.config.tdb_dir, fragment_uris=uris)
+            c = tiledb.Config({"sm.vacuum.mode": "fragments"})
+            tiledb.vacuum(self.config.tdb_dir, c)
+            self.config.log.info(f"Consolidated time slot {proc_num}.")
+        except Exception as e:
+            if retries >= 3:
+                self.config.log.warning("Failed to consolidate time slot "
+                        f"{proc_num} {retries} time(s). Stopping.")
+                raise e
+            self.config.log.warning("Failed to consolidate time slot "
+                    f"{proc_num} {retries+1} time. Retrying...")
+            self.consolidate_shatter(proc_num, retries+1)
 
-            # if self.config.app.debug:
-            #     tiledb.stats_dump()
-            #     tiledb.stats_reset()


Test patch:
diff --git a/tests/conftest.py b/tests/conftest.py
index 306c52a..9849d70 100644
--- a/tests/conftest.py
+++ b/tests/conftest.py
@@ -1,181 +1,136 @@
 import pytest
 import os
-import dask
 import pdal
-from uuid import uuid4
+import copy
+
 from datetime import datetime
+from typing import Generator
 
-from silvimetric.resources import Extents, Bounds, Metrics, Attribute, Storage, Log, Data
-from silvimetric.resources.config import ShatterConfig, StorageConfig, ApplicationConfig
+from silvimetric import Extents, Bounds, Attribute, Storage
+from silvimetric import grid_metrics
+from silvimetric import Log, Metric, ShatterConfig, StorageConfig
+from silvimetric import ApplicationConfig, ExtractConfig
 from silvimetric import __version__ as svversion
 
-
-@pytest.fixture(scope="session", autouse=True)
-def configure_dask():
-    dask.config.set(scheduler="single-threaded")
-
-@pytest.fixture(scope='function')
-def tdb_filepath(tmp_path_factory) -> str:
-    path = tmp_path_factory.mktemp("test_tdb")
-    yield os.path.abspath(path)
-
-@pytest.fixture(scope="function")
-def threaded_dask():
-    dask.config.set(scheduler="threads")
+# pull together fixtures
+pytest_plugins=[
+    'fixtures.shatter_fixtures', 'fixtures.extract_fixtures',
+    'fixtures.command_fixtures', 'fixtures.chunk_fixtures',
+    'fixtures.western_fixtures', 'fixtures.data_fixtures',
+    'fixtures.cli_fixtures', 'fixtures.fusion_fixtures',
+    'fixtures.metric_fixtures', 'fixtures.dask_fixtures'
+]
 
 @pytest.fixture(scope='function')
-def storage_config(tdb_filepath, bounds, resolution, crs, attrs, metrics):
-    log = Log(20)
-    yield StorageConfig(tdb_dir = tdb_filepath,
-                        log = log,
-                        crs = crs,
-                        root = bounds,
-                        resolution = resolution,
-                        attrs = attrs,
-                        metrics = metrics,
-                        version = svversion)
-
-@pytest.fixture(scope="function")
-def storage(storage_config) -> Storage:
-    yield Storage.create(storage_config)
+def tdb_filepath(storage_config) -> Generator[str, None, None]:
+    yield storage_config.tdb_dir
 
 @pytest.fixture(scope='function')
-def app_config(tdb_filepath, debug=True):
+def app_config(tdb_filepath, debug=True) -> Generator[ApplicationConfig, None, None]:
     log = Log(20) # INFO
     app = ApplicationConfig(tdb_dir = tdb_filepath,
                             log = log)
     yield app
 
 @pytest.fixture(scope='function')
-def shatter_config(tdb_filepath, copc_filepath, storage_config, bounds, app_config, storage, date):
-    log = Log(20) # INFO
-    s = ShatterConfig(tdb_dir = tdb_filepath,
-                      log = log,
-                      filename = copc_filepath,
-                      attrs = storage_config.attrs,
-                      metrics = storage_config.metrics,
-                      bounds=bounds,
-                      debug = True,
-                      date=date)
-
-    yield s
+def storage_config(tmp_path_factory, bounds, resolution, crs, attrs, metrics) -> Generator[StorageConfig, None, None]:
+    path = tmp_path_factory.mktemp("test_tdb")
+    p = os.path.abspath(path)
+    log = Log('DEBUG')
 
-@pytest.fixture(scope='function')
-def uneven_storage_config(tdb_filepath, bounds, crs, attrs, metrics):
-    log = Log(20)
-    yield StorageConfig(tdb_dir = tdb_filepath,
+    sc =  StorageConfig(tdb_dir = p,
                         log = log,
                         crs = crs,
                         root = bounds,
-                        resolution = 7,
+                        resolution = resolution,
                         attrs = attrs,
                         metrics = metrics,
                         version = svversion)
+    Storage.create(sc)
+    yield sc
 
 @pytest.fixture(scope='function')
-def uneven_storage(uneven_storage_config):
-    yield Storage.create(uneven_storage_config)
+def storage(storage_config):
+    yield Storage.from_db(storage_config.tdb_dir)
 
 @pytest.fixture(scope='function')
-def uneven_shatter_config(tdb_filepath, copc_filepath, uneven_storage_config, storage, date):
-    log = Log(20) # INFO
-    s = ShatterConfig(tdb_dir = tdb_filepath,
+def shatter_config(copc_filepath, storage_config, bounds, date) -> Generator[ShatterConfig, None, None]:
+    log = Log('INFO') # INFO
+    s = ShatterConfig(tdb_dir = storage_config.tdb_dir,
                       log = log,
                       filename = copc_filepath,
-                      attrs = uneven_storage_config.attrs,
-                      metrics = uneven_storage_config.metrics,
+                      attrs = storage_config.attrs,
+                      metrics = storage_config.metrics,
+                      bounds = bounds,
                       debug = True,
-                      date=date)
-    yield s
+                      date = date, tile_size=10)
 
-@pytest.fixture(scope="function")
-def s3_bucket():
-    yield "silvimetric"
+    yield s
 
 @pytest.fixture(scope='function')
-def s3_uri(s3_bucket):
-    uuid = uuid4()
-    yield f"s3://{s3_bucket}/test_silvimetric/{uuid}"
+def extract_config(tif_filepath, metrics, shatter_config, extract_attrs):
+    from silvimetric.commands import shatter
+    tdb_dir = shatter_config.tdb_dir
+    shatter.shatter(shatter_config)
+    log = Log(20)
+    c =  ExtractConfig(tdb_dir = tdb_dir,
+                       log = log,
+                       out_dir = tif_filepath,
+                       attrs = extract_attrs,
+                       metrics = metrics)
+    yield c
 
-@pytest.fixture(scope="function")
-def s3_storage_config(s3_uri, bounds, resolution, crs, attrs, metrics):
-    yield StorageConfig(bounds, crs, resolution, attrs, metrics,
-                        svversion, tdb_dir=s3_uri)
 
 @pytest.fixture(scope='function')
-def s3_storage(s3_storage_config):
-    import subprocess
-    yield Storage.create(s3_storage_config)
-    subprocess.call(["aws", "s3", "rm", "--recursive", s3_storage_config.tdb_dir])
-
-@pytest.fixture(scope="function")
-def s3_shatter_config(s3_storage, copc_filepath, attrs, metrics, date):
-    config = s3_storage.config
-    yield ShatterConfig(filename=copc_filepath, attrs=attrs, metrics=metrics,
-                        debug=True, tdb_dir=config.tdb_dir, date=date)
-
-@pytest.fixture(scope='session')
-def metrics():
-    yield [Metrics['mean'], Metrics['median']]
-
-@pytest.fixture(scope='session')
-def copc_filepath() -> str:
-    path = os.path.join(os.path.dirname(__file__), "data",
-            "test_data.copc.laz")
-    assert os.path.exists(path)
-    yield os.path.abspath(path)
+def metrics() -> Generator[list[Metric], None, None]:
+    yield [copy.deepcopy(grid_metrics['mean']), copy.deepcopy(grid_metrics['median'])]
 
 @pytest.fixture(scope='function')
-def copc_data(copc_filepath, storage_config) -> Data:
-    d = Data(copc_filepath, storage_config)
-    yield d
-
-@pytest.fixture(scope='class')
-def bounds(minx, maxx, miny, maxy) -> Bounds:
+def bounds(minx, maxx, miny, maxy) -> Generator[Bounds, None, None]:
     b =  Bounds(minx, miny, maxx, maxy)
     yield b
 
-@pytest.fixture(scope='class')
-def extents(resolution, bounds) -> Extents:
+@pytest.fixture(scope='function')
+def extents(resolution, bounds) -> Generator[Extents, None, None]:
     yield Extents(bounds,resolution,bounds)
 
-@pytest.fixture(scope="session")
-def attrs(dims) -> list[str]:
+@pytest.fixture(scope="function")
+def attrs(dims) -> Generator[list[Attribute], None, None]:
     yield [Attribute(a, dims[a]) for a in
            ['Z', 'NumberOfReturns', 'ReturnNumber', 'Intensity']]
 
 @pytest.fixture(scope="session")
-def dims():
+def dims() -> Generator[dict, None, None]:
     yield { d['name']: d['dtype'] for d in pdal.dimensions }
 
-@pytest.fixture(scope='class')
-def resolution() -> int:
+@pytest.fixture(scope='session')
+def resolution() -> Generator[int, None, None]:
     yield 30
 
-@pytest.fixture(scope='class')
-def test_point_count() -> int:
+@pytest.fixture(scope='session')
+def test_point_count() -> Generator[int, None, None]:
     yield 90000
 
-@pytest.fixture(scope='class')
-def minx() -> float:
+@pytest.fixture(scope='session')
+def minx() -> Generator[float, None, None]:
     yield 300
 
-@pytest.fixture(scope='class')
-def miny() -> float:
+@pytest.fixture(scope='session')
+def miny() -> Generator[float, None, None]:
     yield 300
 
-@pytest.fixture(scope='class')
-def maxx() -> float:
+@pytest.fixture(scope='session')
+def maxx() -> Generator[float, None, None]:
     yield 600
 
-@pytest.fixture(scope='class')
-def maxy() -> float:
+@pytest.fixture(scope='session')
+def maxy() -> Generator[float, None, None]:
     yield 600
 
-@pytest.fixture(scope="class")
-def crs():
+@pytest.fixture(scope='session')
+def crs() -> Generator[str, None, None]:
     yield "EPSG:5070"
 
-@pytest.fixture(scope='class')
-def date():
+@pytest.fixture(scope='session')
+def date() -> Generator[datetime, None, None]:
     yield datetime(2011, 1, 1)
\ No newline at end of file
diff --git a/tests/create_test_data.py b/tests/create_test_data.py
index 56b9304..57b4a0a 100644
--- a/tests/create_test_data.py
+++ b/tests/create_test_data.py
@@ -1,6 +1,6 @@
 import pdal
 import numpy as np
-from math import sqrt, ceil
+from math import sqrt, floor
 
 
 # Create a pointcloud dataset that follows expected paramters for data testing
@@ -28,10 +28,10 @@
 diff_maker = 0
 
 x_pos = np.arange(minx, maxx, interval, dtype=np.float32)
-y_pos = np.arange(maxy, miny, -1*interval, dtype=np.float32)
+y_pos = np.arange(miny, maxy, interval, dtype=np.float32)
 # positions = pos[np.where(pos % cell_size != 0)]
 
-alg = lambda y: ceil(y/cell_size) + diff_maker
+alg = lambda y: floor(y/cell_size) + diff_maker
 
 data = np.array([(x, y, alg(y), alg(y), alg(y), alg(y))
                   for x in x_pos for y in y_pos],
diff --git a/tests/data/df_sample b/tests/data/df_sample
new file mode 100644
index 0000000..dad54c6
Binary files /dev/null and b/tests/data/df_sample differ
diff --git a/tests/data/test_data.copc.laz b/tests/data/test_data.copc.laz
index 42eef8a..2323610 100644
Binary files a/tests/data/test_data.copc.laz and b/tests/data/test_data.copc.laz differ
diff --git a/tests/fixtures/chunk_fixtures.py b/tests/fixtures/chunk_fixtures.py
new file mode 100644
index 0000000..e27f052
--- /dev/null
+++ b/tests/fixtures/chunk_fixtures.py
@@ -0,0 +1,26 @@
+import os
+import pytest
+
+from typing import Generator, List
+
+from silvimetric import Extents, Data
+
+@pytest.fixture(scope='function', autouse=True)
+def filtered(copc_data, extents: Extents) -> Generator[List[Extents], None, None]:
+    yield list(extents.chunk(copc_data, 1))
+
+@pytest.fixture(scope='function')
+def unfiltered(extents: Extents) -> Generator[List[Extents], None, None]:
+    yield list(extents.get_leaf_children(30))
+
+@pytest.fixture(scope='session')
+def copc_filepath() -> Generator[str, None, None]:
+    path = os.path.join(os.path.dirname(__file__), "..", "data",
+            "test_data.copc.laz")
+    assert os.path.exists(path)
+    yield os.path.abspath(path)
+
+@pytest.fixture(scope='function')
+def copc_data(copc_filepath, storage_config) -> Generator[Data, None, None]:
+    d = Data(copc_filepath, storage_config)
+    yield d
\ No newline at end of file
diff --git a/tests/fixtures/cli_fixtures.py b/tests/fixtures/cli_fixtures.py
new file mode 100644
index 0000000..d4514e1
--- /dev/null
+++ b/tests/fixtures/cli_fixtures.py
@@ -0,0 +1,15 @@
+import pytest
+from typing import Generator
+from click.testing import CliRunner
+
+from silvimetric.commands import shatter
+
+@pytest.fixture(scope='session')
+def runner():
+    return CliRunner()
+
+@pytest.fixture
+def pre_shatter(shatter_config) -> Generator[int, None, None]:
+    shatter_config.tile_size=1
+    shatter_config.mbr = (((0,4), (0,4)),)
+    yield shatter.shatter(shatter_config)
\ No newline at end of file
diff --git a/tests/fixtures/command_fixtures.py b/tests/fixtures/command_fixtures.py
new file mode 100644
index 0000000..9ef5959
--- /dev/null
+++ b/tests/fixtures/command_fixtures.py
@@ -0,0 +1,36 @@
+import uuid
+import pytest
+
+from time import sleep
+from typing import List, Generator
+from datetime import datetime
+
+from silvimetric.commands import shatter
+from silvimetric import ShatterConfig
+
+@pytest.fixture(scope='function')
+def config_split(shatter_config: ShatterConfig) -> Generator[List[ShatterConfig], None, None]:
+    sc = shatter_config
+    config_split = []
+    day = 1
+    for b in sc.bounds.bisect():
+        day += 1
+        date = datetime(2011, 1, day)
+
+        config_split.append(ShatterConfig(
+            filename=sc.filename,
+            date=date,
+            attrs=sc.attrs,
+            metrics=sc.metrics,
+            bounds=b,
+            name=uuid.uuid4(),
+            tile_size=4,
+            tdb_dir=sc.tdb_dir,
+            log=sc.log
+        ))
+
+    for s in config_split:
+        shatter.shatter(s)
+    sleep(1)
+
+    yield config_split
\ No newline at end of file
diff --git a/tests/fixtures/dask_fixtures.py b/tests/fixtures/dask_fixtures.py
new file mode 100644
index 0000000..0594744
--- /dev/null
+++ b/tests/fixtures/dask_fixtures.py
@@ -0,0 +1,40 @@
+import pytest
+import os
+from typing import Generator
+import dask
+
+@pytest.fixture(scope="session", autouse=True)
+def configure_dask() -> None:
+    dask.config.set(scheduler="single-threaded")
+
+@pytest.fixture(scope="function")
+def threaded_dask() -> Generator[None, None, None]:
+    dask.config.set(scheduler="threads")
+
+    yield None
+    dask.config.set(scheduler="single-threaded")
+
+@pytest.fixture(scope="function")
+def process_dask() -> Generator[None, None, None]:
+    dask.config.set(scheduler="processes")
+
+    yield None
+    dask.config.set(scheduler="single-threaded")
+
+@pytest.fixture(scope='function')
+def dask_proc_client(process_dask) -> Generator[None, None, None]:
+    client = dask.distributed.Client()
+    dask.config.set({'distributed.client': client})
+
+    yield None
+    dask.config.set({'distributed.client': None})
+    client.close()
+
+@pytest.fixture(scope='function')
+def dask_thread_client(threaded_dask) -> Generator[None, None, None]:
+    client = dask.distributed.Client()
+    dask.config.set({'distributed.client': client})
+
+    yield None
+    dask.config.set({'distributed.client': None})
+    client.close()
\ No newline at end of file
diff --git a/tests/fixtures/data_fixtures.py b/tests/fixtures/data_fixtures.py
new file mode 100644
index 0000000..78c0c9a
--- /dev/null
+++ b/tests/fixtures/data_fixtures.py
@@ -0,0 +1,28 @@
+import pytest
+import os
+from typing import Generator
+
+@pytest.fixture(scope='session')
+def no_cell_line_pc() -> Generator[int, None, None]:
+    yield 84100
+
+@pytest.fixture(scope='session')
+def no_cell_line_path() -> Generator[str, None, None]:
+    path = os.path.join(os.path.dirname(__file__), "..", "data",
+            "test_data_2.copc.laz")
+    assert os.path.exists(path)
+    yield os.path.abspath(path)
+
+@pytest.fixture(scope='session')
+def no_cell_line_pipeline() -> Generator[str, None, None]:
+    path = os.path.join(os.path.dirname(__file__), "..", "data",
+            "test_data_2.json")
+    assert os.path.exists(path)
+    yield os.path.abspath(path)
+
+@pytest.fixture(scope='session')
+def autzen_filepath() -> Generator[str, None, None]:
+    path = os.path.join(os.path.dirname(__file__), "..", "data",
+            "autzen-small.copc.laz")
+    assert os.path.exists(path)
+    yield os.path.abspath(path)
\ No newline at end of file
diff --git a/tests/fixtures/extract_fixtures.py b/tests/fixtures/extract_fixtures.py
new file mode 100644
index 0000000..1104d42
--- /dev/null
+++ b/tests/fixtures/extract_fixtures.py
@@ -0,0 +1,39 @@
+import pytest
+import os
+import copy
+import uuid
+from typing import Generator
+
+from silvimetric.commands.shatter import shatter
+from silvimetric import Attribute, ExtractConfig, Log, Bounds
+
+@pytest.fixture(scope='function')
+def tif_filepath(tmp_path_factory) -> Generator[str, None, None]:
+    path = tmp_path_factory.mktemp("test_tifs")
+    yield os.path.abspath(path)
+
+@pytest.fixture(scope='function')
+def extract_attrs(dims)->Generator[list[str], None, None]:
+    yield [Attribute('Z', dtype=dims['Z']), Attribute('Intensity', dtype=dims['Intensity'])]
+
+
+@pytest.fixture(scope='function')
+def multivalue_config(tif_filepath, metric_shatter_config):
+
+    shatter(metric_shatter_config)
+
+    # reset config
+    second_config = copy.deepcopy(metric_shatter_config)
+    second_config.name = uuid.uuid4()
+    second_config.bounds = Bounds(300,300,450,450)
+    second_config.point_count = 0
+    second_config.time_slot = 0
+    second_config.mbr = ()
+
+    shatter(second_config)
+    log = Log(20)
+    tdb_dir = metric_shatter_config.tdb_dir
+    c =  ExtractConfig(tdb_dir = tdb_dir,
+                       log = log,
+                       out_dir = tif_filepath)
+    yield c
diff --git a/tests/fixtures/fake_metrics.py b/tests/fixtures/fake_metrics.py
new file mode 100644
index 0000000..05fd4bb
--- /dev/null
+++ b/tests/fixtures/fake_metrics.py
@@ -0,0 +1,13 @@
+import numpy as np
+from silvimetric import Metric
+
+def count(data):
+    return data.count()
+
+def has_data(data):
+    return data.any()
+
+def metrics() -> list[Metric]:
+    m1 = Metric('count', np.float32, count)
+    m2 = Metric('exists', '?', has_data)
+    return [m1, m2]
\ No newline at end of file
diff --git a/tests/fixtures/fusion_fixtures.py b/tests/fixtures/fusion_fixtures.py
new file mode 100644
index 0000000..5f61c16
--- /dev/null
+++ b/tests/fixtures/fusion_fixtures.py
@@ -0,0 +1,6 @@
+import pytest
+
+@pytest.fixture(scope='function')
+def fusion_data_path():
+
+    pass
\ No newline at end of file
diff --git a/tests/fixtures/metric_fixtures.py b/tests/fixtures/metric_fixtures.py
new file mode 100644
index 0000000..79fbdc5
--- /dev/null
+++ b/tests/fixtures/metric_fixtures.py
@@ -0,0 +1,81 @@
+import os
+import pytest
+from typing import Generator
+import pandas as pd
+import numpy as np
+
+from silvimetric import Log, StorageConfig, ShatterConfig, Storage, Data, Bounds
+from silvimetric import __version__ as svversion
+
+@pytest.fixture(scope='function')
+def autzen_storage(tmp_path_factory: pytest.TempPathFactory) -> Generator[StorageConfig, None, None]:
+    path = tmp_path_factory.mktemp("test_tdb")
+    p = os.path.abspath(path)
+
+    srs = """PROJCS[\"NAD83 / Oregon GIC Lambert (ft)\",
+GEOGCS[\"NAD83\",DATUM[\"North_American_Datum_1983\",SPHEROID[\"GRS 1980\",
+6378137,298.257222101,AUTHORITY[\"EPSG\",\"7019\"]],AUTHORITY[\"EPSG\",\"6269\"]],
+PRIMEM[\"Greenwich\",0,AUTHORITY[\"EPSG\",\"8901\"]],UNIT[\"degree\",
+0.0174532925199433,AUTHORITY[\"EPSG\",\"9122\"]],AUTHORITY[\"EPSG\",\"4269\"]],
+PROJECTION[\"Lambert_Conformal_Conic_2SP\"],PARAMETER[\"latitude_of_origin\",
+41.75],PARAMETER[\"central_meridian\",-120.5],PARAMETER[\"standard_parallel_1\",
+43],PARAMETER[\"standard_parallel_2\",45.5],PARAMETER[\"false_easting\",
+1312335.958],PARAMETER[\"false_northing\",0],UNIT[\"foot\",0.3048,
+AUTHORITY[\"EPSG\",\"9002\"]],AXIS[\"Easting\",EAST],AXIS[\"Northing\",NORTH],
+AUTHORITY[\"EPSG\",\"2992\"]]"""
+    b = Bounds(635579.2,848884.83,639003.73,853536.21)
+    sc = StorageConfig(b, srs, 10, tdb_dir=p)
+    Storage.create(sc)
+    yield sc
+
+@pytest.fixture(scope='function')
+def autzen_data(autzen_filepath: str, autzen_storage: StorageConfig) -> Generator[Data, None, None]:
+    d = Data(autzen_filepath, autzen_storage)
+    yield d
+
+@pytest.fixture(scope='function')
+def metric_data(autzen_data: Data) -> Generator[pd.DataFrame, None, None]:
+    p = autzen_data.pipeline
+    autzen_data.execute()
+    points = p.get_dataframe(0)
+    points.loc[:, 'xi'] = np.floor(points.xi)
+    points.loc[:, 'yi'] = np.ceil(points.yi)
+    points = points.loc[points.xi == 1]
+    points = points.loc[points.yi == 437]
+    yield points[['Z', 'xi', 'yi']]
+
+@pytest.fixture(scope='function')
+def metric_shatter_config(tmp_path_factory, copc_filepath, attrs, metrics, bounds,
+        date, crs, resolution) -> Generator[pd.Series, None, None]:
+
+    path = tmp_path_factory.mktemp("test_tdb")
+    p = os.path.abspath(path)
+    log = Log('DEBUG')
+
+    def dummy_fn(df: pd.DataFrame) -> pd.DataFrame:
+        assert isinstance(df, pd.DataFrame)
+        ndf = df[df['NumberOfReturns'] >= 1]
+        assert isinstance(ndf, pd.DataFrame)
+        return ndf
+
+    metrics[0].add_filter(dummy_fn, 'This is a function.')
+    metrics[0].attributes=attrs
+
+    """Make output"""
+    st_config=StorageConfig(tdb_dir=p,
+                        log=log,
+                        crs=crs,
+                        root=bounds,
+                        resolution=resolution,
+                        attrs=attrs,
+                        metrics=metrics,
+                        version=svversion)
+
+    s = Storage.create(st_config)
+    sh_config = ShatterConfig(tdb_dir=p,
+            log=log,
+            filename=copc_filepath,
+            bounds=bounds,
+            debug=True,
+            date=date)
+    yield sh_config
\ No newline at end of file
diff --git a/tests/fixtures/shatter_fixtures.py b/tests/fixtures/shatter_fixtures.py
new file mode 100644
index 0000000..8f9cff3
--- /dev/null
+++ b/tests/fixtures/shatter_fixtures.py
@@ -0,0 +1,96 @@
+import pytest
+from typing import Generator
+from uuid import uuid4
+import os
+
+from silvimetric import __version__ as svversion
+from silvimetric import StorageConfig, ShatterConfig, Storage, Log, Bounds
+
+@pytest.fixture(scope="function")
+def s3_bucket() -> Generator[str, None, None]:
+    yield "silvimetric"
+
+@pytest.fixture(scope='function')
+def s3_uri(s3_bucket) -> Generator[str, None, None]:
+    uuid = uuid4()
+    yield f"s3://{s3_bucket}/test_silvimetric/{uuid}"
+
+@pytest.fixture(scope="function")
+def s3_storage_config(s3_uri, bounds, resolution, crs, attrs, metrics) -> Generator[StorageConfig, None, None]:
+    yield StorageConfig(bounds, crs, resolution, attrs, metrics,
+                        svversion, tdb_dir=s3_uri)
+
+@pytest.fixture(scope='function')
+def s3_storage(s3_storage_config) -> Generator[Storage, None, None]:
+    import subprocess
+    yield Storage.create(s3_storage_config)
+    subprocess.call(["aws", "s3", "rm", "--recursive", s3_storage_config.tdb_dir])
+
+@pytest.fixture(scope="function")
+def s3_shatter_config(s3_storage, copc_filepath, attrs, metrics, date) -> Generator[ShatterConfig, None, None]:
+    config = s3_storage.config
+    yield ShatterConfig(filename=copc_filepath, attrs=attrs, metrics=metrics,
+                        debug=True, tdb_dir=config.tdb_dir, date=date)
+
+@pytest.fixture(scope='function')
+def uneven_storage_config(tmp_path_factory, bounds, crs, attrs, metrics) -> Generator[StorageConfig, None, None]:
+    log = Log('INFO')
+    path = tmp_path_factory.mktemp("test_tdb")
+    p = os.path.abspath(path)
+
+    log = Log('DEBUG')
+    sc = StorageConfig(tdb_dir = p,
+                        log = log,
+                        crs = crs,
+                        root = bounds,
+                        resolution = 7,
+                        attrs = attrs,
+                        metrics = metrics,
+                        version = svversion)
+    Storage.create(sc)
+    yield sc
+
+@pytest.fixture(scope='function')
+def uneven_shatter_config(copc_filepath, uneven_storage_config, date) -> Generator[ShatterConfig, None, None]:
+    tdb_dir = uneven_storage_config.tdb_dir
+    log = uneven_storage_config.log
+
+    s = ShatterConfig(tdb_dir = tdb_dir,
+                      log = log,
+                      filename = copc_filepath,
+                      attrs = uneven_storage_config.attrs,
+                      metrics = uneven_storage_config.metrics,
+                      debug = True,
+                      date=date)
+    yield s
+
+@pytest.fixture(scope='function')
+def partial_storage_config(tmp_path_factory, crs, attrs, metrics) -> Generator[StorageConfig, None, None]:
+    path = tmp_path_factory.mktemp("test_tdb")
+    p = os.path.abspath(path)
+    log = Log('DEBUG')
+
+    bounds = Bounds(300,300,450,450)
+    sc = StorageConfig(tdb_dir = p,
+                        log = log,
+                        crs = crs,
+                        root = bounds,
+                        resolution = 30,
+                        attrs = attrs,
+                        metrics = metrics,
+                        version = svversion)
+    Storage.create(sc)
+    yield sc
+
+@pytest.fixture(scope='function')
+def partial_shatter_config(copc_filepath, date, partial_storage_config) -> Generator[ShatterConfig, None, None]:
+    tdb_dir = partial_storage_config.tdb_dir
+    psc: StorageConfig = partial_storage_config
+    log = Log('INFO') # INFO
+    yield ShatterConfig(tdb_dir=tdb_dir,
+                      log=log,
+                      attrs=psc.attrs,
+                      metrics=psc.metrics,
+                      filename=copc_filepath,
+                      debug=True,
+                      date=date)
\ No newline at end of file
diff --git a/tests/fixtures/western_fixtures.py b/tests/fixtures/western_fixtures.py
new file mode 100644
index 0000000..8fc991f
--- /dev/null
+++ b/tests/fixtures/western_fixtures.py
@@ -0,0 +1,64 @@
+import tempfile
+from typing import Generator
+
+from silvimetric import Storage, ShatterConfig, Log, Bounds, StorageConfig
+from silvimetric import __version__ as svversion
+import os, shutil
+import pytest
+
+@pytest.fixture(scope="class")
+def WebMercator():
+    yield "EPSG:3857"
+
+@pytest.fixture(scope='class')
+def WesternBounds() -> Generator[Bounds, None, None]:
+    b = Bounds(-14100053.268191, 3058230.975702, -11138180.816218, 6368599.176434)
+    yield b
+
+@pytest.fixture(scope='function')
+def western_filepath(tmp_path_factory) -> Generator[str, None, None]:
+
+    temp_name = next(tempfile._get_candidate_names())
+    path = tmp_path_factory.mktemp(temp_name)
+    yield os.path.abspath(path)
+    shutil.rmtree(path)
+
+
+@pytest.fixture(scope='function')
+def western_config(western_filepath, WesternBounds, resolution, WebMercator,
+        attrs, metrics) -> Generator[StorageConfig, None, None]:
+    log = Log(20)
+    yield StorageConfig(tdb_dir = western_filepath,
+                        log = log,
+                        crs = WebMercator,
+                        root = WesternBounds,
+                        resolution = resolution,
+                        attrs = attrs,
+                        metrics = metrics,
+                        version = svversion)
+
+@pytest.fixture(scope='function')
+def western_storage(western_config) -> Generator[Storage, None, None]:
+    yield Storage.create(western_config)
+
+@pytest.fixture(scope='function')
+def western_pipeline() -> Generator[str, None, None]:
+    path = os.path.join(os.path.dirname(__file__), "data",
+            "western_us.json")
+    assert os.path.exists(path)
+    yield os.path.abspath(path)
+
+@pytest.fixture(scope='function')
+def western_shatter_config(western_pipeline, western_storage, bounds, date) -> Generator[ShatterConfig, None, None]:
+    log = Log(20) # INFO
+    st = western_storage.config
+
+    s = ShatterConfig(tdb_dir = st.tdb_dir,
+        log = log,
+        filename = western_pipeline,
+        attrs = st.attrs,
+        metrics = st.metrics,
+        bounds=bounds,
+        debug = True,
+        date=date)
+    yield s
\ No newline at end of file
diff --git a/tests/test_cli.py b/tests/test_cli.py
new file mode 100644
index 0000000..a6e524b
--- /dev/null
+++ b/tests/test_cli.py
@@ -0,0 +1,148 @@
+from click.testing import CliRunner
+import os
+from datetime import datetime
+
+import numpy as np
+from pytest import TempPathFactory
+
+from silvimetric.cli import cli
+from silvimetric.commands import shatter, info
+from silvimetric import ShatterConfig, Bounds, StorageConfig, Storage, ExtractConfig
+
+class TestCli(object):
+    def test_cli_init(self, tmp_path_factory: TempPathFactory,
+            runner: CliRunner, bounds: Bounds) -> None:
+        path = tmp_path_factory.mktemp("test_tdb")
+        p = os.path.abspath(path)
+        res = runner.invoke(cli.cli, args=["-d", p, "--debug",
+                "--scheduler", "single-threaded", "initialize", "--resolution",
+                '10', '--crs', 'EPSG:3857', '--bounds', str(bounds)],
+                catch_exceptions=False)
+        assert res.exit_code == 0
+
+    def test_cli_metric_groups(self, tmp_path_factory: TempPathFactory,
+            runner: CliRunner, bounds: Bounds) -> None:
+        path = tmp_path_factory.mktemp("test_tdb")
+        p = os.path.abspath(path)
+        res = runner.invoke(cli.cli, args=["-d", p, "--debug", "--scheduler",
+                "single-threaded", "initialize", "--resolution", "10", "--crs", "EPSG:3857",
+                "--bounds", str(bounds), "-m", "stats,p_moments"], catch_exceptions=False)
+        assert res.exit_code == 0
+
+    def test_cli_metric_file(self, tmp_path_factory: TempPathFactory,
+            runner: CliRunner, bounds: Bounds, date: datetime, copc_filepath: str) -> None:
+        # test that we can pass in a file containing metrics and then run shatter over it
+        path = tmp_path_factory.mktemp("test_tdb")
+        p = os.path.abspath(path)
+        fakes = os.path.abspath('tests/fixtures/fake_metrics.py')
+        res = runner.invoke(cli.cli, args=["-d", p, "--debug",
+                "--scheduler", "single-threaded", "initialize",
+                "--resolution", "10", "--crs", "EPSG:3857", "--bounds",
+                str(bounds), "-m", f"{fakes}"],
+                catch_exceptions=False)
+        assert res.exit_code == 0
+        res2 = runner.invoke(cli.cli, args=["-d", p,
+                "--scheduler", "single-threaded",
+                "shatter", copc_filepath,
+                "--date", date.isoformat()+'Z',
+                "--tilesize", '10'], catch_exceptions=False)
+        assert res2.exit_code == 0
+        storage = Storage.from_db(p)
+        with storage.open('r') as a:
+            exists_vals = a[:,:]['m_Z_exists']
+            assert all(exists_vals)
+            count_vals = a[:,:]['m_Z_count']
+            assert all(count_vals == 100)
+
+
+
+    def test_cli_shatter(self, runner: CliRunner, maxy: float, date: datetime,
+            tdb_filepath: str, copc_filepath: str, storage: shatter.Storage) -> None:
+
+        res = runner.invoke(cli.cli, args=["-d", tdb_filepath,
+                "--scheduler", "single-threaded",
+                "shatter", copc_filepath,
+                "--date", date.isoformat()+'Z',
+                "--tilesize", '10'], catch_exceptions=False)
+        assert res.exit_code == 0
+
+        with storage.open('r') as a:
+            assert a[:,:]['Z'].shape[0] == 100
+            xdom = a.schema.domain.dim('X').domain[1]
+            ydom = a.schema.domain.dim('Y').domain[1]
+            assert xdom == 10
+            assert ydom == 10
+
+            for xi in range(xdom):
+                for yi in range(ydom):
+                    a[xi, yi]['Z'].size == 1
+                    a[xi, yi]['Z'][0].size == 900
+                    # this should have all indices from 0 to 9 filled.
+                    # if oob error, it's not this test's fault
+                    assert bool(np.all( a[xi, yi]['Z'][0] == (
+                        (maxy/storage.config.resolution) - (yi + 1)) ))
+
+    def test_cli_scan(self, runner: CliRunner, copc_filepath: str,
+            storage_config: StorageConfig) -> None:
+        tdb_dir = storage_config.tdb_dir
+        res = runner.invoke(cli.cli, args=['-d', tdb_dir, '--scheduler',
+            'single-threaded', 'scan', copc_filepath])
+        print(res.output)
+        assert res.exit_code == 0
+
+    def test_cli_info(self, tdb_filepath: str, runner: CliRunner,
+            shatter_config: ShatterConfig) -> None:
+        shatter.shatter(shatter_config)
+        res = runner.invoke(cli.cli, args=['-d', tdb_filepath,
+            '--scheduler', 'single-threaded', 'info'])
+        assert res.exit_code == 0
+
+    def test_cli_extract(self, runner: CliRunner, extract_config: ExtractConfig,
+            storage: Storage):
+        atts = []
+        for a in extract_config.attrs:
+            atts.append('-a')
+            atts.append(a.name)
+        ms = []
+        for m in extract_config.metrics:
+            ms.append('-m')
+            ms.append(m.name)
+        out_dir = extract_config.out_dir
+        tdb_dir = extract_config.tdb_dir
+
+        res = runner.invoke(cli.cli, args=['-d', tdb_dir, '--scheduler',
+            'single-threaded', 'extract', *atts, *ms, '--outdir' ,out_dir])
+
+        assert res.exit_code == 0
+
+    def test_cli_delete(self, runner: CliRunner, shatter_config: ShatterConfig,
+            pre_shatter: int) -> None:
+        i = shatter_config.name
+        res = runner.invoke(cli.cli, ['-d', shatter_config.tdb_dir,
+                '--scheduler', 'single-threaded', 'delete', '--id', i])
+        assert res.exit_code == 0
+
+        i = info.info(shatter_config.tdb_dir)
+        assert len(i['history']) == 0
+
+    def test_cli_restart(self, runner: CliRunner, shatter_config: ShatterConfig,
+            pre_shatter: int) -> None:
+
+        i = shatter_config.name
+        res = runner.invoke(cli.cli, ['-d', shatter_config.tdb_dir,
+                '--scheduler', 'single-threaded', 'restart', '--id', i])
+        assert res.exit_code == 0
+        i = info.info(shatter_config.tdb_dir)
+        assert len(i['history']) == 1
+        assert ShatterConfig.from_dict(i['history'][0])
+
+    def test_cli_resume(self, runner: CliRunner, shatter_config: ShatterConfig,
+            pre_shatter: int) -> None:
+
+        i = shatter_config.name
+        res = runner.invoke(cli.cli, ['-d', shatter_config.tdb_dir,
+                '--scheduler', 'single-threaded', 'restart', '--id', i])
+        assert res.exit_code == 0
+        i = info.info(shatter_config.tdb_dir)
+        assert len(i['history']) == 1
+        assert ShatterConfig.from_dict(i['history'][0])
\ No newline at end of file
diff --git a/tests/test_commands.py b/tests/test_commands.py
index 4e5fd92..bcabe8a 100644
--- a/tests/test_commands.py
+++ b/tests/test_commands.py
@@ -1,50 +1,17 @@
-import uuid
-from time import sleep
-from typing import List
-from datetime import datetime
-
-from silvimetric.commands import scan, info, shatter
-from silvimetric.resources import ShatterConfig
-
-import pytest
-import conftest
-
-@pytest.fixture(scope='function')
-def config_split(shatter_config: ShatterConfig) -> List[ShatterConfig]:
-    sc = shatter_config
-    config_split = []
-    day = 1
-    for b in sc.bounds.bisect():
-        day += 1
-        date = datetime(2011, 1, day)
-
-        config_split.append(ShatterConfig(
-            sc.filename,
-            date,
-            sc.attrs,
-            sc.metrics,
-            b,
-            uuid.uuid4(),
-            sc.tile_size,
-            tdb_dir=sc.tdb_dir,
-            log=sc.log
-        ))
-
-    for s in config_split:
-        shatter.shatter(s)
-    sleep(1)
-
-    return config_split
+import json
 
+from silvimetric.commands import scan, info, shatter, manage
+from silvimetric import ShatterConfig, Storage
 
 class TestCommands(object):
 
     def test_scan(self, shatter_config):
         s = shatter_config
         res = scan.scan(s.tdb_dir, s.filename, s.bounds, 10, 10, 5)
-        assert res == 1
+        assert res['tile_info']['recommended'] == 1
+
         res = scan.scan(s.tdb_dir, s.filename, s.bounds, depth=5)
-        assert res == 100
+        assert res['tile_info']['recommended'] == 100
 
     def test_info(self, tdb_filepath, config_split):
         i = info.info(tdb_filepath)
@@ -67,4 +34,60 @@ def test_info(self, tdb_filepath, config_split):
         d = config_split[1].date
         i = info.info(tdb_filepath, start_time=d, end_time=d)
         assert len(i['history']) == 1
-        assert i['history'][0] == config_split[1].to_json()
\ No newline at end of file
+        assert i['history'][0] == config_split[1].to_json()
+
+    def test_delete(self, tdb_filepath, config_split):
+        ids = [c.name for c in config_split]
+
+        for i in range(1,len(ids)+1):
+            h = info.info(tdb_dir=tdb_filepath)
+
+            assert bool(h['history'])
+            assert len(h['history']) == 5-i
+            idx = i - 1
+            manage.delete(tdb_filepath, ids[idx])
+
+        h = info.info(tdb_dir=tdb_filepath)
+        assert not bool(h['history'])
+
+    def test_311_failure(self, shatter_config, dask_proc_client):
+        # make sure we handle tiledb contexts correctly within dask
+
+        tdb_dir = shatter_config.tdb_dir
+        shatter.shatter(shatter_config)
+
+        i = info.info(tdb_dir)
+        history = i['history'][-1]
+
+        finished_config = ShatterConfig.from_dict(history)
+        sh_id = finished_config.name
+        manage.delete(tdb_dir, sh_id)
+
+    def test_restart(self, tdb_filepath, config_split):
+        ids = [c.name for c in config_split]
+
+        for i in range(1,len(ids)+1):
+            h = info.info(tdb_dir=tdb_filepath)
+
+            assert bool(h['history'])
+            assert len(h['history']) == len(ids)
+            manage.restart(tdb_filepath, ids[i-1])
+
+        h = info.info(tdb_dir=tdb_filepath)
+        assert bool(h['history'])
+        assert len(h['history']) == 4
+
+    def test_resume(self, shatter_config: ShatterConfig, test_point_count):
+        # test that shatter only operates on cells not touched by nonempty_domain
+        storage = Storage.from_db(shatter_config.tdb_dir)
+
+        # modify shatter config
+        shatter_config.tile_size=1
+        shatter_config.mbr = (((0,4), (0,4)),)
+        # send to shatter
+        pc = shatter.shatter(shatter_config)
+        assert pc == test_point_count - (test_point_count / 4)
+        # check output config to see what the nonempthy_domain is
+        #   - should be only only values outside of that tuple
+        m = json.loads(storage.getMetadata('shatter', 1))
+        assert len(m['mbr']) == 75
diff --git a/tests/test_configuration.py b/tests/test_configuration.py
index 9f81dd0..6c27838 100644
--- a/tests/test_configuration.py
+++ b/tests/test_configuration.py
@@ -1,40 +1,28 @@
-import pytest
-import json
-import os
-import dataclasses
-
-
-from silvimetric.resources import StorageConfig, Bounds, Log, Metric
-
-@pytest.fixture(scope='function')
-def tdb_filepath(tmp_path_factory) -> str:
-    path = tmp_path_factory.mktemp("test_tdb")
-    yield os.path.abspath(path)
-
-@pytest.fixture(scope="function")
-def config(tdb_filepath, resolution, attrs, minx, maxx, miny, maxy,
-           crs) -> StorageConfig:
-
-    b = Bounds(minx, miny, maxx, maxy)
-    log = Log(20)
-    config = StorageConfig(tdb_dir = tdb_filepath,
-                           log = log,
-                           root = b,
-                           resolution = resolution,
-                           crs = crs,
-                           attrs = attrs)
-    yield config
-
+from silvimetric import StorageConfig, ShatterConfig, ExtractConfig
 
 class Test_Configuration(object):
 
-    def test_serialization(self, config: StorageConfig):
-
-        j = str(config)
+    def test_serialization(self, storage_config, shatter_config, extract_config):
+        # storage
+        j = str(storage_config)
         c = StorageConfig.from_string(j)
+
         mean = [ m for m in c.metrics if m.name == 'mean']
         assert len(mean) == 1
 
-        assert int(mean[0]([2,2,2,2])) == 2
-        assert config == c
-
+        assert int(mean[0]._method([2,2,2,2])) == 2
+        cd = c.to_json()
+        scd = storage_config.to_json()
+        cd.pop('log')
+        scd.pop('log')
+        assert scd == cd
+
+        # shatter
+        sh_str = str(shatter_config)
+        sh_cfg = ShatterConfig.from_string(sh_str)
+        assert shatter_config == sh_cfg
+
+        # extract
+        ex_str = str(extract_config)
+        ex_cfg = ExtractConfig.from_string(ex_str)
+        assert extract_config == ex_cfg
\ No newline at end of file
diff --git a/tests/test_data.py b/tests/test_data.py
index e32054e..6a5d106 100644
--- a/tests/test_data.py
+++ b/tests/test_data.py
@@ -1,70 +1,46 @@
-import pytest
-import os
-
-from silvimetric.resources import Data
-from silvimetric.resources import Storage
-from silvimetric.resources import Bounds
-
-import conftest
-
-@pytest.fixture(scope='session')
-def no_cell_line_pc():
-    yield 84100
-
-@pytest.fixture(scope='session')
-def no_cell_line_path():
-    path = os.path.join(os.path.dirname(__file__), "data",
-            "test_data_2.copc.laz")
-    assert os.path.exists(path)
-    yield os.path.abspath(path)
-
-@pytest.fixture(scope='session')
-def no_cell_line_pipeline():
-    path = os.path.join(os.path.dirname(__file__), "data",
-            "test_data_2.json")
-    assert os.path.exists(path)
-    yield os.path.abspath(path)
-
-@pytest.fixture(scope='session')
-def autzen_filepath() -> str:
-    path = os.path.join(os.path.dirname(__file__), "data",
-            "autzen-small.copc.laz")
-    assert os.path.exists(path)
-    yield os.path.abspath(path)
+from silvimetric import Data, Storage, Bounds
 
 class Test_Data(object):
 
-    def test_filepath(self, no_cell_line_path, storage: Storage, no_cell_line_pc, bounds):
+    def test_filepath(self, no_cell_line_path, storage_config, no_cell_line_pc, bounds):
         """Check open a COPC file"""
-        data = Data(no_cell_line_path, storage.config)
+        data = Data(no_cell_line_path, storage_config)
         assert data.is_pipeline() == False
         data.execute()
         assert len(data.array) == no_cell_line_pc
         assert data.estimate_count(bounds) == no_cell_line_pc
 
-    def test_pipeline(self, no_cell_line_pipeline, bounds, storage: Storage, no_cell_line_pc):
+    def test_pipeline(self, no_cell_line_pipeline, bounds, storage_config, no_cell_line_pc):
         """Check open a pipeline"""
-        data = Data(no_cell_line_pipeline, storage.config)
+        data = Data(no_cell_line_pipeline, storage_config)
         assert data.is_pipeline() == True
         data.execute()
         assert len(data.array) == no_cell_line_pc
         assert data.estimate_count(bounds) == no_cell_line_pc
 
-    def test_pipeline_bounds(self, no_cell_line_pipeline, bounds, storage: Storage, no_cell_line_pc):
+    def test_pipeline_bounds(self, no_cell_line_pipeline, bounds, storage_config, no_cell_line_pc):
         """Check open a pipeline with our own bounds"""
         ll = list(bounds.bisect())[0]
-        data = Data(no_cell_line_pipeline, storage.config, bounds = ll)
+
+        #data will be collared upon execution, extra data will be grabbed
+        minx, miny, maxx, maxy = ll.get()
+        collared = Bounds(minx - 30, miny - 30, maxx + 30, maxy + 30)
+
+        data = Data(no_cell_line_pipeline, storage_config, bounds = ll)
         assert data.is_pipeline() == True
         data.execute()
-        assert len(data.array) == no_cell_line_pc / 4
+
+        collared_count = data.count(collared)
+        assert len(data.array) == collared_count
+
         assert data.estimate_count(ll) == no_cell_line_pc
         assert data.count(ll) == no_cell_line_pc / 4
 
 class Test_Autzen(object):
 
-    def test_filepath(self, autzen_filepath, storage: Storage):
+    def test_filepath(self, autzen_filepath, storage_config):
         """Check open Autzen """
-        data = Data(autzen_filepath, storage.config)
+        data = Data(autzen_filepath, storage_config)
         assert data.is_pipeline() == False
         data.execute()
         assert len(data.array) == 577637
diff --git a/tests/test_chunk.py b/tests/test_extents.py
similarity index 89%
rename from tests/test_chunk.py
rename to tests/test_extents.py
index e46da03..5efd829 100644
--- a/tests/test_chunk.py
+++ b/tests/test_extents.py
@@ -1,12 +1,14 @@
 import numpy as np
-import dask
-import pdal
-import pytest
 import datetime
 
-from silvimetric.resources import Extents
+from silvimetric import Extents
 from silvimetric.commands.shatter import run
 
+def check_for_overlap(leaves: list[Extents], chunk: Extents):
+    idx = [i for l in leaves for i in l.get_indices() ]
+    u, c = np.unique(idx, return_counts=True)
+    assert not np.any(np.where(c > 1))
+
 def check_for_holes(leaves: list[Extents], chunk: Extents):
     ind = np.array([], dtype=chunk.get_indices().dtype)
     for leaf in leaves:
@@ -83,19 +85,13 @@ def check_indexing(extents: Extents, leaf_list):
 
 
 class TestExtents(object):
-    @pytest.fixture(scope='function', autouse=True)
-    def filtered(self, copc_data, extents: Extents):
-        return list(extents.chunk(copc_data, 1))
-
-    @pytest.fixture(scope='function')
-    def unfiltered(self, extents: Extents):
-        return list(extents.get_leaf_children(30))
-
     def test_indexing(self, extents, filtered, unfiltered):
         check_indexing(extents, filtered)
         check_for_holes(filtered, extents)
+        check_for_overlap(filtered, extents)
         check_indexing(extents, unfiltered)
         check_for_holes(unfiltered, extents)
+        check_for_overlap(unfiltered, extents)
 
     # def test_cells(self, copc_filepath, unfiltered, resolution):
     #     flag = False
@@ -122,8 +118,8 @@ def test_pointcount(self, filtered, unfiltered, test_point_count, shatter_config
 
         with storage.open('w') as tdb:
             shatter_config.start_time = datetime.datetime.now().timestamp() * 1000
-            fc = run(filtered, shatter_config, storage, tdb)
-            ufc = run(unfiltered, shatter_config, storage, tdb)
+            fc = run(filtered, shatter_config, storage)
+            ufc = run(unfiltered, shatter_config, storage)
 
             assert fc == ufc, f"""
                 Filtered and unfiltered point counts don't match.
diff --git a/tests/test_extract.py b/tests/test_extract.py
index ccad6e9..f16726d 100644
--- a/tests/test_extract.py
+++ b/tests/test_extract.py
@@ -1,80 +1,43 @@
-import pytest
-import os
 from pathlib import Path
 from osgeo import gdal
 from pyproj import CRS
-import copy
-import uuid
 
-from silvimetric.commands.shatter import shatter
-from silvimetric.commands.extract import extract
-from silvimetric.resources import Metrics, Attribute, ExtractConfig, Extents, Storage, Log
-
-@pytest.fixture(scope='function')
-def tif_filepath(tmp_path_factory) -> str:
-    path = tmp_path_factory.mktemp("test_tifs")
-    yield os.path.abspath(path)
-
-@pytest.fixture(scope='function')
-def extract_attrs(dims)->list[str]:
-    yield [Attribute('Z', dtype=dims['Z']), Attribute('Intensity', dtype=dims['Intensity'])]
-
-@pytest.fixture(scope='function')
-def extract_config(tdb_filepath, tif_filepath, metrics, shatter_config, extract_attrs):
-    shatter(shatter_config)
-    log = Log(20)
-    c =  ExtractConfig(tdb_dir = tdb_filepath,
-                       log = log,
-                       out_dir = tif_filepath,
-                       attrs = extract_attrs,
-                       metrics = metrics)
-    yield c
-
-@pytest.fixture(scope='function')
-def multivalue_config(tdb_filepath, tif_filepath, metrics, shatter_config, extract_attrs):
-
-    shatter(shatter_config)
-    e = Extents.from_storage(shatter_config.tdb_dir)
-    b: Extents = e.split()[0]
-
-    second_config = copy.deepcopy(shatter_config)
-    second_config.bounds = b.bounds
-    second_config.name = uuid.uuid4()
-    second_config.point_count = 0
-
-    shatter(second_config)
-    log = Log(20)
-    c =  ExtractConfig(tdb_dir = tdb_filepath,
-                       log = log,
-                       out_dir = tif_filepath,
-                       attrs = extract_attrs,
-                       metrics = metrics)
-    yield c
+from silvimetric import grid_metrics, ExtractConfig, Extents, Log, extract, Storage
 
 def tif_test(extract_config):
     minx, miny, maxx, maxy = extract_config.bounds.get()
     resolution = extract_config.resolution
-    filenames = [Metrics[m.name].entry_name(a.name)
+    filenames = [grid_metrics[m.name].entry_name(a.name)
                     for m in extract_config.metrics
                     for a in extract_config.attrs]
-    e = Extents.from_storage(extract_config.tdb_dir)
-    yimin = e.get_indices()['y'].min()
+    storage = Storage.from_db(extract_config.tdb_dir)
+    e = Extents(extract_config.bounds, extract_config.resolution, storage.config.root)
+    root_maxy = storage.config.root.maxy
 
     for f in filenames:
+        if 'Return' in f:
+            continue
         path = Path(extract_config.out_dir) / f'{f}.tif'
         assert path.exists()
 
         raster: gdal.Dataset = gdal.Open(str(path))
         derived = CRS.from_user_input(raster.GetProjection())
+        rminx, xres, xskew, rmaxy, yskew, yres  = raster.GetGeoTransform()
+        assert rminx == minx
+        assert rmaxy == maxy
+        assert xres == resolution
+        assert -yres == resolution
+
+
         assert derived == extract_config.crs
 
-        xsize = int((maxx - minx) / resolution)
-        ysize = int((maxy - miny) / resolution)
+        xsize = e.x2
+        ysize = e.y2
         assert raster.RasterXSize == xsize
         assert raster.RasterYSize == ysize
 
         r = raster.ReadAsArray()
-        assert all([ r[y,x] == ((maxy/resolution)-(y+yimin))  for y in range(ysize) for x in range(xsize)])
+        assert all([ r[y,x] == ((root_maxy/resolution)-y-1)  for y in range(e.y1, e.y2) for x in range(e.x1, e.x2)])
 
 class Test_Extract(object):
 
diff --git a/tests/test_fusion.py b/tests/test_fusion.py
new file mode 100644
index 0000000..838724f
--- /dev/null
+++ b/tests/test_fusion.py
@@ -0,0 +1,5 @@
+
+
+# Test against Bob's FUSION data
+class TestFusion():
+    pass
\ No newline at end of file
diff --git a/tests/test_metrics.py b/tests/test_metrics.py
index 16be781..d2d273b 100644
--- a/tests/test_metrics.py
+++ b/tests/test_metrics.py
@@ -1,90 +1,97 @@
-import os
-import pytest
 import numpy as np
 import pandas as pd
-import dask
-
-from silvimetric.resources import Metrics, Log, ShatterConfig, StorageConfig
-from silvimetric.resources import Storage, Bounds, Extents
-from silvimetric.commands import shatter
-from silvimetric import __version__ as svversion
-
-@pytest.fixture(scope='function')
-def metric_input():
-    # randomly generated int array
-    yield np.array([54, 64, 54, 18, 76, 82, 85, 12,  7, 32, 89, 79, 78, 58, 17,  5, 77,
-        20, 63, 87, 37, 27,  7, 22, 34, 61, 52, 64, 65, 90, 76, 88,  0, 98,
-        78, 26, 88, 44, 78, 54, 38, 28, 94, 74, 90, 84, 30,  7, 81, 75, 22,
-        12, 67, 10, 46, 62, 79, 52, 24, 88,  4, 73, 75, 35, 75, 16, 66, 43,
-        28, 72, 15, 80, 53,  1, 28, 70, 20, 42, 83, 80, 82,  2, 41, 38, 23,
-        17, 18, 19, 43, 30, 88, 41, 60, 29, 93, 27,  1, 13, 93, 82], np.int32)
-
-@pytest.fixture(scope='session')
-def autzen_filepath() -> str:
-    path = os.path.join(os.path.dirname(__file__), "data",
-            "autzen-small.copc.laz")
-    assert os.path.exists(path)
-    yield os.path.abspath(path)
-
-@pytest.fixture(scope='function')
-def autzen_bounds():
-    yield Bounds(635579.2, 848884.83, 639003.73, 853536.21)
-
-@pytest.fixture(scope='function')
-def metric_st_config(tdb_filepath, autzen_bounds, resolution, crs, attrs, metrics):
-    log = Log(20)
-    yield StorageConfig(tdb_dir = tdb_filepath,
-                        log = log,
-                        crs = crs,
-                        root = autzen_bounds,
-                        resolution = resolution,
-                        attrs = attrs,
-                        version = svversion)
-
-@pytest.fixture(scope='function')
-def metric_shatter(tdb_filepath, autzen_filepath, metric_st_config, autzen_bounds, app_config, date):
-    log = Log(20) # INFO
-    Storage.create(metric_st_config)
-    e = Extents(autzen_bounds, 30, autzen_bounds)
-    e1 = e.split()[1]
-    s = ShatterConfig(tdb_dir = tdb_filepath,
-                      log=log,
-                      filename=autzen_filepath,
-                      attrs=metric_st_config.attrs,
-                      metrics=metric_st_config.metrics,
-                      bounds=e1.bounds,
-                      debug=True,
-                      date=date,
-                      tile_size=700)
-
-    yield s
-
-class Test_Metrics(object):
-
-    def test_input_output(self, metric_input):
-        for m in Metrics:
-            try:
-                out = Metrics[m]._method(metric_input)
-                msg = f"Metric {Metrics[m].name} does not output a valid data type."\
-                        f" Interpretted type: {type(out)}"
-                assert not any([
-                    isinstance(out, np.ndarray),
-                    isinstance(out, list),
-                    isinstance(out, tuple),
-                    isinstance(out, pd.Series),
-                    isinstance(out, pd.DataFrame),
-                    isinstance(out, str)
-                ]), msg
-                assert any((
-                    isinstance(out, np.number),
-                    isinstance(out, int),
-                    isinstance(out, float),
-                    isinstance(out, complex),
-                    isinstance(out, bool)
-                )), msg
-            except:
-                print('yi')
-
-    def test_metric_shatter(self, metric_shatter):
-        dask.config.set(scheduler='processes')
-        shatter.shatter(metric_shatter)
\ No newline at end of file
+
+from silvimetric import shatter, Storage, grid_metrics, run_metrics, Metric
+from silvimetric import all_metrics as s
+from silvimetric import l_moments
+from silvimetric.resources.attribute import Attribute
+
+class TestMetrics():
+
+    def test_metrics(self, metric_data):
+        ms: pd.DataFrame = run_metrics(metric_data,
+            list(grid_metrics.values())).compute()
+        assert isinstance(ms, pd.DataFrame)
+        adjusted = [k.split('_')[-1] for k in ms.keys()]
+
+        assert all(a in grid_metrics.keys() for a in adjusted)
+
+
+    def test_intermediate_metric(self, metric_data):
+        ms = list(l_moments.values())
+        computed = run_metrics(metric_data, ms).compute()
+
+        assert computed.m_Z_l1.any()
+        assert computed.m_Z_l2.any()
+        assert computed.m_Z_l3.any()
+        assert computed.m_Z_l4.any()
+        assert computed.m_Z_lcv.any()
+        assert computed.m_Z_lskewness.any()
+        assert computed.m_Z_lkurtosis.any()
+
+    def test_dependencies(self, metric_data):
+        # should be able to create a dependency graph
+        cv = s['cv']
+        mean = s['mean']
+        stddev = s['stddev']
+        median = s['median']
+
+        mean.dependencies = [ median ]
+        stddev.dependencies = [ median ]
+        cv.dependencies = [ mean, stddev ]
+
+        b = run_metrics(metric_data, [cv, mean]).compute()
+        # cv/mean should be there
+        assert b.m_Z_cv.any()
+        assert b.m_Z_mean.any()
+
+        #and median/stddev should not
+        assert not any(x in b.dtypes for x in  ['m_Z_median', 'm_Z_stddev'])
+
+    def test_filter(self, metric_shatter_config, test_point_count, maxy,
+            resolution):
+
+        m = metric_shatter_config.metrics[0]
+        assert len(m.filters) == 1
+
+        pc = shatter(metric_shatter_config)
+        assert pc == test_point_count
+
+        s = Storage.from_db(metric_shatter_config.tdb_dir)
+        with s.open('r') as a:
+            q = a.query(coords=False, use_arrow=False).df
+
+            nor_mean = q[:]['m_NumberOfReturns_mean']
+            nor = q[:]['NumberOfReturns']
+            assert not nor_mean.isna().any()
+            assert nor.notna().any()
+
+            assert a[:,:]['Z'].shape[0] == 100
+            xdom = a.schema.domain.dim('X').domain[1]
+            ydom = a.schema.domain.dim('Y').domain[1]
+            assert xdom == 10
+            assert ydom == 10
+
+            data = a.query(attrs=['Z'], coords=True, use_arrow=False).df[:]
+            data = data.set_index(['X','Y'])
+
+            for xi in range(xdom):
+                for yi in range(ydom):
+                    curr = data.loc[xi,yi]
+                    curr.size == 1
+                    curr.iloc[0].size == 900
+                    # this should have all indices from 0 to 9 filled.
+                    # if oob error, it's probably not this test's fault
+                    assert bool(np.all( curr.iloc[0] == ((maxy/resolution) - (yi + 1)) ))
+
+    def test_custom(self, metric_data: pd.DataFrame, attrs: list[Attribute]) -> None:
+        def m_over500(data):
+            return data[data >= 500].count()
+        z_att = attrs[0]
+        m_cust = Metric(name='over500', dtype=np.float32, method=m_over500,
+            attributes=[z_att])
+        b = run_metrics(metric_data, [m_cust]).compute()
+
+        assert b.m_Z_over500.any()
+        assert b.m_Z_over500.values[0] == 3
+
diff --git a/tests/test_shatter.py b/tests/test_shatter.py
index 44aff3d..cc10bad 100644
--- a/tests/test_shatter.py
+++ b/tests/test_shatter.py
@@ -4,15 +4,13 @@
 import dask
 import json
 import uuid
-import platform
 
+from silvimetric import Extents, Log, info, shatter, Bounds, Storage
+from silvimetric import StorageConfig, ShatterConfig
 
-from silvimetric.commands.shatter import shatter
-from silvimetric.commands.info import info
-from silvimetric.resources import Storage, Extents, ShatterConfig, Log
 
 @dask.delayed
-def write(x,y,val, s:Storage, attrs, dims, metrics):
+def write(x, y, val, s:Storage, attrs, dims, metrics):
     m_list = [m.entry_name(a.name) for m in metrics for a in attrs]
     data = { a.name: np.array([np.array([val], dims[a.name]), None], object)[:-1]
                 for a in attrs }
@@ -26,35 +24,64 @@ def write(x,y,val, s:Storage, attrs, dims, metrics):
 
 class Test_Shatter(object):
 
-    def test_shatter(self, shatter_config, storage: Storage, maxy):
+    def test_command(self, shatter_config, storage: Storage, maxy):
         shatter(shatter_config)
         with storage.open('r') as a:
-            y = a[:,0]['Z'].shape[0]
-            x = a[0,:]['Z'].shape[0]
-            assert y == 10
-            assert x == 10
-            for xi in range(x):
-                for yi in range(y):
+            assert a[:,:]['Z'].shape[0] == 100
+            xdom = a.schema.domain.dim('X').domain[1]
+            ydom = a.schema.domain.dim('Y').domain[1]
+            assert xdom == 10
+            assert ydom == 10
+
+            for xi in range(xdom):
+                for yi in range(ydom):
                     a[xi, yi]['Z'].size == 1
-                    assert bool(np.all(a[xi, yi]['Z'][0] == ((maxy/storage.config.resolution)-yi)))
+                    a[xi, yi]['Z'][0].size == 900
+                    # this should have all indices from 0 to 9 filled.
+                    # if oob error, it's not this test's fault
+                    assert bool(np.all( a[xi, yi]['Z'][0] == ((maxy/storage.config.resolution) - (yi + 1)) ))
 
+    def test_multiple(self, shatter_config, storage: Storage, maxy):
+        shatter(shatter_config)
+        with storage.open('r') as a:
+            assert a[:,:]['Z'].shape[0] == 100
+            xdom = a.schema.domain.dim('X').domain[1]
+            ydom = a.schema.domain.dim('Y').domain[1]
+            assert xdom == 10
+            assert ydom == 10
+
+            for xi in range(xdom):
+                for yi in range(ydom):
+                    a[xi, yi]['Z'].size == 1
+                    a[xi, yi]['Z'][0].size == 900
+                    # this should have all indices from 0 to 9 filled.
+                    # if oob error, it's not this test's fault
+                    assert bool(np.all( a[xi, yi]['Z'][0] == ((maxy/storage.config.resolution) - (yi + 1)) ))
+
+        # change attributes to make it a new run
         shatter_config.name = uuid.uuid4()
+        shatter_config.mbr = ()
+        shatter_config.time_slot = 2
+
         shatter(shatter_config)
         with storage.open('r') as a:
             # querying flattens to 20, there will 10 pairs of values
+            assert a[:,:]['Z'].shape[0] == 200
             assert a[:,0]['Z'].shape[0] == 20
             assert a[0,:]['Z'].shape[0] == 20
             # now test that the second set is the same as the first set
             # and test that this value is still the same as the original
             # which was set at ((maxy/resolution)-yindex)
-            for xi in range(x):
-                for yi in range(y):
+            for xi in range(xdom):
+                for yi in range(ydom):
                     a[xi, yi]['Z'].size == 2
+                    a[xi, yi]['Z'][0].size == 900
+                    a[xi, yi]['Z'][1].size == 900
                     assert bool(np.all(a[xi, yi]['Z'][1] == a[xi,yi]['Z'][0]))
-                    assert bool(np.all(a[xi, yi]['Z'][1] == ((maxy/storage.config.resolution)-yi)))
+                    assert bool(np.all(a[xi, yi]['Z'][1] == ((maxy/storage.config.resolution) - (yi + 1))))
 
-            m = info(storage.config.tdb_dir)
-            assert len(m['history']) == 2
+        m = info(storage.config.tdb_dir)
+        assert len(m['history']) == 2
 
     def test_parallel(self, storage, attrs, dims, threaded_dask, metrics):
         # test that writing in parallel doesn't affect ordering of values
@@ -75,7 +102,7 @@ def test_parallel(self, storage, attrs, dims, threaded_dask, metrics):
     def test_config(self, shatter_config, storage, test_point_count):
         shatter(shatter_config)
         try:
-            meta = storage.getMetadata('shatter')
+            meta = storage.getMetadata('shatter', shatter_config.time_slot)
         except BaseException as e:
             pytest.fail("Failed to retrieve 'shatter' metadata key." + e.args)
         meta_j = json.loads(meta)
@@ -83,7 +110,7 @@ def test_config(self, shatter_config, storage, test_point_count):
         assert pc == test_point_count
 
     @pytest.mark.parametrize('sh_cfg', ['shatter_config', 'uneven_shatter_config'])
-    def test_sub_bounds(self, sh_cfg, test_point_count, request):
+    def test_sub_bounds(self, sh_cfg, test_point_count, request, maxy, resolution):
         s = request.getfixturevalue(sh_cfg)
         storage = Storage.from_db(s.tdb_dir)
         e = Extents.from_storage(s.tdb_dir)
@@ -91,6 +118,7 @@ def test_sub_bounds(self, sh_cfg, test_point_count, request):
         pc = 0
         for b in e.split():
             log = Log(20)
+            time_slot = storage.reserve_time_slot()
             sc = ShatterConfig(tdb_dir = s.tdb_dir,
                                log = log,
                                filename = s.filename,
@@ -99,30 +127,59 @@ def test_sub_bounds(self, sh_cfg, test_point_count, request):
                                metrics = s.metrics,
                                debug = s.debug,
                                bounds = b.bounds,
-                               date=s.date)
+                               date=s.date,
+                               time_slot=time_slot)
             pc = pc + shatter(sc)
         history = info(s.tdb_dir)['history']
+        assert len(history) == 4
         assert isinstance(history, list)
         pcs = [ h['point_count'] for h in history ]
         assert sum(pcs) == test_point_count
         assert pc == test_point_count
 
+        with storage.open('r') as a:
+            xdom = a.schema.domain.dim('X').domain[1]
+            ydom = a.schema.domain.dim('Y').domain[1]
+
+            data = a.query(attrs=['Z'], coords=True, use_arrow=False).df[:]
+            data = data.set_index(['X','Y'])
+
+            for xi in range(xdom):
+                for yi in range(ydom):
+                    curr = data.loc[xi,yi]
+                    # check that each cell only has one allocation
+                    curr.size == 1
+
+    def test_partial_overlap(self, partial_shatter_config, test_point_count):
+        pc = shatter(partial_shatter_config)
+        assert pc == test_point_count / 4
+
     @pytest.mark.skipif(
         os.environ.get('AWS_SECRET_ACCESS_KEY') is None or
         os.environ.get('AWS_ACCESS_KEY_ID') is None,
         reason='Missing necessary AWS environment variables'
     )
     def test_remote_creation(self, s3_shatter_config, s3_storage):
+        # need processes scheduler to accurately test bug fix
         dask.config.set(scheduler="processes")
         resolution = s3_storage.config.resolution
         maxy = s3_storage.config.root.maxy
         shatter(s3_shatter_config)
         with s3_storage.open('r') as a:
-            y = a[:,0]['Z'].shape[0]
-            x = a[0,:]['Z'].shape[0]
-            assert y == 10
-            assert x == 10
-            for xi in range(x):
-                for yi in range(y):
-                    a[xi, yi]['Z'].size == 1
-                    assert bool(np.all(a[xi, yi]['Z'][0] == ((maxy/resolution)-yi)))
+            assert a[:,:]['Z'].shape[0] == 100
+            xdom = a.schema.domain.dim('X').domain[1]
+            ydom = a.schema.domain.dim('Y').domain[1]
+            assert xdom == 10
+            assert ydom == 10
+            # get all data local so we're not hittin s3 all the time
+            data = a.query(attrs=['Z'], coords=True, use_arrow=False).df[:]
+            data = data.set_index(['X','Y'])
+
+            for xi in range(xdom):
+                for yi in range(ydom):
+                    curr = data.loc[xi,yi]
+                    curr.size == 1
+                    curr.iloc[0].size == 900
+                    # this should have all indices from 0 to 9 filled.
+                    # if oob error, it's not this test's fault
+                    assert bool(np.all( curr.iloc[0] == ((maxy/resolution) - (yi + 1)) ))
\ No newline at end of file
diff --git a/tests/test_storage.py b/tests/test_storage.py
index f358caa..791145e 100644
--- a/tests/test_storage.py
+++ b/tests/test_storage.py
@@ -1,7 +1,10 @@
 import tiledb
 import numpy as np
+import pytest
+import os
+import copy
 
-from silvimetric.resources import Storage, Metrics, Attribute
+from silvimetric import Storage, grid_metrics, Attribute, Attributes, StorageConfig, Log
 from silvimetric import __version__ as svversion
 
 class Test_Storage(object):
@@ -16,6 +19,11 @@ def test_schema(self, storage: Storage, attrs: list[Attribute]):
                 assert s.has_attr(a.name)
                 # assert s.attr(a.name) == a.schema()
 
+    def test_time_reserve(self, storage):
+        for x in range(5):
+            time_slot = storage.reserve_time_slot()
+            assert time_slot == x + 1
+
     def test_local(self, storage: Storage, attrs: list[Attribute]):
         with storage.open('r') as st:
             sc = st.schema
@@ -36,6 +44,27 @@ def test_config(self, storage: Storage):
         assert config.crs == storage.config.crs
         assert storage.config.version == svversion
 
+    def test_metric_dependencies(self, tmp_path_factory, metrics, crs, resolution, attrs, bounds):
+        ms = copy.deepcopy(metrics)
+
+        path = tmp_path_factory.mktemp("test_tdb")
+        p = os.path.abspath(path)
+        log = Log('DEBUG')
+
+        ms[0].dependencies = [Attributes['HeightAboveGround']]
+        sc = StorageConfig(tdb_dir=p, crs=crs, resolution=resolution,
+                attrs=attrs, metrics=ms, root=bounds)
+
+        with pytest.raises(ValueError) as e:
+            Storage.create(sc)
+        assert str(e.value) == 'Missing required dependency, HeightAboveGround.'
+
+        ms[0].dependencies = [Attributes['NumberOfReturns']]
+        s = Storage.create(sc)
+        assert isinstance(s, Storage)
+
+        ms[0].dependencies = []
+
     def test_metrics(self, storage: Storage):
         m_list = storage.getMetrics()
         a_list = storage.getAttributes()
@@ -43,9 +72,9 @@ def test_metrics(self, storage: Storage):
         with storage.open('r') as st:
             s: tiledb.ArraySchema = st.schema
             for m in m_list:
-                assert m.name in Metrics.keys()
+                assert m.name in grid_metrics.keys()
                 def e_name(att):
                     return s.attr(m.entry_name(att.name))
                 def schema(att):
-                    return Metrics[m.name].schema(att)
+                    return grid_metrics[m.name].schema(att)
                 assert all([e_name(a) == schema(a) for a in a_list])
diff --git a/tests/test_western.py b/tests/test_western.py
index 81e305c..ef103a0 100644
--- a/tests/test_western.py
+++ b/tests/test_western.py
@@ -1,69 +1,6 @@
 import tiledb
-import tempfile
 import numpy as np
-from silvimetric.resources import Storage, StorageConfig, Metrics, Attribute, Bounds
-
-from silvimetric.commands.shatter import shatter
-from silvimetric.resources import Storage, Extents, ShatterConfig, Log
-from silvimetric import __version__ as svversion
-import os, shutil
-import pytest
-
-@pytest.fixture(scope="class")
-def WebMercator():
-    yield "EPSG:3857"
-
-@pytest.fixture(scope='class')
-def WesternBounds() -> Bounds:
-    b = Bounds(-14100053.268191, 3058230.975702, -11138180.816218, 6368599.176434)
-    yield b
-
-@pytest.fixture(scope='function')
-def western_filepath(tmp_path_factory):
-
-    temp_name = next(tempfile._get_candidate_names())
-    path = tmp_path_factory.mktemp(temp_name)
-    yield os.path.abspath(path)
-    shutil.rmtree(path)
-
-
-@pytest.fixture(scope='function')
-def western_config(western_filepath, WesternBounds, resolution, WebMercator, attrs, metrics):
-    log = Log(20)
-    yield StorageConfig(tdb_dir = western_filepath,
-                        log = log,
-                        crs = WebMercator,
-                        root = WesternBounds,
-                        resolution = resolution,
-                        attrs = attrs,
-                        metrics = metrics,
-                        version = svversion)
-
-@pytest.fixture(scope='function')
-def western_storage(western_config) -> Storage:
-    yield Storage.create(western_config)
-
-@pytest.fixture(scope='function')
-def western_pipeline():
-    path = os.path.join(os.path.dirname(__file__), "data",
-            "western_us.json")
-    assert os.path.exists(path)
-    yield os.path.abspath(path)
-
-@pytest.fixture(scope='function')
-def western_shatter_config(western_pipeline, western_storage):
-    log = Log(20) # INFO
-    st = western_storage.config
-
-    s = ShatterConfig(tdb_dir = st.tdb_dir,
-        log = log,
-        filename = western_pipeline,
-        attrs = st.attrs,
-        metrics = st.metrics,
-        bounds=bounds,
-        debug = True,
-        date=date)
-    yield s
+from silvimetric import Storage
 
 class Test_Western(object):
 
@@ -72,4 +9,3 @@ def test_schema(self, western_storage: Storage):
             s:tiledb.ArraySchema = st.schema
             assert s.has_attr('count')
             assert s.attr('count').dtype == np.int32
-