([WHQVLRQVRIWKH.

DUSOXV6WURQJ3OXFNHG6WULQJ$OJRULWKP
$XWKRU V 'DYLG$-DIIHDQG-XOLXV26PLWK
5HYLHZHGZRUN V 
6RXUFH&RPSXWHU0XVLF-RXUQDO9RO1R 6XPPHU SS
3XEOLVKHGE\The MIT Press
6WDEOH85/[Link] .
$FFHVVHG

Your use of the JSTOR archive indicates your acceptance of the Terms & Conditions of Use, available at .
[Link]

.
JSTOR is a not-for-profit service that helps scholars, researchers, and students discover, use, and build upon a wide range of
content in a trusted digital archive. We use information technology and tools to increase productivity and facilitate new forms
of scholarship. For more information about JSTOR, please contact support@[Link].

The MIT Press is collaborating with JSTOR to digitize, preserve and extend access to Computer Music
Journal.

[Link]

This content downloaded on Wed, 9 Jan 2013 [Link] PM

All use subject to JSTOR Terms and Conditions
DavidA.JaffeandJulius0. Smith
Center for Computer Research in Music and Extensions of the
Acoustics (CCRMA)
StanfordUniversity Karplus-Strong
Stanford,California 94305
Plucked-String
Algorithm

Introduction The String-SimulationAlgorithm

In 1960, an efficient computational model for vi- The Karplus-Strongplucked-string algorithm is pre-
bratingstrings, based on physical resonating, was sented in this issue of Computer Music Journal.
proposed by McIntyre and Woodhouse (1960). This Fromour point of view, the algorithm consists of
model plays a crucial role in their recent work on a high-order digital filter, which represents the
bowed strings (McIntyre,Schumacher, and Wood- string; and a short noise burst, which represents
house 1981; 1983), and methods for calibratingthe the "pluck." The digital filter is given by the dif-
model to recorded data have been developed (Smith ference equation
1983).
Yn-N ++ yn-(N+1)
YN+
Independently,in 1978, Alex Strong devised an =
= Xn +
Yn Yn Yn-N
2 (1)
efficient special case of the McIntyre-Woodhouse
string model that produces remarkablyrich and re- where x. is the input signal amplitude at sample n,
alistic timbres despite its simplicity (Karplusand y, is the output amplitude at sample n, and N is
Strong 1983). Since then, Strong and Kevin Karplus the (approximate)desired pitch period of the note in
have explored several variations and refinements of samples. The noise burst is defined by
the algorithm, with an emphasis on small-system
Au, n = 0, 1, 2,...,N- 1
implementations. We have found that the Karplus- xn = 0, n N,
Strong algorithm can be used with equally impres-
sive results on fast, high-power equipment. The where A is the desired amplitude, and un, [-1,1] is
availability of multiplies, for example, allows sev- the output of a [Link] out-
eral modifications and extensions that increase its put yn is taken beginning at time n = N in our
usefulness and flexibility. These extensions are de- implementation.
scribed in this paper. The developments were moti-
vated by musical needs that arose during the
composition of May All YourChildren Be Acrobats Analysis of the String Simulator
(1981) for computer-generatedtape, eight guitars,
and voice and Silicon Valley Breakdown (1982) for Beforeproceeding to practical extensions of the al-
four-channel, computer-generatedtape, both writ- gorithm, we will describe the theory on which
ten by David [Link] theoretical approachand many of them are based. Various concepts from
the extensions based on it have also been applied to digital filter theory are employed. For a tutorial in-
the McIntyre-Woodhousealgorithm (Smith 1983). troduction to digital filter theory, see the works by
Smith (1982b) and Steiglitz (1974).
The input-output relation of Eq. (1) may be ex-
David A. Jaffeis also affiliatedwith the Music Departmentat
StanfordUniversity,and Julius O. Smith is also affiliatedwith
the Electrical EngineeringDepartment there.
1. In some situations, the sound more closely resembles a string
ComputerMusic Journal,Vol. 7, No. 2, struck with a hammer or mallet than one plucked with a pick,
Summer 1983, 0148-9267/83/020056-14 $04.00/0, but we will always use the term pluck when referringto the
? 1983 Massachusetts Institute of Technology. excitation.

56 ComputerMusic Journal

This content downloaded on Wed, 9 Jan 2013 [Link] PM

All use subject to JSTOR Terms and Conditions
 Fig.1. Blockdiagramfor
the basic stringsimulator.

pressed differently by means of delay-operatornota-

tion. We define the unit-sample delay operatord by
the relation
dkx A-Xn-k,

where x, is an arbitrarysignal, and k is an integer.

(The symbol A means "is defined as.")Thus, multi-
plying a signal by dk delays the signal in time by k
samples. In these terms, Eq. (1) becomes
dNy + dN+lyn
yn =
tor is the following set of difference equations:

L+d Vn = Yn-N
=x+n 2
V= n + Vn-
Solving for yn yields 2
Xn1 yn = xn + Wn
Yn XI+d' (2)
1- dN The frequency response of a digital filter is de-
2 fined as the transferfunction evaluated at z =
We can convert linear delay-operatorequations e's = cos(oT,) + j sin(o Ts),where T, is the sam-
in seconds (Tsis the inverse of the sam-
immediately to z-transform equations by replacing pling period =
each time signal with its z-transform, and replacing pling rate fj, o 27rfis radian frequency,f is
d with z-'. It is customary to denote a time signal frequency in Hz, and j = \-7. The frequency re-
in lowercase letters (e.g., x) and the corresponding sponse of the string simulator is then
z-transform in uppercase letters (e.g., X(z)). The 1
H(e,Ts) =
transferfunction of a (linear, time-invariant) digital He') 1 - Ha eiw sHb(es)'
filter is the z-transform of the output signal divided
where
by the z-transform of the input. The transferfunc-
tion of the string simulator is then found to be 1 + e- 'T
H|(e
T = ee-' 2cos(coT,/2)
2Ts)
1
H(Z) A Y
Y(z) _=X
H X(z)x 1 l+z - = e-fTscos(rfTs)
2
Hb(es s) = e-iN = e- [Link]
1
1 - Ha(z)Hb(z)' In this paper it is necessary to consider the am-
plitude response and phase delay of the feedback fil-
where ters [Link] amplitude response is defined
as the magnitude of the frequency response, and
1 + z- it gives the gain of the filter as a function of fre-
Ha(z) Q 2
2
quency. The phase delay is defined as minus the
Hb(z) z-N complex angle of the frequency response divided
by radianfrequency, and it gives the time delay
This form of the description is shown in Fig. 1. The (in seconds) experienced by a sinusoid at each
feedback loop consists of a length N delay line frequency.
Hb(z) in series with a two-point average Ha(z). Cor- The amplitude response of each component filter
responding to this breakdown of the string simula- is given by

Jaffeand Smith 57

This content downloaded on Wed, 9 Jan 2013 [Link] PM

All use subject to JSTOR Terms and Conditions
 ___

Ga(f)A IHa(eiTs)I= cos(oT,/2)l = Icos(T7fT)l NTs + Ts/2 sec. Experience shows this to corre-
spond well with perceived pitch.
Gb(f) A IHblee
ws)I 1.
Thus the delay line Hb is lossless, and the two-
point average Ha exhibits a gain that decreases with Decayof "Harmonics"
frequency according to the first quadrantof a co-
sine. We will assume hereafter that all frequencies Since the signal is only quasi-periodic,it does not
are restricted to the Nyquist limit, that is, Ifl consist of discrete sinusoids. Essentially, we have
fs/2. In this range, we have (cos(i7fTj)= cos(1rfTs).
many narrow "bands"of energy decaying to zero at
It is convenient to define phase delay in units differentrates. When these energy bands are cen-
consisting of samples rather than seconds. The tered at frequencies that are an integer multiple of a
phase delays of Ha and Hb in samples are given by lowest frequency, they will be referredto as har-
monics. When the frequency components are not
Pa(f)
"- - / H(e'2)
~oT,
=1
2'
necessarily uniformly spaced, the term partial will
be used to emphasize the possibility of inharmo-
nicity. Consider, then, a partial at frequency f Hz
Pb(f) A - / Hb(e'S) = N.
PfT, circulating in the loop. On each pass through the
loop, it suffers an attenuation equal to the loop-
(Lz denotes the complex angle of z). The two-point
amplitude response, Ga(f)Gb(f) = cos(IrfTs);that is,
averagehas a phase delay equal to half a sample,
and the delay line has a phase delay equal to its one period's attenuation = cos(rfTs).
length. Since the round-trip time in the loop equals N +
Since the total loop consists of Haand Hbin se-
1/2 samples, the number of trips through the loop
ries, the loop gain and effective loop length are after n samples (nTssec) is equal to n/(N + 1/2) =
loop gain = G(f) Gb(f) = cos(rfT,), tfS/(N + 1/2). Thus the attenuation factor at time
t = nT, is given by
and
tf5
loop length = Pa(f) + Pb(f) = N + 1/2 (samples)
af(t) A [cos(7rfTs)]N+ 2. (3)
for each sinusoidal frequency f Hz. Forexample, an initial partial amplitude A at time
In synthesizing a single plucked-string note, we 0 becomes amplitude Aaf(t) at time t seconds,
feed in N samples of white noise at amplitude A where f is the frequency of the partial.
and listen to the output immediately [Link] The time constant of an exponential decay is tra-
is equivalent to initialize the delay line Hb with ditionally defined as the time when the amplitude
sealed random numbers at time 0 and employ no has decayed to l/e - 0.37 times its initial value.
input signal. Since the two-point averageHais con- The time constant at frequency f is found by equat-
stantly changing the contents of the loop, the out- ing Eq. (3) to e-'/f and solving for rf, which gives
put signal is not periodic. It is close to periodic,
however, and we use the term period in this loose Zs
sense. Each period of the synthetic string sound
= -t - ( )
correspondsto the contents of the delay line at a Tf In af(t) = ncos(rfT (seconds). (4)
In cos(rTrfT)
particulartime, and each period equals a somewhat
lowpass version of the previous period. More pre- Foraudio, it is normally more useful to define the
cisely, a running two-point averageof the samples time constant of decay as the time it takes to decay
comprising one period gives the next period in the -60 db, or to 0.001 times the initial value. In this
output waveform. Since the effective loop length is case, we equate Eq. (3) to 0.001 and solve for t. This
N + 1/2 samples, the period is best defined to be value of t is often called t6o. Conversion from f to

58 ComputerMusic Journal

This content downloaded on Wed, 9 Jan 2013 [Link] PM

All use subject to JSTOR Terms and Conditions
 Fig.2. Spectralevolution
duringthe first 16 periods.

t60(/)is accomplished by
t6o(f)= ln(1000),rf 6.91 f. (5)
Forexample, if a sinusoid at frequency f Hz has
amplitude A at time 0, then at time t,o(f)it has am-
plitude Aaf(t6o (f)) = A/1000, or it is 60 db below its
starting level.
The above analysis describes the attenuation due
to "propagation"around the loop. It does not, how-
ever, incorporate the fact that sinusoids that do
not "fit"in the loop are quickly destroyed by self-
interference. This situation is analogous to making
an actual string vibrate. Any signal may be "fed
into" the string, but after the input ceases, the re-
maining energy quickly assumes a quasi-periodic
nature. Thus, even though the loop is initialized
with random numbers, after a very short time the monic amplitudes, since a straight line is drawn
from one harmonic amplitude to the next.
primary frequencies present in the loop are those
that have an integral number of periods in N + 1/2 In certain extensions to the algorithm, HAis
other than a two-point [Link] such a case, the
samples. These frequencies are all multiples of the
attenuation factor of the kth harmonic after t sec-
frequency whose period exactly matches the loop
onds is approximately
length N + 1/2. This lowest frequency provides the
fundamental, or pitch frequency, of the note: tf,

ak(t) = G,(fk )N + Plfk), (10)

A (6)fs
where we require Ga(fl) 1 for stability. The phase
(N+2) 2T~ N4+-
12 delay Pa(fk)of Hamay be used to create inharmonic
Setting f to the harmonic series beginning with f,, spectra (Smith 1983). The spectrum is harmonic
only when Pa,fk)is the same at all harmonic
fs
f, = la, E k k= 2, ... N/2, (7) frequencies fk.

N+ 2 Similarly, when Hais more general, the time con-

2 stant of decay for each harmonic becomes
gives the decay factor at time t for the kthharmonic t
N + P,(fk)
to be Tk= -
f
T
In G(2rfkT)
econ
seconds. (11)
ak(t) = [costlrfkTs)]ft. (8) Having provided an analytic vocabulary,we now
proceed to a detailed examination of our additions
Similarly, the time contrast per harmonic is to the Karplus-Strongalgorithm.
given by
_-t 1
k In ak(t) f In cos(rfkT,) (seconds. (9 Tuning
Figure 2 shows the spectral evaluation during the The fact that the delay-line length N must be an
first 16 periods of a note having a period of 128 integer causes tuning problems. Since the funda-
samples. A 128 length Fast FourierTransform(FFT) mental frequency is f, = f/(N + 1/2), the allowed
was computed every other period. Each curve in the pitches are quantized, especially at high frequency.
figure is interpreted as the envelope of the har- Forlarge values of N (low pitches), the difference

Jaffeand Smith 59

This content downloaded on Wed, 9 Jan 2013 [Link] PM

All use subject to JSTOR Terms and Conditions
Fig. 3. Desired pitch versus
resultingpitch for a 50-
KHz sampling rate.

whole string is now

- 20- H(z) )A 1 - Ha(z)Hb(z)Hc(z)

The filter H, is a first-orderallpass filter, and as

a -/ such it has a constant amplitude response. Indeed,
the amplitude response is simply
4
10- / Gc(f) - IHc(e s)l= 1 + Ce-7T| 1
N

The use of an allpass filter ensures that no modi-

0
CO
O - fication of the decay rate will take place. The loop
gain is G,(f)Gb(f)G,(f) = cos(rTfTs)as before.
We will select the phase delay of H, so as to tune
5 10 15 20 25 fl to the precise desired frequency. This requires
Desiredfrequency(KHz) only the ability to select phase delays between 0
and T, sec, or one sample's worth.
The phase delay of the first-orderallpass H, is
between the pitch at N and N + 1 is very slight.
given by
However, for high pitches, N and N + 1 yield very
differentpitches and tuning becomes crude. Figure
PC(f) A _ /Hc (e's)
3 shows the distortion in frequency for the sam- coTT
pling rate 50 KHz. For lower sampling rates, the -1 C + e-'"T
curve is identical in form, and the distortion occurs
at proportionatelylower frequencies. coTs 1 + Ce-'T
The key to a solution for this problem lies in the L(1 + Ce-"T) _ (C + e-'"T)
expression for loop length in terms of phase delay. coTs coT
The fundamental frequency is given by
= 1tan(- C sinlT)
_
= fs
f N+Pa(fl)'
oT tanI -sin(cT)
tn (13)
where Pa(f) = 1/2 when the two-point average is oT,s C + cos( oTj)
used for Ha. To make up the difference between f
When the arguments to the arc tangent above have
and the desired frequency, we need to introduce
into the feedback loop a filter that can contribute a magnitude less than unity, we can use the power-
series expansion (Abramowitzand Stegun 1966),
small delay without altering the loop gain. The fil-
ter we introduce has the difference equation
tan-'(x) = x --
X3 XS X7
+ + ., Ix< 1.
Yn = CXn + Xn-1 - CYn- (12)
Thus we can approximate the low-frequency phase
and transferfunction
delay by
H,(z)A 1 + Cz-' Pcf) = sin(oT,) C sin(co
T
cT,()C + cos(oT,)) oT,(1 + C cos(oTj))
where C is the only coefficient to be set. Forsta-
1 C 1-C
bility, we must have ICI< 1. It can be shown that C+ 1 +C 1+C' (14)
when the input x, is bounded by 1, the output is
bounded by 21CI+ 1. The transferfunction of the A plot of the exact phase delay is given in Fig. 4 for

60 Computer Music Journal

This content downloaded on Wed, 9 Jan 2013 [Link] PM

All use subject to JSTOR Terms and Conditions
Fig. 4. Phase delay for the
fine-tuning allpass filter
H,(z) = (C + z-')/
(1 +C z-).

2-
and the delay Pc(fL)requiredfrom the allpass filter
become
N A Floor(P, - P(f,) - ) (15)
Pc(fl)A Pi - N - P,(f),
c :
4

where s > 0 is the offset that shifts Pc(fl)into the

' range [e, 1+e], and Pa(f,)is the delay in samples due
cn
-1 _ ...-~
to the filter Ha. In the simple case where Hais a
two-point average, Pa(f,)= 1/2.
C We next solve for the filter coefficient in Eq. (13)
as a function of Pc(f,).Taking the tangent of both
sides, and using an identity for the tangent of a
difference leads to the quadraticequation in C,
n
j C2 sin(co1TJP,(f,)+ cw,Tj + 2C sin(co,TPc(fA)) +
I. I I I I ?
- =
sin(co,TP,(fl) coTj) 0,
0 Frequency F,/2
where co A_27rf. The solution is found, after some
manipulation (Mont-Reynaud 1982), to be
17 values of C equally spaced between -0.999 and
0.999. Note that delays between 0 and 1 sample can -sin(ol TsP(fL))+ sin(olTs)
be provided somewhat uniformly across the fre- sin(wolTsP(f ) + coTj)
quency axis. A delay of 0 samples correspondsto We have introduced an extra root by producing a
C = 1, where the pole and zero of Hc(z)cancel to quadraticequation. The previous approximation
give Hc(z) 1. However, pole-zero cancellation on Eq. (14) indicates that the + sign should be taken.
the unit circle is not a good thing in practice, since Therefore,the final solution is
round-offerrors may yield an unstable filter. There-
fore, it is preferableto shift the range of one-sample - sin(c, Tj - sin(o, TsP,c(f))
delay control to the region s
< Pc < (1 + e) for sin(cl TsPc(fl)+ c Ts)
some small nonnegative e (0 < < 1). It is best not co T - o,
to shift very far, since the phase-delay curves are -sin( TsPc(fl)
less flat in the region beyond one sample's delay. = (16)
Note that the delay curves below the one-sample i( TT,+co, TsPff,
level in Fig. 4 correspond to slightly flattened up-
per partials, while the delay curves above the one- which can be approximated, at low frequencies, by
sample level correspond to slightly sharpenedupper 1 - Pc(f,)
partials. The timbre change due to slight systemat- i + Pc(AM) (17)
ic shifting of the upper partials, of an amount less
than one sample period, was found to be hardly no- Although this technique provides a precise funda-
ticeable. It may even be desirable as a source of sub- mental frequency, it does not guaranteean in-tune
tle timbral variation. In any case, it is important to percept, since the perceived pitch does not always
get the perceived pitch right. coincide with the fundamental [Link] ad-
To tune the instrument precisely to a desired ditional mapping onto a perceptual tuning dimen-
fundamental frequency f, let P1 equal fs/f, the real sion may be needed for the very high notes. In our
value for the period of the first partial, in samples, case, it was found that tuning the octaves slightly
which would give perfect tuning. Then we desire stretched, as is done in piano tuning, gives a more
N + Pa(fi) + P,(fl) = P. The integer buffer length N satisfying in-tune percept.

Jaffe and Smith 61

This content downloaded on Wed, 9 Jan 2013 [Link] PM

All use subject to JSTOR Terms and Conditions
 Alteration
Decay-time rable in duration with notes from a real string. With
the loss factor operative, the decay-time constant
The basic algorithm naturally results in a shorter for the fundamental frequency becomes
decay time for high pitches than for low pitches, 1
reflecting the behavior of real strings. This is due = - (19)
r(p) f, In lp cos(rTf,Tj)I'
to two effects. First, higher frequencies are more
attenuated by the two-point average Ha; second, Decay shortening produces a dampedversion of the
higher pitch means more trips through the atten- algorithm, analogous to substitution of a soft mate-
uating loop in a given time. This may be seen in rial for the bridge of a string instrument.
Eq. (3).
Unfortunately, the range of decay times between
the high and low pitches is too extreme. The high- DecayStretching
pitched notes die away so fast that only a click is
perceived, while the low-pitched notes last for an To stretch the decay, the feedback average(Ha) can
unnaturally long time. In addition, the decay time be changed to a two-point weighted [Link]
of real strings varies with many factors such as ten- reduces the amount of loss at high frequencies.
sion, length, thickness, and material. Consequently, Thus, we replace H(z}) by
we have found it useful to add a means for altering
H,(z, S) = (1 - S) + Sz-, (20)
the note duration. The ability to control decay time
is essential for a realistic simulation as well as for where S, the stretching factor, is between 0 and 1.
musical flexibility. The gain of this filter is
On systems that have separatecontrol of the
G,(f, S) = (1 - S) + Se-'it
sampling rate of each voice, sampling-ratechange
can be used to control decay time. On other sys- = V ((1 - S) + S cos oTJ)2+ (S sin(w)T,))2
tems, however, it is necessary to use other methods = V (1 - S)2 + S2 + 2S(1 - S) cos wTr
to alter decay time.
(21)

DecayShortening With S = 1/2, Ha(z, S) reduces to the previous case

Ha(z). Forstability of the overall string, we must have
To shorten the decay time, a loss factor p can be in- 0 < S < 1. If S = 0 or 1, the frequency-dependent
troduced in the feedback loop. With the loss factor, term disappears,and the gain response is unity for
the difference Eq. (1) for the string becomes all f; in this case, the initial white-noise burst cir-
culates forever in the loop, producingharmonics
(18) that never decay. At intermediate values, 0 < S < 1,
-N
= X + Yn-(N+)18
y,=fx,,+ p Yn
the effective note duration (t60)is finite, and it is
The amplitude envelope of a sinusoid at frequency minimum for S = 1/2. The amplitude trajectory
f, previously given by Eq. (8), is now proportionalto and the decay-time constant for each partial can be
obtained by substituting Eq. (21) into Eqs. (10) and
a,(t, p) = |p cos(TrfTs)jfl= IPIM'taf(t).
(11), respectively.
Thus all partials are affected equally: the relative For the greatest control, both the uniform-loss
decay rates are unchanged. method and the weighted two-point-averagemethod
Note that p cannot be used to lengthen the decay may be used for decay-time alteration. The result-
time, since the amplitude at 0 Hz would increase ing decay time is then a function of loss factor p
In we must have - 1 if and stretch factor S. Karplusand Strong (1983) de-
exponentially. general, Ipl
the string is to be stable. Thus p is used to shorten scribe a method of decay stretching that uses no
the low-pitched notes to make them more compa- multiplies.

62 ComputerMusic Journal

This content downloaded on Wed, 9 Jan 2013 [Link] PM

All use subject to JSTOR Terms and Conditions
Fig. 5. Phase delay for the
decay-stretching, one-zero
filter Ha(z, S) = (1 - S) +
S z-'.

Note that for S > 0 the phase delay is quite flat

over most of the frequency axis. Another point
of interest is that since Ga(z, S) = Ga(z, 1 - S),
2/
we may choose the case that yields the best phase
-s /4 - / delay curve for the fine-tuning allpass H,.
Since the tuning calculation needs to be done
only once per note, the precise form of Eq. (22)
can be used for each new frequency without much
computational expense.

Dynamics
The loudness of the signal output by the algorithm
0 Frequency F,/2 is a function of the amplitude of the input noise
burst. However, this is an unsatisfactory control in
simulating the timbral effect of dynamic level as it
occurs in the case of a real string instrument. The
effect of varying initial amplitude gives the impres-
Effectof DecayStretchingon Tuning sion of a change more in the distance between the
listener and the apparentsource than in dynamics.
ChangingS changesthe effectiveloop lengthas a Since strings plucked hard have more energy in the
functionof frequency,since it changesthe phase
higher partials than strings plucked lightly (due to
delayof the overallloop. We must thereforecom- nonlinearities becoming important), the dynamic
pute Pa(f,)for use in Eq. (15) when fine tuning with simulation is based on modeling this differencein
the allpass filter H,. The phase delay of the weighted
spectral balance. We therefore change the effective
two-point average is given by spectral bandwidth of a note to modulate its appar-
ent intensity.
Pa(f, S) - ?Ha(e'TS, S) The bandwidth is controlled by means of a one-
coT,
pole, lowpass filter applied to the initial noise burst
((L - S) + Se-' T)
(beforeit is fed into the string). This filter will be
coT, referredto as the dynamics filter. The difference
_ - tan-'
tan- -Ssin(oT) ,2
(22) equation of the dynamics filter is
, -
(1 S)+ S cos(ouT) )
yn = (1 - R)xn + Ryn-,
and for low frequencies, relative to the sampling and its transfer function is
rate, we may use the approximation
1-R
S)
P(f, T -S sin(oTs) Hd(z) = 1 - Rz-' (23)
Pa( f T,({1- S) + ScoTs cos()oT,)
where R is a real number between 0 and 1, com-
- S, O < S 1.
puted as a function of fundamental frequency fl and
For S = 1/2, we have the basic string algorithm, the desired dynamic level L. When a series of notes
and the phase delay of Ha is 1/2, as given by the at pitch f, is played while R is moved graduallyto-
above approximation. For other values of S the ward 1, a diminuendo is approximatedin terms of
approximation is always precise at f = 0. Figure5 both decreasing loudness and spectral bandwidth
shows the true phase-delay curves of H,(z, S) as S reduction.
is stepped uniformly through 17 values from 0 to 1. We define the dynamic level L as a bandwidthbe-

Jaffe and Smith 63

This content downloaded on Wed, 9 Jan 2013 [Link] PM

All use subject to JSTOR Terms and Conditions
 Fig. 6. Frequencyresponse
of H(z) = (1 - R)
(1 - R z-l) computed at
dynamic level L = 100 for
six values of f,.

tween 0 and f/2. If L is small, the spectrum is more

lowpass filtered, corresponding to a softer dynamic
level. Conversely, large L gives a bright spectrum
corresponding to louder notes. It is not sufficient to
o -25t1
-- 25
use a fixed lowpass filter for all pitches, since low-
pitched notes would then be louder than high-pitched
notes. Rather, for a given dynamic level, R must be
changed with pitch to yield a uniform perceived -50= 100 Hz
loudness. While this is a difficult problem in gen-
eral, a good approximation is obtained by varying R 50123I ....
so that the amplitude of the fundamental frequency 0 1 2 3 4
is constant. fi fm Frequency (KHz) f.
It remains to be shown how R is computed for a
given pitch fl and dynamic level L. The main steps
are as follows. First, a one-pole, lowpass filter is de-
signed having bandwidth L. Second, the gain of this
filter at a "middle" frequency is computed. Third, Now, for any desired fundamental frequency fL, R
is computed so as to provide gain Gd(fl) = [Link] other
the dynamics filter is computed as a one-pole, low-
words, all fundamental frequencies are made to have
pass filter having this gain at the desired fundamen- the same amplitude. The value of R is found by
tal f,. The remainder of this section gives the equa-
tions needed for these steps. solving
The reference frequency fmis chosen as the loga-
GL =
l -R
R
rithmic middle (geometric mean) of the range to GL - Ts
|I1 Re-i2fI
be used (a function of the particular musical con-
text and the sampling frequency): Squaring both sides of this equation and solving the
resulting quadratic polynomial in R yield
fm= e(l/2)(log(f)+ log(fl)) = /-f

R= 1 - GL cos(2rfr Ts)
R
where f, is the upper pitch limit (< f/2), and fi is 1 - GL
the lower pitch limit.
+ 2GL sinr rf T
V
V' 1 -
T1 Ts)
G2L COS2(rTf, )
The one-pole lowpass filter having bandwidth L 1- G-
is given by
We use whichever value is < 1 in magnitude to en-
I -
HL(Z) A L - RL sure stability.
1 RZ-'
A family of frequency-response curves for Hd is
where shown in Fig. 6 for six fundamental frequencies in
octave steps from f, = 100 Hz to fl = 3,200 Hz. The
RL A e- LTs.
dynamic level in each case is L = 100 Hz. A verti-
cal line is drawn to each curve at the fundamental
The substitution RL = e- LTis a somewhat stan-
dard approximate formula for mapping bandwidth frequency to which it applies. The reference fre-
to pole radius. quency fmis set to 282.84 Hz (the geometric mean
of f, = 20 Hz and fu = f,/2), and the sampling rate is
The gain of the lowpass filter HLat the reference
fs = 8,000 Hz.
frequency is defined as To add to the effect of simulated dynamics, it is
sometimes helpful to do a bit of decay shortening
GL A IHL(e"7sfms)
on the low soft notes, using a loss factor p as de-
1 - RL scribed previously. It is also possible to simulate
11 - RLe-'21rmTs' the spectral characteristics of soft notes by simply

64 Computer Music Journal

This content downloaded on Wed, 9 Jan 2013 [Link] PM

All use subject to JSTOR Terms and Conditions
turning the output of the algorithm on late, after cay shortening, is set to a relatively small value in
some of the high-frequency energy has died away. the last few milliseconds of the note. The duration
This latter method has the effect of diminishing the of the decay is dependent on the loss factor chosen.
attack noise, corresponding to a milder pluck. For A loss factor close to 1 simulates a string being
other decay-softening techniques, see the section dampedwith a soft material such as the flesh of the
entitled VaryingAttack Characteristics. finger, while a smaller loss factor simulates damp-
ing with a hard material, such as a pick. It is useful
to compute p as a function of a desired [Link]-
Restsandthe Endsof Notes tuting Eq. (19) into Eq. (5) and solving for p yields
e- l/fT
In the case of a real string instrument, the player p(t60) = . A t6o
cot6
s(=ICOSfT)I' = ln(1000)'
often plays two notes on the same string, playing
the second note before the first has died away. The Forpitches above about 3 KHz, clicks can appear
basic algorithm handles this without a problem for even at the onset of notes. Onset clicks can also
all but very high pitches, since the discontinuity at occur when a pianissimo setting of the dynamics
the beginning of the new note is perceived as a new filter is used, since the masking effect of the on-
pluck ratherthan as a click. (Witha real instrument, set noise burst is weakened by the dynamic filter.
there is a short time when the pick has muted the To alleviate this problem, it may be necessary to
vibration for the previous note but has not yet set multiply the output of the algorithm by an expo-
the string again into vibration, but this brief silence nential or linear envelope, rather than switching it
was not found to be necessary for realism in the on abruptly.
synthetic model.) A problem arises, however, when
there is a rest after a note. If the output of the algo-
rithm is abruptly turned off (replacedby zeros), a Glissandi
andSlurs
discontinuity in the waveform results, causing a
click. A slur can be simulated by changing the orderof
Even if the note is allowed to decay for a very the filter without reexciting it; that is, by chang-
long time, a note turned off abruptly may cause a ing the delay buffer length N. The result is analo-
click because the feedback loop has unity gain at gous to a performer'srefretting a string without re-
0 Hz. To see this, note that Ga(O)= cos(O)= 1. Thus plucking it (what guitarists call the "hammer-on"
the final value of the waveform is the mean of the and "pull-off"technique). A rapid alternation of as-
initial noise burst. Since the pseudo-random-number cending and descending pitch changes gives a good
sequence used to initially excite the filter has a left-hand-trill effect.
mean of zero only in special cases or when an in- If the buffer length is gradually changed over
finite number of samples is taken, the 0 Hz compo- time, a crude glissando results. For low pitches and
nent can be significant. For uniform pseudorandom high sampling rates, it can be quite smooth, but for
white noise at amplitude A, the statistical variance higher pitches and lower sampling rates, the pitch
is A2/3, which implies a variance in sample mean quantizationresulting from the integer-delaylengths
over N samples equal to A2/(3N). Thus the stan- becomes noticeable. In the upper range, the effect
darddeviation of the mean in a length of N sam- is similar to a glissando on a fretted instrument,
ples is A/\V3N. As N decreases, the probability of where the pitch changes in discrete steps, and it
having a large amount of energy at 0 Hz becomes can be musically useful. Of course, the synthetic
[Link] Strong describe a technique quantization is not in semitone intervals.
they call "dithering"to handle this problem, but A perfectly smooth glissando can be created by
this technique was not economical in our context. ramping C, the tuning coefficient, during the time
Our solution to the discontinuity problem is as between buffer-length changes. This technique can
follows: the loss factor p, as was discussed for de- also be used to create vibrato.

Jaffeand Smith 65

This content downloaded on Wed, 9 Jan 2013 [Link] PM

All use subject to JSTOR Terms and Conditions
 StringSimulation
Sympathetic that a loss factor p, such as was introduced with
referenceto decay shortening, be used to provide
In a real string instrument, sympathetic vibrations energy dissipation.
of other open strings, as well as resonances in the The effect of several sympathetic strings can be
instrument body, give each pitch in the range of the created simply by a bank of parallel sympathetic
instrument an individual characterand thus give strings, as defined above, each tuned to a different
the instrument as a whole a distinctive identifi- [Link] resulting overall string-transferfunc-
able [Link] comparison, the basic form of the tion (omitting fine tuning and decay alteration for
Karplus-Strongalgorithm, like many instrument- clarity) is then
simulation algorithms, has an excessive homogene-
ity of characterthroughout its range. One way to H(z) = l+z
remedy this situation is to feed the output of the 1-
22
string into a body resonator. This technique has
producedimpressive results, but its discussion is M -i
beyond the scope of this paper (see Smith 1983). (1 -
Another way to combat unnatural uniformity is i I "il-p."
p) 2 z-z-)N
to create the effect of an instrument with sym-
pathetic strings, using a modified version of the where y is the fractional part of the plucked-string
basic algorithm. signal sent to the sympathetic strings, M is the
Just as a sympathetic string is set into motion by number of sympathetic strings, pi is the loss factor
the vibration of another string, the illusion of a for the ith sympathetic string, Ni = fJ/fi where f, is
sympathetically vibrating string can be created by the fundamental frequency of the ith sympathetic
exciting one copy of the string simulator by a small string, and fs/N is the pitch of the plucked string
percentageof the output from another (plucked) and hence of the played note.
string, tuned to a different pitch. The sympathetic-string version of the algorithm
In the discussion that follows, the algorithm that is also helpful in creating a stereo or quadraphonic
is excited with the noise burst is referredto as the image. By distributing around the room the outputs
plucked string and the algorithm that is excited of several banks of differently tuned sympathetic
only by the plucked string is referredto as the sym- strings, all fed with the same plucked string, the
pathetic string. All partials of the plucked string effect of being inside a huge guitar can be created.
that do not coincide with those of the sympathetic Attractive musical results have been created by
string will be highly attenuated. Thus the sympa- replacingthe plucked string with another computer
thetic string acts as a bank of very narrow bandpass instrument, so that the bank of sympathetic strings
filters with center frequencies at the partial fre- is used as a "reverberator."One can achieve the
quencies of the sympathetic string. The partials of effect, for example, of a clarinet being played into
the plucked string that will strongly resonate are an open grandpiano with the pedal down.
those for which
fi = fA, Simulation
of a MovingPick
where fkis the frequency of the kthpartial of the
sympathetically resonating string, and f, is the fre-An effective means of simulating pick position is to
quency of the jth partial of the plucked string. introducezeros uniformly distributed over the spec-
A problem can arise after several successive noise trum of the noise burst. This can quite accurately
burstshave excited the plucked string. The repeated simulate the effect of plucking a string at varying
reintroductionof energy into the sympathetic string distances from the bridge. The noise burst is fil-
may cause it to overflow. Therefore, it is essential tered with a comb filter, H,, having the difference

66 ComputerMusic Journal

This content downloaded on Wed, 9 Jan 2013 [Link] PM

All use subject to JSTOR Terms and Conditions
equation ner, the steady state fuses into a single note. In a
real mandolin, the strings are never perfectly tuned,
Yn = x, - Xn-,,N
and the beating effect of the slightly mistuned strings
where /, is the fraction of the string between the is a strong recognition cue. A better mandolin sim-
bridge and pluck point. When /u = 1/2, the even ulation simply uses two parallel forms of the algo-
harmonics are removed, and the effect is that of rithm, differing in pitch by a few cents and excited
plucking a string at its midpoint. Similarly, when at slightly different times.
A = 1/10, every tenth harmonic is suppressed,and
the effect is like plucking a tenth of the way up the
string. With /A= 1/N, the filter approximates a dif- Useof OtherFiltersin FeedbackLoop
ferentiator,creating a sharp sul ponticello sound.
Forthe theory behind the simulation of pick loca- The use of filters other than a one-zero for Hawill
tion, see Smith's paper (1982a). give a different decay characteristic and, in turn, a
different timbre. However, care must be taken that
the amplitude response G, does not reach unity near
the Character
Varying andNumber
of Attacks any partial frequency. Energy at any frequency fk
for which Ga(fk) = 1 will never decay, and if G,(fk)
Since the attack is very important in perception of > 1, the amplitude will grow exponentially until
timbre, it is advantageous to be able to alter its overflow.
[Link] give a noticeably more noisy attack, Forexample, a one-pole filter with the pole be-
approximatingthe sound of a snap or "Bart6kpizzi- tween 0 and 1 gives a tone with the same attack
cato," the duration of the noise burst x, can simply and a more mellow decay than with a one-zero fil-
be increased from tx = NTs to some tx > NTs. Simi- ter; that is, the higher partials decay more rapidly.
larly, the attack can be subdued by making 0 < tx < Placing the pole at z = Q and normalizing the peak-
NT,, though for very small tx the pluck illusion (as amplitude response of the one-pole filter to unity
well as loudness) fades. yield the transferfunction
A variety of other methods can soften the attack.
The string can be excited with a rich harmonic H()
- IQI
= 1- Ql"
spectrum rather than a noise spectrum, or with
some mixture of the two, with the sum of their am-
where IQI< 1 is required for stability of the one-
plitudes not exceeding 1. Another possibility is to
pole filter. The transfer function of the whole string
lowpass filter the noise burst. Yet another way to becomes
soften the attack is to turn on the output of the
algorithm late, after some of the high-frequency en- 1
H(z) = H
ergy has been filtered out. 1z-z 1-IQI
A realistic simulation of the up-and-downpick- 1- Qz-
ing pattern characteristic of a mandolin tremolo 1 - Qz-'
has been created by using a one-pole, lowpass filter, 1 - Qz- - (1 - IQI)z-
H,, to mellow the "up"picks while using the stan-
dardunfiltered noise burst for the "down"picks. Thus the difference equation is
A crude simulation of instruments having multi-
Yn = Xn - Qx,n- + QYn-1 + (1 - IQI)yn-N
ple strings tuned in unison, such as the mandolin
or bazooki, can be created by simply exciting the This version of the algorithm, with 0 < Q < 1 is
string with two successive noise bursts separatedby useful for pitches in the lower half of the [Link]
a short amount of time, on the order of .05 sec. the upper range, the notes die away too fast to be of
While multiple attacks can be achieved in this man- use. This is because the one-pole, lowpass filter,
with Q > 0, filters out the high-frequency energy

Jaffeand Smith 67

This content downloaded on Wed, 9 Jan 2013 [Link] PM

All use subject to JSTOR Terms and Conditions
Fig.7. Blockdiagramcon-
tainingalgorithmexten-
sions discussedin this
paper.

To sympatheticstrings k _ kf,s(k)
fk
Po + Pg(fk) s(l)
0 n)
sl)
_pI
Pg(fk)= f s( )
f1s(k)
X(n) HeHfO 0 Ha
where Pois the length of the loop in the absence of
the allpass Hg(typically P + 1/2). Methods for de-
Hg H,
signing allpass filters with prescribedphase delay
are reviewed by Smith (1983).

in the loop much more drastically than does the Summary

one-zero, lowpass filter. A technique similar to that
used in dynamics simulation could be used to com- Figure 7 shows a block diagramof the string simu-
lator with some of our revisions, where Hais the
pensate for this trend. feedback lowpass filter, Hb is the delay line, H, is
the allpass filter used for tuning, Hd is the lowpass
filter used to simulate dynamics, He is the comb fil-
Simulation of StiffStrings ter that simulates pick position, Hf is the filter that
simulates the difference between "up"and "down"
The spectral components of the basic algorithm
picking, and Hgis the allpass filter used to simulate
that have significant amplitude are almost perfectly
string stiffness.
harmonic after the attack noise has been filtered
The simulator provides a high degree of flexibil-
away,corresponding,in the real world, to a per- that begins to approachthat of a skilled player
ity
fectly flexible string. But since real strings always on a real musical instrument. Many as-
have some degree of stiffness, it is desirable to alter performing
pects of a real string instrument have been simu-
the spectrum of the algorithm [Link] the-
lated. Pitch can be precisely specified, and articula-
ory of stiff strings (Morse 1976) indicates that stiff- tion can be finely tuned. An expressive vocabulary
ness creates a stretching of the partials accordingto
is provided by a wide variety of performancenu-
the approximate formula
ances, including such "left-hand"techniques as
fk = kfo 1 + 1+ k2 Tr A kfos(k), glissandi, slurs, and trills, as well as such "right-
8 282) hand" techniques as variation in dynamic level,
4 pick position, and attack characteristics. These pa-
k = 1, 2,..., k2 < 2 rameterswere found to be sufficient to create shaped
musical phrases. Furthermore,parameterssuch as
where fois near the fundamental frequency,and k sustain time, body resonance, string flexibility,
is the partial number. The parameter8 has been bridge and pick hardness, and degree of sympathetic
called the coefficient of inharmonicity; if 8 = 0, string excitation, which, in the case of real instru-
then perfect harmonicity results. ments, are usually fixed at the time of instrument
This effect can be created, in principle, by intro- construction, are available as performanceparame-
ducing an allpass filter Hg(z)in the string loop (Al- ters. It is important to point out that this variety is
len 1982) much as was done for the fine tuning of at no time achieved at the expense of the integrity
pitch. The phase delay in samples desired for an all- of the basic sound. Rather, as is the case with a real
pass filter inserted in the feedback loop of a harmonic musical instrument, the diversity exists within the
string simulator tuned to fl is given by solving bounds of a clearly defined sound domain.

68 ComputerMusic Journal

This content downloaded on Wed, 9 Jan 2013 [Link] PM

All use subject to JSTOR Terms and Conditions
 I

Conclusion damentals of Bowed String Dynamics." Acustica 43(2):

93-108.
The algorithm originated by Karplus and Strong McIntyre,M. E., R. T. Schumacher, and J. Woodhouse.
and extended by the methods outlined here has 1981. "Aperiodicityin Bowed-StringMotion." Acustica
49(1): 13-32.
proven very useful as a computer instrument. In
the process of composing May All Your Children McIntyre,M. E., R. T. Schumacher, and J. [Link]
press. "On the Oscillations of Musical Instruments."
Be Acrobats and Silicon Valley Breakdown, it was
Journalof the Acoustical Society of America.
found to be sufficiently flexible to allow for a wide
Mont-Reynaud,B. 1982. Private communication.
range of musical expression and sufficiently idio- Morse, P. M. 1976. Vibration and Sound. New York:
syncratic to maintain a characteristic identity. We American Institute of Physics for the Acoustical So-
expect that new refinements of the algorithm will ciety of America. (Originally published in two editions
continue to arise. [1936 and 1948].)
Smith, J. 0. 1982a. "Synthesis of Bowed Strings."Paper
presented at the Acoustical Society of America Con-
References ference, Chicago, Illinois. (Reprintsavailable upon
request.)
Smith, J. 0. 1982b. "Introductionto Digital Filters."
Abramowitz,M., and I. A. Stegun, eds. 1966. Handbook
Typescript.(Copies available upon request.)
of Mathematical Functions. Washington, D.C.: Na-
tional Bureauof Standards. Smith, J. 0. 1983. "Techniquesfor Digital Filter Design
and System Identification with Application to the Vio-
Allen, J. B. 1982. Private communication. lin." Ph.D. Diss., Electrical EngineeringDepartment,
Karplus,K., and A. Strong. 1983. "Digital Synthesis of StanfordUniversity.
Plucked-Stringand Drum Timbres." ComputerMusic
Steiglitz, K. 1974. An Introduction to Discrete Systems.
Journal7(2):43- 55. New York:Wiley.
McIntyre,M. E., and J. Woodhouse. 1960. "On the Fun-

Jaffe and Smith 69

This content downloaded on Wed, 9 Jan 2013 [Link] PM

All use subject to JSTOR Terms and Conditions