<jupyter_start><jupyter_text>Chicago dataThe data can be downloaded from https://catalog.data.gov/dataset/crimes-2001-to-present-398a4 (see the module docstring of `open_cp.sources.chicago`  See also https://data.cityofchicago.org/Public-Safety/Crimes-2001-to-present/ijzp-q8t2In this notebook, we quickly look at the data, check that the data agrees between both sources, and demo some of the library features provided for loading the data.<jupyter_code>import open_cp.sources.chicago as chicago
import geopandas as gpd

import sys, os, csv, lzma
filename = os.path.join("..", "..", "open_cp", "sources", "chicago.csv")
filename_all = os.path.join("..", "..", "open_cp", "sources", "chicago_all.csv.xz")
filename_all1 = os.path.join("..", "..", "open_cp", "sources", "chicago_all1.csv.xz")<jupyter_output><empty_output><jupyter_text>Let us look at the snapshot of the last year, vs the total dataset.  The data appears to be the same, though the exact format changes.<jupyter_code>with open(filename, "rt") as file:
    reader = csv.reader(file)
    print(next(reader))
    print(next(reader))
with lzma.open(filename_all, "rt") as file:
    reader = csv.reader(file)
    print(next(reader))
    print(next(reader))<jupyter_output>['ID', 'Case Number', 'Date', 'Block', 'IUCR', 'Primary Type', 'Description', 'Location Description', 'Arrest', 'Domestic', 'Beat', 'District', 'Ward', 'Community Area', 'FBI Code', 'X Coordinate', 'Y Coordinate', 'Year', 'Updated On', 'Latitude', 'Longitude', 'Location']
['8651563', 'HV322174', '06/05/2012 11:00:00 AM', '022XX N CANNON DR', '0810', 'THEFT', 'OVER $500', 'STREET', 'false', 'false', '1814', '018', '43', '7', '06', '1175057', '1915111', '2012', '02/04/2016 06:33:39 AM', '41.922450893', '-87.632206293', '(41.922450893, -87.632206293)']<jupyter_text>As well as loading data directly into a `TimedPoints` class, we can process a sub-set of the data to GeoJSON, or straight to a geopandas dataframe (if geopandas is installed).<jupyter_code>geo_data = chicago.load_to_GeoJSON()
geo_data[0]
frame = chicago.load_to_geoDataFrame()
frame.head()<jupyter_output><empty_output><jupyter_text>Explore with QGISWe can save the dataframe to a shape-file which can be viewed in e.g. QGIS.To explore the spatial-distribution, I would recommend using an interactive GIS package.  Using QGIS (free and open source) you can easily add a basemap using GoogleMaps or OpenStreetMap, etc.  See http://maps.cga.harvard.edu/qgis/wkshop/basemap.phpI found this to be slightly buggy.  On Windows, QGIS 2.18.7 I found that the following worked:- First open the `chicago.shp` file produced from the line above.- Select the Coordinate reference system "WGS 84 / EPSG:4326"- Now go to the menu "Web" -> "OpenLayers plugin" -> Whatever- The projection should change to EPSG:3857.  The basemap will obscure the point map, so in the "Layers Panel" drag the basemap to the bottom.- Selecting EPSG:3857 at import time doesn't seem to work (which is different to the instructions..!)<jupyter_code># On my Windows install, if I don't do this, I get a GDAL error in
# the Jupyter console, and the resulting ".prj" file is empty.
# This isn't critical, but it confuses QGIS, and you end up having to
# choose a projection when loading the shape-file.
import os
os.environ["GDAL_DATA"] = "C:\\Users\\Matthew\\Anaconda3\\Library\\share\\gdal\\"

frame.to_file("chicago")<jupyter_output><empty_output><jupyter_text>A geoPandas exampleLet's use the "generator of GeoJSON" option shown above to pick out only BURGLARY crimes from the 2001-- dataset (which is too large to easily load into a dataframe in one go).<jupyter_code>with lzma.open(filename_all, "rt") as file:
    features = [ event for event in chicago.generate_GeoJSON_Features(file, type="all")
                if event["properties"]["crime"] == "THEFT" ]
    
frame = gpd.GeoDataFrame.from_features(features)
frame.crs = {"init":"EPSG:4326"} # Lon/Lat native coords
frame.head()
frame.to_file("chicago_all_theft")
with lzma.open(filename_all, "rt") as file:
    features = [ event for event in chicago.generate_GeoJSON_Features(file, type="all")
                if event["properties"]["crime"] == "BURGLARY" ]

frame = gpd.GeoDataFrame.from_features(features)
frame.crs = {"init":"EPSG:4326"} # Lon/Lat native coords
frame.head()
frame.to_file("chicago_all_burglary")
frame["type"].unique()
frame["location"].unique()<jupyter_output><empty_output><jupyter_text>Upon loading into QGIS to visualise, we find that the 2001 data seems to be geocoded in a different way...  The events are not on the road, and the distribution looks less artificial.  Let's extract the 2001 burglary data, and then the all the 2001 data, and save.<jupyter_code>with lzma.open(filename_all, "rt") as file:
    features = [ event for event in chicago.generate_GeoJSON_Features(file, type="all")
                if event["properties"]["timestamp"].startswith("2001") ]

frame = gpd.GeoDataFrame.from_features(features)
frame.crs = {"init":"EPSG:4326"} # Lon/Lat native coords
frame.head()
frame.to_file("chicago_2001")<jupyter_output><empty_output><jupyter_text>Explore rounding errorsWe check the following:- The X and Y COORDINATES fields (which we'll see, in a different notebook, at longitude / latitude coordinates projected in EPSG:3435 in feet) are always whole numbers.- The longitude and latitude data contains at most 9 decimals places of accuracy.In the other notebook, we look at map projections.  The data is most consistent with the longitude / latitude coordinates being the primary source, and the X/Y projected coordinates being computed and rounded to the nearest integer.<jupyter_code>longs, lats = [], []
xcs, ycs = [], []

with open(filename, "rt") as file:
    reader = csv.reader(file)
    header = next(reader)
    print(header)
    for row in reader:
        if len(row[14]) > 0:
            longs.append(row[14])
            lats.append(row[15])
            xcs.append(row[12])
            ycs.append(row[13])
set(len(x) for x in longs), set(len(x) for x in lats)
any(x.find('.') >= 0 for x in xcs), any(y.find('.') >= 0 for y in ycs)<jupyter_output><empty_output><jupyter_text>Repeated dataMostly the "case" assignment is unique, but there are a few exceptions to this.<jupyter_code>import collections

with lzma.open(filename_all, "rt") as file:
    c = collections.Counter( event["properties"]["case"] for event in
                            chicago.generate_GeoJSON_Features(file, type="all") )
    
multiples = set( key for key in c if c[key] > 1 )
len(multiples)
with lzma.open(file_all, "rt") as file:
    data = gpd.GeoDataFrame.from_features(
        event for event in chicago.generate_GeoJSON_Features(file, type="all")
        if event["properties"]["case"] in multiples
        )
    
len(data), len(data.case.uniques())<jupyter_output><empty_output>