 |
 |
| View previous topic :: View next topic |
| Author |
Message |
Marco27293
Joined: 09 May 2020
Posts: 126
|
| PIC18F47J53 FFT optimized code PCH CCS C v5.090 compiler |
 Posted: Wed Mar 15, 2023 5:27 am |
 |
|
Hi,
I'm working with a pic18f47j53, using ccs c pch compiler.
I would like to perform fft of a 128 long audio signal (3khz Fs, each sample is signed int16).
How can I customize this open-source code, found online ?
Do you have any suggestions in order to do this process with very low fixed points computation ?
I tried replacing the short variable declarations in fft.c and fft.h with signed int16 and the long int variable declarations in fft.c and fft.h with signed int32 but it seems not to work.
Unfortunately bit shift division/multiplication seems to not work with signed variables.
I would very appreciate your support.
Regards,
Marco
| Code: |
/************************************************************************
fft.c
FFT Audio Analysis
Copyright (C) 2011 Simon Inns
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 3 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, see <http://www.gnu.org/licenses/>.
Email: simon.inns@gmail.com
************************************************************************/
#ifndef FFT_C
#define FFT_C
#include <htc.h>
#include "fft.h"
// fix_fft.c - Fixed-point in-place Fast Fourier Transform
// All data are fixed-point short integers, in which -32768
// to +32768 represent -1.0 to +1.0 respectively. Integer
// arithmetic is used for speed, instead of the more natural
// floating-point.
//
// For the forward FFT (time -> freq), fixed scaling is
// performed to prevent arithmetic overflow, and to map a 0dB
// sine/cosine wave (i.e. amplitude = 32767) to two -6dB freq
// coefficients.
//
// Written by: Tom Roberts 11/8/89
// Made portable: Malcolm Slaney 12/15/94 malcolm@interval.com
// Enhanced: Dimitrios P. Bouras 14 Jun 2006 dbouras@ieee.org
// Ported to PIC18F: Simon Inns 20110104
/*
fix_fft() - perform forward fast Fourier transform.
fr[n],fi[n] are real and imaginary arrays, both INPUT AND
RESULT (in-place FFT), with 0 <= n < 2**m
*/
void fix_fft(short fr[], short fi[], short m)
{
long int mr = 0, nn, i, j, l, k, istep, n, shift;
short qr, qi, tr, ti, wr, wi;
n = 1 << m;
nn = n - 1;
/* max FFT size = N_WAVE */
//if (n > N_WAVE) return -1;
/* decimation in time - re-order data */
for (m=1; m<=nn; ++m)
{
l = n;
do
{
l >>= 1;
} while (mr+l > nn);
mr = (mr & (l-1)) + l;
if (mr <= m) continue;
tr = fr[m];
fr[m] = fr[mr];
fr[mr] = tr;
ti = fi[m];
fi[m] = fi[mr];
fi[mr] = ti;
}
l = 1;
k = LOG2_N_WAVE-1;
while (l < n)
{
/*
fixed scaling, for proper normalization --
there will be log2(n) passes, so this results
in an overall factor of 1/n, distributed to
maximize arithmetic accuracy.
It may not be obvious, but the shift will be
performed on each data point exactly once,
during this pass.
*/
// Variables for multiplication code
long int c;
short b;
istep = l << 1;
for (m=0; m<l; ++m)
{
j = m << k;
/* 0 <= j < N_WAVE/2 */
wr = Sinewave[j+N_WAVE/4];
wi = -Sinewave[j];
wr >>= 1;
wi >>= 1;
for (i=m; i<n; i+=istep)
{
j = i + l;
// Here I unrolled the multiplications to prevent overhead
// for procedural calls (we don't need to be clever about
// the actual multiplications since the pic has an onboard
// 8x8 multiplier in the ALU):
// tr = FIX_MPY(wr,fr[j]) - FIX_MPY(wi,fi[j]);
c = ((long int)wr * (long int)fr[j]);
c = c >> 14;
b = c & 0x01;
tr = (c >> 1) + b;
c = ((long int)wi * (long int)fi[j]);
c = c >> 14;
b = c & 0x01;
tr = tr - ((c >> 1) + b);
// ti = FIX_MPY(wr,fi[j]) + FIX_MPY(wi,fr[j]);
c = ((long int)wr * (long int)fi[j]);
c = c >> 14;
b = c & 0x01;
ti = (c >> 1) + b;
c = ((long int)wi * (long int)fr[j]);
c = c >> 14;
b = c & 0x01;
ti = ti + ((c >> 1) + b);
qr = fr[i];
qi = fi[i];
qr >>= 1;
qi >>= 1;
fr[j] = qr - tr;
fi[j] = qi - ti;
fr[i] = qr + tr;
fi[i] = qi + ti;
}
}
--k;
l = istep;
}
}
#endif
|
| Code: |
/************************************************************************
fft.h
FFT Audio Analysis
Copyright (C) 2011 Simon Inns
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 3 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, see <http://www.gnu.org/licenses/>.
Email: simon.inns@gmail.com
************************************************************************/
#ifndef _FFT_H
#define _FFT_H
// Definitions
#define N_WAVE 1024 // full length of Sinewave[]
#define LOG2_N_WAVE 10 // log2(N_WAVE)
// Since we only use 3/4 of N_WAVE, we define only
// this many samples, in order to conserve data space.
const short Sinewave[N_WAVE-N_WAVE/4] = {
0, 201, 402, 603, 804, 1005, 1206, 1406,
1607, 1808, 2009, 2209, 2410, 2610, 2811, 3011,
3211, 3411, 3611, 3811, 4011, 4210, 4409, 4608,
4807, 5006, 5205, 5403, 5601, 5799, 5997, 6195,
6392, 6589, 6786, 6982, 7179, 7375, 7571, 7766,
7961, 8156, 8351, 8545, 8739, 8932, 9126, 9319,
9511, 9703, 9895, 10087, 10278, 10469, 10659, 10849,
11038, 11227, 11416, 11604, 11792, 11980, 12166, 12353,
12539, 12724, 12909, 13094, 13278, 13462, 13645, 13827,
14009, 14191, 14372, 14552, 14732, 14911, 15090, 15268,
15446, 15623, 15799, 15975, 16150, 16325, 16499, 16672,
16845, 17017, 17189, 17360, 17530, 17699, 17868, 18036,
18204, 18371, 18537, 18702, 18867, 19031, 19194, 19357,
19519, 19680, 19840, 20000, 20159, 20317, 20474, 20631,
20787, 20942, 21096, 21249, 21402, 21554, 21705, 21855,
22004, 22153, 22301, 22448, 22594, 22739, 22883, 23027,
23169, 23311, 23452, 23592, 23731, 23869, 24006, 24143,
24278, 24413, 24546, 24679, 24811, 24942, 25072, 25201,
25329, 25456, 25582, 25707, 25831, 25954, 26077, 26198,
26318, 26437, 26556, 26673, 26789, 26905, 27019, 27132,
27244, 27355, 27466, 27575, 27683, 27790, 27896, 28001,
28105, 28208, 28309, 28410, 28510, 28608, 28706, 28802,
28897, 28992, 29085, 29177, 29268, 29358, 29446, 29534,
29621, 29706, 29790, 29873, 29955, 30036, 30116, 30195,
30272, 30349, 30424, 30498, 30571, 30643, 30713, 30783,
30851, 30918, 30984, 31049, 31113, 31175, 31236, 31297,
31356, 31413, 31470, 31525, 31580, 31633, 31684, 31735,
31785, 31833, 31880, 31926, 31970, 32014, 32056, 32097,
32137, 32176, 32213, 32249, 32284, 32318, 32350, 32382,
32412, 32441, 32468, 32495, 32520, 32544, 32567, 32588,
32609, 32628, 32646, 32662, 32678, 32692, 32705, 32717,
32727, 32736, 32744, 32751, 32757, 32761, 32764, 32766,
32767, 32766, 32764, 32761, 32757, 32751, 32744, 32736,
32727, 32717, 32705, 32692, 32678, 32662, 32646, 32628,
32609, 32588, 32567, 32544, 32520, 32495, 32468, 32441,
32412, 32382, 32350, 32318, 32284, 32249, 32213, 32176,
32137, 32097, 32056, 32014, 31970, 31926, 31880, 31833,
31785, 31735, 31684, 31633, 31580, 31525, 31470, 31413,
31356, 31297, 31236, 31175, 31113, 31049, 30984, 30918,
30851, 30783, 30713, 30643, 30571, 30498, 30424, 30349,
30272, 30195, 30116, 30036, 29955, 29873, 29790, 29706,
29621, 29534, 29446, 29358, 29268, 29177, 29085, 28992,
28897, 28802, 28706, 28608, 28510, 28410, 28309, 28208,
28105, 28001, 27896, 27790, 27683, 27575, 27466, 27355,
27244, 27132, 27019, 26905, 26789, 26673, 26556, 26437,
26318, 26198, 26077, 25954, 25831, 25707, 25582, 25456,
25329, 25201, 25072, 24942, 24811, 24679, 24546, 24413,
24278, 24143, 24006, 23869, 23731, 23592, 23452, 23311,
23169, 23027, 22883, 22739, 22594, 22448, 22301, 22153,
22004, 21855, 21705, 21554, 21402, 21249, 21096, 20942,
20787, 20631, 20474, 20317, 20159, 20000, 19840, 19680,
19519, 19357, 19194, 19031, 18867, 18702, 18537, 18371,
18204, 18036, 17868, 17699, 17530, 17360, 17189, 17017,
16845, 16672, 16499, 16325, 16150, 15975, 15799, 15623,
15446, 15268, 15090, 14911, 14732, 14552, 14372, 14191,
14009, 13827, 13645, 13462, 13278, 13094, 12909, 12724,
12539, 12353, 12166, 11980, 11792, 11604, 11416, 11227,
11038, 10849, 10659, 10469, 10278, 10087, 9895, 9703,
9511, 9319, 9126, 8932, 8739, 8545, 8351, 8156,
7961, 7766, 7571, 7375, 7179, 6982, 6786, 6589,
6392, 6195, 5997, 5799, 5601, 5403, 5205, 5006,
4807, 4608, 4409, 4210, 4011, 3811, 3611, 3411,
3211, 3011, 2811, 2610, 2410, 2209, 2009, 1808,
1607, 1406, 1206, 1005, 804, 603, 402, 201,
0, -201, -402, -603, -804, -1005, -1206, -1406,
-1607, -1808, -2009, -2209, -2410, -2610, -2811, -3011,
-3211, -3411, -3611, -3811, -4011, -4210, -4409, -4608,
-4807, -5006, -5205, -5403, -5601, -5799, -5997, -6195,
-6392, -6589, -6786, -6982, -7179, -7375, -7571, -7766,
-7961, -8156, -8351, -8545, -8739, -8932, -9126, -9319,
-9511, -9703, -9895, -10087, -10278, -10469, -10659, -10849,
-11038, -11227, -11416, -11604, -11792, -11980, -12166, -12353,
-12539, -12724, -12909, -13094, -13278, -13462, -13645, -13827,
-14009, -14191, -14372, -14552, -14732, -14911, -15090, -15268,
-15446, -15623, -15799, -15975, -16150, -16325, -16499, -16672,
-16845, -17017, -17189, -17360, -17530, -17699, -17868, -18036,
-18204, -18371, -18537, -18702, -18867, -19031, -19194, -19357,
-19519, -19680, -19840, -20000, -20159, -20317, -20474, -20631,
-20787, -20942, -21096, -21249, -21402, -21554, -21705, -21855,
-22004, -22153, -22301, -22448, -22594, -22739, -22883, -23027,
-23169, -23311, -23452, -23592, -23731, -23869, -24006, -24143,
-24278, -24413, -24546, -24679, -24811, -24942, -25072, -25201,
-25329, -25456, -25582, -25707, -25831, -25954, -26077, -26198,
-26318, -26437, -26556, -26673, -26789, -26905, -27019, -27132,
-27244, -27355, -27466, -27575, -27683, -27790, -27896, -28001,
-28105, -28208, -28309, -28410, -28510, -28608, -28706, -28802,
-28897, -28992, -29085, -29177, -29268, -29358, -29446, -29534,
-29621, -29706, -29790, -29873, -29955, -30036, -30116, -30195,
-30272, -30349, -30424, -30498, -30571, -30643, -30713, -30783,
-30851, -30918, -30984, -31049, -31113, -31175, -31236, -31297,
-31356, -31413, -31470, -31525, -31580, -31633, -31684, -31735,
-31785, -31833, -31880, -31926, -31970, -32014, -32056, -32097,
-32137, -32176, -32213, -32249, -32284, -32318, -32350, -32382,
-32412, -32441, -32468, -32495, -32520, -32544, -32567, -32588,
-32609, -32628, -32646, -32662, -32678, -32692, -32705, -32717,
-32727, -32736, -32744, -32751, -32757, -32761, -32764, -32766,
};
// Function prototypes
void fix_fft(short fr[], short fi[], short m);
#endif
|
|
|
 |
