0

I am using this following code for plotting the trajectories...........

library(openair)

load("GDASNDL1000m.Rdata")

trajLevel(traj,method="hexbin",col="jet",xbin=40,parameters=NULL,
orientation=c(90,0,0),projection="mercator")

result https://github.com/adeckmyn/maps/files/2667752/GDASNDL1000m.zip

Here, I would like to change the base world map with my own shape file.

my shape file is follows....

z1=maptools::readShapePoly("/home/sateeshm/shapefiles/ncmrwf/india_map")
library(maps)
map(z1)

https://github.com/adeckmyn/maps/files/2667336/World-India.zip

#

Now, the actual question is how to link z1 to trajLevel?

M. Sateesh
  • 1
  • 2

1 Answers1

0

To avoid the hard-coded call to "world" in openair, you wil have to create a new world database in the same file based format as that of the "maps" package.

Probably the simplest way to do this, is to use the mapMaker package. This package is not on CRAN but can be found on github. It is the package I used to create the standard world map. The documentation is minimal, but if you don't care about polygon names etc, you can create a "quick and dirty" world map as follows:

# get your new map as a simple list of polygons (or lines)
z1=maps::map(maptools::readShapePoly("india_map"), plot=FALSE)
# create internal representation
z2=mapMaker::map.make(z1)
# write binary files:
mapMaker::map.export.bin(z2, "/my/path/to/world")

# To make map() call this new database:
library(maps)
worldMapEnv="MYMAP"
Sys.setenv("MYMAP"="/my/path/to/") # don't add the "world" !

now map("world") will draw your version of the world map.

Alex Deckmyn
  • 1,017
  • 6
  • 11
  • Dear @Alex, Currently, my india_map is 180W to 180E. Is it possible to customize my india_map as 0 to 360 or 90W to 0,180E,270E with mapMaker ?. If so, give some more information. – M. Sateesh Dec 24 '18 at 12:25
  • Yes, this is in fact quite easy. When calling the map() function, just add the option wrap=c(0,360) (or any two latitudes). Note that if fill=T, this will split the polygons that cross the meridian. However, it will not magically blend polygons that were split at the 180 meridian. So there may remain a few superfluous lines in the map. Removing those may require manual editing. – Alex Deckmyn Dec 30 '18 at 01:00
  • > map("world",wrap=c(-60,300)) Warning message: In if (wrap) coord <- map.wrap(coord) : the condition has length > 1 and only the first element will be used \n Here, I used india_map.shp – M. Sateesh Dec 31 '18 at 12:03
  • you probably have an older version of the maps package, without the new wrapping syntax. The current version on CRAN, v3.3, should not give this error. – Alex Deckmyn Jan 01 '19 at 21:50