Accessing every elements of an array except a certain range

Hi all,

I would like to access and modify the contents of an array except for a certain range. For instance, setting all elements of an array to value one, except for those elements whose indices are between 10-20 (like leaving a square within the array untouched). How can I do that in a concise way in Julia?

Thanks in advance

If it’s a 1-dimensional array, a vector, with a single forbidden interval, you could do

julia> a = rand(Int, 20)
20-element Array{Int64,1}:
  9137123603072731695
 -1997914555664844611
 -5259674904209714741
  3157717621575648352
  8597434240321990305
 -4830398965124842278
  8841487105864588126
  1872069802758885653
  7801761719405426744
                    ⋮
  6777382442102949775
  -293044651686397507
  9041568364216053481
   147838476344048725
 -4775248271519179018
  4733728860237532737
 -6728603255363358386
 -5795792524732282326

julia> for i in (union(1:5, 11:20))
   a[i] = 0
end
julia> a
20-element Array{Int64,1}:
                    0
                    0
                    0
                    0
                    0
 -4830398965124842278
  8841487105864588126
  1872069802758885653
  7801761719405426744
                    ⋮
                    0
                    0
                    0
                    0
                    0
                    0
                    0
                    0

julia> 

Perhaps that’s not what you mean?

Thanks! That helped! However, I was mostly referring to 2-d (or n-d) arrays, how to change every element within that array except elements within a certain range (for instance elements that have 10<i<20 and 15<j<25) and if there’s a way to do it without having to use for loop.

Some ways, including a loop:

