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

double get_base_freq(void)
{
	double base_freq;
	char *envvar = getenv("BASEFREQ");
	if (envvar) {
		base_freq = atof(envvar);
	} else {
		base_freq = 440.0;
	}
	return base_freq;
}

int main(void)
{
	int N = 4096;
	int sample_rate = 44100;
	double base_freq;
	fftw_real in[N], out[N], power_spectrum[N/2+1];
	double pitch[N/2+1];
	rfftw_plan p;
	int k;
	base_freq = get_base_freq();
	p = rfftw_create_plan(N, FFTW_REAL_TO_COMPLEX, FFTW_ESTIMATE);
	for (k = 0; k < N/2+1; ++k)
		pitch[k] = (log((k+0.0)/(N/2+1))+log(sample_rate/base_freq/2))/log(2)*12;
	while(1) {
		double max;
		for (k = 0; k < N; ++k) {
			short s;
			int c;
			c = getchar();
			if (c == EOF) exit(0);
			s = c;
			c = getchar();
			if (c == EOF) exit(0);
			s |= (c << 8);
			in[k] = s / 32768.0;

//			if (scanf("%lf", &f) < 0)
//				exit(0);
		}
		rfftw_one(p, in, out);
		power_spectrum[0] = out[0]*out[0];  /* DC component */
		for (k = 1; k < (N+1)/2; ++k)  /* (k < N/2 rounded up) */
			power_spectrum[k] = out[k]*out[k] + out[N-k]*out[N-k];
		if (N % 2 == 0) /* N is even */
			power_spectrum[N/2] = out[N/2]*out[N/2];  /* Nyquist freq. */
		max = -100;
		for (k = 0; k < N/2+1; ++k) {
			power_spectrum[k] = log(power_spectrum[k]/N);
			if (power_spectrum[k] > max)
				max = power_spectrum[k];
		}
		for (k = 1; k < N/2; ++k)
//			if (power_spectrum[k] > max-3 && power_spectrum[k] > power_spectrum[k-1] && power_spectrum[k] > power_spectrum[k+1])
				printf("%g\t%g\n", pitch[k], power_spectrum[k]);
//				printf("%g\t%g\n", (pitch[k-1]*power_spectrum[k-1]+pitch[k]*power_spectrum[k]+pitch[k+1]*power_spectrum[k+1])/(power_spectrum[k-1]+power_spectrum[k]+power_spectrum[k+1]), (power_spectrum[k-1]+power_spectrum[k]+power_spectrum[k+1])/3);
		printf("\n");
	}
	rfftw_destroy_plan(p);
	exit(0);
}