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

#define MAX_R 25

char *pitches[] = {"F2 ", "F#2", "G2 ", "G#2", "A2 ", "A#2", "B2 ", "C3 ", "C#3", "D3 ", "D#3", "E3 ",
	           "F3 ", "F#3", "G3 ", "G#3", "A3 ", "A#3", "B3 ", "C4 ", "C#4", "D4 ", "D#4", "E4 ",
                   "F4 ", "F#4", "G4 ", "G#4", "A4 ", "A#4", "B4 ", "C5 ", "C#5", "D5 ", "D#5", "E5 ",
                   "F5 ", "F#5", "G5 ", "G#5", "A5 ", "A#5", "B5 ", "C6 ", "C#6", "D6 ", "D#6", "E6 "};

double R[MAX_R] = {1e3 , 1.2e3, 1.5e3, 1.8e3, 2.2e3, 2.7e3, 3.3e3, 3.9e3, 4.7e3, 5.6e3, 6.8e3, 8.2e3, 10e3, 12e3, 15e3, 18e3, 22e3, 27e3, 33e3, 39e3, 47e3, 56e3, 68e3, 82e3, INFINITY};  // Standard E12 resistance values 1-82 k Ohms
//double R[MAX_R] = {1.5e3, 2e3, 2.2e3, 2.7e3, 3.3e3, 5.1e3, 5.6e3, 6.2e3, 10e3, 15e3, 20e3, 22e3, 33e3, 39e3, 48.1e3, 62e3, INFINITY};  // This is what we have in metal-film resistors right now!

double filtererror, filtererrormax = 0;

int parts[MAX_R] = {0};  // Serves to count how many of each resistor size is needed to construct all 48 filter networks

double C1 = 0.22e-6;  // Capacitor size used in all 48 RC low-pass filter networks

int filterpick(double);

int main (void)
{
  double xtal_osc = 4e6;
  double idealpitch, actualpitch;
  double fcf = 2.0;  // Filter cutoff frequency factor (how many times lower the f_c is compared to the square-wave tone)
                     // 2.0 produces a sawtooth waveshape with inverse-exponential rising and falling curves
                     // 4.0 produces nearly a triangle wave
  int note, n;
  
  printf("Master oscillator frequency = %.3f MHz \n\n", xtal_osc/1e6);
  printf("Each 12-bit counter's recycle frequency is twice that of its highest-octave pitch due to the 1/2 frequency\n");
  printf("division of the final 4-bit counter used to generate the four octaves from each frequency divider.\n");
  printf("For example, the 12-bit `F' counter takes the crystal oscillator's frequency and divides it by 2863 to \n");
  printf("yield 1397.14 Hz.  This signal clocks a 4-bit binary counter that outputs 698.57 Hz (F5), 349.28 Hz (F4),\n");
  printf("174.64 Hz (F3), and 87.32 Hz (F2) tones.\n\n");
  printf("Pitches and frequency-divider parameters for the highest octave (error in cents same for all octaves): \n");

  for (note = 0 ; note < 12 ; ++note)  // Calculates just a single octave
  {
    printf("%s  ", pitches[note + 36]);  // Calculates only the upper octave

    idealpitch = 440.0 * pow(2,(note + 8)/12.0); // Begins with an F5 pitch, eight semitones above 440 Hz (A4)
    printf("Ideal = %f Hz  ", idealpitch);
    
    printf("Modulus = %i  ", (int) round(xtal_osc / (2.0 * idealpitch)));

    actualpitch = xtal_osc / round(xtal_osc / idealpitch);
    printf("Actual = %f Hz  ", actualpitch);

    printf("Pitch error = %f cents  ", log2(actualpitch / idealpitch) * 1200);

    printf("Beat f = %f Hz  ", (actualpitch - idealpitch));

    printf("PRE = 0x%X ", (int) pow(2,12) - (int) round(xtal_osc / (2.0 * idealpitch)));
    
    printf("\n");
  }


  printf("\n\n");
  printf("Simple first-order low-pass filter RC values for four octaves, with cutoff frequencies 1/%.1f that of the\n", fcf);
  printf("pitch's fundamental to produce a waveshape with less harmonic content than a square wave:\n\n");

  for (note = 0 ; note < 48 ; ++note)  // Calculates four octaves
  {
    printf("%s  ", pitches[note]);

    idealpitch = 55.0 * pow(2,(note + 8)/12.0); // Begins with an F2 pitch, eight semitones above 55 Hz (A1)
    
    actualpitch = xtal_osc / round(xtal_osc / idealpitch);

    printf("Actual pitch = %.5g Hz   Filter cutoff = %.5g Hz   ", actualpitch, actualpitch / fcf);

    filterpick(actualpitch / fcf);

    printf("\n");
  }

  printf("\n");
  printf("Greatest filter error = %f percent\n\n", filtererrormax);


  printf("\n");
  printf("Resistor Bill Of Materials:\n");
  
  for (n = 0 ; n < (MAX_R - 1) ; ++n)  // Print number of resistors needed
  {
    printf("%g k Ohm resistors needed = %i\n", R[n]/1e3, parts[n]);
  }

  return 0;
}


// This function picks R and C values to form a first-order lowpass
// filter network with a cutoff frequency a fraction of the specified
// pitch.  This reduces the 5 Volt peak-peak square-wave signal down 
// to a nearly-triangular waveshape of about 1.5 Volts peak-peak.
int filterpick(double f_c)
{
  double f, rp, error, error_least = 99.9;
  double r1_match, r2_match, r3_match;
  int r1, r2, r3, n;
 
  for (r1 = 0 ; r1 < (MAX_R - 1) ; ++r1)
    for (r2 = 0; r2 < MAX_R; ++r2)
      for (r3 = 0; r3 < MAX_R; ++r3)
      {
        rp = 1 / (1 / R[r1] + 1 / R[r2] + 1 / R[r3]);
        f = 1 / (2 * M_PI * rp * C1);
        error = fabs (f - f_c) / f_c;

        if (error < error_least)
        {
          error_least = error;
          r1_match = R[r1]; 
          r2_match = R[r2]; 
          r3_match = R[r3]; 
        }
      }

  rp = 1 / (1 / r1_match + 1 / r2_match + 1 / r3_match);
  f = 1 / (2 * M_PI * rp * C1);

  error = (f - f_c) / f_c * 100;

  printf ("LP filter error = %f percent \n", error);
  
  printf ("C1 = %.2fu  ", C1 * 1e6);

  printf ("R1 = %gk  ", r1_match / 1e3);

  if (r2_match < 100e3)
    printf ("R2 = %gk  ", r2_match / 1e3);
  
  if (r3_match < 100e3)
    printf ("R3 = %gk  ", r3_match / 1e3);

  printf ("(Rparallel = %.3gk)  \n", rp / 1e3);

  for (n = 0 ; n < MAX_R ; ++n)
  {
    if (r1_match == R[n])
      ++parts[n];
    
    if (r2_match == R[n])
      ++parts[n];
    
    if (r3_match == R[n])
      ++parts[n];
  }

  if (fabs(error) > filtererrormax)
    filtererrormax = fabs(error);

  return 0;
}


