1

My exact problem is that I have a file that report vehicular traffic that has 6 columns and two of these columns are latitude and longitude.for example this is the first row and 45. and 9. are latitude and longitude columns.

77c1b49236f442f37531551e2fc25e51    1430352000000000000 45.450824   9.211725    20  0
                         .
                         .
                         .

and on the other hand I have this shapefile of polygons

              id                                           geometry
0           3939_1_1  POLYGON ((504174.2673271392 5003118.268122713,...
1         3939_1_2_1  POLYGON ((508268.9938896392 5003118.268122713,...
2     3939_1_2_2_3_3  POLYGON ((511851.8796318267 5005064.825739901,...
3     3939_1_2_2_3_0  POLYGON ((511340.0388115142 5005064.825739901,...
4       3939_1_2_2_0  POLYGON ((510316.3571708892 5004415.973200838,...

with {'init': 'epsg:32632'}

As you see each polygon has an id I want to check each pair of latitude and longitude and see each one belongs to which polygon and add the column of polygon ids in my first file. So first I read my data file, then I converted the pair of latitude and longitude into a point and I created a geo data frame and then I read the shapefile and tried to check, This is the piece of code that I have tried

import pandas as pd   
import geopandas as gpd
import geopandas.tools
from shapely.geometry import Point, mapping, shape
data_path = "/home/foroogh/PhD/milano-fcd/test"
filename = "2015-04-30-1430352000000000000"
mapDir = "/home/foroogh/PhD/milano-grid"
mapName = "intersection_Milano_W_GRIDIT_NEW.shp"
map_path = os.path.join(mapDir, mapName)
data = os.path.join(data_path, filename)
places = pd.read_csv(data, sep="\t",header = None, names =   ['id','dateandtime', 'latitude', 'longitude', 'vehiclecode', 'velocity'])
places = places[["id", "latitude", "longitude"]]
places["geometry"] = places.apply(lambda row: Point(row["longitude"],  row["latitude"]), axis = 1)
del(places["latitude"], places["longitude"])
places = gpd.GeoDataFrame(places, geometry = "geometry")
places.crs = {'init': 'epsg:32632'}
milano = gpd.GeoDataFrame.from_file(map_path)
print(milano)
milano = milano[["id", "geometry"]]
print(milano)
result = gpd.tools.sjoin(places, milano, how = "left")
print(result.head())

but I do not get the desired result and the problem is that my shapefile is in utm and my points in geo dataframe are in lat and lon and I did not know which one should be converted and how?

2
  • Welcome to GIS SE. As a new user, please take the Tour. Please Edit the question to include the code you have tried, and what error resulted. You should specify the coordinate system in the .prj file associated with the file.
    – Vince
    Commented Nov 15, 2017 at 20:32
  • It doesn't matter whether you convert the lat/lon values (use 4326) to 32632 or vice versa. geopandas might be smart enough to do the join if the lat/lon data has a crs definition but I'm not familiar with that software.
    – mkennedy
    Commented Nov 16, 2017 at 1:01

1 Answer 1

1

I had the same problem and realized geopandas being built on several other libraries makes it a bit convoluted.

What work for me was:

#import some additional libraries: 
import pandas 
from pyproj import Proj, transform
from shapely.geometry import Point

#create your in and out projections with pyroj:
inProj = Proj({'init': 'epsg:xxxx'}) # use df.crs to get
outProj = Proj({'init': 'epsg:4326'})

#Iterate through your column of points and convert:
empty_list = []
for pt in geometry:
     #need to pull out data form Shapely object
     coords_obj = list(pt.coords)
      # transform points us pyroj
      x,y = transform(inProj,outProj, coords_obj[0][0], coords_obj[0][1])
      # put in your empty list as the desired Shapely object
      new_geo.append(Point(x,y))
 # convert to a pandas series
 empty_list =pandas.Series(empyt_list)
 # replace your geopandas column
 geometry = new_geo

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.