R/check.popdataPLT2.R

In tfrescino/FIESTA: Forest Inventory Estimation and Analysis

check.popdataPLT <- function(dsn, tabs, tabIDs, pltassgn, pltassgnid, 
	pjoinid, module, popType, popevalid, adj, popFilter, 
	nonsamp.pfilter=NULL, unitarea=NULL, areavar, unitvar, unitvar2=NULL, 
	areaunits, unit.action="keep", removetext="unitarea", strata=FALSE, 
	stratalut=NULL, auxlut=NULL, strvar=NULL, stratcombine=TRUE, pivot=FALSE, nonresp=FALSE, 
	strwtvar="strwt", prednames=NULL, predfac=NULL, pvars2keep=NULL, pdoms2keep=NULL, 
	nullcheck=FALSE, gui=FALSE, dsn_driver="SQLite", unitlevels=NULL, defaultVars=FALSE) {

  ###################################################################################
  ## DESCRIPTION: Checks plot data inputs
  ## - Set plt domains to add to cond (pdoms2keep) - STATECD, UNITCD, COUNTYCD,
  ##		INVYR, MEASYEAR, PLOT_STATUS_CD, RDDISTCD, WATERCD, ELEV, ELEV_PUBLIC,
  ##		ECOSUBCD, CONGCD, INTENSITY, DESIGNCD
  ## Check logical parameters: ACI, strata, stratcombine (if strata=TRUE)
  ## - If ACI, add NF_PLOT_STATUS_CD to pvars2keep and NF_COND_STATUS_CD to cvars2keep
  ## - If unit.action='combine', estimation units are combined if less than 10 plots
  ## - If strata, only 1 auxvar allowed, add to pvars2keep
  ## - If module = SA or MA-greg, add prednames to pvars2keep
  ## - If adj="samp", nonsample adjustment factors calculated at strata level
  ## - If adj="plot", nonsample adjustment factors calculated at plot level
  ## Check unit.action ('keep', 'remove', 'combine').
  ## Check predfac, if module = SA, MA-greg
  ## Import and check plt and pltassgn tables
  ## Check corresponding unique identifiers (puniqueid, pltassgnid)
  ## Merge plt, pltassgn tables to check necessary variables
  ## Check plot data
  ## - Check for missing pvars2keep variables
  ## - Check for missing pdoms2keep variables (STATE, INTENSITY, INVYR, DESIGNCD)
  ## - Get state(s) and inventory year(s) (for titles)
  ## - Generate table of sampled/nonsampled plots (if PLOT_STATUS_CD included)
  ## - If ACI, add table of sampled/nonsampled nonforest plots (if NF_PLOT_STATUS_CD included)
  ## - Generate table of plots by strata, including nonsampled plots (P2POINTCNT)
  ## - If nonresp, generate table of nonsampled plots by strata, substrvar
  ## - Generate and apply nonsamp.pfilter (PLOT_STATUS_CD != 3)
  ## - Check for NA values in pvars2keep variables
  ## - If unitvar = NULL, add unitvar=ONEUNIT
  ## Subset variables for pltx and pltassgnx
  ###################################################################################

  ## Define function
  popTabchk <- function(tablst, popnames) {
    ## DESCRIPTION: check name in tabs list
    for (popname in popnames) {
      chk <- findnm(popname, tablst, returnNULL = TRUE)
      if (!is.null(chk)) return(chk)
    }
    return(NULL)
  }
  
  
  ## Set global variables
  STATECD=PLOT_STATUS_CD=PSTATUSCD=plotsampcnt=nfplotsampcnt=INVYR=
	NF_PLOT_STATUS_CD=NF_COND_STATUS_CD=TPA_UNADJ=methodlst=nonresplut=
	plotqry=pfromqry=pltassgnqry=unitareaqry=auxlutqry=pwhereqry=palias=
	pltx=pltassgnx=popwhereqry=projectid <- NULL


  ###################################################################################
  ## Define necessary plot variables
  ###################################################################################
  datindb=unitindb=stratindb=subcycle99 <- FALSE
  dbconn=pstratvars=plotqry <- NULL
  unitvars <- unique(c(unitvar2, unitvar))
  pltassgnvars <- unique(c(projectid, pltassgnid, unitvars)) 
  returndata=pltindb <- FALSE
  
  ## Define plt variables
  #########################################################################
  pvars2keep <- unique(c(pvars2keep, c("SAMPLING_STATUS_CD", "PLOT_STATUS_CD", 
               "PLOT_NONSAMPLE_REASN_CD", "PSTATUSCD", "INTENSITY")))    
  pdoms <- c("STATECD", "UNITCD", "COUNTYCD", "INVYR", "MEASYEAR", "RDDISTCD", 
	  "WATERCD", "ELEV", "ELEV_PUBLIC", "ECOSUBCD", "CONGCD", "DESIGNCD", "EMAP_HEX")

  ## Get tables from tabs
  ########################################################
  plotchk <- popTabchk(names(tabs), c("plt", "plot"))
  if (is.null(plotchk)) {
    plt <- tabs[[plotchk]]
    puniqueid <- tabIDs[[plotchk]]
  }
  
 
  ###################################################################################
  ## Check parameters
  ###################################################################################

  ## Check adj
  ########################################################
  adjlst <- c("samp", "plot", "none")
  adj <- pcheck.varchar(var2check=adj, varnm="adj", gui=gui,
		checklst=adjlst, caption="adj", multiple=FALSE, stopifnull=TRUE)
  if (adj == "plot" && module == "GB") {
    message("adj='plot' is not typical for GA modules")
  }
  if (adj != "none") {
    pvars2keep <- c(pvars2keep, "MACRO_BREAKPOINT_DIA")
  }

  ## Check ACI (if ACI=FALSE, need to filter COND_STATUS_CD == 1)
  ###################################################################################
  ACI <- pcheck.logical(popFilter$ACI, varnm="ACI", title="ACI?", first="NO", gui=gui)
  if (ACI) {
    pvars2keep <- c(pvars2keep, "NF_SAMPLING_STATUS_CD", "NF_PLOT_STATUS_CD")
  }

  ## Check defaultVars
  ###################################################################################
  defaultVars <- pcheck.logical(defaultVars, varnm="defaultVars", 
                       title="Default variables?", first="NO", gui=gui)

  ## Check unit.action
  ########################################################
  unit.actionlst <- c("keep", "remove", "combine")
  unit.action <- pcheck.varchar(var2check=unit.action, varnm="unit.action", gui=gui,
		checklst=unit.actionlst, caption="unit.action", multiple=FALSE, stopifnull=TRUE)


  ## Set additional pvars2keep
  #####################################################################################
  if (popType %in% c("GRM", "CHNG", "LULC")) {
    pvars2keep <- unique(c(pvars2keep, c("PREV_PLT_CN", "REMPER")))
  } else if (popType == "P2VEG") {
    pvars2keep <- c(pvars2keep, "P2VEG_SAMPLING_STATUS_CD", "P2VEG_SAMPLING_LEVEL_DETAIL_CD",
	                              "SAMP_METHOD_CD")
  } else if (popType == "INV") {
    pvars2keep <- c(pvars2keep, "INVASIVE_SAMPLING_STATUS_CD", "INVASIVE_SPECIMEN_RULE_CD")
  }  

  ## Check strata, strvars
  ###################################################################################
  strata <- pcheck.logical(strata, varnm="strata",
		title="Post stratify?", first="YES", gui=gui, stopifnull=TRUE)

  ## Check strata parameters
  ########################################################
  if (strata) {
    ## pivot
    pivot <- pcheck.logical(pivot, varnm="pivot",
		title="Pivot stratalut?", first="NO", gui=gui)
    ## Check nonresp
    nonresp <- pcheck.logical(nonresp, varnm="nonresp",
		title="Post stratify?", first="YES", gui=gui)
    ## Check stratcombine
    stratcombine <- pcheck.logical(stratcombine, varnm="stratcombine",
		title="Combine strata?", first="YES", gui=gui, stopifnull=TRUE)

    ## Check strvar
    if (is.null(strvar)) stop("must include strvar for post-strat estimates")
    if (length(strvar) > 1) stop("invalid strvar... only 1 variable allowed")
    pltassgnvars <- unique(c(pltassgnvars, strvar))

    if (nonresp) {
      pstratvars <- unique(c(pstratvars, c("PLOT_STATUS_CD", "SAMP_METHOD_CD")))
      pltassgnvars <- unique(c(pltassgnvars, prednames))
    } 
  } else {
    strvar <- NULL
  }

  ## Check predfac
  ###################################################################################
  if (!is.null(predfac)) {
    if (!is.character(predfac)) {
      stop("invalid predfac... must be character string")
    }
    notin <- predfac[!predfac %in% prednames]
    if (length(notin) > 0) {
      warning("invalid predfac... not in prednames: ", toString(notin))
      predfac <- predfac[predfac %in% prednames]
      if (length(predfac) == 0) predfac <- NULL
    }
  }

  ###################################################################################
  ## Check if data are in database
  ###################################################################################
  evalid <- popFilter$evalid
  invyrs <- popFilter$invyrs
  intensity <- popFilter$intensity
  measCur <- popFilter$measCur
  measEndyr <- popFilter$measEndyr
  ACI <- popFilter$ACI
  pltx=pltassgnx=plotnm=ppsanm=tablst <- NULL
  
  
  ## Check if plt and pltassgn tables are data.frame objects
  pltx <- pcheck.table(plt, tabnm="plot", 
                      caption="plot table?", returnsf=FALSE)
  if (!is.null(pltx)) {
    plotnm <- "pltx"
    pltflds <- names(pltx)
  }	
  pltassgnx <- pcheck.table(pltassgn, tabnm="pltassgn", 
                            caption="pltassgn table?", returnsf=FALSE)
  if (!is.null(pltassgnx)) {
    ppsanm <- "pltassgnx"
    ppsaflds <- names(pltassgn)
  }	   
  
  ## Check if plt and pltassgn tables are in a database
  if (!is.null(dsn) && getext(dsn) %in% c("sqlite", "db", "db3", "sqlite3", "gpkg")) {
    datindb=pltindb <- TRUE
    dbconn <- DBtestSQLite(dsn, dbconnopen=TRUE, showlist=FALSE)
    tablst <- DBI::dbListTables(dbconn)
	
	  ## Check name of plt
	  if (is.null(pltx) && !is.null(plt)) {
	    plotnm <- popTabchk(tablst, unique(c(plt, "plt", "plot")))
      if (!is.null(plotnm)) {	  
	      pltflds <- DBI::dbListFields(dbconn, plotnm)
	    }
	  }

    ## Check name of pltassgn
	  if (is.null(pltassgnx) && !is.null(pltassgn)) {
	    ppsanm <- findnm(pltassgn, tablst, returnNULL = TRUE)
      if (is.null(pltassgn)) {
	      message("pltassgn does not exist in database")
	   	  if (!is.null(popevalid)) {
		      return(NULL)
		    }
	    } else {
	      ppsaflds <- DBI::dbListFields(dbconn, ppsanm)
	    }
	  }

    ## If pltassgn is a data.frame, we need to import the plot table as data.frame as well.
	  if (!is.null(pltassgnx) && is.data.frame(pltassgn)) {	  
	    if (!is.null(plotnm)) {
        pltx <- pcheck.table(plt, conn = dbconn, 
                       tabnm = plotnm, returnsf = FALSE)
	      plotnm <- "pltx"
	      pltflds <- names(pltx)
		    pltindb <- FALSE
	    }
    }
  }	

  ## Check unique identifiers
  ##################################################################################
  if (!is.null(plotnm)) {
    puniqueid <- pcheck.varchar(var2check=puniqueid, varnm="puniqueid", gui=gui,
		               checklst=pltflds, caption="UniqueID variable of plot",
		               warn=paste(puniqueid, "not in plt table"), stopifnull=TRUE)
  }
  if (!is.null(ppsanm)) {
    pltassgnchk <- unlist(sapply(pltassgnid, findnm, ppsaflds, returnNULL=TRUE))
	  if (length(pltassgnchk) < length(pltassgnid)) {
	    pltassgnidmiss <- pltassgnid[!pltassgnid %in% pltassgnchk]
	    message("pltassgnid variables are missing from dataset: ", toString(pltassgnidmiss))
      return(NULL)
	  }
  }  
##########################
##########################
  
  ## Check pjoinid
  ## If pjoinid is NULL, set to pltassgnid
  if (!is.null(plotnm) && !is.null(ppsanm)) {
    pjoinidchk <- unlist(sapply(pjoinid, findnm, pltflds, returnNULL=TRUE))
	  if (is.null(pjoinidchk)) {
	    pjoinid <- pltassgnid
	  } 
    pjoinidchk <- unlist(sapply(pjoinid, findnm, pltflds, returnNULL=TRUE))
	  if (length(pjoinidchk) < length(pjoinid)) {
	    pjoinmiss <- pjoinid[!pjoinid %in% pjoinidchk]
	    message("pjoinid variables are missing from plt: ", toString(pjoinmiss))
      return(NULL)
	  } else if (length(pjoinidchk) != length(pltassgnid)) {
	    message("pjoinid must be same number of variables as pltassgnid")
	    return(NULL)
	  } else {
      pjoinid <- pjoinidchk
	  }
  }
  joinqry <- getjoinqry(pjoinid, pltassgnid)
  
  
  #######################################################################
  ## Build pfromqry and pselectqry
  #######################################################################
  pwhereqry <- NULL
  pltassgn. <- "ppsa."
  plt. <- "p."
  
  ## Build from qry for all popfilters (epopfromqry)
  if (is.null(plotnm)) {
    pfromqry <- paste0("\nFROM ", ppsanm, " ", pltassgn.)
	  pflds <- ppsaflds
	  pltassgnvars <- pltassgnid	
  } else if (is.null(ppsanm)) {
    pfromqry <- paste0("\nFROM ", plotnm, " ", plt.)
	  pflds <- pltflds
	  pltassgnvars <- puniqueid
  } else {
    joinqry <- getjoinqry(pjoinid, pltassgnid)
	  pfromqry <- paste0("\nFROM ", ppsanm, " ppsa ", 
		               "\nJOIN ", plotnm, " p ", joinqry)
	  pflds <- c(ppsaflds, pltflds[!pltflds %in% ppsaflds])
	  pltassgnvars <- pltassgnid	
  }	

  ## Define pltfromqry to pass through
  pltfromqry <- pfromqry
  if (popType %in% c("GRM", "CHNG", "LULC")) {
	  prev_pltnm <- findnm("PREV_PLT_CN", pflds, returnNULL = TRUE)
	  if (is.null(prev_pltnm)) {
	    message("dataset is missing PREV_PLT_CN for popType = ", popType)
      return(NULL)
	  }
	  pltfromqry <- paste0(pltfromqry, 
		     "\nJOIN ", plotnm, " pplot ON(pplot.", puniqueid, " = p.", prev_pltnm, ")")
  }				  

  ## Get default variables for plot
  if (defaultVars) {
    #pdoms2keep <- DBvars.default()$pltvarlst
	  pdoms2keep <- pdoms
  } else {
    pdoms2keep <- pflds
  }

  ###################################################################################
  ## If data.frames, check pltx and pltassgnx
  ###################################################################################
  if (!datindb && (!is.null(pltx) || !is.null(pltassgnx))) {
    if (!is.null(pltx)) {
      if (any(duplicated(pltx[[puniqueid]]))) {
        dups <- pltx[[puniqueid]][duplicated(pltx[[puniqueid]])]
        warning(paste("plt records are not unique in: plt:", toString(dups)))
      }

      ## Check for NA values in necessary variables in plt table
      pltx.na <- sum(is.na(pltx[[puniqueid]]))
      if (pltx.na > 0) {
        if (length(pltx.na) == 1) {
          message("1 NA value in ", puniqueid)
        } else {
          message(pltx.na, " NA values in ", puniqueid)
        }
        pltx <- pltx[!is.na(pltx[[puniqueid]]), ]
      }
      ## Set key
      setkeyv(pltx, puniqueid)
    }

    if (!is.null(pltassgnx)) {
      popevalid <- unlist(popevalid)
      pltassgnx <- datFilter(pltassgnx, getfilter("EVALID", popevalid, syntax="R"))$xf
      if (nrow(pltassgnx) == 0) {
        message("evalid popFilter removed all records")
		    return(NULL)
      }
      if (any(duplicated(pltassgnx[[pltassgnid]]))) {
        warning("plot records are not unique in: pltassgn")
      }
      ## Check for NA values in necessary variables in plt table
      pltassgnx.na <- sum(is.na(pltassgnx[[pltassgnid]]))
      if (pltassgnx.na > 0) stop("NA values in ", pltassgnid)

      ## Set key
      setkeyv(pltassgnx, pltassgnid)
    }

    ## Merge plot and pltassgn tables
    #########################################################
    if (!is.null(pltx) && !is.null(pltassgnx)) {
      pjoinid <- pcheck.varchar(var2check=pjoinid, varnm="pjoinid", gui=gui,
		            checklst=pltflds, caption="Joinid variable in plot",
		            warn=paste(pjoinid, "not in plt table"))
      if (is.null(pjoinid)) pjoinid <- puniqueid
        setkeyv(pltx, pjoinid)

        pltassgnx <- pltassgnx[, unique(c(pltassgnid,
		    names(pltassgnx)[!names(pltassgnx) %in% names(pltx)])), with=FALSE]
        setkeyv(pltassgnx, pltassgnid)

      } else if (is.null(pltx)) {
        pltx <- pltassgnx
        puniqueid <- pltassgnid
        if (is.null(pjoinid)) pjoinid <- pltassgnid
      }
  
      ## Check for duplicate plots
      locvars <- c("STATECD", "UNITCD", "COUNTYCD", "PLOT")
      locvars <- locvars[locvars %in% pltflds] 
      if (any(pltx[, duplicated(.SD), .SDcols=locvars]) & (!popType %in% c("GRM", "CHNG"))) {
        warning("duplicated plot locations exist... invalid for estimation")
	      return(NULL)
      }
    }
	
  ######################################################################################
  ## Check essential variables in ppsa or plt tables
  ######################################################################################   
  
  ## Check statecd and invyr (to add to titles)
  statecdnm <- findnm("STATECD", pflds, returnNULL = TRUE)
  if (is.null(statecdnm)) {
    message("STATECD not in dataset.. assuming 1 state in dataset")
  } 
  invyrnm <- findnm("INVYR", pflds, returnNULL = TRUE)
  if (is.null(invyrnm)) {
    message("INVYR not in dataset.. assuming inventory years do not span more ",
		        "than 1 cycle of measurements")
  } 
  
  if (!is.null(statecdnm) && !is.null(invyrnm)) {
    ## Create a table of number of plots by statecd, invyr
    statecda. <- ifelse(statecdnm %in% pltflds, plt., pltassgn.)	  
    invyra. <- ifelse(invyrnm %in% pltflds, plt., pltassgn.)	
	
	  stinvyrvars <- c(paste0(statecda., statecdnm), paste0(invyra., invyrnm))
    stinvyrqry <- paste0("SELECT ", toString(stinvyrvars), ", COUNT(*) AS NBRPLOTS", 
                           pltfromqry,
                           pwhereqry,
                         "\nGROUP BY ", toString(stinvyrvars),
                         "\nORDER BY ", toString(stinvyrvars))
    if (pltindb) {      
      invyrtab <- DBI::dbGetQuery(dbconn, stinvyrqry)
    } else {
	    invyrtab <- sqldf::sqldf(stinvyrqry)
    }
	  if (nrow(invyrtab) > 0) {
	    stcds <- sort(invyrtab[[statecdnm]])
	    states <- pcheck.states(stcds)
	    invyrs <- sort(invyrtab[[invyrnm]])
	  }
  } else if (!is.null(statecdnm)) {
  
    ## Get distinct values of statecd
    statecda. <- ifelse(statecdnm %in% pltflds, plt., pltassgn.)	  
    statecdqry <- paste0("SELECT DISTINCT ", statecda., statecdnm, 
                           pltfromqry, 
						   pwhereqry,
						   "\nORDER BY ", statecda., statecdnm)
    if (pltindb) {      
      stcds <- DBI::dbGetQuery(dbconn, statecdqry)[[1]]
    } else {
	    stcds <- sqldf::sqldf(statecdqry)[[1]]
    }
	  states <- pcheck.states(stcds)
	
  } else if (!is.null(invyrnm)) {
  
    ## Get distinct values of invyr
    invyra. <- ifelse(invyrnm %in% pltflds, plt., pltassgn.)	
    invyrqry <- paste0("SELECT DISTINCT ", invyra., invyrnm, 
                           pltfromqry, 
						   pwhereqry,
						   "\nORDER BY ", invyra., invyrnm)
    if (pltindb) {      
      invyrs <- DBI::dbGetQuery(dbconn, invyrqry)[[1]]
    } else {
	    invyrs <- sqldf::sqldf(invyrqry)[[1]]
	  }
  }
   
  ## Check plot_status_cd
  pstatusvars <- c("PLOT_STATUS_CD", "PSTATUSCD")
  pstatuschk <- unlist(sapply(pstatusvars, findnm, pflds, returnNULL=TRUE))
  if (is.null(pstatuschk)) {
    message("PLOT_STATUS_CD not in dataset.. assuming all plots are at least ",
			      "partially sampled")
  } else {  
    ref_plot_status_cd <- ref_codes[ref_codes$VARIABLE == "PLOT_STATUS_CD", ]
    if (length(pstatuschk) > 1) {
      pstatuscdnm <- pstatuschk[1]
    } else {
      pstatuscdnm <- pstatuschk
    }  
 
    ## Generate table of sampled/nonsampled plots (if ACI, nonforest status included)
    plotsampcntqry <- paste0("SELECT p.", pstatuscdnm, ", COUNT(*) AS NBRPLOTS", 
                           pltfromqry, 
						   pwhereqry, 
						   "\nGROUP BY p.", pstatuscdnm,
						   "\nORDER BY p.", pstatuscdnm)
    if (pltindb) {      
      plotsampcnt <- DBI::dbGetQuery(dbconn, plotsampcntqry)
    } else {
	    plotsampcnt <- sqldf::sqldf(plotsampcntqry)
    }
    plotsampcnt <-
	     cbind(PLOT_STATUS_NM = ref_plot_status_cd[match(plotsampcnt$PLOT_STATUS_CD,
	                      ref_plot_status_cd$VALUE), "MEANING"], plotsampcnt)
						  
						  
    ## Generate and apply nonsamp.pfilter to pltx
    #############################################################################
    if (popType != "ALL") {
      if (!is.null(pstatuscdnm) && (is.null(nonsamp.pfilter) || nonsamp.pfilter == "")) {
        nonsamp.pfilter <- paste0(plt., "PLOT_STATUS_CD != 3")
      }
      nbrnonsampled <- plotsampcnt$NBRPLOTS[plotsampcnt$PLOT_STATUS_CD == 3]
	    if (length(nbrnonsampled) > 0) {
        message("removing ", nbrnonsampled, " nonsampled forest plots")
      }
    }						  
  }

  ## If ACI, check nf_plot_status_cd
  if (ACI) {
    nfstatusvars <- c("NF_PLOT_STATUS_CD", "PSTATUSNF")
    nfstatuschk <- unlist(sapply(nfstatusvars, findnm, pflds, returnNULL=TRUE))
    if (is.null(nfstatuschk)) {
      message("NF_PLOT_STATUS_CD not in dataset.. assuming all ACI nonforest plots are at least ",
			"partially sampled")
    } else {  
      ref_nf_plot_status_cd <- ref_codes[ref_codes$VARIABLE == "NF_PLOT_STATUS_CD", ]
      if (length(nfstatuschk) > 1) {
        nfstatuscdnm <- nfstatuschk[1]
      } else {
        nfstatuscdnm <- nfstatuschk
      }  
 
      ## Generate table of sampled/nonsampled plots (if ACI, nonforest status included)
      nfplotsampcntqry <- paste0("SELECT p.", nfstatuscdnm, ", COUNT(*) AS NBRPLOTS", 
                           pltfromqry, 
						   pwhereqry, 
						   "\nGROUP BY p.", nfstatuscdnm,
						   "\nORDER BY p.", nfstatuscdnm)
      if (pltindb) {      
        nfplotsampcnt <- DBI::dbGetQuery(dbconn, nfplotsampcntqry)
      } else {
	      nfplotsampcnt <- sqldf::sqldf(nfplotsampcntqry)
      }
      nfplotsampcnt <- nfplotsampcnt[!is.na(nfplotsampcnt$NF_PLOT_STATUS_CD), ]
	    if (nrow(nfplotsampcnt) > 0) {
        nfplotsampcnt <-
	 	        cbind(NF_PLOT_STATUS_NM = ref_nf_plot_status_cd[match(nfplotsampcnt$NF_PLOT_STATUS_CD,
		                  ref_nf_plot_status_cd$VALUE), "MEANING"], nfplotsampcnt)
	
        if (!is.null(plotsampcnt)) {
          plotsampcnt <- rbindlist(list(plotsampcnt, nfplotsampcnt), use.names=FALSE)
        } else {
          plotsampcnt <- nfplotsampcnt
		    }
		
        ## Generate and apply nonsamp.pfilter to pltx
        #############################################################################
        if (popType != "ALL") {
          nfnonsamp.pfilter <- paste(plt., "NF_PLOT_STATUS_CD != 3")
          nbrnfnonsampled <- nfplotsampcnt$NBRPLOTS[nfplotsampcnt$NF_PLOT_STATUS_CD == 3]
	        if (length(nbrnfnonsampled) > 0) {
            message("removing ", nbrnfnonsampled, " nonsampled nonforest plots")
          }
		      if (!is.null(nonsamp.pfilter)) {
		        nonsamp.pfilter <- paste0(nonsamp.pfilter, " AND ", nfnonsamp.pfilter)
          } else {
		        nonsamp.pfilter <- nfnonsamp.pfilter
      	  }  
		    }		
      }
	  }
  }
    
  ## Check unitvars in ppsa or plt
  if (!is.null(unitvars)) {
    unitvarchk <- unlist(sapply(unitvars, findnm, pflds, returnNULL=TRUE))
    if (is.null(unitvarchk)) {
	    message("unitvars must be included in dataset")
      return(NULL)
	  } else {
	    unitvars <- unitvarchk
	  }	
    unitvarsa. <- ifelse(all(unitvars %in% pltflds), plt., pltassgn.)	  	
	  unitvarqry <- paste0("SELECT DISTINCT ", toString(paste0(unitvarsa., unitvars)), 
                      pltfromqry,
                      pwhereqry,
                      "\nORDER BY ", toString(paste0(unitvarsa., unitvars)))
    if (pltindb) {      
      unitvartab <- data.table(DBI::dbGetQuery(dbconn, unitvarqry))
    } else {
	    unitvartab <- data.table(sqldf::sqldf(unitvarqry))
    }
	  if (nrow(unitvartab) > 0) {
	    punit.vals <- sort(do.call(paste, unitvartab[, unitvars, with=FALSE]))
	    nbrunits <- length(punit.vals)
	  } else {
	    message("unitarea has no rows...")
	    message(unitvarqry)
	  }
	
  } else {
  
    unitvar <- checknm("ONEUNIT", pflds)
    message("no unitvar specified...  adding a variable named ", unitvar)
    unitvar=unitvars <- "ONEUNIT"
	  nbrunits <- 1
	  punit.vals <- 1
  }
  pltassgnvars <- unique(c(pltassgnvars, unitvars))
  
  ## Check strvar in ppsa or plt
  if (strata) {
    if (!is.null(strvar)) {
      strvar <- findnm(strvar, pflds, returnNULL=TRUE)
      if (is.null(strvar)) {
	      message("strata=TRUE, strvar must be included in dataset")
        return(NULL)
	    }	
	  }
	  pltassgnvars <- unique(c(pltassgnvars, strvar))
  }
  
  ## Check prednames
  if (!is.null(prednames)) {
    prednameschk <- unlist(sapply(prednames, findnm, pflds, returnNULL=TRUE))
    if (is.null(prednameschk)) {
	    message("no prednames are not found in dataset")
	    return(NULL)
	  } else if (length(prednameschk) < length(prednames)) {
	    message("prednames are missing from dataset: ", 
	                       toString(prednames[!prednames %in% prednameschk]))
      return(NULL)
	  } else {
	    prednames <- prednameschk
	  }
	  pltassgnvars <- unique(c(pltassgnvars, prednames))
  }

  ######################################################################################
  ## Check unitarea and auxlut
  ######################################################################################
  vars2keep=unitareax <- NULL
  if (module == "SA" && "AOI" %in% names(unitarea)) {
    vars2keep <- "AOI"
  }
  
  if (is.null(unitarea)) {
    message("unitarea is missing... include with population data if total estimates are desired")
  } else {
    if (nbrunits == 1 && !is.data.frame(unitarea)) {
      if (is.vector(unitarea) && length(unitarea) == 1 && is.null(chkdbtab(tablst, unitarea))) {
	      if (is.numeric(unitarea)) {
	        unitarea <- data.table(ONEUNIT = 1, unitarea)
	      } else if (sum(grepl(",", unitarea)) > 0) {
	        unitarea <- data.table(ONEUNIT = 1, as.numeric(gsub(",", "", unitarea)))
	      }
	      unitvar=unitvars <- "ONEUNIT"
		    if (is.null(areavar)) areavar <- "AREA_USED"
	        names(unitarea) <- c(unitvar, areavar)
	      }
	    } else if (is.data.frame(unitarea)) {	
	      ## Check unitarea
        unitareax <- pcheck.table(unitarea, tabnm="unitarea", 
		                   caption="unitarea table?", returnsf=FALSE)
	      if (nbrunits == 1 && unitvar %in% names(unitarea)) {
	        unitarea[[unitvar]] <- 1
	      }
	    } 
	    if (is.null(unitareax) && !is.null(chkdbtab(tablst, unitarea))) {
        unitindb <- TRUE
        unitarea_layer <- chkdbtab(tablst, unitarea)
	      unitareaflds <- DBI::dbListFields(dbconn, unitarea_layer)
        unitareaqry <- paste0("SELECT * \nFROM ", unitarea_layer)
	
	      if (!is.null(popevalid)) {	
	        uevalidnm <- findnm("EVALID", unitareaflds, returnNULL = TRUE) 
		
		      ## Check popevalid values in database	  
	        uevalidqry <- paste0("SELECT DISTINCT ", uevalidnm, 
	                    "\nFROM ", unitarea_layer,
	                    "\nORDER BY ", uevalidnm)
        if (!is.data.frame(unitarea)) {      
          uevalidvals <- DBI::dbGetQuery(dbconn, uevalidqry)[[1]]
	      } else {
	        uevalidvals <- sqldf::sqldf(uevalidqry)[[1]]
	      }
        uevalidmiss <- popevalid[!popevalid %in% uevalidvals]
        if (any(!popevalid %in% uevalidvals)) {
	        message("evalids are missing in unitarea: ", 
		                        toString(popevalid[!popevalid %in% uevalidvals]))
          return(NULL)
        }
		
        if (!is.null(popevalid) && !is.null(uevalidnm)) {
          unitareaqry <- paste0(unitareaqry, 
		                     "\nWHERE ", evalidnm, " IN(", toString(popevalid), ")")
        }
	    } else {
	      unitareaqry <- NULL
	    }
      unitareax <- pcheck.table(unitarea, conn=dbconn,
              tabnm="unitarea", caption="unitarea?", 
              nullcheck=nullcheck, tabqry=unitareaqry, returnsf=FALSE)
			  
	    if (is.null(unitareax)) {
	      message("invalid unitareax")
	      return(NULL)
	    } 
	  
	    unitareax <- unique(unitareax[, c(unitvars, areavar), with=FALSE])	  
	    if (any(duplicated(unitareax[, unitvars, with=FALSE]))) {
	      message("unitarea is invalid... multiple unitvars exist")
		    return(NULL)
	    }
	  }
	
    ## Check areavar 
	  #############################################
    areavar <- pcheck.varchar(var2check=areavar, varnm="areavar", gui=gui,
		   checklst=names(unitareax), caption="Area variable?", stopifnull=TRUE)
		   
    ## Check if areavar column is numeric
    if (!is.numeric(unitareax[[areavar]])) {
      if(sum(grepl(",", unitareax[[areavar]])) > 0)
        unitareax[[areavar]] <- as.numeric(gsub(",", "", unitareax[[areavar]]))
        if (!is.na(unitareax[[areavar]])) {
          stop("invalid areavar in unitarea.. must be a number")
        }
      }

      ## Check for NA values in areavar
      if (any(is.na(unitareax[[areavar]]))) {
        navals <- unitareax[is.na(get(areavar)), ]
        message("there are NA values in area.. removing from table")
        print(navals)
        unitareax <- unitareax[!is.na(get(areavar)), ]
      }
		   
	    ## Check unitvars
	    if (!all(unitvars %in% names(unitareax))) {
        unitvarchk2 <- unlist(sapply(unitvars, findnm, names(unitareax), returnNULL=TRUE))
	    if (!is.null(unitvarchk2)) {
	      setnames(unitareax, unitvarchk2, names(unitvarchk2))
	    }
 	  }
	  if (!all(unitvars %in% names(unitareax))) {
	    message("unitvars are not in unitareax: ", 
	                          toString(unitvars[!unitvars %in% names(unitareax)]))
	    return(NULL)
	  }	
	
	  ## Check unit.vals with punit.vals
	  unit.vals <- sort(do.call(paste, unitareax[, unitvars, with=FALSE]))
	
	  if (any(is.na(match(unit.vals, punit.vals)))) {
	    unit.miss <- unit.vals[is.na(match(unit.vals, punit.vals))]
	    message("unit in unitarea is not in plot data: ", toString(unit.miss))
 	  
	    if (unit.action == "remove") {
        unitareax[, `:=`(MATCH, do.call(paste, .SD)), .SDcols = unitvars]
        unitareax <- unitareax[!MATCH %in% unit.miss, ]
        unitareax[, `:=`(MATCH, NULL)]
      }
		
	  } else if (any(is.na(match(punit.vals, unit.vals)))) {
	    punit.miss <- unit.vals[is.na(match(punit.vals, unit.vals))]
	    message("unit in plot data is not in unitarea: ", toString(punit.miss))
	    return(NULL)
	  }  

	  ## Sum areavar by unitvars to aggregate any duplicate rows
    unitareax <- unitareax[, sum(.SD, na.rm=TRUE), by=c(unitvars, vars2keep), .SDcols=areavar]
    setnames(unitareax, "V1", areavar)
  }	


  ######################################################################################
  ## Check auxlut
  ######################################################################################
  strunitvars <- unitvars
  
  auxvars <- unique(c(strvar, prednames))
  auxtabnm <- ifelse(strata, "stratalut", "auxlut")
  if (!is.null(auxvars)) {
    if (is.null(auxlut)) {
      if (strata) {
        message("strata=TRUE and stratalut is missing... ")
	    } else {
	      message("prednames != NULL and auxlut is missing... ")
	    }
	    return(NULL)
    }	
	  ## Check if auxlut is a data.frame R object
    auxlutx <- pcheck.table(auxlut, tabnm = auxtabnm, 
		                    caption = paste(auxtabnm, " table?"), returnsf=FALSE)

    if (is.null(auxlutx) && !is.null(chkdbtab(tablst, auxlut))) {
      auxindb <- TRUE
      auxlut_layer <- chkdbtab(tablst, auxlut)
	    auxlutflds <- DBI::dbListFields(dbconn, auxlut_layer)
      auxlutqry <- paste0("SELECT * \nFROM ", auxlut_layer)
	  
	    if (!is.null(popevalid)) {
	      evalidnm <- findnm("EVALID", auxlutflds, returnNULL = TRUE) 
        if (!is.null(popevalid) && !is.null(evalidnm)) {
          auxlutqry <- paste0(auxlutqry, 
		                 "\nWHERE ", evalidnm, " IN(", toString(popevalid), ")")
        }
	      auxlutx <- pcheck.table(auxlut, conn = dbconn,
                     tabnm = auxtabnm, caption = paste(auxtabnm, " table?"),
		             tabqry = auxlutqry, returnsf = FALSE)
			
	      if (is.null(auxlutx) || nrow(auxlutx) == 0) {
		      ## Check popevalid values in database	  
	        sevalidqry <- paste0("SELECT DISTINCT ", evalidnm, 
	                    "\nFROM ", auxlut_layer)
          if (!is.data.frame(auxlutx)) {      
            sevalidvals <- DBI::dbGetQuery(dbconn, sevalidqry)[[1]]
	        } else {
	          sevalidvals <- sqldf::sqldf(sevalidqry)[[1]]
	        }
          sevalidmiss <- popevalid[!popevalid %in% sevalidvals]
          if (any(!popevalid %in% sevalidvals)) {
	          message("evalids are missing in stratalut: ", toString(popevalid[!popevalid %in% sevalidvals]))
            return(NULL)
          }
	      }
	    } 
	  }
	
	  if (strata) {
      ## Check strvar
      strvar <- pcheck.varchar(var2check = strvar, varnm = "strvar", gui=gui,
		     checklst = names(auxlutx), caption = "strata variable",
		     warn=paste(strvar, "not in strata table"), stopifnull=TRUE)

 	    ## Check unitlevels for collapsing
	    if (!is.null(unitlevels)) {
	      unitvals <- unique(stratalut[[unitvar]])
	      if (length(unitlevels) != length(unitvals)) {
		      misslevels <- unitvals[!unitvals %in% unitlevels]
	        message("unitlevels does not match unitvals... missing: ", toString(misslevels))
		      return(NULL)
        } else if (any(is.na(match(unitlevels, unitvals)))) {
	        difflevels <- unitvals[is.na(match(unitlevels, unitvals))]
	        message("unitlevels does not match unitvals... missing: ", toString(difflevels))
		      return(NULL)
        }		
      } 	
	
	    if (pivot) {
        ## Pivot auxlut table based on strwtvar
        stratalutx <- strat.pivot(auxlutx, unitvars = unitvars, 
                              strvar = strvar, strwtvar = strwtvar)
      } else {			
        stratalutx <- auxlutx
	    }
	
      strunitvars <- unique(c(unitvars, strvar))
	    strunitvarsqry <- paste0("ppsa.", strunitvars)
	    P2POINTCNTqry <- paste0("SELECT ", toString(strunitvarsqry), ", COUNT(*) AS NBRPLOTS", 
	                        pltfromqry,
							pwhereqry,
							"\nGROUP BY ", toString(strunitvarsqry),
							"\nORDER BY ", toString(strunitvarsqry))
	    if (pltindb) {      
        P2POINTCNT <- data.table(DBI::dbGetQuery(dbconn, P2POINTCNTqry))
      } else {
	      P2POINTCNT <- data.table(sqldf::sqldf(P2POINTCNTqry))
      }
      setkeyv(P2POINTCNT, strunitvars)
 
      ## If nonresp, get Response Homogeneity Groups for WestFest
      #####################################################################
      if (nonresp) {
        nonrespvars <- c("PLOT_STATUS_CD", "SAMP_METHOD_CD")
        nonrespchk <- unlist(sapply(nonrespvars, findnm, pflds, returnNULL=TRUE))
        if (is.null(unitvarchk)) {
	        message("the following vars must be included in dataset: ", toString(nonrespvars))
          return(NULL)
	      }
        pltnrqry <- paste0("SELECT ", toString(c(strunitvarsqry, puniqueid, nonrespvars)),
                           pltfromqry,
                           whereqry,
                           "\nORDER BY ", toString(c(strunitvarsqry, puniqueid, nonrespvars)))						   
	      if (pltindb) {      
          pltnr <- data.table(DBI::dbGetQuery(dbconn, pltnrqry))
        } else {
	        pltnr <- data.table(sqldf::sqldf(pltnrqry))
        }
	
        RHGdat <- getRHG(pltx = pltnr, puniqueid = puniqueid, 
		                 unitvars = unitvars, strvar = strvar)  
        pltnr <- RHGdat$pltx
        RHGlut <- RHGdat$RHGlut
        P2POINTCNT <- RHGdat$P2POINTCNT
        nonresplut <- RHGdat$nonresplut  

        pltassgnvars <- unique(c(pltassgnvars, "RHG"))  
      } 	
    }
  }

   
  ######################################################################################
  ## Query and Import pltassgnx
  ######################################################################################
  ppsaselectvars <- {}
  ppsavars <- pltassgnvars[pltassgnvars %in% ppsaflds]
  if (length(ppsavars) > 0) {
    ppsaselectvars <- paste0(pltassgn., ppsavars)
  }
  pltvars <- pltassgnvars[pltassgnvars %in% pltflds]
  pltvars <- pltvars[!pltvars %in% pltassgnvars]
  if (length(pltvars) > 0) {   
    ppsaselectvars <- c(ppsaselectvars, paste0(plt., pltvars))
  }
  
  ## Add puniqueid to pltassgn table to join to all other tables
  if (!puniqueid %in% ppsaselectvars) {
    ppsaselectvars <- c(paste0(plt., puniqueid), ppsaselectvars)
  }    
  
  ## Build select for pltassgnx
  ppsaselectqry <- paste0("SELECT ", toString(ppsaselectvars))
  if (unitvar == "ONEUNIT" && !"ONEUNIT" %in% pflds) {
    ppsaselectqry <- paste0(ppsaselectqry, ", 1 AS ", unitvar)
  }
  
  ## Build query for pltassgnx
  pltassgnxqry <- paste0(ppsaselectqry, 
                         pltfromqry,
                         pwhereqry)
  if (pltindb) {      
    pltassgnx <- data.table(DBI::dbGetQuery(dbconn, pltassgnxqry))
  } else {
	  pltassgnx <- data.table(sqldf::sqldf(pltassgnxqry))
  }
  setkeyv(pltassgnx, puniqueid)
  
  
  ######################################################################################
  ## Query and Import pltx
  ###################################################################################### 
  pselectvars <- {}
  if (defaultVars) {
    pltxvars <- unique(c(puniqueid, pvars2keep, pdoms2keep))
  } else {
    pltxvars <- unique(c(puniqueid, pvars2keep, pdoms2keep))
  }
  ppsavars <- pltxvars[pltxvars %in% ppsavars]
  if (length(ppsavars) > 0) {
    pselectvars <- paste0(pltassgn., ppsavars)
  }
  pltvars <- pltxvars[pltxvars %in% pltflds]
  if (length(pltvars) > 0) {
    pselectvars <- c(pselectvars, paste0(plt., pltvars))
  }
  
  ## Build select for pltx
  #pselectqry <- paste0("SELECT ", toString(pselectvars))
  #pselectqry <- paste0("SELECT ", paste0(plt., puniqueid))
  pselectqry <- paste0("SELECT ", plt., puniqueid)
  
#  if (popType %in% c("GRM", "CHNG", "LULC")) {
#    pselectqry <- paste0(pselectqry, "\nUNION\n",
#	     gsub("p.", "pplot.", pselectqry))           
#  }
 
  ## Append nonsamp.filte to pqhereqry
  if (!is.null(nonsamp.pfilter)) {
    if (is.null(pwhereqry)) {
	    pwhereqry <- nonsamp.pfilter
	  } else {
	    pwhereqry <- paste0(pwhereqry, " AND ", nonsamp.pfilter)
	  }
  }

  if (returndata) {

    ## Build select for pltx
    pltxqry <- paste0(pselectqry, 
                      pltfromqry,
                      pwhereqry)
    if (pltindb) {      
      pltx <- tryCatch(
	            data.table(DBI::dbGetQuery(dbconn, pltxqry)),
				error=function(e) {
				warning(e)
  			    return(NULL)}
                )
    } else {
	    pltx <- tryCatch(
	            data.table(sqldf::sqldf(pltxqry)),
				error=function(e) {
				warning(e)
  			    return(NULL)}
                )
    }
    setkeyv(pltx, puniqueid)
  } else {
    pltx <- plotnm
  }

  
  #############################################################################  
  ## Build pwithqry, including any population filters in popFilter
  #############################################################################  
  if ((!is.null(measCur) && measCur) || !is.null(intensity)) {
    palias <- "p"
	  if (is.null(plotnm)) {
   	  pnm <- ppsanm
	  } else {
	    pnm <- plotnm
	  }
	
	  if (measCur) {
	    pfromqry <- getpfromqry(plotCur = TRUE, 
	                        varCur = "MEASYEAR",
	                        Endyr = 2021,
	                        intensity = intensity,
	                        plotnm = pnm, pjoinid = pjoinid,
	                        dbconn = dbconn)
	
      pwithqry <- getpwithqry(plotCur = TRUE, 
	                        varCur = "MEASYEAR",
	                        Endyr = 2021,
	                        intensity = intensity,
	                        plotnm = pnm, pjoinid = pjoinid,
	                        dbconn = dbconn)
	  
    } else if (!is.null(invyrs)) { 
      pwithqry <- getpwithqry(invyrs = invyrs, 
		                    intensity = intensity,
		                    plotnm = pnm, pjoinid = pjoinid,
		                    dbconn = dbconn)
	  }
  } else { 
    pwithqry <- paste0(pselectqry, 
                       pltfromqry,
                       pwhereqry)
  }
  

  #############################################################################
  ## Return data
  #############################################################################
  returnlst <- list(pltassgnx=pltassgnx, pltassgnid=pltassgnid, pltx=pltx,
    pwithqry=pwithqry, pselectqry=pselectqry, plotnm=plotnm,
		puniqueid=puniqueid, pjoinid=pjoinid, palias=palias, 
		unitvar=unitvar, unitarea=unitareax, unitvar2=unitvar2, areavar=areavar, 
		areaunits=areaunits, unit.action=unit.action, ACI=ACI, 
 		P2POINTCNT=as.data.frame(P2POINTCNT), unitlevels=unitlevels, 
		plotsampcnt=as.data.frame(plotsampcnt), pdoms2keep=pdoms2keep, 
		states=states, invyrs=lapply(invyrs,I), datindb=datindb, dbconn=dbconn)

  if (module == "GB") {
    returnlst$strata <- strata
    if (strata) {
      returnlst$stratcombine <- stratcombine 
      returnlst$stratalut <- stratalutx
      returnlst$strvar <- strvar
      returnlst$nonresp <- nonresp
    }
    if (nonresp) {
      returnlst$RHGlut <- RHGlut
      returnlst$nonresplut <- nonresplut
    }  
  } 
  if (module %in% c("MA", "SA")) {
    returnlst$auxlut <- auxlutx
    returnlst$prednames <- prednames
    returnlst$predfac <- predfac
  }
  if (ACI) {
    returnlst$nfplotsampcnt <- nfplotsampcnt
  }

  return(returnlst)
}