Conference napalm::commusic_v1

    Sorry, there really isn't any way to do that correctly (the
    undersampling, that is),  unless you go to something like a linear
    predictive code or some other data compression algorithm.
                            
    How much time and effort (and code space) do you want to spend,
    and how much smarts will the reproduction hardware have (as in,
    will the repro hardware be capable of unencoding the LPC, etc.)?
    
    A good (well, it was written by a prof. that I have a lot of respect
    for) is "Communication Systems" by A. Bruce Carlson.  It has several
    good chunks on quantization, data compression, and information
    theoretic implications of the above.  About 25% of the book will
    be directly relevant to what you want.  Another 20% will be fun
    to read just because it's darn good!
    
    	-Bill
    
    (p.s. if you just want some simple algorithms (not code, just
    algorithms) let me know and I will write them up)

    Bill, thanks for your quick reply!  I'll look for the book.  The
    problem I'm having is too much data, I can get about 6.5 seconds
    of high quality digital sound, which takes 1/2 megabyte of storage
    space.  What I have in mind is:
    
    1)	Do some kind of DSP low-pass filter and then throw away every
    other sample to halve the data rate.  Yes, I realise this cuts my
    bandwidth in half, but that's life....
    
    2)	Do some kind of DSP dynamic range compression to cut the data
    down from 12 bits to 8 bits.  The easy way out is to chop off the
    low order nibble.  This again (as long as we are talking bytewide
    EPROMs) cuts the storage requirement in half.
    
    3)	Edit the start and end points to cut the time in half, to pick
    the 3.25 seconds of interest.
    
    Doing all three of these things makes the 'sound' fit in a 27512
    (64K x 8) EPROM.  I want to be able to plug EPROMs into a box in
    my car and 'play' them out (with the appropriate 8-bit DAC, power
    amp, and PA horn) on demand.  Does all this sound do-able from a
    DSP point of view?
    
    Willie
                                        

    I think I have a cheap filter around here, in pascal or basic.
    It averages points the same way smoothing does on the tektronix
    7854.
    I get the feeling you could use a sound effects chip.
    77467 or something, what was that thing called?
    Tom

    I also think that i f you wanted to do speech, you could 
    buy speech ROMs if you buy into the speech decoding chip family
    they're for.
    the commodore amiga is a good potential development tool for this.
    but I don't have a a/d on it, although they're available.
    Tom

    If your sound loops, you could do the standard looping approach
    a la regular old sampler.
    
    Steve

    Don't just "throw away" half of your data points!  What will happen
    is that you'll lose half of your bandwidth and any frequencies in
    the upper half (which you lost) will "alias" down into the lower
    half.  
    
    It will sound _horrible_!
    
    Cutting it down to 4 bits (= 24 dB S/N) and keeping the sample rate up
    would probably sound much better!
    
    How much CPU are you willing to put into this?  If it's a Z-80 class
    machine or better, you can probably use an LPC compression scheme
    and get something that sounds very good. (the Z-80 can use the EEROM
    space that is saved by the compression :-) )
    
    	- Bill
    

    Hey, that might be a good idea for Tycho!  I'll have to build him
    one too...  Nope, this is for my new Golf, kind of a hardware hacker's
    version of those 'terminator' things you can attatch to your dash.
    
    I know I don't want to throw away every other sample, I need to
    low-pass filter it first and then throw away every other sample.
    (or whatever the equivalent is).
    
    I _could_ keep the data rate the same, but it would double my EPROM
    purchasing, and 27512s are around $27 each....
    
    I'm not planning on any CPU power at all, just a counter chain running
    thru all the addresses and clocking them into an 8-bit DAC.  I really
    don't want to get into any kind of compression or anything that
    would distort the sound quality, as I want to be able to recognize
    things like photon torpedos and phasers and blasters and Cylons
    and voices (ie: Darth Vader saying "It is your DESTINY" to clear
    pedestrians out of the way :+)  ).
    
    I don't know that keeping the original sampling rate will get me
    anywhere, as a 20KHz signal isn't going to make it thru a cheap
    PA horn anyway....
    
    I'd love some algorythms!
    
    Thanks,
    
    Willie

