{"cells":[{"metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true},"cell_type":"code","source":"# This Python 3 environment comes with many helpful analytics libraries installed\n# It is defined by the kaggle/python Docker image: https://github.com/kaggle/docker-python\n# For example, here's several helpful packages to load\n\nimport numpy as np # linear algebra\nimport pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)\npd.options.display.max_rows=100\n\n# Input data files are available in the read-only \"../input/\" directory\n# For example, running this (by clicking run or pressing Shift+Enter) will list all files under the input directory\n\nimport os\nfor dirname, _, filenames in os.walk('/kaggle/input'):\n    for filename in filenames:\n        print(os.path.join(dirname, filename))\n\n# You can write up to 20GB to the current directory (/kaggle/working/) that gets preserved as output when you create a version using \"Save & Run All\" \n# You can also write temporary files to /kaggle/temp/, but they won't be saved outside of the current session","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Dataset description and relation to competition\nThis dataset is available here: http://citycarbonfootprints.info/ and derives from a recent publication on carbon footprints of world cities: https://iopscience.iop.org/article/10.1088/1748-9326/aac72a \"Carbon footprints of 13 000 cities\", published in 2018. The data appear to come from earlier than 2018, although it seems much effort went in to producing the output, which includes carbon footprints not only of 13,000 cities around the world, but also carbon footprints on a 250m spanning the globe.\n\nThe study appears to be one of the most extensive efforts yet to characterize carbon footprints with as much spatial coverage as possible. This opens up numerous possibilities for analysis and incorporation into KPIs for the CDP: Unlocking Climate Solutions competition. The global extent should enable actual carbon footprint estimates to be spatially joined to the competition data, which include:\n- City-level information about commitments to improving carbon budgets, as well as\n- Fine-grained spatial analysis within cities at the zip code and census tract level\n\nThe data presented here could be used for both large, multi-city analysis, as well as within-city analysis due to the 250m spatial resolution.\n\nHere I load the data and visualize for Los Angeles County in California, USA."},{"metadata":{},"cell_type":"markdown","source":"### Load Carbon Footprint Data"},{"metadata":{"_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","trusted":true},"cell_type":"code","source":"import rasterio\nimport geopandas as gpd\nco2_df = rasterio.open('../input/global-gridded-model-of-carbon-footprints-ggmcf/GGMCF_v1.0.tif')\nco2_df.crs.wkt","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"## Cities polygons\ncities_poly_df = gpd.read_file('/kaggle/input/geospatial-environmental-and-socioeconomic-data/1_CITIES_landscan/ne_10m_urban_areas_landscan/ne_10m_urban_areas_landscan.shp')\ncities_poly_df.crs","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# create an output directory\n!mkdir /kaggle/working/global-gridded-model-of-carbon-footprints-ggmcf","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"%%time\n\n# The coordinate reference system of these data  need to match other systems, to be spatially joined.\nfrom rasterio.warp import calculate_default_transform, reproject, Resampling\n\n# reprojection of Carbon Footprint data\ndst_crs = 'EPSG:4326'\n\nwith rasterio.open('../input/global-gridded-model-of-carbon-footprints-ggmcf/GGMCF_v1.0.tif') as src:\n    transform, width, height = calculate_default_transform(\n        src.crs, dst_crs, src.width, src.height, *src.bounds)\n    kwargs = src.meta.copy()\n    kwargs.update({\n        'crs': dst_crs,\n        'transform': transform,\n        'width': width,\n        'height': height\n    })\n\n    with rasterio.open('/kaggle/working/global-gridded-model-of-carbon-footprints-ggmcf/GGMCF_v1.0.EPSG4326.tif', 'w', **kwargs) as dst:\n        for i in range(1, src.count + 1):\n            reproject(\n                source=rasterio.band(src, i),\n                destination=rasterio.band(dst, i),\n                src_transform=src.transform,\n                src_crs=src.crs,\n                dst_transform=transform,\n                dst_crs=dst_crs,\n                resampling=Resampling.nearest)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Now we can combine exemplary city mask with carbon footprint"},{"metadata":{"trusted":true},"cell_type":"code","source":"from rasterio.mask import mask\n\ncity_name = 'Warsaw'\n\nwith rasterio.open('/kaggle/working/global-gridded-model-of-carbon-footprints-ggmcf/GGMCF_v1.0.tif') as src:\n    carbon_data, carbon_transform = mask(src, cities_poly_df.loc[cities_poly_df['name_conve']==city_name]['geometry'], crop=True)\n    carbon_meta = src.meta\n\n# let's visualize this:\nimport matplotlib.pyplot as plt\nfrom rasterio.plot import show\n\nfig, axs = plt.subplots(1,2)\ncities_poly_df.loc[cities_poly_df['name_conve']==city_name].plot(ax=axs[0])\nshow(carbon_data, ax=axs[1])","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## What could we calculate?"},{"metadata":{"trusted":true},"cell_type":"code","source":"print(f'Mean co2 footprint of {city_name}: {carbon_data.mean()}')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Combine with CDP's list of cities"},{"metadata":{"trusted":true},"cell_type":"code","source":"import geopandas as gpd\nimport pandas as pd\nimport numpy as np\npd.options.display.max_colwidth=None\npd.options.display.max_rows=200","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"PREFIX = '../input/cdp-unlocking-climate-solutions/Cities/Cities Disclosing/'\n\ncdp_cities_disclosing_df = pd.concat([pd.read_csv(PREFIX+'2018_Cities_Disclosing_to_CDP.csv'),\\\n                                           pd.read_csv(PREFIX+'2019_Cities_Disclosing_to_CDP.csv'),\\\n                                           pd.read_csv(PREFIX+'2020_Cities_Disclosing_to_CDP.csv')])\ncdp_cities_disclosing_df.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"!pip install geotext","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"from geotext import GeoText","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# get valid city name from 'Organization' and then 'City' field of the City Disclosure Form, with the help of GeoText library\n\ncdp_cities_disclosing_df['Organization_clean'] = cdp_cities_disclosing_df.apply(lambda x: GeoText(str(x['Organization'])).cities, axis=1)\\\n                                                        .apply(lambda x: ','.join(x))\\\n                                                        .replace('', np.nan)\n\nprint('cdp_cities_disclosing_df[\\'Organization_clean\\']: ', cdp_cities_disclosing_df['Organization_clean'].dropna().shape)\n\n# check more than one city decoded anywhere in one field\ncdp_cities_disclosing_df.loc[cdp_cities_disclosing_df['Organization_clean'].fillna('').str.split(',').apply(len)>1]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"df_unique_cities = cdp_cities_disclosing_df[['City', 'Country', 'City Location']].drop_duplicates().reset_index(drop=True)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"df_unique_cities.count()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"cdp_cities_disclosing_df['City_clean'] = cdp_cities_disclosing_df.apply(lambda x: GeoText(str(x['City'])).cities, axis=1)\\\n.apply(lambda x: ','.join(x))\\\n.replace('', np.nan)\n# check if over one city decoded anywhere in one field\nprint('cdp_cities_disclosing_df[\\'City_clean\\']: ', cdp_cities_disclosing_df['City_clean'].dropna().unique().shape)\ncdp_cities_disclosing_df.loc[cdp_cities_disclosing_df['City_clean'].fillna('').str.split(',').apply(len)>1]\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"cdp_cities_disclosing_df.loc[cdp_cities_disclosing_df['City_clean']=='Yaoundé']","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"cdp_cities_disclosing_df['City_clean'].value_counts()[:100]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"cdp_cities_disclosing_df['City'].value_counts()[:100]#isna().sum()#.drop_duplicates()#.dropna()#.value_counts()#.dropna().unique().shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"cdp_cities_disclosing_df['City_clean'] = cdp_cities_disclosing_df[['Organization', 'City', 'Organization_clean', 'City_clean']]\\\n                                                .ffill(axis=1)['City_clean']\nprint('cdp_cities_disclosing_df[\\'City_clean\\']: ', cdp_cities_disclosing_df['City_clean'].dropna().unique().shape)\ncdp_cities_disclosing_df.loc[cdp_cities_disclosing_df['City_clean'].fillna('').str.split(',').apply(len)>1].head()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Country"},{"metadata":{"trusted":true},"cell_type":"code","source":"GeoText('Orange County').cities","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"cdp_cities_disclosing_df['Country_clean'] = cdp_cities_disclosing_df['Country']\\\n                                                    .apply(lambda x: GeoText(str(x)).countries)\\\n                                                    .apply(lambda x: ','.join(x))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"## so, if we use the same method for Countries, will that be unambiguous? (there should be zero rows below):\ncdp_cities_disclosing_df.loc[cdp_cities_disclosing_df['Country']\\\n                                  .apply(lambda x: GeoText(str(x)).countries)\\\n                                  .apply(lambda x: len(x))>1]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"cdp_cities_disclosing_df[['Country', 'City', 'Organization', 'Country_clean', 'City_clean']]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"cdp_cities_disclosing_df['City_clean'].dropna().shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"cdp_cities_disclosing_df[['Country_clean', 'City_clean']]\\\n.drop_duplicates()\\\n.sort_values(by=['Country_clean', 'City_clean'])\\\n.shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"## as temp df\n###'Country_clean', 'City_clean', ]]\nunique_cities_disclosing_df = cdp_cities_disclosing_df[['Account Number', 'City', 'Country', 'City Location']]\\\n.drop_duplicates()\\\n.reset_index(drop=True)\nunique_cities_disclosing_df.shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"cdp_cities_disclosing_df.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"(~unique_cities_disclosing_df['City Location'].isnull()).sum()#.dropna() # 687 has geographical point assigned","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Now get the city name from map by Point Coordinates"},{"metadata":{"trusted":true},"cell_type":"code","source":"## as temp df\n###'Country_clean', 'City_clean', ]]\nunique_cities_disclosing_df = cdp_cities_disclosing_df[['Account Number', 'City', 'Country', 'City Location']]\\\n.drop_duplicates()\\\n.reset_index(drop=True)\nunique_cities_disclosing_df.shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"unique_cities_disclosing_df","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"## Cities polygons\ncities_poly_df = gpd.read_file('../input/1-cities-landscan/1_CITIES_landscan/ne_10m_urban_areas_landscan/ne_10m_urban_areas_landscan.shp')\ncities_poly_df.crs","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"import matplotlib.pyplot as plt\n\n# Plot the WGS84\ncities_poly_df.loc[cities_poly_df['name_conve']=='Warsaw'].plot(facecolor='gray');\n# Add title\nplt.title(\"WGS84 (EPSG:4326) projection\");\n# Remove empty white space around the plot\nplt.tight_layout()\n\n# Plot the one with ETRS-LAEA projection (more realistic for this region)\ncities_poly_df.loc[cities_poly_df['name_conve']=='Warsaw'].to_crs(epsg=3035).plot(facecolor='blue');\n# Add title\nplt.title(\"ETRS Lambert Azimuthal Equal Area projection\");\n# Remove empty white space around the plot\nplt.tight_layout()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# We use shapely.wkt sub-module to parse wkt format:\nfrom shapely import wkt\n\n# for WKT to convert correctly, in columns must be valid strings. Execute only once!\nunique_cities_disclosing_df.loc[~unique_cities_disclosing_df['City Location'].isnull(), 'City Location'] = unique_cities_disclosing_df.loc[~unique_cities_disclosing_df['City Location'].isnull(), 'City Location'].apply(wkt.loads)\n\n# cast Pandas DF into Geo Pandas DF\ngdf = gpd.GeoDataFrame(unique_cities_disclosing_df, geometry='City Location', crs=\"EPSG:4326\")","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"poly_df = cities_poly_df.reset_index(drop=True).copy()     #.loc[cities_poly_df['name_conve'].isin(['Warsaw', 'Moscow'])].reset_index(drop=True)\npoly_df.shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"points_df = gdf.reset_index(drop=True).copy()     #.loc[gdf['Country']=='Poland'].reset_index(drop=True)\npoints_df.shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# fig, ax = plt.subplots()\n# poly_df.plot(ax=ax, facecolor='red');\n# points_df.plot(ax=ax, color='blue', markersize=5);\n# plt.tight_layout();","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"assert points_df.crs == poly_df.crs, 'Coordinate Systems do not match!!'","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"points_df.head()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Spatial Join - joining cities with polygons"},{"metadata":{"trusted":true},"cell_type":"code","source":"print('We have coordinates for the following number of cities:')\nprint(unique_cities_disclosing_df.loc[unique_cities_disclosing_df['City Location'].notnull()].shape[0])\n\nprint('We don\\'t have coordinates for the following number of cities:')\nprint(unique_cities_disclosing_df.loc[unique_cities_disclosing_df['City Location'].isnull()].shape[0])","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"## Assigning polygons:\n\nprint(f'Initially we have {points_df.shape} unique cities')\n# this spatial join works well here, but two problems: some polygons in our dataset are wrong like (Yokohama inside Tokyo). \nresult = gpd.sjoin(points_df, poly_df, how=\"left\", op=\"within\")\nprint(f'after left joining: {result.shape}')\n# We lost the geometry so now, we just have to add it back again.\nresult = result.merge(poly_df[['name_conve', 'geometry']], left_on='name_conve', right_on='name_conve', how='left')\nprint(f'and then: {result.shape}')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# But problem: city Tokyo was found in Tokio and in Yokohama polygon \nresult.loc[result['Account Number']==31111].iloc[:,:-1]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# The problem is, that more cities have their polgyons overlapping over each other like Tokyo and Yokohama.\n# I'll assign the polygon with the smaller area in that cases: using 'min_areakm' column (so Yokohama in this case)\n# when duplicated rows (by keys) that should keep only the first, smaller polygon (because sorted ascending)\nresult = result.sort_values(by=['min_areakm']).drop_duplicates(subset=['Account Number', 'City', 'Country', 'City Location'])\n\nprint('Shape: ', result.shape)\n\n# check if Yokohama.\nresult.loc[result['Account Number']==31111].iloc[:,:-3]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"## how many points don't have the pologyon assigned\nprint(result.loc[result['City Location'].notnull() & (result['geometry'].isnull())].shape)\nresult.loc[result['City Location'].notnull() & (result['geometry'].isnull())]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# so we have polygons only for \nprint(f'we have polygons for:', result.loc[(result['geometry'].notnull())].shape[0])","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# test\nfig, ax = plt.subplots()\nresult.loc[result['name_conve']=='Warsaw', 'geometry'].plot(ax=ax, facecolor='gray');\nresult.loc[result['name_conve']=='Warsaw', 'City Location'].plot(ax=ax, facecolor='red');\nplt.tight_layout();","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Saving results"},{"metadata":{"trusted":true},"cell_type":"code","source":"# create an output directory\n!mkdir /kaggle/working/CDP","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"## Saving results\nresult.to_pickle('/kaggle/working/CDP/CDP_cities_with_polygons.pkl')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Geocoding - let's try to geocode more points (the rest of the cities)\nBecause for many cities, coordinates are outside cities"},{"metadata":{"trusted":true},"cell_type":"code","source":"## not null geometries:\nresult.geometry.notnull().sum()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"!pip install geopy","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"import geopy\n\n# create an address if City is not empty (because that would be misleading)\nresult['addr'] = result[['City', 'Country']].fillna('').apply(lambda x: str(x['City']) + ',' + str(x['Country']), axis=1)\nresult.loc[result['City'].isna(), 'addr'] = np.nan\n\n# let's geocode\nfrom geopy.geocoders import Nominatim\ngeolocator = Nominatim(user_agent=\"CDP Kaggle challenge\")\nresult['addr_geocoded'] = None\n## TODO: doesn't work for greater dataset\nresult.loc[result['addr'].notnull(), 'addr_geocoded'] = geocode(result.loc[result['addr'].notnull(), 'addr'][:10])['geometry']\n\n# cast to geometry type\nresult['addr_geocoded'] = gpd.GeoSeries(result['addr_geocoded'], crs='EPSG:4326')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# same Coordinate System?\nresult['City Location'].crs == result['addr_geocoded'].crs\n\n# with pyproj\n# \n# Points are in a lon, lat coordinate system (EPSG:4326 or WGS 84). To calculate a distance in meters, \n# we would need to either use the Great-circle distance or project them in a local coordinate system \n# to approximate the distance with a good precision.\n# compare with http://www.csgnetwork.com/lldistcalc.html\nimport pyproj\ngeod = pyproj.Geod(ellps='WGS84')\ndist_zip = result.loc[(result['addr_geocoded'].notnull()) & (result['City Location'].notnull())]\\\n            .apply(lambda x: geod.inv(x['addr_geocoded'].x, x['addr_geocoded'].y, x['City Location'].x, x['City Location'].y), axis=1)\ndist_zip = dist_zip.apply(pd.Series)\ndist_zip.columns=['angle1','angle2','distance']\ndist_zip.index.name=None\nresult['geocoder_dist_km'] = np.NaN\nresult.loc[(result['addr_geocoded'].notnull()) & (result['City Location'].notnull()), 'geocoder_dist_km'] = (dist_zip['distance']/1000).round(2)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"\nresult['geocoder_dist_km'][:10]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"result.loc[result['addr_geocoded'].notnull(), ['addr_geocoded', 'addr', 'City Location', 'name_conve', 'geocoder_dist_km']]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"result.shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"result.loc[result['geocoded'].notnull()].shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# create an output directory\n!mkdir /kaggle/working/CDP/","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"type(result)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"result.to_pickle('/kaggle/working/CDP/CDP_cities_with_polygons.pkl')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"result_x = pd.read_pickle('/kaggle/working/CDP/CDP_cities_with_polygons.pkl')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"type(result_x)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"## Writing results to Shapefile:\n#result.to_file('/kaggle/working/CDP/CDP_cities_with_polygons.shp')\nresult.to_file(\"countries.geojson\", driver='GeoJSON')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"temp = pd.read_pickle('CDP_cities_with_polygons.pkl')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"temp.columns","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Assign values from different dataset to what we've got (688 for 1032 overall)"},{"metadata":{"trusted":true},"cell_type":"code","source":"# how many points we had originally:\nresult.loc[result['City Location'].notnull()].shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"result.loc[result['geometry'].notnull()].shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"cdp_cities_disclosing_df.loc[cdp_cities_disclosing_df['City']=='Warsaw']","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"## End!","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"fig, ax = plt.subplots()\nresult.loc[result['name_conve']=='Warsaw', 'geometry'].plot(ax=ax, facecolor='gray');\nresult.loc[result['name_conve']=='Warsaw', 'City Location'].plot(ax=ax, facecolor='red');\nplt.tight_layout();","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"\nfig, ax = plt.subplots()\n\nresult.loc[result['name_conve']=='Tokyo', 'geometry'].plot(ax=ax, facecolor='yellow');\nresult.loc[result['Account Number']==31111, 'City Location'].plot(ax=ax, facecolor='red');\nplt.tight_layout();","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"\nfig, ax = plt.subplots()\n\nresult.loc[result['name_conve']=='Yokohama', 'geometry'].plot(ax=ax, facecolor='blue');\nresult.loc[result['Account Number']==31111, 'City Location'].plot(ax=ax, facecolor='red');\nplt.tight_layout();","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"result.loc[result['name_conve']=='Yokohama', 'geometry'].plot()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"import maup\nmaup.resolve_overlaps(poly_df['geometry'].values)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"poly_df['geometry'].values()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"points_df.shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"result.shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"city = points_df.iloc[0]\ncity.crs","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"for city in points_df:\n    for poly in poly_df:\n        city.intersects(poly.unary_union)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"points_df.apply(lambda x: x['City Location'].intersects(poly_df.unary_union))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"points_df[:1].intersects(poly_df.unary_union) ","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"gdf[:1].intersects(cities_poly_df.unary_union)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"gdf[:1].within(cities_poly_df.loc[:10, 'geometry'])","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"cities_poly_df","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Now we can merge by account with Cities Responses"},{"metadata":{"trusted":true},"cell_type":"code","source":"## TODO:....   \n\n\n# import cities response df\ncdp_full_cities_df = pd.read_csv(\"../input/cdp-unlocking-climate-solutions/Cities/Cities Responses/2020_Full_Cities_Dataset.csv\")\n#cdp_full_cities_df.head()\n\ncdp_full_cities_df[['Country', 'Organization']].replace('City of ', '', regex=True)\\\n       .replace('City of ', '', regex=True)\\\n       .replace('Township of ', '', regex=True)\\\n       .drop_duplicates()\\\n       .sort_values(by=['Country', 'Organization'])[:100]","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Combine with raster data for cities"},{"metadata":{"trusted":true},"cell_type":"code","source":"unique_cities_disclosing_df.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"cities_poly_df.shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"cities_poly_df.iloc[:, :-1].head(3).T","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"cities_poly_df['City_clean'] = cities_poly_df['name_conve']\\\n                                .apply(lambda x: GeoText(x).cities)\\\n                                .apply(lambda x: ','.join(x))\\\n                                .replace('', np.nan)\nprint('GeoTexted: ', cities_poly_df['City_clean'].dropna().shape[0], ' out of ', cities_poly_df.shape[0])\n\n# the rest (names not found in GeoText) will be filled with the original names\ncities_poly_df['City_clean'] = cities_poly_df[['name_conve', 'City_clean']].ffill(axis=1)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"unique_cities_disclosing_df.shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# assign polygons to cities\ndf_uniq_cities_polygons_comb = unique_cities_disclosing_df.merge(cities_poly_df, left_on='City_clean', right_on='City_clean', how='left')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"df_uniq_cities_polygons_comb[['geometry']].dropna().shape   ## 383 ","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"## Now we can print on map and calculate Co2 Footprint for those cities ~380 cities .... only... \nfrom rasterio.mask import mask\nimport matplotlib.pyplot as plt\nfrom rasterio.plot import show\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"#df_uniq_cities_polygons_comb['co2_footprint'] = np.nan\n\ni = 0 \nwith rasterio.open('/kaggle/working/global-gridded-model-of-carbon-footprints-ggmcf/GGMCF_v1.0.tif') as src:\n    for city_name in df_uniq_cities_polygons_comb.dropna(subset=['geometry'])['City_clean'].tolist()[:10]:\n        carbon_data, carbon_transform = mask(src, df_uniq_cities_polygons_comb.loc[df_uniq_cities_polygons_comb['City_clean']==city_name]['geometry'], crop=True)\n        carbon_meta = src.meta\n        \n        print(city_name, 'CO2 footprint: ', carbon_data.mean())\n        \n        df_uniq_cities_polygons_comb.loc[df_uniq_cities_polygons_comb['City_clean']==city_name, 'co2_footprint'] = carbon_data.mean()\n\n        #while i < 10:  # show first 10 polygons\n        fig, axs = plt.subplots(1,2)      # (1,2)\n        #df_uniq_cities_polygons_comb.loc[df_uniq_cities_polygons_comb['City_clean']==city_name].plot(ax=axs[0])\n        show(carbon_data, ax=axs[1], title=city_name)\n        i += 1","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"cities_poly_df.loc[cities_poly_df['name_conve']=='Warsaw']['geometry']","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"cities_poly_df.","execution_count":null,"outputs":[]}],"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"pygments_lexer":"ipython3","nbconvert_exporter":"python","version":"3.6.4","file_extension":".py","codemirror_mode":{"name":"ipython","version":3},"name":"python","mimetype":"text/x-python"}},"nbformat":4,"nbformat_minor":4}