[R-sig-Geo] rgdal HFA and spheroid issue

[James Jeuck]

Greetings:

My apologies if I have not posted correctly or research effectively for a
solution to this question...I am brand new to this list and am still
feeling my way around.....however, I did want to post this to see if there
are any suggestions while I dig through the site....

I am working with a large US SE data set derived from NLCD (.img) imagery.
I pull into a "diced" section (10,000X10,000 pixels) into R via:

 map.rast<-readGDAL(fn)
 map.data<-slot(map.rast,"data")

I processed the data to get it ready for, and then run, a randomForest
prediction of land use.

The predicted values were "slot" back into the original map.raster S4:
new.map<-map.rast
slot(new.map,"data")<-data.frame("band1"=prediction.rf$predict,

"band2"=prediction.rf$prob_0,

"band3"=prediction.rf$prob_1)
 outfn<-paste(c("timl_pred",i,j),collapse="_")
 outfn<-paste(c(outfn,".img"),collapse="")


the structure of that S4 object is:

Formal class 'SpatialGridDataFrame' [package "sp"] with 4 slots
  ..@ data       :'data.frame': 6250000 obs. of  3 variables:
  .. ..$ band1: num [1:6250000] 1 1 1 1 1 1 1 1 1 1 ...
  .. ..$ band2: num [1:6250000] 0.99 0.991 0.991 0.988 0.99 ...
  .. ..$ band3: num [1:6250000] 0.0105 0.0095 0.0095 0.012 0.0105 0.011
0.0105 0.013 0.0125 0.0125 ...
  ..@ grid       :Formal class 'GridTopology' [package "sp"] with 3 slots
  .. .. ..@ cellcentre.offset: Named num [1:2] 872640 1751490
  .. .. .. ..- attr(*, "names")= chr [1:2] "x" "y"
  .. .. ..@ cellsize         : num [1:2] 30 30
  .. .. ..@ cells.dim        : int [1:2] 2500 2500
  ..@ bbox       : num [1:2, 1:2] 872625 1751475 947625 1826475
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : chr [1:2] "x" "y"
  .. .. ..$ : chr [1:2] "min" "max"
  ..@ proj4string:Formal class 'CRS' [package "sp"] with 1 slots
  .. .. ..@ projargs: chr "+proj=aea +lat_1=29.5 +lat_2=45.5 +lat_0=23
+lon_0=-96 +x_0=0 +y_0=0 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m
+no_defs"

Notice the "@proj4string at projargs" has GRS80 ellipsoid pulled from the
original raster.

Last, I wrote it back out in HFA format:

writeGDAL(new.map,outfn,drivername="HFA",type="Float32",mvFlag=-1000000,
setStatistics=T)

It all works fine but when I try to pull it into ERDAS Imagine 2015 - I get
an error message saying the Spheriod GRS80 table is missing.  I can pull it
into Arc GIS no problem and use the model builder to reproject to WGS84 -
then pull it into Imagine just fine. However, I am puzzled as to what, if
anything, I am doing wrong.  ERDAS Imagine has GRS80 and I have no trouble
pulling in the original data.  Is there a parameter I am missing in the
"new.map" S4 object, should I edit the @ proj4string, or is there a
parameter in writeGDAL that could fix this?  Many thanks in advance -
jaj
-- 
James Jeuck, Extension Associate
NCSU Extension Forestry
Campus Box 8008
Raleigh, NC 27695
919.515.5574 work
828.734.9795 cell
james_jeuck at ncsu.edu

	[[alternative HTML version deleted]]