Simple question about Array and Subarrays

General Usage
1

Let’s say that for a given dimension d, we need to constrain the matrix to an specific value (or index). In other words, we need the subarray in which that dimension has an specific value (or index). How I can do something like:

B = A[where dim=d has index = i]

Maybe there is an specific method in Julia that I don’t know.

The motivation is that we read a variable A from a netcdf file, and the “position” of an specific coordinate in a A, in such way that depending on input netcdf that dimension (coordinate) can be in arbitrary axis in the array of the variable A. We need to extract a subarray based on the value of an specific coordinate whose position in the numerical array is variable depending on the input netcdf file.

Thank you.

2

Could you provide a MWE?
Also, just in case, there is a NetCDF.jl package.

2 Likes
3

Hello and thank you!

I am using NCDatasets.jl, and I could use GeoData.jl and actually I can do what I want. But I was wondering how to do without GeoData.jl

using NCDatasets
using GeoData

inforule = Dict("zlev" => 3)
diag      = "o3"
finput    = "o3_dataset.nc"

function geomean_atlev(diag, finput, inforule)

    ds_diag = Dataset(finput)[diag]
    data = geoarray(finput)
    subdata = Array( data[ z=inforule["zlev"] ] )
    data = convert(Float64, mean(Array(subdata)))

    return diag, finput, data, inforule

end

For example how to get the subdarray named “subdata” above in the code like this below, where
I know that the axis where I want to constrain the array is “posz”

using NCDatasets

inforule = Dict("zlev" => 3)
diag      = "o3"
finput    = "o3_dataset.nc"

function geomean_atlev(diag, finput, inforule)

    ds_diag = Dataset(finput)[diag]
    dimsizes= dimsize(ds_diag)
    posz=0
    let counter=0
    for (index, value) in collect(pairs(dimsizes))
        counter=counter+1
        if index==:z
           posz=counter
        end
    end
    [...]

    end # let

end #function
4

I don’t fully understand your problem, but this toy example using selectdim, views, and BitMatrix might help. The selectdim function creates a view. Modifying a view modifies the original matrix. The BitMatrix type can be used for logical indexing.

julia> dim = 2; index = 3; oldvalue = 7; newvalue = 10;

julia> A = [5 6 7; 7 4 7; 9 7 6]
3×3 Matrix{Int64}:
 5  6  7
 7  4  7
 9  7  6

julia> Aview = selectdim(A, dim, index)
3-element view(::Matrix{Int64}, :, 3) with eltype Int64:
 7
 7
 6

julia> B = falses(size(A))
3×3 BitMatrix:
 0  0  0
 0  0  0
 0  0  0

julia> Bview = selectdim(B, dim, index)
3-element view(::BitMatrix, :, 3) with eltype Bool:
 0
 0
 0

julia> Bview .= Aview.==oldvalue
3-element view(::BitMatrix, :, 3) with eltype Bool:
 1
 1
 0

julia> B
3×3 BitMatrix:
 0  0  1
 0  0  1
 0  0  0

julia> A[B] .= newvalue
2-element view(::Vector{Int64}, [7, 8]) with eltype Int64:
 10
 10

julia> A
3×3 Matrix{Int64}:
 5  6  10
 7  4  10
 9  7   6
1 Like
5

Thank you @Nathan_Boyer
Yes, I think that selectdim is what I was looking for. If A is a 3D array, then

Aview = selectdim(A, dim, index)

would be the 2D matrix that I was looking for. :grinning:

1 Like
6

If you do have trouble converting to 3D, don’t be afraid to just for loop through the third dimension. Loops are fast in Julia.