>        (p.s. if you just want some simple algorithms (not code, just
>    algorithms) let me know and I will write them up)
 
    Bill,
    
    	is your offer still open?  I'd really appreciate any algorythms
    or other info!
    
    Thanks,
    
    Willie

    Wuppps!  Sorry, I forgot.
    
    OK, here's a simple algorithm to reduce your sample rate, and do the
    low-pass filtering at the same time, in basic. Assume file #1 contains
    the high-rate samples, and file #2 will contain the results. 
                   
    All it does is average five sample points (S1, S2, S3, S4, S5) in a
    weighted average.  The coefficients (c1, c2, c3, c4, c5) are chosen by
    hocus-pocus so that the energy of a wave isn't greatly affected by the
    phase of the sampling program, and the resulting waveform has
    properties we think we'll like.  First job is to figure the 
    coefficients.  We do this by choosing criteria and then writing
    algebraic equations that make the same statement.
                                                           
    Here's an example derivation:
    
    We want our filter to have a 1:1 output for D.C. signals:
    
    c1 + c2 + c3 + c4 + c5 = 1	  	; total response sums to 1 for
    					; a DC signal
                                                             
    
    We also decide we want the filter to look the same in either direction-
    i.e. it's symmetrical:
    
    c1=c5				;filter is symmetrical
    c2=c4                                ;
    
    so we now know
    
    2*c1 + 2*c2 + c3 = 1 			; drop out c4 and c5
    
    We also know that the result for a completely aliasing wave
    (a +1/-1 triangle wave) should be filtered out completely (zero
    response)
    
    2*c1 - 2*c2 + c3  = 0		;a triangle wave has no response
    
    We need one more constraint to decide what the coefficients should
    be.  That can be the "cutoff" of the filter cutoff curve.  We have
    to be careful, however.  If we cutoff too low or too high, then
    we'll leak aliasing through the filter, or cut off signal that 
    we wanted to let through...
    
    Let's make the 3dB point (dreads, now there's three of them :-) )
    be at half the old sample frequency.  That means a square wave
    at input value=�1 should have an output value of .5 (and it will
    end up looking like a triangle wave) 
    
    c1 + c2 - c3 - c2 + c1 = .5
    
    or:
    
    2*c1 - c3 = .5
    
    Three equations, three unknowns.  Solving, we get:
    
    c1= .25, c2 = .25, c3 = 0.0 !
               
    which is rather surprising, but if you work out a few waves, actually
    does seem to generate reasonable looking output waves!
    
    So our program is:
                      
    10 c1 = .25
    20 c2 = .25
    30 c3 = 0.0
    100 read#1 s1 
    110 read#1 s2 
    120 read#1 s3 
    130 read#1 s4
    140 read#1 s5
    150 o1 = c1*s1 + c2*s2 + c3*s3 + c2*s4 + c1*s5
    160 write#2 o1 
    170 s1 = s3
    180 s2 = s4
    190 s3 = s5 
    200 go to 130 
      
    That will average 2 points into one and provide some "tailing" to
    suppress some of the higher aliasing frequencies.  It's not perfect
    (You expect me to remember the Studer intersampling algorithm on
    a Friday night???:-) )  It will also have a slight phasing effect
    but it's small.
    
    -----
    
    You can repeat the above analysis and algebra to "build" any kind
    of filter you want.  You might say you want the response to a triangle
    to be 2x the normal response, zero response to DC, and a 50% response
    to a 2x frequency square wave.  Just do out the algebra, plug in
    the coefficients, and take a listen!
    
    A hint: you will get better fidelity if you have more points on
    the sample curve in your input stream- going to 7 or 9 or even more
    may help.
    
    You also are not restricted to "subsampling"; if you are willing
    to keep the same sample rate you can use the above analysis and
    generate any digital filter you want.  It'll have some interesting
    grunge (because of the way we're choosing coefficients) but it will
    be _interesting_ grunge.
    
    There is also a whole body of theory on how to make it "more accurate",
    things like Laplace and Z-transforms, etc. but that's a senior/grad
    level course in electrical engineering if you're going to learn it
    "right".  (i.e. I took the courses and don't remember them, don't
    make me embarass myself in public!)
    
    Try D.K. Frederick's "Discrete-Time Systems" or A.B. Carlson's
    "Communication Systems" for the hard math approach.  Especially
    the latter, I know A.B.C. and he's an excellent writer. (plug!)
    	
    	-Bill

    
    Many thanks, Bill!
     
    I'll look it over and see what I can come up with.  Prepare yourself
    for the standard silly questions.  See what you get for helping!
    :+}
    
    Willie

    You also should probably either look up a text on linear algebra
    and Cramer's rule (it makes solving 3x3's and 4x4 systems easy,
    almost do-it-in-your-head-easy).
    
    A computer is also a "big win".
    
    -----
    
    If you want to get into "computer EQ's" then you should learn about
    convolution.  The 5� convolution tour follows:
    
    Convolution is the process of taking two waves (call them A and
    B) and munging them together to do something useful.  Convolution
    acts like a digital parametric EQ with some truly fantastic properties.
    
    More specifically, remember the waves are really functions of time.
    So, A is really A(t) and B is really B(t).  Think of the "(t)" as
    being a subscript into the memory array holding your samples.
       
    Of what use is this?  Well, a mathematician named Fourier came up with
    the interesting notion that not only could _any_ infinitely repeating
    wave be built of nothing but sine waves, but _also_ that the amounts of
    each frequency of sine wave could be determined just by multiplying the
    desired wave samples by a sine wave sample and summing the products
    (being careful to keep the � signs straight, of course). 
    
    So, if we have wave A(t), and we want to find out how much Sin(300Hz)
    is in it, we sample A "very finely", generating an array in memory
    called A(t), and we calculate a sine wave at 300 Hz, and term-to-term
    multiply them, and then sum the results of the multiplications.
    The result is how much 300 Hz energy is in the wave. 
    
    Let's assume we have 1024 points in this A(t) array, and 1024 points
    in the sin(300*2pi*t) array.  We do 
    	A(1)*sin(300*2pi*1) + 
    	A(2)*sin(300*2pi*2) +
    	A(3)*sin(300*2pi*3) +
    		.
    		.
    		.
    	A(1024) * sin (300*2pi*1024) = F ( A(t) ) @ 300Hz
    
    We go through the frequencies, multiplying and summing, and we'll
    get another sample curve, but this sample curve isn't in terms of
    time (t), it's in terms of frequency (f).  We can then "sculpt"
    the F(f) curve by hand or algorithm to get whatever filter design
    we like (like absolutely zeroing out everything above and below
    some desired bandpass, etc).
    
    Once we've built our new desired F(f) curve, we can turn the F(f)
    back into A(t) by just summing together a bunch of sine waves, each
    sine wave having the amplitude in that frequency range that the
    F(f) says.  For example, if we wanted a dual bandpass filter with
    no attenuation between 300 and 310 hz, and no attenuation at 700
    to 705 hz, but nothing anywhere else, we would sum up 15 sine waves
    (300, 301, 302, 303...309, 310, 700, 701, 702, 703, 704, 705). 
    with amplitudes that we got from the F(f) array.  The result will
    be the filtered digitized sample B(t).
    
    This really works, but as you can see it's slow.  There are a lot
    of multiplies to be done just to do a meager 1024 points (over a
    million of them).  But don't worry, there are better and faster
    ways.
    
    One way is by FFT (Fast Fourier Transform) which I can't explain
    because I can't remember it.  Look it up if you're interested.
    Another way is by doing "convolution", which is what we started
    out to do in the first place. 
    
    In convolution, you build the desired filter function once, inverse-
    transform it (method of adding sine waves, as above) then use the
    result directly on the A(t) stream.  Yes, this method generates an
    array that, when termwise multiplied and summed over a sample,
    generates the filtered result of that sample.  You then slide the
    convolution array down one sample (or the A(t) array up one sample,
    same difference), and multiply/sum again to generate the next term in
    the output array B(t).  Keep going like this until you have all of the
    output samples you want. (actually, you should flip the convolution
    array end-for-end to be mathematically correct but in your application
    it really makes very little difference!) 
    
    This actually works very well if you're going to burn silicon or build
    circuits. You pre-calculate the factors, and use them to set the
    value of resistors.  Each resistor reads from a tap on a time delay
    line.  The resistors all lead to a summing junction where the
    output value is tapped off.  This is how they can build realtime
    FIR filters with such fantastic specs as they do, as cheaply
    as they do.
                                         
    Lastly, you can do a digital vocoder by taking the A(t) signal,
    transforming it to F(f) by sine wave multiplication or FFT,
    then using the F(f) directly as the convolution filter array.
                                            
    -----
    
    Again, realize that this is just the 5� tour and that you really
    should hit the books to understand this as more than black magic.
    
    	-Bill
    
    

    I wrote a cheap fourier analysis program in basic if anyone wants
    it.
    Tom

    I also wrote a linear algebra package in Ada, a complex arithmetic
    package ada (the fourier can be done complex), and a pascal
    or basic program that synthesizes waves (one cycle)with additive
    harmonics synthesis and graphic output.
    This was for the pdp11 synthesis weekend project.
    Tom

    Bill,
    
    	Thanks for the info.  I did some convolution at N.U. but don't
    remember much about how to make up the filter.  Mostly they just
    gave us a couple of things to convolve and asked for the result.
    How do you pick the number of terms and values of the terms to make
    (say) a low-pass filter using convolution, and how does that compare
    to your other filter (in .9?)?
    
    Tom,
    
    	If you have the stuff up on the network, can you point me at
    it?  I'll probably build some kind of a kludge in Z-80 assembly,
    (I hate waiting!) but it might be useful to have something that's
    known to work.
    
    I thought of trying to do an FFT, carve away all the stuff above
    my 'cuttoff' frequency, and then invert things again, but we are
    talking about 1/2 megabyte of data here, and all those complex math
    operations take a while....
    
    Willie

