Hi

I am computing average values on country /state basis globally. I am 
using raster package in R to compute the zonal values &  I read a shape 
file from ArcGIS in R also. Here is my short process i used. I am not 
able to get the zones to match the polygon in the mask.....

/world<-readShapePoly("world.shp", proj4string=CRS("+proj=longlat"))
     
bdat<-brick("global.obs.Prcp.1950-1999.nc",varname="Prcp")                         
# read the NetCDF file uisng raster package
     ras<-raster(ncol=bdat@ncols,nrow=bdat@nrows)                        
                         #create mask with same attributes with raster
     ras[]<-1
     extent(ras)<-extent(world)
     msk <- rasterize(world, ras)
     bdat <- mask(bdat, msk)

zst<-zonal(bdat*12*30,msk,stat='median')/ /#Compute the zonal statistics /

The shape file has 480 polygons and I expect to get the same number of 
zones but to my surprise I get 337 zones..  As such I cannot map 
back/assign the results of the zonal function to the object in order to 
plot. The mask (msk) has the 480 zones but after the step I get 
different number of zones. What am I missing here.

/world$zonal = zst/

*Error in `[[<-.data.frame`(`*tmp*`, name, value = list(zone = c(3, 4,  :
   replacement has 337 rows, data has 480*

I have not supplied the files but I can do so.

Thank you in advance for suggested solutions

rgds
Byman

	[[alternative HTML version deleted]]

