// TRENDLINE SYNTHESIS
// by athibaul, 2023/05/15

// See Colin Pykett's article on Trendline Synthesis:
// http://www.colinpykett.org.uk/trendline-synthesis-new-music-synthesis-technique.htm
//
// "Trendline synthesis" is an additive synthesis technique for the sound of organ pipes
// based on a small number of parameters describing the spectrum of the sound.
// The sound spectrum is assumed to be a series of harmonic pitches, whose amplitudes
// are described using piecewise linear functions in dB/octave, with a break point at
// a given harmonic number.
// The even-numbered harmonics can be boosted or suppressed independently.
//
// bp: harmonic number of the breakpoint
// s1: slope of group 1 trendline, in dB/octave
// s2: slope of group 2 trendline, in dB/octave
// dl: difference between odd and even pitches, in dB. Use negative numbers to suppress even harmonics.

input.instrument = 6; // min=0, max=6, step=1 (Diapason, Strings, Flute, Stopped Flute, Trumpet, Clarinet, Diapason II)

// Play using the "Pykett 2019" temperament
// http://www.colinpykett.org.uk/more-on-handels-temperament.htm
//const temperament_shift = [0,0,0,0,0,0,0,0,0,0,0,0];
const temperament_shift = [3.11,0.19,0,0.17,0,2.13,0.6,0,0.2,0,1.19,.98];
//const temperament_shift = [0,-24,-7,-31,-14,3,-21,-3,-27,-10,7,-17]; // 1/4 comma meantone
const my_midi_to_hz = (nn) => midi_to_hz(nn + temperament_shift[nn%12]/100) * (415/440);


const diapason = {nn: c4, bp:1, s1: -16, s2: -16, dl: 0};
const strings = {nn: c4, bp:4, s1: 3, s2: -20, dl: 0};
const flute = {nn: c4, bp:4, s1:-25, s2:-30, dl:-6};
const stoppedFlute = {nn: c4, bp: 1, s1:-18, s2:-18, dl:-20};
const trumpet = {nn: c4, bp: 7, s1: -3, s2: -22, dl: 0};
const clarinet = {nn: c4, bp: 6, s1: -6, s2:-30, dl:-15};

// Scaling of a synthetic organ stop
// http://www.colinpykett.org.uk/scaling-synthetic-sample-sets.htm
const diapason8Foot = [
    {nn:c2, bp:2, s1:15, s2:-20, dl:0},
    {nn:c3, bp:4, s1:-10, s2:-24, dl:0},
    {nn:c4, bp:4, s1:-7, s2:-24, dl:0},
    {nn:c5, bp:2, s1:-5, s2:-24, dl:0},
    {nn:c6, bp:2, s1:0, s2:-30, dl:0},
];

const instruments = [[diapason], [strings], [flute], [stoppedFlute], [trumpet], [clarinet], diapason8Foot];

function instr_to_params(instr, nn) {
    // Convert a scaling to a set of parameters for a given note number
    // Assumes the set of pitches is sorted!
    let N = instr.length;
    // If too low, return the params of the lowest note
    if(nn <= instr[0].nn) {
        return instr[0];
    }
    // If too high, return the params of the highest note
    if(nn >= instr[N-1].nn) {
        return instr[N-1];
    }
    // Otherwise, find the corresponding note and interpolate
    let k = 0;
    while(nn > instr[k+1].nn) {
        k++;
        debug.log("k", k);
    }
    /*
    //let t = (nn - instr[k].nn) / (instr[k+1].nn - instr[k].nn);
    return {bp:instr[k].bp, // Math.round(lerp(instr[k].bp, instr[k+1].bp, t))
            s1:instr[k].s1, //lerp(instr[k].s1, instr[k+1].s1, t),
            s2:instr[k].s2, //lerp(instr[k].s2, instr[k+1].s2, t),
            dl:instr[k].dl, //lerp(instr[k].dl, instr[k+1].dl, t),
    };
    */
    let t = (nn - instr[k].nn) / (instr[k+1].nn - instr[k].nn);
    return {bp:Math.round(lerp(instr[k].bp, instr[k+1].bp, t)),
            s1:lerp(instr[k].s1, instr[k+1].s1, t),
            s2:lerp(instr[k].s2, instr[k+1].s2, t),
            dl:lerp(instr[k].dl, instr[k+1].dl, t),
    };
}   

