flat assembler
Message board for the users of flat assembler.

flat assembler > Examples and Tutorials > FIR low pass filter

Author
Thread Post new topic Reply to topic
gens



Joined: 18 Feb 2013
Posts: 159
FIR low pass filter
stereo, 16bit
44100 Hz sampling frequency
4000 Hz cutoff frequency
Blackman window
32 taps
fixed point
attuned to 0dB

forgot to add
amd64 and linux

Code:
format ELF64 executable
entry start

align 16
segment readable writeable

samples equ 1024

align 16
buff_in: rw samples*2
align 16
buff_out: rw samples*2

; 2^16 = 1.0
align 16
coefficients:

dw 5
dw 14
dw 35
dw 57
dw 39
dw -76
dw -328
dw -677
dw -945
dw -828
dw 16
dw 1799
dw 4428
dw 7441
dw 10109
dw 11680
dw 11680
dw 10109
dw 7441
dw 4428
dw 1799
dw 16
dw -828
dw -945
dw -677
dw -328
dw -76
dw 39
dw 57
dw 35
dw 14
dw 5

rq 10 ;just in case
tmp rq 100

sys_read equ 0
sys_write equ 1
sys_exit equ 60


align 16
segment readable executable

start:

xorps xmm7, xmm7
xorps xmm6, xmm6
xorps xmm5, xmm5
xorps xmm4, xmm4
xorps xmm3, xmm3
xorps xmm2, xmm2
xorps xmm1, xmm1
xorps xmm0, xmm0

mov [tmp], rdx
mov [tmp+8], rdx

filter:

mov rdx, samples*2*2
mov rsi, buff_in
mov rdi, 0
mov rax, sys_read
syscall
cmp rax, 0
jz end_fir

xor rax, rax

mov rcx, samples
loopy:
movaps xmm12, [coefficients]
movaps xmm13, [coefficients+16*1]
movaps xmm14, [coefficients+16*2]
movaps xmm15, [coefficients+16*3]
mov edx, [rax+buff_in] ;dx = left channel sample
mov esi, edx
shr esi, 16 ;si = right channel sample

; left channel delay line
pslldq xmm3, 2
movaps xmm10, xmm2
psrldq xmm10, 2*7
orpd xmm3, xmm10

pslldq xmm2, 2
movaps xmm10, xmm1
psrldq xmm10, 2*7
orpd xmm2, xmm10

pslldq xmm1, 2
movaps xmm10, xmm0
psrldq xmm10, 2*7
orpd xmm1, xmm10

pslldq xmm0, 2 ;bytes
movzx rdx, dx
movq xmm10, rdx
orpd xmm0, xmm10

; right channel delay line
pslldq xmm7, 2
movaps xmm10, xmm6
psrldq xmm10, 2*7
orpd xmm7, xmm10

pslldq xmm6, 2
movaps xmm10, xmm5
psrldq xmm10, 2*7
orpd xmm6, xmm10

pslldq xmm5, 2
movaps xmm10, xmm4
psrldq xmm10, 2*7
orpd xmm5, xmm10

pslldq xmm4, 2
movq xmm10, rsi
orpd xmm4, xmm10




movaps xmm8, xmm12
movaps xmm9, xmm13
movaps xmm10, xmm14
movaps xmm11, xmm15


pmaddwd xmm12, xmm0
pmaddwd xmm13, xmm1
pmaddwd xmm14, xmm2
pmaddwd xmm15, xmm3

paddd xmm12, xmm14
paddd xmm13, xmm15

paddd xmm12, xmm13

movhlps xmm13, xmm12
paddd xmm12, xmm13

movss xmm13, xmm12
psrldq xmm12, 32
paddd xmm12, xmm13



pmaddwd xmm8, xmm4
pmaddwd xmm9, xmm5
pmaddwd xmm10, xmm6
pmaddwd xmm11, xmm7

paddd xmm8, xmm10
paddd xmm9, xmm11

paddd xmm8, xmm9

movhlps xmm9, xmm8
paddd xmm8, xmm9

movss xmm9, xmm8
psrldq xmm8, 32
paddd xmm8, xmm9



movd r10d, xmm12
;sal r10d, 3 ;overflow possibility unless compensated
shr r10d, 16
mov [rax+buff_out], r10w

movq r10, xmm8
;sal r10d, 2
shr r10d, 16
mov [rax+buff_out+2], r10w

add rax, 4
dec rcx
jnz loopy

mov rdx, samples*2*2
mov rsi, buff_out
mov rdi, 1
mov rax, sys_write
syscall

jmp filter

