diff options
Diffstat (limited to 'Graphic_Equalizer_v1.0/src/fft.hcc')
| -rw-r--r-- | Graphic_Equalizer_v1.0/src/fft.hcc | 513 |
1 files changed, 513 insertions, 0 deletions
diff --git a/Graphic_Equalizer_v1.0/src/fft.hcc b/Graphic_Equalizer_v1.0/src/fft.hcc new file mode 100644 index 0000000..81098c0 --- /dev/null +++ b/Graphic_Equalizer_v1.0/src/fft.hcc @@ -0,0 +1,513 @@ +/*! \file fft.hcc
+ *
+ * \section generic This modules will take care of the actual FFT calculation
+ * on the samples. Besides the FFT this module also will
+ * equalize the audio signal according to the setting made by the user.
+ *
+ * \section project Project information.
+ * Project Graphic Equalizer\n
+ * \author M. Lauwerijssen
+ * \date 20041110
+ * \version 0.1
+ *
+ * \section copyright Copyright
+ * Copyright ©2004 Koninklijke Philips Electronics N.V. All rights reserved
+ *
+ * \section history Change history
+ * 20041110: M. Lauwerijssen\n Initial version
+ *
+ ********************************************************************/
+#include <stdlib.hch>
+#include "pal_master.hch"
+
+#include "audio.hch"
+#include "weights_256.hch"
+#include "configuration.hch"
+#include "xilinxmult.hch"
+#include "fft.hch"
+
+#if HAVE_DEBUG
+ #include "debug.hch"
+#endif
+
+/* Define two multi-port RAMs for FFT calculation; one for real and one for imaginary values
+ * Extra block RAM settings are defined to make sure read and write actions can be performed
+ * within one clock-cycle.
+ * Left out extra settings on new board the clock changes TODO !!!!
+ */
+#if HARDWARE_MULTIPLY
+mpram
+{
+ ram signed 18 rwrite[256];
+ rom signed 18 read[256];
+} real with {block = "BlockRAM"};
+
+mpram
+{
+ ram signed 18 rwrite[256];
+ rom signed 18 read[256];
+} imaginary with {block = "BlockRAM"};
+#else
+mpram
+{
+ ram signed 24 rwrite[256];
+ rom signed 24 read[256];
+} real with {block = "BlockRAM"};
+
+mpram
+{
+ ram signed 24 rwrite[256];
+ rom signed 24 read[256];
+} imaginary with {block = "BlockRAM"};
+#endif
+
+// multiplication factors for equalizer function
+ram signed 7 eq_settings[16] = {0,2,4,7,10,13,16,19,22,26,30,35,41,48,55,63};
+
+#if HARDWARE_MULTIPLY
+#define DC_COMPONENT 0
+#else
+#define DC_COMPONENT 8470527
+#endif
+
+/*! \fn macro proc multiply(result, op_a, op_b);
+ * \brief Procedure used for multiply-ing
+ *
+ * \param result variable containing the result of the multiply procedure
+ * \param op_a integer value to be multiplied.
+ * \param op_b integer value to be multiplied.
+ *
+ * \return Procedure returns through variable.
+ * \retval signed 36
+ */
+macro proc multiply(result, op_a, op_b)
+{
+#if HARDWARE_MULTIPLY
+ xilinxmult(result, op_a, adjs(op_b,18));
+#else
+ result = (adjs(op_a,38))*(adjs(op_a,38));
+#endif
+}
+
+
+
+
+/*! \fn void calculate_fft(unsigned 1 select_inverse)
+ * \brief This routine performs the Fast Fourier Transform for calculation of the frequency spectrum
+ *
+ * \param select_inverse determines if a FFT or iFFT has to be calculated
+ *
+ * \return nothing
+ * \retval void
+ *
+ * cost 12391 cycles
+ */
+void calculate_fft(unsigned 1 select_inverse)
+{
+ unsigned 4 level;
+ unsigned 8 point1,point2,j,f,k;
+ unsigned 9 e,i;
+ signed 16 weight1,weight2;
+#if HARDWARE_MULTIPLY
+ signed 18 p,q,r,t;
+#else
+ signed 24 p,q,r,t;
+#endif
+ signed a,b;
+
+#if HARDWARE_MULTIPLY
+ // Macro to provide rescaling of 36-bit result of fixed point multiply
+ // down to an 18-bit result. The range of bits selected depends on the
+ // number that represents the value of "1" in the trig function lookup
+ // tables. (Eg. for 16384 == 1, the lowest bit selected should be [14]).
+ macro expr rescale (x) = (x[35] @ x[30:14]);
+#else
+ //Macro to rescale the multiply result down to a 24-bit value.
+ macro expr rescale (x) = ((x>>FRACBITS)<-24);
+#endif
+
+ for(level=1;level<=NUMBER_OF_COLUMNS;level++) // count all the columns
+ {
+ e=1<<(NUMBER_OF_COLUMNS-level+1); // number of points in each block in this column
+ f=(e>>1)<-8; // number of butterflies in each block in this column
+
+ for(j=1;j<=f;j++) // count all the butterflies in each block
+ {
+ par
+ {
+ // Weight factors for real (the same for FFT and iFFT)
+ weight1 = weight_re[((j-1)<<(level-1))<-7];
+
+
+ // Weight factors for imaginary (opposite for FFT and iFFT)
+ weight2 = (!select_inverse) ? (weight_im[((j-1)<<(level-1))<-7]) : -(weight_im[((j-1)<<(level-1))<-7]);
+
+ /* ORIGINAL CODE BELOW, MODIFIED BECAUSE OF MISMATCHING OUTPUT WITH BORLAND TESTAPP
+ weight2 = (!select_inverse) ? -(weight_im[((j-1)<<(level-1))<-7]) : weight_im[((j-1)<<(level-1))<-7];
+ */
+
+
+
+ for(i=0@j;i<=NUMBER_OF_POINTS;i+=e) // count all the blocks in this column
+ {
+ // Butterfly calculation
+ par
+ {
+ point1 = ((i<-8)-1);
+ point2 = (((i<-8)+f)-1);
+ }
+
+ par
+ {
+ p = (real.read[point1] >> 1) + (real.rwrite[point2] >> 1);
+ q = (imaginary.read[point1] >> 1) + (imaginary.rwrite[point2] >> 1);
+ }
+
+ par
+ {
+ r = (real.read[point1] >> 1) - (real.rwrite[point2] >> 1);
+ t = (imaginary.read[point1] >> 1) - (imaginary.rwrite[point2] >> 1);
+ }
+
+ multiply(a,r,weight1);
+ multiply(b,t,weight2);
+
+ par
+ {
+ real.rwrite[point2] = (rescale(a-b));
+ imaginary.rwrite[point1] = q;
+ }
+
+ multiply(a,t,weight1);
+ multiply(b,r,weight2);
+
+ par
+ {
+ real.rwrite[point1] = p;
+ imaginary.rwrite[point2] = (rescale(a+b));
+ }
+
+ }
+ }
+ }
+ }
+
+ j=1;
+ for(i=1;i<NUMBER_OF_POINTS;i++)
+ {
+ if(i<(0@j))
+ {
+ par
+ {
+ point1=j-1;
+ point2=(i-1)<-8;
+ }
+ /*
+ * COPYING ARRAY VALUES FROM ONE PLACE TO ANOTHER IN THE ARRAY MUST BE DONE IN
+ * 2 STEPS. FIRST THE VALUES HAVE TO BE COPIED TO SEPARATE VARIABLES AFTER THAT THEY
+ * ARE COPIED BACK TO THEIR NEW POSITION IN THE ARRAY. THIS MUST BE DONE TO
+ * PREVENT TIMING ISSUES FROM OCCURING.
+ */
+ /* Copy array values to separate variables */
+ par
+ {
+ p = real.read[point1];
+ q = imaginary.read[point1];
+ }
+ par
+ {
+ r = real.read[point2];
+ t = imaginary.read[point2];
+ }
+ /* Copy variables back to their new positions */
+ par
+ {
+ real.rwrite[point1] = r;
+ imaginary.rwrite[point1] = t;
+ }
+ par
+ {
+ real.rwrite[point2] = p;
+ imaginary.rwrite[point2] = q;
+ }
+ }
+
+ k = NUMBER_OF_POINTS>>1;
+
+
+ while(k<j)
+ {
+ j = j-k;
+ k = k>>1;
+ }
+
+ j+=k;
+ }
+ //Bitreversing end
+}
+
+/*! \fn void perform_fft(signed 18 *pcm_audio)
+ * \brief This routine obtains the audio data from the audio I/O component and copies this
+ * data to local arrays for calculating purposes, and calls the FFT algorithm.
+ *
+ * \param *pcm_audio pointer to array containg the audio data
+ *
+ * \return nothing
+ * \retval void
+ *
+ * cost 258 cycles (excl. the calculate FFT function)
+ */
+#if HARDWARE_MULTIPLY
+void perform_fft(signed 18 *pcm_audio)
+#else
+void perform_fft(signed 16 *pcm_audio)
+#endif
+{
+ unsigned 8 k;
+
+ //initialize variables before the copying pipeline
+#if HARDWARE_MULTIPLY
+ signed 18 sample;
+ k=0;
+ sample = adjs(pcm_audio[k],18);
+#else
+ signed 24 sample;
+ k=0;
+ sample = adjs(pcm_audio[k],24);
+#endif
+ // copy audio data to real-array before starting FFT calculation
+ // and set imaginary values to zero
+ do
+ {
+ //Copying the array values has been pipelined to prevent parallel access to the
+ //pcm_audio array. This copying procedure must be finished before another
+ //sample is read from the audio input. The time available for this loop is
+ //determined by the sampling rate of 44,1 Khz
+ par
+ {
+ //COPYING NEEDS TO BE DONE IN 2 STEPS, BECAUSE THE VALUE THAT NEEDS TO WRITTEN
+ //TO THE REAL-RAM NEEDS TO BE AVAILABLE ON THE START OFF THE CLOCKCYCLE.
+#if HARDWARE_MULTIPLY
+ sample = adjs(pcm_audio[k+1],18);
+#else
+ sample = adjs(pcm_audio[k+1],24);
+#endif
+ real.rwrite[k] = sample;
+ imaginary.rwrite[k] = 0;
+ k++;
+ }
+ } while (k);
+
+
+
+#if PERFORM_FFT_CALCULATION
+ calculate_fft(0);
+#endif
+
+
+}
+
+/*! \fn void perform_ifft(signed 18 *modified_audio, unsigned 6 *ifft_info)
+ * \brief This routine calls the ifft algorithm and after completing that it obtains the
+ * modified audio data and copies that to the output arrays of the audio I/O component.
+ * Besides that it also fills the array used by the display routine for displaying the waveform.
+ *
+ * \param *modified_audio pointer to array containg the audio data
+ * \param *ifft_info Pointer to the ifft_info array containing the modified waveform data for display purposes
+ *
+ * \return nothing
+ * \retval void
+ *
+ * cost 258 cycles (excl. the calculate iFFT function)
+ */
+#if HARDWARE_MULTIPLY
+void perform_ifft(signed 18 *modified_audio, unsigned 6 *ifft_info)
+#else
+void perform_ifft(signed 16 *modified_audio, unsigned 6 *ifft_info)
+#endif
+{
+ unsigned 6 k;
+#if HARDWARE_MULTIPLY
+ signed 18 p;
+#else
+ signed 24 p;
+#endif
+#if PERFORM_FFT_CALCULATION
+ calculate_fft(1);
+#endif
+
+ k=0;
+//initialize variables for the copying pipeline
+#if PERFORM_FFT_CALCULATION
+ #if HARDWARE_MULTIPLY
+ p = (real.read[(0@k)+95] << NUMBER_OF_COLUMNS);
+ #else
+ p = (real.read[(0@k)+95] >> NUMBER_OF_COLUMNS);
+ #endif
+#else
+ p = (real.read[(0@k)+95]);
+#endif
+
+ do
+ {
+ //Copying the array values has been pipelined to prevent parallel access to the
+ //pcm_audio array. This copying procedure must be finished before another
+ //sample is read from the audio input. The time available for this loop is
+ //determined by the sampling rate of 44,1 Khz
+ par
+ {
+ /*
+ * Before copying the modified audio from the local real-array
+ * to the output array of the audio I/O component, compensate
+ * for the FFT calculation by shifting the values.
+ * 95 is added to start the output from the middle of the sliding
+ * window, this is done to get a better sound quality.
+ */
+#if PERFORM_FFT_CALCULATION
+ #if HARDWARE_MULTIPLY
+ p = (real.read[(0@k)+95] << NUMBER_OF_COLUMNS);
+ #else
+ p = (real.read[(0@k)+95] >> NUMBER_OF_COLUMNS);
+ #endif
+#else
+ p = (real.read[(0@k)+95]);
+#endif
+ //Copy the modified audio from the local real array to the output array of the audio I/O component.
+#if HARDWARE_MULTIPLY
+ modified_audio[k] = p ;
+#else
+ modified_audio[k] = (p<-16);
+#endif
+ //Fill the array for displaying the waveform, only the 6 MSB are needed.
+ //Because the display routine needs unsigned values a cast is done here, this should be fixed
+ //in such a way a cast is not necessary anymore.
+ ifft_info[k] = (unsigned 6)(32+(p[17:12]));
+ k++;
+ }
+ } while(k);
+}
+
+/*! \fn void equalize_audio(audiodata_t *audiodata)
+ * \brief This routine equalizes the frequencies derived by the FFT calculation,
+ * according to the settings of the equalizer bars.
+ *
+ * \note Cost: 3844 clock cycles (Maximum)
+ *
+ * \param *audiodata pointer to the audiodata struct, containing the eq_info, etc.
+ *
+ * \return void
+ * \retval void
+ *
+ */
+void equalize_audio(audiodata_t *audiodata)
+{
+#if HARDWARE_MULTIPLY
+ signed 18 p,q;
+#else
+ signed 24 p,q;
+#endif
+ signed 18 a;
+ unsigned 8 i, mirror_i, bit, m, n;
+ unsigned 7 old_value;
+ unsigned 9 tmp;
+
+ //macro expr equalize_bar = multiply(q,a)[29:6];
+
+ macro proc equalize_bar(retval)
+ {
+ signed result;
+ multiply(result, q,a);
+#if HARDWARE_MULTIPLY
+ retval = result[23:6]; //drop last 6 bit to compensate the maximum multiplication with 64 from the eq_settings array
+#else
+ retval = result[29:6]; //drop last 6 bit to compensate the maximum multiplication with 64 from the eq_settings array
+#endif
+ }
+
+ p = real.read[0] - DC_COMPONENT; // remove DC component for calculations
+ real.rwrite[0] = p;
+
+ for(i=0;i!=NUMBER_OF_FREQUENCIES;i++)
+ {
+
+ // set multiplication factor (0..64) for current frequency bar, The first frequency band must be equalized at 100% (63) since there is no DC-component taken into account.
+ a = adjs(eq_settings[audiodata->equalizer_levels_ptr[i <- 7]],18);
+
+
+ // multiply frequency with this factor and divide by 64 (drop 6 LSB's)
+ q = real.read[i];
+ equalize_bar(p);
+ real.rwrite[i] = p;
+
+ q = imaginary.read[i];
+ equalize_bar(p);
+ imaginary.rwrite[i] = p;
+
+ // the upper part(128..255) of the spectrum is mirrored to the lower part;
+ // these values need to be adjusted too
+ if ((i<-7)!=0) // if not in DC component bar
+ {
+ mirror_i = (NUMBER_OF_POINTS-1)-i+1;
+ q = real.read[mirror_i];
+ equalize_bar(p);
+ real.rwrite[mirror_i] = p;
+
+ q = imaginary.read[mirror_i];
+ equalize_bar(p);
+ imaginary.rwrite[mirror_i] = p;
+ }
+ }
+
+ //write data to fft_info for display purposes
+ for(i=0;i<NUMBER_OF_FREQUENCIES;i++)
+ {
+ p = real.read[i];
+ q = imaginary.read[i];
+#if HARDWARE_MULTIPLY
+ if (p[17] == 1) p = -p; else delay;
+ if (q[17] == 1) q = -q; else delay;
+#else
+ if (p[23] == 1) p = -p; else delay;
+ if (q[23] == 1) q = -q; else delay;
+#endif
+ p = (p<q) ? q : p; // This is done to get the best visual frequency result
+
+ if (audiodata->display_log)
+ {
+
+ bit=126;
+#if HARDWARE_MULTIPLY
+ while ((p[15] == 0) && (bit != 0))
+#else
+ while ((p[21] == 0) && (bit != 0))
+#endif
+ par
+ {
+ p = p<<1;
+ bit = bit - 18;
+ }
+ old_value = audiodata->fft_info.write[0 @ (i <- 7)];
+ tmp = ((0@old_value) + (0@bit))>>1;
+
+ //Write data to FFT_info array to be displayed as a sound spectrum
+ audiodata->fft_info.write[0 @ (i <- 7)] = (old_value <= (tmp<-7)) ? (tmp<-7) : old_value-1;
+ }
+ else
+ {
+ //Write data to FFT_info array to be displayed as a sound spectrum
+
+ old_value = audiodata->fft_info.write[0 @ (i <- 7)];
+#if HARDWARE_MULTIPLY
+ audiodata->fft_info.write[0 @ (i <- 7)] = (old_value<=(unsigned)(p[15:9])) ? (unsigned)(p[15:9]) : old_value-1;
+#else
+ audiodata->fft_info.write[0 @ (i <- 7)] = (old_value<=(unsigned)(p[21:15])) ? (unsigned)(p[21:15]) : old_value-1;
+#endif
+ }
+ }
+
+ // add DC component again before inverse FFT calculation is performed
+
+ p = real.read[0] + DC_COMPONENT;
+ real.rwrite[0] = p;
+}
+
|