Repository: martinfleis/momepy
Pull request #48 (2019-08-06 16:06:50)

Preparation for release 0.1rc1

codecov[bot]: # [Codecov](https://codecov.io/gh/martinfleis/momepy/pull/48?src=pr&el=h1) Report
> Merging [#48](https://codecov.io/gh/martinfleis/momepy/pull/48?src=pr&el=desc) into [master](https://codecov.io/gh/martinfleis/momepy/commit/0d553b61c449a739a4a6568b69b6ec9394d7ac87?src=pr&el=desc) will **decrease** coverage by `0.2%`.
> The diff coverage is `74.6%`.

[![Impacted file tree graph](https://codecov.io/gh/martinfleis/momepy/pull/48/graphs/tree.svg?width=650&token=i04999UxkK&height=150&src=pr)](https://codecov.io/gh/martinfleis/momepy/pull/48?src=pr&el=tree)

```diff
@@            Coverage Diff             @@
##           master      #48      +/-   ##
==========================================
- Coverage   91.76%   91.55%   -0.21%     
==========================================
  Files           9        9              
  Lines        2125     2203      +78     
==========================================
+ Hits         1950     2017      +67     
- Misses        175      186      +11
```


| [Impacted Files](https://codecov.io/gh/martinfleis/momepy/pull/48?src=pr&el=tree) | Coverage Δ | |
|---|---|---|
| [momepy/graph.py](https://codecov.io/gh/martinfleis/momepy/pull/48/diff?src=pr&el=tree#diff-bW9tZXB5L2dyYXBoLnB5) | `98.36% <ø> (ø)` | :arrow_up: |
| [momepy/distribution.py](https://codecov.io/gh/martinfleis/momepy/pull/48/diff?src=pr&el=tree#diff-bW9tZXB5L2Rpc3RyaWJ1dGlvbi5weQ==) | `94.11% <100%> (ø)` | :arrow_up: |
| [momepy/elements.py](https://codecov.io/gh/martinfleis/momepy/pull/48/diff?src=pr&el=tree#diff-bW9tZXB5L2VsZW1lbnRzLnB5) | `92.54% <100%> (+9.24%)` | :arrow_up: |
| [momepy/dimension.py](https://codecov.io/gh/martinfleis/momepy/pull/48/diff?src=pr&el=tree#diff-bW9tZXB5L2RpbWVuc2lvbi5weQ==) | `98.36% <100%> (ø)` | :arrow_up: |
| [momepy/shape.py](https://codecov.io/gh/martinfleis/momepy/pull/48/diff?src=pr&el=tree#diff-bW9tZXB5L3NoYXBlLnB5) | `95.18% <100%> (ø)` | :arrow_up: |
| [momepy/utils.py](https://codecov.io/gh/martinfleis/momepy/pull/48/diff?src=pr&el=tree#diff-bW9tZXB5L3V0aWxzLnB5) | `77.51% <58.46%> (-11.96%)` | :arrow_down: |
| [momepy/intensity.py](https://codecov.io/gh/martinfleis/momepy/pull/48/diff?src=pr&el=tree#diff-bW9tZXB5L2ludGVuc2l0eS5weQ==) | `88.62% <68.18%> (-3.91%)` | :arrow_down: |
| [momepy/diversity.py](https://codecov.io/gh/martinfleis/momepy/pull/48/diff?src=pr&el=tree#diff-bW9tZXB5L2RpdmVyc2l0eS5weQ==) | `87.38% <91.66%> (+0.96%)` | :arrow_up: |

------

[Continue to review full report at Codecov](https://codecov.io/gh/martinfleis/momepy/pull/48?src=pr&el=continue).
> **Legend** - [Click here to learn more](https://docs.codecov.io/docs/codecov-delta)
> `Δ = absolute <relative> (impact)`, `ø = not affected`, `? = missing data`
> Powered by [Codecov](https://codecov.io/gh/martinfleis/momepy/pull/48?src=pr&el=footer). Last update [0d553b6...d0c390e](https://codecov.io/gh/martinfleis/momepy/pull/48?src=pr&el=lastupdated). Read the [comment docs](https://docs.codecov.io/docs/pull-request-comments).


Patch:
diff --git a/README.md b/README.md
index 13477d2..ecddb7d 100644
--- a/README.md
+++ b/README.md
@@ -8,26 +8,31 @@ Momepy is a project allowing advanced quantitative analysis of urban morphology.
 
 > *momepy* stands for Morphological Measuring in Python
 
-The development of momepy is still in its initial phase. The whole structure of the package can always be changed, until the release of first versions on PIP. Momepy is a result of ongoing research of [Urban Design Studies Unit (UDSU)](http://udsu-strath.com) supported by the Axel and Margaret Ax:son Johnson Foundation as a part of “The Urban Form Resilience Project” in partnership with University of Strathclyde in Glasgow, UK.
+Momepy is a result of ongoing research of [Urban Design Studies Unit (UDSU)](http://udsu-strath.com) supported by the Axel and Margaret Ax:son Johnson Foundation as a part of “The Urban Form Resilience Project” in partnership with University of Strathclyde in Glasgow, UK.
 
 Comments, suggestions, feedback, and contributions, as well as bug reports, are very welcome.
 
 ## Documentation
 Documentation of `momepy` is available at [docs.momepy.org](https://docs.momepy.org/).
 
+## User Guide
+User guide with examples of `momepy` usage is available at [guide.momepy.org](https://guide.momepy.org/).
+
 ## Install
-Momepy is not yet released. At this moment, you can either download the GitHub
-repository to work with `momepy` or install the latest development version
-from GitHub:
+See the [installation instructions](http://docs.momepy.org/en/latest/install.html) for detailed instructions.
+Momepy depends on python geospatial stack, which might cause some dependency issues.
+
+## Contributing to momepy
+Contributions of any kind to momepy are more than welcome. That does not mean new code only, but also improvements of documentation and user guide, additional tests (ideally filling the gaps in existing suite) or bug report or idea what could be added or done better.
 
-    pip install git+git://github.com/martinfleis/momepy.git
+All contributions should go through our GitHub repository. Bug reports, ideas or even questions should be raised by opening an issue on the GitHub tracker. Suggestions for changes in code or documentation should be submitted as a pull request. However, if you are not sure what to do, feel free to open an issue. All discussion will then take place on GitHub to keep the development of momepy transparent.
 
-Please ensure that you have all dependencies installed as `momepy` setup
-currently does not specify dependency requirements.
-Momepy depends on the following packages: `geopandas`, `libpysal`, and `tqdm`. All of the packages are currently required, which will be changed in future.
+If you decide to contribute to the codebase, ensure that you are using an up-to-date `master` branch. The latest development version will always be there, including a significant part of the documentation (powered by `sphinx`). The user guide is located in the `gh-pages` branch and is powered by `Jupyter book`.
 
-## References
-- Dibble J, Prelorendjos A, Romice O, et al. (2017) On the origin of spaces: Morphometric foundations of urban form evolution. Environment and Planning B: Urban Analytics and City Science 24(1): 239980831772507–24.
+Details are available in the [documentation](https://docs.momepy.org/).
 
+## Get in touch
+If you have a question regarding momepy, feel free to open an issue on GitHub. Eventually, you can contact us on [dev@momepy.org](mailto:dev@momepy.org).
 
+---
 Copyright (c) 2018-2019 Martin Fleischmann, University of Strathclyde, Urban Design Studies Unit
diff --git a/docs/api.rst b/docs/api.rst
index ce3b464..a50588c 100644
--- a/docs/api.rst
+++ b/docs/api.rst
@@ -85,11 +85,10 @@ intensity
 
    blocks_count
    courtyards
-   covered_area_ratio
    density
    elements_count
-   floor_area_ratio
    node_density
+   object_area_ratio
    reached
 
 graph
@@ -118,6 +117,7 @@ diversity
    rng
    simpson
    theil
+   gini
 
 utilities
 ---------
@@ -129,5 +129,6 @@ utilities
    network_false_nodes
    nx_to_gdf
    preprocess
+   snap_street_network_edge
    unique_id
    Queen_higher
diff --git a/docs/conf.py b/docs/conf.py
index ee420f3..75162ea 100644
--- a/docs/conf.py
+++ b/docs/conf.py
@@ -18,8 +18,8 @@ import sphinx_bootstrap_theme
 
 sys.path.insert(0, os.path.abspath('../../'))
 
-autodoc_mock_imports = ['numpy', 'scipy', 'pandas', 'shapely', 'rtree',
-                        'geopandas', 'tqdm', 'osmnx',
+autodoc_mock_imports = ['numpy', 'scipy', 'pandas', 'shapely',
+                        'geopandas', 'tqdm', 'rtree', 'osmnx',
                         'libpysal', 'networkx', 'osgeo']
 
 # class Mock(object):
@@ -57,7 +57,7 @@ author = 'Martin Fleischmann'
 # The short X.Y version
 version = '0.1'
 # The full version, including alpha/beta/rc tags
-release = '0.1alpha'
+release = '0.1.0'
 
 
 # -- General configuration ---------------------------------------------------
@@ -176,8 +176,10 @@ html_theme_options = {
     'navbar_links': [
                     # ("Gallery", "auto_examples/index"),
                     ("Installation", "install"),
-                    ("Tutorial", "tutorial"),
+                    ("User Guide", "https://guide.momepy.org/intro"),
                     ("API", "api"),
+                    ("Contributing", "contributing"),
+                    ("References", "references"),
     ],
 }
 
diff --git a/docs/contributing.rst b/docs/contributing.rst
new file mode 100644
index 0000000..a52fed5
--- /dev/null
+++ b/docs/contributing.rst
@@ -0,0 +1,206 @@
+Contributing to momepy
+======================
+
+(Contribution guidelines are largely based on `geopandas`_.)
+
+Contributions of any kind to momepy are more than welcome. That does not mean
+new code only, but also improvements of documentation and user guide, additional
+tests (ideally filling the gaps in existing suite) or bug report or idea what
+could be added or done better.
+
+All contributions should go through our GitHub repository. Bug reports, ideas or
+even questions should be raised by opening an issue on the GitHub tracker.
+Suggestions for changes in code or documentation should be submitted as a pull
+request. However, if you are not sure what to do, feel free to open an issue.
+All discussion will then take place on GitHub to keep the development of
+momepy transparent.
+
+If you decide to contribute to the codebase, ensure that you are using an
+up-to-date `master` branch. The latest development version will always be there,
+including a significant part of the documentation (powered by `sphinx`_). The
+user guide is located in the ``gh-pages`` branch and is powered by
+`Jupyter book`_.
+
+
+Seven Steps for Contributing
+----------------------------
+
+There are seven basic steps to contributing to momepy:
+
+1. Fork the momepy git repository
+2. Create a development environment
+3. Install momepy dependencies
+4. Make a development build of momepy
+5. Make changes to code and add tests
+6. Update the documentation
+7. Submit a Pull Request
+
+Each of the steps is detailed below.
+
+1. Fork the momepy git repository
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+Git can be complicated for new users, but you no longer need to use command line
+to work with git. If you are not familiar with git, we recommend using tools on
+GitHub.org, GitHub Desktop or tools with included git like Atom. However, if you
+want to use command line, you can fork momepy repository using following::
+
+    git clone git@github.com:your-user-name/momepy.git momepy-yourname
+    cd momepy-yourname
+    git remote add upstream git://github.com/martinfleis/momepy.git
+
+This creates the directory momepy-yourname and connects your repository to
+the upstream (main project) momepy repository.
+
+Then simply create a new branch of either master (for codebase and documentation)
+or gh-pages branch (for user guide).
+
+
+2. Create a development environment
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+A development environment is a virtual space where you can keep an independent
+installation of momepy. This makes it easy to keep both a stable version of
+python in one place you use for work, and a development version (which you may
+break while playing with code) in another.
+
+An easy way to create a momepy development environment is as follows:
+
+- Install either `Anaconda <http://docs.continuum.io/anaconda/>`_ or
+  `miniconda <http://conda.pydata.org/miniconda.html>`_
+- Make sure that you have :ref:`cloned the repository <contributing.forking>`
+- ``cd`` to the *momepy* source directory
+
+Tell conda to create a new environment, named ``momepy_dev``, or any other name you would like
+for this environment, by running::
+
+      conda create -n momepy_dev python=3
+
+This will create the new environment, and not touch any of your existing environments,
+nor any existing python installation.
+
+To work in this environment, Windows users should ``activate`` it as follows::
+
+      activate momepy_dev
+
+macOS and Linux users should use::
+
+      conda activate momepy_dev
+
+You will then see a confirmation message to indicate you are in the new development environment.
+
+To view your environments::
+
+      conda info -e
+
+To return to you home root environment::
+
+      deactivate
+
+See the full conda docs `here <http://conda.pydata.org/docs>`__.
+
+At this point you can easily do a *development* install, as detailed in the next sections.
+
+3. Installing Dependencies
+^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+To run *momepy* in an development environment, you must first install
+*momepy*'s dependencies. We suggest doing so using the following commands
+(executed after your development environment has been activated)::
+
+    conda install -c conda-forge geopandas networkx libpysal tqdm pysal mapclassify
+
+This should install all necessary dependencies including optional.
+
+4. Making a development build
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+Once dependencies are in place, make an in-place build by navigating to the git
+clone of the *momepy* repository and running::
+
+    python setup.py develop
+
+This will install momepy into your environment but allows any further changes
+without the need of reinstalling new version.
+
+5. Making changes and writing tests
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+*momepy* is serious about testing and strongly encourages contributors to embrace
+`test-driven development (TDD) <http://en.wikipedia.org/wiki/Test-driven_development>`_.
+This development process "relies on the repetition of a very short development cycle:
+first the developer writes an (initially failing) automated test case that defines a desired
+improvement or new function, then produces the minimum amount of code to pass that test."
+So, before actually writing any code, you should write your tests. Often the test can be
+taken from the original GitHub issue. However, it is always worth considering additional
+use cases and writing corresponding tests.
+
+*momepy* uses the `pytest testing system <http://doc.pytest.org/en/latest/>`_.
+
+Writing tests
+~~~~~~~~~~~~~
+
+All tests should go into the ``tests`` directory. This folder contains many
+current examples of tests, and we suggest looking to these for inspiration.
+
+Running the test suite
+~~~~~~~~~~~~~~~~~~~~~~
+
+The tests can then be run directly inside your Git clone (without having to
+install *momepy*) by typing::
+
+    pytest
+
+6. Updating the Documentation and User Guide
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+*momepy* documentation resides in the `docs` folder. Changes to the docs are
+make by modifying the appropriate file within `doc`.
+*momepy* docs us reStructuredText syntax, `which is explained here <http://www.sphinx-doc.org/en/stable/rest.html#rst-primer>`_
+and the docstrings follow the `Numpy Docstring standard <https://github.com/numpy/numpy/blob/master/doc/HOWTO_DOCUMENT.rst.txt>`_.
+
+Once you have made your changes, you may try if they render correctly by building the docs using sphinx.
+To do so, you can navigate to the doc folder and type::
+
+    make html
+
+The resulting html pages will be located in doc/build/html. In case of any errors,
+you can try to use make html within a new environment based on environment.yml specification in the doc folder.
+Using conda::
+
+    conda env create -f environment.yml
+    conda activate geopandas_docs
+    make html
+
+For minor updates, you can skip whole make html part as reStructuredText syntax is
+usually quite straightforward.
+
+Updating User Guide
+~~~~~~~~~~~~~~~~~~~
+
+Updating user guide might be slightly more complicated as it is powered by `Jupyter book`_
+and consists of collection of markdown files and Jupyter notebooks. Moreover, it sits in
+its own separate branch `gh-pages`. To make changes, you need to fork this branch first,
+not master.
+Changes in markdown are straightforward, changes in notebooks should be done using Jupyter.
+
+Once you are done with changes, you should navigate to parental folder of momepy git folder
+and run::
+
+    jupyter-book build momepy/
+
+Then you can push your fork to GitHub.
+
+7. Submitting a Pull Request
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+Once you've made changes and pushed them to your forked repository, you then
+submit a pull request to have them integrated into the *momepy* code base.
+
+You can find a pull request (or PR) tutorial in the `GitHub's Help Docs <https://help.github.com/articles/using-pull-requests/>`_.
+
+
+.. _geopandas: https://geopandas.org/
+
+.. _Jupyter book: https://jupyter.org/jupyter-book/intro.html
+
+.. _sphinx: https://www.sphinx-doc.org/
diff --git a/docs/index.rst b/docs/index.rst
index 9d9118c..9443b36 100644
--- a/docs/index.rst
+++ b/docs/index.rst
@@ -15,9 +15,7 @@ frameworks for a description of urban structures.
 
     *momepy* stands for Morphological Measuring in Python
 
-The development of momepy is still in its initial phase. The whole
-structure of the package can always be changed, until the release of
-first versions on PIP. Momepy is a result of ongoing research of `Urban
+Momepy is a result of ongoing research of `Urban
 Design Studies Unit (UDSU) <http://udsu-strath.com>`__ supported by the
 Axel and Margaret Ax:son Johnson Foundation as a part of “The Urban Form
 Resilience Project” in partnership with University of Strathclyde in Glasgow,
@@ -31,34 +29,48 @@ https://github.com/martinfleis/momepy
 Install
 -------
 
-Momepy is not yet released. At this moment, you can either download the GitHub
-repository to work with ``momepy`` or install the latest development version
-from GitHub::
+See the :doc:`installation docs <install>` for detailed instructions.
+Momepy depends on python geospatial stack, which might cause some dependency
+issues.
 
-  pip install git+git://github.com/martinfleis/momepy.git
+Contributing to momepy
+----------------------
 
-Please ensure that you have all dependencies installed as ``momepy`` setup
-currently does not specify dependency requirements.
-Momepy depends on the following packages: ``geopandas``, ``libpysal``
-and ``tqdm``. All of the packages are currently required, which will be
-changed in future.
+Contributions of any kind to momepy are more than welcome. That does not mean
+new code only, but also improvements of documentation and user guide, additional
+tests (ideally filling the gaps in existing suite) or bug report or idea what
+could be added or done better.
 
-References
-----------
+All contributions should go through our GitHub repository. Bug reports, ideas or
+even questions should be raised by opening an issue on the GitHub tracker.
+Suggestions for changes in code or documentation should be submitted as a pull
+request. However, if you are not sure what to do, feel free to open an issue.
+All discussion will then take place on GitHub to keep the development of momepy
+transparent.
 
--  Dibble J, Prelorendjos A, Romice O, et al. (2017) On the origin of
-   spaces: Morphometric foundations of urban form evolution. Environment
-   and Planning B: Urban Analytics and City Science 24(1):
-   239980831772507–24.
+If you decide to contribute to the codebase, ensure that you are using an
+up-to-date `master` branch. The latest development version will always be there,
+including a significant part of the documentation (powered by `sphinx`). The
+user guide is located in the `gh-pages` branch and is powered by `Jupyter book`.
 
+Details are available in the :doc:`contributing guide <contributing>`.
 
+Get in touch
+------------
+If you have a question regarding momepy, feel free to open an issue on GitHub.
+Eventually, you can contact us on dev@momepy.org.
+
+
+Documentation contents
+----------------------
 .. toctree::
    :maxdepth: 1
-   :caption: Contents:
 
    install
    tutorial
    api
+   contributing
+   references
 
 
 Indices and tables
diff --git a/docs/install.rst b/docs/install.rst
index 3e9f63b..fa8af8e 100644
--- a/docs/install.rst
+++ b/docs/install.rst
@@ -1,26 +1,101 @@
 Install
 =======
 
-Momepy is not yet released. At this moment, you can either download the GitHub
-repository to work with ``momepy`` or install the latest development version
-from GitHub via PIP::
+.. note::
 
-    pip install git+git://github.com/martinfleis/momepy.git
+    MOMEPY IS NOT YET AVAILABLE AT CONDA-FORGE. INSTRUCTIONS DOES NOT APPLY.
 
-Please ensure that you have all dependencies installed as ``momepy`` setup
-currently does not specify dependency requirements.
-Momepy depends on the following packages: ``geopandas``, ``libpysal``
-and ``tqdm``. All of the packages are currently required, which will be
-changed in future.
+Momepy, similar to GeoPandas, can be a bit complicated to install. However,
+if you follow recommended instructions below, there should be no issue. For
+more details on issues with geospatial python stack, please refer to `GeoPandas
+installation instructions <http://geopandas.org/install.html>`__.
 
 Install via Conda
 -----------------
 
-As `momepy` is dependent on `geopandas` and other spatial packages, we recommend
-to install all dependencies via `conda` from `conda-forge`::
+As ``momepy`` is dependent on `geopandas`_ and other spatial packages, we recommend
+to install all dependencies via `conda`_ from `conda-forge`_::
 
-    conda config --add channels conda-forge  # ensure that conda-forge is in your list of channels on the top
-    conda config --set channel_priority strict  # strict priority to install all from conda-forge
-    conda create --name momepy_env geopandas libpysal tqdm
-    conda activate momepy_env
-    pip install git+git://github.com/martinfleis/momepy.git
+    conda install -c conda-forge momepy
+
+Conda should be able to resolve any dependency conflicts and install momepy
+together with all necessary dependencies.
+
+.. note::
+
+    We strongly recommend to install everything from the *conda-forge* channel.
+    Mixture of conda channels or conda and pip packages can lead to import problems.
+
+
+Creating a new environment for momepy
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+If you want to make sure, that everything will work as it should, you can create
+a new conda environment for momepy. Assuming we want to create a new environment
+called `momepy_env`::
+
+    conda create -n momepy_env -c conda-forge -c defaults --override-channels python=3 momepy
+
+
+Install via pip
+---------------
+
+Momepy is also available on PyPI, but ensure that all dependencies are properly
+installed before installing momepy. Some C dependencies are causing problems with
+installing using pip only::
+
+    pip install momepy
+
+Install from the repository
+---------------------------
+
+If you want to work with the latest development version of momepy, you can do so
+by cloning `GitHub repository <https://github.com/martinfleis/momepy>`__ and
+installing momepy from local directory::
+
+    git clone https://github.com/martinfleis/momepy.git
+    cd momepy
+    pip install .
+
+Alternatively, you can install the latest version directly from GitHub::
+
+    pip install git=git://github.com/martinfleis/momepy.git
+
+Installing directly from repository might face the same dependency issues as
+described above regarding installing using pip. To ensure that environment is
+properly prepared and every dependency will work as intended, you can install
+them using conda before installing development version of momepy::
+
+    conda install -c conda-forge geopandas, networkx, libpysal, tqdm
+
+
+Dependencies
+------------
+
+Required dependencies:
+
+- `geopandas`_
+- `libpysal`_
+- `networkx`_
+- `tqdm`_
+
+Some functions also depend on additional packages, which are optional:
+
+- `pysal`_
+- `mapclassify`_ (>= 2.1.1)
+
+.. _geopandas: https://geopandas.org/
+
+.. _mapclassify: https://mapclassify.readthedocs.io/en/latest/
+
+.. _libpysal: https://libpysal.readthedocs.io
+
+.. _networkx: http://networkx.github.io
+
+.. _tqdm: http://networkx.github.io
+
+.. _pysal: http://pysal.org
+
+.. _conda-forge: https://conda-forge.org/
+
+.. _conda: https://conda.io/en/latest/
diff --git a/docs/references.rst b/docs/references.rst
new file mode 100644
index 0000000..1671a8d
--- /dev/null
+++ b/docs/references.rst
@@ -0,0 +1,7 @@
+References
+==========
+
+-  Dibble J, Prelorendjos A, Romice O, et al. (2017) On the origin of
+   spaces: Morphometric foundations of urban form evolution. Environment
+   and Planning B: Urban Analytics and City Science 24(1):
+   239980831772507–24.
diff --git a/momepy/__init__.py b/momepy/__init__.py
index cc3fa22..5497519 100644
--- a/momepy/__init__.py
+++ b/momepy/__init__.py
@@ -7,3 +7,7 @@ from .distribution import *
 from .graph import *
 from .shape import *
 import momepy.datasets
+
+__author__ = 'Martin Fleischmann'
+__author_email__ = 'martin@martinfleischmann.net'
+__version__ = '0.1rc1'
diff --git a/momepy/dimension.py b/momepy/dimension.py
index 4dd8735..d42d6dd 100644
--- a/momepy/dimension.py
+++ b/momepy/dimension.py
@@ -17,7 +17,7 @@ def area(gdf):
     Calculates area of each object in given shapefile. It can be used for any
     suitable element (building footprint, plot, tessellation, block).
 
-    It is a simple wrapper for geopandas gdf.geometry.area for consistency of momepy.
+    It is a simple wrapper for geopandas `gdf.geometry.area` for consistency of momepy.
 
     Parameters
     ----------
@@ -53,7 +53,7 @@ def perimeter(gdf):
     Calculates perimeter of each object in given shapefile. It can be used for any
     suitable element (building footprint, plot, tessellation, block).
 
-    It is a simple wrapper for geopandas gdf.geometry.length for consistency of momepy.
+    It is a simple wrapper for geopandas `gdf.geometry.length` for consistency of momepy.
 
     Parameters
     ----------
@@ -168,13 +168,13 @@ def floor_area(gdf, heights, areas=None):
 
     Examples
     --------
-    >>> buildings['floor_area'] = momepy.volume(buildings, heights='height_col')
+    >>> buildings['floor_area'] = momepy.floor_area(buildings, heights='height_col')
     Calculating floor areas...
     Floor areas calculated.
     >>> buildings.floor_area[0]
     2185.672484113309
 
-    >>> buildings['floor_area'] = momepy.volume(buildings, heights='height_col', areas='area_col')
+    >>> buildings['floor_area'] = momepy.floor_area(buildings, heights='height_col', areas='area_col')
     Calculating floor areas...
     Floor areas calculated.
     >>> buildings.floor_area[0]
@@ -354,7 +354,7 @@ def mean_character(gdf, values, spatial_weights, unique_id, rng=None):
             values_list = limit_range(values_list.tolist(), rng=rng)
         results_list.append(np.mean(values_list))
 
-    series = pd.Series(results_list)
+    series = pd.Series(results_list, index=gdf.index)
 
     print('Mean character value calculated.')
     return series
@@ -362,7 +362,13 @@ def mean_character(gdf, values, spatial_weights, unique_id, rng=None):
 
 def street_profile(left, right, heights=None, distance=10, tick_length=50):
     """
-    Calculates the street profile widths, standard deviation of width, heights, and ratio height/width
+    Calculates the street profile characters.
+
+    Returns a dictionary with widths, standard deviation of width, openness, heights,
+    standard deviation of height and ratio height/width. Algorithm generates perpendicular
+    lines to `right` dataframe features every `distance` and measures values on intersection
+    with features of `left.` If no feature is reached within
+    `tick_length` its value is set as width (being theoretical maximum).
 
     .. math::
         \\
@@ -372,13 +378,13 @@ def street_profile(left, right, heights=None, distance=10, tick_length=50):
     left : GeoDataFrame
         GeoDataFrame containing streets to analyse
     right : GeoDataFrame
-        GeoDataFrame containing buildings along the streets
+        GeoDataFrame containing buildings along the streets (only Polygon geometry type is supported)
     heights: str, list, np.array, pd.Series (default None)
         the name of the buildings dataframe column, np.array, or pd.Series where is stored building height. If set to None,
         height and ratio height/width will not be calculated.
-    distance : int, optional
+    distance : int (default 10)
         distance between perpendicular ticks
-    tick_length : int, optional
+    tick_length : int (default 50)
         lenght of ticks
 
     Returns
@@ -401,7 +407,11 @@ def street_profile(left, right, heights=None, distance=10, tick_length=50):
     References
     ----------
     Oliveira V (2013) Morpho: a methodology for assessing urban form. Urban Morphology 17(1): 21–33.
-    Araldi and Fusco...
+
+    Araldi A and Fusco G (2017) Decomposing and Recomposing Urban Fabric: The City from the Pedestrian
+    Point of View. In: Gervasi O, Murgante B, Misra S, et al. (eds), Computational Science and Its
+    Applications – ICCSA 2017, Lecture Notes in Computer Science, Cham: Springer International
+    Publishing, pp. 365–376. Available from: http://link.springer.com/10.1007/978-3-319-62407-5.
 
     Examples
     --------
@@ -411,11 +421,6 @@ def street_profile(left, right, heights=None, distance=10, tick_length=50):
     Street profile calculated.
     >>> streets_df['width'] = street_profile['widths']
     >>> streets_df['deviations'] = street_profile['width_deviations']
-
-    Notes
-    -----
-    Add explanation of algorithm
-
     """
 
     print('Calculating street profile...')
@@ -521,7 +526,7 @@ def street_profile(left, right, heights=None, distance=10, tick_length=50):
                 possible_intersections_index = list(sindex.intersection(tick.bounds))
                 possible_intersections = right.iloc[possible_intersections_index]
                 real_intersections = possible_intersections.intersects(tick)
-                get_height = right.iloc[list(real_intersections.index)]
+                get_height = right.loc[list(real_intersections.index)]
                 possible_int = get_height.exterior.intersection(tick)
 
                 if possible_int.any():
@@ -556,7 +561,7 @@ def street_profile(left, right, heights=None, distance=10, tick_length=50):
                             dist = row.geometry.distance(Point(tick.coords[-1]))
                             indices[idx] = dist
                         minim = min(indices, key=indices.get)
-                        m_heights.append(right.iloc[minim][heights])
+                        m_heights.append(right.loc[minim][heights])
 
         openness = (len(lefts) + len(rights)) / len(ticks * 2)
         openness_list.append(1 - openness)
@@ -582,13 +587,13 @@ def street_profile(left, right, heights=None, distance=10, tick_length=50):
                 heights_deviations_list.append(0)
 
     street_profile = {}
-    street_profile['widths'] = pd.Series(results_list)
-    street_profile['width_deviations'] = pd.Series(deviations_list)
-    street_profile['openness'] = pd.Series(openness_list)
+    street_profile['widths'] = pd.Series(results_list, index=left.index)
+    street_profile['width_deviations'] = pd.Series(deviations_list, index=left.index)
+    street_profile['openness'] = pd.Series(openness_list, index=left.index)
 
     if heights is not None:
-        street_profile['heights'] = pd.Series(heights_list)
-        street_profile['heights_deviations'] = pd.Series(heights_deviations_list)
+        street_profile['heights'] = pd.Series(heights_list, index=left.index)
+        street_profile['heights_deviations'] = pd.Series(heights_deviations_list, index=left.index)
         street_profile['profile'] = street_profile['heights'] / street_profile['widths']
 
     print('Street profile calculated.')
@@ -626,19 +631,14 @@ def weighted_character(gdf, values, spatial_weights, unique_id, areas=None):
 
     References
     ----------
-    Jacob
+    Dibble J, Prelorendjos A, Romice O, et al. (2017) On the origin of spaces: Morphometric foundations of urban form evolution.
+    Environment and Planning B: Urban Analytics and City Science 46(4): 707–730.
 
     Examples
     --------
-    >>> buildings_df['w_height'] = momepy.weighted_character(buildings_df, tessellation_df, characters='height', unique_id='uID')
-    Generating weights matrix (Queen) of 3 topological steps...
-    Calculating weighted height...
-    100%|██████████| 144/144 [00:00<00:00, 385.52it/s]
-    Weighted height calculated.
-
     >>> sw = libpysal.weights.DistanceBand.from_dataframe(tessellation_df, threshold=100, silence_warnings=True)
-    >>> buildings_df['w_height_100'] = momepy.weighted_character(buildings_df, tessellation_df, characters='height',
-                                                                 unique_id='uID', spatial_weights=sw)
+    >>> buildings_df['w_height_100'] = momepy.weighted_character(buildings_df, values='height', spatial_weights=sw,
+                                                                 unique_id='uID')
     Calculating weighted height...
     100%|██████████| 144/144 [00:00<00:00, 361.60it/s]
     Weighted height calculated.
@@ -663,18 +663,15 @@ def weighted_character(gdf, values, spatial_weights, unique_id, areas=None):
             building_neighbours = gdf.loc[gdf[unique_id].isin(neighbours)]
 
             if areas is not None:
-                results_list.append((sum(building_neighbours[values]
-                                         * building_neighbours[areas])
-                                    + (row[values] * row[areas]))
-                                    / (sum(building_neighbours[areas]) + row[areas]))
+                results_list.append((sum(building_neighbours[values] * building_neighbours[
+                                    areas]) + (row[values] * row[areas])) / (sum(building_neighbours[areas]) + row[areas]))
             else:
-                results_list.append((sum(building_neighbours[values]
-                                         * building_neighbours.geometry.area)
-                                    + (row[values] * row.geometry.area))
-                                    / (sum(building_neighbours.geometry.area) + row.geometry.area))
+                results_list.append((sum(building_neighbours[values] * building_neighbours.geometry.area
+                                         ) + (row[values] * row.geometry.area)
+                                     ) / (sum(building_neighbours.geometry.area) + row.geometry.area))
         else:
             results_list.append(row[values])
-    series = pd.Series(results_list)
+    series = pd.Series(results_list, index=gdf.index)
 
     print('Weighted {} calculated.'.format(values))
     return series
@@ -693,7 +690,7 @@ def covered_area(gdf, spatial_weights, unique_id):
     ----------
     gdf : GeoDataFrame
         GeoDataFrame containing morphological tessellation
-    spatial_weights : libpysal.weights, optional
+    spatial_weights : libpysal.weights
         spatial weights matrix
     unique_id : str
         name of the column with unique id used as spatial_weights index.
@@ -703,11 +700,13 @@ def covered_area(gdf, spatial_weights, unique_id):
     Series
         Series containing resulting values.
 
-    References
-    ----------
-
     Examples
     --------
+    >>> sw = momepy.Queen_higher(k=3, geodataframe=tessellation_df, ids='uID')
+    >>> tessellation_df['covered3steps'] = mm.covered_area(tessellation_df, sw, 'uID')
+    Calculating covered area...
+    100%|██████████| 144/144 [00:00<00:00, 549.15it/s]
+    Covered area calculated.
 
     """
     # define empty list for results
@@ -721,7 +720,7 @@ def covered_area(gdf, spatial_weights, unique_id):
         areas = gdf.loc[gdf[unique_id].isin(neighbours)].geometry.area
         results_list.append(sum(areas))
 
-    series = pd.Series(results_list)
+    series = pd.Series(results_list, index=gdf.index)
 
     print('Covered area calculated.')
     return series
@@ -744,9 +743,18 @@ def wall(gdf, spatial_weights=None):
     Series
         Series containing resulting values.
 
+    Examples
+    --------
+    >>> buildings_df['wall_length'] = mm.wall(buildings_df)
+    Calculating perimeter wall length...
+    Calculating spatial weights...
+    Spatial weights ready...
+    100%|██████████| 144/144 [00:00<00:00, 4171.39it/s]
+    Perimeter wall length calculated.
+
     Notes
     -----
-    Script is not optimised at all, so it is currently extremely slow.
+    It might take a while to compute this character.
     """
     # define empty list for results
     results_list = []
@@ -791,7 +799,7 @@ def wall(gdf, spatial_weights=None):
     for index, row in tqdm(gdf.iterrows(), total=gdf.shape[0]):
         results_list.append(walls[index])
 
-    series = pd.Series(results_list)
+    series = pd.Series(results_list, index=gdf.index)
     print('Perimeter wall length calculated.')
     return series
 
@@ -810,7 +818,7 @@ def segments_length(gdf, spatial_weights=None, mean=False):
         GeoDataFrame containing streets (edges) to analyse
     spatial_weights : libpysal.weights, optional
         spatial weights matrix - If None, Queen contiguity matrix will be calculated
-        based on streets (note: based on index, not ID).
+        based on streets (note: spatial_weights shoud be based on index, not unique ID).
     mean : boolean, optional
         If mean=False it will return total length, if mean=True it will return mean value
 
@@ -818,6 +826,14 @@ def segments_length(gdf, spatial_weights=None, mean=False):
     -------
     Series
         Series containing resulting values.
+
+    Examples
+    --------
+    >>> streets_df['length_neighbours'] = mm.segments_length(streets_df, mean=True)
+    Calculating segments length...
+    Calculating spatial weights...
+    Spatial weights ready...
+    Segments length calculated.
     """
     results_list = []
 
@@ -838,6 +854,6 @@ def segments_length(gdf, spatial_weights=None, mean=False):
         else:
             results_list.append(sum(dims))
 
-    series = pd.Series(results_list)
+    series = pd.Series(results_list, index=gdf.index)
     print('Segments length calculated.')
     return series
diff --git a/momepy/distribution.py b/momepy/distribution.py
index 68a30bb..6c42f97 100644
--- a/momepy/distribution.py
+++ b/momepy/distribution.py
@@ -13,7 +13,7 @@ import statistics
 
 def orientation(gdf):
     """
-    Calculate orientation (azimuth) of object
+    Calculate the orientation of object
 
     Defined as an orientation of the longext axis of bounding rectangle in range 0 - 45.
     It captures the deviation of orientation from cardinal directions.
@@ -101,7 +101,7 @@ def orientation(gdf):
                 az = az - 2 * diff
             results_list.append(az)
 
-    series = pd.Series(results_list)
+    series = pd.Series(results_list, index=gdf.index)
     print('Orientations calculated.')
     return series
 
@@ -176,7 +176,7 @@ def shared_walls_ratio(gdf, unique_id, perimeters=None):
                 subset = gdf.loc[i]['geometry']
                 length = length + row.geometry.intersection(subset).length
             results_list.append(length / row[perimeters])
-    series = pd.Series(results_list)
+    series = pd.Series(results_list, index=gdf.index)
     print('Shared walls ratio calculated.')
 
     return series
@@ -269,13 +269,13 @@ def street_alignment(left, right, orientations, left_network_id, right_network_i
                 diff = az - 45
                 az = az - 2 * diff
             results_list.append(abs(row[orientations] - az))
-    series = pd.Series(results_list)
+    series = pd.Series(results_list, index=left.index)
 
     print('Street alignments calculated.')
     return series
 
 
-def cell_alignment(left, right, left_orientations, right_orientations, unique_id):
+def cell_alignment(left, right, left_orientations, right_orientations, left_unique_id, right_unique_id):
     """
     Calculate the difference between cell orientation and orientation of object
 
@@ -294,9 +294,10 @@ def cell_alignment(left, right, left_orientations, right_orientations, unique_id
     right_orientations : str, list, np.array, pd.Series
         the name of the dataframe column, np.array, or pd.Series where is stored object orientation value
         (can be calculated using :py:func:`momepy.orientation`)
-    unique_id : str
-        the name of the dataframe column with unique id shared between a cell and a building
-        (must be present in both geodataframes)
+    left_unique_id : str
+        the name of the left dataframe column with unique id shared between left and right gdf
+    right_unique_id : str
+        the name of the right dataframe column with unique id shared between left and right gdf
 
     Returns
     -------
@@ -305,16 +306,13 @@ def cell_alignment(left, right, left_orientations, right_orientations, unique_id
 
     Examples
     --------
-    >>> buildings_df['cell_alignment'] = momepy.cell_alignment(buildings_df, tessellation_df, 'bl_orient', 'tes_orient', 'uID')
+    >>> buildings_df['cell_alignment'] = momepy.cell_alignment(buildings_df, tessellation_df, 'bl_orient', 'tes_orient', 'uID', 'uID')
     Calculating cell alignments...
     100%|██████████| 144/144 [00:00<00:00, 799.09it/s]
     Cell alignments calculated.
     >>> buildings_df['cell_alignment'][0]
     0.8795123936951939
 
-    Notes
-    -----
-    Allow left unique_id and right unique_id.
     """
     print('Calculating cell alignments...')
     left = left.copy()
@@ -331,9 +329,9 @@ def cell_alignment(left, right, left_orientations, right_orientations, unique_id
 
     for index, row in tqdm(left.iterrows(), total=left.shape[0]):
 
-        results_list.append(abs(row[left_orientations] - right[right[unique_id] == row[unique_id]][right_orientations].iloc[0]))
+        results_list.append(abs(row[left_orientations] - right[right[right_unique_id] == row[left_unique_id]][right_orientations].iloc[0]))
 
-    series = pd.Series(results_list)
+    series = pd.Series(results_list, index=left.index)
 
     print('Cell alignments calculated.')
     return series
@@ -396,7 +394,7 @@ def alignment(gdf, spatial_weights, unique_id, orientations):
         else:
             results_list.append(0)
 
-    series = pd.Series(results_list)
+    series = pd.Series(results_list, index=gdf.index)
 
     print('Alignments calculated.')
     return series
@@ -426,7 +424,7 @@ def neighbour_distance(gdf, spatial_weights, unique_id):
     References
     ---------
     Schirmer PM and Axhausen KW (2015) A multiscale classiﬁcation of urban morphology.
-    Journal of Transport and Land Use 9(1): 101–130.
+    Journal of Transport and Land Use 9(1): 101–130. (adapted)
 
     Examples
     --------
@@ -451,7 +449,7 @@ def neighbour_distance(gdf, spatial_weights, unique_id):
         else:
             results_list.append(0)
 
-    series = pd.Series(results_list)
+    series = pd.Series(results_list, index=gdf.index)
 
     print('Distances calculated.')
     return series
@@ -548,7 +546,7 @@ def mean_interbuilding_distance(gdf, spatial_weights, unique_id, spatial_weights
             selection = adj_list[adj_list.focal.isin(neighbours)][adj_list.neighbor.isin(neighbours)]
             results_list.append(np.nanmean(selection.distance))
 
-    series = pd.Series(results_list)
+    series = pd.Series(results_list, index=gdf.index)
     print('Mean interbuilding distances calculated.')
     return series
 
@@ -617,7 +615,7 @@ def neighbouring_street_orientation_deviation(gdf):
             diff = az - 45
             az = az - 2 * diff
         results_list.append(az)
-    series = pd.Series(results_list)
+    series = pd.Series(results_list, index=gdf.index)
 
     gdf['tmporient'] = series
 
@@ -645,7 +643,7 @@ def neighbouring_street_orientation_deviation(gdf):
         else:
             results_list.append(0)
 
-    series = pd.Series(results_list)
+    series = pd.Series(results_list, index=gdf.index)
     print('Street alignments calculated.')
     return series
 
@@ -749,7 +747,7 @@ def building_adjacency(gdf, spatial_weights_higher, unique_id, spatial_weights=N
         else:
             results_list.append(0)
 
-    series = pd.Series(results_list)
+    series = pd.Series(results_list, index=gdf.index)
 
     print('Adjacency calculated.')
     return series
@@ -789,9 +787,8 @@ def neighbours(gdf, spatial_weights, unique_id, weighted=False):
 
     Examples
     --------
-    >>> tessellation_df['neighbours'] = momepy.neighbours(tessellation_df)
-    Calculating spatial weights...
-    Spatial weights ready...
+    >>> sw = libpysal.weights.contiguity.Queen.from_dataframe(tessellation_df, ids='uID')
+    >>> tessellation_df['neighbours'] = momepy.neighbours(tessellation_df, sw, 'uID')
     Calculating neighbours...
     100%|██████████| 144/144 [00:00<00:00, 6909.50it/s]
     Neighbours calculated.
@@ -807,6 +804,6 @@ def neighbours(gdf, spatial_weights, unique_id, weighted=False):
         else:
             neighbours.append(spatial_weights.cardinalities[row[unique_id]])
 
-    series = pd.Series(neighbours)
+    series = pd.Series(neighbours, index=gdf.index)
     print('Neighbours calculated.')
     return series
diff --git a/momepy/diversity.py b/momepy/diversity.py
index 06a5e9e..aebbba0 100644
--- a/momepy/diversity.py
+++ b/momepy/diversity.py
@@ -7,11 +7,12 @@
 from tqdm import tqdm  # progress bar
 import pandas as pd
 import scipy as sp
+import numpy as np
 
 
 def rng(gdf, values, spatial_weights, unique_id, rng=(0, 100), **kwargs):
     """
-    Calculates the range of values within k steps of morphological tessellation
+    Calculates the range of values within neighbours defined in `spatial_weights`.
 
     Uses `scipy.stats.iqr` under the hood.
 
@@ -41,9 +42,17 @@ def rng(gdf, values, spatial_weights, unique_id, rng=(0, 100), **kwargs):
 
     References
     ----------
+    Dibble J, Prelorendjos A, Romice O, et al. (2017) On the origin of spaces: Morphometric foundations of urban form evolution.
+    Environment and Planning B: Urban Analytics and City Science 46(4): 707–730.
 
     Examples
     --------
+    >>> sw = momepy.Queen_higher(k=3, geodataframe=tessellation_df, ids='uID')
+    >>> tessellation_df['area_IQR_3steps'] = mm.rng(tessellation_df, 'area', sw, 'uID', rng=(25, 75))
+    Calculating range...
+    100%|██████████| 144/144 [00:00<00:00, 722.50it/s]
+    Range calculated.
+
 
     """
     # define empty list for results
@@ -66,17 +75,17 @@ def rng(gdf, values, spatial_weights, unique_id, rng=(0, 100), **kwargs):
 
         results_list.append(sp.stats.iqr(values_list, rng=rng, **kwargs))
 
-    series = pd.Series(results_list)
+    series = pd.Series(results_list, index=gdf.index)
 
     print('Range calculated.')
     return series
 
 
-def theil(gdf, values, spatial_weights, unique_id, rng=(0, 100)):
+def theil(gdf, values, spatial_weights, unique_id, rng=None):
     """
-    Calculates the Theil measure of inequality of values within k steps of morphological tessellation
+    Calculates the Theil measure of inequality of values within neighbours defined in `spatial_weights`.
 
-    Uses `pysal.explore.inequality.theil.Theil` under the hood.
+    Uses `pysal.explore.inequality.theil.Theil` under the hood. Requires `pysal` dependency.
 
     .. math::
 
@@ -100,14 +109,20 @@ def theil(gdf, values, spatial_weights, unique_id, rng=(0, 100)):
     Series
         Series containing resulting values.
 
-    References
-    ----------
-
     Examples
     --------
-
+    >>> sw = momepy.Queen_higher(k=3, geodataframe=tessellation_df, ids='uID')
+    >>> tessellation_df['area_Theil'] = mm.theil(tessellation_df, 'area', sw, 'uID')
+    Calculating Theil index...
+    100%|██████████| 144/144 [00:00<00:00, 597.37it/s]
+    Theil index calculated.
     """
-    from pysal.explore.inequality.theil import Theil
+    try:
+        from pysal.explore.inequality.theil import Theil
+    except ImportError:
+        raise ImportError(
+            "The 'pysal' package is required.")
+
     # define empty list for results
     results_list = []
     gdf = gdf.copy()
@@ -132,7 +147,7 @@ def theil(gdf, values, spatial_weights, unique_id, rng=(0, 100)):
             values_list = limit_range(values_list.tolist(), rng=rng)
         results_list.append(Theil(values_list).T)
 
-    series = pd.Series(results_list)
+    series = pd.Series(results_list, index=gdf.index)
 
     print('Theil index calculated.')
     return series
@@ -140,9 +155,9 @@ def theil(gdf, values, spatial_weights, unique_id, rng=(0, 100)):
 
 def simpson(gdf, values, spatial_weights, unique_id, binning='HeadTailBreaks', **classification_kwds):
     """
-    Calculates the Simpson\'s diversity index of values within k steps of morphological tessellation
+    Calculates the Simpson\'s diversity index of values within neighbours defined in `spatial_weights`.
 
-    Uses `mapclassify.classifiers` under the hood for binning.
+    Uses `mapclassify.classifiers` under the hood for binning. Requires `mapclassify>=.2.1.0` dependency.
 
     .. math::
 
@@ -177,10 +192,16 @@ def simpson(gdf, values, spatial_weights, unique_id, binning='HeadTailBreaks', *
 
     References
     ----------
+    Feliciotti A (2018) RESILIENCE AND URBAN DESIGN:A SYSTEMS APPROACH TO THE STUDY OF RESILIENCE
+    IN URBAN FORM. LEARNING FROM THE CASE OF GORBALS. Glasgow.
 
     Examples
     --------
-
+    >>> sw = momepy.Queen_higher(k=3, geodataframe=tessellation_df, ids='uID')
+    >>> tessellation_df['area_Simpson'] = mm.simpson(tessellation_df, 'area', sw, 'uID')
+    Calculating Simpson's diversity index...
+    100%|██████████| 144/144 [00:00<00:00, 455.83it/s]
+    Simpson's diversity index calculated.
     """
     def simpson_di(data):
 
@@ -238,11 +259,92 @@ def simpson(gdf, values, spatial_weights, unique_id, binning='HeadTailBreaks', *
             neighbours = row[unique_id]
         values_list = gdf.loc[gdf[unique_id].isin(neighbours)][values]
 
-        sample_bins = mapclassify.classifiers.User_Defined(values_list, bins)
+        sample_bins = mapclassify.classifiers.UserDefined(values_list, bins)
         counts = dict(zip(bins, sample_bins.counts))
         results_list.append(simpson_di(counts))
 
-    series = pd.Series(results_list)
+    series = pd.Series(results_list, index=gdf.index)
 
     print("Simpson's diversity index calculated.")
     return series
+
+
+def _gini(vals):
+    """Calculate the Gini coefficient of a numpy array."""
+    # based on bottom eq:
+    # http://www.statsdirect.com/help/generatedimages/equations/equation154.svg
+    # from:
+    # http://www.statsdirect.com/help/default.htm#nonparametric_methods/gini.htm
+    # All values are treated equally, arrays must be 1d:
+    vals = vals.flatten()
+    if np.amin(vals) < 0:
+        # Values cannot be negative:
+        vals -= np.amin(vals)
+    # Values cannot be 0:
+    vals += 0.0000001
+    # Values must be sorted:
+    vals = np.sort(vals)
+    # Index per array element:
+    index = np.arange(1, vals.shape[0] + 1)
+    # Number of array elements:
+    n = vals.shape[0]
+    # Gini coefficient:
+    return ((np.sum((2 * index - n - 1) * vals)) / (n * np.sum(vals)))
+
+
+def gini(gdf, values, spatial_weights, unique_id, rng=None):
+    """
+    Calculates the Gini index of values within neighbours defined in `spatial_weights`.
+
+    .. math::
+
+
+    Parameters
+    ----------
+    gdf : GeoDataFrame
+        GeoDataFrame containing morphological tessellation
+    values : str, list, np.array, pd.Series
+        the name of the dataframe column, np.array, or pd.Series where is stored character value.
+    spatial_weights : libpysal.weights
+        spatial weights matrix
+    unique_id : str
+        name of the column with unique id used as spatial_weights index
+    rng : Two-element sequence containing floats in range of [0,100], optional
+        Percentiles over which to compute the range. Each must be
+        between 0 and 100, inclusive. The order of the elements is not important.
+
+    Returns
+    -------
+    Series
+        Series containing resulting values.
+
+    Examples
+    --------
+    >>> sw = momepy.Queen_higher(k=3, geodataframe=tessellation_df, ids='uID')
+    >>> tessellation_df['area_Theil'] = mm.gini(tessellation_df, 'area', sw, 'uID')
+    Calculating Gini index...
+    100%|██████████| 144/144 [00:00<00:00, 597.37it/s]
+    Gini index calculated.
+    """
+    # define empty list for results
+    results_list = []
+    gdf = gdf.copy()
+    print('Calculating Gini index...')
+
+    for index, row in tqdm(gdf.iterrows(), total=gdf.shape[0]):
+        neighbours = spatial_weights.neighbors[row[unique_id]]
+        if neighbours:
+            neighbours.append(row[unique_id])
+
+            values_list = gdf.loc[gdf[unique_id].isin(neighbours)][values].values
+
+            if rng:
+                from momepy import limit_range
+                values_list = np.array(limit_range(values_list, rng=rng))
+            results_list.append(_gini(values_list))
+        else:
+            results_list.append(0)
+    series = pd.Series(results_list, index=gdf.index)
+
+    print('Gini index calculated.')
+    return series
diff --git a/momepy/elements.py b/momepy/elements.py
index edb510e..37050d8 100644
--- a/momepy/elements.py
+++ b/momepy/elements.py
@@ -33,6 +33,12 @@ def buffered_limit(gdf, buffer=100):
     MultiPolygon
         MultiPolygon or Polygon defining the study area
 
+    Examples
+    --------
+    >>> limit = mm.buffered_limit(buildings_df)
+    >>> type(limit)
+    shapely.geometry.polygon.Polygon
+
     """
     study_area = gdf.copy()
     study_area['geometry'] = study_area.buffer(buffer)
@@ -130,27 +136,25 @@ def _regions(voronoi_diagram, unique_id, ids, crs):
     return regions_gdf
 
 
-def _split_lines(polygon, distance, crs):
-    dense = _densify(polygon, distance)
-    boundary = dense.boundary
-
-    def _pair(coords):
-        '''Iterate over pairs in a list -> pair of points '''
-        for i in range(1, len(coords)):
-            yield coords[i - 1], coords[i]
+def _split_lines(polygon, distance):
+    """Split polygon into GeoSeries of lines no longer than `distance`."""
+    list_points = []
+    current_dist = distance  # set the current distance to place the point
 
-    segments = []
+    boundary = polygon.boundary  # make shapely MultiLineString object
     if boundary.type == 'LineString':
-        for seg_start, seg_end in _pair(boundary.coords):
-            segment = LineString([seg_start, seg_end])
-            segments.append(segment)
+        line_length = boundary.length  # get the total length of the line
+        while current_dist < line_length:  # while the current cumulative distance is less than the total length of the line
+            list_points.append(boundary.interpolate(current_dist))  # use interpolate and increase the current distance
+            current_dist += distance
     elif boundary.type == 'MultiLineString':
         for ls in boundary:
-            for seg_start, seg_end in _pair(ls.coords):
-                segment = LineString([seg_start, seg_end])
-                segments.append(segment)
+            line_length = ls.length  # get the total length of the line
+            while current_dist < line_length:  # while the current cumulative distance is less than the total length of the line
+                list_points.append(ls.interpolate(current_dist))  # use interpolate and increase the current distance
+                current_dist += distance
 
-    cutted = gpd.GeoSeries(segments, crs=crs)
+    cutted = shapely.ops.split(boundary, shapely.geometry.MultiPoint(list_points).buffer(0.001))
     return cutted
 
 
@@ -160,15 +164,15 @@ def _cut(tessellation, limit, unique_id):
 
     ADD: add option to delete everything outside of limit. Now it keeps it.
     """
-    # cut infinity of voronoi by set buffer (thanks for script to Geoff Boeing)
-    print('Preparing buffer zone for edge resolving...')
-    geometry_cut = _split_lines(limit, 100, tessellation.crs)
+    print('Preparing limit for edge resolving...')
+    geometry_cut = _split_lines(limit, 100)
 
     print('Building R-tree...')
     sindex = tessellation.sindex
     # find the points that intersect with each subpolygon and add them to points_within_geometry
+    print('Identifying edge cells...')
     to_cut = pd.DataFrame()
-    for poly in geometry_cut:
+    for poly in tqdm(geometry_cut, total=(len(geometry_cut))):
         # find approximate matches with r-tree, then precise matches from those approximate ones
         possible_matches_index = list(sindex.intersection(poly.bounds))
         possible_matches = tessellation.iloc[possible_matches_index]
@@ -324,6 +328,30 @@ def tessellation(gdf, unique_id, limit, shrink=0.4, segment=0.5, queen_corners=F
     GeoDataFrame
         GeoDataFrame of morphological tessellation with the unique id based on original buildings.
 
+    Examples
+    --------
+    >>> tessellation_df = mm.tessellation(buildings_df, 'uID', limit=mm.buffered_limit(buildings_df))
+    Bufferring geometry...
+    Converting multipart geometry to singlepart...
+    Densifying geometry...
+    Generating input point array...
+    100%|██████████| 144/144 [00:00<00:00, 376.15it/s]
+    Generating Voronoi diagram...
+    Generating GeoDataFrame...
+    Vertices to Polygons: 100%|██████████| 33059/33059 [00:01<00:00, 31532.72it/s]
+    Dissolving Voronoi polygons...
+    Preparing buffer zone for edge resolving...
+    Building R-tree...
+    100%|██████████| 42/42 [00:00<00:00, 752.54it/s]
+    Cutting...
+    >>> tessellation_df2.head()
+        uID	geometry
+    0	1	POLYGON ((1603586.677274485 6464344.667944215,...
+    1	2	POLYGON ((1603048.399497852 6464176.180701573,...
+    2	3	POLYGON ((1603071.342637536 6464158.863329805,...
+    3	4	POLYGON ((1603055.834005827 6464093.614718676,...
+    4	5	POLYGON ((1603106.417554705 6464130.215958447,...
+
     Notes
     -------
     queen_corners and sensitivity are currently experimental only and can cause errors.
@@ -376,220 +404,6 @@ def tessellation(gdf, unique_id, limit, shrink=0.4, segment=0.5, queen_corners=F
     return morphological_tessellation
 
 
-def snap_street_network_edge(edges, buildings, tessellation, tolerance_street, tolerance_edge):
-    """
-    Fix street network before performing blocks()
-
-    Extends unjoined ends of street segments to join with other segmets or tessellation boundary.
-
-    Parameters
-    ----------
-    edges : GeoDataFrame
-        GeoDataFrame containing street network
-    buildings : GeoDataFrame
-        GeoDataFrame containing building footprints
-    tessellation : GeoDataFrame
-        GeoDataFrame containing morphological tessellation
-    tolerance_street : float
-        tolerance in snapping to street network (by how much could be street segment extended).
-    tolerance_edge : float
-        tolerance in snapping to edge of tessellated area (by how much could be street segment extended).
-
-    Returns
-    -------
-    GeoDataFrame
-        GeoDataFrame of extended street network.
-
-    """
-    # extrapolating function - makes line as a extrapolation of existing with set length (tolerance)
-    def getExtrapoledLine(p1, p2, tolerance):
-        """
-        Creates a line extrapoled in p1->p2 direction.
-        """
-        EXTRAPOL_RATIO = tolerance  # length of a line
-        a = p2
-
-        # defining new point based on the vector between existing points
-        if p1[0] >= p2[0] and p1[1] >= p2[1]:
-            b = (p2[0] - EXTRAPOL_RATIO * math.cos(math.atan(math.fabs(p1[1] - p2[1] + 0.000001) / math.fabs(p1[0] - p2[0] + 0.000001))),
-                 p2[1] - EXTRAPOL_RATIO * math.sin(math.atan(math.fabs(p1[1] - p2[1] + 0.000001) / math.fabs(p1[0] - p2[0] + 0.000001))))
-        elif p1[0] <= p2[0] and p1[1] >= p2[1]:
-            b = (p2[0] + EXTRAPOL_RATIO * math.cos(math.atan(math.fabs(p1[1] - p2[1] + 0.000001) / math.fabs(p1[0] - p2[0] + 0.000001))),
-                 p2[1] - EXTRAPOL_RATIO * math.sin(math.atan(math.fabs(p1[1] - p2[1] + 0.000001) / math.fabs(p1[0] - p2[0] + 0.000001))))
-        elif p1[0] <= p2[0] and p1[1] <= p2[1]:
-            b = (p2[0] + EXTRAPOL_RATIO * math.cos(math.atan(math.fabs(p1[1] - p2[1] + 0.000001) / math.fabs(p1[0] - p2[0] + 0.000001))),
-                 p2[1] + EXTRAPOL_RATIO * math.sin(math.atan(math.fabs(p1[1] - p2[1] + 0.000001) / math.fabs(p1[0] - p2[0] + 0.000001))))
-        else:
-            b = (p2[0] - EXTRAPOL_RATIO * math.cos(math.atan(math.fabs(p1[1] - p2[1] + 0.000001) / math.fabs(p1[0] - p2[0] + 0.000001))),
-                 p2[1] + EXTRAPOL_RATIO * math.sin(math.atan(math.fabs(p1[1] - p2[1] + 0.000001) / math.fabs(p1[0] - p2[0] + 0.000001))))
-        return LineString([a, b])
-
-    # function extending line to closest object within set distance
-    def extend_line(tolerance, idx):
-        """
-        Extends a line geometry withing GeoDataFrame to snap on itself withing tolerance.
-        """
-        if Point(l_coords[-2]).distance(Point(l_coords[-1])) <= 0.001:
-            if len(l_coords) > 2:
-                extra = l_coords[-3:-1]
-            else:
-                return False
-        else:
-            extra = l_coords[-2:]
-        extrapolation = getExtrapoledLine(*extra, tolerance=tolerance)  # we use the last two points
-
-        possible_intersections_index = list(sindex.intersection(extrapolation.bounds))
-        possible_intersections_lines = network.iloc[possible_intersections_index]
-        possible_intersections_clean = possible_intersections_lines.drop(idx, axis=0)
-        possible_intersections = possible_intersections_clean.intersection(extrapolation)
-
-        if possible_intersections.any():
-
-            true_int = []
-            for one in list(possible_intersections.index):
-                if possible_intersections[one].type == 'Point':
-                    true_int.append(possible_intersections[one])
-                elif possible_intersections[one].type == 'MultiPoint':
-                    true_int.append(possible_intersections[one][0])
-                    true_int.append(possible_intersections[one][1])
-
-            if len(true_int) >= 1:
-                if len(true_int) > 1:
-                    distances = {}
-                    ix = 0
-                    for p in true_int:
-                        distance = p.distance(Point(l_coords[-1]))
-                        distances[ix] = distance
-                        ix = ix + 1
-                    minimal = min(distances.items(), key=operator.itemgetter(1))[0]
-                    new_point_coords = true_int[minimal].coords[0]
-                else:
-                    new_point_coords = true_int[0].coords[0]
-
-                l_coords.append(new_point_coords)
-                new_extended_line = LineString(l_coords)
-
-                # check whether the line goes through buildings. if so, ignore it
-                possible_buildings_index = list(bindex.intersection(new_extended_line.bounds))
-                possible_buildings = buildings.iloc[possible_buildings_index]
-                possible_intersections = possible_buildings.intersection(new_extended_line)
-
-                if possible_intersections.any():
-                    pass
-                else:
-                    network.loc[idx, 'geometry'] = new_extended_line
-        else:
-            return False
-
-    # function extending line to closest object within set distance to edge defined by tessellation
-    def extend_line_edge(tolerance, idx):
-        """
-        Extends a line geometry withing GeoDataFrame to snap on the boundary of tessellation withing tolerance.
-        """
-        if Point(l_coords[-2]).distance(Point(l_coords[-1])) <= 0.001:
-            if len(l_coords) > 2:
-                extra = l_coords[-3:-1]
-            else:
-                return False
-        else:
-            extra = l_coords[-2:]
-        extrapolation = getExtrapoledLine(*extra, tolerance)  # we use the last two points
-
-        # possible_intersections_index = list(qindex.intersection(extrapolation.bounds))
-        # possible_intersections_lines = geometry_cut.iloc[possible_intersections_index]
-        possible_intersections = geometry.intersection(extrapolation)
-
-        if possible_intersections.type != 'GeometryCollection':
-
-            true_int = []
-
-            if possible_intersections.type == 'Point':
-                true_int.append(possible_intersections)
-            elif possible_intersections.type == 'MultiPoint':
-                true_int.append(possible_intersections[0])
-                true_int.append(possible_intersections[1])
-
-            if len(true_int) >= 1:
-                if len(true_int) > 1:
-                    distances = {}
-                    ix = 0
-                    for p in true_int:
-                        distance = p.distance(Point(l_coords[-1]))
-                        distances[ix] = distance
-                        ix = ix + 1
-                    minimal = min(distances.items(), key=operator.itemgetter(1))[0]
-                    new_point_coords = true_int[minimal].coords[0]
-                else:
-                    new_point_coords = true_int[0].coords[0]
-
-                l_coords.append(new_point_coords)
-                new_extended_line = LineString(l_coords)
-
-                # check whether the line goes through buildings. if so, ignore it
-                possible_buildings_index = list(bindex.intersection(new_extended_line.bounds))
-                possible_buildings = buildings.iloc[possible_buildings_index]
-                possible_intersections = possible_buildings.intersection(new_extended_line)
-
-                if possible_intersections.any():
-                    pass
-                else:
-                    network.loc[idx, 'geometry'] = new_extended_line
-
-    network = edges.copy()
-    # generating spatial index (rtree)
-    print('Building R-tree for network...')
-    sindex = network.sindex
-    print('Building R-tree for buildings...')
-    bindex = buildings.sindex
-    print('Dissolving tesselation...')
-    geometry = tessellation.geometry.unary_union.boundary
-
-    print('Snapping...')
-    # iterating over each street segment
-    for idx, row in tqdm(network.iterrows(), total=network.shape[0]):
-
-        line = row['geometry']
-        l_coords = list(line.coords)
-        # network_w = network.drop(idx, axis=0)['geometry']  # ensure that it wont intersect itself
-        start = Point(l_coords[0])
-        end = Point(l_coords[-1])
-
-        # find out whether ends of the line are connected or not
-        possible_first_index = list(sindex.intersection(start.bounds))
-        possible_first_matches = network.iloc[possible_first_index]
-        possible_first_matches_clean = possible_first_matches.drop(idx, axis=0)
-        first = possible_first_matches_clean.intersects(start).any()
-
-        possible_second_index = list(sindex.intersection(end.bounds))
-        possible_second_matches = network.iloc[possible_second_index]
-        possible_second_matches_clean = possible_second_matches.drop(idx, axis=0)
-        second = possible_second_matches_clean.intersects(end).any()
-
-        # both ends connected, do nothing
-        if first and second:
-            continue
-        # start connected, extend  end
-        elif first and not second:
-            if extend_line(tolerance_street, idx) is False:
-                extend_line_edge(tolerance_edge, idx)
-        # end connected, extend start
-        elif not first and second:
-            l_coords.reverse()
-            if extend_line(tolerance_street, idx) is False:
-                extend_line_edge(tolerance_edge, idx)
-        # unconnected, extend both ends
-        elif not first and not second:
-            if extend_line(tolerance_street, idx) is False:
-                extend_line_edge(tolerance_edge, idx)
-            l_coords.reverse()
-            if extend_line(tolerance_street, idx) is False:
-                extend_line_edge(tolerance_edge, idx)
-        else:
-            print('Something went wrong.')
-
-    return network
-
-
 def blocks(tessellation, edges, buildings, id_name, unique_id):
     """
     Generate blocks based on buildings, tesselation and street network
@@ -612,17 +426,43 @@ def blocks(tessellation, edges, buildings, id_name, unique_id):
 
     Returns
     -------
-    buildings, cells, blocks : tuple
+    blocks, buildings_ID, cells_ID : tuple
 
-    buildings : GeoDataFrame
-        GeoDataFrame containing buildings with added block ID
-    cells : GeoDataFrame
-        GeoDataFrame containing morphological tessellation with added block ID
     blocks : GeoDataFrame
         GeoDataFrame containing generated blocks
+    buildings_ID : Series
+        Series derived from buildings with block ID
+    cells_ID : Series
+        Series derived from morphological tessellation with block ID
+
+    Examples
+    --------
+    >>> blocks, buildings_df['blockID'], tessellation_df['blockID'] = mm.blocks(tessellation_df, streets_df, buildings_df, 'bID', 'uID')
+    Buffering streets...
+    Generating spatial index...
+    Difference...
+    Defining adjacency...
+    Defining street-based blocks...
+    Defining block ID...
+    Generating centroids...
+    Spatial join...
+    Attribute join (tesselation)...
+    Generating blocks...
+    Multipart to singlepart...
+    Attribute join (buildings)...
+    Attribute join (tesselation)...
+    >>> blocks.head()
+        bID	geometry
+    0	1.0	POLYGON ((1603560.078648818 6464202.366899694,...
+    1	2.0	POLYGON ((1603457.225976106 6464299.454696888,...
+    2	3.0	POLYGON ((1603056.595487018 6464093.903488506,...
+    3	4.0	POLYGON ((1603260.943782872 6464141.327631323,...
+    4	5.0	POLYGON ((1603183.399594798 6463966.109982309,...
+
     """
 
     cells_copy = tessellation.copy()
+    cells_copy = cells_copy[[unique_id, 'geometry']]
 
     print('Buffering streets...')
     street_buff = edges.copy()
@@ -681,21 +521,16 @@ def blocks(tessellation, edges, buildings, id_name, unique_id):
 
     print('Defining street-based blocks...')
     blocks_single = blocks_gdf.dissolve(by='patch')
+    blocks_single.crs = buildings.crs
 
     blocks_single['geometry'] = blocks_single.buffer(0.1)
 
     print('Defining block ID...')  # street based
-    blocks_single[id_name] = None
-    blocks_single[id_name] = blocks_single[id_name].astype('float')
-    b_id = 1
-    for idx, row in tqdm(blocks_single.iterrows(), total=blocks_single.shape[0]):
-        blocks_single.loc[idx, id_name] = b_id
-        b_id = b_id + 1
+    blocks_single[id_name] = range(len(blocks_single))
 
     print('Generating centroids...')
     buildings_c = buildings.copy()
-    buildings_c['geometry'] = buildings_c.representative_point()  # make centroids
-    blocks_single.crs = buildings.crs
+    buildings_c['geometry'] = buildings_c.representative_point()  # make points
 
     print('Spatial join...')
     centroids_tempID = gpd.sjoin(buildings_c, blocks_single, how='left', op='intersects')
@@ -703,22 +538,19 @@ def blocks(tessellation, edges, buildings, id_name, unique_id):
     tempID_to_uID = centroids_tempID[[unique_id, id_name]]
 
     print('Attribute join (tesselation)...')
-    cells_copy = cells_copy.merge(tempID_to_uID, on=unique_id)
+    cells_copy = cells_copy.merge(tempID_to_uID, on=unique_id, how='left')
 
     print('Generating blocks...')
     blocks = cells_copy.dissolve(by=id_name)
-    cells_copy = cells_copy.drop([id_name], axis=1)
 
     print('Multipart to singlepart...')
     blocks = blocks.explode()
     blocks.reset_index(inplace=True, drop=True)
 
     blocks['geometry'] = blocks.exterior
+    blocks[id_name] = range(len(blocks))
 
-    uid = 1
     for idx, row in tqdm(blocks.iterrows(), total=blocks.shape[0]):
-        blocks.loc[idx, id_name] = uid
-        uid = uid + 1
         blocks.loc[idx, 'geometry'] = Polygon(row['geometry'])
 
     # if polygon is within another one, delete it
@@ -741,13 +573,12 @@ def blocks(tessellation, edges, buildings, id_name, unique_id):
     bl_ID_to_uID = centroids_w_bl_ID2[[unique_id, id_name]]
 
     print('Attribute join (buildings)...')
-    buildings = buildings.merge(bl_ID_to_uID, on=unique_id)
+    buildings_m = buildings[[unique_id]].merge(bl_ID_to_uID, on=unique_id, how='left')
 
     print('Attribute join (tesselation)...')
-    cells = tessellation.merge(bl_ID_to_uID, on=unique_id)
+    cells_m = tessellation[[unique_id]].merge(bl_ID_to_uID, on=unique_id, how='left')
 
-    print('Done')
-    return (buildings, cells, blocks_save)
+    return (blocks_save, buildings_m[id_name], cells_m[id_name])
 
 
 def get_network_id(left, right, unique_id, network_id, min_size=100):
@@ -792,7 +623,9 @@ def get_network_id(left, right, unique_id, network_id, min_size=100):
     # point with edges of size 2*MIN_SIZE, you know a priori that at least one
     # segment is intersected with the box. Otherwise, you could get an inexact
     # solution, there is an exception checking this, though.
-
+    left = left.copy()
+    right = right.copy()
+    
     if not isinstance(unique_id, str):
         left['mm_uid'] = unique_id
         unique_id = 'mm_uid'
@@ -838,8 +671,8 @@ def get_node_id(objects, nodes, edges, node_id, edge_id):
     """
     Snap each building to closest street network node on the closest network edge.
 
-    Adds node ID to objects (preferably buildings). Gets ID of edge, and determines
-    which of its end points is closer.
+    Adds node ID to objects (preferably buildings). Gets ID of edge (:py:func:`momepy.get_node_id`)
+    , and determines which of its end points is closer to building centroid.
 
     Parameters
     ----------
@@ -859,10 +692,8 @@ def get_node_id(objects, nodes, edges, node_id, edge_id):
     node_ids : Series
         Series containing node ID for objects
 
-    Examples
-    --------
-
     """
+    nodes = nodes.copy()
     if not isinstance(node_id, str):
         nodes['mm_noid'] = node_id
         node_id = 'mm_noid'
@@ -886,220 +717,5 @@ def get_node_id(objects, nodes, edges, node_id, edge_id):
             else:
                 results_list.append(startID)
 
-    series = pd.Series(results_list)
+    series = pd.Series(results_list, index=objects.index)
     return series
-
-
-# '''
-# street_edges():
-#
-# Generate street edges based on buildings, blocks, tesselation and street network with street names
-# Adds nID and eID to buildings and tesselation.
-#
-#     buildings = gdf of buildings (with unique id)
-#     streets = gdf of street network (with street names and unique network segment id)
-#     tesselation = gdf of tesselation (with unique id and block id)
-#     street_name_column = column with street names
-#     unique_id_column = column with unique ids
-#     block_id_column = column with block ids
-#     network_id_column = column with network ids
-#     tesselation_to = path to save tesselation with nID, eID
-#     buildings_to = path to save buildings with nID, eID
-#     save_to = path to save street edges
-#
-# Optional:
-# '''
-#
-#
-# def street_edges(buildings, streets, tesselation, street_name_column,
-#                  unique_id_column, block_id_column, network_id_column,
-#                  tesselation_to, buildings_to, save_to):
-#     INFTY = 1000000000000
-#     MIN_SIZE = 100
-#     # MIN_SIZE should be a vaule such that if you build a box centered in each
-#     # point with edges of size 2*MIN_SIZE, you know a priori that at least one
-#     # segment is intersected with the box. Otherwise, you could get an inexact
-#     # solution, there is an exception checking this, though.
-#
-#     def distance(a, b):
-#         return math.sqrt((a[0] - b[0]) ** 2 + (a[1] - b[1]) ** 2)
-#
-#     def get_distance(apoint, segment):
-#         a = apoint
-#         b, c = segment
-#         # t = <a-b, c-b>/|c-b|**2
-#         # because p(a) = t*(c-b)+b is the ortogonal projection of vector a
-#         # over the rectline that includes the points b and c.
-#         t = (a[0] - b[0]) * (c[0] - b[0]) + (a[1] - b[1]) * (c[1] - b[1])
-#         t = t / ((c[0] - b[0]) ** 2 + (c[1] - b[1]) ** 2)
-#         # Only if t 0 <= t <= 1 the projection is in the interior of
-#         # segment b-c, and it is the point that minimize the distance
-#         # (by pythagoras theorem).
-#         if 0 < t < 1:
-#             pcoords = (t * (c[0] - b[0]) + b[0], t * (c[1] - b[1]) + b[1])
-#             dmin = distance(a, pcoords)
-#             return pcoords, dmin
-#         elif t <= 0:
-#             return b, distance(a, b)
-#         elif 1 <= t:
-#             return c, distance(a, c)
-#
-#     def get_rtree(lines):
-#         def generate_items():
-#             sindx = 0
-#             for lid, nid, l in tqdm(lines, total=len(lines)):
-#                 for i in range(len(l) - 1):
-#                     a, b = l[i]
-#                     c, d = l[i + 1]
-#                     segment = ((a, b), (c, d))
-#                     box = (min(a, c), min(b, d), max(a, c), max(b, d))
-#                     # box = left, bottom, right, top
-#                     yield (sindx, box, (lid, segment, nid))
-#                     sindx += 1
-#         return index.Index(generate_items())
-#
-#     def get_solution(idx, points):
-#         result = {}
-#         for p in tqdm(points, total=len(points)):
-#             pbox = (p[0] - MIN_SIZE, p[1] - MIN_SIZE, p[0] + MIN_SIZE, p[1] + MIN_SIZE)
-#             hits = idx.intersection(pbox, objects='raw')
-#             d = INFTY
-#             s = None
-#             for h in hits:
-#                 nearest_p, new_d = get_distance(p, h[1])
-#                 if d >= new_d:
-#                     d = new_d
-#                     # s = (h[0], h[1], nearest_p, new_d)
-#                     s = (h[0], h[1], h[-1])
-#             result[p] = s
-#             if s is None:
-#                 result[p] = (0, 0)
-#
-#             # some checking you could remove after you adjust the constants
-#             # if s is None:
-#             #     raise Warning("It seems INFTY is not big enough. Point was not attached to street. It might be too far.", p)
-#
-#             # pboxpol = ((pbox[0], pbox[1]), (pbox[2], pbox[1]),
-#             #            (pbox[2], pbox[3]), (pbox[0], pbox[3]))
-#             # if not Polygon(pboxpol).intersects(LineString(s[1])):
-#             #     msg = "It seems MIN_SIZE is not big enough. "
-#             #     msg += "You could get inexact solutions if remove this exception."
-#             #     raise Exception(msg)
-#
-#         return result
-#
-#     print('Generating centroids...')
-#     buildings_c = buildings.copy()
-#     buildings_c['geometry'] = buildings_c.centroid  # make centroids
-#
-#     print('Generating list of points...')
-#     # make points list for input
-#     centroid_list = []
-#     for idx, row in tqdm(buildings_c.iterrows(), total=buildings_c.shape[0]):
-#         centroid_list = centroid_list + list(row['geometry'].coords)
-#
-#     print('Generating list of lines...')
-#     # make streets list for input
-#     street_list = []
-#     for idx, row in tqdm(streets.iterrows(), total=streets.shape[0]):
-#         street_list.append((row[street_name_column], row[network_id_column], list(row['geometry'].coords)))
-#     print('Generating rtree...')
-#     idx = get_rtree(street_list)
-#
-#     print('Snapping...')
-#     solutions = get_solution(idx, centroid_list)
-#
-#     print('Forming DataFrame...')
-#     df = pd.DataFrame.from_dict(solutions, orient='index', columns=['street', 'unused', network_id_column])  # solutions dict to df
-#     df['point'] = df.index  # point to column
-#     df = df.reset_index()
-#     df['idx'] = df.index
-#     buildings_c['idx'] = buildings_c.index
-#
-#     print('Joining DataFrames...')
-#     joined = buildings_c.merge(df, on='idx')
-#     print('Cleaning DataFrames...')
-#     cleaned = joined[[unique_id_column, 'street', network_id_column]]
-#
-#     print('Merging with tesselation...')
-#     tesselation = tesselation.merge(cleaned, on=unique_id_column)
-#
-#     print('Defining merge ID...')
-#     for idx, row in tqdm(tesselation.iterrows(), total=tesselation.shape[0]):
-#         tesselation.loc[idx, 'mergeID'] = str(row['street']) + str(row[block_id_column])
-#
-#     print('Dissolving...')
-#     edges = tesselation.dissolve(by='mergeID')
-#
-#     # multipart geometry to singlepart
-#     def multi2single(gpdf):
-#         gpdf_singlepoly = gpdf[gpdf.geometry.type == 'Polygon']
-#         gpdf_multipoly = gpdf[gpdf.geometry.type == 'MultiPolygon']
-#
-#         for i, row in gpdf_multipoly.iterrows():
-#             Series_geometries = pd.Series(row.geometry)
-#             df = pd.concat([gpd.GeoDataFrame(row, crs=gpdf_multipoly.crs).T] * len(Series_geometries), ignore_index=True)
-#             df['geometry'] = Series_geometries
-#             gpdf_singlepoly = pd.concat([gpdf_singlepoly, df])
-#
-#         gpdf_singlepoly.reset_index(inplace=True, drop=True)
-#         return gpdf_singlepoly
-#
-#     edges_single = multi2single(edges)
-#     edges_single['geometry'] = edges_single.exterior
-#     print('Generating unique edge ID...')
-#     id = 1
-#     for idx, row in tqdm(edges_single.iterrows(), total=edges_single.shape[0]):
-#         edges_single.loc[idx, 'eID'] = id
-#         id = id + 1
-#         edges_single.loc[idx, 'geometry'] = Polygon(row['geometry'])
-#
-#     edges_clean = edges_single[['geometry', 'eID', block_id_column]]
-#
-#     print('Isolating islands...')
-#     sindex = edges_clean.sindex
-#     islands = []
-#     for idx, row in edges_clean.iterrows():
-#         possible_matches_index = list(sindex.intersection(row['geometry'].bounds))
-#         possible_matches = edges_clean.iloc[possible_matches_index]
-#         possible_matches = possible_matches.drop([idx], axis=0)
-#         if possible_matches.contains(row['geometry']).any():
-#             islands.append(idx)
-#
-#     edges_clean = edges_clean.drop(islands, axis=0)
-#     print(len(islands), 'islands deleted.')
-#     print('Cleaning edges...')
-#     edges_clean['geometry'] = edges_clean.buffer(0.000000001)
-#
-#     print('Saving street edges to', save_to)
-#     edges_clean.to_file(save_to)
-#
-#     print('Cleaning tesselation...')
-#     tesselation = tesselation.drop(['street', 'mergeID'], axis=1)
-#
-#     print('Tesselation spatial join [1/3]...')
-#     tess_centroid = tesselation.copy()
-#     tess_centroid['geometry'] = tess_centroid.centroid
-#
-#     edg_join = edges_clean.drop(['bID'], axis=1)
-#
-#     print('Tesselation spatial join [2/3]...')
-#     tess_with_eID = gpd.sjoin(tess_centroid, edg_join, how='left', op='intersects')
-#     tess_with_eID = tess_with_eID[['uID', 'eID']]
-#
-#     print('Tesselation spatial join [3/3]...')
-#     tesselation = tesselation.merge(tess_with_eID, on='uID')
-#
-#     print('Saving tesselation to', tesselation_to)
-#     tesselation.to_file(tesselation_to)
-#
-#     print('Buildings attribute join...')
-#     # attribute join cell -> building
-#     tess_nid_eid = tesselation[['uID', 'eID', 'nID']]
-#
-#     buildings = buildings.merge(tess_nid_eid, on='uID')
-#
-#     print('Saving buildings to', buildings_to)
-#     buildings.to_file(buildings_to)
-#
-#     print('Done.')
diff --git a/momepy/graph.py b/momepy/graph.py
index 9dd639d..6782704 100644
--- a/momepy/graph.py
+++ b/momepy/graph.py
@@ -22,7 +22,7 @@ def node_degree(graph, name='degree'):
     graph : networkx.Graph
         Graph representing street network.
         Ideally genereated from GeoDataFrame using :py:func:`momepy.gdf_to_nx`
-    name : str, optional
+    name : str(default 'degree')
         calculated attribute name
 
     Returns
@@ -30,12 +30,9 @@ def node_degree(graph, name='degree'):
     Graph
         networkx.Graph
 
-    References
-    ----------
-
     Examples
     --------
-
+    >>> network_graph = mm.node_degree(network_graph)
     """
     netx = graph
 
@@ -77,11 +74,12 @@ def meshedness(graph, radius=5, name='meshedness', distance=None):
 
     References
     ----------
-    Ale, Song et al 2013
+    Feliciotti A (2018) RESILIENCE AND URBAN DESIGN:A SYSTEMS APPROACH TO THE STUDY OF RESILIENCE IN URBAN FORM.
+    LEARNING FROM THE CASE OF GORBALS. Glasgow.
 
     Examples
     --------
-
+    >>> network_graph = mm.meshedness(network_graph, radius=800, distance='edge_length')
     """
     netx = graph
 
@@ -116,11 +114,9 @@ def mean_node_dist(graph, name='meanlen', length='mm_len'):
     Graph
         networkx.Graph
 
-    References
-    ----------
-
     Examples
     --------
+    >>> network_graph = mm.mean_node_dist(network_graph)
 
     """
     netx = graph
@@ -152,6 +148,8 @@ def cds_length(graph, radius=5, mode='sum', name='cds_len', degree='degree', len
         Ideally genereated from GeoDataFrame using :py:func:`momepy.gdf_to_nx`
     radius : int
         number of topological steps defining the extent of subgraph
+    mode : str (defualt 'sum')
+        if 'sum', calculate total length, if 'mean' calculate mean length
     name : str, optional
         calculated attribute name
     degree : str
@@ -169,12 +167,9 @@ def cds_length(graph, radius=5, mode='sum', name='cds_len', degree='degree', len
     Graph
         networkx.Graph
 
-    References
-    ----------
-
     Examples
     --------
-
+    >>> network_graph = mm.cds_length(network_graph, radius=9, mode='mean')
     """
     # node degree needed beforehand
     netx = graph
@@ -201,7 +196,7 @@ def cds_length(graph, radius=5, mode='sum', name='cds_len', degree='degree', len
 
 def mean_node_degree(graph, radius=5, name='mean_nd', degree='degree', distance=None):
     """
-    Calculates meshedness for subgraph around each node.
+    Calculates mean node degree for subgraph around each node.
 
     Subgraph is generated around each node within set radius. If distance=None,
     radius will define topological distance, otherwise it uses values in distance
@@ -231,12 +226,9 @@ def mean_node_degree(graph, radius=5, name='mean_nd', degree='degree', distance=
     Graph
         networkx.Graph
 
-    References
-    ----------
-
     Examples
     --------
-
+    >>> network_graph = mm.mean_node_degree(network_graph, radius=3)
     """
     netx = graph
 
@@ -283,12 +275,9 @@ def proportion(graph, radius=5, three=None, four=None, dead=None, degree='degree
     Graph
         networkx.Graph
 
-    References
-    ----------
-
     Examples
     --------
-
+    >>> network_graph = mm.proportion(network_graph, three='threeway', four='fourway', dead='deadends')
     """
     if not three and not four and not dead:
         raise ValueError('Nothing to calculate. Define names for at least one proportion to be calculated: three, four, dead.')
@@ -340,11 +329,12 @@ def cyclomatic(graph, radius=5, name='cyclomatic', distance=None):
 
     References
     ----------
-    Bourdic, Salat, Nowacki 2012 Assessing cities
+    Bourdic L, Salat S and Nowacki C (2012) Assessing cities: a new system of cross-scale spatial indicators.
+    Building Research & Information 40(5): 592–605.
 
     Examples
     --------
-
+    >>> network_graph = mm.cyclomatic(network_graph, radius=3)
     """
     netx = graph
 
@@ -389,11 +379,12 @@ def edge_node_ratio(graph, radius=5, name='edge_node_ratio', distance=None):
 
     References
     ----------
-    Dibble
+    Dibble J, Prelorendjos A, Romice O, et al. (2017) On the origin of spaces: Morphometric foundations of urban form evolution.
+    Environment and Planning B: Urban Analytics and City Science 46(4): 707–730.
 
     Examples
     --------
-
+    >>> network_graph = mm.edge_node_ratio(network_graph, radius=3)
     """
     netx = graph
 
@@ -438,10 +429,12 @@ def gamma(graph, radius=5, name='gamma', distance=None):
 
     References
     ----------
-    Dibble
+    Dibble J, Prelorendjos A, Romice O, et al. (2017) On the origin of spaces: Morphometric foundations of urban form evolution.
+    Environment and Planning B: Urban Analytics and City Science 46(4): 707–730.
 
     Examples
     --------
+    >>> network_graph = mm.gamma(network_graph, radius=3)
 
     """
     netx = graph
@@ -461,6 +454,7 @@ def gamma(graph, radius=5, name='gamma', distance=None):
 def _closeness_centrality(G, u=None, distance=None, wf_improved=True, len_graph=None):
     r"""Compute closeness centrality for nodes. Slight adaptation of networkx
     `closeness_centrality` to allow normalisation for local closeness.
+    Adapted script used in networkx.
 
     Closeness centrality [1]_ of a node `u` is the reciprocal of the
     average shortest path distance to `u` over all `n-1` reachable nodes.
@@ -508,27 +502,6 @@ def _closeness_centrality(G, u=None, distance=None, wf_improved=True, len_graph=
     nodes : dictionary
       Dictionary of nodes with closeness centrality as the value.
 
-    See Also
-    --------
-    betweenness_centrality, load_centrality, eigenvector_centrality,
-    degree_centrality
-
-    Notes
-    -----
-    The closeness centrality is normalized to `(n-1)/(|G|-1)` where
-    `n` is the number of nodes in the connected part of graph
-    containing the node.  If the graph is not completely connected,
-    this algorithm computes the closeness centrality for each
-    connected part separately scaled by that parts size.
-
-    If the 'distance' keyword is set to an edge attribute key then the
-    shortest-path length will be computed using Dijkstra's algorithm with
-    that edge attribute as the edge weight.
-
-    In NetworkX 2.2 and earlier a bug caused Dijkstra's algorithm to use the
-    outward distance rather than the inward distance. If you use a 'distance'
-    keyword and a DiGraph, your results will change between v2.2 and v2.3.
-
     References
     ----------
     .. [1] Linton C. Freeman: Centrality in networks: I.
@@ -598,10 +571,12 @@ def local_closeness(graph, radius=5, name='closeness', distance=None, closeness_
 
     References
     ----------
-    Porta
+    Porta S, Crucitti P and Latora V (2006) The network analysis of urban streets: A primal approach.
+    Environment and Planning B: Planning and Design 33(5): 705–725.
 
     Examples
     --------
+    >>> network_graph = mm.local_closeness(network_graph, radius=400, distance='edge_length')
 
     """
     netx = graph
@@ -615,7 +590,7 @@ def local_closeness(graph, radius=5, name='closeness', distance=None, closeness_
 
 def eigenvector(graph, name='eigen', **kwargs):
     """
-    Calculates mean distance to neighbouring nodes.
+    Calculates eigenvector centrality of network.
 
     .. math::
 
@@ -635,11 +610,9 @@ def eigenvector(graph, name='eigen', **kwargs):
     Graph
         networkx.Graph
 
-    References
-    ----------
-
     Examples
     --------
+    >>> network_graph = mm.eigenvector(network_graph)
 
     """
     netx = graph
@@ -670,12 +643,9 @@ def clustering(graph, name='cluster'):
     Graph
         networkx.Graph
 
-    References
-    ----------
-
     Examples
     --------
-
+    >>> network_graph = mm.clustering(network_graph)
     """
     netx = graph
 
diff --git a/momepy/intensity.py b/momepy/intensity.py
index 72b335f..4d28f66 100644
--- a/momepy/intensity.py
+++ b/momepy/intensity.py
@@ -10,9 +10,9 @@ import numpy as np
 import collections
 
 
-def covered_area_ratio(left, right, left_areas, right_areas, unique_id):
+def object_area_ratio(left, right, left_areas, right_areas, unique_id):
     """
-    Calculate covered area ratio of objects.
+    Calculate covered area ratio or floor area ratio of objects.
 
     .. math::
         \\textit{covering object area} \\over \\textit{covered object area}
@@ -23,13 +23,14 @@ def covered_area_ratio(left, right, left_areas, right_areas, unique_id):
         GeoDataFrame containing objects being covered (e.g. land unit)
     right : GeoDataFrame
         GeoDataFrame with covering objects (e.g. building)
-    left_areas : str
-        name of the column of left gdf where is stored area value
-    right_areas : str
-        name of the column of right gdf where is stored area value
+    left_areas : str, list, np.array, pd.Series
+        the name of the left dataframe column, np.array, or pd.Series where is stored area value.
+    right_areas : str, list, np.array, pd.Series
+        the name of the right dataframe column, np.array, or pd.Series where is stored area value
+        representing either projected or floor area.
     unique_id : str
         name of the column with unique id shared amongst left and right gdfs.
-        If there is none, it could be generated by unique_id().
+        If there is none, it could be generated by :py:func:'momepy.unique_id()'.
 
     Returns
     -------
@@ -38,79 +39,36 @@ def covered_area_ratio(left, right, left_areas, right_areas, unique_id):
 
     References
     ---------
+    Schirmer PM and Axhausen KW (2015) A multiscale classiﬁcation of urban morphology.
+    Journal of Transport and Land Use 9(1): 101–130.
 
+    Examples
+    --------
+    >>> tessellation_df['CAR'] = mm.object_area_ratio(tessellation_df, buildings_df, 'area', 'area', 'uID')
+    Calculating object area ratio...
+    Merging DataFrames...
+    Calculating OAR...
+    Object area ratio calculated.
     """
-    print('Calculating covered area ratio...')
-
-    print('Merging DataFrames...')
-    look_for = right[[unique_id, right_areas]]  # keeping only necessary columns
-    look_for.rename(index=str, columns={right_areas: 'lf_area'}, inplace=True)
-    objects_merged = left.merge(look_for, on=unique_id)  # merging dataframes together
-
-    print('Calculating CAR...')
-
-    # define empty list for results
-    results_list = []
-
-    # fill new column with the value of area, iterating over rows one by one
-    for index, row in tqdm(objects_merged.iterrows(), total=objects_merged.shape[0]):
-        results_list.append(row['lf_area'] / row[left_areas])
-
-    series = pd.Series(results_list)
-
-    print('Covered area ratio calculated.')
-    return series
-
-
-def floor_area_ratio(left, right, left_areas, right_areas, unique_id):
-    """
-    Calculate floor area ratio of objects.
-
-    .. math::
-        \\textit{covering object floor area} \\over \\textit{covered object area}
-
-    Parameters
-    ----------
-    left : GeoDataFrame
-        GeoDataFrame containing objects being covered (e.g. land unit)
-    right : GeoDataFrame
-        GeoDataFrame with covering objects (e.g. building)
-    left_areas : str
-        name of the column of left gdf where is stored area value
-    right_areas : str
-        name of the column of right gdf where is stored floor area value
-    unique_id : str
-        name of the column with unique id shared amongst left and right gdfs.
-        If there is none, it could be generated by unique_id().
-
-    Returns
-    -------
-    Series
-        Series containing resulting values.
-
-    References
-    ---------
 
-    """
-    print('Calculating floor area ratio...')
+    print('Calculating object area ratio...')
+    if not isinstance(left_areas, str):
+        left['mm_a'] = left_areas
+        left_areas = 'mm_a'
+    if not isinstance(right_areas, str):
+        right['mm_a'] = right_areas
+        right_areas = 'mm_a'
 
     print('Merging DataFrames...')
     look_for = right[[unique_id, right_areas]]  # keeping only necessary columns
     look_for.rename(index=str, columns={right_areas: 'lf_area'}, inplace=True)
-    objects_merged = left.merge(look_for, on=unique_id)  # merging dataframes together
-
-    print('Calculating FAR...')
+    objects_merged = left[[unique_id, left_areas]].merge(look_for, on=unique_id)  # merging dataframes together
 
-    # define empty list for results
-    results_list = []
+    print('Calculating OAR...')
 
-    # fill new column with the value of area, iterating over rows one by one
-    for index, row in tqdm(objects_merged.iterrows(), total=objects_merged.shape[0]):
-        results_list.append(row['lf_area'] / row[left_areas])
+    series = objects_merged['lf_area'] / objects_merged[left_areas]
 
-    series = pd.Series(results_list)
-
-    print('Floor area ratio calculated.')
+    print('Object area ratio calculated.')
     return series
 
 
@@ -133,9 +91,9 @@ def elements_count(left, right, left_id, right_id, weighted=False):
     right : GeoDataFrame
         GeoDataFrame containing objects to analyse
     left_id : str
-        name of the column where is stored unique ID in aggr
+        name of the column where is stored unique ID in left gdf
     right_id : str
-        name of the column where is stored unique ID of aggregation in elements gdf
+        name of the column where is stored unique ID of aggregation in right gdf
     weighted : bool (default False)
         if weighted=True, count will be divided by the area or length
 
@@ -145,16 +103,20 @@ def elements_count(left, right, left_id, right_id, weighted=False):
         Series containing resulting values.
 
     References
-    ---------
-    Hermosilla T, Ruiz LA, Recio JA, et al. (2012) Assessing contextual descriptive features
+    ----------
+    1. Hermosilla T, Ruiz LA, Recio JA, et al. (2012) Assessing contextual descriptive features
     for plot-based classification of urban areas. Landscape and Urban Planning, Elsevier B.V.
     106(1): 124–137.
-    Feliciotti A (2018) RESILIENCE AND URBAN DESIGN:A SYSTEMS APPROACH TO THE
+    2. Feliciotti A (2018) RESILIENCE AND URBAN DESIGN:A SYSTEMS APPROACH TO THE
     STUDY OF RESILIENCE IN URBAN FORM. LEARNING FROM THE CASE OF GORBALS. Glasgow.
+
+    Examples
+    --------
+    >>> blocks_df['buildings_count'] = mm.elements_count(blocks_df, buildings_df, 'bID', 'bID', weighted=True)
     """
     count = collections.Counter(right[right_id])
     df = pd.DataFrame.from_dict(count, orient='index', columns=['mm_count'])
-    joined = left.join(df['mm_count'], on=left_id)
+    joined = left[[left_id, 'geometry']].join(df['mm_count'], on=left_id)
     joined['mm_count'][np.isnan(joined['mm_count'])] = 0
 
     if weighted:
@@ -176,8 +138,8 @@ def courtyards(gdf, block_id, spatial_weights=None):
     ----------
     gdf : GeoDataFrame
         GeoDataFrame containing objects to analyse
-    block_id : str
-        name of the column where is stored block ID
+    block_id : str, list, np.array, pd.Series
+        the name of the dataframe column, np.array, or pd.Series where is stored block ID
     spatial_weights : libpysal.weights, optional
         spatial weights matrix - If None, Queen contiguity matrix will be calculated
         based on objects. It is to denote adjacent buildings (note: based on index).
@@ -187,14 +149,30 @@ def courtyards(gdf, block_id, spatial_weights=None):
     Series
         Series containing resulting values.
 
+    References
+    ---------
+    Schirmer PM and Axhausen KW (2015) A multiscale classiﬁcation of urban morphology.
+    Journal of Transport and Land Use 9(1): 101–130.
+
+    Examples
+    --------
+    >>> buildings_df['courtyards'] = mm.courtyards(buildings_df, 'bID')
+    Calculating courtyards...
+    Calculating spatial weights...
+    Spatial weights ready...
+    Courtyards calculated.
+
     Notes
     -----
-    Script is not optimised at all, so it is currently extremely slow.
+    Script is not well optimised, so it is currently slow.
     """
     # define empty list for results
     results_list = []
 
     print('Calculating courtyards...')
+    if not isinstance(block_id, str):
+        gdf['mm_bid'] = block_id
+        block_id = 'mm_bid'
 
     # if weights matrix is not passed, generate it from objects
     if spatial_weights is None:
@@ -228,14 +206,14 @@ def courtyards(gdf, block_id, spatial_weights=None):
             try:
                 interiors = len(list(dissolved.interiors))
             except(ValueError):
-                print('Something happened.')
+                print('Something unexpected happened.')
             for b in to_join:
                 courtyards[b] = interiors  # fill dict with values
     # copy values from dict to gdf
     for index, row in tqdm(gdf.iterrows(), total=gdf.shape[0]):
         results_list.append(courtyards[index])
 
-    series = pd.Series(results_list)
+    series = pd.Series(results_list, index=gdf.index)
     print('Courtyards calculated.')
     return series
 
@@ -260,7 +238,6 @@ def blocks_count(gdf, block_id, spatial_weights, unique_id):
     unique_id : str
         name of the column with unique id used as spatial_weights index.
 
-
     Returns
     -------
     Series
@@ -268,15 +245,15 @@ def blocks_count(gdf, block_id, spatial_weights, unique_id):
 
     References
     ----------
-    Jacob
+    Dibble J, Prelorendjos A, Romice O, et al. (2017) On the origin of spaces: Morphometric foundations of urban form evolution.
+    Environment and Planning B: Urban Analytics and City Science 46(4): 707–730.
 
     Examples
     --------
-
-    Notes
-    -----
-    Blocks count or blocks density?
-
+    >>> sw4 = mm.Queen_higher(k=4, geodataframe='tessellation_df', ids='uID')
+    >>> tessellation_df['blocks_within_4'] = mm.blocks(tessellation_df, 'bID', sw4, 'uID')
+    Calculating blocks...
+    Blocks calculated.
     """
     # define empty list for results
     results_list = []
@@ -297,7 +274,7 @@ def blocks_count(gdf, block_id, spatial_weights, unique_id):
 
         results_list.append(len(set(list(vicinity[block_id]))) / sum(vicinity.geometry.area))
 
-    series = pd.Series(results_list)
+    series = pd.Series(results_list, index=gdf.index)
 
     print('Blocks calculated.')
     return series
@@ -340,11 +317,11 @@ def reached(left, right, unique_id, spatial_weights=None, mode='count', values=N
     Series
         Series containing resulting values.
 
-    References
-    ----------
-
     Examples
     --------
+    >>> streets_df['reached_buildings'] = mm.reached(streets_df, buildings_df, 'uID')
+    Calculating reached count...
+    Reached count calculated.
 
     """
     # define empty list for results
@@ -389,7 +366,7 @@ def reached(left, right, unique_id, spatial_weights=None, mode='count', values=N
             else:
                 results_list.append(np.nanstd(right.loc[right[unique_id].isin(ids)].geometry.area))
 
-    series = pd.Series(results_list)
+    series = pd.Series(results_list, index=left.index)
 
     print('Reached {} calculated.'.format(mode))
     return series
@@ -417,9 +394,9 @@ def node_density(left, right, spatial_weights, weighted=False, node_degree=None,
         if True density will take into account node degree as k-1
     node_degree : str
         name of the column of left gdf containing node degree
-    node_start : str
+    node_start : str (default 'node_start')
         name of the column of right gdf containing id of starting node
-    node_end : str
+    node_end : str (default 'node_end')
         name of the column of right gdf containing id of ending node
 
     Returns
@@ -429,10 +406,14 @@ def node_density(left, right, spatial_weights, weighted=False, node_degree=None,
 
     References
     ---------
-    Jacob
+    Dibble J, Prelorendjos A, Romice O, et al. (2017) On the origin of spaces: Morphometric foundations of urban form evolution.
+    Environment and Planning B: Urban Analytics and City Science 46(4): 707–730.
 
-    Notes
-    -----
+    Examples
+    --------
+    >>> nodes['density'] = mm.node_density(nodes, edges, sw)
+    Calculating node density...
+    Node density calculated.
 
     """
     # define empty list for results
@@ -459,14 +440,14 @@ def node_density(left, right, spatial_weights, weighted=False, node_degree=None,
         else:
             results_list.append(0)
 
-    series = pd.Series(results_list)
+    series = pd.Series(results_list, index=left.index)
     print('Node density calculated.')
     return series
 
 
 def density(gdf, values, spatial_weights, unique_id, areas=None):
     """
-    Calculate the density
+    Calculate the gross density
 
     .. math::
 
@@ -492,7 +473,14 @@ def density(gdf, values, spatial_weights, unique_id, areas=None):
 
     References
     ---------
-    Jacob??
+    Dibble J, Prelorendjos A, Romice O, et al. (2017) On the origin of spaces: Morphometric foundations of urban form evolution.
+    Environment and Planning B: Urban Analytics and City Science 46(4): 707–730.
+
+    Examples
+    --------
+    >>> tessellation_df['floor_area_dens'] = mm.density(tessellation_df, 'floor_area', sw, 'uID')
+    Calculating gross density...
+    Gross density calculated.
     """
     # define empty list for results
     results_list = []
@@ -523,6 +511,6 @@ def density(gdf, values, spatial_weights, unique_id, areas=None):
 
         results_list.append(sum(values_list) / sum(areas_list))
 
-    series = pd.Series(results_list)
+    series = pd.Series(results_list, index=gdf.index)
     print('Gross density calculated.')
     return series
diff --git a/momepy/shape.py b/momepy/shape.py
index b9b2034..cb247a4 100644
--- a/momepy/shape.py
+++ b/momepy/shape.py
@@ -25,7 +25,7 @@ def form_factor(gdf, volumes, areas=None):
         GeoDataFrame containing objects
     volumes : str, list, np.array, pd.Series
         the name of the dataframe column, np.array, or pd.Series where is stored volume value.
-        (To calculate volume you can use :py:func:`momepy.dimension.volume`)
+        (To calculate volume you can use :py:func:`momepy.volume`)
     areas : str, list, np.array, pd.Series (default None)
         the name of the dataframe column, np.array, or pd.Series where is stored area value. If set to None, function will calculate areas
         during the process without saving them separately.
@@ -751,7 +751,7 @@ def corners(gdf):
 
         results_list.append(corners)
 
-    series = pd.Series(results_list)
+    series = pd.Series(results_list, index=gdf.index)
 
     print('Corners calculated.')
     return series
@@ -842,7 +842,7 @@ def squareness(gdf):
             deviations.append(dev)
         results_list.append(np.mean(deviations))
 
-    series = pd.Series(results_list)
+    series = pd.Series(results_list, index=gdf.index)
 
     print('Squareness calculated.')
     return series
@@ -912,7 +912,7 @@ def equivalent_rectangular_index(gdf, areas=None, perimeters=None):
         bbox = row['geometry'].minimum_rotated_rectangle
         results_list.append(math.sqrt(row[areas] / bbox.area) * (bbox.length / row[perimeters]))
 
-    series = pd.Series(results_list)
+    series = pd.Series(results_list, index=gdf.index)
 
     print('Equivalent rectangular index calculated.')
     return series
@@ -977,7 +977,7 @@ def elongation(gdf):
 
         results_list.append(elo)
 
-    series = pd.Series(results_list)
+    series = pd.Series(results_list, index=gdf.index)
 
     print('Elongation calculated.')
     return series
@@ -1081,8 +1081,8 @@ def centroid_corners(gdf):
         else:
             results_list.append(np.mean(distances))  # calculate mean
             results_list_sd.append(np.std(distances))  # calculate st.dev
-    means = pd.Series(results_list)
-    st_devs = pd.Series(results_list_sd)
+    means = pd.Series(results_list, index=gdf.index)
+    st_devs = pd.Series(results_list_sd, index=gdf.index)
 
     print('Distances centroid - corner calculated.')
     return means, st_devs
@@ -1129,7 +1129,7 @@ def linearity(gdf):
         euclidean = Point(row['geometry'].coords[0]).distance(Point(row['geometry'].coords[-1]))
         results_list.append(euclidean / row['geometry'].length)
 
-    series = pd.Series(results_list)
+    series = pd.Series(results_list, index=gdf.index)
 
     print('Linearity calculated.')
     return series
@@ -1163,13 +1163,11 @@ def compactness_weighted_axis(gdf, areas=None, perimeters=None, longest_axis=Non
     Series
         Series containing resulting values.
 
-    References
-    ---------
-    Resilience project
-
     Examples
     --------
-
+    >>> blocks_df['cwa'] = mm.compactness_weighted_axis(blocks_df)
+    Calculating compactness-weighted axis...
+    Compactness-weighted axis calculated.
     """
 
     print('Calculating compactness-weighted axis...')
diff --git a/momepy/utils.py b/momepy/utils.py
index 22e9a3e..98f3b75 100644
--- a/momepy/utils.py
+++ b/momepy/utils.py
@@ -3,12 +3,13 @@
 
 import geopandas as gpd
 import libpysal
-from shapely.geometry import Point
+from shapely.geometry import Point, LineString
 import networkx as nx
 import numpy as np
 from tqdm import tqdm
 import operator
 import shapely
+import math
 
 from .shape import circular_compactness
 
@@ -406,3 +407,222 @@ def network_false_nodes(gdf):
     geoms_gdf.crs = streets.crs
     streets = geoms_gdf
     return streets
+
+
+def snap_street_network_edge(edges, buildings, tolerance_street, tessellation=None, tolerance_edge=None):
+    """
+    Fix street network before performing blocks()
+
+    Extends unjoined ends of street segments to join with other segmets or tessellation boundary.
+
+    Parameters
+    ----------
+    edges : GeoDataFrame
+        GeoDataFrame containing street network
+    buildings : GeoDataFrame
+        GeoDataFrame containing building footprints
+    tolerance_street : float
+        tolerance in snapping to street network (by how much could be street segment extended).
+    tessellation : GeoDataFrame (default None)
+        GeoDataFrame containing morphological tessellation
+    tolerance_edge : float (default None)
+        tolerance in snapping to edge of tessellated area (by how much could be street segment extended).
+
+    Returns
+    -------
+    GeoDataFrame
+        GeoDataFrame of extended street network.
+
+    """
+    # extrapolating function - makes line as a extrapolation of existing with set length (tolerance)
+    def getExtrapoledLine(p1, p2, tolerance):
+        """
+        Creates a line extrapoled in p1->p2 direction.
+        """
+        EXTRAPOL_RATIO = tolerance  # length of a line
+        a = p2
+
+        # defining new point based on the vector between existing points
+        if p1[0] >= p2[0] and p1[1] >= p2[1]:
+            b = (p2[0] - EXTRAPOL_RATIO * math.cos(math.atan(math.fabs(p1[1] - p2[1] + 0.000001) / math.fabs(p1[0] - p2[0] + 0.000001))),
+                 p2[1] - EXTRAPOL_RATIO * math.sin(math.atan(math.fabs(p1[1] - p2[1] + 0.000001) / math.fabs(p1[0] - p2[0] + 0.000001))))
+        elif p1[0] <= p2[0] and p1[1] >= p2[1]:
+            b = (p2[0] + EXTRAPOL_RATIO * math.cos(math.atan(math.fabs(p1[1] - p2[1] + 0.000001) / math.fabs(p1[0] - p2[0] + 0.000001))),
+                 p2[1] - EXTRAPOL_RATIO * math.sin(math.atan(math.fabs(p1[1] - p2[1] + 0.000001) / math.fabs(p1[0] - p2[0] + 0.000001))))
+        elif p1[0] <= p2[0] and p1[1] <= p2[1]:
+            b = (p2[0] + EXTRAPOL_RATIO * math.cos(math.atan(math.fabs(p1[1] - p2[1] + 0.000001) / math.fabs(p1[0] - p2[0] + 0.000001))),
+                 p2[1] + EXTRAPOL_RATIO * math.sin(math.atan(math.fabs(p1[1] - p2[1] + 0.000001) / math.fabs(p1[0] - p2[0] + 0.000001))))
+        else:
+            b = (p2[0] - EXTRAPOL_RATIO * math.cos(math.atan(math.fabs(p1[1] - p2[1] + 0.000001) / math.fabs(p1[0] - p2[0] + 0.000001))),
+                 p2[1] + EXTRAPOL_RATIO * math.sin(math.atan(math.fabs(p1[1] - p2[1] + 0.000001) / math.fabs(p1[0] - p2[0] + 0.000001))))
+        return LineString([a, b])
+
+    # function extending line to closest object within set distance
+    def extend_line(tolerance, idx):
+        """
+        Extends a line geometry withing GeoDataFrame to snap on itself withing tolerance.
+        """
+        if Point(l_coords[-2]).distance(Point(l_coords[-1])) <= 0.001:
+            if len(l_coords) > 2:
+                extra = l_coords[-3:-1]
+            else:
+                return False
+        else:
+            extra = l_coords[-2:]
+        extrapolation = getExtrapoledLine(*extra, tolerance=tolerance)  # we use the last two points
+
+        possible_intersections_index = list(sindex.intersection(extrapolation.bounds))
+        possible_intersections_lines = network.iloc[possible_intersections_index]
+        possible_intersections_clean = possible_intersections_lines.drop(idx, axis=0)
+        possible_intersections = possible_intersections_clean.intersection(extrapolation)
+
+        if possible_intersections.any():
+
+            true_int = []
+            for one in list(possible_intersections.index):
+                if possible_intersections[one].type == 'Point':
+                    true_int.append(possible_intersections[one])
+                elif possible_intersections[one].type == 'MultiPoint':
+                    true_int.append(possible_intersections[one][0])
+                    true_int.append(possible_intersections[one][1])
+
+            if len(true_int) >= 1:
+                if len(true_int) > 1:
+                    distances = {}
+                    ix = 0
+                    for p in true_int:
+                        distance = p.distance(Point(l_coords[-1]))
+                        distances[ix] = distance
+                        ix = ix + 1
+                    minimal = min(distances.items(), key=operator.itemgetter(1))[0]
+                    new_point_coords = true_int[minimal].coords[0]
+                else:
+                    new_point_coords = true_int[0].coords[0]
+
+                l_coords.append(new_point_coords)
+                new_extended_line = LineString(l_coords)
+
+                # check whether the line goes through buildings. if so, ignore it
+                possible_buildings_index = list(bindex.intersection(new_extended_line.bounds))
+                possible_buildings = buildings.iloc[possible_buildings_index]
+                possible_intersections = possible_buildings.intersection(new_extended_line)
+
+                if possible_intersections.any():
+                    pass
+                else:
+                    network.loc[idx, 'geometry'] = new_extended_line
+        else:
+            return False
+
+    # function extending line to closest object within set distance to edge defined by tessellation
+    def extend_line_edge(tolerance, idx):
+        """
+        Extends a line geometry withing GeoDataFrame to snap on the boundary of tessellation withing tolerance.
+        """
+        if Point(l_coords[-2]).distance(Point(l_coords[-1])) <= 0.001:
+            if len(l_coords) > 2:
+                extra = l_coords[-3:-1]
+            else:
+                return False
+        else:
+            extra = l_coords[-2:]
+        extrapolation = getExtrapoledLine(*extra, tolerance)  # we use the last two points
+
+        # possible_intersections_index = list(qindex.intersection(extrapolation.bounds))
+        # possible_intersections_lines = geometry_cut.iloc[possible_intersections_index]
+        possible_intersections = geometry.intersection(extrapolation)
+
+        if possible_intersections.type != 'GeometryCollection':
+
+            true_int = []
+
+            if possible_intersections.type == 'Point':
+                true_int.append(possible_intersections)
+            elif possible_intersections.type == 'MultiPoint':
+                true_int.append(possible_intersections[0])
+                true_int.append(possible_intersections[1])
+
+            if len(true_int) >= 1:
+                if len(true_int) > 1:
+                    distances = {}
+                    ix = 0
+                    for p in true_int:
+                        distance = p.distance(Point(l_coords[-1]))
+                        distances[ix] = distance
+                        ix = ix + 1
+                    minimal = min(distances.items(), key=operator.itemgetter(1))[0]
+                    new_point_coords = true_int[minimal].coords[0]
+                else:
+                    new_point_coords = true_int[0].coords[0]
+
+                l_coords.append(new_point_coords)
+                new_extended_line = LineString(l_coords)
+
+                # check whether the line goes through buildings. if so, ignore it
+                possible_buildings_index = list(bindex.intersection(new_extended_line.bounds))
+                possible_buildings = buildings.iloc[possible_buildings_index]
+                possible_intersections = possible_buildings.intersection(new_extended_line)
+
+                if possible_intersections.any():
+                    pass
+                else:
+                    network.loc[idx, 'geometry'] = new_extended_line
+
+    network = edges.copy()
+    # generating spatial index (rtree)
+    print('Building R-tree for network...')
+    sindex = network.sindex
+    print('Building R-tree for buildings...')
+    bindex = buildings.sindex
+    if tessellation is not None:
+        print('Dissolving tesselation...')
+        geometry = tessellation.geometry.unary_union.boundary
+
+    print('Snapping...')
+    # iterating over each street segment
+    for idx, row in tqdm(network.iterrows(), total=network.shape[0]):
+
+        line = row['geometry']
+        l_coords = list(line.coords)
+        # network_w = network.drop(idx, axis=0)['geometry']  # ensure that it wont intersect itself
+        start = Point(l_coords[0])
+        end = Point(l_coords[-1])
+
+        # find out whether ends of the line are connected or not
+        possible_first_index = list(sindex.intersection(start.bounds))
+        possible_first_matches = network.iloc[possible_first_index]
+        possible_first_matches_clean = possible_first_matches.drop(idx, axis=0)
+        first = possible_first_matches_clean.intersects(start).any()
+
+        possible_second_index = list(sindex.intersection(end.bounds))
+        possible_second_matches = network.iloc[possible_second_index]
+        possible_second_matches_clean = possible_second_matches.drop(idx, axis=0)
+        second = possible_second_matches_clean.intersects(end).any()
+
+        # both ends connected, do nothing
+        if first and second:
+            continue
+        # start connected, extend  end
+        elif first and not second:
+            if extend_line(tolerance_street, idx) is False:
+                if tessellation is not None:
+                    extend_line_edge(tolerance_edge, idx)
+        # end connected, extend start
+        elif not first and second:
+            l_coords.reverse()
+            if extend_line(tolerance_street, idx) is False:
+                if tessellation is not None:
+                    extend_line_edge(tolerance_edge, idx)
+        # unconnected, extend both ends
+        elif not first and not second:
+            if extend_line(tolerance_street, idx) is False:
+                if tessellation is not None:
+                    extend_line_edge(tolerance_edge, idx)
+            l_coords.reverse()
+            if extend_line(tolerance_street, idx) is False:
+                if tessellation is not None:
+                    extend_line_edge(tolerance_edge, idx)
+        else:
+            print('Something went wrong.')
+
+    return network
diff --git a/requirements.txt b/requirements.txt
new file mode 100644
index 0000000..6ad21fc
--- /dev/null
+++ b/requirements.txt
@@ -0,0 +1,5 @@
+geopandas>=0.5.1
+networkx>=2.3
+libpysal>=4.0.1
+pandas>=0.24
+tqdm>=4.25
diff --git a/setup.py b/setup.py
index 95ff7e3..035dc62 100644
--- a/setup.py
+++ b/setup.py
@@ -15,13 +15,31 @@ for item in os.listdir("momepy/datasets"):
         elif item.endswith('.gpkg'):
             data_files.append(os.path.join("datasets", item))
 
+with open("README.md", "r") as fh:
+    long_description = fh.read()
+
 setup(name='momepy',
-      version='0.1alpha',
+      version='0.1rc1',
       description='Urban Morphology Measuring Toolkit',
+      long_description=long_description,
+      long_description_content_type="text/markdown",
       license="MIT",
       author='Martin Fleischmann',
       author_email='martin@martinfleischmann.net',
-      url='https://github.com/martinfleis/momepy',
+      keywords=["urban morphology",
+                "urban morphometrics",
+                "tessellation"],
+      url='http://momepy.org',
       packages=['momepy', 'momepy.datasets'],
       package_data={'momepy': data_files},
+      classifiers=["Programming Language :: Python :: 3",
+                   "License :: OSI Approved :: MIT License",
+                   "Operating System :: OS Independent",
+                   "Intended Audience :: Science/Research",
+                   "Topic :: Scientific/Engineering :: GIS"],
+      install_requires=["geopandas>=0.5.1",
+                        "networkx>=2.3",
+                        "libpysal>=4.0.1",
+                        "tqdm>=4.25.0"
+                        ]
       )


Test patch:
diff --git a/tests/test_distribution.py b/tests/test_distribution.py
index 04ecdfc..179fe32 100644
--- a/tests/test_distribution.py
+++ b/tests/test_distribution.py
@@ -43,7 +43,7 @@ class TestDistribution:
     def test_cell_alignment(self):
         self.df_buildings['orient'] = mm.orientation(self.df_buildings)
         self.df_tessellation['orient'] = mm.orientation(self.df_tessellation)
-        self.df_buildings['c_align'] = mm.cell_alignment(self.df_buildings, self.df_tessellation, 'orient', 'orient', 'uID')
+        self.df_buildings['c_align'] = mm.cell_alignment(self.df_buildings, self.df_tessellation, 'orient', 'orient', 'uID', 'uID')
         check = abs(self.df_buildings['orient'][0] - self.df_tessellation[
             self.df_tessellation['uID'] == self.df_buildings['uID'][0]]['orient'].iloc[0])
         assert self.df_buildings['c_align'][0] == check
diff --git a/tests/test_diversity.py b/tests/test_diversity.py
index df5d23c..8d41883 100644
--- a/tests/test_diversity.py
+++ b/tests/test_diversity.py
@@ -37,3 +37,12 @@ class TestDiversity:
         assert ht_sw[0] == 0.385
         quan_sw = mm.simpson(self.df_tessellation, 'area', self.sw, 'uID', binning='quantiles', k=3)
         assert quan_sw[0] == 0.395
+
+    def test_gini(self):
+        full_sw = mm.gini(self.df_tessellation, 'area', self.sw, 'uID')
+        assert full_sw[0] == 0.39453880038576933
+        limit = mm.gini(self.df_tessellation, 'area', self.sw, 'uID', rng=(10, 90))
+        assert limit[0] == 0.28532814171764387
+        self.df_tessellation['negative'] = self.df_tessellation.area - self.df_tessellation.area.mean()
+        neg = mm.gini(self.df_tessellation, 'negative', self.sw, 'uID')
+        assert neg[0] == 0.43135969210184627
diff --git a/tests/test_elements.py b/tests/test_elements.py
index 8e3487a..c10f03d 100644
--- a/tests/test_elements.py
+++ b/tests/test_elements.py
@@ -21,14 +21,10 @@ class TestElements:
         w = libpysal.weights.Queen.from_dataframe(queen_corners)
         assert w.neighbors[14] == [35, 36, 13, 15, 26, 27, 28, 30, 31]
 
-    def test_snap_street_network_edge(self):
-        snapped = mm.snap_street_network_edge(self.df_streets, self.df_buildings, self.df_tessellation, 20, 70)
-        assert sum(snapped.geometry.length) == 5980.041004739525
-
     def test_blocks(self):
-        buildings, cells, blocks = mm.blocks(self.df_tessellation, self.df_streets, self.df_buildings, 'bID', 'uID')
-        assert not buildings.bID.isna().any()
-        assert not cells.bID.isna().any()
+        blocks, buildingsID, cellsID = mm.blocks(self.df_tessellation, self.df_streets, self.df_buildings, 'bID', 'uID')
+        assert not buildingsID.isna().any()
+        assert not cellsID.isna().any()
         assert len(blocks) == 8
 
     def test_get_network_id(self):
@@ -49,8 +45,8 @@ class TestElements:
 
     def test__split_lines(self):
         large = mm.buffered_limit(self.df_buildings, 100)
-        dense = mm.elements._split_lines(large, 100, self.df_buildings.crs)
+        dense = mm.elements._split_lines(large, 100)
         small = mm.buffered_limit(self.df_buildings, 30)
-        dense2 = mm.elements._split_lines(small, 100, self.df_buildings.crs)
-        assert len(dense) == 273
-        assert len(dense2) == 667
+        dense2 = mm.elements._split_lines(small, 100)
+        assert len(dense) == 53
+        assert len(dense2) == 51
diff --git a/tests/test_intensity.py b/tests/test_intensity.py
index 09dce30..617784b 100644
--- a/tests/test_intensity.py
+++ b/tests/test_intensity.py
@@ -18,15 +18,14 @@ class TestIntensity:
         self.df_buildings['area'] = self.df_buildings.geometry.area
         self.df_buildings['fl_area'] = mm.floor_area(self.df_buildings, 'height')
         self.df_buildings['nID'] = mm.get_network_id(self.df_buildings, self.df_streets, 'uID', 'nID')
-        self.df_buildings, self.df_tessellation, self.blocks = mm.blocks(self.df_tessellation, self.df_streets, self.df_buildings, 'bID', 'uID')
+        self.blocks, self.df_buildings['bID'], self.df_tessellation['bID'] = mm.blocks(
+            self.df_tessellation, self.df_streets, self.df_buildings, 'bID', 'uID')
 
-    def test_covered_area_ratio(self):
-        car = mm.covered_area_ratio(self.df_tessellation, self.df_buildings, 'area', 'area', 'uID')
+    def test_object_area_ratio(self):
+        car = mm.object_area_ratio(self.df_tessellation, self.df_buildings, 'area', 'area', 'uID')
         check = 0.3206556897709747
         assert car.mean() == check
-
-    def test_floor_area_ratio(self):
-        far = mm.floor_area_ratio(self.df_tessellation, self.df_buildings, 'area', 'fl_area', 'uID')
+        far = mm.object_area_ratio(self.df_tessellation, self.df_buildings, 'area', 'fl_area', 'uID')
         check = 1.910949846262234
         assert far.mean() == check
 
diff --git a/tests/test_utils.py b/tests/test_utils.py
index 6d3f943..a2f66f2 100644
--- a/tests/test_utils.py
+++ b/tests/test_utils.py
@@ -73,3 +73,9 @@ class TestUtils:
         self.false_network = gpd.read_file(test_file_path2, layer='network')
         fixed = mm.network_false_nodes(self.false_network)
         assert len(fixed) == 55
+
+    def test_snap_street_network_edge(self):
+        snapped = mm.snap_street_network_edge(self.df_streets, self.df_buildings, 20, self.df_tessellation, 70)
+        snapped_nonedge = mm.snap_street_network_edge(self.df_streets, self.df_buildings, 20)
+        assert sum(snapped.geometry.length) == 5980.041004739525
+        assert sum(snapped_nonedge.geometry.length) < 5980.041004739525