对Julia range进行操作,从物理单元到数组索引



我正试图想出一种方法来操作Julia范围,以便它们在某些单位中指定,然后转换为索引。为了更好地解释,请考虑下面的例子:

julia> Δx = 150
150
julia> x = collect(0:Δx:1000)
7-element Array{Int64,1}:
0
150
300
450
600
750
900
julia> y = @. x^2 + 50
7-element Array{Int64,1}:
50
22550
90050
202550
360050
562550
810050

假设我想获得x在100到500之间的每个值。如果我知道相应的索引,我可以通过使用range轻松地根据索引获得数组的一个块:

julia> index_chunk = 2:4
2:4
julia> y[index_chunk]
3-element Array{Int64,1}:
22550
90050
202550

然而,我希望能够得到相同的结果,但创建一个基于单位的索引。就像

julia> phys_chunk = 100:500
100:500
julia> y[phys_chunk/Δx]
ERROR: BoundsError: attempt to access 7-element Array{Int64,1} at index [0.6666666666666666:0.006666666666666667:3.3333333333333335]
Stacktrace:
[1] throw_boundserror(::Array{Int64,1}, ::Tuple{StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}}}) at ./abstractarray.jl:541
[2] checkbounds at ./abstractarray.jl:506 [inlined]
[3] _getindex at ./multidimensional.jl:742 [inlined]
[4] getindex(::Array{Int64,1}, ::StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}}) at ./abstractarray.jl:1060
[5] top-level scope at REPL[24]:1

当然,这行不通。我已经尝试了一些其他的东西沿着相同的路线(如phys_chunk//Δx),但我没有任何运气。有可能做我想做的事吗?

PS:对于那些熟悉Python的人,我试图重现xarray.sel()切片方法的行为:

In [4]: da = xr.DataArray([1, 2, 3, 4, 5], [("x", [0, 1, 2, 3, 4])])
In [5]: da
Out[5]: 
<xarray.DataArray (x: 5)>
array([1, 2, 3, 4, 5])
Coordinates:
* x        (x) int64 0 1 2 3 4
In [6]: da.sel(x=slice(0.9, 3.1))
Out[6]: 
<xarray.DataArray (x: 3)>
array([2, 3, 4])
Coordinates:
* x        (x) int64 1 2 3

逻辑索引呢?

julia> y[100 .< x .< 500]
3-element Vector{Int64}:
22550
90050
202550

编辑:使用DimensionalData添加一点物质。考虑到评论。(注意,我刚刚试用了这个包,所以可能有更好的方法来使用它。)

让我创建一个名为V的2DDimArray,它具有与原始帖子示例相同的x维度,但现在也有一个名为y的第二大维度:

julia> using DimensionalData
julia> Δx = 150
150
julia> x = 0:Δx:1000 # your x (no need to collect it)
0:150:900
julia> y = -500:100:500 # let's make the example 2D (now y is the 2nd dim)
-500:100:500
julia> V = DimArray(x.^2 .+ 2 .* y', (X(x), Y(y))) # V is a 2D DimArray
DimArray (named ) with dimensions:
X: 0:150:900 (Sampled: Ordered Regular Points)
Y: -500:100:500 (Sampled: Ordered Regular Points)
and data: 7×11 Matrix{Int64}
-1000    -800    -600    -400    -200  …     400     600     800    1000
21500   21700   21900   22100   22300      22900   23100   23300   23500
89000   89200   89400   89600   89800      90400   90600   90800   91000
201500  201700  201900  202100  202300     202900  203100  203300  203500
359000  359200  359400  359600  359800     360400  360600  360800  361000
561500  561700  561900  562100  562300  …  562900  563100  563300  563500
809000  809200  809400  809600  809800     810400  810600  810800  811000

请注意,V的值与y == 0的原始数据值相同:

julia> V[Y(At(0))]
DimArray (named ) with dimensions:
X: 0:150:900 (Sampled: Ordered Regular Points)
and referenced dimensions:
Y: 0 (Sampled: Ordered Regular Points)
and data: 7-element Vector{Int64}
[0, 22500, 90000, 202500, 360000, 562500, 810000]

我可以这样选择y=0100<x<500的值:

julia> V[X(Between(100, 500)), Y(At(0))]
DimArray (named ) with dimensions:
X: 150:150:450 (Sampled: Ordered Regular Points)
and referenced dimensions:
Y: 0 (Sampled: Ordered Regular Points)
and data: 3-element Vector{Int64}
[22500, 90000, 202500]

最新更新