In [ ]:
"""
Generate sub-buffers of a buffer
"""

In [1]:
import os
from osgeo import ogr

from senpy.geopnd.flickr import photos_location
from senpy.ogr.analysis.proximity.buffers import dic_buffer_array_to_shp
from senpy.pgsql.db import create_db
from senpy.pgsql.data.w import create_table_by_query
from senpy.topsql import shp_to_pgsql
from senpy.topsql.pnd import geodf_to_pgsql
from senpy.toshp.pnd import pandas_to_shp
from senpy.toshp.psql import psql_to_shp

In [ ]:
# Inputs

inBuffers = {"lisbon"    : {
    'x'      : -89004.994779, # in meters
    'y'      : -102815.866054, # in meters
    'radius' : 10000,
    'epsg'   : 3763
}}


keywords = ['flood', 'accident', 'fire apartment', 'graffiti', 'homeless']

workspace = r'C:\gis\dsn_collect'

conParam = {
    "HOST" : "localhost", "PORT" : "5432", "PASSWORD" : "admin",
    "USER" : "postgres", "DB" : "flickr_db", "TEMPLATE" : "postgis_template"
}

In [2]:
# Inputs

inBuffers = {'london'    : {
    'x'      : -14210.551441, # in meters
    'y'      : 6711542.47559, # in meters
    'radius' : 10000,
    'epsg'   : 3857
}}


keywords = ['flood', 'accident', 'fire apartment', 'graffiti', 'homeless']

workspace = r'C:\gis\dsn_collect'

conParam = {
    "HOST" : "localhost", "PORT" : "5432", "PASSWORD" : "admin",
    "USER" : "postgres", "DB" : "flickr_db", "TEMPLATE" : "postgis_template"
}

In [ ]:
# Inputs

inBuffers = {"rome"      : {
    'x'      : 1388522.02296, # in meters
    'y'      : 5145760.3812, # in meters
    'radius' : 10000,
    'epsg'   : 3857
}}


keywords = ['flood', 'accident', 'fire apartment', 'graffiti', 'homeless']

workspace = r'C:\gis\dsn_collect'

conParam = {
    "HOST" : "localhost", "PORT" : "5432", "PASSWORD" : "admin",
    "USER" : "postgres", "DB" : "flickr_db", "TEMPLATE" : "postgis_template"
}

In [ ]:
# Inputs

inBuffers = {"paris"     : {
    'x'      : 261611.316871, # in meters
    'y'      : 6250363.67084, # in meters
    'radius' : 10000,
    'epsg'   : 3857
}}


keywords = ['flood', 'accident', 'fire apartment', 'graffiti', 'homeless']

workspace = r'C:\gis\dsn_collect'

conParam = {
    "HOST" : "localhost", "PORT" : "5432", "PASSWORD" : "admin",
    "USER" : "postgres", "DB" : "flickr_db", "TEMPLATE" : "postgis_template"
}

In [ ]:
# Inputs

inBuffers = {"frankfurt" : {
    'x'      : 965640.304444, # in meters
    'y'      : 6464913.8033, # in meters
    'radius' : 10000,
    'epsg'   : 3857
}}


keywords = ['flood', 'accident', 'fire apartment', 'graffiti', 'homeless']

workspace = r'C:\gis\dsn_collect'

conParam = {
    "HOST" : "localhost", "PORT" : "5432", "PASSWORD" : "admin",
    "USER" : "postgres", "DB" : "flickr_db", "TEMPLATE" : "postgis_template"
}

In [3]:
def getBuffers(x, y, radius):
    sub_buf = ['north', 'northeast', 'east', 'southeast',
               'south', 'southwest', 'west', 'northwest']
    
    lstSubBuffer = []

    for cardeal in sub_buf:
        if cardeal == 'north':
            _y = y + (radius / 2)
        
        elif cardeal == 'northeast' or cardeal == 'northwest':
            _y =  y + ((radius)**2 / 8.0)**0.5
    
        elif cardeal == 'south':
            _y = y - (radius / 2)
        
        elif cardeal == 'southwest' or cardeal == 'southeast':
            _y = y - ((radius)**2 / 8.0)**0.5
    
        else:
            _y = y
    
    
        if cardeal == 'west':
            _x = x - (radius / 2)
        
        elif cardeal == 'southwest' or cardeal == 'northwest':
            _x = x - ((radius)**2 / 8.0)**0.5
    
        elif cardeal == 'east':
            _x = x + (radius / 2)
        
        elif cardeal == 'southeast' or cardeal == 'northeast':
            _x = x + ((radius)**2 / 8.0)**0.5
    
        else:
            _x = x
    
        lstSubBuffer.append({'X' : _x, 'Y' : _y, 'RADIUS' : radius / 2, 'cardeal' : cardeal})
    
    return lstSubBuffer

