/*
 * C code to generate the ultimate Protracker tuning table reference.
 *
 * Created because of a discussion about the history and imperfections of the
 * original period tables used in the Protracker replay routines.
 *
 * This calculates a reference table for all 16 different 'tunings' using 1024-
 * bits of precision. Just because it's possible, not because it's needed. But
 * useful to suppress discussion about rounding errors.
 *
 * pipe^nature, 2013
 *
 *
 *
 * BUILDING (linux):
 *
 * Download and install the Multiple Precision Floating-Point Reliable Library,
 * GNU MPFR.
 *
 * Compile with: gcc -std=c99 -o pttables pttables.c -lmpfr
 * Run with: ./pttables
 *
 * Watch the glorious numbers.
 *
 * HISTORY
 *
 * V1.1 - Changed the signs on the detune columns.
 *        Fixed the compilation instructions.
 *        Better header before each table.
 * V1.0 - Initial release.
 *
 */

#include <stdlib.h>
#include <stdio.h>
#include <stdarg.h>
#include <mpfr.h>

#define PRECISION 1024

#define CC_NTSC 3579545
#define CC_PAL  3546895



static int mpfr_cents( mpfr_t rop, mpfr_t op1, mpfr_t op2, mpfr_rnd_t rnd );
static const int period_table[16][36];
static const char * const note_name[12];



