Commits

Chris Mutel  committed e50ecda

Removed development cruft and Brightway stuff

  • Participants
  • Parent commits 52b032c

Comments (0)

Files changed (10)

File lcia_autocorrelation/ac_base.py

 import pysal
 from hashlib import md5 as universal_hash
 
+def fake_callback(message, **kwargs):
+    import pprint
+    pprint.pprint(message)
+    pprint.pprint(kwargs)
 
-def fake_callback(message, **kwargs):
-    pass
 
 class AutocorrelationBase(object):
     """

File lcia_autocorrelation/ac_brightway.py

-# -*- coding: utf-8 -*
-import os
-import pysal
-from numpy import zeros
-from django.conf import settings
-from brightway.core.array_proxies import OneDimensionalArrayProxy
-from brightway.core.utils.safe_hash import generate_safe_hash
-from brightway.apps.impact_assessment.models import Weighting
-from brightway.core.caching.cache_interface import CacheInterface
-# from brightway.core.regionalize.amoeba import AMOEBA
-from brightway.core.caching.errors import CacheError
-
-from hashlib import md5 as universal_hash
-from numpy import array
-from amoeba import DjangoAMOEBA as AMOEBA
-from ac_base import AutocorrelationBase
-from brightway.apps.geo import WORLD
-from shapefile import ShapefileWriter
-
-class AutocorrelationBrightway(AutocorrelationBase):
-    def __init__(self, process, method):
-        self.process = process
-        self.method = method
-        self.qs = Weighting.objects.filter(process=process, method=method
-            ).exclude(geography=WORLD()).exclude(amount__lte=0)
-        if not self.qs.exists():
-            raise ValueError("The selected process and method have no" + \
-                " regionalized weightings")
-        self.cache_dir = os.path.join(settings.BRIGHTWAY_ROOT_PATH, 
-            "core/regionalize/cache/")
-        self.mapping, self.variable = self.get_variable()
-        self.reverse = dict([(value, key) for key, value in self.mapping.\
-            iteritems()])
-        self.W = self.weighting_matrix()
-        # self.local = pysal.Moran_Local(self.variable.data, self.W)
-
-    # def get_variable(self, ds):
-    #     return array([x.amount for x in ds])
-
-    def construct_data_source(self, ds):
-        """Parse a Django GeoQuerySet of Brightway Weightings"""
-        ds = ds.exclude(geography=WORLD()).order_by('id')
-        self.method = ds[0].method
-        return ds
-
-    def get_variable(self):
-        mapping = dict([(key, index) for index, key in enumerate( 
-            self.qs.values_list('geography_id', flat=True).distinct())])
-        variable = OneDimensionalArrayProxy(zeros(len(mapping)), mapping)
-        for obj in self.qs:
-            variable[obj.geography_id] = obj.amount
-        return mapping, variable
-
-    def weighting_values(self, qs=None):
-        """Build a pysal spatial weights matrix."""
-        if not qs:
-            qs = self.qs
-        amoeba = AMOEBA(qs).go()
-        # AMOEBA produces results like: {geo1: {geo2: weight, geo3: weight}}. 
-        # Need to transform to: {geo1: [geo2, geo3]}, {geo1: [weight, weight]}
-        # See http://pysal.org/library/weights/weights.html#pysal.weights.W
-        neighbors = {}
-        weights = {}
-        for home in amoeba:
-            home_id = self.mapping[home]
-            neighbors[home_id] = []
-            weights[home_id] = []
-            for target, weight in amoeba[home].iteritems():
-                neighbors[home_id].append(self.mapping[target])
-                weights[home_id].append(weight)
-        return neighbors, weights
-
-    # def shapefile(self):
-    #     """Write Moran's I values and CFs to shapefile for inspection."""
-    #     sf = ShapefileWriter('lisa-output-%s-%s.shp' % (self.method, 
-    #         self.process), self.qs[0].geography.shape, 'lisa', id=int, 
-    #         moran_i=float, cf=float)
-    # 
-    #     for obj in self.qs:
-    #         sf.add_record(obj.geography.shape, id=obj.id, cf=obj.amount, 
-    #             moran_i=float(self.local.Is[self.mapping[obj.geography.id]]))
-    # 
-    #     sf.close()
-

File lcia_autocorrelation/amoeba.py

     def get_weights_dict(self):
         return self.ds.mapping_dict()
 
