1

enter image description here

I am trying to plot a point on top of a shapefile. the Coordinates of the point are accurate. I verified it using the go to XY tool of ArcMap. but when I try to plot overlay using Python the location of both is distorted. following is my code.

loc =gpd.read_file(r'C:\Users\MJ\Desktop\populationFunc\rentalRtes\resapart\gurgaonGISdata\locGISd.shp')
loc = loc.to_crs(epsg = 32643)
loc.plot(figsize = (10,10))
lat = 28.44198794960
long = 77.07818660830
p1 = Point((long, lat))
df = pd.DataFrame({'a':[lat,long]})
po = gpd.GeoDataFrame(geometry = [p1], crs = loc.crs)
po = po.to_crs(epsg = 32643)
po.plot()
fig,ax = plt.subplots(figsize = (15,7))
loc.plot(ax=ax, color = 'none', edgecolor = 'black')
po.plot(ax= ax)`
gene
55.8k3 gold badges115 silver badges196 bronze badges
asked Feb 11, 2021 at 10:54

3 Answers 3

1

This is because you change the loc.crs and then reference it. The following should work:

loc =gpd.read_file(r'C:\Users\MJ\Desktop\populationFunc\rentalRtes\resapart\gurgaonGISdata\locGISd.shp')
loc = loc.to_crs(epsg = 32643)
loc.plot(figsize = (10,10))
lat = 28.44198794960
long = 77.07818660830
p1 = Point((long, lat))
df = pd.DataFrame({'a':[lat,long]})
po = gpd.GeoDataFrame(geometry = [p1], crs = loc.crs)
loc = loc.to_crs(epsg = 32643)
loc.plot(figsize = (10,10))
po = po.to_crs(epsg = 32643)
po.plot()
fig,ax = plt.subplots(figsize = (15,7))
loc.plot(ax=ax, color = 'none', edgecolor = 'black')
po.plot(ax= ax)`

sidenote: Try to not use hard-coded paths in a string. Have a look at Pathlib

answered Feb 11, 2021 at 11:47
1

I dont't know the original projection of loc, nor the content of the GeoDataFrame, but if you fix the new projection of loc, the projection is epsg:32643 (in meters)

 loc = loc.to_crs(epsg = 32643)
 print(loc.crs)
 epsg:32643
 

The unit of EPSG:32643 is meter and not degree (long,lat), therefore

lat = 28.44198794960 #in degree
long = 77.07818660830 #in degree
p1 = Point(long, lat)
po = gpd.GeoDataFrame(geometry = [p1],crs = "EPSG:4326") #in degrees
print(po)
 geometry
0 POINT (77.07819 28.44199)
 print(po.crs)
 EPSG:4326
# conversion to meters
po = po.to_crs(epsg = 32643) # in meters
print(po)
 geometry
0 POINT (703517.314 3147923.581)
print(po.crs)
EPSG:32643

With your solution

po = gpd.GeoDataFrame(geometry = [p1],crs = "EPSG:32643") #in meters
print(po)
 geometry
0 POINT (77.07819 28.44199)
print(po.crs)
EPSG:32643

plot:

fig, ax = plt.subplots(figsize=(6,6))
ax.plot(*po.geometry[0].xy,'o')
ax.axis('equal')
plt.show()

enter image description here

answered Feb 11, 2021 at 12:19
2
  • I have tried using GCS(WGS 84). plotting is accurate in that... but when I try using PCS(UTM Zone 43N) the output is distorted... I have used above bith methods... but not working for me Commented Feb 11, 2021 at 12:59
  • look at my answer correction Commented Feb 11, 2021 at 13:02
-1

Thank you, guys. the issue was with the projection. I changed the projection from GCS to PCS for all files and the issue resolved.

answered Mar 21, 2021 at 10:52

Your Answer

Draft saved
Draft discarded

Sign up or log in

Sign up using Google
Sign up using Email and Password

Post as a guest

Required, but never shown

Post as a guest

Required, but never shown

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

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.