CCS C Software and Maintenance Offers
FAQFAQ   FAQForum Help   FAQOfficial CCS Support   SearchSearch  RegisterRegister 

ProfileProfile   Log in to check your private messagesLog in to check your private messages   Log inLog in 

CCS does not monitor this forum on a regular basis.

Please do not post bug reports on this forum. Send them to support@ccsinfo.com

PIC18F47J53 FFT optimized code PCH CCS C v5.090 compiler

 
Post new topic   Reply to topic    CCS Forum Index -> General CCS C Discussion
View previous topic :: View next topic  
Author Message
Marco27293



Joined: 09 May 2020
Posts: 110

View user's profile Send private message

PIC18F47J53 FFT optimized code PCH CCS C v5.090 compiler
PostPosted: Wed Mar 15, 2023 5:27 am     Reply with quote

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: 19222

View user's profile Send private message

PostPosted: Wed Mar 15, 2023 7:39 am     Reply with quote

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: 110

View user's profile Send private message

PostPosted: Wed Mar 15, 2023 7:44 am     Reply with quote

I tried but it seems to not work.

Maybe signed variables do not work with bit shift operations.
Ttelmah



Joined: 11 Mar 2010
Posts: 19222

View user's profile Send private message

PostPosted: Wed Mar 15, 2023 8:02 am     Reply with quote

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: 110

View user's profile Send private message

PostPosted: Wed Mar 15, 2023 8:03 am     Reply with quote

I use pch ccs c v5.090 compiler
Ttelmah



Joined: 11 Mar 2010
Posts: 19222

View user's profile Send private message

PostPosted: Wed Mar 15, 2023 8:15 am     Reply with quote

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: 110

View user's profile Send private message

PostPosted: Wed Mar 15, 2023 8:17 am     Reply with quote

Thanks!!

What about left shift ?
Ttelmah



Joined: 11 Mar 2010
Posts: 19222

View user's profile Send private message

PostPosted: Wed Mar 15, 2023 8:32 am     Reply with quote

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!...
Display posts from previous:   
Post new topic   Reply to topic    CCS Forum Index -> General CCS C Discussion All times are GMT - 6 Hours
Page 1 of 1

 
Jump to:  
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