Hi,

I am translating some code from Matlab to Julia and can’t find an equvalent to upfirdn. The only working option I have is to use PyCall and use the Python (SciPy) version. Like this:

```
scipy_signal = pyimport("scipy.signal")
s = scipy_signal.upfirdn(filter, s, p, q)
```

I also tried to make my own version using DSP.resample and DSP.Filters.filt. But it doesn’t work…yet .

Have I missed some great library?

1 Like

It seems like you should be able to implement this easily with resample and filter. What is the issue you are having?

I had problems with padding. It works now, but an explanation why it works would be appreciated

```
function upfirdn(s, filter, p, q)
pad = length(s)*p+length(filter)-p
if p != 1
s_padded = [s; zeros(ComplexF64, length(filter))] # need padding for the filter
s = DSP.resample(s_padded, p, [1]) # [1] makes it a zero insertion upsampling
end
s = DSP.Filters.filt(FIRFilter(filter), s)
s = s[1:pad] # Take away padding
if q != 1
#result_length = Int(p * length(s))
s = DSP.resample(s, 1/q, [1]) # downsample
#s = [s; zeros(ComplexF64, result_length - length(s))] # zero padding to get correct length
end
#=
@show size(s), p, q, size(filter)
scipy_signal = pyimport("scipy.signal")
s = scipy_signal.upfirdn(filter, s, p, q)
@show size(s)
=#
return s
end
```

It works but it is slightly slower than using PyCall. So I need to find a way to speed it up.

There are a number of things that you could try, but the first thing would be to

```
using BenchmarkTools
@benchmark upfirdn(...)
```

to figure out how much time it’s taking and how much GC it’s doing.

Your best bets are going to be to reduce allocations. perhaps create a buffer that you reuse by calling filt! and then use a @view to strip padding etc.

2 Likes

Hi @psteffensen,I just needed to point something out here.

What you have implemented is exactly correct for getting the result, but it is a very slow and inefficient algorithm. You should use the FIRFilter stateful filter object with a rational rate, which uses an efficient polyphase resampler.

```
FIRFilter(h[, ratio])
```

Construct a stateful FIRFilter object from the vector of filter taps `h`

. `ratio`

is an optional rational integer which specifies the input to output sample rate relationship (e.g. `147//160`

for converting recorded audio from 48 KHz to 44.1 KHz).

Best regards,

Tom

2 Likes

Thank you Tom (sirtom67) for pointing that out,

I actually ended up using “DSP.Filters.filt!” in my code. I guess it is the same as FIRfilter. I used it because it was much faster, like you say, but I thought I was cheating

My code only needed “upfir” and not the “down” part, so it ended up looking like this:

```
function upfirdn(s, filter, p, q)
## FIR (Filter)
result = zeros(ComplexF64, length(s)*p+length(filter))
DSP.Filters.filt!(result, FIRFilter(filter, p), s)
return result
end
```

Then for the “down” part i guess you can you resample to decimate, like this:

```
if q != 1
s = DSP.resample(s, 1/q, [1])
end
```

Looks good.

The “down” part here suggests part of the problem - for downsampling with factor q > 1, it is inefficient to call filt first, and then throw away (q-1) out of q samples - you’ve done “q”-times as much work as required! So - if it matters (and it may very well not matter if the code runs fast enough due to small signal and / or filter lengths) - it is better to include the “q” downsample factor in the FIRFilter.

1 Like