diff options
Diffstat (limited to 'tests/test6')
| -rw-r--r-- | tests/test6/Makefile | 30 | ||||
| -rwxr-xr-x | tests/test6/compil.sh | 2 | ||||
| -rw-r--r-- | tests/test6/compute.c | 138 | ||||
| -rw-r--r-- | tests/test6/compute.h | 11 | ||||
| -rw-r--r-- | tests/test6/fft.c | 241 | ||||
| -rw-r--r-- | tests/test6/fft.h | 19 | ||||
| l--------- | tests/test6/test.raw | 1 | ||||
| -rw-r--r-- | tests/test6/test6.c | 19 |
8 files changed, 456 insertions, 5 deletions
diff --git a/tests/test6/Makefile b/tests/test6/Makefile new file mode 100644 index 0000000..560ff03 --- /dev/null +++ b/tests/test6/Makefile @@ -0,0 +1,30 @@ +CC=gcc +CFLAGS=-W -Wall -Werror -Wno-error=unused-parameter -g +LDFLAGS=-Werror -g +EXEC=test6 + +#CFLAGS+=$(shell pkg-config --cflags gtk+-2.0 gthread-2.0 libpulse) +#LDFLAGS+=$(shell pkg-config --libs gtk+-2.0 gthread-2.0 libpulse) +LDFLAGS+=-lm + +SRC= $(wildcard *.c) +OBJ= $(SRC:.c=.o) + +all: $(EXEC) + +$(EXEC): $(OBJ) + $(CC) -o $@ $^ $(LDFLAGS) + +#main.o: hello.h + +%.o: %.c + $(CC) -o $@ -c $< $(CFLAGS) + +.PHONY: clean mrproper + +clean: + rm *.o + +mrproper: clean + rm $(EXEC) + diff --git a/tests/test6/compil.sh b/tests/test6/compil.sh deleted file mode 100755 index 2e49576..0000000 --- a/tests/test6/compil.sh +++ /dev/null @@ -1,2 +0,0 @@ -#!/bin/bash -xe -gcc -Wall -g -o test6 test6.c -lm diff --git a/tests/test6/compute.c b/tests/test6/compute.c new file mode 100644 index 0000000..017416b --- /dev/null +++ b/tests/test6/compute.c @@ -0,0 +1,138 @@ +#include "compute.h" + +#include "fft.h" +#include <math.h> + +#define MIN_SAMPLES 256 +#define MAX_SAMPLES 2048 + +//#define MAX(a,b) (a>b?a:b) +//#define MIN(a,b) (a<b?a:b) + +//static inline float todB_a(const float *x); +void compute_spectrum(float *data, int width, float output[PSHalf]); + + +float compute_level(const float *data, size_t nsamples, int rate) { + + size_t i; + float input[MAX_SAMPLES], pwrspec[PSHalf]; + float value; + int f, min_f_index, max_f_index; + + if (nsamples >= MAX_SAMPLES) { + printf("WARN : nsamples >= MAX_SAMPLES : %i >= %i\n", nsamples, MAX_SAMPLES); + nsamples=MAX_SAMPLES; + } + if (nsamples < MIN_SAMPLES) { + printf("WARN : nsamples < MIN_SAMPLES : %i >= %i\n", nsamples, MIN_SAMPLES); + // Replicate with symmetry the sound to obtain an input buffer of the minimal len + for (i=0;i<MIN_SAMPLES;i++) { + if ( (i/nsamples)%2==1 ) + input[i]=data[i]; // First channel only + else + input[i]=data[nsamples-i-1]; + } + nsamples=MIN_SAMPLES; + } else { + for (i=0;i<nsamples;i++) { + input[i]=data[i]; // First channel only + } + } + + compute_spectrum(input, nsamples, pwrspec); + + // Compute the mean power for 200Hz to 2000Hz band + min_f_index=((float)PSHalf)*200.f/(((float)rate)/2.f); + max_f_index=((float)PSHalf)*2000.f/(((float)rate)/2.f); + + value=0.f; + for (f=min_f_index;f<=max_f_index;f++) { + value+=pwrspec[f]; + } + // Mean value + value=value/(max_f_index-min_f_index+1); + + return value; +} +/* +static inline float todB_a(const float *x){ + return (float)((*(int32_t *)x)&0x7fffffff) * 7.17711438e-7f -764.6161886f; +} +*/ +// Adapted from Audacity +void compute_spectrum(float *data, int width, float output[PSHalf]) { + + int i, start, windows; + float temp; + float in[PSNumS]; + float out[PSHalf]; + float processed[PSHalf]={0.0f}; + + start = 0; + windows = 0; + while (start + PSNumS <= width) { + // Windowing : Hanning + for (i=0; i<PSNumS; i++) + in[i] = data[start+i] *(0.50-0.50*cos(2*M_PI*i/(PSNumS-1))); + + // Returns only the real part of the result + PowerSpectrum(in, out); + + for (i=0; i<PSHalf; i++) + processed[i] += out[i]; + + start += PSHalf; + windows++; + } + // Convert to decibels + // But do it safely; -Inf is nobody's friend + for (i = 0; i < PSHalf; i++){ + temp=(processed[i] / PSNumS / windows); + if (temp > 0.0) + output[i] = 10*log10(temp); + else + output[i] = 0; + } +} + +void audio2hsv_1(int audio_level, int *light_h, int *light_s, int *light_v) { + // Dummy code + *light_h=-audio_level; + *light_s=audio_level; + *light_v=65535; +} + + +void hsv2rgb(int h, int s, int v, int *r, int *g, int *b) { + /* + * Purpose: + * Convert HSV values to RGB values + * All values are in the range [0..65535] + */ + float F, M, N, K; + int I; + + if ( s == 0 ) { + /* + * Achromatic case, set level of grey + */ + *r = v; + *g = v; + *b = v; + } else { + I = (int) h/(65535/6); /* should be in the range 0..5 */ + F = h - I; /* fractional part */ + + M = v * (1 - s); + N = v * (1 - s * F); + K = v * (1 - s * (1 - F)); + + if (I == 0) { *r = v; *g = K; *b = M; } + if (I == 1) { *r = N; *g = v; *b = M; } + if (I == 2) { *r = M; *g = v; *b = K; } + if (I == 3) { *r = M; *g = N; *b = v; } + if (I == 4) { *r = K; *g = M; *b = v; } + if (I == 5) { *r = v; *g = M; *b = N; } + } +} diff --git a/tests/test6/compute.h b/tests/test6/compute.h new file mode 100644 index 0000000..d663ba5 --- /dev/null +++ b/tests/test6/compute.h @@ -0,0 +1,11 @@ +#ifndef COMPUTE_H +#define COMPUTE_H +#include <stdlib.h> // for size_t +#include <stdio.h> // for printf + +float compute_level(const float *data, size_t nsamples, int rate); +void audio2hsv_1(int audio_level, int *light_h, int *light_s, int *light_v); +void hsv2rgb(int h, int s, int v, int *r, int *g, int *b); + +#endif + diff --git a/tests/test6/fft.c b/tests/test6/fft.c new file mode 100644 index 0000000..dd91d14 --- /dev/null +++ b/tests/test6/fft.c @@ -0,0 +1,241 @@ +#include "fft.h" +#include <math.h> +#include <stdlib.h> + +/* Everything here comes from Audacity 1.3.13 + (orignally in C++ and with more genericity and functionnality) + Original Author : Dominic Mazzoni + Licenced under GPL 2.0 (see LICENCE) +*/ + +#define MaxFastBits 16 + +int *gFFTBitTable[MaxFastBits]={NULL}; + +void InitFFT(); +int NumberOfBitsNeeded(int PowerOfTwo); +inline int FastReverseBits(int i, int NumBits); +int ReverseBits(int index, int NumBits); + + +void FFT(int NumSamples, + int InverseTransform, + float *RealIn, float *ImagIn, float *RealOut, float *ImagOut) +{ + int NumBits; /* Number of bits needed to store indices */ + int i, j, k, n; + int BlockSize, BlockEnd; + + double angle_numerator = 2.0 * M_PI; + double tr, ti; /* temp real, temp imaginary */ +/* + if (!IsPowerOfTwo(NumSamples)) { + fprintf(stderr, "%d is not a power of two\n", NumSamples); + exit(1); + } +*/ + if (!gFFTBitTable[0]) + InitFFT(); + + if (!InverseTransform) + angle_numerator = -angle_numerator; + + NumBits = NumberOfBitsNeeded(NumSamples); + + /* + ** Do simultaneous data copy and bit-reversal ordering into outputs... + */ + for (i = 0; i < NumSamples; i++) { + j = FastReverseBits(i, NumBits); + RealOut[j] = RealIn[i]; + ImagOut[j] = (ImagIn == NULL) ? 0.0 : ImagIn[i]; + } + + /* + ** Do the FFT itself... + */ + + BlockEnd = 1; + for (BlockSize = 2; BlockSize <= NumSamples; BlockSize <<= 1) { + + double delta_angle = angle_numerator / (double) BlockSize; + + double sm2 = sin(-2 * delta_angle); + double sm1 = sin(-delta_angle); + double cm2 = cos(-2 * delta_angle); + double cm1 = cos(-delta_angle); + double w = 2 * cm1; + double ar0, ar1, ar2, ai0, ai1, ai2; + + for (i = 0; i < NumSamples; i += BlockSize) { + ar2 = cm2; + ar1 = cm1; + + ai2 = sm2; + ai1 = sm1; + + for (j = i, n = 0; n < BlockEnd; j++, n++) { + ar0 = w * ar1 - ar2; + ar2 = ar1; + ar1 = ar0; + + ai0 = w * ai1 - ai2; + ai2 = ai1; + ai1 = ai0; + + k = j + BlockEnd; + tr = ar0 * RealOut[k] - ai0 * ImagOut[k]; + ti = ar0 * ImagOut[k] + ai0 * RealOut[k]; + + RealOut[k] = RealOut[j] - tr; + ImagOut[k] = ImagOut[j] - ti; + + RealOut[j] += tr; + ImagOut[j] += ti; + } + } + BlockEnd = BlockSize; + } + + /* + ** Need to normalize if inverse transform... + */ + + if (InverseTransform) { + float denom = (float) NumSamples; + + for (i = 0; i < NumSamples; i++) { + RealOut[i] /= denom; + ImagOut[i] /= denom; + } + } +} + +void InitFFT() { + int len=2; + int b, i; + for (b=0; b<MaxFastBits; b++) { + gFFTBitTable[b]=malloc(len*sizeof(int)); + + for (i=0; i<len; i++) + gFFTBitTable[b][i] = ReverseBits(i, b+1); + + len <<= 1; + } +} + +int NumberOfBitsNeeded(int PowerOfTwo) +{ + int i; + +/* + if (PowerOfTwo < 2) { + fprintf(stderr, "Error: FFT called with size %d\n", PowerOfTwo); + exit(1); + } +*/ + for (i = 0;; i++) + if (PowerOfTwo & (1 << i)) + return i; +} + +inline int FastReverseBits(int i, int NumBits) +{ + if (NumBits <= MaxFastBits) + return gFFTBitTable[NumBits - 1][i]; + else + return ReverseBits(i, NumBits); +} + +int ReverseBits(int index, int NumBits) +{ + int i, rev; + + for (i = rev = 0; i < NumBits; i++) { + rev = (rev << 1) | (index & 1); + index >>= 1; + } + + return rev; +} + + +/* + * PowerSpectrum + * + * This function computes the same as RealFFT, above, but + * adds the squares of the real and imaginary part of each + * coefficient, extracting the power and throwing away the + * phase. + * + * For speed, it does not call RealFFT, but duplicates some + * of its code. + */ + +void PowerSpectrum(float In[PSNumS], float Out[PSHalf]) { + int i; + + float theta = M_PI / PSHalf; + + float tmpReal[PSHalf]; + float tmpImag[PSHalf]; + float RealOut[PSHalf]; + float ImagOut[PSHalf]; + + for (i=0; i<PSHalf; i++) { + tmpReal[i] = In[2*i]; + tmpImag[i] = In[2*i+1]; + } + + FFT(PSHalf, 0, tmpReal, tmpImag, RealOut, ImagOut); + + float wtemp = sin(0.5 * theta); + + float wpr = -2.0 * wtemp * wtemp; + float wpi = -1.0 * sin(theta); + float wr = 1.0 + wpr; + float wi = wpi; + + int i3; + + float h1r, h1i, h2r, h2i, rt, it; + for (i=1; i < PSHalf/2; i++) { + + i3 = PSHalf - i; + + h1r = 0.5 * (RealOut[i] + RealOut[i3]); + h1i = 0.5 * (ImagOut[i] - ImagOut[i3]); + h2r = 0.5 * (ImagOut[i] + ImagOut[i3]); + h2i = -0.5 * (RealOut[i] - RealOut[i3]); + + rt = h1r + wr * h2r - wi * h2i; + it = h1i + wr * h2i + wi * h2r; + + Out[i] = rt * rt + it * it; + + rt = h1r - wr * h2r + wi * h2i; + it = -h1i + wr * h2i + wi * h2r; + + Out[i3] = rt * rt + it * it; + + wr = (wtemp = wr) * wpr - wi * wpi + wr; + wi = wi * wpr + wtemp * wpi + wi; + } + + rt = (h1r = RealOut[0]) + ImagOut[0]; + it = h1r - ImagOut[0]; + Out[0] = rt * rt + it * it; + rt = RealOut[PSHalf / 2]; + it = ImagOut[PSHalf / 2]; + Out[PSHalf / 2] = rt * rt + it * it; +} + +void DeinitFFT() { + int i; + if (gFFTBitTable[0]) { + for (i=0;i<MaxFastBits;i++) { + free(gFFTBitTable[i]); + } + } +} + diff --git a/tests/test6/fft.h b/tests/test6/fft.h new file mode 100644 index 0000000..6066d52 --- /dev/null +++ b/tests/test6/fft.h @@ -0,0 +1,19 @@ +#ifndef FFT_H +#define FFT_H + +/* Everything here comes from Audacity 1.3.13 + (orignally in C++ and with more genericity and functionnality) + Original Author : Dominic Mazzoni + Licenced under GPL 2.0 (see LICENCE) +*/ + +#define PSNumS 256 +#define PSHalf 128 +void PowerSpectrum(float In[PSNumS], float Out[PSHalf]); + +void FFT(int NumSamples, int InverseTransform, + float *RealIn, float *ImagIn, float *RealOut, float *ImagOut); +void DeinitFFT(); + +#endif + diff --git a/tests/test6/test.raw b/tests/test6/test.raw new file mode 120000 index 0000000..a266d33 --- /dev/null +++ b/tests/test6/test.raw @@ -0,0 +1 @@ +../../../../musique-test/V clean pizz.raw
\ No newline at end of file diff --git a/tests/test6/test6.c b/tests/test6/test6.c index abaa3e4..e5a2a12 100644 --- a/tests/test6/test6.c +++ b/tests/test6/test6.c @@ -2,6 +2,8 @@ #include <stdio.h> #include <stdlib.h> +#include "compute.h" + typedef void (*cb_processdata_t)(int n, float *); @@ -60,10 +62,12 @@ void parse_testfile(cb_processdata_t cb) { FILE *fh=fopen("./test.raw", "r"); if (fh==NULL) return; - n=128; + //n=128; + n=512; while ( (n=fread(f, sizeof(float), n, fh)) > 0 ) { cb(n,f); - n=128+256*(rand()%7); + //n=128+256*(rand()%7); + n=512; } fclose(fh); } @@ -81,10 +85,19 @@ void process_mean_max(int n, float *f) { printf("%+.3f %+.3f %4i\n", mean, max, n); } +void process_level(int n, float *f) { + int rate=24000; + float level; + + level=compute_level(f, n, rate); + printf("%+.3f %4i\n", level, n); +} + int main() { //test_todb_a(); //dump_testfile(); - parse_testfile(process_mean_max); + //parse_testfile(process_mean_max); + parse_testfile(process_level); return 0; } |