ConsonanceScoring.nb

(* Calculate good rational representation of [value] within relative \
error [epsilon], with denominator not bigger than [limit].

Returns List[(approximation, relativeError, absoluteError)]
 *)
RationalReprInside := Function[{limit, epsilon, value}, Module[
   {reprs, real, approx, delta},
   reprs = {};
   (* For each denominator, 
   find the appropriate approximation if possible *)
   For[l = 1, l <= limit, l++,
    real = value*l;
    approx = Round[real];
    delta = Abs[approx - real];
    If[delta < epsilon \[And] GCD[approx, l] == 1,
     reprs = Append[reprs, {approx/l, N[delta, 5], N[delta/l, 5]}]
     ];
    ];
   
   reprs
   ]]


(* Calculate a consonant score of the given [value] *)
ConsonanceScore := Function[{limit, epsilon, value}, Module[
   {rationalRepr, scores},
   (* First calculate the rational approximations *)
   rationalRepr = RationalReprInside[limit, 0.5, value];
   scores = Map[
     (
       (* Each rational representation contribute a score *)
       (* The smaller the denominator and numerator, 
       the more consonant it is *)
       (1/(Denominator[#[[1]]]*Numerator[#[[1]]])^(1/2))
        (* 
        The smaller the disance between it and the rational number, 
        the more consonant it is *)
        *((PDF[NormalDistribution[0, epsilon/3], #] &)[ #[[2]]])
        *10) &, rationalRepr];
   (* Then add up the representations *)
   N[Sum[i, {i, scores}]]
   ]]