src/fht.h - amarok in Amarok

// FHT - Fast Hartley Transform Class
//
// Copyright (C) 2004  Melchior FRANZ - mfranz@kde.org
//
// This program is free software; you can redistribute it and/or
// modify it under the terms of the GNU General Public License as
// published by the Free Software Foundation; either version 2 of the
// License, or (at your option) any later version.
//
// This program is distributed in the hope that it will be useful, but
// WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
// General Public License for more details.
//
// You should have received a copy of the GNU General Public License
// along with this program; if not, write to the Free Software
// Foundation, 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA
//
// $Id$

#ifndef FHT_H
#define FHT_H

/**
 * Implementation of the Hartley Transform after Bracewell's discrete
 * algorithm. The algorithm is subject to US patent No. 4,646,256 (1987)
 * but was put into public domain by the Board of Trustees of Stanford
 * University in 1994 and is now freely available[1].
 *
 * [1] Computer in Physics, Vol. 9, No. 4, Jul/Aug 1995 pp 373-379
 */
class FHT
{
	int	m_exp2;
	int	m_num;
	float	*m_buf;
	float	*m_tab;
	int	*m_log;

	/**
	 * Create a table of "cas" (cosine and sine) values.
	 * Has only to be done in the constructor and saves from
	 * calculating the same values over and over while transforming.
	 */
	void	makeCasTable();

	/**
	 * Recursive in-place Hartley transform. For internal use only!
	 */
	void	_transform(float *, int, int);

   public:
	/**
	* Prepare transform for data sets with @f$2^n@f$ numbers, whereby @f$n@f$
	* should be at least 3. Values of more than 3 need a trigonometry table.
	* @see makeCasTable()
	*/
	FHT(int);

	~FHT();
	inline int sizeExp() const { return m_exp2; }
	inline int size() const { return m_num; }
	float	*copy(float *, float *);
	float	*clear(float *);
	void	scale(float *, float);

	/**
	 * Exponentially Weighted Moving Average (EWMA) filter.
	 * @param d is the filtered data.
	 * @param s is fresh input.
	 * @param w is the weighting factor.
	 */
	void	ewma(float *d, float *s, float w);

	/**
	 * Logarithmic audio spectrum. Maps semi-logarithmic spectrum
	 * to logarithmic frequency scale, interpolates missing values.
	 * A logarithmic index map is calculated at the first run only.
	 * @param p is the input array.
	 * @param out is the spectrum.
	 */
	void	logSpectrum(float *out, float *p);

	/**
	 * Semi-logarithmic audio spectrum.
	 */
	void	semiLogSpectrum(float *);

	/**
	 * Fourier spectrum.
	 */
	void	spectrum(float *);

	/**
	 * Calculates a mathematically correct FFT power spectrum.
	 * If further scaling is applied later, use power2 instead
	 * and factor the 0.5 in the final scaling factor.
	 * @see FHT::power2()
	 */
	void	power(float *);

	/**
	 * Calculates an FFT power spectrum with doubled values as a
	 * result. The values need to be multiplied by 0.5 to be exact.
	 * Note that you only get @f$2^{n-1}@f$ power values for a data set
	 * of @f$2^n@f$ input values. This is the fastest transform.
	 * @see FHT::power()
	 */
	void	power2(float *);

	/**
	 * Discrete Hartley transform of data sets with 8 values.
	 */
	void	transform8(float *);

	void	transform(float *);
};

#endif

1	// FHT - Fast Hartley Transform Class
2	//
3	// Copyright (C) 2004 Melchior FRANZ - mfranz@kde.org
4	//
5	// This program is free software; you can redistribute it and/or
6	// modify it under the terms of the GNU General Public License as
7	// published by the Free Software Foundation; either version 2 of the
8	// License, or (at your option) any later version.
9	//
10	// This program is distributed in the hope that it will be useful, but
11	// WITHOUT ANY WARRANTY; without even the implied warranty of
12	// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
13	// General Public License for more details.
14	//
15	// You should have received a copy of the GNU General Public License
16	// along with this program; if not, write to the Free Software
17	// Foundation, 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA
18	//
19	// $Id$
20
21	#ifndef FHT_H
22	#define FHT_H
23
24	/**
25	* Implementation of the Hartley Transform after Bracewell's discrete
26	* algorithm. The algorithm is subject to US patent No. 4,646,256 (1987)
27	* but was put into public domain by the Board of Trustees of Stanford
28	* University in 1994 and is now freely available[1].
29	*
30	* [1] Computer in Physics, Vol. 9, No. 4, Jul/Aug 1995 pp 373-379
31	*/
32	class FHT
33	{
34	int m_exp2;
35	int m_num;
36	float *m_buf;
37	float *m_tab;
38	int *m_log;
39
40	/**
41	* Create a table of "cas" (cosine and sine) values.
42	* Has only to be done in the constructor and saves from
43	* calculating the same values over and over while transforming.
44	*/
45	void makeCasTable();
46
47	/**
48	* Recursive in-place Hartley transform. For internal use only!
49	*/
50	void _transform(float *, int, int);
51
52	public:
53	/**
54	* Prepare transform for data sets with @f$2^n@f$ numbers, whereby @f$n@f$
55	* should be at least 3. Values of more than 3 need a trigonometry table.
56	* @see makeCasTable()
57	*/
58	FHT(int);
59
60	~FHT();
61	inline int sizeExp() const { return m_exp2; }
62	inline int size() const { return m_num; }
63	float copy(float , float *);
64	float clear(float );
65	void scale(float *, float);
66
67	/**
68	* Exponentially Weighted Moving Average (EWMA) filter.
69	* @param d is the filtered data.
70	* @param s is fresh input.
71	* @param w is the weighting factor.
72	*/
73	void ewma(float d, float s, float w);
74
75	/**
76	* Logarithmic audio spectrum. Maps semi-logarithmic spectrum
77	* to logarithmic frequency scale, interpolates missing values.
78	* A logarithmic index map is calculated at the first run only.
79	* @param p is the input array.
80	* @param out is the spectrum.
81	*/
82	void logSpectrum(float out, float p);
83
84	/**
85	* Semi-logarithmic audio spectrum.
86	*/
87	void semiLogSpectrum(float *);
88
89	/**
90	* Fourier spectrum.
91	*/
92	void spectrum(float *);
93
94	/**
95	* Calculates a mathematically correct FFT power spectrum.
96	* If further scaling is applied later, use power2 instead
97	* and factor the 0.5 in the final scaling factor.
98	* @see FHT::power2()
99	*/
100	void power(float *);
101
102	/**
103	* Calculates an FFT power spectrum with doubled values as a
104	* result. The values need to be multiplied by 0.5 to be exact.
105	* Note that you only get @f$2^{n-1}@f$ power values for a data set
106	* of @f$2^n@f$ input values. This is the fastest transform.
107	* @see FHT::power()
108	*/
109	void power2(float *);
110
111	/**
112	* Discrete Hartley transform of data sets with 8 values.
113	*/
114	void transform8(float *);
115
116	void transform(float *);
117	};
118
119	#endif

Gitorious