Here are some application notes from chip vendors about dsp.  I would
    go for the Intel 2920 applications book,but  I have  it already.
    TRW LSI Products PO Box 2472 La Jolla CA 92038 (714)457-1000
TP7A Hardware Development for a General Purpose Digital Filter Computing Machine
TP4A Digital Signal Processing For Radar Systems 
TP10A An Introduction to Digital Spectrum Analysis Including A High Speed 
	FFT Processor Design 
INTEL
AP-92 Implementation of a Scanning Spectrum Analyzer Using the 2920 Signal 
	Processor 
Intel 2920 Assembly Lanugage manual
2920 Simulator User's Guide
2920 Signal Processing Applications s/w Compiler User's Guide
Using the 2920 Signal Processor in Modem Applications

Texas Instruments 
SPRY015 TI Leading Electronics Press Coverage TMS320 
SPRA002 Floating Point Arithmetic with the TMS32010 

    Tom

    I think I forgot part of the fourier transform bit where you also
    do a cosine array multiplication.  Yeah, that's it; the sine wave
    job only gets you the odd part of the wave (where x(t) = -x(-t)).
    The other half of the wave is the even part (where x(t) = x(-t)),
    you get that by using cosine waves instead of sine waves.
    
    Or if you really want to, you can do it in the complex domain and
    get both parts simultaneously.
    
    	-Bill (rebuilding his system disk)
    

