tuneralgorithm.cpp

Go to the documentation of this file.

/*****************************************************************************
 * Copyright 2016-2017 Karolin Stange, Christoph Wick, and Haye Hinrichsen
 *
 * This file is part of JustIntonation.
 *
 * JustIntonation is free software: you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by the
 * Free Software Foundation, either version 3 of the License, or (at your
 * option) any later version.
 *
 * JustIntonation is distributed in the hope that it will be useful, but
 * WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
 * or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License 
 * for more details.
 *
 * You should have received a copy of the GNU General Public License along 
 * with JustIntonation. If not, see http://www.gnu.org/licenses/.
 *****************************************************************************/

//=============================================================================
//                          Main Tuning Algorithm
//=============================================================================

#include "tuneralgorithm.h"
#include "tuner.h"
#include "tunerdebug.h"

#include <QGenericMatrix>

//-----------------------------------------------------------------------------
//                              Constructor
//-----------------------------------------------------------------------------

TunerAlgorithm::TunerAlgorithm ()
    : mSelectedVariant(128,128)         // selected JI-variant of in interval
    , mProgression(0)                   // pitch progression variable
{
    // Reset just intonation interval sizes
    mSelectedVariant.setZero();
}


//-----------------------------------------------------------------------------
//              Tune in a statically defined unequal temperament
//-----------------------------------------------------------------------------

void TunerAlgorithm::tuneStatically (KeyDataVector &keys,
                                     const QVector<double> pitches,
                                     const int referenceKey,
                                     const int wolfsIntervalShift)
{
    for (KeyData &key : keys) if (key.intensity > 0)
    {
        int distance = (key.key + 240 - referenceKey + wolfsIntervalShift)%12;
        int A4dist   = (69 + 240 - referenceKey + wolfsIntervalShift)%12;
        key.pitch = pitches[distance]-pitches[A4dist];
    }
    mProgression = 0;
}


//-----------------------------------------------------------------------------
//                  Main Tuning Algorithm: Tune dynamically
//-----------------------------------------------------------------------------

