'test program for audio tone detection on ADC
 'uses PLAY SOUND ot generate test signal and noise
 'technology can be used for morse decoder
 'fft(32) takes 3ms, so 10ms time resolution can be possible
 
 Option base 1
 Option default float
 
 'for graphs and statistical data enable debug
 debug=1                   '0=no debug, 1=tone detector, 2=text decoder
 
 's/n threshold
 snt=0.25
 
 'process defines
 n=32                      'samples
 Dim a(n),b(n),c(n),w(n)   'data arrays for FFT and samples and window(=w())
 m=n/2                     'size for n/2 samples (half FFT)
 Dim d(m),scd(m),x(m),ad(m)'data arrays analysis,screen plot,average fft
 dim nb(m)                 'histogram base for frequency bins
 ready=0                   'new ADC data available
 pingpong=0                'fill samples in array a() or b()
 whte=RGB(white)           'frame color
 pkb%=6                    'carrier average default value
 cntl%=0                   'loop counter 10ms tick, 64bits
 
 
 'decoder filter paramters
 flt = 4                   'circular filter depth
 Dim ra(flt)               'circular filter for frequency bin
 
 'create a Hann window array w()
 For i=1 To n:w(i)=(Sin(i*Pi/n))^2:Next
 
 'graph x-axis and frame
 For i=1 To m:x(i)=4*i:Next      'x axis graph increments with 4
 CLS:print @(0,72)"Ti "
 if debug=1 then
   box 0,104,640,36:box 0,158,640,36:box 0,210,640,36
   print @(0,198)"Decode":print @(0,252)"bin"
   print @(0,92)"SD ":print @(0,146)"s/n"
 end if
 
 'hardware parameters
 freq = 3200                     'ADC sampling frequency
 testtone = 400                  'connect audio to adc gp26
 
 'test signal for ADC from PWM audio picomite
 'play tone testtone,testtone:play volume 100,100         'noise free max
 Play sound 1,b,s,testtone,15                            'signal
 Play sound 2,b,n,100,25                                 'noise added
 
 'open the ADC on gp26 and start first conversion
 ADC open freq,1,INT_RDY
 pingpong = 0
 ADC start a()
 
 'the main loop decodes the pre-processed fft sample arrays
 Do
   timer=0
   Do :Loop While ready=0              'wait until new samples available
   print @(24,72)str$(timer,2,3)       'time idle (waiting) in 10ms loop, debug
   
   'average spectrums 2x for better s/n (not more, delays will kill you)
   inc cntl%                           'time reference for decoding                          
   math C_ADD ad(),d(),ad()
   math scale ad(),0.5,ad()            'average 1/2
   
   'low pass filter on highets peak in the spectrum
   pk=Math(MAX ad(),f%)
   Inc i,1:If i>flt Then i=1           'pointer in history array
   ra(i)=f%                            'array with bins
   sdv=Math(sd ra())             'standard deviation bin array (0=homegeneous)
   
   'signal to noise for actual measurement
   sn=pk/((MATH(SUM ad()))-pk)       's/n value
   'sna=(sna*127+sn)/128             'rolling average 128 of sn
   
   'build spectrum histogram to determine morse carrier (bin# and weight)
   if sdv=0 then
     inc cnt,1:inc nb(f%),sn*cnt     'add weight to bin
   else
     cnt=0:pkv=math(max nb(),pkb%)   'find heaviest bin = carrier
     if debug=1 then print @(24,252)str$(pkb%,2,0)+" wght "+str$(pkv,6,0)
   end if
   
   
   if debug=1 then
     
     'diagnosis
     print @(24,92)str$(sdv,2,3)     'standard deviation x bins, 0=perfect
     print @(24,146)str$(sn,2,3)     'signal to noise 1=perfect, 0.12=bad
     
     'trend graphs
     Inc xi,1:If xi>638 Then
       xi=0
       box 0,104,640,36,1,whte,0
       box 0,158,640,36,1,whte,0
       box 0,210,640,36,1,whte,0
     end if
     ysdv=16*min(sdv,1.9)+106:pixel xi,ysdv      'standard deviation graph
     ysn=192-32*sn:pixel xi,ysn                  's/n graph
     
     'try decoding dit-dah when good s/n, correct + stable carrier
     pixel xi,242-30*(sn>snt)*(sdv=0)*(f%=pkb%)  'decoded carrier detect graph
   end if
   
   
   'here comes the actual decoding based on (sn>snt)*(sdv=0)*(f%=pkb%) and cntl%
   if debug=2 then
     
   end if
   
   'scale fft data to graph window
   'Math window d(),60,10,scd()                                 'window y
   Math scale ad(),-min(60/pk,99),scd():Math add scd(),65,scd() 'fixed y
   Box 0,0,70,70,1,whte,0:Line graph x(),scd()                  'show graph
   
   ready=0
   
   select case ASC(inkey$)
       case 32:  math set 1,nb()
       case 128: snt=min(1.1*snt,0.5)
       case 129: snt=max(0.9*snt,0.2)
       case 27: exit do
   end select
   
 Loop
End
 
 'interrupt routine switches between buffers and does pre-processing
Sub INT_RDY
 
 'toggle buffer
 pingpong = 1 - pingpong
 
 'offset removed from input data to avoid FFT element 0 to be dominant
 If pingpong Then
   ADC START b()
   Math add a(),-Math(mean a()),a()  'offset gone
   Math c_mult a(),w(),a()           'window data
   Math fft magnitude a(),c()        'perform fft
 Else
   ADC START a()
   Math add b(),-Math(mean b()),b()  'offset gone
   Math c_mult b(),w(),b()           'window data
   Math fft magnitude b(),c()        'perform fft
 EndIf
 
 'copy lower half of FFT in c() to d()
 Memory copy float Peek(varaddr c()),Peek(varaddr d()),m
 
 ready=1
 
End Sub