100 rem november 1983 fourier analysis program by Thomas E. janzen
120 twopi = 2 * pi
130 rem convert radians to degrees and degrees to radians
135 radians = twopi/360
140 degrees = 360/twopi
160 Print "This is a Fourier analysis program."
161 Print "It receives an even number of points at equal intervals of a function."
162 Print "Data start at 0 degrees.  It returns sine harmonics with phase."
163 Print "for new data enter n, for test data enter captial T"
164 input data$
165 if data$ = 'T' then samples = 20 \ halfsamples = samples / 2
166 if (data$ <> 'T') then gosub 4000
196 DIM F(samples),A(Halfsamples),B(halfsamples),phase(halfsamples),recovr(samples)
197 dim amplitude(halfsamples)
198 if data$ = 'T' then gosub 1000 else goto 205
205 if data$ = 'T' then goto 280
220 print "Enter the ";samples;" sample points (return after each)"
240 for x = 1 to samples
250 input f(x)
260 next x
280 for m = 1 to halfsamples
290 t=0
300 for x = 1 to samples
305 rem find f(x) times sine(m omega t) - (omega = angular velocity)
310 t1 = f(x) * sin(m * (x-1) * (twopi) / 20)
320 t=t1+t
330 next x
340 b(m) = (fix ((100 * (t/10)) + .5))/100
360 next m
370 for m = 1 to halfsamples
380 t = 0
390 for x = 1 to halfsamples
400 t2=f(x)*cos(m*(x-1)*twopi/20)
410 t=t2+t
420 next x
430 a(m)=(fix((100*(t/10))+.5))/100
450 next m
460 t3=0
470 for x = 1 to samples
480 t3=t3+f(x)
490 next x
500 a0=t3/20
510 print
520 print "a0/2";a0
530 print
540 for m = 1 to halfsamples
545 rem get the phase of the angle in radians and convert to degrees
548 gosub 2000
549 phase(m)=theta*degrees
580 amplitude(m)=sqr(a(m)**2+b(m)**2)
595 print "Sine harmonic ";m;"has amplitude ";Amplitude(m);
597 print "and phase ";phase(m);"degrees."
600 next m
610 gosub 900
650 goto 5000
900 for point =1 to samples
905 recovr(point) = 0
910 for harmonic = 1 to halfsamples
915 rem recover(point) is an acumulator for the function at a point
917 rem this is the fourier transform solutionm expression
919 rem sin(m omega t + thetha) plus the accumulator (omega = 2pif)
920 recovr(point)=amplitude(harmonic)*sin(harmonic*(twopi/samples* &
    (point-1))+radians*phase(harmonic))+recovr(point)
