[python #DEM-876252]: Lightning Density

Greetings!

Assuming you have already read your data into Python, the following code 
demonstrates some methods of plotting the data:

    import cartopy.crs as ccrs
    import matplotlib.pyplot as plt
    import numpy as np
    import scipy.stats as ss

    # Random point locations
    lons = -80 + 5 * np.random.randn(3000)
    lats = 40 + 3 * np.random.randn(3000)

    # Perform Gaussian KDE estimate on random data
    kernel = ss.gaussian_kde(np.vstack([lons, lats]))

    # Create a grid of lon/lat values for calculating our density and contouring
    grid_lon = np.linspace(-100, -60, 300)
    grid_lat = np.linspace(20, 60, 299)
    X, Y = np.meshgrid(grid_lon, grid_lat)
    vals = kernel(np.vstack([X.ravel(), Y.ravel()]))

    # Plot raw points
    fig = plt.figure(figsize=(10, 20))
    ax = fig.add_subplot(1, 3, 1, projection=ccrs.PlateCarree())
    ax.plot(lons, lats, 'ro')
    ax.coastlines(zorder=100)

    # Draw a hexbin plot of density
    ax2 = fig.add_subplot(1, 3, 2, projection=ccrs.PlateCarree())
    hb = ax2.hexbin(lons, lats, gridsize=20, mincnt=1)
    ax2.coastlines()

    # Plot contours of our Kernel Density Estimate
    ax3 = fig.add_subplot(1, 3, 3, projection=ccrs.PlateCarree())
    hb = ax3.contour(grid_lon, grid_lat, vals.reshape(X.shape))
    ax3.coastlines()
    ax3.set_extent(ax2.get_extent())

So one method, if you just want to plot the density, is to use `hexbin`, which 
automatically bins the data to hexagons in 2 dimensions, and counts the data 
points that fall in each.

If you want to calculate the density, the only method I'm aware of is Kernel 
Density Estimation, which is built into SciPy:

  
https://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.gaussian_kde.html#scipy.stats.gaussian_kde

The code above demonstrates using this and plotting contours of the estimated 
density.

Hope this helps,

Ryan

> Sir,
> 
> We are trying to read the data from a text file.  We would like to pull all
> of the lat long data from the text file and map it over a map of NY state
> as a density map. We would like to do this in Python if possible.
> Example of a line of the text file:   1996-01-17 03:54:35.853  44.9628
> -78.9399   -37.9
> 
> If you need any more information please let me know.
> 
> Very Respectfully
> Dyllon Main
> 


Ticket Details
===================
Ticket ID: DEM-876252
Department: Support Python
Priority: Low
Status: Closed
===================
NOTE: All email exchanges with Unidata User Support are recorded in the Unidata 
inquiry tracking system and then made publicly available through the web.  If 
you do not want to have your interactions made available in this way, you must 
let us know in each email you send to us.