-
-class DjangoAMOEBA(AMOEBA):
-    NEIGHBOR_QUERY = """select f.id 
-        from geo_geography f
-        where ST_Intersects(
-            f.poly_shape, (
-                select ST_Buffer(ST_Union(poly_shape), 0.00001) 
-                    from geo_geography 
-                    where id in %s -- Existing neighbors
-            )
-        ) and f.id in %s -- All possible geometries
-        and f.id not in %s -- Excluded
-        and f.id not in %s; -- Already present"""
-
-    def expand(self, tiers, excluded):
-        """
-Find Geographies that border beighbors but are not in excluded.
-
-Uses buffer, as border can be inexact.
-        """
-        if tiers.level > 10:
-            # One could continue to expand, but there is not really any 
-            # point; the spatial weights are < 1e-6 anyway.
-            return []
-        cursor = sql(self.NEIGHBOR_QUERY, (tuple(tiers.all_geo), 
-            self.geo_ids, tuple(excluded), tuple(tiers.all_geo)))
-        return [x[0] for x in cursor.fetchall()]
-
-    def get_oid(self, o):
-        return o.geography_id
-
-    def write_to_shapefile(self, tiers, obj):
-        """Write tiers to shapefile for inspection."""
-        sf = ShapefileWriter('amoeba-output-%s-%s.shp' % (obj.process, 
-            obj.id), obj.geography.shape, "amoeba", id=int, weight=float)
-
-        feature_def = layer.GetLayerDefn()
-        for key, weight in tiers.normalized.iteritems():
-            geo = Geography.objects.get(id=key).shape
-            sf.add_record(geo, id=key, weight=weight)
-
-        sf.close()
-
-    def get_weights_dict(self):
-        return dict([(o['geography_id'], o['amount']) for \
-            o in self.ds.values('geography_id', 'amount')])

File lcia_autocorrelation/contiguity_weights.py

