world-ecoregion/data.py
2019-02-03 09:04:28 +03:30

75 lines
2.6 KiB
Python

import geopandas
import os
import rasterio
import pandas as pd
from matplotlib import pyplot
from shapely.geometry import Point
directory = os.path.dirname(os.path.abspath(__file__))
GEODATA = os.path.join(directory, 'geodata')
ECOREGIONS = os.path.join(GEODATA, 'ecoregions', 'Ecoregions2017.shp')
ELEVATION = os.path.join(GEODATA, 'srtm', 'topo30-180.tif')
TEMP = os.path.join(GEODATA, 'air_temp')
PRECIP = os.path.join(GEODATA, 'precipitation')
def read_temp_data(year):
return pd.read_csv(os.path.join(TEMP, 'air_temp.{}'.format(year)), sep='\s+', header=None,
names=['longitude', 'latitude', 'january',
'february', 'march', 'april',
'may', 'june', 'july', 'august',
'september', 'november', 'october',
'december', 'yearly_avg'])
def read_precip_data(year):
return pd.read_csv(os.path.join(PRECIP, 'precip.{}'.format(year)), sep='\s+', header=None,
names=['longitude', 'latitude', 'january',
'february', 'march', 'april',
'may', 'june', 'july', 'august',
'september', 'november', 'october',
'december', 'yearly_avg'])
eco = geopandas.read_file(ECOREGIONS)
elevation = rasterio.open(ELEVATION)
# world = geopandas.read_file(geopandas.datasets.get_path('naturalearth_lowres'))
# world['geometry'] = world['geometry'].unary_union
# world = geopandas.read_file(geopandas.datasets.get_path('naturalearth_lowres'))
# print(world.head())
# world = world[['continent', 'geometry']]
# continents = world.dissolve(level=1)
# continents.plot();
# print(eco.head())
# print(elevation)
# eco['geometry'].boundary.plot()
# eco.dissolve()
# eco.plot()
# # rasterio.plot.show(src)
# # pyplot.imshow(elevation.read(1))
# 51.42
# 35.69
# tehran = eco.geometry.contains()
def get_point_information(longitude, latitude):
p = Point(longitude, latitude)
ecoregion = eco.loc[lambda c: c.geometry.contains(p)]
return {
'biome_num': ecoregion.loc['BIOME_NUM'].iloc[0],
'biome_name': ecoregion.loc['BIOME_NAME'].iloc[0],
}
import time
start_time = time.time()
print('Before call')
print('Tehran', get_point_information(51.42, 35.69))
print("--- %s seconds ---" % (time.time() - start_time))
start_time = time.time()
print('Amazon', get_point_information(-59.78, -5.5))
print("--- %s seconds ---" % (time.time() - start_time))
# print(eco.geometry)
# print(tehran.distance(world.boundary))
# world.boundary.plot()
pyplot.show()