Pascal Source Code for the Kameoka and Kuriyagawa model for sensory dissonance
UNIT FKK;
{----------------------------------------------------------------------}
{ Source of unit that implements the Kameoka and Kuriyagawa model for }
{ sensory dissonance, according to their articles in The Journal of }
{ the Acoustical Society of America, volume 45, number 6, 1969, }
{ pages 1451 to 1469. Pascal text for Borland's Turbo Pascal compi- }
{ ler, version 4.0. This file implements the following: }
{ }
{ CONSTANTS : KK_number_max }
{ KK_version }
{ TYPES : KK_spectrum }
{ KK_parameters }
{ VARIABLES : KK_conditions }
{ PROCEDURES : KK_set_default_conditions }
{ KK_set_conditions }
{ KK_show_conditions }
{ KK_show_stimulus (stimulus,format) }
{ KK_show_info }
{ KK_checks }
{ KK_preprocess (raw, sorted) }
{ FUNCTIONS : KK_dissonance (stimulus) }
{ KK_fb (f1, L) }
{ LOG (a) }
{ POWER (a, b) }
{ ubar_to_dBspl (ubar) }
{ dBspl_to_ubar (dBspl) }
{ }
{ Before using these functions and procedures, first initialize the }
{ variable KK_conditions and call KK_set_conditions, which then up- }
{ dates the derived parameters (in the same global variable). Never }
{ change the derived parameters directly. You may also call the pro- }
{ cedure KK_set_default_conditions to initialize and update. Global }
{ variable KK_conditions influences almost all functions and proce- }
{ dures. So be sure you initialized and updated correctly. }
{ }
{ Pieter Suurmond, 5 march 1994. }
{----------------------------------------------------------------------}
INTERFACE
USES Dos, Crt;
CONST KK_version = '1.06 (13 march 1994)';
KK_number_max = 80; { Maximum number of frequency components. }
TYPE KK_spectrum = RECORD
number : BYTE; { Number partials. }
Hertz : ARRAY[0..KK_number_max] OF REAL; { Frequency. }
ubar : ARRAY[0..KK_number_max] OF REAL; { Sound pressure. }
rad : ARRAY[0..KK_number_max] OF REAL; { Phase. }
END;
KK_parameters = RECORD
{ Primary parameters. These may be set or adjusted by the user. }
{ Always call the procedure KK_set_conditions after changing one of }
{ these (except for Lfb, fl and fh, (the only 3!) which may be }
{ changed without updating the derived parameters). }
{ For the actual model: }
kk0 : REAL; { Noise conditions at auditorium. }
CC0 : REAL;
k0 : REAL;
C0 : REAL;
nh : REAL; { Exponents for sound pressure correction. }
nl : REAL;
ne : REAL;
p57 : REAL; { Sound pressure (microbar) at 57 dB SPL. }
beta : REAL; { Extended power law constant. }
Lm : REAL; { Masking level difference (dB SPL). }
Lfb : REAL; { Fixed level (dBspl) to be used by KK-dis- }
{ sonance at fb calculation. When Lfb = -1, }
{ the average dB level is used. }
{ For the preprocessor: }
match : REAL; { Deviation (cents) identical frequencies. }
fl : REAL; { Bandpass low frequency (Hertz). }
fh : REAL; { Bandpass high frequency (Hertz). }
slope : REAL; { Bandpass cutoff slopes (dBspl/octave). }
treshold : REAL; { Partials below this (dB SPL) are trimmed }
{ (equal freqs will be added first). }
{ Derived parameters. These are calculated from the preceeding ones }
{ by the procedure KK_set_conditions. Never change them directly. }
{ For the actual model: }
p0 : REAL; { Derived from p57. Speeds up dB <-> ubar. }
ne_nh : REAL; { Used at function D2i. }
ne_nl : REAL;
p57_ne : REAL;
beta_inv : REAL;
DDIn : REAL; { Dissonance intensity of noises. }
Dn0_beta_inv : REAL; { Dn0^(1/beta), used by function Dm. }
k0_kk0 : REAL; { Used by KK_dissonance. }
mask_r : REAL; { Pressure ratio perfect masking. }
mask_r_inv : REAL; { Inverse Pressure ratio perfect masking. }
{ For the preprocessor: }
dev : REAL; { Deviation (fraction) equal frequencies. }
expon : REAL; { Exponent derived from bandpass-slope. }
pt : REAL; { Derived from treshold. }
END;
VAR KK_conditions : KK_parameters; { Global var. to store conditions. }
PROCEDURE KK_set_default_conditions;
PROCEDURE KK_set_conditions;
PROCEDURE KK_show_conditions;
PROCEDURE KK_show_stimulus (VAR S : KK_spectrum; format : BYTE);
PROCEDURE KK_show_info;
PROCEDURE KK_checks;
PROCEDURE KK_preprocess (VAR Stim, Stim_sorted : KK_spectrum);
FUNCTION KK_dissonance (VAR stim : KK_spectrum) : REAL;
FUNCTION KK_fb (f1, L : REAL) : REAL;
FUNCTION ubar_to_dBspl (ubar : REAL) : REAL;
FUNCTION dBspl_to_ubar (dBspl : REAL) : REAL;
FUNCTION POWER (X: REAL; Y: REAL): REAL;
FUNCTION LOG (X : REAL) : REAL;
{----------------------------------------------------------------------}
IMPLEMENTATION
PROCEDURE Message (S : STRING);
FORWARD;
{----------------------------- BASIC ARITHMETIC: ----------------------}
FUNCTION POWER (X: REAL; Y: REAL): REAL;
{ Raises X to the Yth power. }
{ X may not be negative. If X is 0, then Y must be positive. }
BEGIN
IF X > 0 THEN
POWER := EXP (Y * LN (X))
ELSE
IF X < 0 THEN
Message ('Error: Cannot raise negative numbers.')
ELSE { Thus X=0. }
IF Y > 0
THEN POWER := 0
ELSE
Message ('Error: Cannot raise 0 to non-positive numbers.')
END;
FUNCTION LOG (X : REAL) : REAL;
{ Calculates the 10 log of X. X must be greater than 0. }
CONST LN10 = 2.302585093;
BEGIN
IF X <= 0 THEN
Message ('Error: Logarithm of a non-positive number.')
ELSE
LOG := LN (X) / LN10
END;
{------------------------ SOUND PRESSURE LEVELS: ----------------------}
FUNCTION ubar_to_dBspl (ubar : REAL): REAL;
BEGIN
ubar_to_dBspl := 20 * LOG (ubar/KK_conditions.p0)
END;
FUNCTION dBspl_to_ubar (dBspl : REAL): REAL;
BEGIN
dBspl_to_ubar := KK_conditions.p0 * POWER (10, dBspl/20)
END;
{------------------------ CALCULATION CONDITIONS: ---------------------}
PROCEDURE KK_set_conditions;
{ Calculates and sets derived parameters in global variable }
{ KK_conditions. Call this procedure after you changed one or }
{ more of the primary parameters in KK_conditions. }
BEGIN
WITH KK_conditions DO
BEGIN
{ For the model: }
p0 := p57 * POWER (10, -57/20); { For level conversions. }
ne_nh := ne - nh;
ne_nl := ne - nl;
p57_ne := POWER (p57, -ne);
beta_inv := 1 / beta;
DDIn := POWER (CC0, beta_inv);
Dn0_beta_inv := POWER (k0*C0, beta_inv);
k0_kk0 := k0 / kk0; { Used by KK_dissonance. }
mask_r := POWER (10, Lm/20); { Pressure ratios at which }
mask_r_inv := 1 / mask_r; { perfect masking occurs. }
{ For the preprocessor: }
pt := dBspl_to_ubar (treshold); { Convert to linear scale. }
dev := -1 + POWER (2, (match/1200)); { cents to fractions. }
expon := slope / (20 * LOG(2)); { dB/octave to frequency- }
END; { ratio-exponent. }
END;
PROCEDURE KK_set_default_conditions;
{ Sets the default calculation conditions of the model in the }
{ KK_conditions (it's thus destructive). You don't need to }
{ call KK_set_conditions, after you've called this procedure. }
BEGIN
WITH KK_conditions DO
BEGIN
kk0 := 1.0; { Noise conditions at auditorium of the }
CC0 := 65.0; { Toshiba Central Research Laboratory. }
k0 := 1.0;
C0 := 65.0;
p57 := 0.141589156; { (microbar) 2*10^-1.15 at 57 dBspl. }
nh := 0.15; { Exponents for sound pressure correc- }
nl := 0.32; { tions, according to K & K. }
ne := 0.20;
beta := 0.25; { Extended power law constant, K & K. }
Lm := 24.00; { (dBspl) Level difference at which }
{ perfect masking occurs. }
Lfb := 57.00; { Use fixed level at fb calculation. }
{ For the preprocessor: }
match := 0.00001; { (cents) Quite strict. }
fl := 20; { (Hertz) Auditive range from. }
fh := 20000;
slope := 60; { (dBspl/octave). }
treshold := 18; { (dBspl) No contribution below this. }
END;
KK_set_conditions; { Calculates derived parameters. }
END; { of KK_set_default_conditions }
{------------------ PREPROCESSOR FOR DISSONANCE MODEL: ----------------}
PROCEDURE KK_preprocess (VAR Stim, Stim_sorted : KK_spectrum);
{ This procedure reads the raw spectral data from Stim and }
{ outputs a sorted copy in Stim_sorted which then can be used }
{ by KK_dissonance. Besides that frequency components will be- }
{ come ordered on frequency (low to high) and some components }
{ may disappear. Dependant on the values in KK_conditions, }
{ "equal" frequencies will be (phase-dependantly) added and }
{ frequencies below fl or above fh will be attenuated. Lastly, }
{ all components below treshold are trimmed from the spectrum. }
VAR p, diff : REAL;
Point : ARRAY[1..KK_number_max] OF BYTE;
i, j : INTEGER;
PROCEDURE Sort_frequencies;
{ Fills the array Point with sorted BYTE-indices to fre- }
{ quencies (from low to high). }
VAR i, l, H : BYTE;
j, a : INTEGER;
BEGIN
FOR i := 1 TO Stim.number DO { Create array of pointers }
Point[i] := i - 1; { to frequency components. }
a := Stim.number DIV 2;
WHILE a > 0 DO
BEGIN
FOR i := a + 1 TO Stim.number DO
BEGIN
j := i - a;
WHILE j > 0 DO
BEGIN
l := j + a;
IF Stim.Hertz[Point[j]] <= Stim.Hertz[Point[l]] THEN
j := 0
ELSE
BEGIN
H := Point[j];
Point[j] := Point[l]; Point[l] := H
END;
j := j - a
END
END;
a := a DIV 2
END
END;
PROCEDURE Add_phase_vectors (VAR amp1, ph1, amp2, ph2 : REAL);
{ Adds second phase vector (with linear magnitude amp2 and }
{ radian phase ph2) to the first vector. The result is thus }
{ written to amp1 and ph1. Note that this is only valid when }
{ the two frequencies are exactly equal. }
VAR i : INTEGER;
Re, Im : REAL;
BEGIN
Re := amp1 * COS(ph1) + amp2 * COS(ph2); { Convert polar form to }
Im := amp1 * SIN(ph1) + amp2 * SIN(ph2); { cartesian and add up. }
amp1 := SQRT (Re*Re + Im*Im); { Resulting magnitude. }
IF Re = 0 THEN { Resulting phase. }
IF Im > 0 THEN
ph1 := PI/2
ELSE
ph1 := -PI/2 { If Re = Im = 0 phase doesn't matter. }
ELSE
BEGIN
ph1 := ARCTAN (Im/Re); { 1st or 4th quadrant. }
IF Re < 0 THEN
IF Im < 0 THEN
ph1 := ph1 - PI { 3rd quadrant. }
ELSE
ph1 := ph1 + PI { 2nd quadrant. }
END
END;
FUNCTION Bandpass (f : REAL) : REAL;
{ Attenuates frequency components below .fl and above .fh }
{ (Hertz), with a slope of .slope (dBspl per octave). To- }
{ gether with .treshold, this function is used for a gradual }
{ removal of components outside the audible range. }
BEGIN
WITH KK_conditions DO
IF (fl <= f) AND (f <= fh) THEN
bandpass := 1
ELSE
IF f > fh THEN
bandpass := POWER (fh/f, expon)
ELSE {thus f < fl}
bandpass := POWER (f/fl, expon)
END;
BEGIN { Begin of procedure KK_preprocess. }
Sort_frequencies;
i := 1; { Reading from pointer array, index 1 to Stim.number. }
j := 0; { Writing to Stim_sorted, index 0 to Stim_sorted.number-1. }
WITH Stim_sorted DO
BEGIN
REPEAT
Hertz[j] := Stim.Hertz[Point[i]]; { Copy frequency compo- }
ubar[j] := Stim.ubar[Point[i]]; { nent to Stim_sorted. }
rad[j] := Stim.rad[Point[i]];
diff := KK_conditions.dev * Hertz[j];
WHILE (i < Stim.number) AND
(ABS (Stim.Hertz[Point[i+1]] - Hertz[j]) < diff) DO
BEGIN { Add up "identical" frequencies. }
i := i + 1;
Add_phase_vectors (ubar[j], rad[j], Stim.ubar[Point[i]],
Stim.rad[Point[i]])
END;
p := ubar[j] * bandpass(Hertz[j]); { Attenuated pressure. }
IF p >= KK_conditions.pt THEN { Skip when too weak. }
BEGIN
ubar[j] := p; { Possible sound pressure adjustment. }
j := j + 1
END;
i := i + 1;
UNTIL i > Stim.number;
number := j;
END;
END;
{------------------ THE KAMEOKA AND KURIYAGAWA MODEL: -----------------}
FUNCTION KK_fb (f1, L: REAL) : REAL;
{ Calculates the most dissonant frequency difference (fb in Hz) }
{ from lowest frequency (f1 in Hz) and the level (L in dB SPL). }
{ See J.A.S.A.: page 1455, equation 4 / page 1462, equation 6. }
BEGIN
IF L <= 17 THEN
BEGIN
KK_fb := 0;
Message
('Error: Most dissonant frequency difference unknown below 17 dBspl.')
END
ELSE
KK_fb := 0.05675 * (L-17) * POWER (f1, 0.477)
END;
FUNCTION KK_dissonance (VAR stim : KK_spectrum) : REAL;
{ Calculates total absolute dissonance of frequency spectrum }
{ stored in stim. The function does NOT destruct stim. Frequen- }
{ cies in stim are assumed to be sorted and unique. So always: }
{ stim.Hertz[a-1] < stim.Hertz[a]. This function skips all dy- }
{ ads with over octave separation (they do not contribute to }
{ total dissonance (independant from Co, Ko, Co' and ko'). }
{ Always uses the calculation conditions written in the global }
{ variable KK_conditions. Furthermore, if, and only if, the }
{ calc-Lfb option is true, the total intensity level will be }
{ left in KK_conditions.Lfb. }
FUNCTION D2i (VAR f1, p1, f2, p2, fb : REAL): REAL;
{ Calculates the real absolute dissonance of a dyad without }
{ internal and external noises. Modified for frequency sepa- }
{ rations smaller than 1%: D'2ei is then (k'o * C'o = 65), so }
{ D2i always becomes 0 (like over octave sep.independant from }
{ k'o, C'o, ko or Co). f1 and f2 are the lower and higher }
{ frequencies (Hertz). p1 and p2 are the sound pressures }
{ (microbars) of respectively the lower and higher frequen- }
{ cies. fb is the most dissonant frequency difference (Hertz).}
{ Over octave separation never occurs here anymore. Perfect }
{ masking occurs when the level difference exceeds Lm }
{ (decibels). The pressure ratio is then greater than mask_r }
{ and D2i becomes 0. }
VAR denom, numer, DD2ei, fD2i, p1_p2 : REAL;
BEGIN
WITH KK_conditions DO
BEGIN { Masking or small separation. }
p1_p2 := p1 / p2; { ATTENTION: p2 must be > 0. }
IF (p1_p2 > mask_r) OR (p1_p2 < mask_r_inv) OR (f2/f1 < 1.01) THEN { SMALL SEPARATION }
D2i := 0
ELSE
BEGIN { Sound pressure correction. }
IF p1_p2 >= 1 THEN
fD2i := POWER (p1,ne_nh) * POWER (p2,nh) * p57_ne
ELSE
fD2i := POWER (p2,ne_nl) * POWER (p1,nl) * p57_ne;
denom := LOG (fb/f1);
numer := LOG ((f2-f1)/f1);
IF f2 - f1 <= fb THEN { Dynamic domain. }
DD2ei := kk0*(CC0 + 100*(2+numer)/(2+denom))
ELSE { Static domain. }
DD2ei := kk0*(CC0 + 10 + 90*numer/denom);
D2i := fD2i*kk0*POWER (POWER (DD2ei/kk0, beta_inv)-DDIn, beta);
END
END
END;
VAR i, j : INTEGER;
octave, accu, L, fb : REAL;
BEGIN { of function KK_dissonance }
WITH Stim, KK_conditions DO
BEGIN
IF number > KK_number_max THEN Message
('Error: Too many frequency components.');
accu := Dn0_beta_inv; { Background dissonance = (k0*c0)^(1/beta). }
IF Lfb <> -1 THEN
L := Lfb { Use specified fixed level at fb calc. }
ELSE
BEGIN
L := 0; { Calculate average dB level. }
FOR i := 0 TO number-1 DO
L := L + ubar_to_dBspl (ubar[i]);
L := L / number;
END;
FOR i := 0 TO number - 2 DO { All possible dyads within }
BEGIN { critical band (1 octave). }
IF Hertz[i] >= Hertz[i+1] THEN { Check frequency rank order. }
Message
('Error: Frequencies in wrong order, preprocessing recommended.');
octave := 2 * Hertz[i];
fb := KK_fb (Hertz[i], L);
j := i + 1;
WHILE (Hertz[j] <= octave) AND (j < Number) DO
BEGIN
accu := accu +
POWER (D2i (Hertz[i],ubar[i],Hertz[j],ubar[j], fb),
beta_inv); { Inverse power Law. }
j := j + 1
END
END;
KK_dissonance := k0_kk0 * POWER (accu, beta) { Power law. }
END
END; { of function KK_dissonance }
{-------------------------- BASIC USER I/O: ---------------------------}
PROCEDURE KK_show_conditions;
BEGIN
WITH KK_conditions DO
BEGIN
WRITELN ('CONDITIONS: nl =',nl:9:6,' ne =',ne:9:6,
' nh =', nh:9:6);
WRITELN (' p57 =', P57:7:4, ' ubar',' ko =', k0:9:6,
' ko',CHR($27),' =',kk0:9:6,' Lm =', Lm:9:4, ' dB');
WRITE (' ',CHR($E1),' =', beta:10:7,' Co =', C0:9:4,
' Co', CHR($27), ' =', CC0:9:4);
IF Lfb <> -1 THEN
WRITELN(' Lfb =', Lfb:9:4, ' dB')
ELSE
WRITELN(' Lfb = total average');
WRITELN (' fl =',fl:7:2,' Hz fh =',fh:6:0,' Hz ',
'tresh =',treshold:6:2,' dB match =', match:9:6, ' cents');
WRITELN ('slope =',slope:7:2,' dB/8');
WRITELN
END
END;
PROCEDURE KK_show_stimulus (VAR S : KK_spectrum; format : BYTE);
VAR i, V, R, A : INTEGER;
BEGIN
WITH S DO
BEGIN
V := 0; R := Number;
WHILE R>0 DO
BEGIN
IF R < 8 THEN
A := R ELSE A := 8;
IF (format AND 1)>0 THEN { format 1 = Print Hertz. }
BEGIN
WRITE ('Hertz:');
FOR i := V TO V+A-1 DO
WRITE (Hertz[i]:9:2)
END;
IF (format AND 2)>0 THEN { format 2 = Print dBspl. }
BEGIN
WRITELN; WRITE ('dBspl:');
FOR i := V TO V+A-1 DO
WRITE (ubar_to_dBspl (ubar[i]):9:2)
END;
IF (format AND 4)>0 THEN { format 4 = Print ubar. }
BEGIN
WRITELN; WRITE (' ubar:');
FOR i := V TO V+A-1 DO
WRITE (ubar[i]:9:5)
END;
IF (format AND 8)>0 THEN { format 8 = Print angle. }
BEGIN
WRITELN; WRITE ('angle:');
FOR i := V TO V+A-1 DO
WRITE (rad[i]:9:4)
END;
WRITELN; WRITELN;
R := R-A; V := V+A;
END;
END
END;
PROCEDURE Message {S: STRING};
BEGIN
IF S<>'' THEN
WRITELN ('> FKK ', S);
WRITE (' Hit any key to continu. ');
WHILE READKEY = '' DO
WRITE;
GOTOXY (1, WHEREY);
WRITE (' ');
GOTOXY (1, WHEREY);
END;
PROCEDURE KK_show_info;
CONST s = ' ';
BEGIN
WRITELN;
WRITELN(s,'旼컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴�');
WRITE(s,'� Using unit FKK.TPU Version ',KK_version);
GotoXY(69,WhereY); WRITELN('�');
WRITELN(s,'� Pascal implementation of the sensory dissonance model ac- �');
WRITELN(s,'� cording to Akio Kameoka and Mamoru Kuriyagawa (published �');
WRITELN(s,'� in Journal of the Acoustical Society of America, 1969). �');
WRITELN(s,'� Pieter Suurmond (1993-95). �');
WRITELN(s,'읕컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴컴�');
END;
{---------------------------- CHECKS: ---------------------------------}
PROCEDURE KK_checks;
{ This procedure provides a primitive user interface to the }
{ dissonance model in which you may perform several checks. }
{ It also contains some examples that show you how to use this }
{ Pascal-unit. }
VAR C: CHAR;
Stimulus, Sorted : KK_SPECTRUM; { Contains input for the model. }
PROCEDURE PRINT_MDFD;
{ Prints most dissonant frequency difference for several }
{ lower frequencies. }
VAR f1, L : REAL;
OCT : INTEGER;
BEGIN
f1 := 110; {Hz}
L := 57; {dB SPL}
WRITELN;
WRITELN ('Most dissonant frequency difference:');
WRITELN;
WRITELN ('f1 [Hz]: level [dB SPL]: fb [Hz]:');
WRITELN;
FOR OCT :=1 TO 7 DO
BEGIN
WRITELN (f1:6:2, ' L=',L:6:2,' --> ', KK_fb(f1,L):6:2 );
f1 := 2 * f1;
END
END;
PROCEDURE KK_Expt_1;
VAR i: INTEGER;
PROCEDURE TAKE_KK_Expt_1 (stim_numb: INTEGER);
{ Fills Stimulus with one of the stimuli (1 to 11) used at }
{ Expt.1 of Kameoka and Kuriyagawa (page 1464 in J.A.S.A.). }
CONST stim : ARRAY[1..11, 1..2, 0..5] OF REAL =
{1} (((440, 0, 0, 0, 0, 0), { freq. in Hertz (0=end). }
(60, 0, 0, 0, 0, 0)), { levels in dB SPL. }
{2} ((440, 455, 484, 581, 0, 0),
(57.5, 49.5, 44, 45, 0, 0)),
{3} ((440, 455, 484, 880, 0, 0),
(56.5, 54, 53.5, 44.5, 0, 0)),
{4} ((440, 484, 0, 0, 0, 0),
(57, 57.5, 0, 0, 0, 0)),
{5} ((440, 455, 484, 581, 880, 0),
(55, 57, 54, 54.5, 46, 0)),
{6} ((440, 455, 484, 581, 880, 0),
(50, 54, 58.5, 54.5, 50, 0)),
{7} ((440, 455, 484, 581, 880, 0),
(46, 50.5, 55.5, 57, 54.5, 0)),
{8} ((440, 455, 484, 581, 880, 0),
(42, 46, 50, 55, 58, 0)),
{9} ((440, 484, 581, 880, 0, 0),
(51, 54, 50.5, 58, 0, 0)),
{10} ((440, 455, 484, 581, 880, 0),
(54, 43.5, 58, 50.5, 58.5, 0)),
{11} ((440, 455, 484, 880, 0, 0),
(57, 50.5, 57.5, 50.5, 0, 0)));
BEGIN
WITH Stimulus DO
BEGIN
number := 0; { Clear Stimulus. }
WHILE stim[stim_numb, 1, number] <> 0 DO
BEGIN
Hertz[number] := stim[stim_numb, 1, number];
ubar[number] := dBspl_to_ubar (stim[stim_numb,2,number]);
rad[number] := 0;
number := number + 1
END
END
END;
BEGIN
WRITELN;
WRITELN ('Kameoka and Kuriyagawa experiment 1. ',
'See J.A.S.A. page 1464.');
WRITELN;
KK_set_default_conditions;
KK_show_conditions;
FOR i := 1 TO 11 DO
BEGIN
WRITE (' Stimulus ', i, ': ');
TAKE_KK_Expt_1 (i);
WRITELN ('AD =', KK_dissonance(Stimulus):7:2);
KK_show_stimulus (Stimulus, 3); { 3 means Hertz and dBspl. }
Message('');
END
END;
PROCEDURE KK_Expt_2;
VAR i: INTEGER;
PROCEDURE TAKE_KK_Expt_2 (stim_numb: INTEGER);
{ Fills Stimulus with one of the stimuli (1 to 14) used at }
{ Expt.2 of Kameoka and Kuriyagawa (page 1464 in JASA). }
CONST stim : ARRAY[1..14, 1..2, 0..12] OF REAL =
{1} (((440, 461, 484, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0),
(56, 49.5, 42, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0)),
{2} ((440, 461, 484, 573, 880, 0, 0, 0, 0, 0,
0, 0, 0),
(55.5, 50.5, 48.5, 44.5, 42, 0, 0, 0, 0, 0,
0, 0, 0)),
{3} ((440, 459, 484, 572, 880, 0, 0, 0, 0, 0,
0, 0, 0),
(56, 52, 49.5, 44.5, 42.5, 0, 0, 0, 0, 0,
0, 0, 0)),
{4} ((440, 460, 482, 570, 880, 0, 0, 0, 0, 0,
0, 0, 0),
(55, 56, 52, 49.5, 46, 0, 0, 0, 0, 0,
0, 0, 0)),
{5} ((440, 484, 880, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0),
(57.5, 56.5, 50, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0)),
{6} ((440, 482, 880, 975, 1300, 1402, 1860, 0, 0, 0,
0, 0, 0),
(56, 56.5, 49.5, 55, 54.5, 49.5, 42, 0, 0, 0,
0, 0, 0)),
{7} ((440, 497, 880, 934, 1306, 1450, 1716, 1890, 2135, 2372,
0, 0, 0),
(55, 55, 52.5, 49.5, 52.5, 51.5, 48.5, 44.5, 47, 42.5,
0, 0, 0)),
{8} ((440, 497, 880, 924, 1305, 1452, 1718, 1886, 2131, 2347,
2580, 2880, 0),
(56, 56, 54.5, 54.5, 52, 51.5, 49.5, 49, 49, 49.5,
45.5, 44, 0)),
{9} ((440, 496, 880, 912, 1305, 1334, 1698, 1803, 2108, 2180,
2543, 2642, 0),
(56.5, 56, 49, 49.5, 51, 51, 50, 50, 53, 53.5,
44, 44, 0)),
{10} ((440, 496, 880, 912, 1296, 1337, 1698, 1822, 2100, 2200,
2555, 2642, 0),
(55.5, 55, 49.5, 50, 50.5, 50.5, 49.5, 50, 51.5, 51.5,
42.5, 43, 0)),
{11} ((440, 704, 880, 1275, 1406, 1712, 2031, 2133, 2580, 2706,
0, 0, 0),
(56, 56, 55, 52, 55, 50, 52, 46, 44, 50,
0, 0, 0)),
{12} ((440, 704, 880, 1275, 1406, 1712, 2031, 2133, 2580, 2706,
3355, 3920, 0),
(56, 56, 54.5, 52, 53.5, 49.5, 52, 46, 45, 50.5,
47.5, 43, 0)),
{13} ((440, 880, 1321, 1760, 2153, 2620, 0, 0, 0, 0,
0, 0, 0),
(55.5, 55.5, 55, 55, 55.5, 55, 0, 0, 0, 0,
0, 0, 0)),
{14} ((440, 880, 1321, 1760, 2153, 2620, 0, 0, 0, 0,
0, 0, 0),
(57.5, 57.5, 55.5, 51.5, 47, 45, 0, 0, 0, 0,
0, 0, 0)));
VAR i : INTEGER;
BEGIN
WITH Stimulus DO
BEGIN
number := 0;
WHILE stim[stim_numb, 1, number] <> 0 DO
BEGIN
Hertz[number] := stim[stim_numb, 1, number];
ubar[number] := dBspl_to_ubar (stim[stim_numb,2,number]);
rad[number] := 0;
number := number + 1;
END
END
END;
BEGIN
WRITELN;
WRITELN ('Kameoka and Kuriyagawa experiment 2. ',
'See J.A.S.A. page 1464.');
WRITELN;
KK_set_default_conditions;
KK_show_conditions;
FOR i := 1 TO 14 DO
BEGIN
WRITE (' Stimulus no. ', i, ': ');
TAKE_KK_Expt_2 (i);
WRITELN ('AD =', KK_dissonance(Stimulus):7:2);
KK_show_stimulus(Stimulus, 3); { 3 means Hertz and dBspl. }
Message('');
END
END;
PROCEDURE KK_Expt_3;
VAR r : REAL;
PROCEDURE TAKE_KK_Expt_3 (ratio: REAL);
{ Fills Stimulus with two complex harmonic tones. Fundamental }
{ frequency of the first is fixed at 440 Hz, the second de- }
{ pends on ratio. Spectral structures accord to experiment 3 }
{ of Kameoka and Kuriyagawa. All phases are assumed to be 0. }
CONST levels : ARRAY[1..8] OF REAL = { Levels of the eight }
(-12, 0, 0, -4, -8, -12, -14, -16); { harmonics (0=57 dB SPL). }
VAR i, j : INTEGER;
BEGIN
WITH Stimulus DO
BEGIN
j := 0;
number := 16;
FOR i := 1 TO 8 DO
BEGIN
Hertz[j] := i * 440; { Fixed tone. }
ubar[j] := dBspl_to_ubar (57 + levels[i]);
rad[j] := 0;
j := j + 1;
Hertz[j] := Hertz[j-1] * ratio; { Variable tone. }
ubar[j] := ubar[j-1];
rad[j] := 0;
j := j + 1
END
END
END;
BEGIN
WRITELN (' K & K expt.3. See fig.7 on page 1465 in J.A.S.A..');
KK_set_default_conditions; { match := 0.01; Strict. }
WRITE (' Input parameter match (in cents) ');
READLN (KK_conditions.match);
KK_set_conditions;
KK_show_conditions;
WRITELN;
WRITE ('> Input frequency ratio ');
READLN (r);
TAKE_KK_Expt_3 (r);
KK_preprocess (Stimulus, Sorted);
WRITELN (' AD =', KK_dissonance (Sorted):7:2);
KK_show_stimulus (Sorted, 11); { 11 means Hertz + dBspl + angle. }
END;
PROCEDURE KK_Fig_9;
VAR r : REAL;
sp : INTEGER;
PROCEDURE TAKE_KK_Fig_9 (spec : INTEGER; ratio : REAL);
{ Fills Stimulus with two complex (harmonic) tones. With spec }
{ you may choose the spectrum (7, 8 or 9) as given at fig. 9 }
{ on page 1466 in J.A.S.A.. All phases are assumed to be 0. }
CONST freqs : ARRAY[7..9, 1..4] OF REAL = { Frequencies of }
((440, 880, 1760, 2640), { harmonics. }
(440, 1320, 2200, 3080),
(440, 880, 1760, 3520));
VAR i, j : INTEGER;
BEGIN
WITH Stimulus, KK_conditions DO
BEGIN
j := 0;
number := 8;
FOR i := 1 TO 4 DO
BEGIN
Hertz[j] := freqs[spec,i]; { Fixed tone. }
ubar[j] := p57;
rad[j] := 0;
j := j + 1;
Hertz[j] := ratio*freqs[spec,i]; { Variable tone. }
ubar[j] := p57;
rad[j] := 0;
j := j + 1;
END
END;
END;
BEGIN
WRITELN (' K & K fig. 9 in J.A.S.A..');
KK_set_default_conditions; { match := 0.01; Strict. }
WRITE (' Input parameter frequency match (in cents) ');
READLN (KK_conditions.match);
KK_set_conditions;
KK_show_conditions;
WRITELN;
WRITE ('> Input frequency ratio ');
READLN (r);
FOR sp := 7 TO 9 DO
BEGIN
TAKE_KK_Fig_9 (sp, r);
KK_preprocess (Stimulus, Sorted);
WRITE ('Spectrum ', sp);
WRITELN (' AD =', KK_dissonance (Sorted):7:2);
KK_show_stimulus (Sorted, 11); { 11 = Hertz + dBspl + angle. }
Message ('')
END
END;
CONST s = '--------------------------------------';
BEGIN { KK_checks. }
KK_set_default_conditions;
C := '0';
WRITELN (' ',s,s);
REPEAT
WRITELN ('MENU (CHECKING THE FKK UNIT):');
WRITELN ('M = Most dissonant freq. diff. ',
'1 = Kameoka and Kuriyagawa expt 1.');
WRITELN ('2 = Kameoka and Kuriyagawa expt 2. ',
'3 = Kameoka and Kuriyagawa expt 3.');
WRITELN ('9 = Kameoka and Kuriyagawa fig 9. ',
'Q = Quit FKK-checks.');
REPEAT
C := READKEY;
IF (CHR($60) < C) AND (C < CHR($7B)) THEN { Up case. }
C := CHR(ORD(C) - $20);
UNTIL C IN ['M','1','2','3','9','Q'];
CASE C OF
'M': PRINT_MDFD;
'1': KK_Expt_1;
'2': KK_Expt_2;
'3': KK_Expt_3;
'9': KK_Fig_9
END;{Case}
WRITELN;
UNTIL C = 'Q';
WRITELN (' ',s,s)
END;
{--------------- INITIALIZATION SECTION OF THE FKK UNIT: --------------}
BEGIN
KK_show_info; { Only once. }
KK_set_default_conditions;
END.
{----------------------------------------------------------------------}