Sorry, this should have been:

Alternatively, you could create three RasterLayer objects
etopo1 <- raster("Etopo1.tif", band=1)
etopo2 <- raster("Etopo1.tif", band=2)
etopo3 <- raster("Etopo1.tif", band=3)


On Fri, Jul 23, 2010 at 9:53 AM, Robert J. Hijmans <r.hijm...@gmail.com> wrote:
> If you want a multi-layer object you can create a RasterBrick in stead
> of a RasterLayer. I.e. 'brick()' in stead of  'raster()'
>
> library(raster)
> etopo1 <- brick("Etopo1.tif") # change here
> e <- extent(-90, -32, -60, 15)
> r <- crop(etopo1, e)
>
>
> Alternatively, you could create three RasterLayer objects
> etopo1 <- brick("Etopo1.tif", band=1)
> etopo2 <- brick("Etopo1.tif", band=2)
> etopo3 <- brick("Etopo1.tif", band=3)
>
>
> Robert
>
> On Fri, Jul 23, 2010 at 7:54 AM, Michael Sumner <mdsum...@gmail.com> wrote:
>> I know this is rather simplistic, and has bugs that I just don't have
>> the impetus to track down at the moment, but it should work for your
>> problem. I didn't realize that raster doesn't support (AFAICS)
>> multiband rasters, such as 3-band images.
>>
>> Please use with caution - there are bugs! and very bad practices . . .
>>
>> readpartGDAL <- function(x, xlim = NULL, ylim = NULL, ...) {
>>        require(rgdal)
>>        info <- GDALinfo(x)
>>        offs <- info[c("ll.x", "ll.y")]
>>        scl <- info[c("res.x", "res.y")]
>>        dimn <- info[c("columns", "rows")]
>>        ## brain dead, but easy
>>        xs <- seq(offs[1], by = scl[1], length = dimn[1]) + scl[1]/2
>>        ys <- seq(offs[2], by = scl[2], length = dimn[2]) + scl[2]/2
>>
>>        if (!is.null(xlim)) {
>>                if (!is.numeric(xlim)) stop("xlim must be numeric")
>>                if (!length(xlim) == 2) stop("xlim must be of length 2")
>>                if (!diff(xlim) > 0) stop("xlim[1] must be less than xlim[2]")
>>                xind <- which(xs >= xlim[1] & xs <= xlim[2])
>>        }
>>
>>        if (!is.null(ylim)) {
>>                if (!is.numeric(ylim)) stop("ylim must be numeric")
>>                if (!length(ylim) == 2) stop("ylim must be of length 2")
>>                if (!diff(ylim) > 0) stop("ylim[1] must be less than ylim[2]")
>>                yind <- which(ys >= ylim[1] & ys <= ylim[2])
>>        }
>>        ## probably need a sign check for info["ysign"]
>>        ## reverse for y/x order in readGDAL
>>        rgdal.offset <- rev(c(min(xind), dimn[2] - max(yind)))
>>        rgdal.dim <- rev(c(length(xind), length(yind)))
>>        readGDAL(x, offset = rgdal.offset, region.dim = rgdal.dim, ...)
>> }
>>
>> f <- "Etopo1.tif"
>>
>> d <- readpartGDAL(f, xlim = c(-90, -32), ylim = c(-60, 15))
>>
>> ## assuming you have 3-bands describing RGB
>> ## (see ?SGDF2PCT
>> cols <- SGDF2PCT(d)
>> d$idx <- cols$idx
>> image(d, "idx", col = cols$ct)
>>
>>
>>
>> On Sat, Jul 24, 2010 at 12:09 AM, Rodrigo Aluizio <r.alui...@gmail.com> 
>> wrote:
>>> Well Etopo1 have the colored image as well 
>>> (ftp://ftp.ngdc.noaa.gov/mgg/global/relief/ETOPO1/image/color_etopo1_ice_full.tif.gz).
>>>  I'm not able to import it using readGDAL (because of the vector size), but 
>>> normally when I do so with georeferenced images, the object created have, 
>>> in the data slot, the three band (RGB) color codes as columns. Using 
>>> raster{raster} and then as(x,'SpatialGridDataFrame') these data seems to 
>>> get truncated and unrecognizable to image{sp}.
>>>
>>> Regards.
>>>
>>> Rodrigo
>>>
>>> -----Mensagem original-----
>>> De: Michael Sumner [mailto:mdsum...@gmail.com]
>>> Enviada em: sexta-feira, 23 de julho de 2010 10:53
>>> Para: Rodrigo Aluizio
>>> Cc: R Help
>>> Assunto: Re: [R-sig-Geo] RES: Import Part of a GeoTiff using readGDAL{rgdal}
>>>
>>> I don't know what you mean by "band data structure" and "original colors". 
>>> I just have the generic (ice/cell registered) Etopo1 which is elevation 
>>> values as 16-bit integers. Perhaps you could describe the actual file you 
>>> have more, where you got it, and how you have used it before.
>>>
>>> Cheers, Mike.
>>>
>>> On Fri, Jul 23, 2010 at 8:43 PM, Rodrigo Aluizio <r.alui...@gmail.com> 
>>> wrote:
>>>> Excellent suggestion, but I'm facing another little problem. When 
>>>> importing through raster{raster}, I lost the band data structure (but the 
>>>> data seems to still be there) which turns impossible to reproduce the 
>>>> original colors. Is there a way to bypass this?
>>>>
>>>> Thanks.
>>>>
>>>> Rodrigo.
>>>>
>>>> -----Mensagem original-----
>>>> De: Michael Sumner [mailto:mdsum...@gmail.com] Enviada em:
>>>> quinta-feira, 22 de julho de 2010 12:50
>>>> Para: Rodrigo Aluizio
>>>> Cc: R Help
>>>> Assunto: Re: [R-sig-Geo] Import Part of a GeoTiff using
>>>> readGDAL{rgdal}
>>>>
>>>> With readGDAL you must figure out the impled coordinates, referenced from 
>>>> the top left - it's confusing, but with some experimenting you will get it.
>>>>
>>>> But, far easier is to use the raster package, which will let you deal with 
>>>> the entire Etopo1 if you want by accessing the file as needed, and there 
>>>> are simpler methods to crop the big raster - and coerce to 
>>>> SpatialGridDataFrame if need be.
>>>>
>>>> E.g.
>>>> ## all this works pretty fast as raster takes care of the memory
>>>> library(raster)
>>>> etopo1 <- raster("Etopo1.tif")
>>>>  e <- extent(-90, -32, -60, 15)
>>>> r <- crop(etopo1, e)
>>>> plot(r)
>>>>
>>>> ## if we want to enforce the full memory use as an sp object x <-
>>>> as(r, "SpatialGridDataFrame")
>>>> image(x)
>>>> summary(x)
>>>> Object of class SpatialGridDataFrame
>>>> Coordinates:
>>>>   min max
>>>> s1 -90 -32
>>>> s2 -60  15
>>>> Is projected: FALSE
>>>> proj4string :
>>>> [+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs +towgs84=0,0,0] Number 
>>>> of points: 2 Grid attributes:
>>>>   cellcentre.offset   cellsize cells.dim
>>>> s1         -89.99167 0.01666667      3480
>>>> s2         -59.99167 0.01666667      4500 Data attributes:
>>>>   Min. 1st Qu.  Median    Mean 3rd Qu.    Max.
>>>>  -8086   -4188   -2583   -2025     148    6494
>>>>
>>>>
>>>>
>>>> sessionInfo()
>>>> R version 2.11.1 (2010-05-31)
>>>> x86_64-pc-mingw32
>>>>
>>>> locale:
>>>> [1] LC_COLLATE=English_Australia.1252  LC_CTYPE=English_Australia.1252
>>>> [3] LC_MONETARY=English_Australia.1252 LC_NUMERIC=C [5]
>>>> LC_TIME=English_Australia.1252
>>>>
>>>> attached base packages:
>>>> [1] stats     graphics  grDevices utils     datasets  methods   base
>>>>
>>>> other attached packages:
>>>> [1] rgdal_0.6-27 raster_1.2-6 sp_0.9-65
>>>>
>>>> loaded via a namespace (and not attached):
>>>> [1] grid_2.11.1    lattice_0.18-8 tools_2.11.1
>>>>
>>>> On Fri, Jul 23, 2010 at 1:14 AM, Rodrigo Aluizio <r.alui...@gmail.com> 
>>>> wrote:
>>>>> Hi List.
>>>>>
>>>>> I’m trying to import only a part of interest from a raster image map
>>>>> (let say, ETOPO 1 Global Relief GeoTiff). I realize that, PROBABLY,
>>>>> the readGDAL function can deal with my issue through the arguments:
>>>>>
>>>>>
>>>>>
>>>>> offset - Number of rows and columns from the origin (usually the
>>>>> upper left
>>>>> corner) to begin reading from; presently ordered (y,x) - this may
>>>>> change
>>>>>
>>>>> region.dim - The number of rows and columns to read from the dataset;
>>>>> presently ordered (y,x) - this may change
>>>>>
>>>>>
>>>>>
>>>>> But, I don’t know exactly how to specify that I want the rows (and
>>>>> columns?) which show only South America, I do know the SA coordinate
>>>>> limits, but not the corresponding rows.
>>>>>
>>>>> The real problem is that I don’t know how to use the function
>>>>> resources properly. I’m trying to limit the region during import to
>>>>> avoid memory limit issues, once importing the ETOPO 1 Geotiff tries to 
>>>>> generate a huge vector.
>>>>>
>>>>> Really Sorry About It, Thanks for Your Patience and Attention.
>>>>>
>>>>>
>>>>>
>>>>> -------------------------------------------------------------
>>>>>
>>>>> MSc.  <mailto:r.alui...@gmail.com> Rodrigo Aluizio
>>>>>
>>>>> Centro de Estudos do Mar/UFPR
>>>>> Laboratório de Micropaleontologia
>>>>> Avenida Beira Mar s/n - CEP 83255-000 Pontal do Paraná - PR - Brasil
>>>>>
>>>>>
>>>>>
>>>>>
>>>>>        [[alternative HTML version deleted]]
>>>>>
>>>>>
>>>>> _______________________________________________
>>>>> R-sig-Geo mailing list
>>>>> R-sig-Geo@stat.math.ethz.ch
>>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>>>>
>>>>>
>>>>
>>>>
>>>>
>>>> --
>>>> Michael Sumner
>>>> Institute for Marine and Antarctic Studies, University of Tasmania Hobart, 
>>>> Australia
>>>> e-mail: mdsum...@gmail.com
>>>>
>>>> _______________________________________________
>>>> R-sig-Geo mailing list
>>>> R-sig-Geo@stat.math.ethz.ch
>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>>>
>>>
>>>
>>>
>>> --
>>> Michael Sumner
>>> Institute for Marine and Antarctic Studies, University of Tasmania
>>> Hobart, Australia
>>> e-mail: mdsum...@gmail.com
>>>
>>>
>>
>>
>>
>> --
>> Michael Sumner
>> Institute for Marine and Antarctic Studies, University of Tasmania
>> Hobart, Australia
>> e-mail: mdsum...@gmail.com
>>
>> _______________________________________________
>> R-sig-Geo mailing list
>> R-sig-Geo@stat.math.ethz.ch
>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>
>

_______________________________________________
R-sig-Geo mailing list
R-sig-Geo@stat.math.ethz.ch
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Reply via email to