From 3e8f0301d3a09d504493a6a2306bbe30b3b2f2b9 Mon Sep 17 00:00:00 2001 From: Nathanael Wong Date: Mon, 9 Sep 2024 01:31:49 -0400 Subject: [PATCH] Splitting up functions between rectilinear and generalized grids --- src/extract/generalized.jl | 65 ++++++++++++++++++ src/{extract.jl => extract/rectilinear.jl} | 8 +-- src/grid/generalized.jl | 79 +++++++++++++++++++++ src/{grid.jl => grid/rectilinear.jl} | 80 ---------------------- 4 files changed, 148 insertions(+), 84 deletions(-) create mode 100644 src/extract/generalized.jl rename src/{extract.jl => extract/rectilinear.jl} (97%) create mode 100644 src/grid/generalized.jl rename src/{grid.jl => grid/rectilinear.jl} (82%) diff --git a/src/extract/generalized.jl b/src/extract/generalized.jl new file mode 100644 index 0000000..ef95e43 --- /dev/null +++ b/src/extract/generalized.jl @@ -0,0 +1,65 @@ +""" + extract( + odata :: AbstractArray{<:Real}, + ggrd :: RegionGrid + ) -> 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. + +!!! 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. + +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 +""" +function extract( + odata :: AbstractArray{<:Real}, + ggrd :: GeneralizedGrid +) + + mask = ggrd.mask + ndata = zeros(size(mask)) + for ii in eachindex(mask) + ndata[ii] = odata[ii] * mask[ii] + end + + return ndata + +end + +""" + extract!( + odata :: AbstractArray{<:Real}, + ndata :: AbstractArray{<:Real}, + ggrd :: RegionGrid + ) -> 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. + +This allows for iterable in-place modification to save memory space and reduce allocations if the dimensions are fixed. + +!!! 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. + +Arguments +========= +- `odata` : An array of dimension N containing gridded data in the region we are interesting in extracting from +- `ndata` : An array of dimension N meant as a placeholder for the extracted data to minimize allocations +- `ggrd` : A `RegionGrid` containing detailed information on what to extract +""" +function extract!( + ndata :: AbstractArray{<:Real}, + odata :: AbstractArray{<:Real}, + ggrd :: GeneralizedGrid +) + + mask = ggrd.mask + for ii in eachindex(mask) + ndata[ii] = odata[ii] * mask[ii] + end + + return + +end \ No newline at end of file diff --git a/src/extract.jl b/src/extract/rectilinear.jl similarity index 97% rename from src/extract.jl rename to src/extract/rectilinear.jl index 65e659b..b8ed5b7 100644 --- a/src/extract.jl +++ b/src/extract/rectilinear.jl @@ -16,7 +16,7 @@ Arguments """ function extract( odata :: AbstractArray{<:Real,2}, - ggrd :: RegionGrid + ggrd :: RectilinearGrid ) ilon = ggrd.ilon; nlon = length(ggrd.ilon) @@ -54,7 +54,7 @@ Arguments function extract!( ndata :: AbstractArray{<:Real,2}, odata :: AbstractArray{<:Real,2}, - ggrd :: RegionGrid + ggrd :: RectilinearGrid ) ilon = ggrd.ilon; nlon = length(ggrd.ilon) @@ -89,7 +89,7 @@ end function extract!( ndata :: AbstractArray{<:Real,3}, odata :: AbstractArray{<:Real,3}, - ggrd :: RegionGrid + ggrd :: RectilinearGrid ) ilon = ggrd.ilon; nlon = length(ggrd.ilon) @@ -126,7 +126,7 @@ end function extract!( ndata :: AbstractArray{<:Real,4}, odata :: AbstractArray{<:Real,4}, - ggrd :: RegionGrid + ggrd :: RectilinearGrid ) ilon = ggrd.ilon; nlon = length(ggrd.ilon) diff --git a/src/grid/generalized.jl b/src/grid/generalized.jl new file mode 100644 index 0000000..2e3e9cd --- /dev/null +++ b/src/grid/generalized.jl @@ -0,0 +1,79 @@ +""" + RegionGrid( + geo :: GeoRegion, + lon :: Array{<:Real,2}, + lat :: Array{<:Real,2} + ) -> RegionGrid + +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 +""" +function RegionGrid( + geo :: GeoRegion, + lon :: Array{<:Real,2}, + lat :: Array{<:Real,2} +) + + @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)) + + 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 + end + end + + return RegionMask{FT}(lon,lat,mask,wgts) + +end + +function VectorGrid( + geo :: GeoRegion, + lon :: Vector{<:Real}, + lat :: Vector{<: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 length(lon) != length(lat) + error("$(modulelog()) - The size of the longitude and latitude arrays are not the same.") + end + + mask = zeros(length(lon)) + + for ii in eachindex(lon) + ipnt = Point2(lon[ii],lat[ii]) + if in(ipnt,geo) + mask[ii] = 1 + else; mask[ii] = NaN + end + end + jj = .!isnan.(mask) + + return VectorMask{FT}(lon,lat,mask,lon[jj],lat[jj]) + +end \ No newline at end of file diff --git a/src/grid.jl b/src/grid/rectilinear.jl similarity index 82% rename from src/grid.jl rename to src/grid/rectilinear.jl index a3e8205..32124d1 100644 --- a/src/grid.jl +++ b/src/grid/rectilinear.jl @@ -34,86 +34,6 @@ RegionGrid( geo::PolyRegion, lon::AbstractRange{<:Real}, lat::AbstractRange{<:Real} ) = PolyGrid(geo,collect(lon),collect(lat)) -""" - RegionGrid( - geo :: GeoRegion, - lon :: Array{<:Real,2}, - lat :: Array{<:Real,2} - ) -> RegionGrid - -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 -""" -function RegionGrid( - geo :: GeoRegion, - lon :: Array{<:Real,2}, - lat :: Array{<:Real,2} -) - - @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)) - - 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 - end - end - - return RegionMask{FT}(lon,lat,mask,wgts) - -end - -function VectorGrid( - geo :: GeoRegion, - lon :: Vector{<:Real}, - lat :: Vector{<: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 length(lon) != length(lat) - error("$(modulelog()) - The size of the longitude and latitude arrays are not the same.") - end - - mask = zeros(length(lon)) - - for ii in eachindex(lon) - ipnt = Point2(lon[ii],lat[ii]) - if in(ipnt,geo) - mask[ii] = 1 - else; mask[ii] = NaN - end - end - jj = .!isnan.(mask) - - return VectorMask{FT}(lon,lat,mask,lon[jj],lat[jj]) - -end - function RectGrid( geo :: RectRegion, lon :: Vector{<:Real},