Create geoprocessing API django app

src/mmw/apps/geocode/tests.py

@@ -12,7 +12,7 @@ class GeocodeTestCase(TestCase):
     def assert_candidate_exists_for(self, address):
         try:
             response = requests.get(
-                'http://localhost:80/api/geocode/?search=%s' % address).json()
+                'http://localhost:80/mmw/geocode/?search=%s' % address).json()
         except requests.RequestException:
             response = {}
 

src/mmw/apps/geoprocessing_api/calcs.py

@@ -0,0 +1,140 @@
+# -*- coding: utf-8 -*-
+from __future__ import print_function
+from __future__ import unicode_literals
+from __future__ import absolute_import
+
+from django.contrib.gis.geos import GEOSGeometry
+
+from django.conf import settings
+from django.db import connection
+
+from apps.modeling.mapshed.calcs import (animal_energy_units,
+                                         get_point_source_table)
+
+DRB = settings.DRB_PERIMETER
+
+ANIMAL_DISPLAY_NAMES = {
+    'sheep': 'Sheep',
+    'horses': 'Horses',
+    'turkeys': 'Turkeys',
+    'layers': 'Chickens, Layers',
+    'beef_cows': 'Cows, Beef',
+    'hogs': 'Pigs/Hogs/Swine',
+    'dairy_cows': 'Cows, Dairy',
+    'broilers': 'Chickens, Broilers',
+}
+
+
+def animal_population(geojson):
+    """
+    Given a GeoJSON shape, call MapShed's `animal_energy_units` method
+    to calculate the area-weighted county animal population. Returns a
+    dictionary to append to the outgoing JSON for analysis results.
+    """
+    geom = GEOSGeometry(geojson, srid=4326)
+    aeu_for_geom = animal_energy_units(geom)[2]
+    aeu_return_values = []
+
+    for animal, aeu_value in aeu_for_geom.iteritems():
+        aeu_return_values.append({
+            'type': ANIMAL_DISPLAY_NAMES[animal],
+            'aeu': int(aeu_value),
+        })
+
+    return {
+        'displayName': 'Animals',
+        'name': 'animals',
+        'categories': aeu_return_values
+    }
+
+
+def point_source_pollution(geojson):
+    """
+    Given a GeoJSON shape, retrieve point source pollution data
+    from the `ms_pointsource` or `ms_pointsource_drb` table to display
+    in the Analyze tab.
+
+    Returns a dictionary to append to the outgoing JSON for analysis
+    results.
+    """
+    geom = GEOSGeometry(geojson, srid=4326)
+    drb = geom.within(DRB)
+    table_name = get_point_source_table(drb)
+    sql = '''
+          SELECT city, state, npdes_id, mgd, kgn_yr, kgp_yr, latitude,
+                 longitude, {facilityname}
+          FROM {table_name}
+          WHERE ST_Intersects(geom, ST_SetSRID(ST_GeomFromText(%s), 4326))
+          '''.format(facilityname='facilityname' if drb else 'null',
+                     table_name=table_name)
+
+    with connection.cursor() as cursor:
+        cursor.execute(sql, [geom.wkt])
+
+        if cursor.rowcount != 0:
+            columns = [col[0] for col in cursor.description]
+            point_source_results = [
+                dict(zip(columns,
+                         [row[0], row[1], row[2],
+                          float(row[3]) if row[3] else None,
+                          float(row[4]) if row[4] else None,
+                          float(row[5]) if row[5] else None,
+                          float(row[6]) if row[6] else None,
+                          float(row[7]) if row[7] else None,
+                          row[8]]))
+                for row in cursor.fetchall()
+            ]
+        else:
+            point_source_results = []
+
+    return {
+        'displayName': 'Point Source',
+        'name': 'pointsource',
+        'categories': point_source_results
+    }
+
+
+def catchment_water_quality(geojson):
+    """
+    Given a GeoJSON shape, retrieve Catchment Water Quality data
+    from the `drb_catchment_water_quality` table to display
+    in the Analyze tab.
+
+    Returns a dictionary to append to the outgoing JSON for analysis
+    result
+    """
+    geom = GEOSGeometry(geojson, srid=4326)
+    table_name = 'drb_catchment_water_quality'
+    sql = '''
+          SELECT nord, areaha, tn_tot_kgy, tp_tot_kgy, tss_tot_kg,
+          tn_urban_k, tn_riparia, tn_ag_kgyr, tn_natural, tn_pt_kgyr,
+          tp_urban_k, tp_riparia, tp_ag_kgyr, tp_natural, tp_pt_kgyr,
+          tss_urban_, tss_rip_kg, tss_ag_kgy, tss_natura,
+          tn_yr_avg_, tp_yr_avg_, tss_concmg,
+          ST_AsGeoJSON(ST_Simplify(geom, 0.0003)) as geom
+          FROM {table_name}
+          WHERE ST_Intersects(geom, ST_SetSRID(ST_GeomFromText(%s), 4326))
+          '''.format(table_name=table_name)
+
+    with connection.cursor() as cursor:
+        cursor.execute(sql, [geom.wkt])
+
+        if cursor.rowcount != 0:
+            columns = [col[0] for col in cursor.description]
+            catchment_water_quality_results = [
+                # The TN, TP, and TSS values return as type Decimal,
+                # but we want floats.
+                dict(zip(columns,
+                         ([int(row[0]) if row[0] else None] +
+                          [float(value) if value else None
+                           for value in row[1:22]] +
+                          [row[22]])))
+                for row in cursor.fetchall()
+            ]
+        else:
+            catchment_water_quality_results = []
+    return {
+        'displayName': 'Water Quality',
+        'name': 'catchment_water_quality',
+        'categories': catchment_water_quality_results
+    }