x = rand(6,10);
@. x = ifelse((2 < $(axes(x,1)) < 5) & (3 < $(axes(x,2)') < 8), x, 0)

using Tullio
x = rand(6,10);
@tullio x[i,j] = ifelse(2<i<5 && 3<j<8, x[i,j], 0)

x = rand(6,10);
for I in CartesianIndices(x)
  2<I[1]<5 && 3<I[2]<8 && continue
  x[I] = 0
end
x

I wondered if InvertedIndices.jl had an easy way, but x[Not(CartesianIndices((3:4, 4:7)))] .= 0; doesn’t appear to do this.

2 Likes

What is the reason you don’t want to write a for-loop? I’ve heard of a couple reasons why people may avoid explicit loops, but those don’t seem to match your case;

  1. functional styles avoid loops and use functions like map and filter. However, they also make copies instead of mutating objects, which you want here.

  2. Some languages have slow loops, so they promote a vectorized style to use faster loops coded in a lower-level language. However, Julia’s loops are already fast, so whether you write loops or vectorize is more a matter of convenience. In your case, a loop could be easier to read and write.

A simple approach:

A = fill(rand(),5,6)
idx, idy = 2:4, 3:4;
B = fill(1.0,size(A));  # Apply mask value everywhere
B[idx,idy] = A[idx, idy];   # Assign inner array values

Output B:

julia> B
5×6 Matrix{Float64}:
 1.0  1.0  1.0       1.0       1.0  1.0
 1.0  1.0  0.378626  0.378626  1.0  1.0
 1.0  1.0  0.378626  0.378626  1.0  1.0
 1.0  1.0  0.378626  0.378626  1.0  1.0
 1.0  1.0  1.0       1.0       1.0  1.0
1 Like

Here’s a nice way (I think) with Cartesian indices:

x = rand(6,10)
for I in CartesianIndices(x)
   if I ∉ CartesianIndex(2,5):CartesianIndex(3,8)
       x[I] = 0
   end
end

I think this has good performance, and it generalizes to n dimensions :slight_smile:

Shorter but probably less efficient:

I = setdiff(CartesianIndices(x), CartesianIndex(2,5):CartesianIndex(3,8))
x[I] .= 0
2 Likes

@mcabbott, those are all amazing solutions! How may we see how Julia interprets the first one, which uses ifelse and interpolation? It is mind-blowing the mix of simplicity & complexity at the same time.

So the $ is just protecting axes from broadcasting, to avoid axes.(x,1). The rest is not so bad if you take it apart, although still probably much harder to read a month later than just writing the loop!

julia> x = rand(6,10);

julia> ax1 = axes(x,1); ax2 = axes(x,2)'
1×10 adjoint(::Base.OneTo{Int64}) with eltype Int64:
 1  2  3  4  5  6  7  8  9  10

julia> 3 .< ax2 .< 8
1×10 BitMatrix:
 0  0  0  1  1  1  1  0  0  0

julia> @. (2 < ax1 < 5) & (3 < ax2 < 8)
6×10 BitMatrix:
 0  0  0  0  0  0  0  0  0  0
 0  0  0  0  0  0  0  0  0  0
 0  0  0  1  1  1  1  0  0  0
 0  0  0  1  1  1  1  0  0  0
 0  0  0  0  0  0  0  0  0  0
 0  0  0  0  0  0  0  0  0  0

julia> @. x = ifelse((2 < ax1 < 5) & (3 < ax2 < 8), x, 0)
6×10 Matrix{Float64}:
 0.0  0.0  0.0  0.0       0.0       0.0       0.0        0.0  0.0  0.0
 0.0  0.0  0.0  0.0       0.0       0.0       0.0        0.0  0.0  0.0
 0.0  0.0  0.0  0.310199  0.734067  0.504184  0.0244042  0.0  0.0  0.0
 0.0  0.0  0.0  0.962183  0.364773  0.62817   0.105928   0.0  0.0  0.0
 0.0  0.0  0.0  0.0       0.0       0.0       0.0        0.0  0.0  0.0
 0.0  0.0  0.0  0.0       0.0       0.0       0.0        0.0  0.0  0.0
2 Likes

@mcabbott, :clap:

Thanks a lot! Regarding your second solution, I tried using Pkg; Pkg.add("Tullio") but it failed, which is a shame because Tullio seems so much capable. Would you mind telling me what went wrong? Here’s the error:

Resolving package versions...
Unsatisfiable requirements detected for package GeometryTypes [4d00f742]:
 GeometryTypes [4d00f742] log:
 ├─possible versions are: [0.6.0-0.6.2, 0.7.0-0.7.10, 0.8.0-0.8.3] or uninstalled
 ├─restricted by compatibility requirements with Plots [91a5bcdd] to versions: 0.7.0-0.7.10
 │ └─Plots [91a5bcdd] log:
 │   ├─possible versions are: [0.12.1-0.12.4, 0.13.0-0.13.1, 0.14.0-0.14.2, 0.15.0-0.15.1, 0.16.0, 0.17.0-0.17.4, 0.18.0, 0.19.0-0.19.3, 0.20.0-0.20.6, 0.21.0, 0.22.0-0.22.5, 0.23.0-0.23.2, 0.24.0, 0.25.0-0.25.3, 0.26.0-0.26.3, 0.27.0-0.27.1, 0.28.0-0.28.4, 0.29.0-0.29.9, 1.0.0-1.0.14, 1.1.0-1.1.4, 1.2.0-1.2.6, 1.3.0-1.3.7, 1.4.0-1.4.4, 1.5.0-1.5.9, 1.6.0-1.6.12] or uninstalled
 │   └─restricted to versions 0.28.4 by an explicit requirement, leaving only versions 0.28.4
 ├─restricted by compatibility requirements with FixedPointNumbers [53c48c17] to versions: [0.7.4-0.7.10, 0.8.0-0.8.3] or uninstalled, leaving only versions: 0.7.4-0.7.10
 │ └─FixedPointNumbers [53c48c17] log:
 │   ├─possible versions are: [0.5.0-0.5.3, 0.6.0-0.6.1, 0.7.0-0.7.1, 0.8.0-0.8.4] or uninstalled
 │   └─restricted by compatibility requirements with Plots [91a5bcdd] to versions: 0.6.0-0.6.1
 │     └─Plots [91a5bcdd] log: see above
 └─restricted by compatibility requirements with ColorTypes [3da002f7] to versions: 0.8.2-0.8.3 or uninstalled — no versions left
   └─ColorTypes [3da002f7] log:
     ├─possible versions are: [0.7.0-0.7.5, 0.8.0-0.8.1, 0.9.0-0.9.1, 0.10.0-0.10.9] or uninstalled
     └─restricted by compatibility requirements with Colors [5ae59095] to versions: 0.10.0-0.10.9
       └─Colors [5ae59095] log:
         ├─possible versions are: [0.9.0, 0.9.2-0.9.6, 0.10.0-0.10.2, 0.11.0-0.11.2, 0.12.0-0.12.4] or uninstalled
         └─restricted to versions 0.12.3 by an explicit requirement, leaving only versions 0.12.3

Stacktrace:
 [1] #propagate_constraints!#61(::Bool, ::typeof(Pkg.GraphType.propagate_constraints!), ::Pkg.GraphType.Graph, ::Set{Int64}) at D:\buildbot\worker\package_win64\build\usr\share\julia\stdlib\v1.3\Pkg\src\GraphType.jl:1007
 [2] propagate_constraints! at D:\buildbot\worker\package_win64\build\usr\share\julia\stdlib\v1.3\Pkg\src\GraphType.jl:948 [inlined]
 [3] #simplify_graph!#121(::Bool, ::typeof(Pkg.GraphType.simplify_graph!), ::Pkg.GraphType.Graph, ::Set{Int64}) at D:\buildbot\worker\package_win64\build\usr\share\julia\stdlib\v1.3\Pkg\src\GraphType.jl:1462
 [4] simplify_graph! at D:\buildbot\worker\package_win64\build\usr\share\julia\stdlib\v1.3\Pkg\src\GraphType.jl:1462 [inlined] (repeats 2 times)
 [5] resolve_versions!(::Pkg.Types.Context, ::Array{Pkg.Types.PackageSpec,1}) at D:\buildbot\worker\package_win64\build\usr\share\julia\stdlib\v1.3\Pkg\src\Operations.jl:321
 [6] #add#112(::Bool, ::Pkg.BinaryPlatforms.Windows, ::typeof(Pkg.Operations.add), ::Pkg.Types.Context, ::Array{Pkg.Types.PackageSpec,1}, ::Array{Base.UUID,1}) at D:\buildbot\worker\package_win64\build\usr\share\julia\stdlib\v1.3\Pkg\src\Operations.jl:1010
 [7] #add at .\none:0 [inlined]
 [8] #add#25(::Bool, ::Pkg.BinaryPlatforms.Windows, ::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}, ::typeof(Pkg.API.add), ::Pkg.Types.Context, ::Array{Pkg.Types.PackageSpec,1}) at D:\buildbot\worker\package_win64\build\usr\share\julia\stdlib\v1.3\Pkg\src\API.jl:102
 [9] add(::Pkg.Types.Context, ::Array{Pkg.Types.PackageSpec,1}) at D:\buildbot\worker\package_win64\build\usr\share\julia\stdlib\v1.3\Pkg\src\API.jl:72
 [10] #add#24 at D:\buildbot\worker\package_win64\build\usr\share\julia\stdlib\v1.3\Pkg\src\API.jl:69 [inlined]
 [11] add at D:\buildbot\worker\package_win64\build\usr\share\julia\stdlib\v1.3\Pkg\src\API.jl:69 [inlined]
 [12] #add#21 at D:\buildbot\worker\package_win64\build\usr\share\julia\stdlib\v1.3\Pkg\src\API.jl:67 [inlined]
 [13] add at D:\buildbot\worker\package_win64\build\usr\share\julia\stdlib\v1.3\Pkg\src\API.jl:67 [inlined]
 [14] #add#20(::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}, ::typeof(Pkg.API.add), ::String) at D:\buildbot\worker\package_win64\build\usr\share\julia\stdlib\v1.3\Pkg\src\API.jl:66
 [15] add(::String) at D:\buildbot\worker\package_win64\build\usr\share\julia\stdlib\v1.3\Pkg\src\API.jl:66
 [16] top-level scope at In[109]:1