Functions for just intonation, consonance, roughness, dissonance, etc. in Python

consonance.py

"""
Collection of functions related to musical consonance, Just intonation, etc.
"""

from fractions import Fraction
from math import log
from functools import reduce


def product(*iterable):
    p = 1
    for n in iterable:
        p *= n
    return p


def gcd(*numbers):
    """
    Return the greatest common divisor of the given integers

    The result will have the same sign as the last number given (so that
    when the last number is divided by the result, the result comes out
    positive).
    """
    def gcd(a, b):
        while b:
            a, b = b, a % b
        return a

    return reduce(gcd, numbers)


def lcm(*numbers):
    """
    Return the least common multiple of the given integers
    """
    def lcm(a, b):
        if a == b == 0:
            return 0
        return (a * b) // gcd(a, b)

    # LCM(a, b, c, d) = LCM(a, LCM(b, LCM(c, d)))
    return reduce(lcm, numbers)


def reduced_form(*numbers):
    """
    Return a tuple of numbers which is the reduced form of the input,
    which is a list of integers
    """
    return tuple(int(a // gcd(*numbers)) for a in numbers)


def prime_factors(n):
    """
    Return a list of the prime factors of the integer n.

    Don't use this for big numbers; it's a dumb brute-force method.
    """
    factors = []
    lastresult = n

    while lastresult > 1:
        c = 2
        while lastresult % c > 0:
            c += 1
        factors.append(c)
        lastresult /= c

    return factors


def euler(*numbers):
    """
    Euler's "gradus suavitatis" (degree of sweetness) function

    Return the "degree of sweetness" of a musical interval or chord expressed
    as a ratio of frequencies a:b:c, according to Euler's formula

    Greater values indicate more dissonance
    """
    factors = prime_factors(lcm(*reduced_form(*numbers)))
    return 1 + sum(p - 1 for p in factors)


def wiseman(a, b):
    """
    Return a value corresponding to the dissonance of a musical interval
    of ratio a:b

    This is taken from Gus Wiseman's paper at
    http://www.nafindix.com/math/sensory.pdf

    It apparently is derived from the total period of two tones, so a tone
    of period 2 and a tone of period 3 when summed will give a total period
    of 6 (which is LCM(2,3))

    Since a ratio of frequencies is just the inverse of a ratio of periods,
    and inverting the terms has no effect on this formula, I *think* it's
    valid to just plug in frequency ratios, too.

    Same as "harmonic value"?

    TODO: a/gcd(all) and b/gcd(all) is just the reduced form, right?
    and if a,b are in reduced form, then lcm(a,b) =a*b, right?  So this is
    just a*b when a,b are in reduced form.  So this should be benedetti().
    https://en.wikipedia.org/wiki/Giambattista_Benedetti#Music

    product(*reduced_form(a,b))

    http://xenharmonic.wikispaces.com/Benedetti+height
    """
    return lcm(a / gcd(a, b), b / gcd(a, b))


def tenney(a, b):
    """
    Return Tenney's harmonic distance of a musical interval of ratio a:b

    "Tenney's HD function, in one of its simplest forms, is perhaps the most
    direct. It is HD(a/b) = log(ab)"

    Also called "Tenney height" with log2:
        http://xenharmonic.wikispaces.com/Tenney+Height
    Tenney's paper says "harmonic distance" and "log", but relates it to a
    lattice so maybe the log2 is implied?

    http://lumma.org/tuning/faq/#heights says
    Of these, the product of numerator and denominator, often
    written n*d, for a ratio n/d in lowest terms, has been found to
    best agree with informal rankings by listeners as well as
    published psychoacoustic data.

    The n*d rule is also called Tenney height, after James Tenney^1,

    Tenney height generalizes to geomean(a*b*c...) for a chord
    a:b:c..., which gives sqrt(n*d) for dyads.
    """
    return log(product(*reduced_form(a, b)))

See also Just intonation classes for music theory experiments in Python
(moved to a proper repo: https://github.com/endolith/just_intonation/)