Dear All,

I'm trying to improve a land-use classification based on MODIS product from
2000-2010 (one map per year). Each pixel was classified into 0 (no class), 1
(LU1), 2 (LU2), or 3 (LU3). The overall result of this classification is
quite good, but there are some cases that do not make much sense. For
example, some pixels have LU1 in year one, LU3 in year 2, and LU1 again in
year 3, followed by LU1 during the remaining of the time-series.

        I wrote a  simple function to deal with this problem, but haven't
figured out why it is not working with calc.

#function (very basic) -- sure can be improved...
f = function(x)
{
    if(any(x%in%1) & any(x%in%3))#if  land uses 1 and 3 are present in that
pixel over time...
    {
        r.0 =  grepl("3",x)#find the postion of LU3
        r.M1 = grepl("1",x[which(r.0)-1])#pixel yr-1 of LU3 with LU1
        r.P1 = grepl("1",x[which(r.0)+1])#pixel yr+1 of LU3 with LU1

        if(any(any(r.0)&any(r.M1&r.P1)))#if there is a LU3 somewhere and LU1
before and after that -- problematic when there are two 131 patterns...
        {
            x[r.0]=1#replace the LU3 by LU1
            return(x)}
        else{return(x)}
    }
    else{return(x)}
}

#testing the function:
x = c(2,2,2,0,1,1,3,1,1,1)
f(x)
#[1] 2 2 2 0 1 1 1 1 1 1

#it seems to work (at least when there is only one series of 1-3-1), but
when I include it in the calc #function, something unexpected happens.

#Example:
library(raster)

#creating rasters...
r  <- raster(ncol=5, nrow=5)
r2 <- raster(ncol=5, nrow=5)
r3 <- raster(ncol=5, nrow=5)
r4 <- raster(ncol=5, nrow=5)

set.seed(23)
r[]  <- as.numeric(sample(rep(1:3,1000),25))
r2[] <- as.numeric(sample(rep(1:3,1000),25))
r3[] <- as.numeric(sample(rep(1:3,1000),25))
r4[] <- as.numeric(sample(rep(1:3,1000),25))

rr = stack(r,r2,r3,r4)

#############
f = function(x)
{
    if(any(x%in%1) & any(x%in%3))#if  land uses 1 and 3 are present in that
pixel over time...
    {
        r.0 =  grepl("3",x)#find the postion of LU3
        r.M1 = grepl("1",x[which(r.0)-1])#pixel yr-1 of LU3 with LU1
        r.P1 = grepl("1",x[which(r.0)+1])#pixel yr+1 of LU3 with LU1

        if(any(any(r.0)&any(r.M1&r.P1)))#if there is a LU3 somewhere and LU1
before and after that
        {
            x[r.0]=1#replace the LU3 by LU1
            return(x)}
        else{return(x)}
    }
    else{return(x)}
}

x = c(2,2,2,0,1,1,3,1,1,1)#example of land use over time for a single pixel
(each value represents a yr)

f(x)

############################
ee = calc(rr, fun=f)

ss = stack(rr,ee)
plot(ss)

#it removes all LU2?

Thanks, Paulo

	[[alternative HTML version deleted]]