function amplitudes(params, nharm) {
    let amps = new Float32Array(nharm);
    amps[0] = 0;
    amps[1] = 1;
    let s1_exponent = params.s1 * 0.166096404744; // s1 * log(10) / (20 * log(2))
    let s2_exponent = params.s2 * 0.166096404744; // s2 * log(10) / (20 * log(2))
    for(let i=2; i<nharm; i++) {
        if(i <= params.bp) {
            amps[i] = amps[i-1] * (i / (i-1)) ** s1_exponent;
        } else {
            amps[i] = amps[i-1] * (i / (i-1)) ** s2_exponent;
        }
    }
    // Apply the amplitude reduction to the even harmonics
    let factor = 10**(params.dl/20);
    for(let i=2; i<nharm; i+=2) {
        amps[i] *= factor;
    }
    // Normalize the total RMS amplitude
    let sum_of_squares = amps.reduce((s,v) => s + v**2);
    amps = amps.map((v) => v / Math.sqrt(sum_of_squares));
    return amps;
}

function wtable(amps, nsamp) {
    let wt = new Float32Array(nsamp);
    for(let i=1; i<amps.length; i++) {
        let phase = 0; //Math.random() * 2 * Math.PI;
        for(let k=0; k<nsamp; k++) {
            wt[k] += amps[i] * Math.sin(2 * Math.PI * i * k / nsamp + phase);
        }
    }
    return wt;
}

// Define a modified adsr envelope
/*
const myAdsr = env.def(class extends envBase {
    constructor(options) {
        super();
        const duration = options.duration;
        const curve = options.curve;

        this.init({
            levels: [0, 0, 1, options.sustain, 0],
            times: [options.dly, options.attack, options.decay, options.release],
            duration: duration,
            curves: Array.isArray(curve) ? curve : [0, 0, curve, curve],
            releaseNode: 3
        });
    }
}, {
    dly: 0,
    attack: 0,
    decay: 0,
    duration: 0,
    sustain: 1,
    release: .5,
    curve: -2,
    name: 'myAdsr'
});
*/

const dlyAdsr = env.def(
    class {
        constructor(options) {
            this.v = 0;
        }

        duration(options) {
            return options.duration + options.dly + 3 * options.release;
        }

        value(tick, options) {
            if(tick < options.dly) {
                return 0;
            }
            if(tick < options.duration) {
                let a0 = clamp01(3 * ditty.dt / options.attack);
                this.v += a0 * (1 - this.v);
                return this.v;
            }
            let a0 = clamp01(3 * ditty.dt / options.release);
            this.v += a0 * (0 - this.v);
            return this.v;
        }
    },
    { attack: 60e-3, release:0.5, dly: 0, duration: 1, name: 'dlyAdsr' }
);

const tsSynth = synth.def(class {
    constructor(o) {}
    process(n,e,t,o) {
        if(!this.initialized) {
            this.freq = my_midi_to_hz(o.note);
            let nsamp = Math.ceil(1 / (this.freq * ditty.dt));
            debug.log("nsamp", nsamp);
            let nharm = Math.floor(20000 / this.freq);
            let instr = instruments[input.instrument];
            let params = instr_to_params(instr, o.note);
            let amps = amplitudes(params, nharm);
            this.wt = wtable(amps, nsamp);
            this.pos = 0;
            this.dpos = ditty.dt * this.freq * nsamp;
            this.initialized = true;
        }
        
        let N = this.wt.length;
        this.pos = (this.pos + this.dpos) % N;
        let k = ~~this.pos;
        let kf = this.pos - k;
        return lerp(this.wt[k], this.wt[(k+1)%N], kf) * e.value;
        //return this.pos / N;
        //return this.wt[k];
    }
}, {env:dlyAdsr, amp:0.1});

//tsSynth.play(c4, {params: clarinet});




// ======================================================================================
// Play a score using this synth
// Adapted from my own ditty Rameau | Allemande
// https://dittytoy.net/ditty/bd3257fcfd
// ======================================================================================



