From 67d5f2ad90a71fcb838ebcd7688ec01645fd6b1c Mon Sep 17 00:00:00 2001 From: George Datseris Date: Thu, 4 Nov 2021 12:00:26 +0100 Subject: [PATCH] New `sametimespan` function and some bugfixes. (#76) * fix ambiguity bug in zonal mean * impelemnt tropics_extratropics for Coord data * correct sametimespan function * aDD TESTS for sametimespan * bump version * fix rtests --- Project.toml | 2 +- docs/src/index.md | 1 + src/ClimateBase.jl | 2 + src/exports.jl | 12 ++++++ src/physical_dimensions/spatial_equalarea.jl | 32 ++++++++++++--- src/physical_dimensions/temporal.jl | 42 ++++++++++++++------ test/runtests.jl | 12 ++++++ 7 files changed, 84 insertions(+), 19 deletions(-) create mode 100644 src/exports.jl diff --git a/Project.toml b/Project.toml index e30060f7..e331db1b 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "ClimateBase" uuid = "35604d93-0fb8-4872-9436-495b01d137e2" authors = ["Datseris ", "Philippe Roy "] -version = "0.13.6" +version = "0.13.7" [deps] Dates = "ade2ca70-3891-5945-98fb-dc099432e06a" diff --git a/docs/src/index.md b/docs/src/index.md index 771901f9..624d9916 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -162,6 +162,7 @@ temporal_sampling realtime_days realtime_milliseconds seasonality +sametimespan ``` ## Spatial diff --git a/src/ClimateBase.jl b/src/ClimateBase.jl index 61fcbfe3..8b95eb8c 100644 --- a/src/ClimateBase.jl +++ b/src/ClimateBase.jl @@ -17,4 +17,6 @@ include("climate/albedo.jl") include("tsa/continuation.jl") include("tsa/decomposition.jl") +include("exports.jl") + end # module diff --git a/src/exports.jl b/src/exports.jl new file mode 100644 index 00000000..42bb4235 --- /dev/null +++ b/src/exports.jl @@ -0,0 +1,12 @@ +# This file contains the exports for ClimateBase module. +# Having all exports at a single file is overall cleaner +# and leads to easier-to-identify public API. + +# TODO: Actually move all exports here. + +# Temporal +export monthday_indices, maxyearspan, daymonth, realtime_days, realtime_milliseconds, +temporal_sampling, timemean, timeagg, monthlyagg, yearlyagg, temporalrange, seasonalyagg, +season, DAYS_IN_ORBIT, HOURS_IN_ORBIT, seasonality, sametimespan + +# diff --git a/src/physical_dimensions/spatial_equalarea.jl b/src/physical_dimensions/spatial_equalarea.jl index 33dfb8c9..16e87863 100644 --- a/src/physical_dimensions/spatial_equalarea.jl +++ b/src/physical_dimensions/spatial_equalarea.jl @@ -27,7 +27,7 @@ function zonalmean(::UnstructuredGrid, A::AbDimArray{T, 1}, ::Nothing) where {T} end # zonal mean with weights -function zonalmean(::UnstructuredGrid, A::ClimArray, W) +function zonalmean(::UnstructuredGrid, A::ClimArray, W::AbstractArray) @assert size(A) == size(W) idxs, lats = uniquelats(A) other = otherdims(A, Coord()) @@ -40,7 +40,7 @@ function zonalmean(::UnstructuredGrid, A::ClimArray, W) end return R end -function zonalmean(::UnstructuredGrid, A::ClimArray{T, 1}, W) where {T} +function zonalmean(::UnstructuredGrid, A::ClimArray{T, 1}, W::AbstractArray) where {T} @assert size(A) == size(W) idxs, lats = uniquelats(A) res = zeros(T, length(lats)) @@ -52,7 +52,7 @@ end """ - uniquelats(A::AbDimArray) → idxs, lats + uniquelats(A::ClimArray) → idxs, lats uniquelats(c::Vector{<:AbstractVector}) → idxs, lats Find the unique latitudes of `A`. Return the indices (vector of ranges) that each latitude in `lats` covers, as well as the latitudes themselves. @@ -121,13 +121,33 @@ Return the indices of coordinates belonging to the north and south hemispheres. function hemisphere_indices(c) idxs, lats = uniquelats(c) i = findfirst(x -> x > 0, lats) - shi = idxs[1][1]:idxs[i][end] - nhi = idxs[i+1][1]:idxs[end][end] + shi = idxs[1][1]:idxs[i-1][end] + nhi = idxs[i][1]:idxs[end][end] return nhi, shi end latitudes(::UnstructuredGrid, A) = unique!([x[2] for x in dims(A, Coord)]) +function tropics_extratropics(::UnstructuredGrid, A; lower_lat=30) + # TODO: Support `higher_lat` keyword + c = dims(A, Coord).val + idxs, lats = uniquelats(c) + i1 = findlast(x -> x < -lower_lat, lats) + i2 = findfirst(x -> x > lower_lat, lats) + # tropics indices (accounting for hemispheric only data as well) + t1 = isnothing(i1) ? 0 : i1 + t2 = isnothing(i2) ? length(idxs)+1 : i2 + i_tropics = idxs[t1+1][1]:idxs[t2-1][end] + # extratropics indices (accounting for hemispheric only data as well) + i_sh_extra = isnothing(i1) ? (1:0) : idxs[1][1]:idxs[i1][end] + i_nh_extra = isnothing(i2) ? (1:0) : idxs[i2][1]:idxs[end][end] + i_extra = vcat(i_sh_extra, i_nh_extra) + + tropics = A[Coord(i_tropics)] + extratropics = A[Coord(i_extra)] + return tropics, extratropics +end + ######################################################################### # Extention of convenience indexing of `Coord` ######################################################################### @@ -145,4 +165,4 @@ end function DimensionalData.sel2indices(c::Coord, sel::Tuple{<:Lat{ <: Between}}) l1, l2 = sel[1].val.val return coord_latitudes_between(c.val, l1, l2) # this is Vector{Int} -end +end \ No newline at end of file diff --git a/src/physical_dimensions/temporal.jl b/src/physical_dimensions/temporal.jl index 910db7e3..3a7dc648 100644 --- a/src/physical_dimensions/temporal.jl +++ b/src/physical_dimensions/temporal.jl @@ -2,12 +2,6 @@ Handling of time in data as a physical quantity, and time-related data processing =# using Statistics, StatsBase -export monthday_indices, maxyearspan, daymonth, realtime_days, realtime_milliseconds -export temporal_sampling -export timemean, timeagg -export monthlyagg, yearlyagg, temporalrange, seasonalyagg, season -export DAYS_IN_ORBIT, HOURS_IN_ORBIT -export seasonality ######################################################################### # Datetime related ######################################################################### @@ -62,12 +56,12 @@ function temporal_sampling(t::AbstractVector{<:TimeType}) return :other end end -temporal_sampling(t::AbstractVector) = :other -temporal_sampling(t::StepRange{<:Any,Month}) = :monthly -temporal_sampling(t::StepRange{<:Any,Year}) = :yearly -temporal_sampling(t::StepRange{<:Any,Day}) = :daily -temporal_sampling(t::StepRange{<:Any,Hour}) = :hourly -temporal_sampling(t::StepRange{<:Any,<:Any}) = :other +temporal_sampling(::AbstractVector) = :other +temporal_sampling(::StepRange{<:Any,Month}) = :monthly +temporal_sampling(::StepRange{<:Any,Year}) = :yearly +temporal_sampling(::StepRange{<:Any,Day}) = :daily +temporal_sampling(::StepRange{<:Any,Hour}) = :hourly +temporal_sampling(::StepRange{<:Any,<:Any}) = :other "return the appropriate subtype of `Dates.Period` or `nothing`." function tsamp2period(tsamp) @@ -214,6 +208,30 @@ end realtime_milliseconds(A) = realtime_milliseconds(dims(A, Ti).val) +""" + sametimespan(Xs...) → Ys +Given several `ClimArray`s, return the same `ClimArray`s but now accessed in the `Time` +dimension so that they all have span the same time interval. + +`sametimespan` also has more intelligent handling of monthly or yearly sampled data. +""" +function sametimespan(Xs...) + mint = maximum(minimum(dims(X, Time).val) for X in Xs) + maxt = minimum(maximum(dims(X, Time).val) for X in Xs) + # Make an intelligent decision for monthly/yearly sampled data + tsamps = temporal_sampling.(Xs) + if all(isequal(:monthly), tsamps) + mint = Date(year(mint), month(mint), 1) + d = daysinmonth(maxt) + maxt = Date(year(maxt), month(maxt), d) + elseif all(isequal(:yearly), tsamps) + mint = Date(year(mint), 1, 1) + maxt = Date(year(maxt), 12, 31) + end + map(X -> X[Time(Between(mint, maxt))], Xs) +end +sametimespan(Xs::Tuple) = sametimespan(Xs...) + ######################################################################### # temporal statistics ######################################################################### diff --git a/test/runtests.jl b/test/runtests.jl index 7cbcf935..9c398ea0 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -131,6 +131,18 @@ end @test temporal_sampling(tother) == :other end + @testset "sametimespan" begin + tm1 = Date(2000, 3, 15):Month(1):Date(2010, 3, 15) + tm2 = Date(2001, 1, 1):Month(1):Date(2011, 1, 1) + A1 = ClimArray(rand(length(tm1)), (Time(tm1),)) + A2 = ClimArray(rand(length(tm2)), (Time(tm2),)) + B1, B2 = sametimespan(A1, A2) + @test size(B1) == size(B2) + @test dims(B1, Time)[1] == Date(2001, 1, 15) + @test dims(B1, Time)[end] == Date(2010, 3, 15) + @test dims(B2, Time)[end] == Date(2010, 3, 1) + end + @testset "monthlyagg and co." begin # test yearly temporal weights (i.e., no weighting) X = ClimArray(rand(3,3), (Lon(1:3), Time(tyearly[1:3])))