Reputation: 431
I have a Hindu Kush Himalayan region shapefile that can be found in http://rds.icimod.org/Home/DataDetail?metadataId=3924 , and I have lat long of Nepal that have taken from this website (http://rds.icimod.org/Home/DataDetail?metadataId=19590&searchlist=True).
The geographical extent of Nepal is Geographic Extent East: 88.19456,Geographic Extent West:80.0522,Geographic Extent North:30.42472,Geographic Extent South: 26.36836
Now I am trying to subset Nepal from Hindu Kush Himalayan shapefile. This is my code:
mountains<-readOGR("outline.shp") #hindukushhimalayanshapefile
sub <- crop(mountains, extent( 80.0522, 88.18456, 26.36836, 30.42472))
plot(sub)
But the plot sub (that is of Nepal) is not shown properly with a proper outline. There are straight lines showing in the top. How can I get a proper subset of Nepal with the right outline. Am I putting the extent wrong? Help would be appreciated
Upvotes: 1
Views: 220
Reputation: 1233
You cannot subset country boundaries from that shapefile as it does not contain information on country boundaries. You would have to use a shapefile such as this one which contains information on country boundaries: http://rds.icimod.org/Home/DataDetail?metadataId=1218, or subset using a separate Nepal boundary shapefile instead of an extent:
require(maptools)
require(rgdal)
hkh_shp=readOGR("/Downloads/data/outline.shp")
data(wrld_simpl) ##this shapefile is quite coarse, you could substitute another one
nepal_shp=wrld_simpl[which(wrld_simpl$NAME=="Nepal"),]
##CRS are similar but not identical, so need to transform
nepal_shp=spTransform(nepal_shp,crs(hkh_shp))
plot(hkh_shp)
lines(nepal_shp,col="red")
##crop
hkh_sub_shp=crop(hkh_shp,nepal_shp)
plot(hkh_sub_shp) ##note, will look better with higher resolution shapefile
Upvotes: 1