// Adapted from "Allemande" by Jean-Philippe Rameau
// in Pièces de clavecin avec une méthode, RCT 2-4 (Rameau, Jean-Philippe)
// https://imslp.org/wiki/Pi%C3%A8ces_de_clavecin_avec_une_m%C3%A9thode,_RCT_2-4_(Rameau,_Jean-Philippe)
const rh = ("'4 g ,5b>'3+d5>e>>g 3 g 3.V>a3 >b 4.>e 3+f ^4.+>>f 3e +^4.>d 3+c ,>b a >@g +f >4.g3 b '4e3 ,g >5+f3 b >4a3 '+d 4>.+d4 e,>b>>g "
+ "3 b 'v>c d 5>c4 ,>e >5+f3 'c, >^b a >2a 3Vb2 3>a >g a b 'c d ,b >>4.'e 2+f g >>4.@c 3,b >5+f@>a "
+ "3 'd >,a 'd ,b g b g >Va 'd ,4d3 'd ,>b >g >b >g >>Va >'d >,4d3 >'d >>,4.^b 3a 4.>+f>@a 3g >5.g>d>>,b");
const lh = ("5 ,e4 5b4 'e d @c ,b 5a'4>a g 5+f,4b '+c +d ,b '5e4 ,g Va b 5.e4. 3b '5e4 ,e "
+ "5a3 ,a b 'c ,5d4 'd ,5g '3g a b g '4c3 5g3 ,4b a3 '4+f3 ,4g3 'g 5d4 3,d e +f g a +f 5,g4 'g ^+f d ,5g4 'g ^+f d v5g4 c 5d4 ,d 5.g'>>>>>g");



ditty.bpm = 40;
const dly_ticks = 1/20;
const dur_to_duration = (dur, dot) => 2**(dur-5) * (1 + 0.5*dot); // Convert Musescore numbers to duration
const nn_to_pan = (nn) => clamp((nn - 64) / 40, -1, 1);

function playScore(score) {
    var alt = 0; // Accidentals: +1 for sharp, -1 for flat (entering a note resets the accidentals)
    var dur = 5; // Duration in Musescore numbers (5 is quarter note, lower is shorter)
    var dot = 0; // Does it have a dot? (changing duration resets the dot)
    var orn = 0; // Ornamentation: 'v' is "pincé" (semitone), 'V' is pincé (whole tone), '^' is trill (semitone), '@' is trill (whole tone)
    // (Sleeping resets the ornamentation)
    var dly = 0; // Delay the onset of the note : ">" is delayed, ">>" is more delayed, etc.
    // (Sleeping resets the delay)
    var octave = 4; // The score starts in octave 4. All notes entered are set in the current octave.
    // An octave shift is performed with "'" and ",".
    // Letters a-g play the notes, but do NOT sleep for the designated duration.
    // Space " " sleeps for the duration of the note (so "g " is one note followed by its duration, but "g  " is one note followed by a rest).
    // This allows a "hack" to make notes sustain longer than their rythmic value.
    // e.g. "5g4 b " plays two eighth notes, a G followed by a B, but sustains the G for the full duration of the beat.
    
    for(let char of score) {
        switch(char) {
            case "+":
                alt++; break;
            case "-":
                alt--; break;
            case ".":
                dot = 1 - dot; break;
            case "v": 
            case "V":
            case "^":
            case "@":
                orn = char; break;
            case ">":
                dly+=dly_ticks; break;
            case "'":
                octave++; break;
            case ",":
                octave--; break;
            case "2":
                dur = 2; dot = 0; break;
            case "3":
                dur = 3; dot = 0; break;
            case "4":
                dur = 4; dot = 0; break;
            case "5":
                dur = 5; dot = 0; break;
            case "a":
            case "b":
            case "c":
            case "d":
            case "e":
            case "f":
            case "g":
                // Play the note
                let nn = notes[char+octave] + alt;
                playNote(nn, dur_to_duration(dur, dot), orn, dly);
                alt = 0;
                break;
            case " ":
                sleep(dur_to_duration(dur, dot));
                orn = 0;
                dly = 0;
                break;
            default:
                break;
        }
    }
}

function playNote(nn, duration, orn, dly) {
    var gap = 0;
    if(orn) {
        switch(orn) {
            case "V":
                gap = -2; break;
            case "v":
                gap = -1; break;
            case "^":
                gap = 1; break;
            case "@":
                gap = 2; break;
        }
        // Play note with trill/pincé
        var t = dly;
        var dt = 1/24 * (1 + 0.1*Math.random());
        var count = 0;
        tsSynth.play(nn,     {dly:t, duration: t+=2*dt, pan:nn_to_pan(nn)});
        tsSynth.play(nn+gap, {dly:t, duration: t+=dt, pan:nn_to_pan(nn)});
        while(duration - t > 4*dt && ++count < 3) {
            tsSynth.play(nn,     {dly:t, duration: t+=dt, pan:nn_to_pan(nn)});
            tsSynth.play(nn+gap, {dly:t, duration: t+=dt, pan:nn_to_pan(nn)});
        }
        tsSynth.play(nn,     {dly:t, duration: duration, pan:nn_to_pan(nn)});
    } else {
        // Simply play the note
        tsSynth.play(nn, {duration: duration, dly:dly + 0.01*Math.random(), pan:nn_to_pan(nn)});
    }
}



