add precompile for anomaly

Project.toml

@@ -10,7 +10,6 @@ DocStringExtensions = "ffbed154-4ef7-542d-bbb7-c09d3a79fcae"
 Glob = "c27321d9-0574-5035-807b-f59d2c89b15c"
 LambdaFn = "fcefd6c4-8de0-11e9-0d15-03bc7c8ec982"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
-LoopVectorization = "bdcacae8-1622-11e9-2a5c-532679323890"
 MPI = "da04e1cc-30fd-572f-bb4f-1f8673147195"
 MutableNamedTuples = "af6c499f-54b4-48cc-bbd2-094bba7533c7"
 NaNStatistics = "b946abbf-3ea7-4610-9019-9858bfdeaf2d"
@@ -35,7 +34,6 @@ CategoricalArrays = "0.10"
 DocStringExtensions = "0.8, 0.9"
 Glob = "1.3"
 LambdaFn = "0.3"
-LoopVectorization = "0.12"
 MPI = "0.20"
 MutableNamedTuples = "0.1"
 NaNStatistics = "0.6.27"

src/Climate/warming_level.jl

@@ -2,13 +2,16 @@
 Calculate yearly air temperature.
 
 # Description
-we use the fixed thresholds and add the seasonal warming signal. 
-Thus, thresholds are defined as a fixed baseline (such as for the fixed threshold) plus
-seasonally moving mean warming of the corresponding future climate based on the
-31-year moving mean of the warmest three months.
+we use the fixed thresholds and add the seasonal warming signal. Thus,
+thresholds are defined as a fixed baseline (such as for the fixed threshold)
+plus seasonally moving mean warming of the corresponding future climate based on
+the 31-year moving mean of the warmest three months.
 
 # Details
-This function calculates the yearly air temperature based on the input temperature data and dates. If `only_summer` is true, it only calculates the temperature for summer months. The function applies the calculation along the specified dimensions.
+This function calculates the yearly air temperature based on the input
+temperature data and dates. If `only_summer` is true, it only calculates the
+temperature for summer months. The function applies the calculation along the
+specified dimensions.
 
 # Arguments
 - `A::AbstractArray{T,N}`: input array of temperature data.
@@ -57,6 +60,6 @@ function cal_warming_level(A::AbstractArray{T,N}, dates;
 
   T_year_clim = selectdim(T_year, dims, inds_clim)
   T_clim = apply(T_year_clim, dims; fun=nanmean)
-  
+
   T_year .- T_clim
 end

src/Statistics/apply.jl

@@ -56,7 +56,7 @@ function apply(x::AbstractArray, dims=3, args...; by=nothing, fun::Function=mean
         # ans = fun(data, args...; kw...)
         ans = mapslices(fun2, data, dims=dims)
         inds_bad = findall(size(ans) .== 1)
-        length(inds_bad) >= 1 && (ans = dropdims(ans; dims=inds_bad[1]))
+        length(inds_bad) >= 1 && ndims(ans) > 1 && (ans = dropdims(ans; dims=inds_bad[1]))
         ans
       end, grps)
     # permutedims(A, perm)

src/dates.jl

@@ -28,6 +28,6 @@ date_ym(dates) = make_date.(year.(dates), month.(dates))
 
 
 export dates_miss, dates_nmiss,
-    DateTime, year, month, day, Year, Month, Day, format, 
+    DateTime, Date, year, month, day, Year, Month, Day, format, 
     make_datetime, make_date, 
     date_year, date_ym

src/precompile.jl

@@ -27,7 +27,8 @@ precompile(str_extract_strip, (String, String))
 @setup_workload begin
   str = "hello world! hello world!"
   strs = [str, "hello", "world"]
-
+  dates = make_date(1960, 1, 1):Day(1):make_date(1961, 12, 31) |> collect
+
   @compile_workload begin
     str_extract(str, "hello")
     str_extract(strs, "hello")
@@ -61,5 +62,29 @@ precompile(str_extract_strip, (String, String))
     end
 
     @pipe str |> _
+
+    ## precompile for cal_anomaly_quantile
+    function test_anomaly2(; T=Float32, dims=(4,))
+      ntime = length(dates)
+
+      A = rand(T, dims..., ntime)
+      kw = (; parallel=true, p1=1960, p2=1960, na_rm=false, probs=[0.5, 0.9])
+
+      @time anom_season = cal_anomaly_quantile(A, dates; kw..., method="season")
+      @time anom_base = cal_anomaly_quantile(A, dates; kw..., method="base")
+      @time anom_full = cal_anomaly_quantile(A, dates; kw..., method="full")
+
+      @assert size(anom_base) == (dims..., ntime, length(kw.probs))
+      @assert size(anom_base) == size(anom_full)
+      @assert size(anom_base) == size(anom_season)
+    end
+    # set_seed(1)
+    l_dims = [(), (4,), (4, 4), (4, 4, 4)]
+    for T in (Float32, Float64)
+      for dims = l_dims
+        test_anomaly2(; T, dims)
+      end
+    end
+
   end
 end

test/test-stat_anomaly.jl

@@ -12,6 +12,29 @@
 end
 
 
+@testset "_cal_anomaly_3d" begin
+  set_seed(1)
+  ny = 10
+  dates = Date(2010):Day(1):Date(2010 + ny - 1, 12, 31)
+  ntime = length(dates)
+
+  dims = (100, 100)
+  set_seed(1)
+  nprob = ()
+
+  arr = rand(dims..., ntime)
+  TRS = rand(dims..., 366, nprob...)
+  T_wl = rand(dims..., ny)
+
+  r1 = _cal_anomaly(arr, TRS, dates; T_wl)
+  r3 = _cal_anomaly_3d(arr, TRS, dates; T_wl)
+
+  @time r1 = _cal_anomaly(arr, TRS, dates; T_wl)
+  @time r3 = _cal_anomaly_3d(arr, TRS, dates; T_wl)
+  @test r1 == r3
+end
+
+
 @testset "cal_anomaly 3d" begin
   dates = make_date(1961, 1, 1):Day(1):make_date(2000, 12, 31) |> collect
   n = length(dates)