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?

| improve this question | | | | |
  • 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 Nov 15 '17 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 Nov 16 '17 at 1:01
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
| improve this answer | | | | |

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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