930 next harmonic
935 RECOVR(POINT)=recovr(point)+a0
940 print "recovr";point;"is";recovr(point),"Function was:";f(point)
950 next point
960 return
1000 Print "enter s for sine, c for cosine, t for sawtooth"
1002 input wave$
1004 for x = 1 to 20
1010 if (wave$ = 'S') then f(x)=sin((x-1)*twopi/20)
1011 if (wave$ = 'C') then f(x) = cos((x-1)*twopi/20)
1012 if (wave$ = 'T') then f(x) = 20-x
1020 next x
1030 return
2000 rem this is a subroutine for calculating a full-circle arctangent
2010 IF ((b(m)=0) AND (a(m)=0)) OR ((b(m)>0) AND (a(m)=0)) THEN theta = 0
2020 IF ((b(m)>0) AND (a(m)>0)) OR ((b(m)>0) AND (a(m)<0)) THEN theta=atn(a(m)/b(m))
2030 IF (b(m)=0)  AND (a(m)>0)     THEN theta = pi/2
2040 IF (b(m)=0)  AND (a(m)<0)     THEN theta = 3/twopi
2050 IF (b(m)<0)  AND (a(m)=0)     THEN theta=pi
2060 IF (b(m)<0)  AND (a(m)>0)     THEN theta=pi-atn(a(m)/abs(b(m)))
2070 IF (b(m)<0)  AND (a(m)<0)     THEN theta=atn(abs(a(m))/abs(b(m)))+pi
2080 return
4000 print "enter the number of samples:"
4010 input samples
4020 halfsamples = samples/2
4030 return
5000 end