In [4]:
for city in inBuffers:
    inBuffers[city]["list_buffer"] = [{
        'X' : inBuffers[city]["x"], 'Y' : inBuffers[city]["y"], 'RADIUS' : inBuffers[city]['radius'],
        'cardeal' : 'major'
    }] + getBuffers(inBuffers[city]["x"], inBuffers[city]["y"], inBuffers[city]["radius"])
    
    dic_buffer_array_to_shp(
        inBuffers[city]["list_buffer"],
        os.path.join(workspace, 'buffers_{}.shp'.format(city)),
        inBuffers[city]['epsg'], fields={'cardeal' : ogr.OFTString}
    )

In [5]:
# Retrieve data for each keyword and buffer
# Record these elements in the dataframe

for city in inBuffers:
    c = None
    tblData = None
    for bf in inBuffers[city]["list_buffer"]:
        for k in keywords:
            tmpData = photos_location(
                bf, inBuffers[city]["epsg"], keyword=k,
                epsg_out=inBuffers[city]["epsg"],
                onlySearchAreaContained=False
            )
        
            if type(tmpData) == int:
                print "NoData finded for buffer '{}' and keyword '{}'".format(
                    bf['cardeal'], k
                )
                continue
        
            tmpData["keyword"] = k
            tmpData["buffer_or"] = bf['cardeal']
        
            if not c:
                tblData = tmpData
                c = 1
            else:
                tblData = tblData.append(tmpData, ignore_index=True)
        
    inBuffers[city]["data"] = tblData

In [ ]:
# Intersect all retrieved data with all buffers
# See if a point was obtained with all buffers which him intersects

In [6]:
create_db(conParam, conParam["DB"])
conParam["DATABASE"] = conParam["DB"]

In [ ]:
conParam["DATABASE"] = conParam["DB"]

In [7]:
_ct = inBuffers.keys()[0]
dataColumns = [
    c for c in inBuffers[_ct]["data"].columns.values
    if c != 'geom' and c != 'keyword' and c != 'buffer_or' and c != 'geometry'
]
print dataColumns


[u'datetaken', u'dateupload', u'title', 'fid', 'url', 'description']

In [ ]:
print inBuffers["paris"]["data"]

In [8]:
for city in inBuffers:
    cols = inBuffers[city]["data"].columns.values
    
    if 'geometry' in cols:
        cgeom = "geometry"
    else:
        cgeom = "geom"
    
    print cgeom
    
    geodf_to_pgsql(
        conParam, inBuffers[city]["data"],
        'flickr_{}'.format(city),
        inBuffers[city]["epsg"], 'POINT', colGeom=cgeom
    )
    inBuffers[city]["table"] = 'flickr_{}'.format(city)


geometry

In [9]:
# Send buffers data to postgresql
for city in inBuffers:
    shp_to_pgsql(
        os.path.join(workspace, 'buffers_{}.shp'.format(city)),
        inBuffers[city]["epsg"], conParam, 'buffers_{}'.format(city)
    )
    
    inBuffers[city]["pg_buffer"] = 'buffers_{}'.format(city)

In [10]:
q = (
    "SELECT srcdata.*, "
    "array_agg(buffersg.cardeal ORDER BY buffersg.cardeal) "
    "AS intersect_buffer FROM ("
        "SELECT {cols}, keyword, geom, "
        "array_agg(buffer_or ORDER BY buffer_or) AS extracted_buffer "
        "FROM {pgtable} "
        "GROUP BY {cols}, keyword, geom"
    ") AS srcdata, ("
        "SELECT cardeal, geom AS bfg FROM {bftable}"
    ") AS buffersg "
    "WHERE ST_Intersects(srcdata.geom, buffersg.bfg) IS TRUE "
    "GROUP BY {cols}, keyword, geom, extracted_buffer"
)

for city in inBuffers:
    inBuffers[city]["table"] = create_table_by_query(
        conParam, 
        "filter_{}".format(inBuffers[city]["table"]),
        q.format(
            cols=", ".join(dataColumns),
            pgtable=inBuffers[city]["table"],
            bftable=inBuffers[city]["pg_buffer"]
        )
    )

In [11]:
# Get Buffers table with info related
# pnt_obtidos = nr pontos obtidos usando esse buffer
# pnt_obtidos_fora = nt pontos obtidos fora desse buffer, mas obtidos com ele
# pnt_intersect = nt pontos que se intersectam com o buffer
# pnt_intersect_non_obtain = nr pontos que se intersectam mas nao foram obtidos como buffer