src/mmw/apps/geoprocessing_api/tasks.py

@@ -0,0 +1,151 @@
+# -*- coding: utf-8 -*-
+from __future__ import print_function
+from __future__ import unicode_literals
+from __future__ import absolute_import
+
+import os
+import logging
+
+from ast import literal_eval as make_tuple
+
+from celery import shared_task
+
+import requests
+
+from django.conf import settings
+
+from apps.modeling.tr55.utils import aoi_resolution
+
+from apps.geoprocessing_api.calcs import (animal_population,
+                                          point_source_pollution,
+                                          catchment_water_quality,
+                                          )
+
+logger = logging.getLogger(__name__)
+
+DRB = 'drb'
+
+RWD_HOST = os.environ.get('RWD_HOST', 'localhost')
+RWD_PORT = os.environ.get('RWD_PORT', '5000')
+
+ACRES_PER_SQM = 0.000247105
+
+
+@shared_task
+def start_rwd_job(location, snapping, data_source):
+    """
+    Calls the Rapid Watershed Delineation endpoint
+    that is running in the Docker container, and returns
+    the response unless there is an out-of-watershed error
+    which raises an exception.
+    """
+    lat, lng = location
+    end_point = 'rwd' if data_source == DRB else 'rwd-nhd'
+    rwd_url = 'http://%s:%s/%s/%f/%f' % (RWD_HOST, RWD_PORT, end_point,
+                                         lat, lng)
+
+    # The Webserver defaults to enable snapping, uses 1 (true) 0 (false)
+    if not snapping:
+        rwd_url += '?snapping=0'
+
+    logger.debug('rwd request: %s' % rwd_url)
+
+    response_json = requests.get(rwd_url).json()
+    if 'error' in response_json:
+        raise Exception(response_json['error'])
+
+    return response_json
+
+
+@shared_task
+def analyze_animals(area_of_interest):
+    """
+    Given an area of interest, returns the animal population within it.
+    """
+    return {'survey': animal_population(area_of_interest)}
+
+
+@shared_task
+def analyze_pointsource(area_of_interest):
+    """
+    Given an area of interest, returns point sources of pollution within it.
+    """
+    return {'survey': point_source_pollution(area_of_interest)}
+
+
+@shared_task
+def analyze_catchment_water_quality(area_of_interest):
+    """
+    Given an area of interest in the DRB, returns catchment water quality data
+    within it.
+    """
+    return {'survey': catchment_water_quality(area_of_interest)}
+
+
+@shared_task(throws=Exception)
+def analyze_nlcd(result, area_of_interest=None):
+    if 'error' in result:
+        raise Exception('[analyze_nlcd] {}'.format(result['error']))
+
+    pixel_width = aoi_resolution(area_of_interest) if area_of_interest else 1
+
+    histogram = {}
+    total_count = 0
+    categories = []
+
+    # Convert results to histogram, calculate total
+    for key, count in result.iteritems():
+        total_count += count
+        histogram[make_tuple(key[4:])] = count  # Change {"List(1)":5} to {1:5}
+
+    for nlcd, (code, name) in settings.NLCD_MAPPING.iteritems():
+        categories.append({
+            'area': histogram.get(nlcd, 0) * pixel_width * pixel_width,
+            'code': code,
+            'coverage': float(histogram.get(nlcd, 0)) / total_count,
+            'nlcd': nlcd,
+            'type': name,
+        })
+
+    return {
+        'survey': {
+            'name': 'land',
+            'displayName': 'Land',
+            'categories': categories,
+        }
+    }
+
+
+@shared_task(throws=Exception)
+def analyze_soil(result, area_of_interest=None):
+    if 'error' in result:
+        raise Exception('[analyze_soil] {}'.format(result['error']))
+
+    pixel_width = aoi_resolution(area_of_interest) if area_of_interest else 1
+
+    histogram = {}
+    total_count = 0
+    categories = []
+
+    # Convert results to histogram, calculate total
+    for key, count in result.iteritems():
+        total_count += count
+        s = make_tuple(key[4:])  # Change {"List(1)":5} to {1:5}
+        s = s if s != settings.NODATA else 3  # Map NODATA to 3
+        histogram[s] = count + histogram.get(s, 0)
+
+    for soil, (code, name) in settings.SOIL_MAPPING.iteritems():
+        categories.append({
+            'area': histogram.get(soil, 0) * pixel_width * pixel_width,
+            'code': code,
+            'coverage': float(histogram.get(soil, 0)) / total_count,
+            'type': name,
+        })
+
+    return {
+        'survey': {
+            'name': 'soil',
+            'displayName': 'Soil',
+            'categories': categories,
+        }
+    }