1 //
2 // firdecim_crcf_example.c
3 //
4 // This example demonstrates the interface to the firdecim (finite
5 // impulse response decimator) family of objects.
6 //
7 // SEE ALSO: firinterp_crcf_example.c
8 //
9 
10 #include <stdio.h>
11 #include <stdlib.h>
12 #include <math.h>
13 #include <getopt.h>
14 
15 #include "liquid.h"
16 
17 #define OUTPUT_FILENAME "firdecim_crcf_example.m"
18 
19 // print usage/help message
usage()20 void usage()
21 {
22     printf("firdecim_crcf_example:\n");
23     printf("  -h         : print usage/help\n");
24     printf("  -M <decim> : decimation factor, M > 1           default: 2\n");
25     printf("  -m <delay> : filter delay (symbols), m > 0,     default: 2\n");
26     printf("  -s <atten> : filter stop-band attenuation [dB], default: 60\n");
27     printf("  -n <num>   : number of samples (after decim),   default: 8\n");
28 }
29 
30 
main(int argc,char * argv[])31 int main(int argc, char*argv[]) {
32     // options
33     unsigned int M           = 6;       // decimation factor
34     unsigned int m           = 8;       // filter delay
35     float        As          = 60.0f;   // filter stop-band attenuation
36     unsigned int num_samples = 120;     // number of samples (after decim)
37 
38     int dopt;
39     while ((dopt = getopt(argc,argv,"hM:m:s:n:")) != EOF) {
40         switch (dopt) {
41         case 'h': usage();                    return 0;
42         case 'M': M           = atoi(optarg); break;
43         case 'm': m           = atoi(optarg); break;
44         case 's': As          = atof(optarg); break;
45         case 'n': num_samples = atoi(optarg); break;
46         default:
47             usage();
48             return 1;
49         }
50     }
51 
52     // validate options
53     if (M < 2) {
54         fprintf(stderr,"error: %s, decim factor must be greater than 1\n", argv[0]);
55         return 1;
56     } else if (m < 1) {
57         fprintf(stderr,"error: %s, filter delay must be greater than 0\n", argv[0]);
58         return 1;
59     } else if (As <= 0.0) {
60         fprintf(stderr,"error: %s, stop-band attenuation must be greater than zero\n", argv[0]);
61         return 1;
62     } else if (num_samples < 1) {
63         fprintf(stderr,"error: %s, must have at least one sample\n", argv[0]);
64         return 1;
65     }
66 
67     // data arrays
68     float complex x[M*num_samples]; // number of samples before decimation
69     float complex y[  num_samples]; // number of samples after decimation
70 
71     // initialize input array
72     unsigned int i;
73     unsigned int w_len = (unsigned int)(0.9*M*num_samples);
74     float f0 = 0.017f;
75     float f1 = 0.021f;
76     for (i=0; i<M*num_samples; i++) {
77         x[i]  = 0.6f*cexpf(_Complex_I*2*M_PI*f0*i);
78         x[i] += 0.4f*cexpf(_Complex_I*2*M_PI*f1*i);
79         x[i] *= (i < w_len) ? hamming(i,w_len) : 0;
80     }
81 
82     // create deimator object adn set scale
83     firdecim_crcf decim = firdecim_crcf_create_kaiser(M, m, As);
84     firdecim_crcf_set_scale(decim, 1.0f/(float)M);
85 
86     // execute decimator
87     firdecim_crcf_execute_block(decim, x, num_samples, y);
88 
89     // destroy decimator object
90     firdecim_crcf_destroy(decim);
91 
92     //
93     // export results to file
94     //
95     FILE*fid = fopen(OUTPUT_FILENAME,"w");
96     fprintf(fid,"%% %s : auto-generated file\n", OUTPUT_FILENAME);
97     fprintf(fid,"clear all;\n");
98     fprintf(fid,"close all;\n");
99     fprintf(fid,"M  = %u;\n", M);
100     fprintf(fid,"m  = %u;\n", m);
101     fprintf(fid,"num_samples=%u;\n", num_samples);
102 
103     // inputs
104     for (i=0; i<M*num_samples; i++)
105         fprintf(fid,"x(%3u) = %12.4e + j*%12.4e;\n", i+1, crealf(x[i]), cimagf(x[i]));
106 
107     // outputs
108     for (i=0; i<num_samples; i++)
109         fprintf(fid,"y(%3u) = %12.4e + j*%12.4e;\n", i+1, crealf(y[i]), cimagf(y[i]));
110 
111     // plot results
112     fprintf(fid,"figure('position',[100 100 600 800]);\n");
113     fprintf(fid,"tx = [0:(M*num_samples-1)];\n");
114     fprintf(fid,"ty = [0:(  num_samples-1)]*M - M*m;\n");
115     fprintf(fid,"nfft=3*2^nextpow2(M*num_samples);\n");
116     fprintf(fid,"fx = [0:(nfft-1)]/nfft-0.5;\n");
117     fprintf(fid,"fy = fx/M;\n");
118     fprintf(fid,"X=20*log10(abs(fftshift(fft(  x,nfft))));\n");
119     fprintf(fid,"Y=20*log10(abs(fftshift(fft(M*y,nfft))));\n");
120     fprintf(fid,"subplot(3,1,1);\n");
121     fprintf(fid,"  plot(tx,real(x),'-',ty,real(y),'s');\n");
122     fprintf(fid,"  grid on;\n");
123     fprintf(fid,"  axis([-M*m M*num_samples -1.2 1.2]);\n");
124     fprintf(fid,"  xlabel('Input sample index');\n");
125     fprintf(fid,"  ylabel('Real');\n");
126     fprintf(fid,"  legend('original','decimated','location','northeast');");
127     fprintf(fid,"subplot(3,1,2);\n");
128     fprintf(fid,"  plot(tx,imag(x),'-',ty,imag(y),'s');\n");
129     fprintf(fid,"  grid on;\n");
130     fprintf(fid,"  axis([-M*m M*num_samples -1.2 1.2]);\n");
131     fprintf(fid,"  xlabel('Input sample index');\n");
132     fprintf(fid,"  ylabel('Imag');\n");
133     fprintf(fid,"  legend('original','decimated','location','northeast');");
134     fprintf(fid,"subplot(3,1,3);\n");
135     fprintf(fid,"  plot(fx,X,'Color',[0.5 0.5 0.5],fy,Y,'LineWidth',2);\n");
136     fprintf(fid,"  grid on;\n");
137     fprintf(fid,"  axis([-0.5 0.5 -40 60]);\n");
138     fprintf(fid,"  xlabel('normalized frequency');\n");
139     fprintf(fid,"  ylabel('PSD [dB]');\n");
140     fprintf(fid,"  legend('original/real','transformed/decimated','location','northeast');");
141 
142     fclose(fid);
143     printf("results written to %s\n", OUTPUT_FILENAME);
144 
145     printf("done.\n");
146     return 0;
147 }
148