GCC Code Coverage Report

Directory: ../src/
File: /home/joels/Current/lispbm/src/extensions/dsp_extensions.c
Date: 2025-10-28 15:15:18
Exec Total Coverage
Lines: 115 118 97.5%
Functions: 4 4 100.0%
Branches: 54 64 84.4%
Line Branch Exec Source
1 /*
2 Copyright 2025 Joel Svensson svenssonjoel@yahoo.se
3
4 This program is free software: you can redistribute it and/or modify
5 it under the terms of the GNU General Public License as published by
6 the Free Software Foundation, either version 3 of the License, or
7 (at your option) any later version.
8
9 This program is distributed in the hope that it will be useful,
10 but WITHOUT ANY WARRANTY; without even the implied warranty of
11 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12 GNU General Public License for more details.
13
14 You should have received a copy of the GNU General Public License
15 along with this program. If not, see <http://www.gnu.org/licenses/>.
16 */
17
18 #include "extensions.h"
19 #include "lbm_utils.h"
20
21 #include <math.h>
22 #include <string.h>
23 #include <stdio.h>
24
25 #ifdef LBM_OPT_DSP_EXTENSIONS_SIZE
26 #pragma GCC optimize ("-Os")
27 #endif
28 #ifdef LBM_OPT_DSP_EXTENSIONS_SIZE_AGGRESSIVE
29 #pragma GCC optimize ("-Oz")
30 #endif
31
32 #if defined(__BYTE_ORDER__) && __BYTE_ORDER__ == __ORDER_LITTLE_ENDIAN__
33 #define LBM_SYSTEM_LITTLE_ENDIAN true
34 #else
35 #define LBM_SYSTEM_LITTLE_ENDIAN false
36 #endif
37
38
39 static lbm_uint sym_inverse = 0;
40 static lbm_uint sym_big_endian = 0;
41 static lbm_uint sym_little_endian = 0;
42
43 426 static inline uint32_t byte_order_swap(uint32_t w) {
44 426 uint32_t r = w;
45 426 uint8_t *bytes = (uint8_t*)&r;
46 426 uint8_t tmp = bytes[0];
47 426 bytes[0] = bytes[3];
48 426 bytes[3] = tmp;
49 426 tmp = bytes[1];
50 426 bytes[1] = bytes[2];
51 426 bytes[2] = tmp;
52
53 426 return r;
54 }
55
56
57 23 void lbm_fft(float *real, float *imag, int n, int inverse) {
58 23 int k = 0;
59
2/2
✓ Branch 0 taken 107 times.
✓ Branch 1 taken 23 times.
 130 for (int i = 1; i < n; i++) {
60 107 int bit = n >> 1;
61
2/2
✓ Branch 0 taken 54 times.
✓ Branch 1 taken 107 times.
 161 while (k & bit) {
62 54 k ^= bit;
63 54 bit >>= 1;
64 }
65 107 k ^= bit;
66
67
2/2
✓ Branch 0 taken 36 times.
✓ Branch 1 taken 71 times.
 107 if (i < k) {
68 36 float temp = real[i];
69 36 real[i] = real[k];
70 36 real[k] = temp;
71
72 36 temp = imag[i];
73 36 imag[i] = imag[k];
74 36 imag[k] = temp;
75 }
76 }
77
78
2/2
✓ Branch 0 taken 53 times.
✓ Branch 1 taken 23 times.
 76 for (int len = 2; len <= n; len <<= 1) {
79
2/2
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 43 times.
 53 float angle = (inverse ? 2.0f : -2.0f) * (float)M_PI / (float)len;
80 53 float wlen_r = cosf(angle);
81 53 float wlen_i = sinf(angle);
82
83
2/2
✓ Branch 0 taken 107 times.
✓ Branch 1 taken 53 times.
 160 for (int i = 0; i < n; i += len) {
84 107 float w_r = 1.0f;
85 107 float w_i = 0.0f;
86
87
2/2
✓ Branch 0 taken 177 times.
✓ Branch 1 taken 107 times.
 284 for (int j = 0; j < len / 2; j++) {
88 177 float u_r = real[i + j];
89 177 float u_i = imag[i + j];
90 177 float v_r = real[i + j + len / 2];
91 177 float v_i = imag[i + j + len / 2];
92
93 177 float t_r = w_r * v_r - w_i * v_i;
94 177 float t_i = w_r * v_i + w_i * v_r;
95
96 177 real[i + j] = u_r + t_r;
97 177 imag[i + j] = u_i + t_i;
98 177 real[i + j + len / 2] = u_r - t_r;
99 177 imag[i + j + len / 2] = u_i - t_i;
100
101 177 float next_w_r = w_r * wlen_r - w_i * wlen_i;
102 177 float next_w_i = w_r * wlen_i + w_i * wlen_r;
103 177 w_r = next_w_r;
104 177 w_i = next_w_i;
105 }
106 }
107 }
108
109
2/2
✓ Branch 0 taken 5 times.
✓ Branch 1 taken 18 times.
 23 if (inverse) {
110
2/2
✓ Branch 0 taken 20 times.
✓ Branch 1 taken 5 times.
 25 for (int i = 0; i < n; i++) {
111 20 real[i] /= (float)n;
112 20 imag[i] /= (float)n;
113 }
114 }
115 23 }
116
117 27 static lbm_value ext_fft_f32(lbm_value *args, lbm_uint argn) {
118 27 lbm_value r = ENC_SYM_TERROR;
119
4/4
✓ Branch 0 taken 26 times.
✓ Branch 1 taken 1 times.
✓ Branch 2 taken 24 times.
✓ Branch 3 taken 2 times.
 53 if (argn >= 2 &&
120
2/2
✓ Branch 0 taken 23 times.
✓ Branch 1 taken 1 times.
 50 lbm_is_array_r(args[0]) &&
121 24 lbm_is_array_r(args[1])) {
122
123 23 int inverse = 0;
124 23 bool be = true;
125
126
3/4
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 15 times.
✓ Branch 2 taken 8 times.
✗ Branch 3 not taken.
 31 if (argn > 2 &&
127 8 lbm_is_symbol(args[2])) {
128 8 lbm_uint sym = lbm_dec_sym(args[2]);
129
2/2
✓ Branch 0 taken 5 times.
✓ Branch 1 taken 3 times.
 8 if (sym == sym_inverse) {
130 5 inverse = 1;
131
2/2
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 1 times.
 3 } else if (sym == sym_little_endian) {
132 2 be = false;
133 }
134 }
135
136
3/4
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 21 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
 25 if (argn > 3 &&
137 2 lbm_is_symbol(args[3])) {
138
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
 2 if (lbm_dec_sym(args[3]) == sym_little_endian) {
139 2 be = false;
140 }
141 }
142
143 23 bool swap_byte_order =
144 (be && LBM_SYSTEM_LITTLE_ENDIAN) ||
145 (!be && !LBM_SYSTEM_LITTLE_ENDIAN);
146
147 23 lbm_array_header_t *real_arr = lbm_dec_array_r(args[0]);
148 23 lbm_array_header_t *imag_arr = lbm_dec_array_r(args[1]);
149
150
2/4
✓ Branch 0 taken 23 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 23 times.
 23 if (!real_arr || !imag_arr) {
151 ✗ return ENC_SYM_TERROR;
152 }
153
154
1/2
✓ Branch 0 taken 23 times.
✗ Branch 1 not taken.
 23 if (real_arr->size == imag_arr->size &&
155
1/2
✓ Branch 0 taken 23 times.
✗ Branch 1 not taken.
 23 real_arr->size > sizeof(float) && // need more than 4 bytes for a float
156
1/2
✓ Branch 0 taken 23 times.
✗ Branch 1 not taken.
 46 (real_arr->size % sizeof(float) == 0)) {
157 23 lbm_uint n = real_arr->size / sizeof(float);
158 23 lbm_uint padded_size = n;
159
2/2
✓ Branch 0 taken 6 times.
✓ Branch 1 taken 17 times.
 23 if ((n & (n - 1)) != 0) { // not a power of two. needs padding.
160 6 padded_size = 1;
161
2/2
✓ Branch 0 taken 17 times.
✓ Branch 1 taken 6 times.
 23 while (padded_size < n) padded_size <<= 1;
162 }
163
164 lbm_value real_copy;
165 lbm_value imag_copy;
166
2/4
✓ Branch 0 taken 23 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 23 times.
✗ Branch 3 not taken.
 46 if (lbm_heap_allocate_array(&real_copy, padded_size * sizeof(float)) &&
167 46 lbm_heap_allocate_array(&imag_copy, padded_size * sizeof(float))) {
168 // cppcheck-suppress invalidPointerCast
169 23 float *real_data = (float*)real_arr->data;
170 // cppcheck-suppress invalidPointerCast
171 23 float *imag_data = (float*)imag_arr->data;
172
173 23 lbm_array_header_t *real_copy_arr = lbm_dec_array_r(real_copy);
174 23 lbm_array_header_t *imag_copy_arr = lbm_dec_array_r(imag_copy);
175 // cppcheck-suppress invalidPointerCast
176 23 float *real_copy_data = (float*)real_copy_arr->data;
177 // cppcheck-suppress invalidPointerCast
178 23 float *imag_copy_data = (float*)imag_copy_arr->data;
179 23 memset(real_copy_data, 0, padded_size * sizeof(float));
180 23 memset(imag_copy_data, 0, padded_size * sizeof(float));
181
182
2/2
✓ Branch 0 taken 19 times.
✓ Branch 1 taken 4 times.
 23 if (swap_byte_order) {
183
2/2
✓ Branch 0 taken 99 times.
✓ Branch 1 taken 19 times.
 118 for (lbm_uint w = 0; w < n; w ++) {
184 // cppcheck-suppress invalidPointerCast
185 99 uint32_t swapped = byte_order_swap(((uint32_t*)real_data)[w]);
186 // cppcheck-suppress invalidPointerCast
187 99 ((uint32_t*)real_copy_data)[w] = swapped;
188 // cppcheck-suppress invalidPointerCast
189 99 swapped = byte_order_swap(((uint32_t*)imag_data)[w]);
190 // cppcheck-suppress invalidPointerCast
191 99 ((uint32_t*)imag_copy_data)[w] = swapped;
192 }
193 } else {
194 4 memcpy(real_copy_data, real_data, real_arr->size);
195 4 memcpy(imag_copy_data, imag_data, real_arr->size);
196 }
197
198 23 lbm_fft(real_copy_data, imag_copy_data, (int)padded_size, inverse);
199
200
2/2
✓ Branch 0 taken 19 times.
✓ Branch 1 taken 4 times.
 23 if (swap_byte_order) {
201
2/2
✓ Branch 0 taken 114 times.
✓ Branch 1 taken 19 times.
 133 for (lbm_uint w = 0; w < padded_size; w ++) {
202 // cppcheck-suppress invalidPointerCast
203 114 uint32_t swapped = byte_order_swap(((uint32_t*)real_copy_data)[w]);
204 // cppcheck-suppress invalidPointerCast
205 114 ((uint32_t*)real_copy_data)[w] = swapped;
206 // cppcheck-suppress invalidPointerCast
207 114 swapped = byte_order_swap(((uint32_t*)imag_copy_data)[w]);
208 // cppcheck-suppress invalidPointerCast
209 114 ((uint32_t*)imag_copy_data)[w] = swapped;
210 }
211 }
212 23 r = lbm_cons(real_copy, imag_copy);
213 } else {
214 ✗ r = ENC_SYM_MERROR;
215 }
216 } else {
217 ✗ r = ENC_SYM_NIL;
218 }
219 }
220 27 return r;
221 }
222
223 286 void lbm_dsp_extensions_init(void) {
224 286 lbm_add_symbol("inverse", &sym_inverse);
225 286 lbm_add_symbol("little-endian", &sym_little_endian);
226 286 lbm_add_symbol("big-endian", &sym_big_endian);
227
228 286 lbm_add_extension("fft", ext_fft_f32);
229 286 }
230