////////////////////////////////// REVERB ////////////////////////////////////
class Delayline {
    constructor(n) {
        this.n = n;
        this.p = 0;
        this.data = new Float32Array(n);
    }
    current() {
        return this.data[this.p]; // using lastOut results in 1 sample excess delay
    }
    clock(input) {
        this.data[this.p] = input;
        if(++this.p >= this.n) {
            this.p = 0;
        }
    }
}
const ISQRT2 = Math.sqrt(0.5);
const SQRT2 = Math.sqrt(2);
const SQRT8 = Math.sqrt(8);
const ISQRT8 = 1/SQRT8;
const fibodelays = [1410,1662,1872,1993,2049,2114,2280,2610]; // X^8 = X+1
//const fibodelays = [1467,1691,1932,2138,2286,2567,3141,3897]; // X^8 = X+2
const fiboshort = [465, 537, 617, 698, 742, 802, 963, 1215]; // X^8 = X+2
let meandelay = fibodelays[3] * ditty.dt;

const reverb = filter.def(class {
    constructor(options) {
        this.outgain = 0.3;
        
        this.dls = [];
        this.s0 = new Float32Array(8); // Lowpass filter memory
        for(let i=0; i<8; i++) {
            this.dls.push(new Delayline(fibodelays[i]));
            //this.dls.push(new Delayline(fiboshort[i]));
            this.s0[i] = 0;
        }
        this.a0 = clamp01(2*Math.PI*options.cutoff*ditty.dt);
        
        this.dls_nr = []; // Non-recirculating delay lines
        for(let i=0; i<4; i++) {
            this.dls_nr.push(new Delayline(fiboshort[i*2]));
        }
    }
    process(input, options) {
        // First stage: non-recirculating delay lines
        let s0 = input[0], s1 = input[1];
        for(let i=0; i<4; i++) {
            let u0 = 0.6*s0 + 0.8*s1, u1 = 0.8*s0 - 0.6*s1;
            let u1p = this.dls_nr[i].current();
            this.dls_nr[i].clock(u1);
            s0 = u0;
            s1 = u1p;
        }
        let v0 = (s0 + s1);
        let v1 = (s0 - s1);
        //let s0 = -0.6*input[0] + 0.8*input[1];
        //let s1 = 0.8*input[0] + 0.6*input[1];
        
        // Second stage: recirculating delay lines
        let rt60 = options.rt60;
        let loopgain = 10 ** (-3*meandelay / rt60) * ISQRT8;
        let higain = 10 ** (-3*meandelay / options.rtHi) * ISQRT8;
        let v = this.dls.map( (dl) => dl.current());
        // Fast Walsh-Hadamard transform
        // https://formulasearchengine.com/wiki/Fast_Walsh%E2%80%93Hadamard_transform
        let w = [v[0]+v[4], v[1]+v[5], v[2]+v[6], v[3]+v[7],
                 v[0]-v[4], v[1]-v[5], v[2]-v[6], v[3]-v[7]];
        let x = [w[0]+w[2], w[1]+w[3], w[0]-w[2], w[1]-w[3],
                 w[4]+w[6], w[5]+w[7], w[4]-w[6], w[5]-w[7]];
        let y = [x[0]+x[1], x[0]-x[1], x[2]+x[3], x[2]-x[3],
                 x[4]+x[5], x[4]-x[5], x[6]+x[7], x[6]-x[7]];
        y[0] += v0;
        y[2] += v1;
        for(let i=0; i<8; i++) {
            let hipass = y[i] - this.s0[i];
            this.dls[i].clock(this.s0[i] * loopgain + hipass * higain); // Mix lowpass and hipass in different amounts
            this.s0[i] += this.a0 * hipass;
        }
        return [lerp(input[0], v[0], options.mix),
                lerp(input[1], v[2], options.mix)];
    }
}, {mix:0.8, rt60:2, cutoff:8000, rtHi:0.7});
///////////////////////////////////////////////////////////////////////////////////////////////




loop( () => {
    playScore(rh);
    sleep(3);
}, {name: "Right hand"}).connect(reverb.create());

loop( () => {
    playScore(lh);
    sleep(3);
}, {name: "Left hand"}).connect(reverb.create());