1 //
2 // firpfbch_crcf_analysis_example.c
3 //
4 // Example of the analysis channelizer filterbank. The input signal is
5 // comprised of several signals spanning different frequency bands. The
6 // channelizer downconverts each to baseband (maximally decimated), and
7 // the resulting spectrum of each is plotted.
8 //
9
10 #include <stdio.h>
11 #include <math.h>
12 #include <complex.h>
13
14 #include "liquid.h"
15
16 #define OUTPUT_FILENAME "firpfbch_crcf_analysis_example.m"
17
main()18 int main() {
19 // options
20 unsigned int num_channels = 8; // number of channels
21 unsigned int m = 4; // filter delay
22 float As = 60; // stop-band attenuation
23 unsigned int num_frames = 25; // number of frames
24
25 //
26 unsigned int i;
27 unsigned int k;
28
29 // derived values
30 unsigned int num_samples = num_frames * num_channels;
31
32 // data arrays
33 float complex x[num_samples]; // time-domain input [size: num_samples x 1 ]
34 float complex y[num_samples]; // channelized output [size: num_channels x num_frames]
35
36 // initialize input with zeros
37 for (i=0; i<num_samples; i++)
38 x[i] = 0.0f;
39
40 // generate input signal(s)
41 unsigned int num_signals = 4;
42 float fc[4] = {0.0f, 0.25f, 0.375f, -0.375f}; // center frequencies
43 float bw[4] = {0.035f, 0.035f, 0.035f, 0.035f}; // bandwidths
44 unsigned int pulse_len = 137;
45 float pulse[pulse_len];
46 for (i=0; i<num_signals; i++) {
47 // create pulse
48 liquid_firdes_kaiser(pulse_len, bw[i], 50.0f, 0.0f, pulse);
49
50 // add pulse to input signal with carrier offset
51 for (k=0; k<pulse_len; k++)
52 x[k] += pulse[k] * cexpf(_Complex_I*2*M_PI*fc[i]*k) * bw[i];
53 }
54
55 // create prototype filter
56 unsigned int h_len = 2*num_channels*m + 1;
57 float h[h_len];
58 liquid_firdes_kaiser(h_len, 0.5f/(float)num_channels, As, 0.0f, h);
59
60 #if 0
61 // create filterbank channelizer object using internal method for filter
62 firpfbch_crcf q = firpfbch_crcf_create_kaiser(LIQUID_ANALYZER, num_channels, m, As);
63 #else
64 // create filterbank channelizer object using external filter coefficients
65 firpfbch_crcf q = firpfbch_crcf_create(LIQUID_ANALYZER, num_channels, 2*m, h);
66 #endif
67
68 // channelize input data
69 for (i=0; i<num_frames; i++) {
70 // execute analysis filter bank
71 firpfbch_crcf_analyzer_execute(q, &x[i*num_channels], &y[i*num_channels]);
72 }
73
74 // destroy channelizer object
75 firpfbch_crcf_destroy(q);
76
77 //
78 // export results to file
79 //
80 FILE * fid = fopen(OUTPUT_FILENAME,"w");
81 fprintf(fid,"%% %s: auto-generated file\n\n", OUTPUT_FILENAME);
82 fprintf(fid,"clear all;\n");
83 fprintf(fid,"close all;\n");
84 fprintf(fid,"num_channels = %u;\n", num_channels);
85 fprintf(fid,"m = %u;\n", m);
86 fprintf(fid,"num_frames = %u;\n", num_frames);
87 fprintf(fid,"h_len = 2*num_channels*m+1;\n");
88 fprintf(fid,"num_samples = num_frames*num_channels;\n");
89
90 fprintf(fid,"h = zeros(1,h_len);\n");
91 fprintf(fid,"x = zeros(1,num_samples);\n");
92 fprintf(fid,"y = zeros(num_channels, num_frames);\n");
93
94 // save prototype filter
95 for (i=0; i<h_len; i++)
96 fprintf(fid," h(%6u) = %12.4e;\n", i+1, h[i]);
97
98 // save input signal
99 for (i=0; i<num_samples; i++)
100 fprintf(fid," x(%6u) = %12.4e + 1i*%12.4e;\n", i+1, crealf(x[i]), cimagf(x[i]));
101
102 // save channelized output signals
103 for (i=0; i<num_frames; i++) {
104 for (k=0; k<num_channels; k++) {
105 float complex v = y[i*num_channels + k];
106 fprintf(fid," y(%3u,%6u) = %12.4e + 1i*%12.4e;\n", k+1, i+1, crealf(v), cimagf(v));
107 }
108 }
109
110 // plot results
111 fprintf(fid,"\n");
112 fprintf(fid,"nfft = 1024;\n"); // TODO: use nextpow2
113 fprintf(fid,"f = [0:(nfft-1)]/nfft - 0.5;\n");
114 fprintf(fid,"H = 20*log10(abs(fftshift(fft(h/num_channels,nfft))));\n");
115 fprintf(fid,"X = 20*log10(abs(fftshift(fft(x,nfft))));\n");
116 fprintf(fid,"figure;\n");
117 fprintf(fid,"subplot(2,1,1);\n");
118 fprintf(fid," plot(f, H, 'Color', [0 0.5 0.25], 'LineWidth', 2);\n");
119 fprintf(fid," axis([-0.5 0.5 -100 10]);\n");
120 fprintf(fid," grid on;\n");
121 fprintf(fid," xlabel('Normalized Frequency [f/F_s]');\n");
122 fprintf(fid," ylabel('Prototype Filter PSD');\n");
123 fprintf(fid,"subplot(2,1,2);\n");
124 fprintf(fid," plot(f, X, 'Color', [0 0.25 0.5], 'LineWidth', 2);\n");
125 fprintf(fid," axis([-0.5 0.5 -100 0]);\n");
126 fprintf(fid," grid on;\n");
127 fprintf(fid," xlabel('Normalized Frequency [f/F_s]');\n");
128 fprintf(fid," ylabel('Input PSD');\n");
129
130 // compute the PSD of each output and plot results on a square grid
131 fprintf(fid,"n = ceil(sqrt(num_channels));\n");
132 fprintf(fid,"figure;\n");
133 fprintf(fid,"for i=1:num_channels,\n");
134 fprintf(fid," Y = 20*log10(abs(fftshift(fft(y(i,:),nfft))));\n");
135 fprintf(fid," subplot(n,n,i);\n");
136 fprintf(fid," plot(f, Y, 'Color', [0.25 0 0.25], 'LineWidth', 1.5);\n");
137 fprintf(fid," axis([-0.5 0.5 -120 20]);\n");
138 fprintf(fid," grid on;\n");
139 fprintf(fid," title(num2str(i-1));\n");
140 fprintf(fid,"end;\n");
141
142 fclose(fid);
143 printf("results written to %s\n", OUTPUT_FILENAME);
144
145 printf("done.\n");
146 return 0;
147 }
148
149