Signals can be more than 2d

[danmackinlay]

for SampledSignals, all audio is always 2d arrays (even 1d signals are effectively 2d in SampledBuffer)

However, signals might be higher dimension - for example a sliding window 256-point windowed DTFT of a stereo signal is most naturally 256x2xT for some length T. Is it wise to enforce the 2d thing?

Systems like pytorch, for example, choose a default dimension (IIRC the first dimension) for time series is chosen to be time, and allow the signal to be some arbitrary series of tensors.

[ssfrr]

There's a balance here between trying to make something that's super general vs. something simpler and more focused. It sounds like you might be looking for something more like AxisArrays.jl. We've actually talked in the past about replacing the built-in SampleBuf and SpectrumBuf types with AxisArrays - see issue #38 for more details.

I think that if we were to try to generalize to higher-dimensional objects, that be the thing that convinces me to just use AxisArrays, given that there's substantial overlap in the goals of the packages.

[ssfrr]

(also BTW I don't mean to shut down discussion by closing these issues, they're just not super actionable at the moment, we can definitely re-open in the future if there turns out to be something to be done).

[danmackinlay]

Thanks for the heads-up about AxisArrays - that does look handy indeed. Maybe even a little TOO featureful for me 😉 but I'l give it a try.

While I sort out my own needs, I wrote my own temporary version of SampledBuf/Spectrum that is much simpler, by dropping support for Julia<1.0 and by having only one type, with a time axis, and no opinions about the other axes (i.e. it doesn't treat 2d specially). I'm not necessarily advocating that as the way forward for SampledSignals, but it is startling that SampledSignals offers methods to take FFTs, types to store FTT spectra, and multichannel audio, but not a data structure than can store FFTs of multichannel audio. And it is suggestive that there are some simplifications maybe available.

Just for reference, I'll include my temporary code here

export AbstractSampledArray, SampledArray, domain, nchannels, nframes
"""
This is stripped down from SampledSignals which I can't use because
in SampledSignals everything is basically a matrix, but I want n-d arryas

Eventually I would like perhaps to extend this so that every axis could have a
sampling period in a different unit (not even neccessarily time). But I don't need
that right now (although it would be nice to have a choice of which axis was
time.)
"""
abstract type AbstractSampledArray{T, N} <: AbstractArray{T, N} end

export AbstractSampledArray

"""
Represents a regularly-sampled buffer that stores its own sample
rate (in samples/second). The wrapped data is an N-dimensional array.
The first axis is time, to which the sample-rate pertains.
The next dimensions are whatever you want.
So a 1-second stereo audio buffer sampled at
44100Hz with 32-bit floating-point samples in the time domain would have the
type SampledArray{Float32, 2}.
"""
struct SampledArray{T, N} <: AbstractSampledArray{T, N}
    data::Array{T, N}
    samplerate::Float64
end

# define constructor so conversion is applied to `sr`
SampledArray(arr::Array{T, N}, sr::Real) where {T, N} = SampledArray{T, N}(arr, sr)

SampledArray(T::Type, sr, dims...) = SampledArray(Array{T}(undef, dims...), sr)
SampledArray(T::Type, sr, len::Quantity) = SampledArray(T, sr, inframes(Int,len,sr))
SampledArray(T::Type, sr, len::Quantity, ch) =
    SampledArray(T, sr, inframes(Int,len,sr), ch)
function SampledArray(buf::AbstractSampleBuf)
    if nchannels(buf)==1
        SampledArray(vec(buf.data), samplerate(buf))
    else
        SampledArray(buf.data, samplerate(buf))
    end
end
SampleBuf(arr::AbstractSampledArray) = SampleBuf(arr.data, samplerate(arr))

# terminology:
# sample - a single value representing the value of 1 channel at some point in time (or frequency)
# frame - a collection of samples from each channel that were sampled simultaneously

# audio methods
samplerate(buf::AbstractSampledArray) = buf.samplerate
nframes(buf::AbstractSampledArray) = size(buf.data, 1)


# it's important to define Base.similar so that range-indexing returns the
# right type, instead of just a bare array
Base.similar(buf::SampledArray, ::Type{T}, dims::Dims) where {T} = SampledArray(Array{T}(undef, dims), samplerate(buf))
"""
Sample times for a given sample rate.
We don't use SampledSignals domain because it returns a
StepRangeLen, which we don't like because TwicePrecision poisons ForwardDiff.
"""
function domain(s::AbstractSampledArray)
    LinRange(
        0.0,
        (nframes(s)-1)/samplerate(s),
        nframes(s)
    )
end

import Base.broadcast

const ArrayIsh = Union{Array, SubArray, Base.LinRange, StepRangeLen}

# Broadcasting in Julia 0.7
# `find_buf` has borrowed from https://docs.julialang.org/en/latest/manual/interfaces/#Selecting-an-appropriate-output-array-1
find_buf(args::Tuple) = find_buf(find_buf(args[1]), Base.tail(args))
find_buf(x) = x

find_buf(bc::Base.Broadcast.Broadcasted{Broadcast.ArrayStyle{SampledArray{T,N}}}) where {T,N} = find_buf(bc.args)
find_buf(::SampledArray, args::Tuple{SampledArray}) = args[1]
find_buf(::Any, args::Tuple{SampledArray}) = args[1]
find_buf(a::SampledArray, rest) = a

Base.BroadcastStyle(::Type{SampledArray{T,N}}) where {T,N} = Broadcast.ArrayStyle{SampledArray{T,N}}()
function Base.similar(bc::Broadcast.Broadcasted{Broadcast.ArrayStyle{SampledArray{T,N}}}, ::Type{ElType}) where {T,N,ElType}
    A = find_buf(bc)
    SampledArray(Array{ElType}(undef, length.(axes(bc))), samplerate(A))
end

typename(::SampledArray{T, N}) where {T, N} = "SampledArray{$T, $N}"
unitname(::SampledArray) = "s"
srname(::SampledArray) = "Hz"

# the index types that Base knows how to handle. Separate out those that index
# multiple results
const BuiltinMultiIdx = Union{Colon,
                        Vector{Int},
                        Vector{Bool},
                        AbstractRange{Int}}
const BuiltinIdx = Union{Int, BuiltinMultiIdx}
# the index types that will need conversion to built-in index types. Each of
# these needs a `toindex` method defined for it
const ConvertIdx{T1 <: Quantity, T2 <: Int} = Union{T1,
                                                # Vector{T1}, # not supporting vectors of Quantities (yet?)
                                                # Range{T1}, # not supporting ranges (yet?)
                                                ClosedInterval{T2},
                                                ClosedInterval{T1}}

"""
    toindex(buf::SampledArray, I)

Convert the given index value to one that Base knows how to use natively for
indexing
"""
function toindex end

toindex(buf::SampledArray, t::Number) = t
toindex(buf::SampledArray, t::FrameQuant) = inframes(Int, t) + 1
toindex(buf::SampledArray, t::Unitful.Time) = inframes(Int, t, samplerate(buf)) + 1
toindex(buf::SampledArray, t::Quantity) = throw(Unitful.DimensionError(t, s))

# indexing by vectors of Quantities not yet supported
toindex(buf::AbstractSampledArray, I::ClosedInterval{Int}) =
    toindex(buf, minimum(I)*frames):toindex(buf, maximum(I)*frames)
toindex(buf::AbstractSampledArray, I::ClosedInterval{T}) where {T <: Quantity} =
    toindex(buf, minimum(I)):toindex(buf, maximum(I))

# AbstractArray interface methods
Base.size(buf::AbstractSampledArray) = size(buf.data)
Base.IndexStyle(::Type{T}) where {T <: AbstractSampledArray} = Base.IndexLinear()
# this is the fundamental indexing operation needed for the AbstractArray interface
Base.getindex(buf::AbstractSampledArray, i::Int) = buf.data[i];

Base.getindex(buf::AbstractSampledArray, I::ConvertIdx) = buf[toindex(buf, I)]
Base.getindex(buf::AbstractSampledArray, I1::ConvertIdx, I2::BuiltinIdx) =
    buf[toindex(buf, I1), I2]

function Base.setindex!(buf::AbstractSampledArray, val, i::Int)
    buf.data[i] = val
end

# equality
import Base.==
==(buf1::AbstractSampledArray, buf2::AbstractSampledArray) =
    samplerate(buf1) == samplerate(buf2) &&
    buf1.data == buf2.data

FFTW.fft(buf::SampledArray) = SampledArray(FFTW.fft(buf.data), nframes(buf)/samplerate(buf))
FFTW.ifft(buf::SampledArray) = SampledArray(FFTW.ifft(buf.data), nframes(buf)/samplerate(buf))