Ttelmah
Joined: 11 Mar 2010
Posts: 19539
|
|
| Posted: Wed Mar 15, 2023 7:39 am |
|
|
Change every reference to 'short', to 'signed int16', and every 'long int'
to 'signed int32'. That is all you should need to do.
I hate people using 'short' and 'long', since these historically differ on
just about every compiler in existence. However the use above is the
defined size for an ANSI compiler.
Everything else is standard C and should run as is. |
|
|
Marco27293
Joined: 09 May 2020
Posts: 126
|
|
| Posted: Wed Mar 15, 2023 7:44 am |
|
|
I tried but it seems to not work.
Maybe signed variables do not work with bit shift operations. |
|
|
Ttelmah
Joined: 11 Mar 2010
Posts: 19539
|
|
| Posted: Wed Mar 15, 2023 8:02 am |
|
|
Should do. However there are differences in 'how' the sign is handled.
On negative numbers this is 'implementation dependant'...
What compiler version are you on?. CCS changed their handling of
shifts on signed values a while ago, to make it more compatible.
They do now sign extend, which is what I'd expect the code to want.
Are you sure you changed all the casts correctly?. The shifts need to
be done on signed int32 values or they won't work correctly. So not just
the variable declarations, but these as well.
Last edited by Ttelmah on Wed Mar 15, 2023 8:03 am; edited 1 time in total |
|
|
Marco27293
Joined: 09 May 2020
Posts: 126
|
|
| Posted: Wed Mar 15, 2023 8:03 am |
|
|
| I use pch ccs c v5.090 compiler |
|
|
Ttelmah
Joined: 11 Mar 2010
Posts: 19539
|
|
| Posted: Wed Mar 15, 2023 8:15 am |
|
|
OK. 5.090 doesn't sign extend on a shift right. Instead of the >> operators
use:
| Code: |
#define arithright_shift(val, amount) (val<0) ? -((signed int32)(-val)>>amount):(signed int32)val>>amount
|
|
|
|
Marco27293
Joined: 09 May 2020
Posts: 126
|
|
| Posted: Wed Mar 15, 2023 8:17 am |
|
|
Thanks!!
What about left shift ? |
|
|
Ttelmah
Joined: 11 Mar 2010
Posts: 19539
|
|
| Posted: Wed Mar 15, 2023 8:32 am |
|
|
Given your values have been cast to signed int32, all the left bits will
be filled with ones, so there will be a one at the top after the left shift,
so it'll still be negative. It is right shift where the result is undefined.
Even ANSI, says that right shifting a signed value gives an undefined
result!... |
|
|
|
|
You cannot post new topics in this forum
You cannot reply to topics in this forum
You cannot edit your posts in this forum
You cannot delete your posts in this forum
You cannot vote in polls in this forum
|
Powered by phpBB © 2001, 2005 phpBB Group
|
FAQ 
Search
Register
Profile 
Log in to check your private messages 
Log in
