R/findPointXYindex.R

In SantanderMetGroup/loadeR: A climate4R package for data access <http://meteo.unican.es/climate4r>

#' @title Define user-defined geolocation parameters
#' @description Finds the XY index position of single point selections
#' 
#' @param lonLim x domain definition
#' @param latLim y domain definition
#' @param gcs the grid coordinate system, formally a \sQuote{GridCoordSys} java object
#' @return a list with the new lonLim and latLim parameters (if changed) and the pointXYindex
#'  definition, consisting of a integer vector of length two.
#' @details The function applies the java method \sQuote{findXYindexFromCoord}, being the input
#' data the point coordinates in the native units. Thus, if longitudes are coded in the form [0-360] in the dataset,
#' the input lonLim coordinate [-180 - 180] is converted prior to passing it to the method.
#'  The resulting pointXYindex is ultimately passed to the java method \sQuote{readDataSlice} after subsetting.
#' The function also takes into account the possibility of selection of \dQuote{strips}, i.e.,
#' a selection across a whole or part of a meridian/parallel.
#' @author J Bedia 
#' @references \url{http://www.unidata.ucar.edu/software/thredds/v4.3/netcdf-java/v4.3/javadocAll/ucar/nc2/dt/grid/GridCoordSys.html#findXYindexFromLatLon\%28double,\%20double,\%20int[]\%29}
#' @keywords internal
#' @export
#' @importFrom rJava .jnull


findPointXYindex <- function(lonLim, latLim, gcs, spatialTolerance = NULL)  {
  pointXYindex <- c(-1L, -1L)
  bboxDataset <- gcs$getLatLonBoundingBox()
  if (any(lonLim < 0) & bboxDataset$getLonMin() >= 0) {
    lon.aux <- sort(lonLim[which(lonLim < 0)] + 360)
  } else {
    if (!is.null(lonLim)) {
      lon.aux <- lonLim
    } else {
      lon.aux <- bboxDataset$getCenterLon()
    }
  }
  if (length(lonLim) == 1) {
    if (is.null(latLim)) {
      lat.aux <- bboxDataset$getLatMin()
    } else {
      lat.aux <- latLim[1]
    }
    if (!is.null(spatialTolerance)){
      if ((lat.aux[1] > bboxDataset$getLatMax()) & (lat.aux[1] - spatialTolerance <= bboxDataset$getLatMax())){
        lat.aux[1] <- bboxDataset$getLatMax()
      } else if ((lat.aux[1] < bboxDataset$getLatMin()) & (lat.aux[1] + spatialTolerance >= bboxDataset$getLatMin())){
        lat.aux[1] <- bboxDataset$getLatMin()
      } 
      if ((lon.aux[1] > bboxDataset$getLonMax()) & (lon.aux[1] - spatialTolerance <= bboxDataset$getLonMax())){
        lon.aux[1] <- bboxDataset$getLonMax()
      } else if ((lon.aux[1] < bboxDataset$getLonMin()) & (lon.aux[1] + spatialTolerance >= bboxDataset$getLonMin())){
        lon.aux[1] <- bboxDataset$getLonMin()
      } 
    }
    pointXYindex[1] <- gcs$findXYindexFromCoord(lon.aux, lat.aux, .jnull())[1]
    if (pointXYindex[1] < 0){
      proj <- gcs$getProjection()
      if (!proj$isLatLon() | proj$getName() == "LambertConformal") {
        pointXYindex[1] <- gcs$findXYindexFromLatLon(lat.aux, lon.aux, .jnull())[1]
      }
    }
    if ((pointXYindex[1] < 0) & (bboxDataset$getLonMin() < bboxDataset$getLonMax())) {
      stop("Selected X point coordinate is out of range")
    }
    lonLim <- NULL   
  }
  if (length(latLim) == 1) {
    if (!is.null(spatialTolerance)){
      if ((latLim[1] > bboxDataset$getLatMax()) & (latLim[1] - spatialTolerance <= bboxDataset$getLatMax())){
        latLim[1] <- bboxDataset$getLatMax()
      } else if ((latLim[1] < bboxDataset$getLatMin()) & (latLim[1] + spatialTolerance >= bboxDataset$getLatMin())){
        latLim[1] <- bboxDataset$getLatMin()
      } 
      if ((lon.aux[1] > bboxDataset$getLonMax()) & (lon.aux[1] - spatialTolerance <= bboxDataset$getLonMax())){
        lon.aux[1] <- bboxDataset$getLonMax()
      } else if ((lon.aux[1] < bboxDataset$getLonMin()) & (lon.aux[1] + spatialTolerance >= bboxDataset$getLonMin())){
        lon.aux[1] <- bboxDataset$getLonMin()
      } 
    }
    pointXYindex[2] <- gcs$findXYindexFromCoord(lon.aux[1], latLim, .jnull())[2]
    if (pointXYindex[2] < 0){
      proj <- gcs$getProjection()
      if (!proj$isLatLon() | proj$getName() == "LambertConformal") {
        pointXYindex[2] <- gcs$findXYindexFromLatLon(latLim, lon.aux[1], .jnull())[2]
      }
    }
    if ((pointXYindex[2] < 0) & (bboxDataset$getLatMin() < bboxDataset$getLatMax())) {
      stop("Selected Y point coordinate is out of range")
    }
    latLim <- NULL
  }
  return(list("lonLim" = lonLim, "latLim" = latLim, "pointXYindex" = pointXYindex))
}
# End