R/check.popdataPLT.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
  ###################################################################################

  ## 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
  ########################################################
  for (tabnm in names(tabs)) {
    assign(tabnm, tabs[[tabnm]])
  }
  puniqueid <- tabIDs[["plt"]]
 
  ###################################################################################
  ## 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)
  }	   
  
  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 <- chkdbtab(tablst, plt)
	  if (is.null(plotnm)) {
	    if (!grepl(plt, "plot", ignore.case=TRUE) && !is.null(chkdbtab(tablst, "plt"))) {
		  plotnm <- chkdbtab(tablst, "plt")
		} else if (!grepl(plt, "plt", ignore.case=TRUE) && !is.null(chkdbtab(tablst, "plot"))) {
		  plotnm <- chkdbtab(tablst, "plot")
        }
      }	
      if (!is.null(plotnm)) {	  
	    pltflds <- DBI::dbListFields(dbconn, plotnm)
	  }
	}
	
    ## Check name of pltassgn
	if (is.null(pltassgnx) && !is.null(pltassgn)) {
	  ppsanm <- chkdbtab(tablst, pltassgn)
      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
	}
  }

  #######################################################################
  ## 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
  }
   
  #######################################################################
  ## Check popevalid pop filter in ppsa and plt
  #######################################################################
  if (!is.null(popevalid)) {	
    evalidnm <- findnm("EVALID", pflds, returnNULL = TRUE)	  
    if (is.null(evalidnm)) {
      message("the EVALID field does not exist in data set...")
	  return(NULL)
	}
		
    ## Check popevalid values in database
	evalida. <- ifelse(evalidnm %in% ppsaflds, pltassgn., plt.)	  
	evalidqry <- paste0("SELECT DISTINCT ", evalida., evalidnm, 
	                    pltfromqry,
						"\nORDER BY ", evalida., evalidnm)
    if (pltindb) {      
      evalidvals <- DBI::dbGetQuery(dbconn, evalidqry)[[1]]
	} else {
	  evalidvals <- sqldf::sqldf(evalidqry)[[1]]
	}
    evalidmiss <- popevalid[!popevalid %in% evalidvals]
    if (any(!popevalid %in% evalidvals)) {
	  message("evalids are missing: ", toString(popevalid[!popevalid %in% evalidvals]))
      return(NULL)
    }
		
	## Build where query to include popevalid popfilter 
	ewhereqry <- paste0(evalida., evalidnm, " IN(", toString(popevalid), ")")
    if (is.null(pwhereqry)) {
      pwhereqry <- paste0("\nWHERE ", ewhereqry)
    } else {
      pwhereqry <- paste0(pwhereqry, 
	                       "\n  AND ", ewhereqry)
    }
  } 

  #######################################################################
  ## Check intensity pop filter in ppsa and plt
  #######################################################################
  if (!is.null(intensity)) { 	   
    intensitynm <- findnm("INTENSITY", pflds, returnNULL = TRUE)
    if (is.null(intensitynm)) {
      message("the INTENSITY field does not exist in data set...")
	  return(NULL)
    }
		
    ## Check intensity values in database
    intensitya. <- ifelse(intensitynm %in% ppsaflds, pltassgn., plt.)	  
    intensityqry <- paste0("SELECT DISTINCT ", intensitya., intensitynm, 
                            pltfromqry, 
						    pwhereqry,
						    "\nORDER BY ", intensitya., intensitynm)
							
    if (pltindb) {      
      intensityvals <- DBI::dbGetQuery(dbconn, intensityqry)[[1]]
    } else {
	  intensityvals <- sqldf::sqldf(intensityqry)[[1]]
    }
    intensitymiss <- intensity[!intensity %in% intensityvals]
    if (any(!intensity %in% intensityvals)) {
	  message("intensity are missing: ", toString(intensity[!intensity %in% intensityvals]))
      return(NULL)
    }
		
    ## Build where query to include invyr popfilter 
    iwhereqry <- paste0(intensitya., intensitynm, " IN(", toString(intensity), ")")
    if (is.null(pwhereqry)) {
      pwhereqry <- paste0("\nWHERE ", iwhereqry)
    } else {
      pwhereqry <- paste0(pwhereqry, 
	                       "\n  AND ", iwhereqry)
    }
  }
  
  ## If Change Plots, remove plots that have no remeasurement data
  ######################################################################################
  if (popType %in% c("GRM", "CHNG", "LULC")) {
    rempernm <- findnm("REMPER", pflds, returnNULL = TRUE)
	if (is.null(rempernm)) {
      message("REMPER is not in dataset... assuming all remeasured plots")
	} else {
	
      ## Build where query to include remper NA removing
	  rempera. <- ifelse(rempernm %in% ppsaflds, pltassgn., plt.)	  

      remper.filter <- paste0(rempera., rempernm, " > 0")
      if (is.null(pwhereqry)) {
        pwhereqry <- paste0("\nWHERE ", remper.filter)
      } else {
        pwhereqry <- paste0(pwhereqry, 
	                       "\n AND ", remper.filter)
      }	
    }
  }

  if (is.null(popevalid)) {

    ## Add subcycle to where statement
    if (!is.null(subcycle99) && !subcycle99) {
	  subcyclenm <- findnm("SUBCYCLE", pflds, returnNULL=TRUE)
	  if (is.null(subcyclenm)) {
	    message("SUBCYCLE not in data... assuming all SUBCYCLE <> 99")
	  } else {
        subcycle.filter <- paste0("p.", subcyclenm, " <> 99")
        if (is.null(pwhereqry)) {
          pwhereqry <- paste0("\nWHERE ", subcycle.filter)
        } else {
          pwhereqry <- paste0(pwhereqry, 
	                       "\n  AND ", subcycle.filter)
        }
	  }
    }

    #######################################################################
    ## Check designcd in ppsa and plt
    #######################################################################
    designcdnm <- findnm("DESIGNCD", pflds, returnNULL = TRUE)
    if (is.null(designcdnm)) {
      message("DESIGNCD is not in dataset... assuming one plot design")
    } else {
		
      ## Check designcd values in database
      designcda. <- ifelse(designcdnm %in% ppsaflds, pltassgn., plt.)	  
      designcdqry <- paste0("SELECT DISTINCT ", designcda., designcdnm, 
                           pltfromqry, 
						   pwhereqry,
						   "\nORDER BY ", designcda., designcdnm)
      if (pltindb) {      
        designcdvals <- DBI::dbGetQuery(dbconn, designcdqry)[[1]]
      } else {
	    designcdvals <- sqldf::sqldf(designcdqry)[[1]]
      }
	
      if (length(designcdvals) > 1) {
        if (any(!designcdvals %in% c(1, 501:505, 230:242, 311:323, 328))) {
          if (adj == "samp") {
            message("designcds include: ", toString(designcdvals))
            message("samp adjustment for trees is only for annual inventory designs... see FIA database manual")
          } else {
            warning("more than 1 plot design, calculate separate estimates by design")
          }
		}
      }
    }
  }  
  
  ###################################################################################
  ## 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... mutliple 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)
}