Nothing
'[.ursaRaster' <- function(x,i,j,...,drop=FALSE)
{
verbose <- isTRUE(getOption("ursaDevel"))
if (verbose) {
cat("----------\n")
on.exit(cat("==========\n"))
print(match.call())
##~ a <- as.list(match.call())$i
##~ print(class(a))
##~ print(is.symbol(a))
##~ print(is.language(a))
##~ print(is.call(a))
##~ q()
}
dimx <- dim(x$value)
clValue <- class(x$value)
res <- x
missingJ <- missing(j) ## read all lines and (selected) bands
missingI <- missing(i) ## read all bands and (selected) lines
if ((missingJ)&&(missingI))
{
if (length(dimx)==2L)
return(x) ## return(x$value) identiical(x$value,as.matrix(x))==TRUE
}
if (!is.null(dimx))
{
if ((!missingI)&&(is.ursa(i)))
{
# return(res*i) ## removed 20160805
if (nband(i)==1)
res$value[which(is.na(i$value)),] <- NA
else
res$value[which(is.na(i$value))] <- NA
return(res) ## added 20160805
}
if (missingJ)
j <- seq(dimx[1])
if (missingI)
i <- seq(dimx[2])
if ((is.numeric(i))&&(any(abs(i)>dimx[2])))
i <- as.character(i)
if (is.character(i))
{
args <- list(...)
regexp <- if ("regexp" %in% names(args)) args[["regexp"]] else FALSE
i1 <- .getBand(res,i,regexp=regexp)
if (is.null(i1)) {
ind <- .grep(i,ursa(res,"category"))
if (length(ind) == -10000+1) {
return(res[res==ind-1L])
}
else {
if (length(ind)) {
res$value[!(x$value %in% (ind-1L))] <- NA
return(res)
}
else
return(res[integer()])
}
}
else
i <- i1
j <- seq(dimx[1])
missingJ <- TRUE
missingI <- FALSE
}
# if (is.con(res$con))
# res$con$indexZ <- res$con$indexZ[i]
# print(c(j=range(j),i=i))
# res$name <- make.unique(res$name[i],"_")
if (all(i<0))
{
if (!is.na(res$con$posZ[1]))
i <- seq(res$con$posZ)[i]
else
i <- seq(res$dim[2])[i]
}
# res$dim[2] <- length(i)
if (!missingJ)
{
if (verbose)
.elapsedTime("expandIndex:start")
if (is.na(x$con$indexC)[1])
s <- x$grid$columns
else
s <- length(x$con$indexC)
j2 <- j
j <- as.integer(c(t(col(array(NA,dim=c(length(j2),s)))+(j2-1)*s)))
if (max(j)>x$dim[1])
{
opW <- options(warn=0)
warning(paste("Expected",length(j),"bytes, but found",max(j),"bytes"))
options(opW)
j <- j[j<=x$dim[1]]
}
if (verbose)
.elapsedTime("expandIndex:finish")
}
if (.is.con(res$con))
{
if (!missingJ) {
if (is.na(res$con$posR[1]))
res$con$posR <- j2
else
res$con$posR <- res$con$posR[j2]
}
if (!missingI)
{
if (is.na(res$con$posZ[1]))
{
if (is.logical(i))
res$con$posZ <- which(i)
else
res$con$posZ <- i
}
else
res$con$posZ <- res$con$posZ[i]
}
}
else
res$name <- res$name[i]
cl <- class(res$value)
sp <- attr(res$value,"sparse")
res$value <- res$value[j,i,drop=drop]
class(res$value) <- cl
attr(res$value,"sparse") <- sp
res$dim <- dim(res$value)
return(res)
}
res$con$compress <- 0L
grid <- res$grid
con <- res$con
indF <- length(con$fname)
opW <- options(warn=-1)
intOverflow <- is.na(with(con,samples*lines*bands*sizeof))
options(opW)
if (verbose & intOverflow) {
print(c(intOverflow=with(con,as.double(samples)*as.double(lines)*
as.double(bands)*as.double(sizeof))))
# intOverflow <- FALSE
}
intOverflow <- FALSE
internalReading <- intOverflow | con$connection %in% c("gzfile")
externalReading <- !internalReading
if ((1)&&(!missingJ)&&(is.character(j)))
stop("TODO: is.character(j)")
if ((1)&&(!missingI)&&(is.character(i)))
{
args <- list(...)
regexp <- if ("regexp" %in% names(args)) args[["regexp"]] else FALSE
i <- .getBand(res,i,regexp=regexp)
if (is.null(i))
{
op <- options(warn=1)
warning(paste("Speciefied name",paste0("'",i,"'"),"is not in bandname"))
options(op)
return(NULL)
}
missingJ <- TRUE
missingI <- FALSE
}
if ((!missingI)&&(is.logical(i))) {
i <- which(i)
}
if ((!missingI)&&(all(i<0)))
{
if (!is.na(res$con$posZ[1]))
i <- seq(res$con$posZ)[i]
else
i <- seq(res$dim[2])[i]
}
if ((missingJ)&&(missingI)) ## read full
{
if (con$driver %in% c("ENVI","EGDAL")) {
n <- prod(with(con,samples*lines*bands))
xdim <- with(con,c(lines,samples,bands))
if ((con$seek)&&(con$interleave %in% c("bsq","bil"))&&
(externalReading)&&(TRUE)) {
seek(con,where=0L,origin="start",rw="r")
if (con$interleave=="bsq") {
if (con$mode=="integer") {
res$value <- with(con,.Cursa(C_readBsqBandInteger
,fname=con$fname[indF],dim=xdim,index=seq(bands),n=bands
,datatype=datatype,swap=swap
,res=integer(bands*samples*lines)))$res
}
else
res$value <- with(con,.Cursa(C_readBsqBandDouble
,fname=con$fname[indF],dim=xdim,index=seq(bands),n=bands
,datatype=datatype,swap=swap
,res=double(bands*samples*lines)))$res
}
else if (con$interleave=="bil") {
if (con$mode=="integer") {
res$value <- with(con,.Cursa(C_readBilLineInteger2
,fname=con$fname[indF],dim=xdim,index=seq(lines),n=lines
,datatype=datatype,swap=swap
,res=integer(bands*samples*lines)))$res
}
else {
res$value <- with(con,.Cursa(C_readBilLineDouble2
,con$fname[indF],dim=xdim
,lines=seq(lines)
,nline=lines,datatype=datatype,swap=swap
,res=double(with(con,lines*samples*bands))))$res
}
}
dim(res$value) <- with(con,c(samples,lines,bands))
}
else {
res$value <- with(con,.readline(handle,datatype,n,endian))
if (con$interleave=="bil") ##bil[col,band,row] -> R[col,row,band]
{
dim(res$value) <- with(con,c(samples,bands,lines))
res$value <- aperm(res$value,c(1,3,2)) ##(3,1,2)
}
else if (con$interleave=="bip") ##bip[band,col,row] -> R[col,row,band]
{
dim(res$value) <- with(con,c(bands,samples,lines))
res$value <- aperm(res$value,c(2,3,1)) ##c(3,2,1)
}
else if (con$interleave=="bsq") ##bsq[col,row,band] -> R[col,row,band]
{
dim(res$value) <- with(con,c(samples,lines,bands))
# res$value <- aperm(res$value,c(1,2,3)) ##c(2,1,3)
}
else
stop("read ENVI: Error in input header file ",con$interleave
," incorrect interleave type")
}
}
else if (con$driver=="RGDAL") { ## read full
res$value <- .rgdal_getRasterData(con$handle)
dim(res$value) <- with(con,c(samples,lines,bands))
}
else if (con$driver=="GDALRASTER") { ## read full
if (con$datatype %in% c(4,5))
res$value <- array(NA_real_,dim=res$dim)
else
res$value <- array(NA_integer_,dim=res$dim)
if (verbose)
cat("read")
for (b2 in seq_len(con$bands)) {
if (verbose)
cat(".")
res$value[,b2] <- con$handle$read(band=b2
,xoff=0,yoff=0
,xsize=con$samples,ysize=con$lines
,out_xsize=con$samples,out_ysize=con$lines
)
}
if (verbose)
cat(" done!\n")
dim(res$value) <- with(con,c(samples,lines,bands))
}
else if (con$driver=="SF") { ## read full
if (con$datatype %in% c(4,5)) {
res$value <- attr(sf::gdal_read(con$handle$filename,read_data=TRUE),"data")
attr(res$value,"units") <- NULL
}
else {
res$value <- as.integer(attr(sf::gdal_read(con$handle$filename
,read_data=TRUE),"data"))
}
dim(res$value) <- with(con,c(samples,lines,bands))
}
else if (con$driver=="NCDF") { ## read full
# stop("NCDF -- read full")
nc <- ncdf4::nc_open(con$fname[indF])
varName <- con$handle
flip <- attr(varName,"flip")
permute <- attr(varName,"permute")
indT <- attr(varName,"temporal")
indS <- attr(varName,"spatial")
indV <- c(indS,indT)
indL <- seq_along(con$offset)[-indV]
w <- attr(varName,"weight")
# level <- attr(varName,"level")
attributes(varName) <- NULL
# print(data.frame(var=varName,flip=flip,permute=permute,time="???"))
nc.start <- rep(1,length(con$offset))
nc.count <- con$offset
if (length(w)) {
ind <- which(w>0)
isW <- all(diff(ind))==1
if (isW) {
nc.start[indL] <- ind[1]
nc.count[indL] <- length(ind)
w <- w[ind]
}
}
else
isW <- FALSE
res$value <- ncdf4::ncvar_get(nc,varName,start=nc.start,count=nc.count
,collapse_degen=!FALSE)
if (length(dim(res$value))==1)
return(res$value)
dim(res$value) <- nc.count
if (length(indL)==1) {
if (nc.count[indL]>1) { #(con$offset[indL]>1)
if (FALSE)
res$value <- apply(res$value,indV,function(x) sum(x*w))
else {
val <- aperm(res$value,c(indV,indL))
dima <- dim(val)
dim(val) <- c(prod(dima[1:3]),dima[4])
val <- .average(val,weight=w)
dim(val) <- nc.count[indV]
res$value <- val
rm(val)
}
}
else
dim(res$value) <- con$offset[indV]
}
if (permute) {
dima <- length(dim(res$value))
if (dima==2)
res$value <- aperm(res$value,c(2,1))
else if (dima==3)
res$value <- aperm(res$value,c(2,1,3))
else if (dima==4)
res$value <- aperm(res$value,c(2,1,3,4))
else
stop(paste("NCDF: transpose is not implemented for dimension:",dima))
}
if (flip) {
if (length(dim(res$value))>2)
res$value <- res$value[,rev(seq(dim(res$value)[2])),,drop=FALSE]
else
res$value <- res$value[,rev(seq(dim(res$value)[2])),drop=FALSE]
}
dim(res$value) <- with(con,c(samples,lines,bands))
ncdf4::nc_close(nc)
}
else
stop(paste("unknown driver:",con$driver))
if ((con$samples!=res$grid$columns)||(con$lines!=res$grid$rows))
res$value <- res$value[con$indexC,con$indexR,,drop=FALSE]
dimy <- dim(res$value)
dim(res$value) <- c(dimy[1]*dimy[2],dimy[3]) ## t()
res$dim <- dim(res$value)
if (!is.na(con$nodata)) {
if (abs(con$nodata)<1)
res$value[abs(res$value-con$nodata)<1e-27] <- NA
else
res$value[abs(res$value/con$nodata-1)<1e-6] <- NA
}
class(res$value) <- clValue
return(res)
}
else if ((missingJ)&&(!missingI)) ## read band
{
if (verbose)
cat("read bands ",min(i),":",max(i),"\n")
if (is.list(i))
i <- unlist(i)
nb <- length(i)
i <- as.integer(i)
nline <- if (!is.na(con$indexR[1L])) length(con$indexR) else con$lines
if (is.na(con$indexR)[1])
minJ <- 0L
else
minJ <- min(con$indexR)-1L
minI <- min(i)
toWarp <- with(con,(!is.na(indexR)[1])&&(length(indexR)!=lines)||
(!is.na(indexC)[1])&&(length(indexC)!=samples))
if (con$driver %in% c("ENVI","EGDAL")) {
if (con$interleave=="bil")
{
if (externalReading)
{
if (con$seek)
seek(con,where=0L,origin="start",rw="r")
xdim <- with(con,c(lines,samples,bands))
if (con$mode=="integer")
val <- .Cursa(C_readBilBandInteger,con$fname[indF],dim=xdim,index=i
,n=nb,datatype=con$datatype,swap=con$swap
,res=integer(with(con,nb*samples*lines)))$res
else {
val <- .Cursa(C_readBilBandDouble,con$fname[indF],dim=xdim,index=i
,n=nb,datatype=con$datatype,swap=con$swap
,res=double(with(con,nb*samples*lines)))$res
}
}
else
{
n <- nb*con$samples
val <- array(NA,dim=c(n,nline))
conseq <- all(diff(sort(i))==1)
if (conseq) {
for (r in seq(nline)) {
# print(with(con,c(minJ=minJ,rsi=bands*(r-1)+minI-1
# ,samples=samples,sizeof=sizeof)))
pos <- with(con,(minJ+bands*(r-1)+minI-1)*samples*sizeof)
# print(c(r=r,pos=pos,n=n))
if (con$seek)
seek(con,where=pos,origin="start",rw="r")
val[,r] <- with(con,.readline(handle,datatype,n,endian))
}
}
else {
for (r in seq(nline))
{
for (s in seq_along(i)) {
pos <- with(con,(minJ+bands*(r-1)+i[s]-1)*samples*sizeof)
if (con$seek)
seek(con,where=pos,origin="start",rw="r")
s2 <- seq((s-1)*con$samples+1,s*con$samples)
val[s2,r] <- with(con,.readline(handle,datatype,con$samples,endian))
}
}
}
}
if (toWarp) ## added 2013-06-14
{
dim(val) <- with(con,c(samples,nb,lines))
val <- val[,,con$indexR,drop=FALSE]
}
dim(val) <- with(con,c(samples,nb,nline))
val <- aperm(val,c(1,3,2))
}
else if (con$interleave=="bip")
{
n <- with(con,samples*bands)
val <- array(NA,dim=c(nb*con$samples,nline))
ind <- which(with(con,(seq(bands*samples)-1)%%bands+1) %in% i)
for (r in seq(nline))
{
pos <- with(con,(minJ+(r-1)*samples*bands)*sizeof)
if (con$seek)
seek(con,where=pos,origin="start",rw="r")
val[,r] <- with(con,.readline(handle,datatype,n,endian))[ind]
}
dim(val) <- with(con,c(nb,samples,nline))
val <- aperm(val,c(2,3,1))
}
else if (con$interleave=="bsq")
{
isSeq <- identical(i,min(i):max(i))
# val <- array(NA,dim=c(con$samples*nline,nb))
if ((externalReading)&&(TRUE))
{
if (con$seek)
seek(con,where=0L,origin="start",rw="r")
xdim <- with(con,c(lines,samples,bands))
# str(list(i=i,dim=xdim,nb=nb,fname=con$fname))
if (con$mode=="integer")
{
val <- .Cursa(C_readBsqBandInteger,fname=con$fname[indF],dim=xdim,index=i
,n=nb,datatype=con$datatype,swap=con$swap
,res=integer(with(con,nb*samples*lines)))$res
}
else
{
val <- .Cursa(C_readBsqBandDouble,fname=con$fname[indF],dim=xdim,index=i
,n=nb,datatype=con$datatype,swap=con$swap
,res=double(with(con,nb*samples*lines)))$res
}
}
else
{
val <- array(NA,dim=c(con$samples*nline,nb))
if ((1)&&(isSeq))
{
n <- with(con,nline*samples*nb)
pos <- with(con,(minJ+(minI-1)*lines*samples)*sizeof)
if (con$seek)
seek(con,where=pos,origin="start",rw="r")
val[] <- with(con,.readline(handle,datatype,n,endian))
}
else
{
n <- with(con,nline*samples)
for (r in seq(along=i))
{
pos <- with(con,(minJ+(i[r]-1)*lines*samples)*sizeof)
if (con$seek)
seek(con,where=pos,origin="start",rw="r")
val[,r] <- with(con,.readline(handle,datatype,n,endian))
}
}
}
if (toWarp)
{
dim(val) <- with(con,c(samples,lines,nb))
val <- val[,con$indexR,,drop=FALSE] ## before 2012-12-23
}
dim(val) <- with(con,c(samples,nline,nb))
# val <- aperm(val,c(1,2,3))
}
else
stop("Error in input header file ",con$interleave
," incorrect interleave type")
}
else if (con$driver=="RGDAL") { ## read band
val <- .rgdal_getRasterData(con$handle,band=i)
dim(val) <- with(con,c(samples,nline,nb))
}
else if (con$driver=="GDALRASTER") { ## read band
if (con$datatype %in% c(4,5))
val <- array(NA_real_,dim=c(res$dim[1],length(i)))
else
val <- array(NA_integer_,dim=c(res$dim[1],length(i)))
if (verbose)
cat("read chunk band")
for (b2 in seq_along(i)) {
if (verbose)
cat(".")
val[,b2] <- con$handle$read(band=i[b2]
,xoff=0,yoff=0
,xsize=con$samples,ysize=con$lines
,out_xsize=con$samples,out_ysize=con$lines
)
}
if (verbose)
cat(" done!\n")
dim(val) <- with(con,c(samples,nline,nb))
}
else if (con$driver=="SF") { ## read band
rasterio <- list(bands=i)
if (con$datatype %in% c(4,5)) {
val <- attr(sf::gdal_read(con$handle$filename,read_data=TRUE
,RasterIO_parameters=rasterio),"data")
attr(val,"units") <- NULL
}
else {
val <- as.integer(attr(sf::gdal_read(con$handle$filename,read_data=TRUE
,RasterIO_parameters=rasterio),"data"))
}
dim(val) <- with(con,c(samples,nline,nb))
}
else if (con$driver=="NCDF") { ## read band
# stop("NCDF -- read band")
# str(con)
# print(c(i=i))
di <- diff(i)
isC <- ((!length(di))||(abs(unique(di))==1))
nc <- ncdf4::nc_open(con$fname)
varName <- con$handle
flip <- attr(varName,"flip")
permute <- attr(varName,"permute")
indT <- attr(varName,"temporal")
indS <- attr(varName,"spatial")
indV <- c(indS,indT)
indL <- seq_along(con$offset)[-indV]
# print(c(indS=indS,indT=indT,indL=indL))
w <- attr(varName,"weight")
attributes(varName) <- NULL
nc.start <- rep(1,length(con$offset))
nc.count <- con$offset
if (length(w)) {
ind <- which(w>0)
isW <- all(diff(ind))==1
if (isW) {
nc.start[indL] <- ind[1]
nc.count[indL] <- length(ind)
w <- w[ind]
}
}
else
isW <- FALSE
# level <- attr(varName,"level")
# print(data.frame(var=varName,flip=flip,permute=permute,time=time))
if (isC) { #(isC) {
nc.start[indT] <- min(i)
nc.count[indT] <- length(i)
val <- ncdf4::ncvar_get(nc,varName,start=nc.start,count=nc.count
,collapse_degen=FALSE)
# if (length(indL)==1)
# val <- apply(val,c(1,2,indT),mean)
}
else {
nc.count[indT] <- 1
dima <- con$offset
dima[indT] <- length(i)
val <- array(NA_real_,dim=dima)
for (i2 in seq_along(i)) {
nc.start[indT] <- i[i2]
if (length(indL)==1) {
if (indT==4)
val[,,,i2] <- ncdf4::ncvar_get(nc,varName,start=nc.start
,count=nc.count,collapse_degen=FALSE)
else
stop("NCDF column 'time' index?")
}
else {
if (indT==3)
val[,,i2] <- ncdf4::ncvar_get(nc,varName,start=nc.start
,count=nc.count,collapse_degen=FALSE)
else if (indT==1)
val[i2,,] <- ncdf4::ncvar_get(nc,varName,start=nc.start
,count=nc.count,collapse_degen=FALSE)
}
}
}
if (length(indL)==1) {
if (nc.count[indL]>1) { ## con$offset[indL]>1
if (FALSE)
val <- apply(val,indV,function(x) sum(x*w))
else {
val <- aperm(val,c(indV,indL))
dima <- dim(val)
dim(val) <- c(prod(dima[1:3]),dima[4])
val <- .average(val,weight=w)
dim(val) <- nc.count[indV]
}
}
else {
dim(val) <- nc.count[indV]
}
}
if (permute) {
dima <- length(dim(val))
if (dima==2)
val <- aperm(val,c(2,1))
else if (dima==3)
val <- aperm(val,c(2,1,3))
else if (dima==4)
val <- aperm(val,c(2,1,3,4))
else
stop(paste("NCDF: transpose is not implemented for dimension:",dima))
}
if (flip) {
if (length(dim(val))>2)
val <- val[,rev(seq(dim(val)[2])),,drop=FALSE]
else
val <- val[,rev(seq(dim(val)[2])),drop=FALSE]
}
dim(val) <- with(con,c(samples,nline,nb))
ncdf4::nc_close(nc)
}
if (!is.na(con$indexC[1]))
res$value <- val[con$indexC,,,drop=FALSE]
else
res$value <- val
dimy <- dim(res$value)
dim(res$value) <- c(dimy[1]*dimy[2],dimy[3])
if (!is.na(con$nodata)) {
if (abs(con$nodata)<1)
res$value[abs(res$value-con$nodata)<1e-27] <- NA
else {
# print(abs(res$value/con$nodata-1))
res$value[abs(res$value/con$nodata-1)<1e-6] <- NA
}
}
res$con$posZ <- i
class(res$value) <- clValue
return(res)
}
else if ((!missingJ)&&(missingI)) ## read line
{
if (is.list(j))
j <- unlist(j)
j <- as.integer(seq(min(j),max(j)))
nline <- length(j)
minJ <- (min(j)-1L)+min(con$indexR-1L)
if (con$driver %in% c("ENVI","EGDAL")) {
if (con$interleave=="bil") ##bil[col,band,row] -> R[col,row,band]
{
if ((externalReading)&&(TRUE))
{
xdim <- with(con,c(lines,samples,bands))
if (con$mode=="integer")
val <- .Cursa(C_readBilLineInteger,con$fname[indF],dim=xdim
,lines=j+as.integer(min(con$indexR-1L))
,nline=nline,datatype=con$datatype,swap=con$swap
,res=integer(with(con,nline*samples*bands)))$res
else
val <- .Cursa(C_readBilLineDouble,con$fname[indF],dim=xdim
,lines=j+as.integer(min(con$indexR-1L))
,nline=nline,datatype=con$datatype,swap=con$swap
,res=double(with(con,nline*samples*bands)))$res
# print(summary(val))
}
else
{
n <- with(con,nline*bands*samples)
# print(sprintf("seek=%d",with(con,minJ*bands*sizeof*samples)))
if (con$seek)
seek(con,where=with(con,minJ*bands*sizeof*samples)
,origin="start",rw="r")
val <- with(con,.readline(handle,datatype,n,endian))
# print(summary(val))
}
dim(val) <- with(con,c(samples,bands,nline))
val <- aperm(val,c(1,3,2))
}
else if (con$interleave=="bip") ##bip[band,col,row] -> R[col,row,band]
{
n <- with(con,nline*bands*samples)
seek(con,where=with(con,minJ*bands*sizeof*samples),origin="start",rw="r")
val <- with(con,.readline(handle,datatype,n,endian))
dim(val) <- with(con,c(bands,samples,nline))
val <- aperm(val,c(2,3,1))
}
else if (con$interleave=="bsq") ##bsq[col,row,band] -> R[col,row,band]
{
if ((externalReading)&&(TRUE))
{
xdim <- with(con,c(lines,samples,bands))
if (con$mode=="integer")
val <- .Cursa(C_readBsqLineInteger,con$fname[indF],dim=xdim,lines=j
,nline=nline,datatype=con$datatype,swap=con$swap
,res=integer(with(con,nline*samples*bands)))$res
else
val <- .Cursa(C_readBsqLineDouble,con$fname[indF],dim=xdim,lines=j
,nline=nline,datatype=con$datatype,swap=con$swap
,res=double(with(con,nline*samples*bands)))$res
}
else
{
n <- nline*con$samples
val <- array(NA,dim=c(n,con$bands))
for (nb in seq(con$bands))
{
pos <- minJ+(nb-1)*con$lines
if (con$seek)
seek(con,where=with(con,pos*sizeof*samples),origin="start",rw="r")
val[,nb] <- with(con,.readline(handle,datatype,n,endian))
}
}
dim(val) <- with(con,c(samples,nline,bands))
# val <- aperm(val,c(1,2,3))
}
}
else if (con$driver=="RGDAL") { ## read line
nline <- length(j)
minJ <- (min(j)-1L)+min(con$indexR-1L)
val <- .rgdal_getRasterData(con$handle,offset=c(minJ,0)
,region.dim=c(nline,con$samples))
dim(val) <- with(con,c(samples,nline,bands))
}
else if (con$driver=="GDALRASTER") { ## read line
nline <- length(j)
minJ <- (min(j)-1L)+min(con$indexR-1L)
if (con$datatype %in% c(4,5))
val <- array(NA_real_,dim=c(con$samples*nline,res$dim[2]))
else
val <- array(NA_integer_,dim=c(con$samples*nline,res$dim[2]))
if (verbose)
cat("read chunk line")
for (b2 in seq_len(con$bands)) {
if (verbose)
cat(".")
val[,b2] <- con$handle$read(band=b2
,xoff=0,yoff=minJ
,xsize=con$samples,ysize=nline
,out_xsize=con$samples,out_ysize=nline
)
}
if (verbose)
cat(" done!\n")
dim(val) <- with(con,c(samples,nline,bands))
}
else if (con$driver=="SF") { ## read line
nline <- length(j)
minJ <- (min(j)-1L)+min(con$indexR-1L)
rasterio <- list(nXOff=1,nYOff=minJ+1L,nXSize=con$samples,nYSize=nline
# ,nBufXSize=2,nBufYSize=2
)
if (con$datatype %in% c(4,5)) {
val <- attr(sf::gdal_read(con$handle$filename,read_data=TRUE
,RasterIO_parameters=rasterio),"data")
attr(val,"units") <- NULL
}
else {
val <- as.integer(attr(sf::gdal_read(con$handle$filename,read_data=TRUE
,RasterIO_parameters=rasterio),"data"))
}
dim(val) <- with(con,c(samples,nline,bands))
}
else if (con$driver=="NCDF") { ## read line
nc <- ncdf4::nc_open(con$fname[indF])
varName <- con$handle
# str(con$offset)
flip <- attr(varName,"flip")
permute <- attr(varName,"permute")
indT <- attr(varName,"temporal")
indS <- attr(varName,"spatial")
indV <- c(indS,indT)
indL <- seq_along(con$offset)[-indV]
indC <- indS[ifelse(permute,2,1)]
indR <- indS[ifelse(permute,1,2)]
# indR <- ifelse(permute,1L,2L)
# indC <- 3L-indR
w <- attr(varName,"weight")
if (flip)
j <- rev(rev(seq(con$offset[[indR]]))[j])
nline <- length(j)
dj <- diff(j)
isC <- ((!length(dj))||(abs(unique(dj))==1))
# level <- attr(varName,"level")
attributes(varName) <- NULL
nc.start <- rep(1,length(con$offset))
nc.count <- con$offset
##~ print(nc.start)
##~ print(nc.count)
if (isC) { #(isC) {
nc.start[indR] <- min(j)
nc.count[indR] <- length(j)
# print(c(start=nc.start))
# print(c(count=nc.count))
val <- ncdf4::ncvar_get(nc,varName,start=nc.start,count=nc.count
,collapse_degen=FALSE)
}
else { ## slow!
nc.count[indR] <- 1
dima <- con$offset
dima[indR] <- length(j)
val <- array(NA_real_,dim=dima)
for (j2 in seq_along(j)) {
nc.start[indR] <- j[j2]
# print(nc.start)
# print(nc.count)
if (length(indL)==1) {
if (indR==2)
val[,j2,,] <- ncdf4::ncvar_get(nc,varName,start=nc.start
,count=nc.count,collapse_degen=FALSE)
else
stop("NCDF column 'line' index? (multi-level)")
}
else {
# print(indV)
if (indR==2)
val[,j2,] <- ncdf4::ncvar_get(nc,varName,start=nc.start
,count=nc.count,collapse_degen=FALSE)
else
stop("NCDF column 'line' index? (single-level)")
}
}
}
if (length(indL)==1) {
if (con$offset[indL]>1) {
# val <- apply(val,indV,mean)
if (FALSE)
val <- apply(val,indV,function(x) sum(x*w))
else {
val <- aperm(val,c(indV,indL))
dima <- dim(val)
dim(val) <- c(prod(dima[1:3]),dima[4])
val <- .average(val,weight=w)
dim(val) <- nc.count[indV]
}
}
else {
dim(val) <- nc.count[indV] # con$offset[indV]
}
}
if (permute) {
dima <- length(dim(val))
if (dima==2)
val <- aperm(val,c(2,1))
else if (dima==3)
val <- aperm(val,c(2,1,3))
else if (dima==4)
val <- aperm(val,c(2,1,3,4))
else
stop(paste("NCDF: transpose is not implemented for dimension:",dima))
}
if (flip) {
if (length(dim(val))>2)
val <- val[,rev(seq(dim(val)[2])),,drop=FALSE]
else
val <- val[,rev(seq(dim(val)[2])),drop=FALSE]
}
dim(val) <- with(con,c(samples,nline,bands))
ncdf4::nc_close(nc)
}
else
stop("read ENVI: Error in input header file ",con$interleave
," incorrect interleave type")
res$value <- val[con$indexC,,,drop=FALSE]
dimy <- dim(res$value)
dim(res$value) <- c(dimy[1]*dimy[2],dimy[3]) ## t()
if (!is.na(con$nodata)) {
if (abs(con$nodata)<1)
res$value[abs(res$value-con$nodata)<1e-27] <- NA
else
res$value[abs(res$value/con$nodata-1)<1e-6] <- NA
}
res$con$posR <- j
res$dim <- dim(res$value)
class(res$value) <- clValue
return(res)
}
else
stop("TODO")
stop("UNREACHABLE CODE")
}
'.readline' <- function(Fin,datatype,size,endian)
{
if (FALSE) ## benchmark is needed
datatype <- datatype+1000L
if (datatype==1001) ## data.type==1 -> 1-byte short
val <- readBin(Fin,integer(),n=size,size=1,signed=FALSE)
else if (datatype==1) ## data.type==1 -> 1-byte short, unsigned
val <- readBin(readBin(Fin,raw(),n=size*1)
,integer(),n=size,size=1,signed=FALSE)
else if (datatype==1002) ## data.type==2 -> 2-byte short, signed
val <- readBin(Fin,integer(),n=size,size=2,endian=endian,signed=TRUE)
else if (datatype==2) ## data.type==2 -> 2-byte short, unsigned
val <- readBin(readBin(Fin,raw(),n=size*2)
,integer(),n=size,size=2,endian=endian,signed=TRUE)
else if (datatype==1003) ## data.type==3 -> 4-byte int
val <- readBin(Fin,integer(),n=size,endian=endian,signed=TRUE)
else if (datatype==3) ## data.type==3 -> 4-byte int
val <- readBin(readBin(Fin,raw(),n=size*4)
,integer(),n=size,endian=endian,signed=TRUE)
else if (datatype==1004) ## data.type==4 -> 4-byte float
val <- readBin(Fin,double(),n=size,size=4,endian=endian)
else if (datatype==4) ## data.type==4 -> 4-byte float
val <- readBin(readBin(Fin,raw(),n=size*4)
,double(),n=size,size=4,endian=endian)
else if (datatype==5) ## data.type==5 -> 8-byte double
val <- readBin(Fin,double(),n=size,endian=endian)
else if (datatype %in% c(11,1011)) ## data.type==11 -> 1-byte signed short integer
val <- readBin(Fin,integer(),n=size,size=1,endian=endian,signed=TRUE)
else if (datatype %in% c(12,1012)) ## data.type==12 -> 2-byte unsigned short integer
val <- readBin(Fin,integer(),n=size,size=2,endian=endian,signed=FALSE)
else if (datatype %in% c(13,1013)) ## data.type==13 -> 4-byte unsigned integer
{
val <- readBin(Fin,integer(),n=size,size=4,endian=endian,signed=!FALSE) ## need signed=FALSE!
ind <- which(val>(2^31-1))
if (length(ind)) {
opW <- options(warn=1)
warning("Datatype 'UInt32' is supported with restriction. See '?readBin'.")
options(opW)
val[ind] <- val[ind]-2^32
}
}
else
stop("Unsupported data data type: ",datatype)
gc(reset=TRUE)
val
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.