-# -*- coding: utf-8 -*
-from __future__ import division
-from numpy import *
-from math import sqrt
-try:
-    from brightway.core.sql import sql
-    from brightway.apps.geo.models import Geography
-except ImportError:
-    brightway = False
-import progressbar as pb
-from tier_dictionary import TierDictionary
-from shapefile import ShapefileWriter, ShapefileReader
-from utils import buffer_gdal_geometry
-
-
-class ContiguityWeights(object):
-    def __init__(self, qs):
-        self.qs = qs
-        self.geo_ids = list(qs.values_list('geography_id', flat=True))
-
-    def iterate(self, obj):
-        first_level_query = None
-
-        qs = Geography.objects.filter(id__in=self.geo_ids, 
-            poly_shape__touches=Geography.objects.get(id=obj).poly_shape 
-            ).exclude(obj)
-        first_level = None
-        # ke: {geo1: {geo2: weight, geo3: weight}}. 
-
-
-
-
-# class NoNeighborsError(StandardError):
-#     pass
-# 
-# 
-# class AMOEBA(object):
-#     """
-# Calculates a spatial weights matrix using AMOEBA algorithm.
-# 
-# The basic algorithm is::
-# 
-#     for each geographic unit:
-#         neighbors = [original unit]
-#         calculate g_star of neigbors
-#         while True
-#             find all neighbors, excluding excluded neighbors
-#             find set of new neighbors that maximizes g_star
-#             if new g_star < old g_star
-#                 break
-#             neighbors += new neighbors
-#             excluded += rejected neighbors
-# 
-# .. rubric:: Inputs
-# 
-# * obj: A shapefile or Django queryset of geography objects and weights.
-# 
-# Aldstadt, J. & Getis, A. (2006). Using AMOEBA to Create a Spatial Weights Matrix and Identify Spatial Clusters. Geographical Analysis, 38(4), 327--343. 
-#     """
-#     def __init__(self, ds, export_shp=False, pb=False):
-#         self.ds = ds
-#         self.export_shp = export_shp
-#         self.pb = pb
-#         self.weights_dict = self.get_weights_dict()
-#         self.geo_ids = tuple(self.weights_dict.keys())
-#         self.number = float(len(self.geo_ids))
-#         values = array(self.weights_dict.values())
-#         self.average = float(average(values))
-#         self.pseudo_std = float(sqrt(
-#             (values**2).sum() / self.number - self.average**2)
-#             )
-# 
-#     def go(self, objs=None):
-#         # Used for testing, with new QS
-#         objs = objs or self.ds
-# 
-#         results = {}
-#         try:
-#             count = objs.count()
-#         except TypeError:
-#             count = len(objs)
-# 
-#         if self.pb:
-#             widgets = ['Weightings: ', pb.Percentage(), ' ', 
-#                 pb.Bar(marker=pb.RotatingMarker()), ' ', pb.ETA()]
-#             pbar = pb.ProgressBar(widgets=widgets, maxval=count).start()
-# 
-#         for index, obj in enumerate(objs):
-#             oid = self.get_oid(obj)
-#             results[oid] = self.iterate(obj)
-# 
-#             if self.pb:
-#                 pbar.update(index)
-# 
-#         if self.pb:
-#             pbar.finish()
-# 
-#         return results
-# 
-#     def g_star(self, values_sum, values_num):
-#         """
-# Calculate g-star statistic for spatial unit in relation to neighbors. 
-# 
-# The Gi* star statistic can be considered a standard normal variate.  So, the the cumulative probability of G*i(k_max) would be the cdf value for the standard normal distribution of the value G*i(k_max).  This value could be close to 1, but could also be much closer to zero in some cases.  So, using equation 2 of the AMOEBA paper, elements in the k_max level of contiguity would have zero weights with unit i.  These units are at the edge of the ecotope centered on i.  The units closer to i would have positive weights that increase as proximity to i increases. The value Gi*(0) is the Gi* value for the unit itself.
-#         """
-#         return (values_sum - self.average * values_num) / (
-#             self.pseudo_std * sqrt(
-#                 (self.number - values_num) * values_num / (self.number - 1)
-#                 )
-#             )
-# 
-#     def get_trend(self, obj, tiers):
-#         """Do a first round of expansion, and determine the trend direction"""
-#         neighbors = self.expand(tiers, set([0]))
-#         if not neighbors:
-#             raise NoNeighborsError
-# 
-#         # Create array with col 1: geo ids, col 2: weights
-#         data = array([self.weights_dict[key] for key in neighbors])
-#         # Sort by weight ascending
-#         data = data[argsort(data)]
-# 
-#         orig_value = float(sum(tiers.all_weights()))
-#         if tiers.score < 0:
-#             operator = max
-#         else:
-#             operator = min
-# 
-#         up_score = operator([self.g_star(float(data[:index].sum()) + \
-#             orig_value, index + 1) for index in range(1, data.shape[0] + 1)])
-#         down_score = operator([self.g_star(float(data[::-1][:index].sum()) + \
-#             orig_value, index + 1) for index in range(1, data.shape[0] + 1)])
-#         if tiers.score < 0:
-#             trend = "downwards" if (down_score < up_score) else "upwards"
-#         else:
-#             trend = "downwards" if (down_score > up_score) else "upwards"
-#         return trend
-# 
-#     def iterate(self, obj):
-#         oid = self.get_oid(obj)
-#         weight = self.weights_dict[oid]
-#         tiers = TierDictionary(oid, weight, self.g_star(weight, 1))
-#         excluded = set([0])
-# 
-#         try:
-#             trend = self.get_trend(obj, tiers)
-#         except NoNeighborsError:
-#             return tiers.normalized
-# 
-#         while True:
-#             neighbors = self.expand(tiers, excluded)
-#             if not neighbors:
-#                 # No new neighbors found; stop iteration
-#                 break
-# 
-#             # Create array with col 1: geo ids, col 2: weights
-#             data = hstack((array(neighbors).reshape((-1,1)), array([ 
-#                 self.weights_dict[key] for key in neighbors]).reshape((-1,1) 
-#                 )))
-#             # Sort by weight
-#             if trend == "upwards":
-#                 data = data[argsort(data[:,1])]
-#             else:
-#                 data = data[argsort(data[:,1])[::-1]]
-# 
-#             already_values = tiers.all_weights()
-#             already_number = float(already_values.shape[0])
-#             already_sum = float(already_values.sum())
-#             base_score = loop_score = tiers.score
-#             new_neighbors = set()
-# 
-#             for index in range(1, data.shape[0]):
-#                 values = data[:index, 1]
-# 
-#                 new_score = self.g_star(float(values.sum()) + already_sum,
-#                     float(values.shape[0]) + already_number)
-# 
-#                 if not (new_score > loop_score > 0 or new_score < \
-#                         loop_score < 0):
-#                     break
-#                 else:
-#                     new_neighbors = set([int(x) for x in data[:index,0]])
-#                     loop_score = new_score
-# 
-#             if loop_score == base_score:
-#                 # No better set of neighbors was found
-#                 break
-#             else:
-#                 # Add neighbors and set new high score
-#                 excluded = excluded.union(set(neighbors).difference( 
-#                     set(new_neighbors)))
-#                 tiers.add(new_neighbors, self.weights_dict, loop_score)
-# 
-#         if self.export_shp:
-#             self.write_to_shapefile(tiers, obj)
-# 
-#         return tiers.normalized
-# 
-# 
-# class ShapefileAMOEBA(AMOEBA):
-#     """AMOEBA algorithm adapted to the objects defined by the Python GDAL shapefile interface."""
-#     def expand(self, tiers, excluded):
-#         if tiers.level > 10:
-#             # One could continue to expand, but there is not really any 
-#             # point; the spatial weights are < 1e-6 anyway.
-#             return []
-#         return list(
-#             reduce(
-#                 set.union, 
-#                 [set(self.ds.intersects( 
-#                     buffer_gdal_geometry(
-#                         self.ds.layer[index].geom, 0.00001), 
-#                     set.union(excluded, set([index]), tiers.all_geo)
-#                 )) for index in tiers.geo[tiers.level]]
-#             ))
-# 
-#     def get_oid(self, o):
-#         return o.fid
-# 
-#     def write_to_shapefile(self, tiers, obj):
-#         """Write tiers to shapefile for inspection."""
-#         filename = self.ds.filepath.split("/")[-1].split(".")[0]
-#         sf = ShapefileWriter('amoeba-output-%s-%s' % (filename, 
-#             self.get_oid(obj)), self.ds.layer[0].geom, "amoeba", id=int, 
-#             weight=float, cf=float)
-# 
-#         for key, weight in tiers.normalized.iteritems():
-#             obj = self.ds.layer[key]
-#             sf.add_record(obj.geom, id=key, weight=weight, 
-#                 cf=obj.get(self.ds.amount_field))
-# 
-#         sf.close()
-# 
-#     def get_weights_dict(self):
-#         return self.ds.mapping_dict()
-# 
-# 
-# class DjangoAMOEBA(AMOEBA):
-#     NEIGHBOR_QUERY = """select f.id 
-#         from geo_geography f
-#         where ST_Intersects(
-#             f.poly_shape, (
-#                 select ST_Buffer(ST_Union(poly_shape), 0.00001) 
-#                     from geo_geography 
-#                     where id in %s -- Existing neighbors
-#             )
-#         ) and f.id in %s -- All possible geometries
-#         and f.id not in %s -- Excluded
-#         and f.id not in %s; -- Already present"""
-# 
-#     def expand(self, tiers, excluded):
-#         """
-# Find Geographies that border beighbors but are not in excluded.
-# 
-# Uses buffer, as border can be inexact.
-#         """
-#         if tiers.level > 10:
-#             # One could continue to expand, but there is not really any 
-#             # point; the spatial weights are < 1e-6 anyway.
-#             return []
-#         cursor = sql(self.NEIGHBOR_QUERY, (tuple(tiers.all_geo), 
-#             self.geo_ids, tuple(excluded), tuple(tiers.all_geo)))
-#         return [x[0] for x in cursor.fetchall()]
-# 
-#     def get_oid(self, o):
-#         return o.geography_id
-# 
-#     def write_to_shapefile(self, tiers, obj):
-#         """Write tiers to shapefile for inspection."""
-#         sf = ShapefileWriter('amoeba-output-%s-%s.shp' % (obj.process, 
-#             obj.id), obj.geography.shape, "amoeba", id=int, weight=float)
-# 
-#         feature_def = layer.GetLayerDefn()
-#         for key, weight in tiers.normalized.iteritems():
-#             geo = Geography.objects.get(id=key).shape
-#             sf.add_record(geo, id=key, weight=weight)
-# 
-#         sf.close()
-# 
-#     def get_weights_dict(self):
-#         return dict([(o['geography_id'], o['amount']) for \
-#             o in self.ds.values('geography_id', 'amount')])

File lcia_autocorrelation/export_brightway.py

-import os
-from shapefile import ShapefileWriter
-
-def export_weighting_qs_to_shapefile(qs):
-    filename = "bw-export-%s-%s" % (qs[0].method, qs[0].process)
-    filepath = os.path.join(os.path.dirname(__file__), "django-export", 
-        filename)
-    sf = ShapefileWriter(filepath, qs[0].geography.poly_shape, 'cfs', 
-        cf=float)
-    for obj in qs:
-        sf.add_record(obj.geography.poly_shape, cf=obj.amount)
-    sf.close()
-
-if __name__ == "__main__":
-    from brightway.apps.geo import WORLD
-    from brightway.apps.impact_assessment.models import Weighting, Method
-    for x in (142, 146, 294, 295, 296, 297):
-        print Method.objects.get(id=x)
-        if x == 146:
-            # Nitrogen oxide
-            qs = Weighting.objects.filter(method__id=x, process__id=249114 
-                ).exclude(geography=WORLD())
-        elif x == 142:
-            qs = Weighting.objects.filter(method__id=x, process__id=249186 
-                ).exclude(geography=WORLD())
-            # Sulfur dioxide
-        else:
-            qs = Weighting.objects.filter(method__id=x).exclude( 
-                geography=WORLD())
-        export_weighting_qs_to_shapefile(qs)

