Mapping

Show Only ...
Maps - Photos - Videos

Another Example of GeoPandas and Python Spatial Joins

Another Example of GeoPandas and Python Spatial Joins πŸ—Ί

Yesterday, I posted a much more complicated piece of code that pulled addresses from the SAM (State Address Management) database and did a spatial join to add a column to the file with Assembly District and Municipality. This was a bit too complex, so I made a simpler one for other purposes that doesn’t require the coordinates to be obtained from SAM.

This python script takes two parameters:

  • Path to a CSV file that contains an X and Y coordinate
  • Path to a Shapefile or Geopackage to Join Against

Then the code will create a new CSV file with the spatially joined attributes pulled from the Shapefile. I have only run it on a few large data sets, but I found it took roughly 1 second to join 1,000 records from call to end of end of script.

#!/usr/bin/python

import requests,sys,json,os,csv

import pandas as pd
import geopandas as gpd

lines=[]

# read list of addresses from parameter 1
with open(sys.argv[-2], newline='') as csvfile:	
	for line in csv.DictReader(csvfile):
		lines.append(line)

# convert to pandas
locPd = pd.DataFrame(lines,columns=lines[0].keys())
locPd.convert_dtypes()

locPd = gpd.GeoDataFrame(locPd,  geometry=gpd.points_from_xy(locPd.x.astype('float32'), locPd.y.astype('float32')))

# run spatial joins against parameter 2
ad = gpd.read_file(sys.argv[-1])
locPd = gpd.sjoin(locPd, ad, op="within")

# remove added geometery and index columns
del locPd['geometry']
del locPd['index_right']

# write pandas back to out csv
locPd.to_csv (os.path.splitext(sys.argv[-2])[0]+'-output.csv', index = False, header=True)

I was reading the state GIS’ email list, and noticed that their Civil Boundaries ArcMap has linked to it Census population counts for 1990, 2000, 2010 and 2020 for those interested in population changes for counties and localities without having to download and link the Census data

I was reading the state GIS’ email list, and noticed that their Civil Boundaries ArcMap has linked to it Census population counts for 1990, 2000, 2010 and 2020 for those interested in population changes for counties and localities without having to download and link the Census data. Get it here:

https://gisservices.its.ny.gov/arcgis/rest/services/NYS_Civil_Boundaries/FeatureServer/

NOAA’s GEOS-East Latest Image Retrieval

Earlier today I posted a Leaflet Map showing the current Satellite Image from NOAA’s GEOS-East derived product, the GeoColor. It’s nice but to look really attractive you need to increase the brightness and saturation. If you want a 1920×1080 image that looks nice that you can download, you can use this python code. This ARC GIS Image Server is relatively straightforward to work with — supply a bounding box in Web Mercator, size and a few other parameters and it will spit back the image you want.

import requests, io
from PIL import Image, ImageEnhance, ImageDraw, ImageFont
from datetime import datetime

# new york - new england
url = "https://satellitemaps.nesdis.noaa.gov/arcgis/rest/services/MERGED_GeoColor/ImageServer/exportImage/exportImage?bbox=-9011008.3905,4747656.7008,-7345292.6701,6195679.7647&size=1920%2C1080&format=jpgpng&transparent=false&bboxSR=3857&imageSR=3857&f=image"

# full usa
#url = "https://satellitemaps.nesdis.noaa.gov/arcgis/rest/services/MERGED_GeoColor/ImageServer/exportImage/exportImage?bbox=-14401959.1214%2C2974317.6446%2C-7367306.5342%2C6447616.2099&size=1920%2C1080&format=jpgpng&transparent=false&bboxSR=3857&imageSR=3857&f=image&fbclid=IwAR1HvlVi55Bd933eOhM-Cvou2JgPdKAlcp2go4f8B5pqcioAPvcJyfYT0BA"

response = requests.get(url)
image_bytes = io.BytesIO(response.content)
img = PIL.Image.open(image_bytes)

enhancer = ImageEnhance.Color(img)
img = enhancer.enhance(2.5)

enhancer = ImageEnhance.Brightness(img)
img = enhancer.enhance(1.2)


img.save('/tmp/latest_image.jpg')