[R-sig-Geo] changing the data type of a gdal dataset

[Alan Swanson]

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
>>
>>     
>
>
>