changing the data type of a gdal dataset
Thanks Robert for your advice. Due to the details of my situation I decided it was easier to tough it out in the rgdal package, but fully intend to use the raster package for future projects. In case anybody is interested, I eventually blundered across the answer to my original question. Rather than changing the data type of a copy of a file, I am creating a new transient dataset of the desired data type and dimensions, then adding spatial reference info gleaned from an existing file. Alan # get spatial reference info from existing image file gi <- GDALinfo(infile.name) dims <- as.vector(gi)[1:2] ps <- as.vector(gi)[6:7] ll <- as.vector(gi)[4:5] pref<-attr(gi,"projection") # calculate position of upper left corner and get geotransform ala http://www.gdal.org/gdal_datamodel.html ul <- c(ll[1]-ps[1]/2,ll[2]+(dims[1]+.5)*ps[2]) gt<-c(ul[1],ps[1],0,ul[2],0,ps[2]) # create transient dataset to hold output maps, and add spatial reference info tds <- new("GDALTransientDataset",new('GDALDriver','GTiff'),rows=dims[1],cols=dims[2],type="Float32") .Call("RGDAL_SetProject", tds, pref, PACKAGE = "rgdal") .Call("RGDAL_SetGeoTransform", tds, gt, PACKAGE = "rgdal") # use putRasterData() in a loop to build the output file # save the dataset saveDataset(tds,outfile.name) GDAL.close(tds)
Robert J. Hijmans wrote:
Alan, You can have a look at the code in the raster package on R-forge (see the writeGDAL function). Or perhaps use this package: one of the primary reasons for developing it was to automate (and hide) block by block raster processing --- and quite successfully, it is being used to process raster files that are too large for the leading commercial gis software. See raster::predict to apply a model to a set or raster predictors. Robert On Thu, Nov 19, 2009 at 8:54 PM, Alan Swanson <alan.swanson at umontana.edu> wrote:
Dear R gurus,
I have a function that applies various model prediction functions over a set
of large image files, producing a single output file with the same spatial
extent. Due to memory issues, I'm breaking the input and output files into
tiles. I have this working except for one small issue regarding data types.
I create a new gdal transient dataset by copying an existing one using:
handle <- GDAL.open(fullnames[1],read.only=T)
tds <-
copyDataset(handle,driver=new('GDALDriver','GTiff'),strict=F,options=NULL)
...
putRasterData(tds,t(preds), offset= c(strt[1], 0))
...
saveDataset(tds,outfile.p)
GDAL.close(tds)
Which works great, except that my output always needs to be floating
point, but the input may be byte or integer, in which case the output
dataset retains the format of the input file. So I either need to change
the data type of the new file, or create the new file using:
tds <- new("GDALTransientDataset",driver,dims[1],dims[2], type="Float32")
and then copy the spatial reference information from an existing dataset. I
can't figure out how to do either of these. Your help would be much
appreciated.
Cheers,
Alan
_______________________________________________ R-sig-Geo mailing list R-sig-Geo at stat.math.ethz.ch https://stat.ethz.ch/mailman/listinfo/r-sig-geo