1using LazyStack
 2using PaddedViews
 3using SparseArrays
 4
 5@views function eigvals_binning(λs; digits::Int=2, δ=10.0^(-digits))
 6    # Minimum number of eigenvalues; to be consistent for the all grid.
 7    N_max = maximum([maximum(length.(λd)) for λd in λs])
 8    λs_restr = lazystack(PaddedView.(NaN, lazystack(λs), Ref((N_max,))))
 9
10    λs_min, λs_max = extrema(x->isnan(x) ? λs_restr[1] : x, λs_restr)
11
12    scale = ceil(λs_min-δ; digits):δ:ceil(λs_max+δ; digits)
13    bin = spzeros(Int, length(scale))
14    for λ in λs_restr
15        index = findfirst(x -> x ≥ λ, scale)
16        bin[index] += 1
17    end
18    (; bin1d=bin, scale, λs=λs_restr, δ)
19end
20
21# Used to consume too much memory before views.
22@views function eigvals_binning_3d(λs; digits::Int=2, δ=10.0^(-digits))
23    #N_max = maximum([maximum(length.(λd)) for λd in λs])
24    # Minimum number of eigenvalues; to be consistent for the all grid.
25    N_max = minimum([minimum(length.(λd)) for λd in λs])
26    λs_restr = lazystack(PaddedView.(NaN, lazystack(λs), Ref((N_max,))))
27
28    λs_min, λs_max = extrema(x->isnan(x) ? λs_restr[1] : x, λs_restr)
29
30    scale = ceil(λs_min-δ; digits):δ:ceil(λs_max+δ; digits)
31    bin = spzeros(Int, length(scale)*size(λs_restr, 2)*size(λs_restr, 3))
32    bin = reshape(bin, length(scale), size(λs_restr, 2), size(λs_restr, 3))
33    bin = zeros(Int, length(scale),size(λs_restr, 2),size(λs_restr, 3))
34
35    for ids in 1:size(λs_restr, 3)
36        for idk in 1:size(λs_restr, 2)
37            for λ in λs_restr[:, idk, ids]
38                index = findfirst(x -> x ≥ λ, scale)
39                bin[index, idk, ids] += 1
40            end
41        end
42    end
43    (; bin, scale, λs=λs_restr, δ, bin1d=vec(sum(bin; dims=(2,3))))
44end
45
46function binning(x, bins)
47    ivalue = findfirst(e -> e ≥ x, bins.scale)
48    return bins.bin1d[ivalue] ./ bins.δ
49end