From 09c3a8419304463381057d4263e77b7c3ef7ab7d Mon Sep 17 00:00:00 2001 From: Nathanael Wong Date: Thu, 3 Oct 2024 13:10:32 -0400 Subject: [PATCH] Updates to functionality and to docstrings --- src/RegionGrids.jl | 44 +++++++++------- src/extract/generalized.jl | 8 +-- src/extract/rectilinear.jl | 102 +++++-------------------------------- src/grid/generalized.jl | 101 ++++++++++++++++++++---------------- src/grid/rectilinear.jl | 59 ++++++++------------- 5 files changed, 122 insertions(+), 192 deletions(-) diff --git a/src/RegionGrids.jl b/src/RegionGrids.jl index 7e420dd..ac8ce61 100644 --- a/src/RegionGrids.jl +++ b/src/RegionGrids.jl @@ -2,10 +2,8 @@ module RegionGrids ## Modules Used using Dates -using GeometryBasics using GeoRegions using Logging -using PolygonOps import Base: show import GeoRegions: isgridinregion @@ -45,21 +43,10 @@ All `RectilinearGrid` types contain the following fields: """ abstract type RectilinearGrid <: RegionGrid end -""" - GeneralizedGrid - -A `GeneralizedGrid` is a `RegionGrid` that is created based on longitude/latitude grids that are **not** rectilinear - this can range from curvilinear grids to unstructured grids. It has its own subtypes: `RegionMask` and `VectorMask`. - -All `GeneralizedGrid` type will contain the following fields: -* `lon` - An array of `Float`s, defining longitude points -* `lat` - An array of `Float`s, defining latitude points -* `mask` - An array of NaNs and 1s, defining the region within the original field in which data points are valid -* `weights` - An array of `Float`s, defining the latitude-weights of each valid point in the grid. -""" -abstract type GeneralizedGrid <: RegionGrid end - """ RectGrid <: RegionGrid + +Information on a `RectilinearGrid` type that is extracted based on a `RectRegion` type. """ struct RectGrid{FT<:Real} <: RectilinearGrid grid :: Vector{Int} @@ -73,6 +60,8 @@ end """ PolyGrid <: RegionGrid + +Information on a `RectilinearGrid` type that is extracted based on a `PolyRegion` type. """ struct PolyGrid{FT<:Real} <: RectilinearGrid grid :: Vector{Int} @@ -87,7 +76,9 @@ end """ TiltGrid <: RegionGrid -A `TiltGrid` type will also contain the following fields: +Information on a `RectilinearGrid` type that is extracted based on a `TiltRegion` type. + +In addition to all the fields common to the `RegionGrid` abstract type, `TiltGrid`s type will also contain the following fields: * `rotX` - A vector of `Float`s, defining indices of the parent longitude vector describing the region * `rotY` - A vector of `Float`s, defining indices of the parent latitude vector describing the region """ @@ -103,8 +94,23 @@ struct TiltGrid{FT<:Real} <: RectilinearGrid rotY :: Array{FT,2} end +""" + GeneralizedGrid + +A `GeneralizedGrid` is a `RegionGrid` that is created based on longitude/latitude grids that are **not** rectilinear - this can range from curvilinear grids to unstructured grids. It has its own subtypes: `RegionMask` and `VectorMask`. + +All `GeneralizedGrid` type will contain the following fields: +* `lon` - An array of `Float`s, defining longitude points +* `lat` - An array of `Float`s, defining latitude points +* `mask` - An array of NaNs and 1s, defining the region within the original field in which data points are valid +* `weights` - An array of `Float`s, defining the latitude-weights of each valid point in the grid. +""" +abstract type GeneralizedGrid <: RegionGrid end + """ RegionMask <: GeneralizedGrid + +Information on a `GeneralizedGrid` type that is extracted based on arrays of longitude/latitude points. """ struct RegionMask{FT<:Real} <: GeneralizedGrid lon :: Array{FT,2} @@ -116,9 +122,11 @@ end """ VectorMask <: GeneralizedGrid +Information on a `GeneralizedGrid` type that is extracted based on vectors of longitude and latitude points. + A `VectorMask` type will also contain the following fields: -* `olon` - A vector of `Float`s, defining indices of the parent longitude vector describing the region -* `olat` - A vector of `Float`s, defining indices of the parent latitude vector describing the region +* `olon` - A vector of `Float`s, defining indices of the original longitude vector describing the region +* `olat` - A vector of `Float`s, defining indices of the original latitude vector describing the region """ struct VectorMask{FT<:Real} <: GeneralizedGrid lon :: Vector{FT} diff --git a/src/extract/generalized.jl b/src/extract/generalized.jl index a31c8b8..85448f2 100644 --- a/src/extract/generalized.jl +++ b/src/extract/generalized.jl @@ -5,19 +5,19 @@ crop :: Bool = false ) -> Array{<:Real} -Extracts data from odata, an Array of dimension N (where N > 2) that contains data of a Parent `GeoRegion`, into another Array of dimension N. +Using a `RegionMask` for a region of interest, retrieve valid data - all other data outside the region of interest !!! warning "Dimension Order" Please ensure that the 1st dimension is longitude and 2nd dimension is latitude. Arguments ========= -- `odata` : An array of dimension N containing gridded data in the region we are interesting in extracting from -- `ggrd` : A `RegionGrid` containing detailed information on what to extract +- `odata` : An array of containing the original gridded data. +- `ggrd` : A `RegionMask` type containing detailed information on the ranges of valid longitude/latitude data to be extracted. Keyword Arguments ================= -- `crop` : If `true`, crop new array to fit extracted data to save size +- `crop` : If `true`, crop new array to valid data to save size. """ function extract( odata :: AbstractArray{<:Real}, diff --git a/src/extract/rectilinear.jl b/src/extract/rectilinear.jl index 58def8b..16e4409 100644 --- a/src/extract/rectilinear.jl +++ b/src/extract/rectilinear.jl @@ -1,10 +1,10 @@ """ - extractGrid( - odata :: AbstractArray{<:Real,N}, + extract( + odata :: AbstractArray{<:Real}, ggrd :: RegionGrid - ) where N <: Int -> Array{<:Real,N} + ) -> Array{<:Real} -Extracts data from odata, an Array of dimension N (where N ∈ 2,3,4) that contains data of a Parent `GeoRegion`, into another Array of dimension N, containing _**only**_ within a sub `GeoRegion` we are interested in. +Extracts data from odata, an Array that contains data of a Parent `GeoRegion`, into another Array of dimension N, containing _**only**_ within a sub `GeoRegion` we are interested in. !!! warning Please ensure that the 1st dimension is longitude and 2nd dimension is latitude before proceeding. The order of the 3rd and 4th dimensions (when used), however, is not significant. @@ -15,9 +15,9 @@ Arguments - `ggrd` : A `RegionGrid` containing detailed information on what to extract """ function extract( - odata :: AbstractArray{<:Real}, - ggrd :: RectilinearGrid -) + odata :: AbstractArray{FT}, + ggrd :: RegionGrid +) where FT <: Real ilon = ggrd.ilon; nlon = length(ggrd.ilon) ilat = ggrd.ilat; nlat = length(ggrd.ilat) @@ -42,11 +42,11 @@ function extract( end """ - extractGrid!( - odata :: AbstractArray{<:Real,N}, - ndata :: AbstractArray{<:Real,N}, + extract!( + odata :: AbstractArray{<:Real}, + ndata :: AbstractArray{<:Real}, ggrd :: RegionGrid - ) where N <: Int -> nothing + ) -> nothing Extracts data from odata, an Array of dimension N (where N ∈ 2,3,4) that contains data of a Parent `GeoRegion`, into ndata, another Array of dimension N, containing _**only**_ within a sub `GeoRegion` we are interested in. @@ -62,9 +62,9 @@ Arguments - `ggrd` : A `RegionGrid` containing detailed information on what to extract """ function extract!( - ndata :: AbstractArray{<:Real}, odata :: AbstractArray{<:Real}, - ggrd :: RectilinearGrid + ndata :: AbstractArray{<:Real}, + ggrd :: RegionGrid ) ilon = ggrd.ilon; nlon = length(ggrd.ilon) @@ -85,78 +85,4 @@ function extract!( return -end - -# function extract( -# odata :: AbstractArray{<:Real,3}, -# ggrd :: RegionGrid -# ) - -# ilon = ggrd.ilon; nlon = length(ggrd.ilon) -# ilat = ggrd.ilat; nlat = length(ggrd.ilat) -# n3D = size(odata,3) -# mask = ggrd.mask -# ndata = zeros(nlon,nlat,n3D) -# for i3D = 1 : n3D, glat in 1 : nlat, glon in 1 : nlon -# ndata[glon,glat,i3D] = odata[ilon[glon],ilat[glat],i3D] * mask[glon,glat] -# end - -# return ndata - -# end - -# function extract!( -# ndata :: AbstractArray{<:Real,3}, -# odata :: AbstractArray{<:Real,3}, -# ggrd :: RectilinearGrid -# ) - -# ilon = ggrd.ilon; nlon = length(ggrd.ilon) -# ilat = ggrd.ilat; nlat = length(ggrd.ilat) -# n3D = size(odata,3) -# mask = ggrd.mask -# for i3D = 1 : n3D, glat in 1 : nlat, glon in 1 : nlon -# ndata[glon,glat,i3D] = odata[ilon[glon],ilat[glat],i3D] * mask[glon,glat] -# end - -# return - -# end - -# function extract( -# odata :: AbstractArray{<:Real,4}, -# ggrd :: RegionGrid -# ) - -# ilon = ggrd.ilon; nlon = length(ggrd.ilon) -# ilat = ggrd.ilat; nlat = length(ggrd.ilat) -# n3D = size(odata,3) -# n4D = size(odata,4) -# mask = ggrd.mask -# ndata = zeros(nlon,nlat,n3D,n4D) -# for i4D = 1 : n4D, i3D = 1 : n3D, glat in 1 : nlat, glon in 1 : nlon -# ndata[glon,glat,i3D,i4D] = odata[ilon[glon],ilat[glat],i3D,i4D] * mask[glon,glat] -# end - -# return ndata - -# end - -# function extract!( -# ndata :: AbstractArray{<:Real,4}, -# odata :: AbstractArray{<:Real,4}, -# ggrd :: RectilinearGrid -# ) - -# ilon = ggrd.ilon; nlon = length(ggrd.ilon) -# ilat = ggrd.ilat; nlat = length(ggrd.ilat) -# n3D = size(odata,3) -# n4D = size(odata,4) -# mask = ggrd.mask -# for i4D = 1 : n4D, i3D = 1 : n3D, glat in 1 : nlat, glon in 1 : nlon -# ndata[glon,glat,i3D,i4D] = odata[ilon[glon],ilat[glat],i3D,i4D] * mask[glon,glat] -# end - -# return - -# end \ No newline at end of file +end \ No newline at end of file diff --git a/src/grid/generalized.jl b/src/grid/generalized.jl index 2e3e9cd..a4f5216 100644 --- a/src/grid/generalized.jl +++ b/src/grid/generalized.jl @@ -1,45 +1,42 @@ """ RegionGrid( - geo :: GeoRegion, - lon :: Array{<:Real,2}, - lat :: Array{<:Real,2} - ) -> RegionGrid + geo :: GeoRegion, + pnts :: Array{Point2{FT}} + ) where FT <: Real -> gmsk :: RegionMask Creates a `RegionMask` type based on the following arguments. This method is more suitable for non-rectilinear grids of longitude and latitude points, such as in the output of WRF or CESM. Arguments ========= - - `geo` : A GeoRegion of interest -- `lon` : An array containing the longitude points -- `lat` : An array containing the latitude points +- `pnts` : An array containing the longitude points + +Returns +======= +- `gmsk` : A `RegionMask` """ function RegionGrid( - geo :: GeoRegion, - lon :: Array{<:Real,2}, - lat :: Array{<:Real,2} -) + geo :: GeoRegion, + pnts :: Array{Point2{FT}} +) where FT <: Real @info "$(modulelog()) - Creating a RegionMask for the $(geo.name) GeoRegion based on an array of longitude and latitude points" - if eltype(lon) <: AbstractFloat - FT = eltype(lon) - end - - if size(lon) != size(lat) - error("$(modulelog()) - The size of the longitude and latitude arrays are not the same.") - end - - mask = zeros(size(lon)) - wgts = zeros(size(lon)) + npnt = size(pnts) + mask = zeros(npnt) + wgts = zeros(npnt) + lon = zeros(npnt) + lat = zeros(npnt) for ii in eachindex(lon) - ipnt = Point2(lon[ii],lat[ii]) - if in(ipnt,geo) - mask[ii] = 1 - wgts[ii] = cosd.(lat[ii]) - else; mask[ii] = NaN - wgts[ii] = 0 + lon[ii] = pnts[ii][1] + lat[ii] = pnts[ii][2] + if in(pnts[ii],geo) + mask[ii] = 1 + wgts[ii] = cosd.(pnts[ii][2]) + else + mask[ii] = NaN + wgts[ii] = 0 end end @@ -47,33 +44,47 @@ function RegionGrid( end -function VectorGrid( - geo :: GeoRegion, - lon :: Vector{<:Real}, - lat :: Vector{<:Real} -) +""" + VectorGrid( + geo :: GeoRegion, + pnts :: Vector{Point2{FT}} + ) where FT <: Real -> gmsk :: VectorMask - @info "$(modulelog()) - Creating a RegionMask for the $(geo.name) GeoRegion based on an array of longitude and latitude points" +Creates a `VectorMask` type based on a vector of - if eltype(lon) <: AbstractFloat - FT = eltype(lon) - end +Arguments +========= +- `geo` : A GeoRegion of interest +- `pnts` : A `Vector` of `Float` Types, containing the longitude points - if length(lon) != length(lat) - error("$(modulelog()) - The size of the longitude and latitude arrays are not the same.") - end +Returns +======= +- `gmsk` : A `VectorMask` +""" +function VectorGrid( + geo :: GeoRegion, + pnts :: Vector{Point2{FT}} +) where FT <: Real - mask = zeros(length(lon)) + @info "$(modulelog()) - Creating a RegionMask for the $(geo.name) GeoRegion based on an array of longitude and latitude points" - for ii in eachindex(lon) - ipnt = Point2(lon[ii],lat[ii]) - if in(ipnt,geo) + npnt = length(pnts) + mask = zeros(npnt) + wgts = zeros(npnt) + lon = zeros(npnt) + lat = zeros(npnt) + + for ii in 1 : length(pnts) + lon[ii] = pnts[ii][1] + lat[ii] = pnts[ii][2] + if in(pnts[ii],geo) mask[ii] = 1 + wgts[ii] = cosd.(pnts[ii][2]) else; mask[ii] = NaN + wgts[ii] = 0 end end - jj = .!isnan.(mask) - return VectorMask{FT}(lon,lat,mask,lon[jj],lat[jj]) + return VectorMask{FT}(lon,lat,mask,wgts) end \ No newline at end of file diff --git a/src/grid/rectilinear.jl b/src/grid/rectilinear.jl index 32124d1..2a042af 100644 --- a/src/grid/rectilinear.jl +++ b/src/grid/rectilinear.jl @@ -3,16 +3,19 @@ geo :: GeoRegion, lon :: Union{Vector{<:Real},AbstractRange{<:Real}, lat :: Union{Vector{<:Real},AbstractRange{<:Real} - ) -> RectGrid, PolyGrid + ) -> ggrd :: RectilinearGrid Creates a `RectGrid` or `PolyGrid` type based on the following arguments. This method is suitable for rectilinear grids of longitude/latitude output such as from Isca, or from satellite and reanalysis gridded datasets. Arguments ========= - - `geo` : A GeoRegion of interest - `lon` : A vector or `AbstractRange` containing the longitude points - `lat` : A vector or `AbstractRange` containing the latitude points + +Returns +======= +- `ggrd` : A `RectilinearGrid` """ RegionGrid( geo::RectRegion, lon::Vector{<:Real}, lat::Vector{<:Real} @@ -36,23 +39,15 @@ RegionGrid( function RectGrid( geo :: RectRegion, - lon :: Vector{<:Real}, - lat :: Vector{<:Real} -) + lon :: Vector{FT}, + lat :: Vector{FT} +) where FT <: Real @info "$(modulelog()) - Creating a RegionGrid for the $(geo.name) GeoRegion" @debug "$(modulelog()) - Determining indices of longitude and latitude boundaries in the given dataset ..." - N = geo.N - S = geo.S - E = geo.E - W = geo.W - - if eltype(lon) <: AbstractFloat - FT = eltype(lon) - end - igrid = regiongrid([N,S,E,W],lon,lat); + igrid = regiongrid([N(geo),S(geo),E(geo),W(geo)],lon,lat); iN = igrid[1]; iS = igrid[2]; iE = igrid[3]; iW = igrid[4] nlon = deepcopy(lon) @@ -64,7 +59,7 @@ function RectGrid( @info "$(modulelog()) - Creating vector of longitude indices to extract ..." if iW < iE; iWE = vcat(iW:iE) - elseif iW > iE || (iW == iE && E != W) + elseif iW > iE || (iW == iE && E(geo) != W(geo)) iWE = vcat(iW:length(lon),1:iE); nlon[1:(iW-1)] .+= 360 else; iWE = [iW]; end @@ -79,8 +74,8 @@ function RectGrid( wgts[:,ilat] .= cosd(nlat[ilat]) end - while maximum(nlon) > geo.E; nlon .-= 360 end - while minimum(nlon) < geo.W; nlon .+= 360 end + while maximum(nlon) > E(geo); nlon .-= 360 end + while minimum(nlon) < W(geo); nlon .+= 360 end return RectGrid{FT}(igrid,nlon,nlat,iWE,iNS,mask,wgts) @@ -88,21 +83,15 @@ end function TiltGrid( geo :: TiltRegion, - lon :: Vector{<:Real}, - lat :: Vector{<:Real} -) + lon :: Vector{FT}, + lat :: Vector{FT} +) where FT <: Real @info "$(modulelog()) - Creating a RegionGrid for the $(geo.name) GeoRegion" @debug "$(modulelog()) - Determining indices of longitude and latitude boundaries in the given dataset ..." - N,S,E,W = getTiltBounds(geo) - - if eltype(lon) <: AbstractFloat - FT = eltype(lon) - end - - igrid = regiongrid([N,S,E,W],lon,lat); + igrid = regiongrid([N(geo),S(geo),E(geo),W(geo)],lon,lat) iN = igrid[1]; iS = igrid[2]; iE = igrid[3]; iW = igrid[4] nlon = deepcopy(lon) @@ -114,7 +103,7 @@ function TiltGrid( @info "$(modulelog()) - Creating vector of longitude indices to extract ..." if iW < iE; iWE = vcat(iW:iE) - elseif iW > iE || (iW == iE && E != W) + elseif iW > iE || (iW == iE && E(geo) != W(geo)) iWE = vcat(iW:length(lon),1:iE); nlon[1:(iW-1)] .+= 360 else; iWE = [iW]; end @@ -154,8 +143,8 @@ function TiltGrid( end end - while maximum(nlon) > geo.E; nlon .-= 360 end - while minimum(nlon) < geo.W; nlon .+= 360 end + while maximum(nlon) > E(geo); nlon .-= 360 end + while minimum(nlon) < W(geo); nlon .+= 360 end return TiltGrid{FT}(igrid,nlon,nlat,iWE,iNS,mask,wgts,rotX,rotY) @@ -163,9 +152,9 @@ end function PolyGrid( geo :: PolyRegion, - lon :: Vector{<:Real}, - lat :: Vector{<:Real} -) + lon :: Vector{FT}, + lat :: Vector{FT} +) where FT <: Real @info "$(modulelog()) - Creating a RegionGrid for the $(geo.name) GeoRegion" @@ -175,10 +164,6 @@ function PolyGrid( E = geo.E W = geo.W - if eltype(lon) <: AbstractFloat - FT = eltype(lon) - end - igrid = regiongrid([N,S,E,W],lon,lat); iN = igrid[1]; iS = igrid[2]; iE = igrid[3]; iW = igrid[4] nlon = deepcopy(lon)