end_fir:
mov rax, sys_exit
syscall    


Image
frequency response, needs some work Smile (but not much)

Code:
define int(x) {
    auto old_scale   /* variables global by default */
    old_scale=scale  /* scale is global */
    scale=0; ret=x/1
    scale=old_scale
    return ret
}
define round(x) {
    if (x<0) x-=.5 else x+=.5
    return int(x)
}

scale=50

#pi=3.14159265358979323846264338327950288
pi=4*a(1)
print "pi is "
pi
print "\n"

fs = 44100
fc = 4000
omega = (2*pi*fc)/fs
omega/pi

nf=31
nt=nf+1
m=nf/2

print "ideal filter coefficients\n"
for (n=0; n<nt; n++) {
    if (n == m) {
    hd[n]=(omega/pi)
    hd[n]
    continue }
    
    hd[n]=( s(omega*(n-m)) ) / (pi*(n-m))
    hd[n]
}


#print "Blackman-Harris window coefficients\n"
#for (n=0; n<nt; n++) {
#    w[n] = 0.35875 - 0.48829*c( (2*pi*n)/nf ) + 0.14128*c( (4*pi*n)/nf ) - 0.01168*c( (6*pi*n)/nf )
#    w[n]
#}

print "Blackman window coefficients\n"
for (n=0; n<nt; n++) {
    w[n] = 7938/18608 -  (9240/18608)*c( (2*n*pi)/nf ) + (1430/18608)*c( (4*n*pi)/nf )
    w[n]
}

#print "Blackman window coefficients\n"
#for (n=0; n<nt; n++) {
#    w[n] = 0.42 -  0.5*c( (2*n*pi)/nf ) + 0.08*c( (4*n*pi)/nf )
#    w[n]
#}

for (n=0; n<nt; n++) {
    h[n] = w[n]*hd[n]
    h[n]
}

print "coefficients sum\n"
for (n=0; n<nt; n++) {
    coeff_sum+=h[n]
}
coeff_sum
print "one over coefficients sum\n"
1/coeff_sum

print "final filter coefficients\n"
for (n=0; n<nt; n++) {
    h[n] = h[n] / coeff_sum
    #h[n]
}

print "normalised to 2^16\n"
for (n=0; n<nt; n++) {
    hn[n] = round( h[n]*(2^16) )
    print "dw "
    hn[n]
    total = total + hn[n]
}

print "new sum\n"
total
print "aim\n"
2^16    

filter calculator in bc language
needs bc ofc, run it with bc -l filter_calc.bc


Last edited by gens on 28 Mar 2014, 02:58; edited 1 time in total
Post 26 Mar 2014, 13:31
View user's profile Send private message Reply with quote
neville



Joined: 13 Jul 2008
Posts: 507
Location: New Zealand
Very nice gens Cool

About 30 years ago I demonstrated an FIR low-pass filter (also an IIR version) using a pair of 6800 CPU's each clocked at about 600kHz! Yes, kilohertz Smile

One of the CPU's was dedicated to producing the variable-frequency digital sine sequences because I didn't have access to suitably fast enough A/D converters at the time.

Needless to say, the specs of my filter weren't quite as impressive as yours. The useable bandwidth was about 100 Hz. I remember using cut-off frequencies starting at 20 Hz! And of course a decent analogue anti-alising filter!

_________________
FAMOS - the first memory operating system
Post 27 Mar 2014, 03:38
View user's profile Send private message Visit poster's website Reply with quote
gens



Joined: 18 Feb 2013
Posts: 159
ty

seems you knew more back then then i do now Smile


this is a step towards a resampling filter for a full sound server
the code also has some deficiencies (3 numbers are added then discarded, etc) but that will be easy to fix when its done

oh ye, and it is very fast since it is simple
in fact the bottleneck is the read/write calls in the kernel
shows how fast cpu's are
Post 27 Mar 2014, 04:27
View user's profile Send private message Reply with quote
gens



Joined: 18 Feb 2013
Posts: 159
was an error in the code
or'ed with the wrong register in the last quarter of the delay lines

fixed now, spectrogram remade to match
Post 28 Mar 2014, 02:59
View user's profile Send private message Reply with quote
Display posts from previous:
Post new topic Reply to topic

Jump to:  


< Last Thread | Next Thread >
Forum Rules:
You cannot post new topics in this forum
You cannot reply to topics in this forum
You cannot edit your posts in this forum
You cannot delete your posts in this forum
You cannot vote in polls in this forum
You cannot attach files in this forum
You can download files in this forum


Copyright © 1999-2019, Tomasz Grysztar.

Powered by rwasa.