[GRASS-user] Obtain land use within 50 km buffers around vector points

[Johannes Radinger]

Hi Moritz,
Hi All...

here my python code snippets that worked for me. I think there might be
many ways to simplify the code...

import grass.script as grass
import sqlite3
import pandas

# Def variable names
respondents_coord="respondents_coord"
respondents_points="respondents_points"
respondents_buffer="respondents_buffer"
respondents_buffer_cat_2_layer="respondents_buffer_cat_2_layer"
respondents_buffer_cat_2_raster="respondents_buffer_cat_2_raster"
correspondance_table="correspondance_table"
clc_summary_table="clc_summary_table"
clc_summary_pivot="clc_summary_pivot"

# Import points from X-Y table
grass.run_command("v.in.db",
overwrite=True,
table=respondents_coord,
output=respondents_points)

#### Calculate buffers around points
grass.run_command("v.buffer",
overwrite=True,
flags="t",
input=respondents_points,
type="point",
output=respondents_buffer,
distance=50000)
grass.run_command("v.category",
overwrite=True,
input=respondents_buffer,
option="add",
layer=2,
out=respondents_buffer_cat_2_layer)

# Get CORINE land use special classes per buffer
correspondance_table_values = []
for line in grass.read_command('v.category',
input_=respondents_buffer_cat_2_layer,
layer='1,2',
option='print').splitlines():
layers=line.split('|')
l1 = layers[0].split('/')
l2 = layers[1]
for cat in l1:
correspondance_table_values.append((cat, l2))

grass.run_command("db.execute",
sql="CREATE TABLE correspondance_table (buffer_cat_1 INTEGER, buffer_cat_2
INTEGER)")

conn = sqlite3.connect('/path/to/sqlite/sqlite.db')
cur = conn.cursor()
cur.executemany('INSERT INTO correspondance_table (buffer_cat_1,
buffer_cat_2) values (?,?)',correspondance_table_values)
conn.commit()
conn.close()

grass.run_command("v.to.rast",
overwrite=True,
input=respondents_buffer_cat_2_layer,
layer=2,
output=respondents_buffer_cat_2_raster,
use="cat")

if clc_summary_pivot in
grass.read_command("db.tables",flags="p").split("\n"):
grass.run_command("db.droptable",
flags="f",
table=clc_summary_pivot)

grass.run_command("db.execute",
sql="CREATE TABLE {} (buffer_cat_2 INTEGER, CLC_built_up INTEGER,
CLC_arable INTEGER, CLC_permanent_crops INTEGER, CLC_grassland INTEGER,
CLC_forest INTEGER, CLC_others INTEGER, CLC_intertidal_coastal INTEGER,
CLC_water_bodies INTEGER, CLC_sea INTEGER)".format(clc_summary_pivot))

clc_values = []
for line in grass.read_command("r.stats",
flags="cn",
input="{},{}".format(respondents_buffer_cat_2_raster,"CLC_reclass at Corine_LandCover
")).splitlines():
clc_values.append((line.split(' ')[0], line.split(' ')[1],line.split('
')[2]))

df = pandas.DataFrame(clc_values, columns=['buffer_cat_2', 'CLC_class',
'count'])
CLC_class_cols = ["1","2","3","4","5","6","7","8","9"]
clc_summary_pivot_values = [tuple(x) for x in
df.pivot(index='buffer_cat_2', columns='CLC_class',
values='count').reindex(columns=CLC_class_cols).fillna(0).astype(int).to_records(index=True)]

conn = sqlite3.connect('/path/to/sqlite/sqlite.db')
cur = conn.cursor()
cur.executemany('INSERT INTO {} (buffer_cat_2, CLC_built_up, CLC_arable,
CLC_permanent_crops, CLC_grassland, CLC_forest, CLC_others,
CLC_intertidal_coastal, CLC_water_bodies, CLC_sea) values
(?,?,?,?,?,?,?,?,?,?)'.format(clc_summary_pivot),clc_summary_pivot_values)
conn.commit()
conn.close()

grass.run_command("db.execute",
sql="CREATE TABLE {} AS SELECT buffer_cat_1,SUM(CLC_built_up) AS
CLC_built_up, SUM(CLC_arable) AS CLC_arable, SUM(CLC_permanent_crops) AS
CLC_permanent_crops, SUM(CLC_grassland) AS CLC_grassland, SUM(CLC_forest)
AS CLC_forest, SUM(CLC_others) AS CLC_others, SUM(CLC_intertidal_coastal)
AS CLC_intertidal_coastal, SUM(CLC_water_bodies) AS CLC_water_bodies,
SUM(CLC_sea) AS CLC_sea  FROM {} AS A LEFT JOIN {} AS B ON
A.buffer_cat_2=B.buffer_cat_2 GROUP BY
buffer_cat_1".format(clc_summary_table,correspondance_table,clc_summary_pivot))

# Join information on land use back to original center points of buffers
grass.run_command("v.db.join",
map=respondents_points,
column="id_INT",
other_table=clc_summary_table,
other_column="resp_id")




cheers, Johannes

On Wed, Jun 14, 2017 at 11:10 AM, Moritz Lennert <
mlennert at club.worldonline.be> wrote:

> Hi Johannes,
>
> On 07/06/17 10:20, Johannes Radinger wrote:
>
>> Thank you Moritz,
>>
>> your suggestion using r.stats and the rasterized layer 2 of buffers
>> works really nice. It took me just a while to summarize and join all
>> data and get them back into the db in the right format. However, I think
>> I managed this task now. Thank you for your help!
>>
>
> Would you be willing to share the final version of your approach ? This
> might be a nice seed for a module that would offer this function.
>
> Moritz
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/grass-user/attachments/20170616/bb2dcc08/attachment-0001.html>