q = (
    "SELECT main.*, get_obtidos.pnt_obtidos, obtidos_fora.pnt_obtidos_fora, "
    "intersecting.pnt_intersect, int_not_obtained.pnt_intersect_non_obtain "
    "FROM {bf_table} AS main "
    "LEFT JOIN ("
        "SELECT gid, cardeal, COUNT(gid) AS pnt_obtidos "
        "FROM {bf_table} AS bf "
        "INNER JOIN {dt_table} AS dt "
        "ON bf.cardeal = ANY(dt.extracted_buffer) "
        "GROUP BY gid, cardeal"
    ") AS get_obtidos ON main.gid = get_obtidos.gid "
    "LEFT JOIN ("
        "SELECT gid, cardeal, COUNT(gid) AS pnt_obtidos_fora "
        "FROM {bf_table} AS bf "
        "INNER JOIN {dt_table} AS dt "
        "ON bf.cardeal = ANY(dt.extracted_buffer) "
        "WHERE ST_Intersects(bf.geom, dt.geom) IS NOT TRUE "
        "GROUP BY gid, cardeal"
    ") AS obtidos_fora ON main.gid = obtidos_fora.gid "
    "LEFT JOIN ("
        "SELECT gid, cardeal, COUNT(gid) AS pnt_intersect "
        "FROM {bf_table} AS bf "
        "INNER JOIN {dt_table} AS dt "
        "ON bf.cardeal = ANY(dt.intersect_buffer) "
        "GROUP BY gid, cardeal"
    ") AS intersecting ON main.gid = intersecting.gid "
    "LEFT JOIN ("
        "SELECT gid, cardeal, COUNT(gid) AS pnt_intersect_non_obtain "
        "FROM {bf_table} AS bf "
        "INNER JOIN {dt_table} AS dt "
        "ON bf.cardeal = ANY(dt.intersect_buffer) "
        "WHERE NOT (bf.cardeal = ANY(dt.extracted_buffer)) "
        "GROUP BY gid, cardeal"
    ") AS int_not_obtained "
    "ON main.gid = int_not_obtained.gid "
    "ORDER BY main.gid"
)

for city in inBuffers:
    inBuffers[city]["pg_buffer"] = create_table_by_query(
        conParam, 
        "dt_{}".format(inBuffers[city]["pg_buffer"]),
        q.format(
            bf_table=inBuffers[city]["pg_buffer"],
            dt_table=inBuffers[city]["table"],
        )
    )

In [12]:
# Get Points table with info related
# nobtido = n vezes um ponto foi obtido
# obtido_e_intersect = n vezes um ponto foi obtido usando um buffer com o qual se intersecta
# obtido_sem_intersect = n vezes um ponto foi obtido usando um buffer com o qual não se intersecta
# nintersect = n vezes que um ponto se intersecta com um buffer
# intersect_sem_obtido = n vezes que um ponto nao foi obtido apesar de se intersectar com o buffer

q = (
    "SELECT {cols}, dt.keyword, dt.geom, "
    "CAST(dt.extracted_buffer AS text) AS extracted_buffer, "
    "CAST(dt.intersect_buffer AS text) AS intersect_buffer, "
    "array_length(extracted_buffer, 1) AS nobtido, "
    "COUNT(fid) AS obtido_e_intersect, "
    "(array_length(extracted_buffer, 1) - COUNT(fid)) AS obtido_sem_intersect, "
    "array_length(intersect_buffer, 1) AS nintersect, "
    "(array_length(intersect_buffer, 1) - COUNT(fid)) AS intersect_sem_obtido "
    "FROM {dt_table} AS dt "
    "INNER JOIN {bf_table} AS bf "
    "ON bf.cardeal = ANY(dt.extracted_buffer) "
    "WHERE ST_Intersects(bf.geom, dt.geom) IS TRUE "
    "GROUP BY {cols}, dt.keyword, dt.geom, "
    "dt.extracted_buffer, dt.intersect_buffer"
)

for city in inBuffers:
    inBuffers[city]["table"] = create_table_by_query(
        conParam, "info_{}".format(city),
        q.format(
            dt_table=inBuffers[city]["table"],
            bf_table=inBuffers[city]["pg_buffer"],
            cols=", ".join(["dt.{}".format(x) for x in dataColumns])
        )
    )

In [13]:
# Export data to ESRI Shapefile
for city in inBuffers:
    psql_to_shp(
        conParam,
        inBuffers[city]["table"],
        os.path.join(workspace, "{}.shp".format(inBuffers[city]["table"])),
        epsg=inBuffers[city]["epsg"],
        geom_col='geom'
    )
    
    psql_to_shp(
        conParam,
        inBuffers[city]["pg_buffer"],
        os.path.join(workspace, "{}.shp".format(inBuffers[city]["pg_buffer"])),
        epsg=inBuffers[city]["epsg"], geom_col="geom"
    )

In [ ]: