MODIS (version 1.2.11)

extractBits: Extract Bit-Encoded Information and create Weights Raster

Description

This function applies bitops::bitAnd() and bitops::bitShiftR() to convert bit-encoded information. It is also possible to convert this information to a scale from 0 to 1 in order to use it as weighting information in functions like whittaker.raster() or smooth.spline.raster().

Usage

extractBits(
  x,
  bitShift = 2,
  bitMask = 15,
  filename = "",
  datatype = "INT1U",
  NAflag = 255,
  ...
)
makeWeights(
  x,
  bitShift = 2,
  bitMask = 15,
  threshold = NULL,
  decodeOnly = FALSE,
  filename = "",
  datatype = "INT1U",
  NAflag = 255,
  ...
)
maskWater(
  X,
  bitShift = NULL,
  bitMask = NULL,
  keep = NULL,
  datatype = "INT1U",
  NAflag = 255,
  ...
)

Value

A Raster* object.

Arguments

x

matrix, vector or Raster* object.

bitShift

integer. Bit starting point, see Examples and detectBitInfo().

bitMask

integer. Bit mask size, see Examples and detectBitInfo().

filename

character passed to raster::writeRaster(). If not specified, output is written to a temporary file.

datatype

character. Default "INT1U" used for Raster* object. Output datatype, see raster::writeRaster().

NAflag

integer. Default 255 used for Raster* object. Set specific NA value, see raster::writeRaster().

...

Other arguments passed to raster::writeRaster().

threshold

integer. Threshold for valid quality.

decodeOnly

logical. If FALSE (default), convert bits to weights from 0 (not used) to 1 (best quality). If TRUE, only extract selected bits and convert to decimal system.

X

Raster* object.

keep

If NULL (default), bits are only encoded, else an integer vector of values you want to keep (becomes TRUE), the rest becomes NA. See Examples.

Functions

• extractBits(): Extract bit-encoded information from Raster* file

• maskWater(): Masks water (additional information required)

Author

Matteo Mattiuzzi

See Also

detectBitInfo().

Examples

if (FALSE) {

# example MOD13Q1 see MODIS Vegetation Index User's Guide (MOD13 Series; 
# https://lpdaac.usgs.gov/sites/default/files/public/product_documentation/mod13_user_guide.pdf)
# enter in Layers
# See in TABLE 5: Descriptions of the VI Quality Assessment Science Data Sets (QA SDS).
# column 1 (bit) row 2 VI usefulness 
bitShift = 2
bitMask = 15 # ('15' is the decimal of the binary '1111')  
# or try to use
detectBitInfo("MOD13Q1") # not all products are available!
viu <- detectBitInfo("MOD13Q1","VI usefulness") # not all products are available!
viu

# simulate bit info
bit <- round(runif(10*10,1,65536))

# extract from vector
makeWeights(bit,bitShift,bitMask,decodeOnly=TRUE)
# the same as
extractBits(bit,bitShift,bitMask)

# create a Raster object
VIqual <- raster(ncol=10,nrow=10)
VIqual[] <- bit

# extract from Raster object
a <- makeWeights(VIqual,bitShift,bitMask,decodeOnly=TRUE)

# linear conversion of 0 (0000) to 15 (1111) to 1 fo 0 
b <- makeWeights(VIqual,bitShift,bitMask,decodeOnly=FALSE)

threshold=6 # every thing < threshold becomes a weight = 0
c <- makeWeights(VIqual,bitShift,bitMask,threshold,decodeOnly=FALSE)

res <- round(cbind(a[],b[],c[]),2)
colnames(res) <- c("ORIG","Weight","WeightThreshold") 
res

#####
# water mask
tif = runGdal(product="MOD13A2",begin="2009001",end="2009001", extent=extent(c(-9,-3 ,54,58)),
              SDSstring="001",job="delme") # 6.4 MB
x <- raster(unlist(tif))

res1 <- maskWater(x)
plot(res1)

res2 <- maskWater(x,keep=1) # 1 = Land (nothing else)
x11()
plot(res2)

# Land (Nothing else but land) + Ocean coastlines and lake shorelines + shallow inland Water,
# the rest becomes NA
x11()
res3 <- maskWater(x,keep=c(1,2,3)) 
plot(res3)

###############

# as on Linux you can read HDF4 directly you can also do:
if(.Platform$OS.type=="unix")
{
  x <- getHdf(HdfName="MOD13A2.A2009001.h17v03.005.2009020044109.hdf", wait=0) # 6.4 MB
  
  detectBitInfo(x) # just info
  getSds(x) # just info
  
  x <- getSds(x)$SDS4gdal[3] # you need 'VI Quality'
  x  <- raster(x) 
  # plot(x)
  # ex <- drawExtent()
  ex <- extent(c(-580779,-200911,5974929,6529959))
  x  <- crop(x,ex) # just for speed-up
  
  res1 <- maskWater(x)
  plot(res1)
  
  res2 <- maskWater(x,keep=1) # 1 = Land (Nothing else but land), the rest becomes NA
  x11()
  plot(res2)
  
  # Land (Nothing else but land) + Ocean coastlines and lake shorelines + shallow inland Water,
  # the rest becomes NA
  res3 <- maskWater(x,keep=c(1,2,3))
  x11()
  plot(res3)
}
}