int main ( void )
{
    /*
     * Declare and initialize every GMP variable that we need.
     */
    mpfr_t sr_ideal, per_ideal, sr_rounded, sr_pt, diff_rounded, diff_pt;
    mpfr_t diff_rounded_pt;
    mpfr_inits2( PRECISION, sr_ideal, per_ideal, sr_rounded, (mpfr_ptr)NULL );
    mpfr_inits2( PRECISION, sr_pt, diff_rounded, diff_pt, (mpfr_ptr)NULL );
    mpfr_inits2( PRECISION, diff_rounded_pt, (mpfr_ptr)NULL );

    /*
     * There are 16 different tuning steps, where 0 is the nominal.
     * A tuning of -8 is one full semitone lower. A tuning step of '1' equals
     * 1/8 of a semitone, or 12.5 cents.
     */
    for ( int tuning = -8; tuning < 8; tuning++ )
    {

        /*
         * Print a header before each tuning section.
         *
         * Sorry guys, breaking my 79 column limit here.
         */
        mpfr_printf( "# Semitone            NTSC PERIOD                        SAMPLE RATE                          DETUNE (IN CENTS)        \n" );
        mpfr_printf( "#   Tuning                                                                            Rounded    Protracker  Protracker\n" );
        mpfr_printf( "#       Note     Ideal  Rounded PT Diff      Ideal        Rounded      Protracker    vs. ideal   vs. ideal   vs. round.\n" );
        mpfr_printf( "#--  --  ---   ----------  ---  ---  --   ------------  ------------  ------------   ----------  ----------  ----------\n" );

        /*
         * For each tuning step, there are 36 semitones, which fills 3 octaves.
         * Semitone 0 is C-0, semitone 35 is B-3.
         */
        for ( int semitone = 0; semitone < 36; semitone++ )
        {

            /*
             * For each semitone and tuning combination we calculate the ideal
             * sample rate. The ideal sample rate is taken from the AMIGA®
             * Hardware Reference Manual, Third Edition, page 158, Table 5-7.
             *
             * It's derived from a wish to replay a 16 byte long sample at a
             * pitch of 880 Hz. This gives a sample rate of 880*16 = 14080 Hz.
             *
             * 14080 Hz is the ideal sample rate for an unspecified semitone
             * 'A', and converting this unspecified semitone to an offset in
             * the original Protracker tables, we see that this semitone is
             * represented with the period at index 21.
             *
             * From this we get that the formula to calculate the ideal
             * frequency from a semitone number and tuning is;
             *
             *    ideal = 880 * 16 * 2^( ( (semitone-21) + tuning/8 ) / 12 )
             *
             *    ideal = 14080 * 2 ^ ( (semitone*8+tuning-168) / 96 )
             */
            mpfr_set_si( sr_ideal, semitone*8 + tuning - 21*8, MPFR_RNDN );
            mpfr_div_ui( sr_ideal, sr_ideal, 12*8, MPFR_RNDN );
            mpfr_ui_pow( sr_ideal, 2, sr_ideal, MPFR_RNDN );
            mpfr_mul_ui( sr_ideal, sr_ideal, 880*16, MPFR_RNDN );

            /*
             * From the ideal frequency we can calculate the ideal period.
             * The original tracker usees the NTSC frequency for this, which
             * is a little odd.
             *
             * The period value is how many clock cycles Paula will hold each
             * sample until changing the output. It is calculated from the
             * clock constant (also known as "CC") and the wanted sample rate.
             *
             * CC is 3579545 in the NTSC Amiga, and 3546895 in the PAL Amiga.
             *
             *    period = CC / sample rate
             */
            mpfr_ui_div( per_ideal, CC_NTSC, sr_ideal, MPFR_RNDN );

            /*
             * The ideal period value is now rounded to the nearest integer.
             *
             * This value is what should go into the tables.
             */
            long per_rounded = mpfr_get_si( per_ideal, MPFR_RNDN );

            /*
             * For our amusement, we can now calculate a few more values.
             *
             * Let's start with the actual sample rate after rounding the
             * period value. This is calculated with this formula:
             *
             *    sample rate = CC / period value
             */
            mpfr_set_ui( sr_rounded, CC_NTSC, MPFR_RNDN );
            mpfr_div_si( sr_rounded, sr_rounded, per_rounded, MPFR_RNDN );

            /*
             * Let's compare our calculated period value with the one used in
             * Protracker.
             */
            long per_pt = period_table[ tuning+8 ][ semitone ];
            long per_diff = per_pt - per_rounded;

            /*
             * Let's calculate the sample rate from the Protracker period
             * value.
             */
            mpfr_set_ui( sr_pt, CC_NTSC, MPFR_RNDN );
            mpfr_div_si( sr_pt, sr_pt, per_pt, MPFR_RNDN );

            /*
             * Finally we can calculate the difference between the ideal and
             * rounded sample rate, and between the ideal and the Protracker
             * sample rate. We describe this difference in cents, and the just-
             * noticeable difference for human hearing is around 10 cents.
             */
            mpfr_cents( diff_rounded, sr_ideal, sr_rounded, MPFR_RNDN );
            mpfr_cents( diff_pt, sr_ideal, sr_pt, MPFR_RNDN );
            mpfr_cents( diff_rounded_pt, sr_rounded, sr_pt, MPFR_RNDN );

            /*
             * Present this information as one huge line, for the user to cut
             * up as he or she wishes.
             */
            mpfr_printf(
                    " %2d  "     // Semitone (0 to 35)
                    "%+2d  "     // Tuning (-8 to 7)
                    "%2s%1d  "   // Note name
                    " "
                    "%10.6Rf  "  // Ideal period value
                    "%3ld  "     // Rounded period value
                    "%3ld  "     // Protracker period value
                    "% 2ld  "    // Period difference
                    " "
                    "%12.6Rf  "  // Ideal sample rate
                    "%12.6Rf  "  // Rounded sample rate
                    "%12.6Rf  "  // Protracker sample rate
                    " "
                    "%+10.6Rf  " // Rounded sample rate difference in cents
                    "%+10.6Rf  " // Protracker sample rate difference in cents
                    "%+10.6Rf"   // Protracker vs. Rounded sample rate diff
                    "\n",
                    semitone,
                    tuning,
                    note_name[ semitone % 12 ], 1 + (int)( semitone / 12 ),
                    per_ideal,
                    per_rounded,
                    per_pt,
                    per_diff,
                    sr_ideal,
                    sr_rounded,
                    sr_pt,
                    diff_rounded,
                    diff_pt,
                    diff_rounded_pt
                    );

        }

        /*
         * Print two empty lines after each tuning section to make it easier
         * to parse the output.
         */
        mpfr_printf( "\n\n" );

    }


    /*
     * Clean up after MPFR, which is not strictly necessary.
     */
    mpfr_clears( sr_ideal, per_ideal, sr_rounded, (mpfr_ptr)NULL );
    mpfr_clears( sr_pt, diff_rounded, diff_pt, (mpfr_ptr)NULL );
    mpfr_clears( diff_rounded_pt, (mpfr_ptr)NULL );

    return 0;
}



