[R-sig-Geo] Writing an ENVI classification file?

Roger Bivand Roger.Bivand at nhh.no
Sun Jan 6 12:25:53 CET 2013


On Sun, 6 Jan 2013, Jonathan Greenberg wrote:

> Hi Roger:
>
> Not sure if I'm following what you are asking about -- ENVI classification
> files have four characteristics that differentiate themselves from "ENVI
> Standard" files, all of which are stored in the ENVI header:

Hi Jonathan,

What ENVI does is irrelevant, it is what the GDAL ENVI driver does that is 
crucial, since nobody is going to re-do all the drivers separately. It 
turns out that with a small patch to handle the fact that your input file 
has more colors and names than unique values (SVN revision 423 on 
R-Forge) - earlier versions required the numbers to agree - this works:

# I've taken the liberty of shipping your file with rgdal for reference
#fn <- system.file("pictures/test_envi_class.envi", package = "rgdal")[1]
fn <- "test_envi_class.envi"
Gi <- GDALinfo(fn, returnColorTable=TRUE, returnCategoryNames=TRUE)
CT <- attr(Gi, "ColorTable")[[1]]
CT
attr(Gi, "CATlist")[[1]]
with <- readGDAL(fn)
# reports input values and names, because factor creates a 1-based 
# sequence, turn off the console output with:
with <- readGDAL(fn, silent=TRUE)
table(with$band1)
table(as.numeric(with$band1))
# read input as-is to retrieve color table 0-based look-up
with1 <- readGDAL(fn, as.is=TRUE)
table(with1$band1)
# beware 0-base table
ind <- sort(unique(na.omit(with1$band1)))+1
# is this as ENVI?
spplot(with, col.regions=CT[ind])
# write out same data (but coded as factor)
tf <- tempfile()
cN <- levels(with$band1)
CTi <- t(col2rgb(CT[ind]))
# beware 0-base table, convert from factor to integer
with$band1 <- as.integer(with$band1)-1
writeGDAL(with, tf, drivername="ENVI", type="Int16", colorTable=list(CTi),
  catNames=list(cN), mvFlag=7L)
# header OK?
cat(paste(readLines(paste(tf, "hdr", sep=".")), "\n", sep=""), "\n")
wGi <- GDALinfo(tf, returnColorTable=TRUE, returnCategoryNames=TRUE)
CTN <- attr(wGi, "ColorTable")[[1]]
CTN
attr(wGi, "CATlist")[[1]]
# re-read
withN <- readGDAL(tf)
table(withN$band1)
withN1 <- readGDAL(tf, as.is=TRUE)
table(withN1$band1)
# beware 0-base table
ind <- sort(unique(na.omit(withN1$band1)))+1
spplot(withN, col.regions=CTN[ind])

with the same data moved into R, out of R, and back in again. The changes 
induced are in the factor codes given to the values in the raster, but 
this isn't easy to handle sensibly if the input codes are arbitrary 
integers.

As you can see, had you investigated the ENVI driver in GDAL, you would 
have found that it does provide access to the metadata you need, and that 
what the GDAL driver sees is also visible in rgdal. That is why GDAL is:

Geospatial Data Abstraction Library

to avoid messing around with file-format specifics.

In the complete thread, what we don't have so far is a mechanism for 
quantization, that is packaging rules for break points in a numeric raster 
to quantize it into discrete color table entries and category names (class 
intervals). I think that it might be done in the Raster Attribute Table, 
and would welcome input from others that applies generally through GDAL.

Please report back as soon as you've tried rgdal from SVN.

Roger

