Search code examples
pythoncsvmatplotlibcartopy

Using a for loop to make a cartopy plot of lines from different latitude and longitude coordinates from a csv file

The goal for my code is to make a rough roadmap using the latitude and longitude of the exits on the pennsylvania turnpike drawing a line between each exit.

I am using a for loop to plot a line on the map every time it loops. This works if i hard code the latitude and longitude but as soon as i plug in my variables nothing gets plotted. Since the coordinates are in order I am just increasing the index every time it loops to get the next coordinates. I have printed the variables inside the loop and verified they have the desired value. I have tried putting the values in ordered pairs but the plot function didn't like me using nparrays. I'm not sure if there is something simple i am missing, but I appreciate any input.

import netCDF4 as nc
import numpy as np
import matplotlib.pyplot as plt
import cartopy
import cartopy.crs as ccrs
from datetime import datetime, timedelta

# Open the file for highway metadata to read csv data
highway_metadata = open('milestone3data.csv', 'r')
metafile = csv.reader(highway_metadata, delimiter = ',')

# Create empty lists with highway data
highway_loc = []
highway_name = []
highway_lat = []
highway_lon = []
highway_dist = []

# Loop to transfer the csv file's data into the lists
for i in metafile:
    highway_loc.append(i[0])
    highway_name.append(i[1])
    highway_lat.append(float(i[2]))
    highway_lon.append(float(i[3]))
    highway_dist.append(i[4])

def road_map():

    enhighway_lat = enumerate(highway_lat)
    enhighway_lon = enumerate(highway_lon)
    orthographic = ccrs.Orthographic()
    platecarree = ccrs.PlateCarree()
    proj = ccrs.Orthographic(central_longitude = -75, central_latitude = 41)
    ax = plt.axes(projection=proj)
  # Set up the background
    ax.add_feature(cartopy.feature.COASTLINE)
    ax.add_feature(cartopy.feature.STATES)
    ax.set_extent((-85,-70,36,45),crs=ccrs.PlateCarree())

    for i,j in enhighway_lat:
        for k,l in enhighway_lon:
            if i or k <= 30:
                plt.plot([highway_lon[k], highway_lon[k+1]], [highway_lat[i], highway_lat[i+1]], color='black', linewidth=1, marker='o', markersize=3, transform=ccrs.PlateCarree())

    plt.savefig('cartopytest7.png')
    plt.show
road_map()

[This is my most recent output from the program][1]

  [1]: https://i.sstatic.net/lgFrN.png

CSV file contents: (mile marker, name of exit, latitude, longitude, miles from beginning of turnpike)

2,Gateway (Ohio Connection),40.90419167,-80.47158333,1.43
10,New Castle,40.83018056,-80.34196111,10.7
13,Beaver Valley,40.8143,-80.307925,12.87
28,Cranberry,40.67983889,-80.09537778,28.47
30,Warrendale,40.65533889,-80.06116667,31
39,Butler Valley,40.60913611,-79.91924444,39.1
48,Allegheny Valley,40.542025,-79.81022222,47.73
57,Pittsburgh,40.43808889,-79.74956944,56.44
67,Irwin,40.31342778,-79.65476111,67.22
75,New Stanton,40.22173333,-79.59573333,75.39
91,Donegal,40.10915,-79.35231944,90.69
110,Somerset,40.02033056,-79.05208056,109.91
146,Bedford,40.05013889,-78.48615,145.5
161,Breezewood,39.98721667,-78.24472778,161.5
180,Fort Littleton,40.05010556,-77.93954444,179.44
189,Willow Hill,40.09674167,-77.78441389,188.59
201,Blue Mountain,40.15755278,-77.58403333,201.29
226,Carlisle,40.22814722,-77.14782222,226.54
236,Gettysburg Pike,40.19569444,-76.95665556,236.22
242,Harrisburg West Shore,40.21216667,-76.85765278,241.87
247,Harrisburg East,40.21501111,-76.78060278,247.38
266,Lebanon-Lancaster,40.22974444,-76.43095,266.45
286,Reading,40.21805,-76.05189167,286.09
298,Morgantown,40.15990278,-75.88311667,298.33
312,Downingtown,40.06838611,-75.66450278,311.93
320,SR29,40.07641667,-75.52881944,319.33
326,Valley Forge,40.09296667,-75.39591111,326.62
333,Norristown,40.11101111,-75.27921389,333.28
339,Fort Washington,40.13231944,-75.17092222,338.36
340,Virginia Dr,40.13854444,-75.16268611,339.8
343,Willow Grove,40.16166111,-75.11271111,342.91
351,Bensalem,40.13200278,-74.96229444,351.49
352,Street Rd,40.13150833,-74.96445,351.89
353,Neshaminy Falls,40.12916667,-74.94150278,352.67
Solution

  • Okay, based on the discussion above, see below for a solution.

    Notes:

    • Am using pandas DataFames to easily work with the .csv file. the names field is the column names.
    • Am not using orthographic projection at all.
    • Am iterating through the list of highway exits one exit at a time; at each index, am extracting the current and next exits' data - am sure there's a more 'pythonic' way to do this, but this is readable at least.
    • edit: the final index in the loop is length-1
    • Update: Thanks to @SimonWillerton, I've removed the loop.
    import netCDF4 as nc
import numpy as np
import matplotlib.pyplot as plt
import cartopy
import cartopy.crs as ccrs
import pandas as pd

def road_map():    
    # Open the file for highway metadata to read csv data
    highway_metadata = pd.read_csv('milestone3data.csv', names=["loc", "name", "lat", "lon", "dist"])

    proj = ccrs.PlateCarree(central_longitude = -75)
    ax   = plt.axes(projection=proj)
    ax.add_feature(cartopy.feature.COASTLINE)
    ax.add_feature(cartopy.feature.STATES)
    ax.set_extent((-85,-70,36,45),crs=ccrs.PlateCarree())
    
    plt.plot(highway_metadata['lon'], highway_metadata['lat'], \
             color='black', linewidth=1, marker='o', markersize=3, transform=ccrs.PlateCarree())
                       
    plt.savefig('cartopytest7.png')

if __name__ == '__main__':
    road_map()

    This produces the following image: highway_connected

    And, based on this image of the Pennsylvania Turnpike from Wikipedia (source=https://commons.wikimedia.org/wiki/File:Pennsylvania_Turnpike_map.svg#file) I think we have success wiki_pennsylvania_turnpike