/* Spectrum_to_Ltas.c */
/* Paul Boersma, January 5, 1998 */

#include "Spectrum_to_Ltas.h"

Ltas Spectrum_to_Ltas (Spectrum me, double bandWidth) {
	Ltas thee = NULL;
	long numberOfBands = ceil ((my xmax - my xmin) / bandWidth);
	long i, j;
	if (bandWidth <= my dx) {
		Melder_error ("Bandwidth must be greater than %f.", my dx);
		goto error;
	}
	if (! (thee = new (Ltas)) ||
		! Matrix_init (thee, my xmin, my xmax, numberOfBands,
			bandWidth, my xmin + 0.5 * bandWidth, 1, 1, 1, 1, 1))
		goto error;
	for (i = 1; i <= numberOfBands; i ++) {
		double fmin = my xmin + (i - 1) * bandWidth;
		double fmax = my xmin + i * bandWidth;
		long jmin = Sampled_xToHighIndex (me, fmin);
		long jmax = Sampled_xToLowIndex (me, fmax);
		if (jmin < 1) jmin = 1; else if (jmin > my nx) jmin = my nx;
		if (jmax < 1) jmax = 1; else if (jmax > my nx) jmax = my nx;
		for (j = jmin; j <= jmax; j ++)
			thy z [1] [i] += my z [1] [j] * my z [1] [j] + my z [2] [j] * my z [2] [j];
		if (thy z [1] [i] != 0.0) {
			double oneSidedPowerSpectralDensity = 2 * thy z [1] [i] / (jmax - jmin + 1) * my dx;
			thy z [1] [i] = 10 * log10 (oneSidedPowerSpectralDensity / 4.0e-10);
		}
	}
	return thee;
error:
	forget (thee);
	return Melder_errorp ("Spectrum_to_Ltas: not performed.");
}

Ltas Spectrum_to_Ltas_1to1 (Spectrum me) {
	long i;
	Ltas thee = new (Ltas); cherror
	Matrix_init (thee, my xmin, my xmax, my nx, my dx, my x1, 1, 1, 1, 1, 1); cherror
	for (i = 1; i <= my nx; i ++) {
		double value = my z [1] [i] * my z [1] [i] + my z [2] [i] * my z [2] [i];
		if (value != 0.0) {
			double oneSidedPowerSpectralDensity = 2 * value * my dx;
			value = 10 * log10 (oneSidedPowerSpectralDensity / 4.0e-10);
		}
		thy z [1] [i] = value;
	}
end:
	iferror forget (thee);
	return thee;
}

/* End of file Spectrum_to_Ltas.c */

