I'm having trouble with splitting a shapefile of all the regions in Mongolia into a smaller shapefile based off some input coordinates.
The big shapefile I'm trying to split can be found here: https://marine.rutgers.edu/~cfree/gis-data/mongolia-gis-data/ under Administrative areas - soum boundaries.
I currently have input coordinates of upper left corner:
uly, ulx = 47.82259512683848,108.3279267594774
and bottom right corner:
lry, lrx = 46.85958308337844,109.81664846692111
I've tried importing the shapefile with pyshp and the ogr, and have been having the most success with pyshp, using:
import shapefile sf = shapefile.Reader('file path') sf.bbox()
The results of sf.bbox() return
[-325212.4151875854, 4620806.467562388, 2096670.9522781507, 5777542.4030963285], and I've been trying to find out what projection this is in, I haven't found a clear way to do this with pyshp.
Once I've worked out the current projection of the shapefile and reprojected it into lat/long coordinates, I then need to split it to get rid of everything outside the lat/long coordinates described above.
I've tried to split it so far using the following code, but this produced a shapefile containing nothing.
w = shapefile.Writer(sf.shapeType) for i in sf.iterShapeRecords(): geom = i.shape rec = i.record sxmin, symin, sxmax, symax = geom.bbox if sxmin < ulx: continue elif sxmax > lrx: continue elif symin < uly: continue elif symax > lry: continue w._shapes.append(geom) w.records.append(rec) w.save('Local File')
I can try to provide more details if required.