# Density per grid cell

I have a shapefile with the following format. These are storm tracks and each unique "Serial" represents one storm track.

``````import matplotlib.pyplot as plt
import geopandas as gpd

input = "lines.shp"
print(linesplot)
``````
``````Serial  Year    Month   Day  Hour   Lat     Lon         geometry
2003002 2003    2       12   0      6.98    125.42      LINESTRING (124.61000 7.87000, 125.42000 6.98000)
2003002 2003    2       12   6      7.87    124.61      LINESTRING (122.98000 8.60000, 124.61000 7.87000
...     ...     ...     ...  ...    ...     ...         ...
2003155 2003    12      23   0      7.48    141.56      LINESTRING (140.62000 7.93000, 141.56000 7.48000)
2003155 2003    12      23   6      7.93    140.62      LINESTRING (138.91000 8.33000, 140.62000 7.93000)
``````

I would like to create a track density map something similar to this figure here. In the density map, I would like to count the number of storm tracks passing through within each 1x1 degree grid cell over the entire domain.

--- Update ---

As suggested, I tried using a spatial join. I have also found this post where he did it with points.

I first created 1x1 grid cells

``````# create the cells in a loop
grid_cells = []
for x0 in np.arange(xmin, xmax+cell_size, cell_size ):
for y0 in np.arange(ymin, ymax+cell_size, cell_size):
# bounds
x1 = x0-cell_size
y1 = y0+cell_size
grid_cells.append( shapely.geometry.box(x0, y0, x1, y1)  )
cell = geopandas.GeoDataFrame(grid_cells, columns=['geometry'],
crs=crs)
``````

and used `merged = geopandas.sjoin(linesplot, cell, how='left', op='within') ` to merge the dataframes

as suggested in the post,

``````# make a simple count variable that we can sum
merged['lines']=1
# Compute stats per grid cell -- aggregate lines to grid cells with dissolve
dissolve = merged.dissolve(by="index_right", aggfunc="count")
# put this into cell
cell.loc[dissolve.index, 'lines']=dissolve.lines.values
del dissolve
``````

However, when plotting, I don't seem to get the right number of counts for each grid cell.

``````ax = cell.plot(column='lines', vmin=0, vmax=10, edgecolor="none")

``````

the highlighted grid cells should have the same color as they have the same number of line count • Wouldn't it be better to test the length of a given storm path within a given cell, sort of like the way anti-aliasing works. This post might be interesting stackoverflow.com/questions/36231289/… Jul 28, 2020 at 5:59
• Intersect or spatial join storms to grid then Count unique storm ids per grid.
– BERA
Jul 28, 2020 at 6:01

Try using Spatial Join like this:

``````import geopandas as gpd

storms = r"C:\GIS\data\jl_riks.shp"
storm_id_field = "stormid"
grid = r"C:\GIS\data\rutnat_no_duplicates.shp"
grid_id_field = "RUTA" 