double TunerAlgorithm::tuneDynamically (KeyDataVector &keyDataVector,
                                        const QVector<double> intervals,
                                        const QVector<double> weights,
                                        bool optimizedJI)
{
    using namespace Eigen;

    // optimizedJI=false; // turn off optimazation by force

    // First we set up maps <intensity,keynumber> of the
    // playing keys (the ones which are audible) and the memorized
    // key (the playing ones AND those that were audible short time ago).

    QMultiMap<double,int> audibleKeys, memorizedKeys;
    for (auto &key : keyDataVector)
    {
        if (key.intensity>0) audibleKeys.insert  (-key.intensity, key.key);
        if (key.memory>0)    memorizedKeys.insert(-key.memory, key.key);
    }

    // The entries in the maps are automatically sorted by intensity.
    // In order to limit the complexity of the tuning process, we
    // take at most 8 keys and copy them to the vectors
    // playingKeys and memorizedKeys:

    const int Pmax=8, Nmax=16;
    QVector<int> keys = audibleKeys.values().mid(0,Pmax).toVector();
    const int P = keys.size();
    keys += memorizedKeys.values().mid(0,Nmax-P).toVector();
    const int N = keys.size();

    if (P==0) return 0;

    // Copy pitch, define the significance as sum of the actual volume
    // and the level of memorization. Find out whether we have newly
    // pressed keys.
    VectorXd pitch(N), significance(N);
    bool newKeys = false;
    for (int i=0; i<N; ++i)
    {
        KeyData& key = keyDataVector[keys[i]];
        pitch(i) = key.pitch;
        significance(i) = key.intensity + key.memory;
        if (key.newlyPressed) newKeys = true;
    }

    // Determine the properties of the intervals
    MatrixXi variants(N,N), semitones(N,N), direction(N,N);
    MatrixXd weight(N,N); weight.setZero();
    for (int i=0; i<N; ++i) for (int j=0; j<N; ++j) if (i<P or j<P)
    {
        semitones(i,j) = abs(keys[j] - keys[i]);
        direction(i,j) = (keys[j] >= keys[i] ? 1:-1);
        variants(i,j)  = cJustIntonation[semitones(i,j)%12].size();
        if (i!=j)
        {
            weight(i,j) = weights[semitones(i,j)%12] *
                          std::pow(octaveWeightFactor, (semitones(i,j)-1)/12) *
                          sqrt(significance(i) * significance(j));
            if (i>P or j>P) if (keys[i]!=keys[j]) weight(i,j) *= memoryWeight;
        }
    }

    // Construct the matrix A and its inverse AI (independent of optimization)
    VectorXd diagonal = VectorXd(weight.array().rowwise().sum());
    MatrixXd A = -weight.block(0,0,P,P) + MatrixXd::Identity(P,P)*epsilon
               + MatrixXd(diagonal.head(P).asDiagonal());
    MatrixXd AI = A.inverse();

    // Determine the interval sizes (optimized or without optimization)
    MatrixXd interval(N,N); interval.setZero();

    if (optimizedJI and newKeys)
    {
        // ####################################################################
        // Determine interval sizes for optimized JI according to the hardcoded
        // in cJustIntonation in the header file, iterating over all possible
        // combinations and searching for the lowest degree of tempering.

        MatrixXi m(N,N); m.setZero();
        for (int i=0; i<N; i++) for (int j=0; j<N; j++) if (i<P or j<P)
        {
            if (keyDataVector[keys[i]].newlyPressed or
                keyDataVector[keys[j]].newlyPressed) m(j,i)=m(i,j)=0;
            else m(i,j) = mSelectedVariant(keys[i],keys[j]);
        }
        bool searching = true;              // Flag for searching
        MatrixXi mmin = m;                  // Optimal indices
        double Pmin = 1e100;                // Min value of dissipated power
        QTimer timer; timer.start(20);      // Timeout 20ms

        while (searching and timer.remainingTime()>0)
        {
            MatrixXd optimalInterval(N,N); optimalInterval.setZero();
            for (int i=0; i<N; ++i) for (int j=0; j<N; ++j) if (i<P or j<P)
            {
                optimalInterval(i,j) = direction(i,j) *
                        cJustIntonation[semitones(i,j)%12][m(i,j)];
                if (j>=P) optimalInterval(i,j) -= pitch(j);
                if (i>=P) optimalInterval(i,j) += pitch(i);
            }
            VectorXd B = (optimalInterval * weight).diagonal().head(P)
                         - epsilon*pitch.head(P);
            VectorXd V = -AI * B;
            double C = (optimalInterval.array() * optimalInterval.array()
                        * weight.array()).sum() / 4;
            double P = V.dot(A*V)/2 + B.dot(V) + C;
            if (P<Pmin-1E-7) { Pmin=P; mmin=m; }

            // interate over all combinations:
            searching = false;
            for (int i=0; i<N; ++i) for (int j=0; j<N; ++j) if (i<P or j<P) if (i<j)
                if (keyDataVector[keys[i]].newlyPressed
                    or keyDataVector[keys[j]].newlyPressed)
                    if (variants(i,j)>1 and not searching)
            {
                m(j,i) = m(i,j) = (m(i,j)+1) % variants(i,j);
                if (m(i,j)>0) searching = true;
            }
        }

        // After finding optimum (or timeout) copy the optimal result
        for (int i=0; i<N; i++) for (int j=0; j<N; j++) if (i<P or j<P)
        {
            mSelectedVariant(keys[i],keys[j]) = mmin(i,j);
            interval(i,j) = direction(i,j) *
                            cJustIntonation[semitones(i,j)%12][mmin(i,j)];
        }
        // ####################################################################
    }

    else // if not optimized
    {
        // ####################################################################
        // Standard procedure: Non-optimized tuning according to a given table
        // of interval sizes passed in the second parameter named 'intervals'
        for (int i=0; i<N; ++i) for (int j=0; j<N; ++j) if (i<P or j<P)
            interval(i,j) = direction(i,j) * intervals[semitones(i,j)%12];
    }

    // Incorporate the -lambda terms on the right hand side (see paper)
    for (int i=0; i<P; ++i) for (int j=P; j<N; ++j)
    {
        interval(i,j) -= pitch(j);
        interval(j,i) += pitch(j);
    }

    // Tune the keys microtonally:
    VectorXd B = (interval * weight).diagonal().head(P) - epsilon*pitch.head(P);
    VectorXd U = - AI * B;
    for (int i=0; i<P; i++) keyDataVector[keys[i]].pitch = U(i);

    // Cancel newly-pressed-key flags
    for (int i=0; i<P; ++i) keyDataVector[keys[i]].newlyPressed = false;

    // Compute potential as a measure of how much the chord is tempered
    double C = (interval.array() * interval.array() * weight.array()).sum() / 4;
    return C - B.dot(AI*B)/2;
}