3D Arrays : Selecting values of array B around some bit index based on array A

I have a 3d array of some values of varA (time, lat, lon). I want to access some varB based on the values of varA.

index = varA.>threshold

I know I can access these values in varB with the syntax

selected_B = varB[index]

But what if I want to select the past 30-time sequence before the varA is greater than the threshold.

For example if at an instance varA[t,lat, lon] > threshold. Then I want to access varB[t-30:t, lat, lon]

What is the most efficient way of doing it, or is it just suggested to run a loop to get these values.

solved on slack. The solution was to use a different data layout and a loop.

1 Like

In case you were the one who solved it, could you also post the solution here, for non slack-users? Thanks!

1 Like

The solution was to store a Vector that stores (time, (lat,lon)) and is sorted by time, at which point this is just a simple loop (with a binary search to find the right entry)

2 Likes

Keyed arrays are often convenient in such scenarios:

julia> using AxisKeys, IntervalSets

# your data array: some measurements on a grid of time, lat, and lon:
julia> A = KeyedArray(rand(101, 360, 181), time=0:10:1000, lat=0:359, lon=-90:90)
3-dimensional KeyedArray(NamedDimsArray(...)) with keys:
↓   time ∈ 101-element StepRange{Int64,...}
→   lat ∈ 360-element UnitRange{Int64}
◪   lon ∈ 181-element UnitRange{Int64}
And data, 101×360×181 Array{Float64, 3}:
...

# compute t where your threshold is exceded, here I just set it explicitly
julia> t0 = 500

# select all values for a given lat/lon, for a time interval:
julia> A(time=(t0 - 30)..t0, lat=250, lon=-45)
1-dimensional KeyedArray(NamedDimsArray(...)) with keys:
↓   time ∈ 4-element StepRange{Int64,...}
And data, 4-element view(::Array{Float64, 3}, 51:54, 251, 46) with eltype Float64:
 (500)  0.7403947749562407
 (510)  0.8438983921119444
 (520)  0.5526915123954188
 (530)  0.24840834241774157

Oscar is right in that the for loop is likely the simplest, fastest, and easiest to read method. Also aplvain is right in that there are nice packages that can help.

That said you could use an array comprehension to do this.

julia> A = rand(Int8, 3,4,5); B = rand(Int8, 3, 4,5); C = similar(B, 1);

julia> threshold = 1; lookback = 1
1

julia> findall(>(threshold), A)
33-element Vector{CartesianIndex{3}}:
 CartesianIndex(1, 1, 1)
 CartesianIndex(2, 1, 1)
 CartesianIndex(3, 1, 1)
 CartesianIndex(1, 2, 1)
 CartesianIndex(2, 2, 1)
 CartesianIndex(1, 3, 1)
 CartesianIndex(2, 3, 1)
 CartesianIndex(3, 3, 1)
 CartesianIndex(1, 1, 2)
 CartesianIndex(3, 1, 2)
 CartesianIndex(2, 3, 2)
 CartesianIndex(2, 4, 2)
 CartesianIndex(3, 4, 2)
 CartesianIndex(1, 1, 3)
 CartesianIndex(3, 1, 3)
 CartesianIndex(2, 2, 3)
 CartesianIndex(1, 3, 3)
 CartesianIndex(2, 3, 3)
 CartesianIndex(1, 4, 3)
 CartesianIndex(3, 4, 3)
 CartesianIndex(3, 1, 4)
 CartesianIndex(2, 2, 4)
 CartesianIndex(3, 2, 4)
 CartesianIndex(1, 3, 4)
 CartesianIndex(2, 3, 4)
 CartesianIndex(2, 4, 4)
 CartesianIndex(3, 4, 4)
 CartesianIndex(1, 2, 5)
 CartesianIndex(3, 2, 5)
 CartesianIndex(1, 3, 5)
 CartesianIndex(2, 3, 5)
 CartesianIndex(1, 4, 5)
 CartesianIndex(3, 4, 5)

julia> [view(B,ind[1]-lookback:ind[1], ind[2], ind[3]) for ind in findall(>(threshold), A) if ind[1] > lookback]
22-element Vector{SubArray{Int8, 1, Array{Int8, 3}, Tuple{UnitRange{Int64}, Int64, Int64}, true}}:
 [-33, -36]
 [-36, -29]
 [-68, 99]
 [107, 82]
 [82, -44]
 [-119, -122]
 [127, -16]
 [62, -113]
 [-113, 104]
 [7, 103]
 [-4, -94]
 [-91, -26]
 [-85, -99]
 [-17, -12]
 [100, 40]
 [40, -50]
 [-39, -2]
 [49, 4]
 [4, 0]
 [38, -57]
 [82, 122]
 [14, 0]
Value of `A` and `B`
julia> A
3×4×5 Array{Int8, 3}:
[:, :, 1] =
 71    89   47  -112
 50    24  111   -35
 38  -102   41  -120

[:, :, 2] =
   35  -17  -31  -65
 -115  -55   19  114
   57    0  -71   69

[:, :, 3] =
 114   -33   103    91
 -82    57   114  -125
  39  -103  -109   117

[:, :, 4] =
 -124  -107    63  -119
  -44    86    39    38
   35    63  -112    60

[:, :, 5] =
 -111   21    90    66
  -75  -30    25  -103
  -79  102  -116    53

julia> B
3×4×5 Array{Int8, 3}:
[:, :, 1] =
 -33  -68  107  -48
 -36   99   82  108
 -29   55  -44  -91

[:, :, 2] =
  -57   91  127    62
 -119  -93  -16  -113
 -122  -85   13   104

[:, :, 3] =
 -21    -4  -91  -125
   7   -94  -26   -85
 103  -124   16   -99

[:, :, 4] =
  38  100  -39  49
 -17   40   -2   4
 -12  -50   87   0

[:, :, 5] =
    9  -20   82  111
 -103   38  122   14
  -19  -57  -35    0

The idea is to filter the array elements that are below the threshold and then take the n largest of these.
For this latter task, I make use of a homemade function Nlargest() that came up in an old discussion.
The equivalent nlargest function exported by the module could obviously be used

using DataStructures


function Nlargest(v,N)
  maxn = heapify!(v[1:N])
  maxn1=maxn[1]
     @inbounds for i in N+1:length(v)
      e=v[i]    
      if maxn1 < e
          heappop!(maxn)
          heappush!(maxn,e)
          maxn1=maxn[1]
          end
      end
  sort!(maxn,rev=true)
end


A=[(i, (rand(0:359), rand(0:180))) for i in 1:100]


struct  Fatto
  coordinates :: Tuple{Int,Int}
  time :: Int
end

fatti=[Fatto(p, t) for (t,p) in A]

Base.isless(a::Fatto,b::Fatto) = isless(a.time,b.time)

thr=47

fA=filter(f->f.time <=(thr),fatti)

Nlargest(fA,10)