/*
 * Will calculate the difference between two frequencies, and return the value
 * in cents.
 */
static int mpfr_cents( mpfr_t rop, mpfr_t op1, mpfr_t op2, mpfr_rnd_t rnd )
{
    mpfr_div( rop, op2, op1, rnd );
    mpfr_log2( rop, rop, rnd );
    return mpfr_mul_ui( rop, rop, 12*100, rnd );
}



/*
 * This is just for the output routine, so that one can more easily find the
 * data for a given note.
 */
static const char * const note_name[12] =
{
    "C-", "C#", "D-", "D#", "E-", "F-", "F#", "G-", "G#", "A-", "A#", "B-"
};



/*
 * The original Protracker period values, as taken from the PT2.3a replay
 * source code at http://16-bits.org/pt_src/replayer/PT2.3a_replay_vbl.s
 */
static const int period_table[16][36] =
{
    {   907,856,808,762,720,678,640,604,570,538,508,480,
        453,428,404,381,360,339,320,302,285,269,254,240,
        226,214,202,190,180,170,160,151,143,135,127,120 },
    {   900,850,802,757,715,675,636,601,567,535,505,477,
        450,425,401,379,357,337,318,300,284,268,253,238,
        225,212,200,189,179,169,159,150,142,134,126,119 },
    {   894,844,796,752,709,670,632,597,563,532,502,474,
        447,422,398,376,355,335,316,298,282,266,251,237,
        223,211,199,188,177,167,158,149,141,133,125,118 },
    {   887,838,791,746,704,665,628,592,559,528,498,470,
        444,419,395,373,352,332,314,296,280,264,249,235,
        222,209,198,187,176,166,157,148,140,132,125,118 },
    {   881,832,785,741,699,660,623,588,555,524,494,467,
        441,416,392,370,350,330,312,294,278,262,247,233,
        220,208,196,185,175,165,156,147,139,131,123,117 },
    {   875,826,779,736,694,655,619,584,551,520,491,463,
        437,413,390,368,347,328,309,292,276,260,245,232,
        219,206,195,184,174,164,155,146,138,130,123,116 },
    {   868,820,774,730,689,651,614,580,547,516,487,460,
        434,410,387,365,345,325,307,290,274,258,244,230,
        217,205,193,183,172,163,154,145,137,129,122,115 },
    {   862,814,768,725,684,646,610,575,543,513,484,457,
        431,407,384,363,342,323,305,288,272,256,242,228,
        216,203,192,181,171,161,152,144,136,128,121,114 },
    {   856,808,762,720,678,640,604,570,538,508,480,453,
        428,404,381,360,339,320,302,285,269,254,240,226,
        214,202,190,180,170,160,151,143,135,127,120,113 },
    {   850,802,757,715,674,637,601,567,535,505,477,450,
        425,401,379,357,337,318,300,284,268,253,239,225,
        213,201,189,179,169,159,150,142,134,126,119,113 },
    {   844,796,752,709,670,632,597,563,532,502,474,447,
        422,398,376,355,335,316,298,282,266,251,237,224,
        211,199,188,177,167,158,149,141,133,125,118,112 },
    {   838,791,746,704,665,628,592,559,528,498,470,444,
        419,395,373,352,332,314,296,280,264,249,235,222,
        209,198,187,176,166,157,148,140,132,125,118,111 },
    {   832,785,741,699,660,623,588,555,524,495,467,441,
        416,392,370,350,330,312,294,278,262,247,233,220,
        208,196,185,175,165,156,147,139,131,124,117,110 },
    {   826,779,736,694,655,619,584,551,520,491,463,437,
        413,390,368,347,328,309,292,276,260,245,232,219,
        206,195,184,174,164,155,146,138,130,123,116,109 },
    {   820,774,730,689,651,614,580,547,516,487,460,434,
        410,387,365,345,325,307,290,274,258,244,230,217,
        205,193,183,172,163,154,145,137,129,122,115,109 },
    {   814,768,725,684,646,610,575,543,513,484,457,431,
        407,384,363,342,323,305,288,272,256,242,228,216,
        204,192,181,171,161,152,144,136,128,121,114,108 },
};