> file type = ENVI Classification
> --> as opposed to "ENVI Standard"
> and 3 class parameters relating to 1) the number of classes, 2) the class
> colors in RGB format, and 3) the class names.
>
> I'm attaching two small example files to this email, one is a "raw" file
> (test_envi_noclass.envi + .hdr) generated from a
> writeRaster(...,.format="ENVI",datatype="INT2S",NAflag=0); and one is what
> the file should look like, given the variables below (test_envi_class.envi
> + .hdr) -- note the hdr is the only thing that differs.
>
> I have written a script that allows me to simply append this info after
> I've written the file to ENVI format with writeRaster, but I wanted to know
> there was some way to do this from within writeRaster.
>
> Incidentally, here are the R variables:
> class_names=c("Unclassified",
> "temperature_min",
> "temperature_max",
> "temperature_var",
> "precipitation_mean",
> "precipitation_var",
> "wnd_mean",
> "rad_mean",
> "eto_mean",
> "aet_mean",
> "def_mean"
> )
> class_number=length(class_names)
>
> # I don't know if the formatting matters for ENVI
> class_colors=format(
> as.character(c(0,0,0,as.vector(col2rgb(rainbow(n=(class_number-1)))))),
> width=4,justify="right")
>
> # The string to append to the header file should be:
> if(length(var_value)==1)
> {
> header_string=paste(var_name,var_value,sep=" = ")
> } else
> {
> header_string=paste(var_name," = {",(paste(var_value,collapse=",
> ")),"}",sep="")
> # Use brackets
> }
>
> # where var_name would be
> var_name="classes"
> var_value=class_number
> # and..
> var_name="class lookup"
> var_value=class_colors
> # and...
> var_name="class names"
> var_value=class_names
>
> # I haven't written code to swap in "file type = ENVI Classification" yet,
> but that also would need to be modded.
>
> Thanks!
>
> --j
>
> On Sat, Jan 5, 2013 at 2:55 PM, Roger Bivand <Roger.Bivand at nhh.no> wrote:
>
>> This is a continuation of thread: gdal color tables:
>>
>> https://stat.ethz.ch/**pipermail/r-sig-geo/2013-**January/017120.html<https://stat.ethz.ch/pipermail/r-sig-geo/2013-January/017120.html>
>>
>> continued from:
>>
>> https://stat.ethz.ch/**pipermail/r-sig-geo/2012-**December/017071.html<https://stat.ethz.ch/pipermail/r-sig-geo/2012-December/017071.html>
>>
>> Could you please try to use the results there with the ENVI driver? If you
>> have a small example file showing what you want to achieve, that might help.
>>
>> Roger
>>
>>
>>
>> On Sat, 5 Jan 2013, Jonathan Greenberg wrote:
>>
>>  Folks:
>>>
>>> I was wondering if there is any support for writing an ENVI classification
>>> file (which contains additional header info to assign the number of
>>> classes, the color as RGB values, and class names to each class).  Here's
>>> an example from the ENVI header:
>>>
>>> classes = 11
>>> class lookup = {
>>>   0,   0,   0, 255,   0,   0, 139,   0,   0, 205,   0,   0,   0,   0, 255,
>>>   0, 255, 255, 238, 210, 238, 255, 255,   0, 160,  32, 240,   0,   0, 139,
>>> 255, 165,   0}
>>> class names = {
>>> Unclassified, temperature_min, temperature_max, temperature_var,
>>> precipitation_mean, precipitation_var, wnd_mean, rad_mean, eto_mean,
>>> aet_mean, deficit_mean}
>>>
>>> Specifically I'd like to use writeRaster to do it, but I assume this is an
>>> issue with passing rgdal the right info...
>>>
>>> --j
>>>
>>>
>>>
>> --
>> Roger Bivand
>> Department of Economics, NHH Norwegian School of Economics,
>> Helleveien 30, N-5045 Bergen, Norway.
>> voice: +47 55 95 93 55; fax +47 55 95 95 43
>> e-mail: Roger.Bivand at nhh.no
>>
>>
>
>
>

-- 
Roger Bivand
Department of Economics, NHH Norwegian School of Economics,
Helleveien 30, N-5045 Bergen, Norway.
voice: +47 55 95 93 55; fax +47 55 95 95 43
e-mail: Roger.Bivand at nhh.no



More information about the R-sig-Geo mailing list