File lcia_autocorrelation/multi.py

-import os
-import cPickle as pickle
-from ac_shapefile import AutocorrelationShapefile
-from multiprocessing import Pool
-
-def ac_function(filename):
-    print "Starting:", filename
-    a = AutocorrelationShapefile("django-export/" + filename)
-    a.global_autocorrelation()
-    print "Finishing: %s, %s" % (filename, a.score)
-    return (filename, a.score)
-    # ls = ShapefileLISA("data2/" + filename)
-    # ls.export_shapefile()
-    # with open(filename[:-4] + ".pickle", "wb") as f:
-    #     pickle.dump(ls.local, f)
-    # print "Finishing:", filename, ls.local.Is.sum() / len(ls.local.Is)
-    # return (filename, ls.local.Is.sum() / len(ls.local.Is))
-
-if __name__ == '__main__':
-    # files = [x for x in os.listdir("data2") if x[:3] == "agg" and x[-3:
-    #     ] == "shp"]
-    # data = [(10, 0.338), 
-    #     (12, 0.234), 
-    #     (14, 0.232), 
-    #     (16, 0.36), 
-    #     (20, 0),
-    #     (30, 0),
-    #     (59, 0.391),
-    #     (109, 0.477), 
-    #     (158, 0.433), 
-    #     (208, 0.387), 
-    #     (257, 0.540),
-    #     (307, 0.465),
-    #     (356, 0.375),
-    #     (455, 0.396),
-    #     (554, 0.425),
-    #     (653, 0.614),
-    # ]
-    # for d in data:
-    #     files = [f for f in files if ".%s." % d[0] not in f]
-    files = [x for x in os.listdir("django-export") if x[-3:] == "shp"]
-    pool = Pool(processes=2)
-    print pool.map(ac_function, files)
-

File lcia_autocorrelation/num_objects.py

-import os
-from shapefile import ShapefileReader
-
-def get_length(filename):
-    sf = ShapefileReader(filename)
-    print len(sf.layer), filename
-
-def process_dir(dirname):
-    files = filter(lambda x: x[-4:] == ".shp", os.listdir(dirname))
-    for o in files:
-        get_length(os.path.join(os.path.dirname(__file__), dirname, o))
-
-if __name__ == '__main__':
-    process_dir("data")
-    process_dir("data2")
-    process_dir("django-export")
-

File lcia_autocorrelation/single.py

-import os
-import cPickle as pickle
-import json
-from autocorrelation_shapefile import AutocorrelationShapefile
-from lisa_shapefile import ShapefileLISA
-
-def f(filename):
-    print "Starting:", filename
-    a = AutocorrelationShapefile("data2/" + filename)
-    a.global_autocorrelation()
-    return (filename, a.score)
-    # ls = ShapefileLISA("data/" + filename)
-    # ls.export_shapefile()
-    # with open(filename[:-4] + ".pickle", "wb") as f:
-    #     pickle.dump(ls.local, f)
-    # return filename, ls.local.Is.sum() / len(ls.local.Is)
-
-if __name__ == '__main__':
-    files = [x for x in os.listdir("data2") if x[:3] == "agg" and x[-3:
-        ] == "shp"]
-    results = []
-    for filename in files:
-        results.append((filename, f(filename)))
-    with open("data2.json", "w") as f:
-        json.dump(results, f)

File lcia_autocorrelation/tier_dictionary.py

 from __future__ import division
 from numpy import array
 from scipy.stats import norm
-try:
-    from brightway.core.utils.flatten import flatten
-except ImportError:
-    from utils import flatten
+from utils import flatten
 
 class TierDictionary(object):
     """Data structure for storing AMOEBA ecotopes; keeps level and objects."""

File lcia_autocorrelation/write_to_shapefile.py

-from brightway.apps.inventory.models import *
-from brightway.apps.impact_assessment.models import *
-from brightway.apps.geo import WORLD
-from shapefile import ShapefileWriter
-
-
-class MethodShapefileWriter(object):
-    def __init__(self, process, method):
-        self.process = process
-        self.method = method
-        self.qs = Weighting.objects.filter(process=process, method=method, 
-            amount__gt=0).exclude(geography=WORLD())
-        filename = str(self.method.compoundname.replace(" ", "-").strip())
-        self.writer = ShapefileWriter(filename, 
-            self.qs[0].geography.poly_shape, cf=float)
-        for o in self.qs:
-            self.writer.add_record(o.geography.poly_shape, cf=o.amount)
-        self.writer.close()
-