here are some more digital signal processing references 
From VaxLab signal processing manual:

Antoniou, A. Digital Filters 1975 McGraw Hill

Hamming , R.W. Digital Filters 1977 Prentice-Hall

Horlick, Betty, Analytical Chemistry vol 49, 1977

Madden, Analytical Chemistry vol 50 1978

Oppenheim, Shafer, Digital Signal Processing Prentice-Hall 1975
Savitsky, Analytical chemistry Vol 36, 1964

Steiner, Termonia, Deltour, Analytical Chemistry Vol 44 1972
Tom

    OK, well here we are 2� years later.  I've gotten a little DSP OJT, and
    I have also noticed the proliferation of modern signal processors on
    the market.  Could it be time for a little experimentation?
    
    I have something of an understanding of the concepts presented in
    previous notes.  One thing I noticed is that many (perhaps most 
    or all) of the discussion centered around after-the-fact processing of
    sampled waveforms.  So the first question is, what types of DSP
    algorithms and concepts lend themselves to the real-time signal
    processors in popular use today?
    
    The next question has to do with the types of effects commonly found on
    these new-fangled boxes.  Can somebody explain mathematically or
    intuitively the following:
    
    	-	Reverb (in all it's many flavors)
    	-	Chorus
    	-	Flange
    	-	Phase Shifter
    	-	Various distortions
    
    Perhaps there's a book out by now... "DSP techniques for electronic
    music?"
    
            Marc_thinking_of_connecting_some_wires_and_writing_some_code
    
           

Hello Marc-

Two places to look right off are


(1) "The Element of Computer Music," by F. Benjamin Moore, gives a great
overview of the mathematical basis for the full range of current synthetic
sound generation techniques, as well as for signal processing. It's pretty
heavily oriented towards math and the "c" language, as well the language
"cmusic," invented by Moore (and probably zillions of grad students).

(2) EM (Electronic Music) ran an article called something like "State of the
Art in Signal Processing" sometime this last year. It focused primarily on
non-run-of-the-mill stuff, as I recall, so it may not be quite what you're
looking for. EM has been regularly running "educational" articles for some time
now, and (if you don't already know it) also addresses bitheads on a regular
basis. Sounds like you need to give yourself an early Christmas present.


Hope this is a start.

Cheers,
	Bob

	Reverberation and echo are done with convolution.
	Pitch shifting come to think of it could be done in a similar way
	while skipping or duplicating indices.
	Convolution is key.
	Once you know the unit impulse response of Carnegie Hall when it's
	full of people who brought winter coats, as heard from the middle
	seat in the tenth row, you're in!
Tom

    
    
    
    re .21 - In the real (as opposed to mathematical) world, it's easier to
    do echo by simple delay. ;^)
    
    re .19 - These are all delay based effects.  The basic idea is to delay
    the input signal by some amount, add it to the input, and feed some of
    it back.  I.e.,
    
               
                                               +---+
    IN  >------+------------------------------>| m |
               |   +---+      +-------+        | i |-----------> OUT
               +-->| m |  +-->| delay |---+--->| x |     
                   | i |--+   +-------+   |    +---+     
               +-->| x |                  |               
               |   +---+                  |               
               +--------------------------+
    
    The delay time may be modulated (e.g., by an LFO) to make things more
    interesting.
     
    Some people claim that flanging and phasing differ only in the value of
    the delay time.  Others claim that phasing uses a different circuit.
    
    This simple delay scheme does not do chorusing or reverb very well.  Both
    require many different delays working in parallel with modulation and
    often other signal processing (e.g., EQ).
    
    len.