Bytebeat implementation of a classical piece of music, using a single sinusoid and comparing the naive, discontinuous signal with an improved, phase-accumulation based version.

classical_star_phase.c

/*
 * Date: 2018-11-26
 * Author of this C code: Wiebe-Marten Wijnja (Qqwy).
 * Author of original JavaScript-version: Unknown.
 *
 * This version uses phase-accumulation to ensure that the signal remains continuous throughout,
 * and that there are no nasty clicks when we change to the next note.
 * You can add the directive `DISCONTINUOUS VERSION` (using e.g. the `-D` flag of `gcc` and `clang`)
 * to compile it as the original version, that bases its samples only directly on the current sample,
 * and therefore includes clicks whenever the note changes.
 *
 * Compilation:
 * `gcc classical_star_phase.c -lm -o classical_star_phase` for the normal version based on phase-offset.
 * `gcc classical_star_phase.c -lm -o classical_star_phase -DDISCONTINUOUS_VERSION` for the version with the time-based implementation that results in a discontinous signal.
 *
 * Example usage:
 * `./classical_star_phase | aplay`  plays at 8 kHz.
 * `./classical_star_phase 128000 | aplay -r128000` plays at 128 kHz.
 *
 * If you do not have `aplay` but do have `play` (which is part of Sox), you can do:
 * `./classical_star_phase | play  -t raw -b 8 -e unsigned -r 8000 -
 *
 *
 * To create a spectrogram,
 * it was first stored as a .wav file using `./classic_star_phase 8000 | sox  -t raw -b 8 -e unsigned -r 8000 - classic_star_phase.wav`
 * and then the call: `sox classic_star_phase.wav -n spectrogram -x2000 -Y2000 -z80 -oclassic_star_spectrogram.png`  was used.
 *
 * Based on:
 * http://wurstcaptures.untergrund.net/music/?oneliner=30*cos(t*Math.pow(2%2C(%22B*918%2F916-918%2F91B*918%2F916-918%2F91%3E*%3B2%3A1%3B26%2F%3B2%3A1%3B2%3E*%3B2%3A1%3B26%2F%22%2B%20%22%3B2%3A1%3B2A*%3B291%3B28%2F%3B291%3B2A*%3B291%3B28%2F%3B291%3B2B*%3D-%3B%2C%3D-91%3D-%3B%2C%3D-B*%3D-%3B%2C%3D-91%3D-%3B%2C%3D-E*%3E6%3D4%3E69%22%2B%20%222%3E6%3D4%3E6E*%3E6%3D4%3E692%3E6%3D4%3E6D*%3C3%3A1%3C380%3C3%3A1%3C3D*%3C3%3A1%3C380%3C3%3A1%3C3D(%3D4%3C3%3D481%3D4%3C3%3D4D(%3D4%3C3%3D4%22%2B%20%2281%3D4%3C3%3D4B(%3A18%2F%3A16.%3A18%2F%3A1B(%3A18%2F%3A16.%3A18%2F%3A1B%26%3B2%3A1%3B26%2F%3B2%3A1%3B2B%26%3B2%3A1%3B26%2F%3B2%3A1%3B2%40%26%3B%2C9*%3B%22%2B%20%22%2C8%2F%3B%2C9*%3B%2C%40%26%3B%2C9*%3B%2C8%2F%3B%2C9*%3B%2C%40%25%3D-%3B%2C%3D-91%3D-%3B%2C%3D-%40%25%3D-%3B%2C%3D-91%3D-%3B%2C%3D-%3E*%3D-%3B%2C%3D-92%3D-%3B%2C%3D-%3E*%3D-%3B%2C%22%2B%20%22%3D-92%3D-%3B%2C%3D-%3E%2C8%2F6-8%2F428%2F6-8%2F%3E%2C8%2F6-8%2F428%2F6-8%2F%3D-412%2F4192412141%3D-412%2F4192412141%3B-634%22%2B%20%221613%2F634163%3B-6341613%2F634163%3B%2C8%2F6-8%2F528%2F6-8%2F%3B%2C8%2F6-8%2F528%2F6%22).charCodeAt(t%3E%3E9)%2F12-7))&rate=8000
 * Thanks to the people on the Pouet Forums, especially those in the ByteBeat topic! http://www.pouet.net/topic.php?which=8357
 * */


#include <stdio.h>
#include <math.h>
#include <stdlib.h>

#define SAMPLE_RATE 8000

const char notes[] = "B*918/916-918/91B*918/916-918/91>*;2:1;26/;2:1;2>*;2:1;26/;2:1;2A*;291;28/;291;2A*;291;28/;291;2B*=-;,=-91=-;,=-B*=-;,=-91=-;,=-E*>6=4>692>6=4>6E*>6=4>692>6=4>6D*<3:1<380<3:1<3D*<3:1<380<3:1<3D(=4<3=481=4<3=4D(=4<3=481=4<3=4B(:18/:16.:18/:1B(:18/:16.:18/:1B&;2:1;26/;2:1;2B&;2:1;26/;2:1;2@&;,9*;,8/;,9*;,@&;,9*;,8/;,9*;,@%=-;,=-91=-;,=-@%=-;,=-91=-;,=->*=-;,=-92=-;,=->*=-;,=-92=-;,=->,8/6-8/428/6-8/>,8/6-8/428/6-8/=-412/4192412141=-412/4192412141;-6341613/634163;-6341613/634163;,8/6-8/528/6-8/;,8/6-8/528/6-8/";

double freq(int t, int sample_rate) {
  size_t index = (t * 16) / (sample_rate);
  if(index >= sizeof(notes)) {
    exit(EXIT_SUCCESS);
  }
  index = index % 2 ? index - 1 : index;
  char note = notes[index % sizeof(notes)];
  return pow(2.0L,(double)note/12.0L - 9.0) * 8000 / sample_rate;
}

int double2int(double sample) {
  return trunc(sample * 128.0 + 128.0);
}


int main(int argc, char** argv) {
  int sample_rate = 8000;
  if(argc >= 2){
    sample_rate = atoi(argv[1]);
  }
  double curr_phase = 0;
	for(int t = 0;;++t) {
    double sample;
    double curr_freq = freq(t, sample_rate);
#ifdef DISCONTINUOUS_VERSION
    sample = cos(2 * M_PI * t * curr_freq);
#else
    // Note: We don't use irrational numbers like PI here directly to reduce accumulated rounding errors somewhat:
    curr_phase += curr_freq;
    curr_phase = fmod(curr_phase, 1);
    sample = cos(2 * M_PI * curr_phase);
#endif
    putchar(double2int(0.5*sample));
  }
}

Here is a comparison of the two variants as seen as spectgrograms:

Continuous (Using phase-accumulation)

Discontinuous (Using time-in-samples multiplication)