1 //
2 // firpfbch_crcf_synthesis_example.c
3 //
4 // Example of the synthesis 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_synthesis_example.m"
17
main()18 int main() {
19 // options
20 unsigned int num_channels = 16; // number of channels
21 unsigned int m = 5; // 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_channels][num_frames]; // channelized input
34 float complex y[num_samples]; // time-domain output [size: num_samples x 1]
35
36 // create narrow-band pulse
37 unsigned int pulse_len = 17; // pulse length [samples]
38 float bw = 0.30f; // pulse width (bandwidth)
39 float pulse[pulse_len]; // buffer
40 liquid_firdes_kaiser(pulse_len, bw, 50.0f, 0.0f, pulse);
41
42 // generate input signal(s)
43 int enabled[num_channels]; // signal enabled?
44 for (i=0; i<num_channels; i++) {
45 // pseudo-random channel enabled flag
46 enabled[i] = ((i*37)%101)%2;
47
48 if (enabled[i]) {
49 // move pulse to input buffer
50 for (k=0; k<num_frames; k++)
51 x[i][k] = k < pulse_len ? bw*pulse[k] : 0.0f;
52 } else {
53 // channel disabled; set as nearly zero (-100 dB impulse)
54 for (k=0; k<num_frames; k++)
55 x[i][k] = (k == 0) ? 1e-5f : 0.0f;
56 }
57 }
58
59 // create prototype filter
60 unsigned int h_len = 2*num_channels*m + 1;
61 float h[h_len];
62 liquid_firdes_kaiser(h_len, 0.5f/(float)num_channels, As, 0.0f, h);
63
64 #if 0
65 // create filterbank channelizer object using internal method for filter
66 firpfbch_crcf q = firpfbch_crcf_create_kaiser(LIQUID_SYNTHESIZER, num_channels, m, As);
67 #else
68 // create filterbank channelizer object using external filter coefficients
69 firpfbch_crcf q = firpfbch_crcf_create(LIQUID_SYNTHESIZER, num_channels, 2*m, h);
70 #endif
71
72 // channelize input data
73 float complex v[num_channels];
74 for (i=0; i<num_frames; i++) {
75 // assemble input vector
76 for (k=0; k<num_channels; k++)
77 v[k] = x[k][i];
78
79 // execute synthesis filter bank
80 firpfbch_crcf_synthesizer_execute(q, v, &y[i*num_channels]);
81 }
82
83 // destroy channelizer object
84 firpfbch_crcf_destroy(q);
85
86 //
87 // export results to file
88 //
89 FILE * fid = fopen(OUTPUT_FILENAME,"w");
90 fprintf(fid,"%% %s: auto-generated file\n\n", OUTPUT_FILENAME);
91 fprintf(fid,"clear all;\n");
92 fprintf(fid,"close all;\n");
93 fprintf(fid,"num_channels = %u;\n", num_channels);
94 fprintf(fid,"m = %u;\n", m);
95 fprintf(fid,"num_frames = %u;\n", num_frames);
96 fprintf(fid,"h_len = 2*num_channels*m+1;\n");
97 fprintf(fid,"num_samples = num_frames*num_channels;\n");
98
99 fprintf(fid,"h = zeros(1,h_len);\n");
100 fprintf(fid,"x = zeros(num_channels, num_frames);\n");
101 fprintf(fid,"y = zeros(1,num_samples);\n");
102
103 // save prototype filter
104 for (i=0; i<h_len; i++)
105 fprintf(fid," h(%6u) = %12.4e;\n", i+1, h[i]);
106
107 // save channelized input signals
108 for (i=0; i<num_frames; i++) {
109 for (k=0; k<num_channels; k++) {
110 float complex v = x[k][i];
111 fprintf(fid," x(%3u,%6u) = %12.4e + 1i*%12.4e;\n", k+1, i+1, crealf(v), cimagf(v));
112 }
113 }
114
115 // save output time signal
116 for (i=0; i<num_samples; i++)
117 fprintf(fid," y(%6u) = %12.4e + 1i*%12.4e;\n", i+1, crealf(y[i]), cimagf(y[i]));
118
119 // compute the PSD of each input and plot results on a square grid
120 fprintf(fid,"nfft = 1024;\n"); // TODO: use nextpow2
121 fprintf(fid,"f = [0:(nfft-1)]/nfft - 0.5;\n");
122 fprintf(fid,"n = ceil(sqrt(num_channels));\n");
123 fprintf(fid,"figure;\n");
124 fprintf(fid,"for i=1:num_channels,\n");
125 fprintf(fid," X = 20*log10(abs(fftshift(fft(x(i,:),nfft))));\n");
126 fprintf(fid," subplot(n,n,i);\n");
127 fprintf(fid," plot(f, X, 'Color', [0.25 0 0.25], 'LineWidth', 1.5);\n");
128 fprintf(fid," axis([-0.5 0.5 -120 20]);\n");
129 fprintf(fid," grid on;\n");
130 fprintf(fid," title(num2str(i-1));\n");
131 fprintf(fid,"end;\n");
132
133 // plot results
134 fprintf(fid,"\n");
135 fprintf(fid,"H = 20*log10(abs(fftshift(fft(h/num_channels,nfft))));\n");
136 fprintf(fid,"Y = 20*log10(abs(fftshift(fft(y/num_channels,nfft))));\n");
137 fprintf(fid,"figure;\n");
138 fprintf(fid,"subplot(2,1,1);\n");
139 fprintf(fid," plot(f, H, 'Color', [0 0.5 0.25], 'LineWidth', 2);\n");
140 fprintf(fid," axis([-0.5 0.5 -100 10]);\n");
141 fprintf(fid," grid on;\n");
142 fprintf(fid," xlabel('Normalized Frequency [f/F_s]');\n");
143 fprintf(fid," ylabel('Prototype Filter PSD');\n");
144 fprintf(fid,"subplot(2,1,2);\n");
145 fprintf(fid," plot(f, Y, 'Color', [0 0.25 0.5], 'LineWidth', 2);\n");
146 fprintf(fid," axis([-0.5 0.5 -100 0]);\n");
147 fprintf(fid," grid on;\n");
148 fprintf(fid," xlabel('Normalized Frequency [f/F_s]');\n");
149 fprintf(fid," ylabel('Output PSD');\n");
150
151 fclose(fid);
152 printf("results written to %s\n", OUTPUT_FILENAME);
153
154 printf("done.\n");
155 return 0;
156 }
157
158