[R-sig-Geo] Error in if (dxsd > tolerance) { : missing value where TRUE/FALSE needed

[David Wahlund]

Hi
I got the same error trying to do something similar. For me the
problem was with the polygon data. One of the shapes was corrupted. My
tip is to try and identify what polygon renders the error and confirm
that it is OK.

-David Wahlund

On Wed, Mar 18, 2009 at 16:01, Wesley Roberts <wroberts at csir.co.za> wrote:
> Dear R-sig-Geo'ers
>
> I am currently trying to extract some NDVI values from a group of NDVI images for specific shapefile polygons. I have 17 individual polygons and 23 NDVI images. I have written the script below to select individual polygons and to then extract the average NDVI within each polygon and store this in a list. Once all 23 images have been processed the script moves to the next polygon and so on and so on, until all polygons have been average for all 23 images. The code essentially extracts the NDVI data as a time series for each polygon.
>
> Unfortunately, as everyone knows, Landsat 7 ETM has a broken scan line corrector meaning that all post 2003 data have large areas of missing data. In my pre-processing I convert these areas to nan using gdal_merge.py (I have to stich and subset two scenes). The code below works well until a large portion of a polygon lies on the nan area of the image (some polygons that only have a small area in the nan process fine). When this happens I get the following error message
>
> Error in if (dxsd > tolerance) { : missing value where TRUE/FALSE needed
>
> The error is related to the following command
>
> x.clip <- temp[!is.na(overlay(temp, brede_add)),]
>
> I have no idea how to solve this error, as a backup plan I am using r.fillnulls in grass to get rid of the stripes but this takes a very long time and is not ideal.
>
> Any additional help or pointers would be greatly appreciated
>
> Many thanks and kind regards
>
> ****************************************************************************************************************************************************************************************************************************************************
>
> library(maptools)
> library(rgdal)
> library(spatstat)
>
> input <- "/media/win-drive/wes2006/Projects/Evapotranspiration/Data/Landsat_NDVI"
> maps <- data.frame((list.files(input, pattern=".img",full.names=FALSE)),(list.files(input, pattern=".img", full.names=TRUE)))
> names(maps) <- c("File", "Path")
> fn=data.frame(maps$File)
>
> rastx<-284055.000
> rasty<-6302475.000
>
> Lst <- list()
>
> brede <- readShapePoly("BredeCatchmentSites.shp", IDvar="COUNT", proj4string=CRS("+proj=utm +zone=34 +south +ellps=WGS84 +datum=WGS84 +units=m +no_defs"))
>
> x=0
> while (x <= 17)
> {
>        x<-x+1
>                nfiles <- length(maps$Path)
>                for (n in 1:nfiles)
>                {
>                   brede_add <- brede[brede$COUNT[x], ]
>                   box <- as.list(as.numeric(bbox(brede_add)))
>                   xstart <- (as.numeric(box[1])-rastx)/30
>                   ystart <- abs((as.numeric(box[4])-rasty)/30)
>
>                   row <- round((as.numeric(box[3])-as.numeric(box[1]))/30)
>                   col <- round((as.numeric(box[4])-as.numeric(box[2]))/30)
>
>                   temp <- readGDAL(as.character(maps$Path[n]),silent=TRUE, offset=c(ystart,xstart), region.dim=c(col,row))
>                   fullgrid(temp)=FALSE
>
>                   x.clip <- temp[!is.na(overlay(temp, brede_add)),]   <<<<<<<<<<<-------------HERE IS THE PROBLEM-------------
>
>                   a <- mean(as.data.frame(x.clip$band1))
>
>                   Lst[n] <- list(as.numeric(a))
>                   print(n)
>                   x.clip<-0
>                 }
>
>        data <- as.numeric(Lst)
>        data2 <- cbind(data,fn)
>        write.table(data, paste(brede$NBALID[x]))
>        print(paste(brede$NBALID[x]))
> }
>
>
>
> Wesley Roberts MSc.
> Researcher: Earth Observation (Ecosystems)
> Natural Resources and the Environment
> CSIR
> Tel: +27 (21) 888-2490
> Fax: +27 (21) 888-2693
>
> "To know the road ahead, ask those coming back."
> - Chinese proverb
>
>
> --
> This message is subject to the CSIR's copyright terms and conditions, e-mail legal notice, and implemented Open Document Format (ODF) standard.
> The full disclaimer details can be found at http://www.csir.co.za/disclaimer.html.
>
>
> and is believed to be clean.  MailScanner thanks Transtec Computers for their support.
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at stat.math.ethz.ch
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>