Read single band from raster

[KrOstir]

When working with large rasters, e.g. with satellite images that can be GB in size, it is useful to be able to read only one band (or a selection of them) to GeoArray. Now the best option is to read the whole dataset and subset by band.

As suggested by @visr ArchGDAL has the ArchGDAL.getband() that could be exported to GeoArrays. GeoArrays.jl could add a method getband(::GeoArray, i::Int)::GeoArray to enable reading by band.

See also the conversation at Julia Discourse: GeoArrays band math.

[evetion]

We use getband from ArchGDAL in the read function. I'd propose to extend it:

function read(fn::AbstractString; masked=true, band=1/"name"/nothing)

Would that satisfy your usecase?

[KrOstir]

function read(fn::AbstractString; masked=true, band=1/"name"/nothing)

Yes, that would be perfect. getband(dataset::AbstractDataset, i::Integer) just has the integer for band number, and that can be used in the read. Default is to read all bands.

[evetion]

Btw, @KrOstir if you have memory issues, you can lazily load the dataset with ga = read(fn::AbstractString; masked=false) and only access the data that fits into memory, e.g. ga[100:200,100:200,1].

[KrOstir]

This works as expected. I will update the documentation to reflect the changes.

[KrOstir]

When reading only one band, I get an error when I omit masked=false.

This throws and error.

julia> fn = download("https://github.com/yeesian/ArchGDALDatasets/blob/master/pyrasterio/RGB.byte.tif?raw=true")
julia> ga_band = GeoArrays.read(fn, band=2)
ERROR: BoundsError: attempt to access 791×718×1 BitArray{3} at index [1:791, 1:718, 2]

And this not.

julia> ga_band = GeoArrays.read(fn, masked=false, band=2)